科研費新学術領域研究「コンピューティクス」研究会

オーダーN法第一原理計算プログラム
CONQUESTにおける局在基底の最適化
物質・材料研究機構(NIMS)

理論計算科学ユニット

中田彩子、Sergiu Arapan、宮崎
University College London (UCL)

D. R. Bowler

剛Outline

オーダーN法第一原理計算手法プログラムCONQUEST

局在軌道の最適化
PAO基底関数を用いた時の精度検証
PAO係数の最適化について

京速コンピュータでの並列化効率

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

オーダーN法第一原理計算プログラム CONQUEST
密度行列最適化法 (DMM法)(ρ(r,rʼ) = Σν fν ψKSν∗ (r) ψKSν (rʼ) )密度行列の局所性

ρ(r,rʼ) = Σiα, jβ φiα∗ (r) Kiα, jβ φjβ (rʼ)（ν は状態（バンド）、iは原子、αは軌道。 φiα(r)は局在軌道)

密度行列の局所性-> Lにcutoﬀ RLを導入し、DFT(Kohn-Sham)の全エネルギーを補助行列Liα, jβ と φiα (r) について最小化(Liα, jβ=0 for |Ri-Rj| > RL )(LNV method： K=3LSL-2LSLSL)

密度行列のカットオフ半径を長くすると1. 全エネルギーが正確になる（変分法）
2. 計算量は増える。

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTでの最適化構造
Movement of atoms (MD or Structure Opt.)

Optimisation of φiα (r) (support func.)

1. Full DFT
Optimisation of density n(r)(achieving SCF)

2. SCF-AITB

Optimisation ofdensity matrix Kiα, jβ (Liα, jβ )orDiagonalisation of Hiα, jβ

3. NSC-AITBK=3LSL-2LSLSL
L: auxiliary density matrix
S: overlap matrix

2012.3.16 新学術領域「コンピューティクス」研究会@東大工応用例：イオンチャネル gramicidin A
有田通朗（東京理科大D1）

通常のSCF計算

2012.3.16 新学術領域「コンピューティクス」研究会@東大工応用例：イオンチャネル gramicidin A
有田通朗（東京理科大D1）

密度行列＋SCFの同時最適化

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTでの最適化構造
Movement of atoms (MD or Structure Opt.)

Optimisation of φiα (r) (support func.)

1. Full DFT
Optimisation of density n(r)(achieving SCF)

2. SCF-AITB

Optimisation ofdensity matrix Kiα, jβ (Liα, jβ )orDiagonalisation of Hiα, jβ

3. NSC-AITB
K=3LSL-2LSLSL
L: auxiliary density matrix
S: overlap matrix

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道✦ 密度行列は、局在軌道 φiα (r) で表される。
ρ(r,rʼ) = Σiα, jβ φiα∗ (r) Kiα, jβ φjβ (rʼ)
φiα(r) : centred on atom i with the orbital α.
✦ 局在軌道 φiα (r) は、基底関数で表される。（２種類準備）
blip functions

高精度(B-splines on regular grids around each atom i
Pseudo Atomic Orbitals (PAOs)

高効率

blip functions

PAOs

Localized Atomic Orbital Basis Sets inthe Projector AugmentedWave Methodiα

φiα(r) = Σlmτ clmτ Rlτ (r)Ylm(θ,φ)

Rreg
grid-spacing h <-> cutoﬀ energy in PW calc.
Ask Hjorth Larsen, s021864

Version 1.1. June 24, 20082012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道
blip 基底を用いた場合：
systematicに精度の制御が可能（自由度大）最小基底で平面波精度を持つ事が可能

PAO基底を用いた場合： φiα(r) = Σlmτ ciα lmτ Rlτ (r)Ylm(θ,φ)
PAO基底＝局在軌道（support func.)として計算する。（最適化の必要なし。全ての係数ciα lmτ =1。reference systemについて最適化。）
計算量は局在軌道（support func.)の軌道の数の３乗で増加する。
（SZ:43, DZ:83, SZP:93, DZP:133, TZP:173, TZDP:223, ....)
PAO係数の最適化の回数 <-> 軌道の数の増加による計算量の増加

PAO係数の最適化
効率的なソルバー
良い初期値
2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道
基底関数blipの時 (bulk Ge)

Etot- Region Radius（局在基底の半径依存性）

Etot- Grid Spacing（平面波計算でのcutoﬀ energyに対する収束）

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道
基底関数がPAOの場合
E-V curve
例：bulk GaAs

例：bulk SioSZP	 :	 局在軌道(support	 func.)の数が９(s:1,	 p:3,	 d:	 5)
2012.3.16 新学術領域「コンピューティクス」研究会@東大工

オーダーN法第一原理計算プログラム CONQUEST
例：bulk GaAsCQ blip vs VASP
phonon disp.

E-V curve
CQ blip vs CQ DZP

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道
例：bulk GaAs

CQ blip vs CQ DZPphonon disp.

E-V curve
CQ DZP vs CQ oSZP

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道
Raysonʼs method (PRB 80, 205104 (2009))
各原子から距離がRの範囲内の部分系を作成する
部分系のKohn-Sham 軌道 |i〉(energy λi）を求める|tkʼ〉 = Σi |i 〉f(λi)〈i | tk〉 (f: Fermi-Dirac 関数）| tk〉の選び方によって基底の精度が決まる

CONQUESTへの導入
PAO係数の初期値として用いる
可能であれば、one-shot 計算
multi-site PAOで一つの局在軌道を作るO(N)でない場合にも有効。
2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道
Total energy of H2O [hartree]
total energy

!" fromprimitive

primitive

-76.409366

6-31G**

-76.398461 0.01091

primitive:
7 gaussians for H
28 gaussians for O(GAMESS)

Raysonrcut = 4.0 -76.409360 0.00001
(SZ)

2.0 -76.408078 0.00129
0.5 -76.218438 0.19093

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの局在軌道• !"#$%&'(')*+&",&-./0.1
SZ
primitive:
7 gaussians for H28 gaussians for C(GAMESS)

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの並列化効率について
京速コンピュータ試験
利用期間中の測定データ

CONQUESTの並列化
疎行列の行列積

Ciα, jβ = ΣkγAiα,kγBkγ,jβ• Sends all B data
MPI_isend• Loop over remote nodes• Loop over sets of atoms, K•receive B data for next K
MPI_irecv (or local copy)•check MPI_irecv for present K
MPI_wait for MPI_irecv•Perform multiply-and-add(m_kern_max, m_kern_min)•End loop over sets of atoms, K“ノンブロッキング通信を使い演算と通信をオーバーラップさせる”

call MPI_test for MPI_isend
do k ∈ K
perform matrix_multiplication
call MPI_test for MPI_irecv
enddo• End loop over remote nodes

MPI_waitall (for MPI_isend)2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの並列性、計算規模

2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの並列性、計算規模
京速コンピュータ試験
運用期間中の測定データ
bulk Si (3万ー50万原子, 64-1024 nodes）
LDA, 最小基底、Grid Cutoff= 36 Ry、RL=14 bohr
昨年末： SZP （基底の数：９/原子）245,760原子/1920 nodes: 387.7 sec
491,520原子/3840 nodes: 398.6 sec
先月：DZP（基底の数：13/原子）245,760原子/1920 nodes: 846.7 sec
491,520原子/3840 nodes: 878.9 sec
2012.3.16 新学術領域「コンピューティクス」研究会@東大工

CONQUESTの並列性、計算規模

並列化効率 88.5%

393,216-atom system (DZP)
3072 node : 1116.29 sec
12288 node : 349.14 sec (79.9 %)2012.3.16 新学術領域「コンピューティクス」研究会@東大工

- まとめ 局在軌道（support function）をPAO基底で表した時のPAO係数の最適化について最適化を行えば、軌道の数はSZPで十分な精度を持つ。
現在Raysonの方法などを取り入れて、効率的なPAO係数の最適化を模索している。
「reference systemに対する最適化」の検討要

CONQUESTの並列化効率に対して
MPI̲testallを呼ぶ事によって計算時間（通信時間）を大きく削減できる
DZPの計算であれば、4atoms/coreくらいまで十分なstrong scalingの並列化効率が得られる。

2012.3.16 新学術領域「コンピューティクス」研究会@東大工