カー・パリネロ法

カー・パリネロ法（カーパリネロほう、英: Car–Parrinello method、CP法）は、1985年、カー（R. Car）とパリネロ（M. Parrinello）によって考案されたバンド計算の手法である^[1]。従来用いられていた行列要素の対角化を行わずに固有値（及び固有ベクトル）を求めることにより、計算を大幅に高速化した。これにより、系の電子状態と共に、その構造の最適化（この部分は古典的分子動力学法を用いる）も可能とした。

名称に関して編集

本手法は、第一原理分子動力学法（first-principles molecular dynamics, FPMD, またはab initio molecular dynamics, AIMD）などとも呼ばれる。それぞれ微妙に異なるものを意味している場合もある。現在では、カー、パリネロによるオリジナルな手法が用いられることはほとんどなく、より効率化、高速化を図った手法に置き換わっている。ただし、これら現在主流となっている手法も広い意味でのカー・パリネロ法の範疇にあると言える。

カー・パリネロ法使用の利点編集

計算量の減少 : 対角化を使用しないため、計算量のオーダーをN³から最大 N log N 程度まで減らすことができる。N は使用する基底関数の展開数。
メモリの節約 : 同様にして必要なメモリも大体N²から、NMのオーダーにすることが出来る。M はバンドの数。ここで基底関数の数は、バンド数より一桁以上大きいことが前提（←平面波基底の場合、N ≫ M）
系の構造の分子動力学計算や最適化が電子状態計算と同時に行える。

ラグランジアンが出発点編集

L=\sum _{i,{\vec {k}}}\int m|{\dot {\psi }}_{i,{\vec {k}}}({\vec {r}})|^{2}d{\vec {r}}+\sum _{I}{1 \over 2}M_{I}{\dot {\vec {R_{I}}}}^{2}-E_{\text{tot}}[{\psi _{i,{\vec {k}}}},{{\vec {R}}_{I}}]+\sum _{i,j,{\vec {k}}}\Lambda _{i,j,{\vec {k}}}\left(\int {\psi ^{*}}_{i,{\vec {k}}}({\vec {r}})\psi _{j,{\vec {k}}}({\vec {r}})d{\vec {r}}-\delta _{i,j}\right)

i

: バンドの指標

{\vec {k}}

: k点

m

: 波動関数に対する仮想質量

\psi

: 波動関数

M_{I}

: イオン芯（原子核）部分の質量（Iは指標）

{\vec {R}}_{I}

: イオン芯の座標の位置ベクトル

E_{\text{tot}}[\cdot ]

: 系の全エネルギー

\Lambda _{i,j,{\vec {k}}}

: ラグランジュの未定係数

ドットは時間微分を表し、上付きの"*"は複素共役を表す。

\int {\psi ^{*}}_{i,{\vec {k}}}({\vec {r}})\psi _{j,{\vec {k}}}({\vec {r}})d{\vec {r}}=\delta _{i,j}

規格直交性（基底関数に対する制約）

{\frac {d}{dt}}\left({\frac {\delta L}{\delta {\dot {\psi ^{*}}}_{i,{\vec {k}}}}}\right)={\frac {\delta L}{\delta {\psi ^{*}}_{i,{\vec {k}}}}}

　波動関数部分

{\frac {d}{dt}}\left({\frac {\partial L}{\partial {\dot {{\vec {R}}_{I}}}}}\right)={\frac {\partial L}{\partial {\vec {R}}_{I}}}

←　イオン芯部分

以上から、以下の二つの運動方程式が得られる。t は時間（波動関数に対する時間tは仮想的な時間で実時間ではない）、δは変分（汎関数微分）、H はハミルトニアンを意味する。

{\begin{aligned}m{\ddot {\psi }}_{i,{\vec {k}}}({\vec {r}},t)=-{\frac {\delta E_{tot}}{\delta {\psi ^{*}}_{i,{\vec {k}}}({\vec {r}},t)}}+\sum _{j}\Lambda _{i,j,{\vec {k}}}\psi _{j,{\vec {k}}}({\vec {r}},t)=-H\psi _{i,{\vec {k}}}({\vec {r}},t)+\sum _{j}\Lambda _{i,j,{\vec {k}}}\psi _{j,{\vec {k}}}({\vec {r}},t)\\M_{I}{\ddot {{\vec {R}}_{I}}}=-\nabla _{{\vec {R}}_{I}}E\end{aligned}}

参照 : ラグランジュ力学

上記、最後の2式から系の電子状態及び系の構造の最適化を行う。波動関数に関しての仮想的な運動方程式は、時間の2階微分を1階微分に置き換えると最急降下法となる。他にもいくつかの手法（例 : 共役勾配法など）への発展形がある。

過去、現在、未来編集

カー・パリネロ法が出てくる以前は、電子状態の計算（バンド計算）を行いつつ、構造の最適化（〔準〕安定構造の探索）を同時に行うことは、計算（対角化）で求められる原子間に働く力から、手動でユニットセル内の原子を動かして次のステップに回すか（手動で動かす段階で一旦計算が終了してしまう）、計算量が膨大な対角化を用いて、計算から得られる力をもとにユニットセル内の原子を分子動力学を用いて動かすことが行われていた。いずれにしても、大変効率が悪く、扱える原子数はせいぜい数個のオーダーであった^{[注 1]}。

カー・パリネロ法の出現は、この扱える原子数を一気に数十個のオーダーに引き上げた。数十個という規模ならスラブ近似を用いれば比較的扱い易いシリコンの表面系の安定構造を求める計算や、構造最適化の動力学的な過程を追っていくことが可能となった。

初期のカー・パリネロ法では、系の原子（より正確にはイオン芯または原子核）を分子動力学によって解くのと同時に、電子の波動関数を仮想的に時間発展していくものとして、波動関数に関する仮想的な運動方程式を分子動力学的手法と同じ形式を用いて原子系と連動して解いていた（上記、ラグランジアンの式変形参照）。これは後により効率の良い方法へと発展していくこととなる。

当初は、専ら擬ポテンシャル＋平面波基底による方法で、カー・パリネロ的手法が実現されたが、その後、全電子手法であるAPW法に対し、カー・パリネロ的手法を取り入れたものが出現している^[2]^[3]。更に、タイトバインディング法との結合 (TBMD) や、混合基底を使った手法でもカー・パリネロ的手法を導入したものが登場している。

電子状態を解くための手法面でも、初めは電子の波動関数を最急降下法やベレの方法で逐次的に解く方法が使われたが、その後、共役勾配法（いくつかの改良版がある）や、より洗練された方法（例 : RMM-DIIS法、ダビッドソン法）が使われるようになっている。

波動関数を計算する部分（電子状態計算部分）も、原子の動力学と同時に解くのではなく、原子を動かす時間幅を大きくとるようにして、原子を動かす毎に電子状態部分を常にボルン‐オッペンハイマー面まで収束させるようにする手法が主流となっていった（ペインのアルゴリズム^[4]）。また、この手法では原子の移動（構造の最適化）も同時に行われるが、この時原子の移動量による波動関数の変化を外挿により予想する方法^[5]など、高速化、効率化を図る数多くの提案がなされている。

最近では、オーダーN法や、ハイブリッド法内での利用などの拡張も行われている。

温度に関して編集

ユニットセル内の原子（イオン芯）を分子動力学手法で扱う場合は有限温度での計算が可能であるが、一方、同時に行う電子状態計算部分は密度汎関数法の枠内での計算なので絶対零度での計算となる。ペインのアルゴリズムのように、原子を動かす毎に電子状態をボルン‐オッペンハイマー面に収束させる場合、計算上の支障はないが、電子状態計算はあくまで絶対零度での計算であることに注意が必要である。電子状態計算部分も有限温度へ拡張するアプローチも存在する。