改良テトラへドロン-Brillouin領域積分法による金属中のフォノンの第一原理計算 B 東大理、東大物性研河村光晶A 、合田義弘A 、常行真司A,B A 新学術領域研究「コンピューティクスによる物質デザイン:複合相関と非平衡ダイナミクス」平成 25 年度第 2 回研究会 2014年3月10日　東京大学本郷キャンパス工学部 6 号館 1/18 目次 Motivation ● 超伝導体のフォノン計算 Brillouin領域積分の計算コストと精度通常の線形テトラへドロン法 ● テトラへドロン法の原理線形テトラへドロン法によるフォノン計算結果改良テトラへドロン法 ● 原理結果 2/18 まとめ ● フォノンの第一原理計算フォノン振動数の計算=動力学行列の計算 i q R˙ 2 R l s =R l +¿ s + u q s e 1 ∂ E ® ¯ C̃ s t (q )= 原子位置の周期的変化 N C ∂ u* ∂ u qs® qt¯ l [( ) * 2 ∂ V ion (r) ∂ ½ (r) ∂ V ion (r) 1 ion ̃ ® ¯ 3 = Cst + d r + ½ (r) ∫ NC ∂ uq s® ∂ uq t ¯ ∂ u* ∂ u 電子による応答 q s® Electrons I qt¯ I ] I µ ("F −"n k )−µ("F −" n k +q ) * ∂ ½ (r) V 3 = d k∑ ' nk (r)' n k +q (r) ∂ uq s ® (2 ¼)3 ∫BZ " −" nk n k +q nn SCF ¢ ∂V (r ) 3 * ¢ q ×∫ d r ' n k +q (r ) ' n k (r¢ ) ∂ uq s ® 密度汎関数摂動理論(DFPT) ¢ ¢ ¢ ¢ ¢ S. Baroni et al., RMP 73, 515 (2001). 4/18 Brillouin領域積分の計算コストと精度 V 1 3 3 ∫ d k⋯≈ Nk (2 ¼) finite ⋯ ∑ k µ("F–"nk)をそのままつかうと、分割数(k点数)に対して激しく振動する。 k点数依存性を穏やかにするためには、µ("F–"nk)をなめらかに変化する関数で置き換える。 N 2 x x (−1)n y y µ (x)=∫−∞ d y ± (y )≈∫−∞ d y ∑ H 2n exp − n ¾ ¾ n=0 ¾ n!4 √ ¼ M. Methfessel and A. T. Paxton, PRB 40, 3616 (1989). k点を多くとると・・・・ ( ) [ ( )] 積分の精度が良くなる。計算コストが増える。例/DFPT:O(Nk),hybrid-DFT,GW:O(Nk2) etc. 5/18 いかに精度を犠牲にせず計算コストを下げるか。事例2/FCC-Li(20GPa) K点横波モードウルトラソフト擬ポテンシャル(1Sを価電子に含む) GGA-PBE汎関数波動関数cutoff : 80 Ry Smearing : Methfessel-Paxton(1次) ¡ K 7/18 テトラへドロン法 "1 VT "2 O. Jepsen and O. K. Andersen, Solid State Commun. 9, 1763 (1971) "4 "3 V BZ V T= 6Nk 四面体内でKohn-Shamエネルギー等を線形だと見なすと・・・例 " 1≤"<" 2≤" 3 ≤"4 ならば、 "−"1 "−" 1 "−"1 ∫T d k µ ("−"k )≈V " ≤"=V T "2−"1 "3−"1 "4−"1 解析的に計算可能 3 k これが出来るのは一次補間であるときのみ。 2次補間、3次補間・・・では占有されている領域の面積を解析的に計算できない！ 8/18 通常のテトラへドロン法による結果 Quantum ESPRESSO のDFPT計算部分にテトラへドロン法を実装。ブロードニング幅収束性を調べる手間は無くなったが、k点数についてもっと早く収束しないだろうか。 9/18 k点が少ないときの誤差の原因 Kohn-Shamエネルギー等を線形補間することによる誤差 T T 補間される関数が下に凸ならば、線形補間ではつねに過大評価する。「補間」ではなく「フィット」 10/18 改良テトラへドロン法通常のテトラへドロン法改良テトラへドロン法 T T k 1, KS固有値等をk点間で一次補間する。 k 1, 2次以上の多項式でKS固有値等を補間する。 2, 補間した多項式を四面体内部で一 2, 補間した一次式を用いてテトラへ次式で最小2乗フィットする。ドロン法を行う。 3, 多項式にフィットした一次式を用いて通常のテトラへドロン法を行う。 11/18 改良テトラへドロン法による結果1 MgB2 全エネルギー線形テトラへドロン法の系統的な誤差を減らしている。 ※Blöchl補正は全エネルギー、電荷密度の計算にしか使えない。 12/18 P. E. Blöchl, O. Jepsen, and O. K. Andersen, PRB 49, 16223 (1994). MgB2 ボンド伸縮モード本研究の手法 0.01 Ry 線形テトラへドロン 0.03 Ry 13/18 k点数についての収束を大幅に改善している。 Li K点横波モード 14/18 McMillan公式によるTC見積り V 3 ∗ ¸= d q ¸ ¹ =0.1 qº 3 ∫BZ (2 ¼ ) V 3 ! log =exp d q log (! q º )¸ q º 3 ∫BZ (2 ¼) W. L. McMillan, PR 167, 331 (1968). [ ! log 1.04(1+¸) T C= exp − * 1.20 ¸−¹ (1+0.62 ¸) ¸q º = ] ( 1 V 3 qº 2 d k ∣g ∣ ± (" n k −"F )± ("n k+q −" F −! q º ) ∑ 3∫ nk n k +q N ("F )! q º (2 ¼) nn ¢ ¢ ¢ q点(フォノン): 6×6×4(MgB2) 8×8×8(Li) k点(KS軌道): 24×24×24(MgB2) 32×32×32(Li) !qやgqのDFPT計算に用いるk点数だけを変化させる。 15/18 ) MgB2のTC 16/18 LiのTC 17/18 Appendix 1:Phonon Periodic distortion Dynamical matrix C̃ 2 ®¯ st 1 ∂ E (q )= N C ∂ u*q s ® ∂ u q t ¯ 1,Frozen phonon 2 1 ∂ E 1 ∂ F ls® ®¯ →Fourier transform C s t (R l−R l )= = N C ∂ Rl s® ∂ Rl t ¯ N C ∂ Rl t ¯ ¢ ¢ 2, Linear response C̃ ®s t¯ (q)=elC̃ s®t¯ + ion C̃ ®s t¯ 19/18 ¢ Density response & SCF cycle Density linear response ∂ n(r) =2 ∑ 〈n k∣r 〉 〈r∣ ∂ ∣n k 〉 ∂ uq s ® ∂ uq s ® nk ̂ SCF ∂ V ̂ k + q ) ∂ ∣n k 〉 =−[ µ(" F −" n k )−P̂ n ,k ,k +q ] ̂ SCF −"n k + Q ∣n k 〉 (H ∂ uq s ® ∂ uq s ® ¢ ¢ ∂ V̂ SCF ∂ V̂ ion ± v (r) ∂ n(r ) ∂ n(r ) 1 3 ¢ 3 ¢ xc = + d r +∫ d r ¢ ∂u ¢ ∂ u q s ® ∂ uq s ® ∫ ∣r−r ∣ q s ® ± n(r ) ∂ uq s ® 20/18 3rd order polynomial −3 f 1 +6f 2 −2 f 5 −f 14 −3f 1 +6 f 3 −2f 9 −f 11 −3f 1 +6 f 4 −f 8 −2 f 13 x+ y+ z 6 6 6 −2 f 1 +f 2 +f 5 2 −2f 1 +f 3 +f 9 2 −2f 1 +f 4 +f 13 2 + x+ y + z 2 2 2 −6 f 1 +6f 3 +5 f 5 −f 6 +f 9 −f 11 +f 14 +f 15 −6f 18 + xy 6 −6 f 1 +6f 3 +f 7 +f 8 +f 9 −f 11+5f 13 −f 16 −6 f 20 + yz 6 −4 f 1 +2 f 2 +2 f 4 +f 5 +f 8 −f 10−f 12 +f 13 +f 14 −2f 17 + xz 2 3 f 1 −3f 2−f 5 +f 14 3 3 f 1 −3f 3 −f 9 +f 11 3 3 f 1−3 f 4 +f 8 −f 13 3 + x+ y + z 6 6 6 3 f 1 −6f 3 −4 f 5 −f 6 +f 9 +2f 11 +f 14−2 f 15 +6 f 18 2 f 1 −2f 3 −f 5 +f 11−f 15 +2f 18 2 + x y+ xy 6 2 9 f 1 −6f 2 −6 f 4 −2f 5−2f 8 +f 10 +2 f 12−f 13 −f 14 +6 f 17 2 + x z 6 9 f 1 −6f 2 −6 f 4 −f 5 −f 8 +2 f 10 +f 12 −2f 13 −2 f 14 +6 f 17 2 + xz 6 f 1 −2f 3 −f 7 +f 11 −f 13 +2 f 20 2 3 f 1 −6f 3 −2 f 7 +f 8 +f 9 +2f 11−4 f 13 −f 16 +6f 20 2 + y z+ yz 2 6 3 f 1 −6f 3 −2 f 5 +f 6 −f 7 −f 8 +f 9 +2f 11−2 f 13−f 14 −f 15 +f 16 +3 f 17 +3 f 18−3 f 19 +3 f 20 + xyz 3 f 3T (x )=f 1 + 22/18