構造最適化のための理論

 ここでは構造最適化のための理論について纏める。
-----------------------------------------------------------------------------
■ 一般対角化問題の式→ハウスホルダー変換→三重対角化→OR分解
  DVXa法 SCAT がこの方法で計算している。DVXa法 SCATは対角化で実数のみの計算となる。Gaussianのように結晶軌道法を導入しようとすると虚数が不可避になるため、結晶軌道法の導入の場合には対角化のプログラム(ソルバー)を変更する必要がある。

■ 最急降下法 [1]
 この極小化問題に対する原始的な解法は、まず試行点x0 で式の(-∇f(x) = b - A・x )を評価し勾配を求め、次にその勾配の方向に一次元最小化を行う。その線上の最小点x1 が定まったら、次にその点で再び勾配を求め、またその方向に最小化する。という過程を繰り返すことである。これが最急勾配法(steepest descent method)である。これは直感的に分かりやすい方法であるが、実は図1.3 に示されるようにポテンシャルの谷間では何度もステップを繰り返すことになり非常に効率が悪い。

■ CG法(Conjugate Gradient Method:共役勾配法)[2]
 この最急勾配法の欠点を良く解析すると、同じ探索方向を何度も繰り返している点が問題である。過去探索した方向は二度と繰り返さない方法がある。それが共役勾配法である。
 
■ ダヴィッドソン法 [3]
 ハウスホルダー法による行列全体の対角化を行うとN^3(基底関数^3)のオーダーが計算量が必要になり、Nが数万以上の大規模計算になると計算が非常に困難になる。しかしながら、通常すべての固有状態が必要になることはなく、エネルギー固有値の低い方から電子数程度Nbの固有状態のみ求めれば充分であることが多い。Nb << Nであれば、繰り返し対角化法を用いることによりハウスホルダー法による対角化よりも計算量をかなり減じることができる。(中略)
 ダヴィッドソン法は初期の試行波動関数が真の解に直交していなければ、大抵の場合正しい解に安定に収束する。そのため、自己無撞着計算の初期のイタレーションにはダヴィッドソン法を用いるとよい。しかし、ダヴィッドソン法では波動関数の直交化および部分空間での行列Aijの対角化はイタレーションごとに実行する必要がある。直交化および対角化のステップは系の大きさが大きくなってくると系のサイズの3乗で計算量が増えるため、負担が大きくなってくる。
※ 僅かに原子位置を変えても(=微分して)0になるようにして、エネルギーの変化がなくなるようにする。
※ 計算を早くするために、エネルギー固有値の低い方から電子数程度の固有状態のみを求める。加えて、補正ベクトルの計算には対角成分のみを計算している。
※ 残差ベクトル R ≒勾配と見なして計算していると考えてみると良い。

■ RMM-DIIS法 [3]
 RMM-DIIS(Residual Minimization Method-Direct Inversion in the Iterative Subspace)法は、系のサイズが大きくなると相対的な負担が大きくなる波動関数の直交化および部分空間での行列Aijの対角化を原理的には省略可能となるため、効率的な計算方法となる。(省略)RMM-DIIS法は、原理的には波動関数の直交化および部分空間での行列Aijの対角化を行わずとも、試行関数を最適化していくことにより残差ベクトルがゼロの真の固有ベクトルに収束していくはずであり、大きな系に対して計算量をかなり減らすことが可能である。しかしながら、試行関数は最も近い固有ベクトルに収束してく性質があり、初期に用意する試行関数は真の固有ベクトルに近い状態である必要がある。そのため、RMM-DIIS法を自己無撞着計算の最初から適応するとうまく収束せず、初期の計算はダビッドソン法を用いる必要がある。
※ 補正ベクトルの計算はダヴィッドソン法と同じ。それにxを掛けて波動関数に混ぜる。
※ 文献[3]のP.71の式(41)に記載されているxは収束を速くするためのパラメータだと考えると良い。

■ CP法(Car-Parrinello Method:カー・パリネロ法)[3,4]
 カー・バリネロ法の大きな特徴は、これまで絶対に守らねばならないと考えられてきた断熱近時を大胆にも破棄したことである。こうすることで、計算効率を飛躍的に発展させることができた。従来の解法では、原子を動かす前に波動関数をセルフコンシステントに解かなければいけない。すべてが一つ一つ逐次的に行われていた。CP法では、波動関数も原子位置もどちらが従属という関係ではなく、どちらも同等な扱いを受け、どちらに対しても時間発展を同時に計算する。いわば並列処理とでも呼べる。各時間ステップでは波動関数はセルフコンシステントになっていない。それゆえ計算時間は早くなる。
 ハミルトニアンの波動関数への作用(波動関数の真の値からの誤差)にならない成分は、すべてその波動関数に対する「力」として作用する。この場合、重要な点は、波動関数に関する時間発展が、通常のニュートンの運動方程式のように2次であることである。その結果、波動関数の真の値からの誤差はその波動関数の時間発展に関して「復元力」として働く。波動関数の真の値からの誤差が大きい場合、その時間発展はそれを解消する方向に作用する。このおかげで、波動関数の真の値からの誤差は常に断熱線の周りを振動し、「系統的誤差」は時間平均するとキャンセルすることになる。このような巧妙なカラクリのため、カー・パリネロ法はうまく機能するのである。
 実際のところ、カー・パリネロ法では、波動関数の直交規格化条件の課し方など詳細で違いが生じて、上で述べたメリットはいつでも成り立つわけではない。カー・バリネロ法がうまく動作する条件として「断熱条件からあまり離れていない」ことが挙げられる。このとき波動関数は断熱条件になっていないが、その周りをふらつきながら平均的には断熱曲線を追従することになる。図5.5はその様子を模式的に示している(よく講演などで目にするので、文献[3]を見ておくと良いだろう)。図5.5における断熱曲線からのズレは、だいたい電子の仮想的運動エネルギーKeくらいで与えられる。したがって、このエネルギーがイオンの実際の運動エネルギーKiに比べて小さい Ke << Ki が条件になる。

■ PORT
mixer MSEC1, MSR1, MSR1a

[1] 最急降下法, wikipedia 及び http://www.cmp.sanken.osaka-u.ac.jp/~koun/o2knano/docs/user.pdf
[2] 共役勾配法, wikipedia 及び http://www.cmp.sanken.osaka-u.ac.jp/~koun/o2knano/docs/user.pdf
[3] 笠井秀明ら著、計算機マテリアルデザイン入門、大阪大学出版会
[4] カー・パリネロ法, wikipedia
[5] アドバンスシミュレーション, 20 (2014) 75-85.(情報技術誌アドバンスシミュレーションは、アドバンスソフト株式会社ホームページのシミュレーション図書館からPDF形式がダウンロードできます)

◇ 構造最適化に関する情報[5]
 平面波基底を用いることは、周期性を持った材料の計算では高速フーリエ変換が利用可能であるため有効である。孤立分子系を計算する場合には、ガウス関数基底や原子軌道基底などの局在基底を用いる方が便利である。しかし、局在基底は平面波基底と異なり、基底関数が原子核の位置に依存する。そのため、原子に働く力を計算するとき、基底関数の原子位置による微分が計算として必要となる。このことから、原子位置の最適化を行うような計算は、平面波基底を用いる方が適している。表面・界面や孤立分子系の計算には、大きな真空層を挟んで計算対象が周期的に並んでいるスラブモデルを用いることによって、平面波基底を用いた計算が可能である。この点から、平面波基底を用いた計算は、孤立分子系(0次元)、チューブや原子チェーン(1次元)、表面・界面(2次元)、結晶(3次元)を同じスキームで行うことができる手法であるといえる。
-----------------------------------------------------------------------------
QRコード
携帯用QRコード
アクセス数
ページビュー数
[無料でホームページを作成] [通報・削除依頼]