地域密着型エリア広告配信リクルートの地域密着型広告ドコイク?アドネットワーク Ads by ドコイク?[無料でホームページを作成] [通報・削除依頼]
[無料でホームページを作成] [通報・削除依頼]

第一原理計算 理論の学習の仕方

 ここでは、第一原理計算における理論の学習の仕方を説明する。
------------------------------------------------------------------------------
※相関ポテンシャルについては下記のJSTのWebラーニングプラザも参照すると良い(下記の手順が終わった後に復習として見てみると良い)。 http://weblearningplaza.jst.go.jp/
------------------------------------------------------------------------------
1) まず始めに 『電子構造論による化学の探究』での「付録Aの理論的背景」を読む。(昔は一研究室、一冊に限り千円で『電子構造論による化学の探究』が購入できたそうだが、現在は6千円だそうです。後進の研究者には厳しい条件になり残念でもあります)
 これで分子軌道法は十分である。もう少しプログラミングまで踏み込みたい場合は、京都大名誉教授の足立裕彦先生が書かれた書籍を読み込むと良い。

2) 続いて、『固体の中の電子』で密度汎関数の記述の部分を読み込む。
 これで密度汎関数の理解はLDAからGGA、PBEまで具体的に理解が出来る。この本は、第一原理計算をしていて、何か気になった場合に、その都度、読み込んでいけば良いだろう。

3)  分子軌道法の場合はDVXα研究会の「夏の学校」に参加する。KKR法及びバンド計算法の場合はCMDに参加する。
  分子軌道法やKKR法とバンド計算法をFortranでプログラミングする場合の詳細がよく分かるようになる。更に分子軌道法を理解したい場合は『DV-Xα法による電子状態計算 -そのプログラムと解説-』を読むと良い。

4) AbinitやWIEN2kのマニュアル類や説明を読む。

5) 上記でまだ分からないことがあれば下記のHPを参照する。
A. 中山先生のhttp://homepage3.nifty.com/mnakayama/research/tips/qm.htm の概要を参照。その他は、実際に計算する場合に気になったところを読んでいけば良い。
B. ノーベル賞を受賞したホーエンベルンとコーンの理論を理解したければ、尾形先生のhttp://www-tph.cheme.kyoto-u.ac.jp/c/ms2007/note/ogata2/chap2_DFT.pdf を参照する。
C. スピン軌道相互作用: http://202.11.2.126/iwaya/ryoushi/h19/12-2.pdf
D. 『バンド理論―物質科学の基礎として 』でさらに詳細に密度汎関数法を理解すると・・・・・・。(これはかなり理解が難しい)
E. 『実践 量子化学計算プログラミング―アドバンスソフトシミュレーションシリーズ[2]』を読めば、具体的な励起状態計算のプログラミング方法が理解できるようになる。
F. 擬ポテンシャル法:平面波の係数とそのエネルギーを得るために対角化してそれらを計算するのではなく、共役勾配法( Conjugate Gradient Method : CG法 )を用いて最もエネルギーが低くなる平面波の係数を探す問題に変えて高速化している。通常の対角化計算→最急降下法→共役勾配法 の順に計算速度が速くなる。下記では擬ポテンシャル作成コード atomの記述もある。
http://www.cmp.sanken.osaka-u.ac.jp/~koun/o2knano/docs/user.pdf
http://www.cmp.sanken.osaka-u.ac.jp/~koun/o2knano/docs/user2.pdf 
http://www.cmp.sanken.osaka-u.ac.jp/~koun/o2knano/docs/user3.pdf 
http://rmp.aps.org/abstract/RMP/v64/i4/p1045_1
・ http://flex.phys.tohoku.ac.jp/riron/ronbun/m04mesa.pdf
G. KKR法: 下記を参照するのが良いだろう。英語版の方が記述が豊富な場合もあるのでチェックしておくと良い(KKRは全電子計算法なのでMCDなどの計算も可能である)。
http://pmt.sakura.ne.jp/wiki/index.php?title=Main_Page での Doument for AkaiKKR
http://kkr.phys.sci.osaka-u.ac.jp/jp/document.cgi, http://www.dyn.ap.eng.osaka-u.ac.jp/CMD/CMD4/CMD4_TEXT/CMD4_MACHIKANEYAMA2000_1.pdf
http://olymp.cup.uni-muenchen.de/index.php?option=com_remository&Itemid=20&func=startdown&id=51&lang=en
http://repo.lib.hosei.ac.jp/bitstream/10114/3486/1/07R6108_%E6%B5%A6%E5%B1%B1%E6%81%92%E5%A4%AA%E9%83%8E.pdf
H. 笠井秀明ら著、「計算機マテリアルデザイン入門」、大阪大学出版会 はかなりレベルが高いががんばって読み込んで欲しい。
I. 白井先生 http://www.cmp.sanken.osaka-u.ac.jp/~koun/ にある講義資料を読んでみる。
・ 密度汎関数理論: http://www.cmp.sanken.osaka-u.ac.jp/~koun/Lecs/dft.pdf
http://www.cmp.sanken.osaka-u.ac.jp/~koun/Lecs/var_main.pdf
・ 講義資料: http://www.cmp.sanken.osaka-u.ac.jp/~koun/Lecs/Material_Design13.pdf
J. 赤井久純ら著、「マテリアルデザインへの応用」、シュプリンがー・ジャパン はかなりレベルが高いががんばって読み込んで欲しい。
K. http://rare-events.org/tsuneda/Kyoto-text.pdf (必要な情報が纏められている)
L.  上記で理解が出来なかった場合は『物質の電子状態 上』と『物質の電子状態 下』や『Electronic Structure - Basic Theory and Practical Methods』(Richard M. Martin Cambridge University Press)にすがるしかない。しかし、論文をまとめたような本なので、私はほとんど理解できなかった。

その他:
全体の理解に役立つ: http://ci.nii.ac.jp/els/110007161500.pdf?id=ART0009115362&type=pdf&lang=jp&host=cinii&order_no=&ppv_type=0&lang_sw=&no=1378024711&cp=
個別の計算問題: https://www.jstage.jst.go.jp/browse/jsssj/28/3/_contents/-char/ja/
https://www.jstage.jst.go.jp/article/jsssj/28/3/28_3_129/_pdf
https://www.jstage.jst.go.jp/article/jsssj/28/3/28_3_135/_pdf
https://www.jstage.jst.go.jp/article/jsssj/28/3/28_3_144/_pdf
https://www.jstage.jst.go.jp/article/jsssj/28/3/28_3_150/_pdf
https://www.jstage.jst.go.jp/article/jsssj/28/3/28_3_160/_pdf
-----------------------------------------------------------------------
■ 圧力による相転移
 このような議論で基本的に重要な量として、体積Ωと温度Tを独立変数とする自由エネルギーF(Ω,T)=E(Ω,T) - TS(Ω,T), あるいは圧力Pと温度Tを独立変数とするGibbsの自由エネルギーG(P,T)=H(P,T) - TS(P,T)がある。
 エンタルピーHは、
 H = E + PΩ
で与えられる。温度T=0では、圧力一定で構造が安定な条件はエンタルピーが最小であることである。
 E(Ω)を計算し、2つの相それぞれのE(Ω)曲線の共通接線をGibbsの方法で引けば、その傾きがこの2つの相の間の相転移を引き起こす圧力を与える。
R.M.マーチンら「物質の電子状態(上)」丸善
-----------------------------------------------------------------------
■ 弾性:応力とひずみの関係
 例えば、縦軸を応力/ひずみ=t11/η11、横軸をひずみη11とした図で、線の傾きから弾性定数が求められる。非線形な変化分は非線形定数を決める。
 非線形ひずみに対してはLagrangeの応力ひずみを使えば線形の普通の表式に還元されるので、最も便利である。
 線形の弾性定数は多くの材料について計算されていて、その値は一般に実験値と非常によく一致しており、ずれは大体5-10%である。非線形定数については測定がずっと難しくなるので、理論値は予測値に過ぎない場合が多い。
R.M.マーチンら「物質の電子状態(上)」丸善
-----------------------------------------------------------------------
■ フォノンと変位型相転移
 理論と実験の周波数の一致の程度は通常は約5%以内である。

◇ 凍結フォノン(フローズンフォノン)法
 原子の位置の関数として全エネルギーを直接計算する方法。これは「凍結フォノン」法と呼ばれることも多い。
「凍結フォノン」という言葉は、核が位置{RI}に凍結されているという仮定の下で全エネルギーや力を計算する直接的な方法を意味している。
 この方法は他の問題の場合とまったく同じ計算手法が使えるという大きい利点がある。例えば、同じプログラムで(ほんの少し異なった入力をするだけで)フォノンの分散曲線、表面と界面の構造などその他の多くの特性を計算できる。

◇ 応答関数法(またはDFPT)
 エネルギーの任意の次数の導関数を陽に計算する。これは「応答関数」または「Green関数」法と呼ばれている。
 応答関数を使った方法とは、平衡位置からの変位のべきで展開することにより力の定数を計算するという方法を意味する。この方法の大きい利点は、実験で測定可能であり、1960年代に定式化された応答関数の理論に基づいていることにある。最近の発展により式の表現が計算に非常に使いやすい形になり、現在ではフォノンの分散曲線の計算はごく普通にできるようになった。
 遷移金属については実験との一致はそれほど良くない。同様の結果が多くの材料について見られるが、実験の周波数に関して5%以内で一致している。
 応答関数法は、誘電関数、有効電荷、電子-フォノン行列要素、等々の計算には特に効率が良い。
 全ての調和振動子の力の定数、弾性定数、等々はエネルギーの2次微分のみを含んでおり、2次の摂動論を使えばこれら全てが求められるからである。さらに、このことは、多くの場合、全エネルギーそのものよりも小さい変化の方がより正確に計算できるという事実にも基づいている。
R.M.マーチンら「物質の電子状態(上)」丸善
-----------------------------------------------------------------------
■ 電気伝導度と光学的性質
 誘電関数と電気伝導度は凝縮物質物理では最も重要な応答関数である。なぜならこの2つは物質の光学的性質、電気伝導、および多くの技術的応用を決めるものだからである。また、光学スペクトルは電子励起を研究する際に最も広く使われている測定手段であろう。
 分極性あるいは導電性の物質が存在する場合のMaxellの方程式の現象論的定式化は周波数に依存する複素誘電関数または複素電気伝導度を使って行われる。半導体GaAsの光学スペクトルでは、厳密交換であるEXXや電子-正孔相互作用を取り入れた「GW」準粒子で実験とよい一致を示す。
 さらに、2粒子を扱うBethe-Salpeter方程式(BSE)を解いて電子-正孔相互作用を取り込んだものは、全体としてずっとよく実験と一致するようになっている。
 電子-正孔相互作用の効果は、GaF2のような広いバンドギャップを持つ物質ではさらに大きくなる。この場合には定性的な変化が生じ、吸収の大部分は束縛された励起子状態の方に移る。同様の結果がLiFや他の広いバンドギャップを持つ絶縁体についても見られている。
 電子数が変わらない場合の励起スペクトルの計算には、時間依存密度汎関数(TDDFT)と呼ばれる方法がある。
 これは時間依存シュレーディンガー方程式に従って、原理的にはn(t)の厳密な解を与える。この方法は近似的交換相関汎関数と共に使われ、分子やクラスターのような小さな系の光学スペクトルや固体中の磁気励起に対してはかなりの成功を収めている。
R.M.マーチンら「物質の電子状態(上)」丸善
-----------------------------------------------------------------------
■ 結晶における分極の計算
 分極の計算例としては、有効電荷を考えるのが適切である。というのは、実験では光学フォノンの縦モードと横モードの振動数の分離から有効電荷を測定できるので、理論の結果と比較できるからである。
 例えば、いくつかのグループが、強誘電体への転移を示すBaTiO3などのペロブスカイトで異常有効電荷を算出している。このような物質には、いくつかの赤外(IR)活性モードがあり、すべてのIRモードはお互いに相互作用をして混ざっているので、実験から直接原子の個々の有効電荷を決めることは不可能である。
 これらの結果は、理論的モデルから格子動力学についての情報を使ってのみいろいろな原子からの寄与に分解できる。これとは逆に、非経験的計算では、原子の有効電荷と振動モードが決められ、固有モードの有効電荷が直接計算され、実験と比較できる。計算結果によれば、これらのIRモードは強く混ざっており、その結果最低の横光学(TO)モードは最高の縦光学(LO)モードと最も密接に結びついており、そのモードの有効電荷を非常に大きくしている。すなわち、強誘電相転移において結晶をソフトにするフォノンである。ペロブスカイト構造のBサイトの原子と、その原子に向かう線に沿って動いているO原子の異常に大きい有効電荷は共有結合の結果であると解釈されている。
 線形応答法を使っても、本質的には同じ結果が得られている。しかし、Berry位相法は強誘電状態でも発現する有限の分極に直接適用できるという利点がある。
 線形と非線形の感受率の計算例では、Wannier関数とBloch軌道間の変換及び[2n+1]定理を使って行われる。
※[2n+1]定理とは、0からn次までのすべての次数に対する波動関数が分かっていれば2n+1次までのエネルギーが求められるというものである。
R. M. マーチンら「物質の電子状態(下)」丸善
-----------------------------------------------------------------------
■ 自発分極
自発分極は反転対称の中心がない結晶ではどの結晶にでも起きる。
R. M. マーチンら「物質の電子状態(下)」丸善
-----------------------------------------------------------------------
■ 実空間での計算法
 密度の計算において既に指摘したように、多くの演算は実空間で行った方が簡単である。例えば、もしψ(R)がグリッド上で具体的に分かっていれば、n(r)=Σ|ψ(R)|^2の計算には平面波法で必要であったFFTなどは必要がない。HartreeポテンシャルはFFT法でも実空間多重グリッド法でも求めることができ、後者は高度に最適化されている。ポアソン方程式とシュレーディンガー方程式の解を求めるには実空間での方法は有限な系に対しては特に有利である。有限な系では波動関数は境界の外では消え、クーロンポテンシャルは一般には周期的境界条件には従わないからである。
 実空間での計算によるKohn-Sham方程式の解はフーリエ変換が不必要であるということ以外は平面波の場合と同じ形式になる。
R.M.マーチンら「物質の電子状態(下)」丸善
-----------------------------------------------------------------------
■ DFPT(勉強不足なので誤っているところがあるかもしれません)
 調和振動子近似の範囲内では、周波数ωの振動モードは変位
 ui(t) = Ri(t) - R0i(t) ≡ ui*exp(iωt)
によって記述でき、それを使えば下記の古典的連立運動方程式
 Mi*∂2Ri/∂t2 = Fi(R) = -∂E(R)/∂Ri
に対して、
 -ω^2*Mi*ui = -Σ(∂E(R)/(∂Ri∂Rj))
となる。全振動状態に対する解は、各々が周波数ωを持つ独立な振動子の集合であり、古典的方程式
 det|(∂E(R)/(∂Ri∂Rj))/√(MiMj)-ω^2|=0
によって決められる。ここで、質量MiとMjに対する依存性は対称的な形で入っている。

 ∂E(R)/(∂Ri∂Rj) = ∂Eii(R)/(∂Ri∂Rj) + ∫((∂2Vext(r)/(∂Ri∂Rj))*n(r))dr + ∫((n(r)/∂Ri)*(∂Vext(r)/∂Rj))dr
右辺の最初の2項は無摂動の密度のみを含んでいる。しかし、最後の項では(n(r)/∂Ri)の情報が必要である。

 任意の周期的外部摂動ΔVext(r)に対する応答の計算用のDFPTアルゴリズムは、(n(r)/∂Ri)で用いられるΔψを得るために1次の摂動論より得られる式、さらに、その式の中で用いられるΔV、(n(r)/∂Ri)の3つからなる連立方程式を解くことである。
 DFPT方程式は、与えられた波数ベクトル、例えば、フォノンの波数ベクトルで、ブロッホの定理に従う原子変位(の固有ベクトル)を持つ摂動のある結晶の場合には、驚くほど簡単になる。

 これらの計算は、効率の良い平面波法では標準的な高速フーリエ変換を使って行うことができる。
R.M.マーチンら「物質の電子状態(上)」丸善
-----------------------------------------------------------------------
■ 開殻原子:非球対称ポテンシャル
「開殻」という用語はスピンが対をなさず、あるいは、与えられた角運動量lの状態m=-l, ...,lが完全には詰まっていない場合を指す。従って分類をするのであれば与えられた全角運動量J={j,mj}を持つ多重項で行うのが良い。ここでJは空間(L={l,mj})とスピン(S={s,ms})の一次結合である。一般には複数の行列式でできた波動関数を扱わねばならない。たとえ外部ポテンシャル(原子では核のポテンシャル)が球対称であっても、有効独立粒子ポテンシャルVeff(r)は軌道の占有状態に依存するので球対称ではない。簡単化できることといえば量子化軸を選んでポテンシャルVeff(r)=Veff(r,Θ)が軸対称性を持つようにすることぐらいである。幸いなことに、Slaterの方法は対称性を使って適切な多重項を選ぶことにより、必要な計算すべてを動径方向だけの計算に簡約する仕方を示している。このような対処の仕方に一般的な法則はないが、いろいろな場合を集めたものが文献の付録に載っている。原子の中ではVeffは完全に動径のみの関数Veff(r)が、あるいは円筒座標の関数Veff(r,Θ)の場合のみを考えればよい。開殻問題に対しては、以下のような精度の異なるいろいろな近似がある。
・制限付き:問題を球対称(球対称ではない項にについては角度平均をとってVeff(r)とする)、
かつスピンによらない(スピン状態について平均をとり、軌道が各スピンについて同じになるようにする)として扱う。
これは閉殻、スピン=0の系、に対しては正しい形式であり、注意深く扱えば開殻に対する近似と見なすこともできる。
・スピン非制限:球対称であるがポテンシャルと軌道はスピンに依存するとして扱う。これは殻が半分占有され、
スピンが最大値をとり、各スピンが別々に閉殻構造をとっているときに正しい式である。
・非制限:問題をそのまま扱う。そこでは全mlとmsのみが良い量子数である。この場合にはSlaterの方法を使って問題を簡単化することができる。

◇開殻構造に対する方程式
 完全に非制限の場合には、電子-電子相互作用項のために一段と複雑になる。軸を選んで密度n(r,Θ)とポテンシャルVxc(r,Θ)が軸対称となるようにする。Kohn-Sham法においては、波動関数は球面調和関数で展開される。また、Vxc(r)とn(r)の関数が非線形であるため、Vxc(r,Θ)の球面調和関数による展開においてLにはカットオフの最大値が存在しないので、Vxc(r,Θ)の角度積分は数値的に行わなければならない。さらに、Coulombポテンシャルは多重極モーメントを持ち、Poisson方程式の解は球対称の場合のようには簡単ではない。開殻構造の場合には、Hartree-Fock方程式の方が、交換項が球面調和関数の有限な和で展開できるため、実際には簡単である。電子-電子相互作用の行列要素を計算するためには球面調和関数を使ったよく知られた展開を使うことができ、それはr1とr2を含む項に因数分解される。
 角度の積分は解析的に行うことができ、Gaunt係数が出てくる。これは文献に具体的な式が与えられている。
 動径方向のHartree-Fock方程式は動径関数に関するエネルギーの汎関数微分から出てくる。
 Hartreeポテンシャルのlimi成分をGaunt係数とGlebsch-Gordan係数の関係式を使って表すことができる。同様に、状態ni,li,mi,σiに作用する交換ポテンシャルも表すことができる。

◇相対論的Dirac方程式とスピン軌道相互作用
 重い原子では相対論的効果を無視することはできない。幸い、この効果は内殻の深いところで生じるので、相対論的方程式は球対称と仮定して解くことができる。その結果は分子や固体において本質的な修正なしに使うことができる。固体や分子に関する実際の計算では、相対論的効果は補強法(補強法ではこの計算は1個の原子の場合と等価である)の中に直接に、
あるいは別の方法、例えば擬ポテンシャルの中に間接的に、取り込むことができる。相対論的効果は原子の相対論の計算結果を使って擬ポテンシャルを作ることにより、擬ポテンシャルの中に組み込むことができる。そして、非相対論のシュレーディンガー方程式の中でこの擬ポテンシャルを使えば、相対論的効果を含んだ価電子状態を求めることができる。
R. M. マーチンら「物質の電子状態(上)」丸善
-----------------------------------------------------------------------
■ Chebyshev多項式
Taylor級数は、変数のべきによる展開である。
 f(x) → c0 + c1*x + c2*x2 + ... + cN*xM
フェルミ演算子の展開の式で必要とされているような演算子の展開を使い利点は、引き続いて出てくる各項がxn+1=x*xnを順に行えば簡単に得られることである。しかし、高次のべきでは不安定性の問題が出てくることがあり、xが大きくなるにつれ収束は悪くなる。一方、第1種のChebyshev多項式Tn(x)は、領域[-1,+1]で直交するように定義されているため、この領域上では、どんな関数でもTn(x)の線形結合として一意に展開することができる。さらに、この展開をすると、その全領域で関数f(x)の近似展開の誤差の絶対値の最大値が最小になるという性質があり、さらにその多項式は再帰的に計算することができるという特性がある。多項式は、最初の2項と他の項を表す漸化式
 T0(x)=1, T1(x)=x, Tn+1(x)=2*x*Tn(x)-Tn-1(x)
を定義することで表現できる。最終的な展開式は
 f(x) → c0/2 + ∑(Cn*Tn(x))
と書ける。最初の数項を導出し、その直交性を示すのは簡単な演習問題となる。
 Chebyshev展開法は、もし基底の集合が小さければ、例えばNbasisがNより2桁大きい程度の強束縛モデルであれば、非常に効率の良い方法である。しかし、限界のないスペクトルに対してはうまくいかない。というのは、多項式展開は、スペクトルの幅が広ければそれだけ高次までとらなければならないからである。典型的な平面波計算では、エネルギー領域が広く、かつNbasis≫Nなので、高次のべきが必要になる。

◇密度行列の逆べき展開
 高エネルギーで適切に収束する演算子を使い、限界のないスペクトルを扱う方法がいくつかある。その1つは逆べき法でありそれは本質的にはGreen関数法であり、リカージョン法と密接に関係している。まずFermi関数を次式のように展開する。
 F[H] = 1/(1+exp[β*[H-μ]]) → ∑(ωi/[H-zi])
経路積分について良く知られた関係式を使えば、実軸上の極についての和は、複素平面上の極を囲む積分に変換できる。複素平面上のziに極を持つ各項に対して、逆演算子Gi(zi)は、1次方程式(H-zi)Gi(zi)=Iを基底を使って解くことにより求められる。実空間の基底関数を使えば、演算子Gi(zi)は大きい複素数ziに対してはより局在化しており、このことは計算には都合が良い。その最大の範囲は経路が実軸と公差するところである。絶縁体における経路積分を精確に記述するために必要な極の数はMpole∝(μ-Emin)/Egapである。金属ではギャップがなく、実軸付近の積分の正確な計算には間隔がTに比例するより密な極の分布が必要である。この方法の利点は、べき展開とは異なって、高いエネルギースペクトルとは無関係に常に収束することである。多重散乱Green関数計算のなどに使われている。
※ Chebyshev多項式展開は、全エネルギーや状態密度のような示量性の量の統計的な見積もりを求めるためにランダムなベクトルと共に使うこともできる。これは「最大エントロピー」法あるいはこれに関連した方法である。この方法は対角和をとることができないような非常に大きい行列に対しては有効である。
※ Fermi分布関数を、超幾何関数を使って収束の早い連分数展開で表すことが可能である。
R. M. マーチンら「物質の電子状態(上)」丸善
-----------------------------------------------------------------------
■ 数値計算法
 下記の方法は、波動関数の反復法による改善(反復的対角化)、Kohn-Sham自己無撞着ループにおける電荷密度の更新、および構造緩和における原子の変位の計算に主として使われている。問題の規模や性質は多様であるので、それぞれに応じてより適した方法がある。

◇ 数値積分とNumerov法
 物理学において特別重要な役割を演じる2次微分方程式の例であり、この他の例としてはKohn-Sham方程式やHartree-Fock方程式の自己無撞着解を求める際に解かねばならないPoisson方程式などがある。これらの方程式はグリッド上で離散化し、2次微分に対する数値的近似を使って積分できる。効率の良い方法は、原点からは外向きに、無限大からは内向きに、接合点まで方程式を積分するNumerovのアルゴリズムである。その解は波動関数とその微分がある選ばれた半径Rcで一致するという条件によって与えられる。波動関数の振幅は一致するように調整できるので(全体の振幅は規格化によって決められる)、実際にはφの対数微分を表す比が一致すればよい。
 誤差がh^6オーダーの方法がNumerov法として知られている。
 Numerov法の基礎にある考え方はどの次元にも拡張可能であり、そのグリッド点の作る図形は「Mehrstellen」と呼ばれている。要点は、微分方程式そのものを使い、もとの差分方程式より高い次数で成り立つ運動エネルギーとポテンシャルエネルギーに対する表式を見つけることである。

◇ 最急降下(SD)法
 最急降下方向を目的とする関数Fに対する微分で得るために、SD法では各ステップの始点では直線探索の方向は関数の等高線と直交する。この方法は「(Zeno)ゼノンのパラドックス」の具体例であり、厳密に最小値に到達することはない。
 SD法は、もし関数Fが変数ごとに異なる依存性を示し、その結果、最小点の周囲で細長く長い谷を作るような場合には、特に悪い方法である。

◇ 共役勾配法
 驚くべきことに思えるかもしれないが、いつも「下り坂」の最急降下方向に従って進むよりも、最小点に速く到達する方法がある。最初のステップの後では、その点での勾配だけでなく、その前の点での値と勾配も得ている。この付加情報を使えば、直線探索によってより低いエネルギーに導くためのより適切な方向を選ぶことができる。事実、N次元の2次汎関数に対しては、共役勾配(CG)法を使えばNステップで最小値に到達することが保証される。
 CG法はもっと複雑な汎関数(例えば、Kohn-Sham汎関数など)にも適用でき、汎関数は最小値の近傍で2次形式であるので、多くの利点が期待される。2次形式ではない問題にCG法を使うには、共役方向を定義し、該当の非線形汎関数に直線探索を行うことである。新しい勾配が現在の方向に直交するためには直線探索が済んでいなければならず、その結果、汎関数は各方向について順に最小化されることになるので、CGアルゴリズムにとっては直線探索は本質的な事柄である。

◇ 準Newton-Raphson法
 方程式 F(x)=x を解く問題を考える。ここでxは多次元のベクトルである。例えば、出力密度nout(r)(これは入力密度nin(r)の関数である)が入力密度nin(r)に等しいようなKohn-Sham方程式の解を求める問題である。
 もし密度がM個の関数の集合で、nin(r)=ΣkMxkhk(r)、ただし、x={xk}のように展開されていれば、この問題は厳密にF(x)=xの形となる。これは、残差のノルム|R[x]|の最小化問題となり、ここで、
 R[x]≡F(x)-x
である。もしRxの線形関数であり、そのJacobi行列
 J≡∂R/∂x
がわかっている領域を扱うのであれば、残差を最小化する準Newton-Raphson法が使用できる。
 ステップiでのxiを扱うと、次の反復でRi+1=0となる値は、
 xi+1=xi-J-1Ri
である。
 問題は、一般的にはJacobi行列は未知であり(あるいは逆をとることが困難であり)、関数空間(Krylovの部分空間)において解に至るまで反復するという他の方法に頼らねばならないことである。

◇ PulayのDIIS法
 DIIS(discrete inversion in the iterative subspace)法の根底にある考え方は、それまでに生成した全てのベクトルの最適な結合を使って、すなわち、全Krylov部分空間を使って、どのステップiにおいても残差を最小にすることである。それ以前の全ステップの結果から新しいベクトルの係数を求め、それを残差の式に代入すれば、各ステップiでの新しい最適ベクトルが得られる。極度に大規模な問題、例えばシュレーディンガー方程式の多くの固有ベクトルを扱うなど、に対してはベクトルの多数の集合を保存するのは不可能である。しかし、密度混合では、特に密度の少数の厄介な成分についてだけであれば、数段階前の密度は保存できる。
 KresseとFurthmullerは、上記のPulayのDIIS法が修正Broyden法に密接に関係したJacobi行列の更新と等価であることを示した。さら、van LentheとPulayはDavidsonとDIISのアルゴリズムを各ステップでたった3個のベクトルを使って実行できることを示した。

◇ BroydenのJacobi行列更新法
 Broyden法は、反復法を使う過程で逆Jacobi行列を順次生成する方法である。修正Broyden法は1度に2つの状態のみを陽に扱うということを除けばDIIS法と似ている。この方法は広く使われている。平面波法では密度のほんの少しの厄介な成分だけをこの方法で扱えばよい。しかし、多くの計算で必要となる大きいグリッド上の電荷密度の更新などの、Jacobi行列全部の保存ができないときには役に立たない。
 修正Broyden法は、VanderbiltとLouieにより提案され、Johnsonにより記憶容量が少なくて済むSrivistavaの方法の利点を取り入れるように改良された。その考え方は、引き続く2つのステップが同じ結果を出すという条件はあまりに厳しすぎるというもので、改良されたアルゴリズムでは、それ以前の反復時の情報を取り入れることができる。

◇伝統的なハウスホルダー3重対角化法
 伝統的なハウスホルダー3重対角化法を使ってエルミート行列のスペクトル特性を直接決定するためには、N^3でスケールされる計算量がいる。
R. M. マーチンら「物質の電子状態(下)」丸善
-----------------------------------------------------------------------
■ 電子状態の反復計算法
 エネルギーの最小化対残差の最小化、あるいは単一ベクトルの更新対部分空間全体の反復、等々。これらの方法は1つの共通の枠組みに入れることができ、そこでの重要な特性は:
・行列の対角化を、反復(Krylov)部分空間におけうr波動関数ψiの反復方程式で置き換える。
・波動関数(多分さらに以前の波動関数)と勾配(dE/dψ=Hψ,このステップのアルゴリズムは方法により異なる)を使い新しい波動関数を求める。
・平面波の場合は、密行列の乗算を高速フーリエ変換(FFT)で置き換える。

◇ なぜ反復法を使うのか?
 電子状態の計算は基底関数の形によって2つのグループに分けられる。LCAOやLMTOなどの方法は、サイズNbという最小基底の構築を目標にして構築されている。基底の構築に努力が払われ、そのれは問題の種類に応じて高度に最適化されることになる。非常に大きい系を除いては、ハミルトニアンは小さく、サイズNb*Nbの密行列で表され、伝統的な密行列対角化という手法に適している。その場合の計算労力は、求める固有ベクトルの数をNeとすれば、Nb^3またはNe*Nb^2とスケールされる。 一方、平面波のような一般手kな基底とグリッドを使う方法では、しばしば欲しい固有状態の数よりずっと多数(Ne>>Nb)の基底関数が絡んでくる。ハミルトニアンは非常に簡単に作成され、疎にすることができる。すなわち、
行列要素は主としてゼロであり、ゼロ以外の要素のみを計算したり保存したりすればよい。小規模の問題を除けば、反復法が非常に効率が良く、この場合にはNb*Nbのハミルトニアンは明確に作られることはなく、計算量は、Ne^2*NbまたはNe^2*Nb*ln(Nb)でスケールされる。このような手法は、平面波法、実空間法、有限差分法、有限要素法、多重グリッド法、ウェーブレット法に適用されたとき最も成功を収めてきた。
 APW, LAPW, およびPAW法はある程度この中間にあり、両者の利点をとることができるかもしれない。

◇ 前処理
「前処理」についての基本的な考え方は、関数の変数依存性を修正し、より等方的にすることである。つまり、異なる変数に関する曲率をより均一にすることであり、それは収束性改善の考え方そのものである。電子状態計算で遭遇する問題には、最初の定式化が非常に悪い条件になってしまっていることがしばしばある。しかし、物理的な観点から、そのような条件の改善の仕方を見つけるのは容易である。一般に、どの形式を選ぶかは問題によって異なる。

◇ 反復法での(Krylovの)部分空間
 反復法はある演算子Aを繰り返し適用し新しいベクトルを作ることに基づいている。
(省略)
 各繰り返しにおいて具体的に取り扱われるKrylov部分空間の範囲は、方法によって異なっている。Jacobiアルゴリズムのような簡単な緩和法では新しい近似固有ベクトルを前ステップの近似固有ベクトルをのみを使って求めている。これは最急降下最小化法とにている。この他のLanczos法, Davidson法, RMM-DIIS法などでは、それまでの繰り返しで発生した全部分空間を使っている。一般には非常に大きい問題では全Krylov部分空間を保持しておくのは膨大な負荷がかかる。しかし広く使われている
 Lanczos法や共役勾配(CG)最小化法は全部分空間法であり、近似固有ベクトルは直前の2つの近似固有ベクトルをのみを使って求めているとしても、それ以前の全ベクトルに直交する(あるいは共役な)新しいベクトルを生成できる。従って、これらの方法は、繰り返しの各ステップで必要事項が多少増してはいるが、簡単な緩和法よりはるかに強力である。 Davidson法は、もともとは全部分空間を保持する必要があったが、CG法を使って3個のベクトルを必要とする形式に変えることもできる。

◇ Lanczos法
 Lanczos法は、現代の計算機を使って固有値問題を解く最初の反復法の1つであった。それは物理的な解釈や類似性を引き出す道具としては非常に単純であり、驚くほど強力である。そのアルゴリズムは自動的に直交基底(Krylovの、あるいは反復部分空間)を生成し、そこでは与えられた演算子Aは3重対角行列になる。(電子状態の問題ではA=Hであり、Hはハミルトニアンであることが多い)。大きい行列の、小さい(あるいは大きい)方からいくつかの固有エネルギーとそれらの固有ベクトルの生成にはとりわけ強力である。
 最も簡単なバージョンでは、欲しい固有ベクトルの数が増加すると、数値的な丸め誤差により発生する偽りの解という「Lanczos病」に悩まされるが、これは何回かの反復ステップの後で、直交化を行えば容易に制御できる。さらに、これは連分数として定式化でき、そこからスペクトル分布のモーメントを求める強力な方法が導かれる。
 各ベクトルψnは帰納法で示されるように、他の全てのベクトルに直交している。
 ステップMまで進めば、3重対角行列が得られ、計算時間がO(M)で、次数Mの3重対角行列の固有値を求めるための標準的なルーチンがある。そのアイディアは簡単で、多項式の再帰的集合で定義された関数の根の形で表すことができる。固有ベクトルは固有値が分かれば逆反復によって簡単に求めることができる。もし最初の基底にN個のベクトルがあれば、この3重対角行列形でハミルトニアンを生成するためにはNステップが必要である。
 この方法で重要な点は、スペクトルの最小値と最大値の近くの固有値と固有ベクトルは、行列全体を生成しなくても非常に正確に決められることである。なぜ最大値と最小値のスペクトルが一番最初に生成されるかは、次のことに注意すればわかる。もし固有ベクトルの線形結合である試行ベクトルψから始めたとすれば、3重対角化のための各ステップはHを演算することであり、それは試行状態の平均エネルギーからエネルギー的に最も離れた固有ベクトルの重みを増す操作である。平均エネルギーはスペクトルの中央あたりになければならないので、両端の状態は最も効率良くはじき出される。数百万の状態が係わる典型的な場合には、アルゴリズムの数ステップで、十分に高い精度で基底状態を求めることができる。
 Lanczosアルゴリズムが使える方法には多くのヴァリエーションがある。例えば、再帰的関係式はサイズM*Mの行列を生成するために使用でき、それを対角化すれば欲しい状態、例えば最低の固有状態などを求めることができる。このようにして求められた状態は、次の反復のための初期ベクトルとして使用でき、望みの精度が達成されるまでこれを繰り返す。この方法はとりわけ安定であり、いくつかの状態を求めるためにまとめた形で使えば便利である。この手法はM=2のような小さいMに対しても使うことができ、多くのベクトルを保存する必要がないので、基底状態を求めるためには特に効率のよいアルゴリズムとなる。
 多体問題のように極端に大きい問題では、少数の数値的な固有値を求めるためにLanczos法が広く使われており、この場合にはそれは「厳密な対角化」と同意義である。
 大規模電子状態計算における最も完璧な応用例は、WangとZungerが示したもので、彼らは偽状態の発生を避ける手段を組み込んでいる。
 Arnoldi法は、Gram-Schmidt法によって各ベクトルψnをそれ以前のすべてのベクトルに明確に直交させるLanczos法の変形である。この方法は不安定性を除き、望みの通りの多数の固有ベクトルの計算を可能にする。この方法は非エルミート行列にも適用できる。
 Lanczos法の3重対角行列形式の見事な結果は、動的相関関数のスペクトルの複素周波数zの関数として表した連分数表式である。これからリカージョン法が導かれる。

◇ Davidson法
 Davidsonは、現在、電子状態問題に広く適用されている方法を考案した。その変形版は沢山ある。要点は、Davidson法はLanczos法に密接に関連しているが、演算子が対角優位な問題に対してさらに効率を高めていることである。対角優位な問題は、例えば平面波アルゴリズムにおけるように、電子状態問題ではよくあることである。
 表式は摂動理論によく似ており、ハミルトニアンの対角部が優位であれば急速に収束する反復法に向いていることがわかる。

◇ RMM-DIIS (部分空間での残差最小化法)
 これまでに述べてきた方法(および以下で述べる最小化法)では、基底状態は大局的な最小値であるから、まったく問題なく最低状態に収束する。
 より高い状態を求めるためには、Lanczos法のように暗に、または明確な直交化によって、直交性を確かなものにしなければならない。Pulayが提案した残差最小化法(RMM)ではこの条件を除き、エネルギーではなく「残差ベクトル」のノルムを最小化するので、スペクトル中で試行固有値εに最も近い固有値を持つ状態に収束する。Pulay法はそれまでの反復で生成された全Krylov反復空間で残差を最小化するので「反復部分空間で直接に逆をとる残差最小化法(residual minimization method by direct inversion inthe iterative subspace)」、略して、RMM-DIISとして知られている。(省略)
 時間のかかるステップはHψの演算であり、一般にNbを基底のサイズとすれば、各固有ベクトルψにつきO(Nb^2)の演算を必要とする行列演算である。しかし疎行列の演算に対しては、これは大きい基底の場合にはO(Nb)まで減り、FFTを使えばO(Nb*ln(Nb))まで減る。
 実際には、大きい問題では、多くのベクトルを記憶することはないので、この反復を再度進める前の数ステップのnに対する小さな行列のみを対角化している。電子状態計算におけるこの方法の運用は、WoodとZungerが初めて行い、
続いて多くの研究者が種々の修正を加えて使用した。
 DIIS法は部分空間サイズの行列の全要素を求めることが必要であるが、行列が小さく、その空間を張る全ベクトルが保存可能であれば効率が良い。これは2つの領域でKohn-Sham方程式を解くことで達成される。もし必要な固有状態の数が小さければ、RMM-DIISを使って即座に全状態が求められる。もし必要な固有状態の数が大きければ、問題をいくつかのエネルギー領域に分割でき、選んだエネルギーに最も近い固有値を持つ少数の状態を小さい行列方程式を解いて求めることができる。この場合には、固有状態を見逃したり、重複して数えたりしないように注意しなければならない。
 共役勾配法に比べてこの方法が非常に優れているところは、他のベクトルとの直交性を明確に要請することなく、スペクトルの中央であってもどのような固有値でも求められることである。
R. M. マーチンら「物質の電子状態(下)」丸善
-----------------------------------------------------------------------
■ FFT
 逆格子空間で計算される第一原理計算コードにおいて、FFTが必要になるときとはどのようなときであろうか?
 それは、逆格子空間での波動関数をフーリエ変換で実空間にし、それに有効ポテンシャルを掛け、その後に実空間から逆格子空間に戻すために逆フーリエ変換が用いられる。
 こうすることによって、交換相関ポテンシャルで使われる電子密度を、実空間の電子密度Σ|ψi|^2で取り扱うことで計算を容易にしている。
 逆格子空間でそれを行おうとすると無限大の集合を扱うことになり、厳密な取り扱いが出来なくなる。
R. M. マーチンら「物質の電子状態(下)」丸善
-----------------------------------------------------------------------
■ 擬ポテンシャル

◇ 転用可能性と硬さ
「硬さ」という言葉には2通りの意味がある。
 1つ目の意味は実空間での変化の尺度があり、それはFourier空間におけるポテンシャルの広がりで定量化される。一般に「硬い」ポテンシャルは局在した硬いイオン芯の性質を表しており、1つの物質から他の物質へ転用しやすい。ポテンシャルを「柔らかく」(すなわち、滑らかに)しようとすると転用性を悪くする傾向がある。しかし、正確で転用性があり、しかもFourier空間であまり広がらない、つまり、「最適化された」擬ポテンシャルを作るべく、現在相当な努力が払われている。
 2つ目の意味は、系の環境の変化に対する応答を、擬ポテンシャルを使って適切に記述する能力の尺度である。ノルム保存は原子の電子状態がエネルギーの変化に対して正しい1次の微分を持つことを保証している。「硬さ」のこのひょうな意味はポテンシャルの変化に対する応答の忠実さの尺度である。ポテンシャルは球対称の摂動(電荷、状態、動径ポテンシャルの変化)に対して通常の球対称原子のコードを使ってテストすることができる。GoedeckerとMaschkeは内殻領域の電荷密度の応答を使って洞察力豊かな解析を行っている。これは密度が密度汎関数論の中心的な量であり、積分された密度はノルム保存の条件と深く係わっているので、当を得たことである。球対称ではない摂動を使ったテストでもまた、特に電界中の分極率に関する性能を確認している。

◇ 最良の擬ポテンシャル
 どのような元素に対しても唯一の最良の擬ポテンシャルというものは存在しない。あるのは沢山の「最良の選択」であって、それぞれは擬ポテンシャルを特定の使用に最適化したものである。
 一般には、以下のような2つの競合する要素がある。
・精度と転用性の良さを求める場合には、一般に小さな切断半径Rcとなり、「硬い」ポテンシャルを選ぶことになる。というのは、原子の近傍でできる限り波動関数を正確に記述したいからである。
・擬ポテンシャルの滑らかさを重視すると、一般に大きい切断半径Rcとなり、「柔らかい」ポテンシャルを選ぶことになる。というのは、波動関数をできる限り少数の基底関数(例えば平面波)で記述したいからです。

◇ ウルトラソフト擬ポテンシャル
 擬ポテンシャルの1つの目標は、できる限り滑らかでしかも正確な擬波動関数を作ることである。
 例えば、平面波の計算では、価電子関数がFourier成分で展開され、計算の負荷はその計算で必要なFourier成分の数のべきで増大する。従って滑らかさを最大にするということの1つの意味のある定義は、決められた精度で価電子の特性を記述するために必要なFourier空間の範囲を最小化するということである。
「ノルム保存」擬ポテンシャルは、通常「滑らかさ」をある程度犠牲にして精度という目的を果たしている。これまでとは違う「ウルトラソフト擬ポテンシャル」と呼ばれる方法は、正確な計算という目標を、滑らかな関数と各イオン芯の周りの補助関数を使って問題を再表現するという変換によって達成している。この補助関数は密度が激しく変化する部分を表す。

◇ 射影演算子補強波(PAW)法:全電子波動関数の使用
 射影演算子補強波(PAW)法は、電子構造問題の解を求める一般的な方法であり、OPW法を再定式化して、全エネルギー、力、応力の計算の最新技術に適応させたものである。「ウルトラソフト」擬ポテンシャル法と同様に、これは射影演算子と補助局在関数を導入している。PAW法はまた全エネルギーに対する補助関数を含む汎関数を定義し、一般的な固有値方程式の解を効率良く求めるアルゴリズムの進歩を取り入れている。しかし、PAW法は、一般的なOPW表示と類似の形式で、全電子波動関数を保持しているという違いがある。全電子波動関数は核の近くで急激に変化するので、すべての積分は全空間に広がる滑らかな関数の積分と補強平面波(APW)法と同様にマフィンティン球における動径方向の積分で求められる局在した量の寄与とを組み合わせたものとして計算される。

◇ OPW
 価電子に対する内殻の影響を取り入れている。局在関数を上手に選べば、価電子状態に対する基底関数を滑らかな部分と局在部分に分けることができる。結晶中では滑らかな関数は平面波で表現するのが都合が良い。
 内殻の状態は、分子や固体中でも原子におけるものと同じであると仮定できることも多く、これがOPW法の実際の計算における基礎になっている。
 擬ポテンシャル変換の最大の利点は、擬ポテンシャルの形式的な特性と、同じ散乱特性が、違うポテンシャルからでも得られるという両方の事実を利用できることである。従って、擬ポテンシャルの一義性がないことを利用して、もとのポテンシャルより滑らかでかつ弱く選ぶことができる。たとえポテンシャルの演算子が簡単な局所ポテンシャルより複雑であったとしても、それが弱くて滑らか(すなわち、より少ないFourier成分で表現できる)であれば、論理的にも計算する上からも大変便利である。特に、波動関数は内殻の波動関数と直交していなければならないので、自由電子とは程遠いものであるにも関わらず、多くの物質で価電子バンドがほとんど自由電子に近いバンドになっているという見かけ上の矛盾
を簡単に解くことができる。その理由は、バンドは弱いポテンシャルに関連する滑らかなほとんど自由電子に近い波動関するに対する永年方程式によって決まるためである。

◇ 全電子波動関数の再構成
 擬ポテンシャルを使った計算では、直接求められるのは擬波動関数のみである。しかし、重要な物理的性質、例えば、核磁気共鳴実験におけるKnightシフトや化学シフトなどを求めるためには全電子波動関数が必要である。これらのシフトは核の周囲の状況や価電子状態を知る非常に感度の高いプローブであるが、得られる情報は内殻状態の摂動の影響を大きく受ける。その他の実験、例えば内殻準位の光電子放出と光吸収の実験では、内殻状態が直接関わってくる。OPW法とPAW法では内殻の波動関数が求められる。ところで通常の擬ポテンシャル計算から内殻の波動関数を再現することは可能であろうか? 答えは、ある近似の範囲内であれば可である。その手続きはPAW変換に密接に関係している。滑らかな擬関数を作るどのような「非経験的手法に対しても、分子や固体において滑らかな擬関数が計算されていれば、全電子波動関数を再構成する具体的な方法を定式化することができる。Mauriとその共同研究者たちはこのような再構成法を使って核の化学シフトを計算している。

◇ 非局所Exc汎関数
 Exc汎関数が、Hartree-Fockや厳密交換(EXX)のように本質的に非局所的な場合には価電子からの寄与の除去は複雑になる。一般に、非局所効果は全半径に及ぶので、擬波動関数が内殻半径の外では、ともの全電子問題の波動関数と同じになるように保つポテンシャルを作るのは不可能である。

◇ 球対称境界条件
 固体のいくつかの様相は原子やイオンに異なる球対称の境界条件を課すモデルでよく表現できる。その結果、価電子波動関数は原子のときよりさらに核の近くに集中するという傾向が分かっている。
R. M. マーチンら「物質の電子状態(上)」丸善
------------------------------------------------------------------------------
アクセス数
ページビュー数