Transcript 密度汎関数法の基礎
第一原理電子状態計算の方法の開発
鳥取大工 小谷岳生(Takao Kotani)
July28-30, 2009@osaka-u
Chap.0.
Introduction
目標、現状と問題点を概観する。
~0.5コマ
Chap.1. 理論的な基礎
BeyondLDAを考えるための基礎知識(の一部)
~3.5コマ
Chap.2.GW近似、QSGW近似
~1.5コマ
Chap.3.一体問題の解法
Linearized augmentation methods. 現代の方法。PMT法
~1.5コマ
Chap.4.数値計算技術の実際
我々のプロジェクト「ecalj」の説明(googleする)。
Live CDを用いた実習。 木曜0.5コマ、金曜0.5コマ
1
Chap.0 Introduction
Intro.
• 第一原理電子状態計算とは何か?
(時間的制約で、本講義では原子核の位置を固定した話のみに限定)
以下の問題のユニバーサルな算法を開発すること
=1 、 v ext (r)は原子核から のポテンシャ ルの和
さらには、相対論的補正、ベクトルポテンシャルとの相互作用。
2
何のために?
「物理学」「物質科学」
Intro.
究極の目標:
*物質の特性の高精度な定量的計算
(物質に依存するパラメーターをもちいない)。
*物質を自由自在にあやつる
原料
物質
物性
プロ セス
測定
特性の改良、物質設計
注:現在できることは限定的。
まだかなり幼いレベル(米国の西部開拓時代ぐらい?)
*モデル計算でできることは原則的に第一原理計算でも可能であるはず
第一原理計算は(日本では偉い先生にすら)誤解されている。
少なくとも(数学でないなら)「明瞭な関連性」を確立していく必要あり。
モデル手法の限界:定量性。適用限界が理論内部で判断できない。
昔の理論:e.g.5つの内部パラメーターを含むモデルを、5つの実験で決める。
3
計算したい物理量
Intro.
• 基本量
– 独立粒子近似の描像を与える一体ポテンシャル
『 H H0 ( H H0 )と 分割する こ と 』
– 準粒子描像、独立粒子間の有効相互作用
– 全エネルギー
– 複雑な系(異種物質の接続、表面+分子など)
• 複合的な量
– 光学的、磁気的、電気的な線形応答。T=0 or 有限温度
– 構造に関するエネルギー、フォノン、熱力学的な量、ダイ
ナミクス。
– これらの組み合わせ、 非線形応答、超電導転移温度。モ
デルとの組み合わせ。
4
現状の概観
Intro.
• 「密度汎関数法(DF)の局所密度近似(LDA、
GGA)」を軸に発展してきた。状況(物質のタイプや物理
量)によってはかなりの予言力を持つ(?)。
「DF-LDA(GGA)は局所ポテンシャルをもつ一体問題での平均場近似。」
• 信頼性のない場合も多い(Publication biasに注意!)
– 自己相互作用の問題:水素、鏡像ポテンシャル
– 局在性の強い電子
(次ページで追加説明)
―方法論における最近の新しい進展の方向性―
*新しい物理量の計算法
*一体問題の解法(高速、信頼性、安定性)
*電子相関を取り入れていくこと(Beyond LDA).
* 低エネルギー励起を取り入れた計算(をしたい)。
magnon, phonon-GW
5
現在できてることは限定的
• 個人の能力だけでは不能、いろいろな技術の発明と
インテグレーションしていくことが必要。
「なぜトランジスタが発明できたか?」
「なぜ月まで行って帰ってくるロケットがつく
れたのか?」
工学も理学もない。究極の目標は明瞭。研究の
方法の方法から考えていかないといけない。
いままでの学問の枠を超えていく必要もある
。
6
「鳥取大学乾燥地研究センター」が世界的に著名。
鳥取大学教授で中国で~300万本の木を植えた男「
遠山正瑛」がいた。いまだに世界中から留学生がやっ
てくる。
ちょっとした手品みたいな「物性理論、実験」を組み
立ててみてもあまりおもしろくない。(悪い意味で)
趣味の研究は趣味でしかない。
注:ただ、いろんな研究者がいます。自分にとっての、「趣味
の研究」で大成果を上げるひともいます。
7
Intro.
LDAやGGAの結果
結論: 「DF+Kohn&Sham方程式」の枠組みに限界あり。よい非
8
局所的平均場ハミルトニアンH0を決める方法。
If,
Bonding
HOMO
Anti-Bonding, LUMO
これらを二乗した電子密度は同じ。それゆえ
bonding V (r ) bonding anti-bonding V (r ) anti-bonding .
Local potential shows no difference between LUMO and HOMO.
Or Non-local potential can add the energy difference
without deforming eigenfunctions.
This is a reason why local potential is problematic.
9
パッケージ化されていて利用者の数は急増している
Intro.
• Vasp,Gaussian
メジャー。利用が簡単らしい。ブラックボックスとして利用する。
• 諸派(Wien2K, Abinit, HiLAPW,AkaiKKR,QMAS…)
小谷,木野ら はecaljを主催live CDで実習。
(+元基礎工の長柄先生に参画してもらっている)。
*分業化せざるをえない
方法論レベルの開発をしていかねばならない。
自分たちで改造・改良できるコードが必要。
そして、それを理論家or実験家に渡していかねばならない。
*「物理学」と「数値計算技術(計算機のスキルも)」が必要
*私の戦略:
安易な近似手法を開発しても何がどうな
ってるのか不明になる。まずは真面目に解く。
10
Chap.1. 理論的な基礎
1.1 汎関数法を用いた多体摂動論の基本
1.2 密度汎関数法の基礎づけ (ルジャンドル変換)
1.3 時間依存密度汎関数法 (虚時間形式)
1.4 断熱接続、結合定数積分
1.5 Random Phase Approximation(RPA)
1.6 Luttinger-Ward 汎関数と自己エネルギー
1.7 その他のコメント
11
1.1 汎関数法を用いた多体摂動論の基本
12
多体論での汎関数の方法 1
Ref:Negele&Orland等
• グリーン関数の生成母関数 J, J*をprobe として
Z ( J *, J )
Generating functional method
*
†
= Tr T exp H Nˆ J (1) (r1 ) (r1 ) J (1) d 1
0
(r), † (r)に加えJ (r, ), J * (r, )を導入。r1,1 1 dr1の積分記号、 スピン略。
– 時間順序積
dU ( )
X ( )U ( ),U (0) 1の解は
dt
U ( ) 1
X ( M ) 1
X ( M 1) .. 1
X (2) 1
X (1)
M
M
M
M
M> M 1> M 2>... こ れを U ( ) T exp X ( )d と 書く 。
0
– J, J* (反交換演算子変数)での 『汎関数微分記号の意味』の取り決め
ˆ ˆ ˆ ˆ ... 0 と し て , J 1を 代入し た 結果を Z と 書く 。
Z 0 J ABCD
J
Z
た と え ば
J
と 書いた 方がわかり よ いかも 。 と 思う 。
J 1
J,J*は演算子。内積の外に(形式的に)とり出すときの取り決めがいる。
13
多体論での汎関数の方法2
汎関数微分
Z ( J *, J ) exp(-W ( J *, J ))
*
†
= Tr T exp H Nˆ J (1) (r1 ) (r1 ) J (1) d
0
をJ*, Jで汎関数微分していくと、グリーン関数が得られる。
(cf.キュムラント展開)
W ( J *, J )
1
Z ( J *, J )
J
Z ( J *, J )
J
相互作用がな い場合を H 0と し て Z 0 Z 0 (0, 0) ,
d1d 2
Z 0 ( J *, J ) Z 0 exp
H H 0 +Hi nt であ れば、
J * (1)g 0 (1 2) J (2) と な る 。
次ページ
Z ( J *, J ) exp d Hi nt
, * Z 0 ( J *, J )
J J
0
14
多体論での汎関数の方法3
• 摂動展開。
Z 0 exp d Hi nt
0
Z ( J *, J ) exp(-W ( J *, J ))
, * exp
J J
*
d
1
d
2
J
(1) g 0 (1 2) J (2)
g 0 (1 2) Tr T e( Hˆ Nˆ )( ) ˆ (1) e( Hˆ Nˆ )( ) ˆ † (2) e( Hˆ Nˆ )
0
1
0
1
2
e ( )(1 2 ) (1 2 )(1 n ) ( 2 1 )n
0
2
Tr e( H0 N )
ˆ
ˆ
大雑把にいえば、 Hi nt , * を コ ネク タ ーで、 g 0 (1 2)と いう ラ イ ン
J J
を 結線し ていく と Feynman diagramを 得る 。 次ページ
W ( J * , J )を 汎関数微分し て いく と 、 connect ed gr aphが
得ら れる 。 ルジ ャ ン ド ル変換し て *, を 考え る と
「 足を 落と し た既約グ リ ーン 関数が得ら れる 」 。 。 。
Ref .Negele & Orlandの教科書な ど
15
温度T=0の場合
[i ,0] [T (1 i ), T (1 i )]にし て、 T
G0
2
1
Z ( J *, J ) exp dt W ( J *, J )
H i nt
,
J J*
Z 0 exp dtHi nt
, * exp d1d 2 J * (1)G 0 (1 2) J (2)
J J
G 0 (1 2)
0 T H (r1 , t1 ) H † (r2 , t2 ) 0
Negle&Orland
こ こ で Hはハイ ゼンベルグ演算子
0 H (r1 , t1 ) H † (r2 , t2 ) 0 (t1 t2 ) 0 H † (r2 , t2 ) H (r1 , t1 ) 0 (t2 t1 )
electron
unocc.
(r ) (r ) e
*
1
2
i ( )( t1 t2 )
occ.
hole
(t1 t2 ) (r1 ) * (r2 ) ei ( )(t1 t2 ) (t2 t1 )
(r1 ) * (r2 ) occ. (r1 ) * (r2 )
G 0 (r1, r2 , )
( ) i
( ) i
unocc.
G0
プロパゲーターの図、摂動のイメージ図
E W (0,0)が基底状態のエネルギー
宿題:Z(J,J*)からどういうダイアグラムが得られるか書いてみる。
G0
16
汎関数法のメ リ ッ ト :
( a) 系の対称性( 保存則) が明瞭( e. g. ゲージ対称性 WT 恒等式)
形式論など に適し ている 。
( b) グリ ーン 関数が, エネルギー母関数の微分を
取っ ていく こ と で系統的に求ま る 。 系統的な近似法。
デメリット:
(物性ではありがたみが少ないのであまり好まれない?)
*通常の摂動論(中間状態など)との対応がよくない。
* J,J*が物理的なものでないし物理的な描像
がわかりにくい。(密度汎関数などではまし)。
「物理的な描像(仮定)にもとづく近似法」が重要。それをチェックする。
*独立粒子近似、準粒子
*断熱接続での全エネルギー
*密度汎関数法
*RPAによる密度ゆらぎ
*2体問題
17
1.2 密度汎関数法の基礎づけ
(ルジャンドル変換の立場から)
18
密汎
密度汎関数法の基礎
(DFの構成にはHohenberg-Kohnの定理、Levyの方法などは不要)
符号が変。
19
BとMの1対1対応はconvexity(凸性)で保障される
密汎
[ B (1 ) B '] [ B] (1 )[ B ']
が証明できる(有限系、 e.g. M.Valiev cond-mat/9702247v1)。これにより、
B M=-
の一対一対応が保障さ れる 。
B
それでLegendre変換が可能になる(F[M]もconvex)。多変数の場合も。
• 有限温度であっても無限大の系。 「Legendre変換」と「系∞」の非可換性。
強磁性の場合の ( B)を 考える 。 複素解析関数と し ての ( B)は,
多く の対数的特異点を も つが、 有限系においては実軸を 含む開領域
において解析的であ る 。 し かし 、 系の大き さ 無限大の極限では、
B 0を 含む線( 虚軸?) で解析的に分断さ れ、 微分に不連続が出る
( B 0と B 0で別の解析関数と なる ) 。 と 思う 。 (議論してほしいです)
20
密度汎関数
Hˆ J (r )nˆ (r )
Z [ J (r )] exp( W [ J (r )]) Tr exp
kBT
W [ J ]
n(r )
J (r )
W [ J ]
E[n] W [ J ] J
: ユニバーサルな汎関数
J (r )
E[n]
を も ちいれば、
J (r )と なる 。
n(r )
E [n] E[n]+ ( vext (r) )n(r)
total
21
Etotal [n] E[n]+ ( v ext (r ) )n(r )
Etotal [n]
を も ちいれば、
0を 解けばよ い。
n(r )
原子核から の外場の入っ て る のは二項目だけ!
E[ n ]=E0 [ n] EH [ n] E残り [ n]と いう 近似で評価-->L D A
他の変数を 加え る 場合には、 +Qˆ (r ) P (r )を 加え てルジ ャ ン ド ル変換。
Etotal [n, Q] E[n, Q]+ ( v ext (r ) )n(r )
Etotal [n]
など と なる 。
0も 付け加わる 。
Q(r )
Q(r )はス ピ ン 密度、 電流、 様々な秩序変数など 。
22
J (r)と n(r) の1対1対応はconvexity(凸性)で保障される
[ J1 (1 ) J 2 ] [ J1 ] (1 )[ J 2 ]
が有限系で証明できる。(e.g. M.Valiev cond-mat/9702247v1)
それで、ルジャンドル変換可能となる。一般的な変数で証明できる。
• T=0の極限(基底状態)では、1対1対応が保障されない縮
退を生じるような 外場vext (r)が存在する 。
孤立性の高い特異点であって、問題にはならない。迂回するこ
とができる(それを避けて解析接続できる)。と思う。
n(r)
det
0と なる 特異点。
J (r ')
• 任意のn(r)に対してJ(r)が存在する必要があるか?
凸性があるので、たまたま 「微分=0」となるn(r)が
見つかればよい。ただし無限系では、E[n]は
解析関数として分断されうるので注意が必要。と思う。
23
• n(r,r' ) を用いるような汎関数はつくれない。
n(r,r' )
(r) * (r ')であり 、 すべての固有値が1 かゼロ である 。
Occupied
Vext (r,r' ) n(r,r' )の対応はT 0で保てない。 異常になる 。
• 近年はルジャンドル変換をもちいて密度汎関数法を基礎づける文書が
多くみられるようになったが、まだ日本ではポピュラーではない。
•Density Functional Theory -- an introduction
http://arxiv.org/abs/physics/9806013v2, N.Argaman and G. Makov
• Electronic structure calculations with dynamical mean-field
theory G. Kotliar et al. Rev.Mod.Phys.78 (2006) 865
•高田康民「多体問題特論」朝倉書店。 p.20で言及している。
ただし引用192)に小谷らの論文がないのは不自然。
•高田康民「多体問題」朝倉書店も癖もあるがおおむね丁寧でそれなりに勧め
れる。
24
1.3 時間依存密度汎関数法
(虚時間形式)
25
密度汎関数の虚時間形式への拡張
Ω[J]はユニバーサルな関数。
26
27
ルジャ ン ド ル変換
F[n] [ J ] d1J (1)
J (1)
既約なグリーン関数を与える。
28
• 「Runge-Grossの定理」を使わなくても
Time-dependent DFが作れるのがわかる。
• Kohn-Sham方程式に直す時、DFにおけるlocal
potentialの弊害はそのまま引きずる。
• Time-dependent DFではRPA+αというような近似に
なる。バンドギャップが相当な過小評価をうけると
きに、+αの部分が本質的に事態を改善するか?
(当然、エキシトニックな効果は入りえない)。
•
n(r, t )
det
0と なる 特異点は問題ないか?
J (r ', t ')
29
1.4.断熱接続、結合定数積分
30
1.4 断熱接続
Hˆ Hˆ 0 Hˆ int
結合定数積分はいちばん簡単な場合
注:一般にはHintは、2次の項も含む。
(くりこみ理論のcounter termのようなものQSGW)
Nˆ (r)は面倒なので略。
ˆ
ˆ
ˆ
Z [Q, ] exp W [Q, ] =Tr T exp H 0 H int PQ(1) d
0
W [Q, ]
ˆ (1) d Z [Q, ]
d Tr T Hˆ int ( ) exp Hˆ 0 Hˆ int PQ
0
0
d Hˆ int ( )
0
これらに注意
Q ,
W [Q, ]
1
W [Q,1] W [Q,0]
d W [Q,0] d Hˆ int ( )
0 0
0
1
1
Q ,
d
1
Q 0のと き は、 W [Q 0,1] W [Q 0,0] Hˆ int
0
Q ,
d
31
Hˆ Veff nˆ (Vˆee Vext nˆ Veff nˆ ) に対し て適用し てみる 。
2m
温度ゼロでは,WをEと書くことにして、
1
E i
Veff nˆ i d Vext nˆ Veff nˆ Vˆee
2m
occ.
0
1
i i (Vext Veff )nˆ i d Vˆee
occ.
occ.
occ.
0
こ こ では「 を 変化さ せていっ ても nˆ
e
i
Vext nˆ i d
2m
2 0
2 1
二項目は、 EH EX EC
2
e n(r1 ) n(r2 )
2 | r1 r2 |
e
2
2
があ ま り 変化し ない場合」 を 考えた 。
† (r2 , t ) † (r1 , t ) (r1 , t ) (r2 , t )
| r1 r2 |
+
e
2
2
| r1 r2 |
n(r1 , r2 ) n(r2 , r1 )
1
d
0
こ こ でn(r1, t) nˆ(r1, t) n(r1)
n(r2 , t ) n(r1 , t )
n(r2 , t ) n(r1 , t )
0
| r1 r2 |
32
次ページのRPAで評価する
1.5 Random Phase Approximation(RPA)
ま ず外場摂動 + Vext( 1) nˆ( r1 ) に対する 密度の線形応答
を 考える 。 式で書けば、
n( 2) DR (2,1) Vext( 1)
であ る 。 線形応答理論から
DR (2,1) i n(r2 , t2 ), n(r1 , t1 ) (t2 t1 )と 書ける 。
こ こ でn(r1, t) nˆ(r1, t) n(r1)
*こ れがわかればただちに
D (r2 , t2 ) i Tn(r2 , t 2 ) n(r1 , t1 ) が評価でき る 。
33
Time-dependent Hartee近似(RPA)
2
Veff Vext (t ) k (r, t ) 0
i
t 2m
Veff Vext VH Vxc
k (x, t )が変化する 電子密度変化 VHが変化(
特に長距離ではVHが支配的)
2
Veff Veff (t ) k (r, t ) 0
i
t 2m
e2 n(r ', t )
Veff (t ) VH (t ) Vext (t ), VH (r, t ) dr '
|r r'|
t
n(t ) dt 'dr ' P(rt , r ' t ') Veff (r ', t ')
P(rt , r ' t ')( VH (t ) Vext (t )) P(rt , r ' t ')(v(r, r ') n(t ) Vext (t ))
よ っ て、 n(t )
( 1- Pv)1 P Vext (t )
分極関数P (rt , r ' t ')について は以下のよ う に求ま る 。
2
Veff (r, t ) G0 (r, r ', t , t ') (r r ') (t t ')
i
t 2m
書き 直すと 境界条件は n(r, t ) -iG0 (r, r ', t , t 0) ,
1
G0 (r, r ', t , t ')
2
Veff
i
t 2m
,
P(rt , r ' t ')
G0
G0
n(r, t ) - iG0 (r, r ', t , t 0)
こ れを Veff で微分し て 、
P (rt , r ' t ')
n(r, t )
Veff (r ', t ')
t
iG0 (r, r ', t t ')G0 (r ', r, t ' t )
と 書ける 。 ( 最終的に得た も のは t - t 'と し て よ い) 。
全エネルギーの表式
E i
occ.
EH
e
2
2
i Vext (r)n(r) EH EX EC
2m
n(r1 ) n(r2 )
e
i
1
d
0
2
, EX
| r1 r2 |
2
Ecについては Tn(r1 , t1 )n(r2 , t2 )
e2
Ec = 2
Hˆ Hˆ 0 Hˆ intに対し て
n(r2 , t )n(r1 , t )
n(r1 , r2 ) n(r2 , r1 )
| r1 r2 |
iP(1 vP) 1を も ちいる と 、
n(r2 , t )n(r1 , t )
0
| r1 r2 |
2
(t1 t2 )
e
1
d
vP
vP
(1
vP
)
i
vP
log
1
vP
,
こ
こ
で
v
0
2 | r1 r2 |
1
*クーロン力はテスト電荷の周りを分極させる。
これをとりこむにはRPAは必須。
*原理的にはvan der Wallsに対応するダイアグラムなど
を含む。
36
Hˆ Hˆ 0 Hˆ int
どのH0からスタートするのか?が重要。
結局のところ、self-consistentに決定するしかない。(何らかの
最小化or最適化)。
たとえば、「重い電子系」では「Z、くりこみ因子で重くなる」?
QSGW法
37
Ec i vP log 1 vP i vPvP vPvPvP vPvPvPvP ...
EH EX EC
…
*RPAでの自己エネルギー
EX
ECを G0の汎関数と みる 。
それでこ れを G0で微分する と 『 自己エネルギー』 を 得る 。
*全エネルギーのG0での微分= 2m Vext (r ) VH (r) (r , r ', )
これはJanakの定理に似たものを与えるChap.2で。
(実はGW近似は「グリーン関数」を求めるものではない)。
38
*密度汎関数F[n]の断熱接続はすこし違う。
違い:「与えた密度n」を維持するようにΛに対して外場を調整する。
違いに注意
39
40
41
1.6 Luttinger-Ward 汎関数
(Gを求めるための厳密なFormalism)
42
Luttinger-Ward汎関数
ˆ † 1 (r1 )ˆ 2 (r2 )
(注:ベクトルポテンシャル入れる)
?
43
高田康民「多体問題」より
44
前ページ
45
46
*RPAの全エネルギーと同様の近似式を得ることができる。
それをself-consistent に解く。Full self-consistent GW近似
ダイアグラムとしては通常のGW近似とおなじ。
理論的にきれいな気がする。
(停留問題であってエネルギー最小化ではない)。
しかしながら、非常に問題がある。この方向での進展には
望みがうすいと思う。
次ページ
47
Full self-consistent GW too problematic
E[G]
1
0 G
G
[G]
i GW
W
W and Γare given as a functional of G.
Difficulty 1. Z-factor cancellation
Zi
G
(incoherent part)
i
G
1
at q 0, 0
Z
Thus, you can not set 1 if we use G
Difficulty 2. If you use RPA like formula, D iG G,
This only contains QP weights by ZxZ.
Appendix in PRB76、165106(2007)
48
1.7 その他のコメント
•断熱接続(adiabatic connection)が重要。
これにより、H0Hの接続をおこなった。
このとき、相互作用のない粒子QPにつながる。
ˆ
ˆ
ˆ
ˆ ...
Z
QP for =0 1 =1では
QP
Λをswitch onしていく操作はくり込み操作。くりこみのcounter term (ポテ
ンシャルの差の項)を無視するわけにはいかない。
最終的には、 H0+Λ(H-H0)の分割におけるH0は「QP
(準粒子)」のものである。
49
cf.ゼロ音波
ハバードモデルでは、q0での電荷揺らぎとスピンの
縦揺らぎ、が抑えられない。
*第一にクーロンの長距離性を考慮する必要がある。
RPAが有効(q0のmost divergent termをあつめる形)。
これでスクリーンされた相互作用が得られる。
すくなくとも長距離極限がほぼ定量的によさそう。
50
*RPAでのd電子間相互作用の虚部は、Feなどで~0.6Ryぐらいに
大きい虚部をもつ。いろいろな局所的な電荷揺らぎモード。
d電子間相互作用はd電子自身の揺らぎモードでスクリーンされる。
T.Kotani J.Phys.C(2000)
51
*近距離のフェルミ相関。
Exchange-pair diagram:
クーロン相互作用を点に置き換えたときに同じになる
ダイアグラムのペア。
点に置き換えたとき、Exchange-pairは符号
(フェルミオンループの数)が違うだけで同じであり、
相殺をおこす。これを同時にとるのがパウリ排他律を
保つのに必要(図を描く)。RPAではこのことは無視される。
(ハバードモデルではHugenholz型相互作用をとるので自動的に入っている)。
*近距離のクーロン相関。
electron-holeの二体問題をマルチバンドで解く
(Bethe-Salpeter eqとよぶ)。これを取り込んだ分極関数の計算。
しかしクーロン(をスクリーンしたもの)は長距離なので大変。
*FLEX近似に対応することをするには、RPAレベルのすべての
Exchange pairをとる必要あり。しかし、RPAでは最初から
バランスを崩している(スピン横揺らぎなどだけ入れた方がよい?) 。52
•LDA+U+DMFT法など。 LW汎関数E[G]を基本にする。
木に竹をつぐような理論構成(?)。
•Double counging問題
•スペクトル関数の正定値性
•RPAを超える方法。第一原理的な部分的解決としては、magnonGW, phonon-GWやそれの組み合わせ。
53