#はじめに
とりあえず「新しい量子化学(上)」を読んでわかったことをまとめておきます。
実装はPython3.6を用いて行っています。
#制限Hartree-fock
量子化学ではShorodinger方程式を解き、波動関数を得ることを目指します。
ℋ|\Phi\rangle = ℰ|\Phi\rangle
Hartree-fock法では波動関数が一つのスレーター行列式で表せると仮定します。
|\Phi\rangle=\frac{1}{ \sqrt{N!} }
\begin{vmatrix}
\chi_{1} (1)& \cdots & \chi_{N}(1) \\
\vdots & \ddots & \vdots \\
\chi_{1} (N) & \cdots & \chi_{N} (N)
\end{vmatrix}\\
\\
\chi_{i}(j) はi番目のスピン軌道にj番目の電子が入っているということ\\
スピン軌道は規格直交しています。\\
\langle\chi_{a}|\chi_{b}\rangle= \int dx_{1}\chi *_{a}(1)\chi_{b}(1) = \delta_{ab}
スピン軌道の規格直交条件のもとEについて変分原理を用い、スピン軌道をユニタリー変換するとHartree-fock方程式が得られます。
f(1)|\chi_{a}\rangle = \epsilon_{a}|\chi_{a}\rangle \quad(a=1,...,N)\\
|\chi_{a}\rangleはユニタリ―変換後のスピン軌道\\
fはFock演算子と呼ばれ以下のように定義されます。\\
f(1)=h(1)+v^{HF}(1)\\
h(1)=-\frac{1}{2}\nabla^{2}_{1} - \sum_{A=1} \frac{Z_{A} }{r_{1A} }\\
Z_{A}は核電荷、r_{1A}は核Aと電子1の距離\\
v^{HF}(1)=\sum_{b} C_{b}(1)-K_{b}(1)\\
総和はすべての電子についてとります。(bは電子についてのインデックス)\\
C_{b}(1),K_{b}(1)は以下のように定義されます。\\
C_{b}(1)\chi_{a}(1)=\biggl(\int dx_{2}\chi *_{b}(2)r^{-1}_{12}\chi_{b}(2)\biggr)\chi_{a}(1)\\
K_{b}(1)\chi_{a}(1)=\biggl(\int dx_{2}\chi *_{b}(2)r^{-1}_{12}\chi_{a}(2)\biggr)\chi_{b}(1)\\
C_{b}はクーロン演算子、K_{b}は交換演算子と呼ばれます。\\
スピン軌道は空間関数とスピン関数の積であらわされます。
\chi_{i}(x) = \psi_{i}(r)\alpha(\omega) \quad or \quad \chi_{i}(x) = \psi_{i}(r)\beta(\omega)\\
\psi_{i}(r)はi番目の空間軌道\\
\alpha(\omega),\beta(\omega)はそれぞれ上向きスピン関数、下向きスピン関数\\
閉殻系(電子が2個ずつですべての軌道を埋めている系)では空間関数は同じであると考え、以下のように制限を付けます。
\chi_{i}(x) = \left\{
\begin{array}{ll}
\psi_{j}(r)\alpha(\omega) & (i=1,3,...,N-1)\\
\psi_{j}(r)\beta(\omega) & (i=2,4,...,N)\\
\end{array}
\right.\\
(j=1,2,...n) \quad n=\frac{N}{2}\\
このときHartree-fock方程式は形を変えません。\\
f(1)\psi_{a}(r_1) = \epsilon_{a}\psi_{a}(r_1) \quad (a=1,...,n)\\
次に基底関数を用いて空間関数を展開します。
\{ \phi_{\mu}(r)|\mu = 1,2,...K \}\\
\psi_{i}(r)=\sum_{\mu}^{K}C_{\mu i}\phi_{\mu}(r) \quad (i=1,...,K)\\
ここでiが1からKまでなので占有軌道だけでなく非占有軌道も含んでいることに注意してください。
これをHartree-Fock方程式に代入すると
添え字を\mu からvに変えました。\\
f(1) \sum_{v}^{K}C_{v i}\phi_{v}(r_1)= \epsilon_{i}\sum_{v}^{K}C_{v i}\phi_{v}(r_1) \quad (i=1,...,n)\\
両辺、左から\phi*_{\mu}(r_1)をかけて積分すると\\
\sum_{v}^{K}C_{v i}\int dr_1 \phi*_{\mu}(r_1)f(1)\phi_{v}(r_1)=\epsilon_{i}\sum_{v}^{K}C_{v i} \int dr_1 \phi*_{\mu}(r_1)\phi_{v}(r_1)
ここで以下の要素を持つK×Kの重なり行列Sを定義します。
重なり行列Sは対称行列です。
S_{\mu v}=\int dr_1\phi*_{\mu}(r_1)\phi_{v}(r_1)\\
次にK×KのFock行列を定義します。
これはfock演算子の基底\{\phi_{\mu}\}による表現行列です。\\
F_{\mu v}=\int dr_1 \phi*_{\mu}(r_1)f(1)\phi_{v}(r_1)\\
重なり行列SとFock行列Fを用いると
\sum_{v}^{K}F_{\mu v}C_{v i}=\epsilon_{i}\sum_{v}^{K}S_{\mu v}C_{v i} \quad (i=1,...K)
これをRoothaan方程式と呼びます。
さらに係数CについてのK×K行列を定義します。
C=
\left(
\begin{array}{ccc}
C_{11}& \cdots & C_{1K} \\
\vdots & \ddots & \vdots \\
C_{K1} & \cdots & C_{KK}
\end{array}\right)\\
またεについての対角行列を定義します。
\epsilon=
\left(
\begin{array}{ccccc}
\epsilon_1& 0 & 0 \\
0 & \ddots & 0 \\
0 & 0 & \epsilon_K
\end{array}\right)\\
\epsilon_{i j}=\epsilon_{i}\delta_{i j}
これらの行列を用いるとK個のRoothaanの方程式が以下の1つの行列方程式で書けます。
FC=SC\epsilon
これを解きやすい形にするために、対称直交化または正準直交化で以下の条件を満たすK×K行列Xを見つけます。
X^{†}SX=1 \quad (1はK×Kの単位行列)\\
†はエルミート共役(複素共役をとって転置)
これを使って
C=XC'
となるC'を考えます。行列方程式に代入すると
FXC'=SXC'\epsilon\\
左からX^{†}をかけて\\
X^{†}FXC'=X^{†}SXC'\epsilon\\
F'=X^{†}FXとおくと\\
F'C'=C'\epsilon\\
ここまで簡単になります。
C'=(C'_1 ,...,C'_K) \quad (C'_i は列ベクトル)とすると\\
よってF'C'=(F'C'_1,...,F'C'_K)\\
またC'\epsilon=(\epsilon_1C'_1,...,\epsilon_KC'_K)\\
よって列で比較して
F'C'_i=\epsilon_iC'_i \quad (i=1,...,K)\\
つまり変換後のFock行列の固有値、固有ベクトルをならべて行列をつくった物が行列方程式の解C'、εです。
C=XC'
これからCが求められ波動関数が求められます。
しかし実際にはFock行列が係数行列Cに依存するので(fock演算子自体がスピン軌道に依存するため)繰り返して解きます。(解くべき行列方程式は以下のようになっています。)
F(C)C=SC\epsilon
係数行列Cの初期値を決定し(実際には密度関数Pに関して初期値を定めます)、これに対してFock行列を定めます。上記の手続きで行列方程式を解き、新たに係数行列Cを得ます。この操作を繰り返し、1つ前の係数行列Cと行列方程式から得たCが十分に収束したら計算は終了で波動関数が求められたことになります。(実際にはエネルギーEが10^(-6)オーダーまで一致したら収束したことにしています。)
予想より大分長くなってしまったので制限Hartree-Fockの概要をPart1としてPart2に具体的なコードを紹介します。
#まとめ
大分長くなったので、仮定だけまとめておきます。
- 波動関数は一つのスレーター行列式で表せる。(Post-Hartree-Fockの配置間相互作用につながる)
- 閉殻系のスピン軌道には同じ空間軌道を持ち、スピンの異なる軌道がペアで存在する。
- 空間関数は有限個の基底関数で展開できる。(基底関数を無限個にすると実質的に誤差がなくなる。これをHartree-Fock極限という)
#参考文献
[1] A.ザボ, N.S.オストランド『新しい量子化学―電子構造の理論入門〈上〉』東京大学出版会(1987)