【PySCF】CASSCF・CASSCIで高度な電子構造解析【Pythonで始める量子化学計算】

本記事では、量子化学計算フレームワークPySCFを使用して、CASSCFおよびCASCIという高度な手法を用いて酸素分子（O₂）の電子構造を解析します。記事では、PySCFを用いたCASSCFおよびCASCI計算の実行方法と結果の解析を具体的に紹介し、さらにTD-DFTとの比較を行います。これにより、各計算手法の利点と限界を理解し、研究において適切な電子構造計算手法を選択するための知識を得ることができます。

Gaussianを使った量子化学計算の初心者向け技術書を販売中

しばしば出くわすエラーへの対処法をはじめ
Gaussianと無料ソフトウェア Avogadro を組み合わせた物性解析手法が学べます！

技術書ページへ

量子化学計算用のソフトウェア: PySCFとは

参考文献[1]から引用

PySCF（Python Simulations of Chemistry Framework）は、Pythonベースのオープンソースの量子化学計算フレームワークです。PySCFは、高度な電子構造計算をシンプルかつ効率的に実行するためのツールを提供し、研究者や学生が量子化学の計算を容易に行えるように設計されています。

PySCFの主な特徴は以下の通りです：

• モジュラー設計：PySCFは、計算の各ステップ（分子の構造定義、基底状態の計算、励起状態の計算など）をモジュールとして提供しており、必要な部分だけを組み合わせて使用することができます。
• 高性能計算：内部に高度なアルゴリズムを使用しており、大規模な計算でも高速に処理を行うことができます。
• 柔軟性: Pythonベースであるため、スクリプトを通じて簡単に操作や拡張が可能です。ユーザーは、計算手法やパラメータをカスタマイズして特定のニーズに合わせることができます。
• 拡張性: オープンソースであるため、研究者や開発者が自分の研究に合わせてPySCFを拡張することが容易です。新しいアルゴリズムや方法論を追加することで、独自の研究を進めることができます。また、拡張機能が開発され続けている事も特徴です。

本記事では、PySCFを用いて酸素分子のCASSCFおよびCASCI計算を実行し、その結果を詳細に解析します。

CASSCFとCASCIとは？

CASSCFの基本概念

CASSCF（Complete Active Space Self-Consistent Field）法は、電子相関を厳密に取り扱うための高精度な量子化学計算手法の一つです。特に、基底状態と励起状態が密接に絡み合う系や、電子相関が強い系に対して有効です。

CASSCFでは、以下のような手順で計算を行います：

1. 活性空間の定義：全軌道の中から、化学的に重要な電子遷移が起こる軌道を選びます。これらの軌道と電子で構成される空間を「活性空間」と呼びます。
2. 完全配置相互作用（CI）計算：活性空間内で、可能な限りの電子配置を考慮したCI計算を行います。これにより、活性空間内の電子相関を厳密に取り扱います。
3. 軌道最適化：活性空間内の軌道と占有数を変化させながら、エネルギーを最小化します。これにより、軌道と電子配置が自己無撞着に最適化されます。

CASSCF法は、活性空間内の電子相関を完全に取り入れることで、高精度なエネルギー計算が可能となります。各軌道の役割を理解することが、正確な計算のためには重要です。

図: CASSCF計算における活性空間、非占有軌道、非活性軌道の概念

CASCIの基本概念

CASCI（Complete Active Space Configuration Interaction）法は、CASSCFと似た手法ですが、軌道の最適化を行わない点が異なります。つまり、事前に定義された軌道を用いて、活性空間内での完全CI計算のみを行います。

CASCIは以下の特徴を持ちます：

• 計算コストが低い：軌道最適化を行わないため、CASSCFよりも計算コストが低くなります。
• 軌道依存性：軌道の選択によって結果が大きく変わる可能性があります。

CASCIは、既知の軌道を用いて電子相関を解析する場合や、初期的な解析を行う際に有用です。

CASSCFとCASCIの違い

• 軌道最適化：CASSCFは軌道の最適化を行いますが、CASCIは行いません。
• 計算精度：CASSCFは軌道と電子配置を同時に最適化するため、より高い精度が得られます。
• 計算コスト：CASCIは計算コストが低く、計算時間が短いです。

理論背景

CASSCFとCASCIは、多配置自己無撞着場法（MCSCF）の一種であり、電子相関を適切に取り扱うための手法です。これらの手法は、以下のような系に対して適用されます：

• 開殻系：不対電子を含む系。
• 励起状態：基底状態と励起状態が混ざり合う系。
• 化学反応の遷移状態：結合の切断や生成が起こる状態。

これらの手法を用いることで、電子相関が強い系に対して高精度なエネルギー計算が可能となり、分子の電子構造や反応機構の理解が深まります。

CASSCF/CASCIを行う時の注意点

CASSCFおよびCASCI計算を行う際には、以下のポイントに注意する必要があります。

活性空間の選択

活性空間の選択は、CASSCF/CASCI計算の結果に大きな影響を与えます。活性空間には、化学的に重要な電子遷移や反応に関与する軌道を含める必要があります。一般的な指針として：

• 価電子軌道を含める：結合や反応に関与する軌道。
• 非占有軌道も含める：励起状態や反応過程で重要となる軌道。

活性空間が小さすぎると、重要な電子相関を見落とす可能性があります。一方で、大きすぎると計算コストが増大します。

基底関数の選択

基底関数セットの選択も重要です。CASSCF/CASCI計算では、電子相関を正確に取り扱うために、高品質な基底関数を使用することが推奨されます。具体的には：

• 拡張基底関数：6-31G以上の基底関数を使用。
• 分散関数や偏極関数を含める：電子の広がりや偏極を正確に表現。6-31＋G**の＋や＊がこれに当たる。

初期軌道の選択

CASCIでは軌道の最適化を行わないため、初期軌道の選択が結果に直結します。適切な初期軌道を選択することで、信頼性の高い結果を得ることができます。

計算コスト

CASSCFは計算コストが高いため、分子のサイズや活性空間の大きさに注意が必要です。計算資源や時間を考慮し、適切な妥協点を見つけることが重要です。

多重度とスピン状態

CASSCF/CASCI計算では、多重度やスピン状態を正確に設定する必要があります。特に開殻系や励起状態を扱う場合、スピン汎関数や多重度の設定が結果に大きく影響します。

主な特徴

TD-DFTとの比較

TD-DFT（Time-Dependent Density Functional Theory）は、励起状態のエネルギーや性質を計算する方法です。計算コストが低く大規模分子にも適用可能ですが、強相関系や多重度の変化を伴う遷移の扱いが難しい場合があります。

CASSCF/CASCIは強相関を持つ系や多重度の変化を伴う励起状態の解析に信頼性の高い結果を提供します。一方、TD-DFTは計算コストが低いため、大規模な分子系や高次の励起状態の計算に適していますが、強相関系では精度が低下する可能性があります。特に正確なHOMO-LUMOを計算するにはCASSCF/CASCIが重要です。

PySCFを用いた酸素分子のCASSCF・CASCI計算

以下では、PySCFを使用して酸素分子のCASSCFおよびCASCI計算を行い、その結果を解説します。以下実行コードです。

import pyscf  # PySCFライブラリをインポート

# 酸素分子（O₂）を定義
mol = pyscf.M(
    atom = 'O 0 0 0; O 0 0 1.2',  # 酸素原子の座標（オングストローム単位）
    basis = 'ccpvdz',              # 使用する基底関数系（cc-pVDZ）
    spin = 2)                       # スピン多重度（2S=2、S=1の三重項状態）

# 制限ハートリー・フォック（RHF）計算を実行
myhf = mol.RHF().run()

# Complete Active Space Self-Consistent Field（CASSCF）計算を実行
# 活性空間として6つの軌道と8つの電子を選択
mycas = myhf.CASSCF(6, 8).run()

# 注意: 以下の方法でもmycasオブジェクトを作成可能
#
# from pyscf import mcscf
# mycas = mcscf.CASSCF(myhf, 6, 8).run()

# 分析の詳細度を設定（4は詳細な出力）
mycas.verbose = 4

# CASSCF計算結果の解析を実行（自然軌道占有数やMulliken atomic chargesなど）
mycas.analyze()

必要なライブラリのインポート

import pyscf
from pyscf import gto, scf, mcscf

解説:

• pyscf：量子化学計算のためのライブラリ。
• gto：分子の定義。
• scf：ハートリー・フォック計算。
• mcscf：CASSCFやCASCI計算。

分子の定義

酸素分子（O₂）を定義します。酸素原子間の距離は1.2オングストローム、基底関数系にはcc-pVDZを使用し、スピン多重度は3（三重項、スピン数2）です。

mol = gto.M(
    atom = 'O 0 0 0; O 0 0 1.2',
    basis = 'ccpvdz',
    spin = 2  # 2Sの値（S=1）
)

ハートリー・フォック計算の実行

RHF（制限ハートリー・フォック）計算を行います。

myhf = scf.RHF(mol)
myhf.kernel()

CASCI計算の実行

活性空間として6つの軌道と8つの電子を選び、CASCI計算を行います。

ncas = 6
nelecas = 8

mycasci = mcscf.CASCI(myhf, ncas, nelecas)
mycasci.kernel()

CASSCF計算の実行

同じ活性空間を用いて、CASSCF計算を行います。

mycasscf = mcscf.CASSCF(myhf, ncas, nelecas)
mycasscf.kernel()

計算結果の解析

以上を実行すると、酸素分子のCASSCFとCASCI計算の設定と実行が完了します。次に、出てくる計算結果の詳細な解析と、CASSCFとCASCIの比較、そしてTD-DFTとの比較を行います。

ハートリー・フォック（HF）エネルギー

converged SCF energy = -149.609461121647

解説:

• HF計算により、酸素分子の基底状態エネルギーが-149.609461121647 ハートリーで収束しました。
• このエネルギーは、電子相関を考慮しない近似的なエネルギーです。

CASSCF・CASCIエネルギー

CASSCF energy = -149.708367180420
CASCI E = -149.708367180420  E(CI) = -21.8104798827154  S^2 = 2.0000000

• CASSCFエネルギーは-149.708367180420 ハートリーで、HFエネルギーよりも低くなっています。
• CASCI Eも同じ値であり、活性空間内の電子相関が反映されています。
• E(CI)は活性空間内の電子相関エネルギーを示しています。
• S² = 2.0000000はスピン多重度S=1（三重項状態）を示します。

自然軌道占有数

Natural occ [1.96299816 1.96299816 1.96098518 1.03637641 1.03637641 0.04026568]

解説:

• 最初の3つの軌道はほぼ完全に占有されています。
• 次の2つの軌道は部分的に占有されており、不対電子の存在を示しています。
• 最後の軌道はほとんど占有されていません。

自然軌道の展開

活性空間内の自然軌道が原子軌道（AO）の線形結合で表されている様子を示します。

Natural orbital (expansion on meta-Lowdin AOs) in CAS space
                   #1        #2        #3        #4        #5
0 O 2px         -0.63632   0.30560   0.00000   0.57258  -0.41480
0 O 2py          0.30560   0.63632  -0.00000   0.41480   0.57258
0 O 2pz         -0.00000  -0.00000  -0.69915   0.00000  -0.00000
...

解説:

• 軌道#1と#2は主にOの2pxおよび2py軌道から構成され、結合性のπ軌道を示しています。
• 軌道#3はOの2pz軌道からの寄与が大きく、σ結合軌道を表します。
• 軌道#4と#5は反結合性のπ*軌道で、不対電子が占有しています。
• 軌道#6は2s軌道や3s軌道からの寄与が見られますが、占有数が低いです。

MCSCFとHF軌道の重なり

<mo_coeff-mcscf|mo_coeff-hf>  1      1        0.99999759
<mo_coeff-mcscf|mo_coeff-hf>  2      2        0.99999999
...

解説:

• CASSCF計算で得られた分子軌道（MO）が、HF計算でのMOとどれだけ重なっているかを示します。
• 重なりが1に近い軌道は、CASSCF計算によってあまり変化していないことを示します。
• 重なりが変化している軌道は、電子相関の効果で軌道が再配置されていることを示します。

最大CI成分

 Largest CI components 
  [alpha occ-orbitals] [beta occ-orbitals]            CI coefficient
  [0 1 2 3 4]          [0 1 2]                         0.971982386835
  [0 1 2 3 4]          [0 3 4]                        -0.155027419535

解説:

• 最も寄与の大きい電子配置とその寄与率を示しています。
• 主要な配置はアルファ電子が軌道0～4を占有し、ベータ電子が軌道0～2を占有するもので、CI係数は約0.972です。
• 次に寄与する配置はベータ電子が軌道0, 3, 4を占有し、CI係数は約-0.155です。
• これらの結果は、酸素分子の基底状態が不対電子を持つことを反映しています。

Mulliken atomic charges

 Mulliken pop  
pop of  0 O 1s            1.99999
pop of  0 O 2s            1.79496
pop of  0 O 2px           1.49624
pop of  0 O 2py           1.49624
pop of  0 O 2pz           1.18121
...
 Mulliken atomic charges  
charge of    0O =      0.00000
charge of    1O =     -0.00000

解説:

• Mulliken atomic chargesでは、各原子軌道に存在する電子数を示します。
• 酸素原子の1s軌道はほぼ完全に占有されており、2sや2p軌道も高い占有数を持っています。
• 2p軌道は部分的に占有され、結合性および反結合性軌道の形成が見られます。
• 原子電荷はほぼゼロであり、電子が酸素原子間で対称的に分布していることを示します。

Mulliken atomic spins

Mulliken atomic spins:
spin of  0O =      1.00000
spin of  1O =      1.00000

解説:

• Mulliken atomic spinsでは、各原子における不対電子の分布を示します。
• 各酸素原子が1.00000のスピンを持ち、合計で2.00000となっています。
• これは、酸素分子が三重項状態（スピン多重度が3）であることを確認しています。

酸素分子のCASSCFの計算結果から分かる事

• 酸素分子は強相関系であり、二つの不対電子を持つ。
• CASSCF/CASCIはこのような系で正確な電子構造を得るのに適している。一方でTD-DFTでは、強相関や多重度の変化を正確に扱えない可能性があります。

おわりに

本記事では、PySCFを用いて酸素分子の電子構造をCASSCFとCASCIで解析し、その結果を詳細に解説しました。,また強相関系の取り扱いが必要な場合、CASSCFやCASCIが有効であることを示しました。これらの手法を理解し、自身の研究に活用してみてください！

参考文献

Gaussianを使った量子化学計算の初心者向け技術書を販売中

しばしば出くわすエラーへの対処法をはじめ
Gaussianと無料ソフトウェア Avogadro を組み合わせた物性解析手法が学べます！

技術書ページへ