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

#0. はじめに
第一原理計算ソフトQuantunm ESPRESSO https://www.quantum-espresso.org の入力ファイルの書き方について、複数回にわたって解説していきます。今回は第１６回目です。前回は擬ポテンシャルファイルの作成方法を説明しましたが、今回は作成した擬ポテンシャルファイルのテスト方法について解説します。

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

#1. 概要
擬ポテンシャルファイルは、孤立原子に対するDFT計算の結果を用いて作成されます。このため、作成した擬ポテンシャルは孤立原子系に対しては厳密に正しい結果を導きます。しかしながら、現実の結晶系や分子系においては、各々の原子は孤立系とは異なる電子配置を持ちます。ノルム保存型、ウルトラソフト型、PAW型などの擬ポテンシャル法では、孤立原子系からの電子配置のズレを補外する機構が備わっています。ただし、当該機構にも限界があり、電子配置が大きく変化する場合(大きな価数を持つイオンなど)については事前に擬ポテンシャルファイルの精度を検証する必要があります。

Quantum ESPRESSOには、擬ポテンシャルファイルをテストする機能が実装されています。当該機能の実行には、擬ポテンシャルファイル作成に用いた実行ファイルld1.xを使います。擬ポテンシャルファイルの作成またはテストのいずれを実行するかは、入力ファイルにて指示します。

#2. 入力ファイル
擬ポテンシャルファイルのテストで使用する入力ファイルのサンプルを、以下に示します。Al原子の例です。

&INPUT
  zed     = 13.0
  rel     = 1
  config  = '[Ne] 3s2.0 3p1.0 3d0.0'
  iswitch = 2
  dft     = 'PBE'
/

&TEST
  nconf       = 4
  file_pseudo = 'Al.UPF'
/
2
3S  1  0  2.00  0.00  2.60  2.60  1
3P  2  1  1.00  0.00  2.90  2.90  1
2
3S  1  0  2.00  0.00  2.60  2.60  1
3P  2  1  0.00  0.00  2.90  2.90  1
2
3S  1  0  1.00  0.00  2.60  2.60  1
3P  2  1  0.00  0.00  2.90  2.90  1
2
3S  1  0  0.00  0.00  2.60  2.60  1
3P  2  1  0.00  0.00  2.90  2.90  1

&INPUT の内容は、擬ポテンシャルファイルを作成する場合と同じ設定ですhttps://qiita.com/xa_member/items/5feb8f432170141c89cd。ただし、iswitchは2に設定する必要があります。

擬ポテンシャルファイルの作成では &INPUTP というブロックがありました。擬ポテンシャルファイルのテストを行う際には、代わりに &TEST というブロックを設定します。file_pseudo という項目には、テスト対象の擬ポテンシャルファイルの名前を設定します。nconf には、テストする電子配置の数を設定します。&TEST ブロックより下の行に各電子配置を設定します。この場合は nconf = 4 なので、４つの電子配置が設定されています。電子配置の記述方法は、擬ポテンシャルファイル作成時(iswitch = 3)と同様です。

上記の例では、中性のAl原子、＋１価のカチオン、＋２価のカチオン、＋３価のカチオンの電子配置が設定されています。

#3. 出力ファイル
擬ポテンシャルファイルのテスト結果は、標準出力に表示されます。「Testing the pseudopotential」という行の以降に、データが記述されます。

     ---------------------- Testing the pseudopotential ----------------------

                                                                                
     scalar relativistic calculation

     atomic number is 13.00   valence charge is  3.00
     dft =PBE   lsd =0 sic =0 latt =0  beta=0.20 tr2=1.0E-14
     mesh =1135 r(mesh) = 100.48924 xmin = -7.00 dx = 0.01250

     n l     nl             e AE (Ry)        e PS (Ry)    De AE-PS (Ry) 
     1 0     3S   1( 2.00)       -0.56982       -0.56982        0.00000
     2 1     3P   1( 1.00)       -0.19934       -0.19934        0.00000

上記は、中性のAl原子に対する計算結果です。下から2番目の行がs軌道に対するエネルギー、最後の行がp軌道のエネルギーです。「e AE (Ry)」および「e PS (Ry)」が、それぞれ、Kohn-Sham方程式を厳密に解いて得られたエネルギーおよび擬ポテンシャルファイルで計算されたエネルギーです。「De AE-PS (Ry)」は両エネルギーの差分です。中性原子の場合には、AE と PS の結果が完全に合致しています。

カチオン状態での計算結果は、以下の通りです。

＋１価イオン　：

     n l     nl             e AE (Ry)        e PS (Ry)    De AE-PS (Ry) 
     1 0     3S   1( 2.00)       -1.08854       -1.08850       -0.00004
     2 1     3P   1( 0.00)       -0.66543       -0.66566        0.00023

＋２価イオン　：

     n l     nl             e AE (Ry)        e PS (Ry)    De AE-PS (Ry) 
     1 0     3S   1( 1.00)       -1.74543       -1.73731       -0.00812  !
     2 1     3P   1( 0.00)       -1.26440       -1.26049       -0.00391

＋３価イオン　：

     n l     nl             e AE (Ry)        e PS (Ry)    De AE-PS (Ry) 
     1 0     3S   1( 0.00)       -2.40341       -2.07542       -0.32800  !
     2 1     3P   1( 0.00)       -1.85342       -1.59191       -0.26151  !

＋１価イオンの場合は比較的に精度が良く、＋２および＋３価イオンでは AE と PS の差異が大きくなっています。特に誤差が大きい場合には、軌道エネルギーの行の右側に「!」が自動的に付きます。誤差が過大の場合には、イオン化された電子配置での擬ポテンシャルファイルの作成を検討する必要があります。

#4. ゴースト準位
また、ld1.xには、作成した擬ポテンシャルでゴースト準位が発生するか否かをチェックする機能も実装されています。具体的には、カットオフエネルギーを設定して、Bessel関数を基底として擬ポテンシャル系のKohn-Sham方程式を解きます。このときのエネルギー固有値が、前章にて計算された「e AE (Ry)」および「e PS (Ry)」に合致していればゴースト準位が存在しないことが保証されます(数値誤差程度の差異は可)。エネルギー固有値が大きく異なっている場合にはゴースト準位が存在するので、当該の擬ポテンシャルファイルはSCF計算にて利用できません。ゴースト準位チェック用の入力ファイルは、以下の通りです。

&INPUT
  zed     = 13.0
  rel     = 1
  config  = '[Ne] 3s2.0 3p1.0 3d0.0'
  iswitch = 2
  dft     = 'PBE'
/

&TEST
  nconf       = 1
  file_pseudo = 'Al.UPF'

  ecutmin = 50
  ecutmax = 200
  decut   = 50
/
2
3S  1  0  2.00  0.00  2.60  2.60  1
3P  2  1  1.00  0.00  2.90  2.90  1

基本的な内容は、先述の入力ファイルと同じです。ただし、ecutmin, ecutmax, decut が追加されています。それぞれ、Bessel関数のカットオフエネルギーの最小値、最大値、刻み幅です。この場合、50Ry, 100Ry, 150Ry, 200Ry の4パターンのカットオフエネルギーにてテスト計算が実施されることを意味します。この入力ファイルで計算した結果は、以下の通りです。

     -------------- Test with a basis set of Bessel functions ----------

     Box size (a.u.) :   30.0

     Cutoff (Ry) :   50.0
                           N = 1       N = 2       N = 3
     E(L=0) =        -0.5695 Ry   -0.0252 Ry    0.0261 Ry
     E(L=1) =        -0.1991 Ry    0.0110 Ry    0.0495 Ry
     E(L=2) =         0.0282 Ry    0.0611 Ry    0.1193 Ry

     Cutoff (Ry) :  100.0
                           N = 1       N = 2       N = 3
     E(L=0) =        -0.5695 Ry   -0.0252 Ry    0.0261 Ry
     E(L=1) =        -0.1991 Ry    0.0110 Ry    0.0495 Ry
     E(L=2) =         0.0282 Ry    0.0611 Ry    0.1193 Ry

     Cutoff (Ry) :  150.0
                           N = 1       N = 2       N = 3
     E(L=0) =        -0.5695 Ry   -0.0252 Ry    0.0261 Ry
     E(L=1) =        -0.1991 Ry    0.0110 Ry    0.0495 Ry
     E(L=2) =         0.0282 Ry    0.0611 Ry    0.1193 Ry

     Cutoff (Ry) :  200.0
                           N = 1       N = 2       N = 3
     E(L=0) =        -0.5695 Ry   -0.0252 Ry    0.0261 Ry
     E(L=1) =        -0.1991 Ry    0.0110 Ry    0.0495 Ry
     E(L=2) =         0.0282 Ry    0.0611 Ry    0.1193 Ry

     -------------- End of Bessel function test ------------------------

3Sおよび3P軌道ともに、数値誤算の範囲内で「e AE (Ry)」および「e PS (Ry)」に合致しており、ゴースト順位が存在しないことが確認された。

なお、Winmostarを利用するとQuautum ESPRESSOをGUI上から簡単かつ高度に活用することができます。学生は無料で利用でき、学生以外も無料トライアルを入手することができますので、是非 WinmostarのWebサイトをご覧ください。

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