Quantum ESPRESSO入力ファイル作成手順２．SCF計算の設定

クロスアビリティ Winmostarサポートチームです。
WinmostarはQuautum ESPRESSOを簡単かつ高度に利用するための統合GUI環境となっています。学生は無料で利用でき、学生以外も無料トライアルを入手することができますので、是非 WinmostarのWebサイト をご覧ください。

#0. はじめに
第一原理計算ソフトQuantunm ESPRESSO https://www.quantum-espresso.org の入力ファイルの書き方について、複数回にわたって解説していきます。前回は、結晶構造の設定方法を紹介しました。今回は第２回、SCF計算の設定です。

この記事では、初歩的な内容しか紹介しません。
入力ファイルの詳細について知りたい方は、マニュアル http://www.quantum-espresso.org/Doc/INPUT_PW.html をご参照下さい。

#1. 入力ファイル
SCF計算を実施する場合、入力ファイルをこのように書きます。

&CONTROL
    calculation = "scf"
    pseudo_dir  = 擬ポテンシャルの格納場所
/

&SYSTEM
    ibrav       = 4
    a           = 3.09201
    c           = 5.07335
    nat         = 4
    ntyp        = 2

    ecutwfc     = 25.0
    ecutrho     = 225.0

    occupations = "smearing"
    smearing    = "gaussian"
    degauss     = 0.01
/

&ELECTRONS
    conv_thr    = 1.0e-06
    mixing_beta = 0.4
/

K_POINTS {automatic}
 4  4  2  0 0 0

ATOMIC_SPECIES
Si  -1.0  Si.pbe-rrkj.UPF
C   -1.0  C.pbe-rrkjus.UPF

ATOMIC_POSITIONS {angstrom}
Si     -0.000002   1.785172   2.534588
Si      1.546005   0.892585   5.071262
C      -0.000002   1.785172   4.441264
C       1.546005   0.892585   1.904590

上記は、SiCの例です。必要最小限の内容を記したものです。入力ファイルの各項目について、以降の章で解説します。

#2. 基本設定
先ずは、calculationという変数に"scf"と設定します。これを設定すると、SCF計算が実施されます。必須の内容ですので、忘れずに設定して下さい。

&CONTROL
    calculation = "scf"
/

また、計算の実施には、入力ファイルの他に擬ポテンシャルファイルが必要となります。擬ポテンシャルファイルの格納場所を、以下のように設定します。

&CONTROL
    pseudo_dir  = 擬ポテンシャルの格納場所
/

入力ファイルのこの設定は省略可能です。その場合、環境変数ESPRESSO_PSEUDOの値が参照されます。入力ファイルに都度都度pseudo_dirを記入するのが面倒な場合には、環境変数ESPRESSO_PSEUDOを事前に設定しておくのが便利です。

#3. 平面波カットオフ
Quantum ESPRESSOでは、波動関数と電荷密度を平面波で展開します。それぞれのカットオフエネルギーを以下のように設定します。

&SYSTEM
    ecutwfc = 25.0
    ecutrho = 225.0
/

ecutwfcが波動関数のカットオフエネルギー、ecutrhoが電荷密度のカットオフエネルギーです。単位はRyです。使用する平面波の運動エネルギーの最大値が、カットオフエネルギーです。

カットオフエネルギーの具体的な値は、使用する擬ポテンシャルに依存します。擬ポテンシャルファイルの中にカットオフエネルギーの推奨値が記載されている場合には、その値を使ってください。推奨値の記載ばない場合には、ユーザーが自分で適当な値を探す必要があります。しかしながら、大抵の場合には ecutwfc = 25.0 Ry 程度を設定しておけば、問題ない場合が多いように思えます。（定量的に厳密な計算をしたい場合には、必ず、ecutwfcを少しずつ大きくして 全エネルギーの変化などを確認して下さい。）

使用する擬ポテンシャルファイルが全てノルム保存型の場合には、ecutrhoの設定は不要です。ecutrho = 4 * ecutwfc が自動的に設定されます。少なくとも一つのウルトラソフト型擬ポテンシャルまたはPAWが含まれる場合には、ecutrhoの設定が必須となります。擬ポテンシャルファイルの中に推奨値の記載がない場合には、ecutrho = 9 * ecutwfc 程度を設定しておけば問題ないと思います。

#4. 占有数
各バンドに対する電子の占有数を設定します。フェルミの排他律によれば、占有数は０または２のいずれかです(２はスピンの自由度)。系が絶縁体または半導体の場合には、下のバンドから順番に電子数分だけ占有数２として、それより上のバンドの占有数は０とします。このような場合には、入力ファイルを以下のように設定します。

&SYSTEM
    occupations = "fixed"
/

系が金属の場合には、単純に占有数を０，２とすることは出来ません。正確に言うと、有限個のk点(後述)を用いているためフェルミ面がガタガタになります。フェルミ面がバンドギャップ内にある絶縁体と半導体の場合には このガタガタは問題とならないのですが、フェルミ面がバンド内部にある金属の場合には重大な問題となります。そこで、金属の場合には、このガタガタに対してアンチエイリアス処理を施して、精度の向上を図ります。これがsmearingです。以下のように設定します。

&SYSTEM
    occupations = "smearing"
    smearing    = "gaussian"
    degauss     = 0.01
/

smearingの方法にはいくつかのバリエーションがあるのですが、ガウス関数法(smearing = "gaussian")が最もポピュラーです。ガウス関数を使って、フェルミ面近傍のバンドに対する占有数を０～２の中間値とします。使用するガウス関数の幅は、変数degaussにて設定します。単位は、Ryです。degaussの値としては、0.0数Ry程度に設定するのが一般的です。degaussを有限値とすることは、電子に微小な温度分布を持たせることに等価です。このため、金属系の全エネルギーには、当該の温度分布に対応したエンタルピー項が加算されます。たいていの場合、k点サンプリングが十分であれば、degaussの値の任意性はエンタルピー項に吸収されてしまいます。つまり、degaussが少々変化しても、エンタルピー項を含めた全エネルギーはほとんど変化しないわけです。ただし、磁性などのフェルミ面の構造に敏感な物性を扱う場合に限っては、degaussの値を微調整(最適化)する必要があります。

#5. 収束条件
SCF計算の収束に関する設定をします。収束したか否かの判定は、全エネルギーの推定誤差を用います。推定誤差が指定された閾値以下となった場合に、SCF計算が収束したと判定されます。閾値は以下の変数で設定します。単位は、Ryです。

&ELECTRONS
    conv_thr = 1.0e-06
/

SCF計算が収束しない場合、特に全エネルギーが振動する場合には、電荷密度の混合率を小さめに設定すると良いです。以下のように設定します。

&ELECTRONS
    mixing_beta = 0.4
    mixing_mode = "plain"
/

0.0 < mixing_beta < 1.0 の範囲の値が使用できます。スピン分極がある場合や強相関系では特に有効です。これらの系では、0.2 ～ 0.4 程度を設定するのがお勧めです。また、mixing_modeにて、電荷密度混合法のアルゴリズムを選択可能です。デフォルトでは、"plain"です。"plain"でSCF計算が収束しない場合には、"TF"または"local-TF"に変更すると良いです。計算対象が金属や半導体の場合には、"TF"が良いです。金属酸化物などのバンドギャップの大きい絶縁体や分子系では、"local-TF"が適しています。

#6. k点サンプリング
Brillouinゾーン内のk点サンプリングを設定します。k点サンプリングは、入力ファイルのK_POINTSという項目で設定します。SCF計算の場合には、一様メッシュでサンプリングしますので、{automatic}オプションを適用します。以下のように書きます。

K_POINTS {automatic}
 4  4  2  0 0 0

２列目にて、k点のメッシュを指定しています。１列目、２列目、３列目は、それぞれ、X方向のメッシュ数、Y方向のメッシュ数、Z方向のメッシュ数を意味します。４列目、５列目、６列目は、k点のシフトを意味するのですが、特段の意図が無ければ３つとも0として問題ないです。k点メッシュ数は、格子定数に反比例するように設定します。具体的な値は目的の如何ではありますが、概ね、「格子定数」x 「k点メッシュ数」＝ 10 ~ 12 Å となるように選べば十分な精度が得られます。

---

次回は、擬ポテンシャルについて解説したいと思います。
なお、Winmostarを利用するとQuautum ESPRESSOをGUI上から簡単かつ高度に活用することができます。学生は無料で利用でき、学生以外も無料トライアルを入手することができますので、是非 WinmostarのWebサイト をご覧ください。

Quantum ESPRESSO入力ファイル作成手順シリーズ
１．結晶構造の作成
２．SCF計算の設定
３．擬ポテンシャルファイルの選択方法
４．原子座標のみの最適化
５．格子ベクトル及び格子内部の原子座標の最適化
６．状態密度の計算
７．局所状態密度の計算
８．バンド構造の計算
９．バンド数の設定
１０．van der Waals相互作用
１１．LDA+U法
１２．NEB法
１３．Phonon計算（特定q点）
１４．Phonon計算（バンド構造）
１５．擬ポテンシャルの作成
１６．擬ポテンシャルのテスト
１７．SCF計算の入力ファイルの実例