第一原理分子動力学プログラム STATE Senri Wiki

STATE クイックリファレンス †

入力ファイルの例 (CO分子) †

 0  0  0  0  0  0                      : dummy line (6 integers)
 5.50 20.00  2  2  2                   : GMAX, GMAXP, NTYP, NATM, NATM2
 1  0                                  : space group number, bravis lattice type
 6.00  4.00  4.00  90.00  90.00  90.00 : a, b, c, alpha, beta, gamma
 1  1  1  1  1  1                      : knx, kny, knz, k-point shift
 1  0                                  : NCORD, NINV
 0.0000  0.0000  0.0000  1  1  1       : cps, iwei, imdtyp, ityp
 2.2000  0.0000  0.0000  1  1  2       : cps, iwei, imdtyp, ityp
 6  0.1500  51577.50 3 1 0.d0          : IATOMN, ALFA, AMION, ILOC, IVAN, ZETA1
 8  0.1500  51577.50 3 1 0.d0          : IATOMN, ALFA, AMION, ILOC, IVAN, ZETA1
 0  0  0  0  0                         : ICOND, INIPOS, INIVEL, ININOSE, INIACC
 0  1                                  : IPRE, IPRI
 200  200   0    57200.00  0           : NMD1, NMD2, iter_last, CPUMAX, ifstop
 3   1                                 : way_mix, mix_what
 0    8  0.8                           : starting mixing, kbxmix,alpha
 0.60  0.50  0.60  0.70  1.00          : DTIM1, DTIM2, DTIM3, DTIM4, dtim_last
 30.00    2     1  0.10D-08 1.d-06     : DTIO, IMDALG, IEXPL, EDELTA
  0.0010  0.10D+02    0                : WIDTH, FORCCR, ISTRESS
ggapbe          1                      : XCTYPE, nspin
  1.00     3                           : destm, n_stm
102                                    : NBZTYP
   0   0   0                           : NKX,  NKY,  NKZ  (dummy)
   0   0   0                           : NKX2, NKY2, NKZ2 (dummy)
   8                                   : NEG (# of bands)
       1                               : NEXTST (1: G-space, 0: R-space)
       0                               : 0; random numbers, 1; matrix diagon
       2                               : imsd (2: Davidson, 1: RMM)
       0                               : eval. eko diff.: .0 = no ,1 = yes
       0                               : npdosao
       0    0.0                        : SM_N doping

ダミー行 †

 0  0  0  0  0  0                      : dummy line (6 integers)

6個の整数を指定。歴史的経緯で存在するダミー行。
※過去のバージョンでは1カラム目の変数が1の場合は基本並進ベクトルをカーテシアン座標で与え、0の場合はa, b, c, alpha, beta, gammaで与えていたようである。
※現状STATEでは無視されているが、ユーティリティプログラムrepeat.fでxyzファイルを作る際に参照されている。

カットオフと原子数 †

 5.50 20.00  2  2  2                   : GMAX, GMAXP, NTYP, NATM, NATM2

• GMAX: 波動関数を展開する平面波のカットオフエネルギーを指定。E_{cut}^{wf} = GMAX**2
  
• GMAXP: 補強電荷を展開する平面波のカットオフエネルギー。E_{cut}^{dens} = GMAXP**2
  

GMAX=5.00 GMAXP=15.0 の条件で計算を行う場合は、一般的にpbe3 pbe1sのpseudo potentialのセットを用いるとよい
GMAX=5.50-6.0 GMAXP=20.0 の条件で計算を行う場合はpbe1のpseudo potentialセットを用いるとよい
※これは一般的なルールなのでエネルギーと力のカットオフエネルギーについての収束性は必ずチェックすべきである。

• NTYP: 原子の種類
  
• NATM: 空間反転対称性（NINV=1）がある場合の非等価な原子の数
  
• NATM2:全原子数、空間反転対称性が無い場合はNATM=NATM2である（現在、空間反転対称性はメンテナンスされていないため常にNATM=NATM2とすること）。
  

格子の対称性と種類 †

 1  0                                  : space group number, bravais lattice type

• NUM_SPACE_GROUP: 空間群の番号
  
• TYPE: ブラベー格子の種類を番号で指定:

  #	bravais lattice type
  0	simple
  1	body-centered
  2	face-centered
  3	a-face-centered
  4	b-face-centered
  5	c-face-centered
  6	rhombohedral

†

 6.00  4.00  4.00  90.00  90.00  90.00 : a, b, c, alpha, beta, gamma

• a, b, c: 格子ベクトルの長さ
  
• alpha, beta, gamma: 格子ベクトルのなす角度

これらの代わりに以下のように"Cartesian"と記述するとそれ以下の行で格子ベクトルを直接カーテシアン座標で指定することができる。

Cartesian
 6.00  0.00  0.00
 0.00  4.00  0.00
 0.00  0.00  4.00

†

 1  1  1  1  1  1                      : knx, kny, knz, k-point shift

Brillouin zone積分で用いるk点のメッシュ（knx, kny, knz）を最初の３カラムで、k点のオフセットを次の3カラムで指定する。 4-6カラムでのk点のオフセットの指定方法は、ガンマ点を含む場合には1を、ガンマ点を含まない場合（Monkhorst-Pack）は2を指定する。 六方晶の場合には対称性の破れが起こるためk点はオフセットしない方が良い。

原子位置の指定方法と空間反転対称性 †

 1  0                                  : NCORD, NINV

• NCORD: 1: カーテシアン座標で原子座標を入力（単位はBohr）、0: 基本並進ベクトルを単位として入力、2: conventional unit cellの格子ベクトルを単位とする
  
• NINV: 0: 系に空間反転対称性がない、1: 空間反転対称性がある。
  このオプションを使用すると、空間反転対称性が存在する場合には明示的に考慮する原子数をおよそ半分にすることが可能で、計算時間の短縮をすることが可能である。

原子位置と原子種 †

 0.0000  0.0000  0.0000  1  1  1       : cps, iwei, imdtyp, ityp
 2.2000  0.0000  0.0000  1  1  2       : cps, iwei, imdtyp, ityp

• 1-3カラム: cps(katm,3) (pos(katm,3)): 原子座標。NCORDが1の場合カーテシアン座標で、NCORDが0の場合、基本並進ベクトルを単位として、NCORDが2の場合、コンベンショナルセルの格子ベクトルを単位として指定する。
• 4カラム: iwei(katm): 空間反転対称性について等価な原子数
• 5カラム: imdtyp(katm): Hellmann-Feynman力に従って原子を動かす場合には1、そうではない場合には0を指定。有限温度分子動力学(FTMD)の場合、1000以上の数を指定するとその原子は熱浴に接することを指定。
• 6カラム: ityp(katm): 原子の種類を番号で指定。原子番号等は後に続く列に記述する。

原子（擬イオン）に関する情報 †

 6  0.1500  51577.50 3 1 0.d0          : IATOMN, ALFA, AMION, ILOC, IVAN, ZETA1
 8  0.1500  51577.50 3 1 0.d0          : IATOMN, ALFA, AMION, ILOC, IVAN, ZETA1

• 1カラム: 原子番号
• 2カラム: 初期電荷（ダミー）
• 3カラム: 原子の質量(Da(旧のa.m.u.)(＝水素原子が約1))。後述するが、上の値は現実のものとは異なる。Quenched MD(IMDALG=2)を行うときには原子質量は同じ値、例えばCの12.0107やSiの28.0855などにしておくと計算を安定して行うことができる。 有限温度の分子動力学を実行するときには現実の質量に設定するべきである。
• 4, 5カラム: Local potentialの角運動量とウルトラソフト擬ポテンシャルの指定（ダミー）
• 6カラム: zeta=(rho(up)-rho(down))/(rho(up)+rho(down))の初期値を与える。スピン分極の計算を行う際はnspin=2としてzetaとしては0でない有限の値を設定する。zeta=0とすると非磁性の解に収束するだろう。
• 6カラム: 磁化率の初期値(一般的に分極してそうなものは0.5、そうでないものは0.0を入力する)
• その他: ilocはローカルポテンシャルの軌道量子数l+1、ivanはVanderbiltのultrasoft pseudopotentialを使用するかどうか、といったパラメータであるが、入力ファイルからではなく擬ポテンシャルのファイルから読むので実際にはダミーパラメータである。ちなみにAlFAという変数は現在は使用されていない。電荷密度を乱数から発生させるときに使用していたようである。

計算の継続（初期化）に関するオプション †

 0  0  0  0  0                         : ICOND, INIPOS, INIVEL, ININOSE, INIACC

• INIVEL: 原子の速度
  • INIPOSと同様、1の場合はrestart.dataから、2の場合はGEOMETRYから速度を読み込む。INIPOSの時と同様にrestart.dataが必要。
• ININOS: 能勢の温度制御
• INIACC: accumulator:物理量の統計平均の計算の継続に関するオプション

計算の継続する例 (1): SCF計算の継続

• SCF計算のみを継続したいとき（geometryは入力ファイルから読み込む）

   1  0  0  0  0                         : ICOND,INIPOS,INIVEL,ININOS, INIACC

• 構造最適化を継続 (restart.dataを参照)

   1  1  0  0  0                         : ICOND,INIPOS,INIVEL,ININOS, INIACC

• 構造最適化を継続 (GEOMETRYを参照 (restart.dataは構造最適を継続する際必ず必要))

   1  2 0  0  0                         : ICOND,INIPOS,INIVEL,ININOS, INIACC

• 構造最適化を継続 (restart.dataを参照、波動関数を初期化)

   0  1  0  0  0                         : ICOND,INIPOS,INIVEL,ININOS, INIACC

• MDの継続（IMDALG <=0）

   1  1  1  0  0                         : ICOND,INIPOS,INIVEL,ININOS, INIACC

ストレス計算、出力レベル †

 0  1                                  : IPRE, IPRI

• 圧力（ストレス）の計算（現在は実装されていないのでダミー変数である）
• 標準出力のレベル（通常は1、1よりも大きいと多くの情報が出力される。デバッグ用のオプションではあるが、出力される情報はあまり参考にならないようである）

繰り返しの回数と計算時間に関するオプション †

 200  200   0    57200.00  0           : NMD1, NMD2, iter_last, CPUMAX, ifstop

• NMD1: 電子系のSCF計算の回数の上限を指定する。
• NMD2: 分子動力学 (構造最適化) のステップ数の上限を指定する。NMD2+1が実際のMD（構造最適化）の回数の上限である。
• iter_last:
• CPUMAX: CPU時間の上限を秒単位で指定
• ifstop:

ミキシング †

 3   1                                 : way_mix, mix_what

ミキシングに関するパラメータ †

 0    8  0.8                           : starting mixing, kbxmix,alpha

• start mixing: dummy
• kbmix: ミキシングに用いる履歴の回数
• alpha=mix_alpha: ミキシングパラメータ

RMMのミキシングに関するパラメータ †

 0.60  0.50  0.60  0.70  1.00          : DTIM1, DTIM2, DTIM3, DTIM4, dtim_last 

• DTIM2: RMMで使用、その他はダミー変数.
• dtim_last: ダミー変数

構造最適化または分子動力学(MD)に関するオプション †

 30.00    2     1  0.10D-08 1.d-06     : DTIO,IMDALG,IEXPL,EDELTA

• DTIO: イオン系のMDで用いる時間刻み（タイムステップ）(a.u.)。Quenched MD (IMDALG=2)またはMDで有効。　1a.u.=0.024188fs, 41.3428a.u.=1fs
• IMDALG: MDまたは構造最適化の方法を指定。正の値は構造最適化、0または負の場合は有限温度MDを実行
• Molecular dynamicsも参考にしてください。

• IEXPL: dummy
• EDELTA: 電子系の収束の閾値を規定。収束は一原子あたりの全エネルギーで評価している。

†

  0.0010  0.10D+02    0                : WIDTH, FORCCR, ISTRESS

• WIDTH: Fermi面のぼかしの幅. 正の場合はデルタ関数を多項式で、負の場合にはHermite-Gaussianの一次式で近似する。WIDTH<-10.0のときにはTetrahedron法によるBrillouine zone積分を行う。この時はNBZTYPが有効になる。
• FORCCR: Forceの閾値。Forceがこの値よりも小さくなると構造最適化を終了する。IMDALG>0のみ有効である。
• ISTRESS: ストレスの計算用オプション。現在は使用されていない。

†

ggapbe          1                      : XCTYPE, nspin

• XCTYPE: 交換相関汎関数を指定.
  • ggapbe: Perdew, Burke, Ernzerhof GGA (1996)
  • ggapw91: Perdew-Wang GGA (1991)
  • ldapw91: Perdew-Wang L(S)DA (1991)
  • RPBE: revised PBE of Hammer et al.
  • revPBE: revised PBE of Zhang and Yang
  • WC: Wu-Cohen GGA
  • vdW-DF1/vdW-DF2/optB88-vdW/optB86b-vdW/rev-vdW-DF2
• NSPIN: 1ならば非磁性、2でcollinearな磁性体の計算を行う。

STM用オプション †

  1.00     3                           : destm, n_stm

†

102                                    : NBZTYP

• WIDTH < -10.0の時のみ有効。101が推奨値。

ダミー行 †

   0   0   0                           : NKX,  NKY,  NKZ  (dummy)
   0   0   0                           : NKX2, NKY2, NKZ2 (dummy)

計算に用いるバンドの数 †

   8                                   : NEG (# of bands)

STATEの計算では、常に非占有状態を含めて計算を行う必要がある。非磁性のSiならば、例えば4 (電子数/2)+4などとするのが良いだろう。

nonlocal pseudopotentialの計算方法 †

       1                               : NEXTST (1: G-space, 0: R-space)

• 0: R-space (Davidsonの時は使用しない）
• 1: G-space

ダミー行 †

       0                               : 0; random numbers, 1; matrix diagon

対角化法の選択 †

       2                               : imsd (2: Davidson, 1: RMM)

• 1: RMM-DIIS
• 2: Davidson

大規模な計算をするときにはRMM-DIISとReal space projectionの組み合わせ（NEXTST=0、IMSD=1）が推奨されますが、初期波動関数は乱数から発生させるものではなく、Davidsonである程度収束させたものを使うのが良いでしょう。そうしなければ最安定な電子状態に収束しない可能性がある。

†

       0                               : eval. eko diff.: .0 = no ,1 = yes

前回のiterationにおける擬固有値との比較を行うかどうかを指定。ほぼダミー変数

PDOSの計算オプション †

       0                                                   : npdosao

0でない場合にPDOSの計算を行う。但しnpdosaoは原子の数より小さくなければならない。具体的な計算方法は、柳澤氏による以下の文書を参照のこと。

f電子用の半経験的パラメータ用のオプション †

       0    0.0                                          : SM_N doping

1コラム目は整数を、2コラム目は実数を入力。現在使用している人はいないため全くテストされていません。

PDOS/AOLDOSの設定 †

AOLDOS。Layer-resolved DOSの計算の指定は以下を参照。

原子軌道に射影したPDOS(Partial Density Of States、もしくはAO_LDOS)の計算方法 †

(以下、pdos用のオプション。仮にnpdosao = 5とし、通し番号が1,3,4,9,10の原子の軌道に投影したDOSを計算、原子の種類は2種類。ただし、lda+uでない場合)

   5                                   : NPDOSAO
   1                                   
   3
   4
   9
  10
-15.00 5.00 0.20 501                   :EPDOS(1), EPDOS(2), EPDOS(3), NPDOSE
2.2   0.3                              :RAD, WIDTH FOR TYPE 1
1.0   0.3                              :RAD, WIDTH FOR TYPE 2
0.2   12                               :DR, NR
&OTHERS
 GAUSSDOS
&END

• npdosao: DOSを投影させる原子数。その次の行にnpdosao個の原子の通し番号を続ける
• EPDOS(1,2): pdosを表示するエネルギー領域(Fermi準位との相対値、eV)
• EPDOS(3): gaussianを用いる場合の標準偏差、tetrahedronを用いる場合には未使用
  • GaussianBroadeningを使うためには次のオプション指定が必要なので注意。

    &OTHERS
     GAUSSDOS
    &END

• NPDOSE: EPDOS(1)～EPDOS(2)の間のエネルギーメッシュ数
• rad: DOSを投影させる原子中心の波動関数(擬波動関数)の、動径方向のカットオフ値を指定
  (bohr単位。通常、原子の価電子半径程度を指定。e.g. H:0.6 程度, C,N,O:1.4 程度)
• width: 上記のrad前後で、原子の擬波動関数をFermi分布関数にしたがってぼかす際のぼかし幅 (0.3程度でよい模様)
• DR,NR: (おそらく)未使用
  (元々は、擬波動関数同士やプロジェクターbetaとの重なり積分を実空間で計算する際の刻み幅・積分範囲のパラメータを意図していた模様)

Layer-resolved PDOSの計算 †

xy面に平行な面でスライスしてできた領域内に含まれる電子のDOSを算出する。 いったん電子状態を求めてzaj.data,potential.dataを出力した後に、ICOND=33として実行する。下記のようなnfaldos.dataというファイルを用意する。フェルミエネルギーはHartree単位で与える（直前のSCF計算の最後に得られるFermiEnergyをnfoutから読み取って書く）、Layerの境界はBohr単位のz座標である。Layerの数より一つ大きい数のz座標を与える。DOS計算のためのエネルギー範囲と分割数、ぼかし幅(Gaussian broadening factor)を与える。

&ALDOS
PRINT_WEIGHT
FERMI_LEVEL  -0.14395389
NUMBER_OF_LAYER 7
LAYER
-7.7938
-5.3974
-2.9894
-0.5968
 1.41385
 3.832933
 6.64160692845
 10.60270297025
 EMIN -15.D0
 EMAX   5.D0
 EWIDTH 0.2D0
 NPDOSE  401
&END