Atsushi M. Ito web

DFTの計算ではカットオフエネルギーとかカットオフ波数というものを設定します。入力ファイルにおいてOpenMXではカットオフエネルギーを設定し、xTAPPではカットオフ波数を設定することになります。これらはバンド理論における逆格子ベクトルと関連していますが、数値計算の観点からラフに見ると、実空間グリッドのグリッド数を指定していると言えます。その様に見た方が、数値計算に慣れている人にとっては親しみやすいかと思いますので、ここで少し解説します。

カットオフエネルギーとカットオフ波数

カットオフエネルギーをE_c、カットオフ波数をk_cと書いたとき、両者の間には

k_c^2=E_c \textrm{[Ry]}

という関係があります。ここで1 Ry=27.21138 eVです。

xTAPPの場合には波数を指定します。k_c=7.0を指定した場合、カットオフエネルギーとしては49Ryを指定したことに相当します。

グリッド数への変換 (xTAPPの場合)

ここで、xTAPPの場合にはカットオフ波数k_cとして設定したい値は入力ファイルのTAPPINPUT.CUTOFF_WAVE_FUNCTIONにセットします。加えて、xTAPPにはTAPPINPUT.CUTOFF_LOCAL_POTENTIALというもう一つのカットオフ波数が設定できます。この局所ポテンシャルカットオフ波数をk_lとします。TAPPINPUT.CUTOFF_LOCAL_POTENTIALの値を設定する場合は、

k_l > 2k_c

を満たすように設定する必要があります。また、TAPPINPUT.CUTOFF_LOCAL_POTENTIALの値を設定しなかった場合は、デフォルトでk_l=2k_c+0.1と設定されるようです。xTAPPでは局所ポテンシャルカットオフ波数k_lを利用してのグリッド数を算出します。

まず、計算するシミュレーションボックス(システムサイズ)の3つの格子ベクトルを\vec{a}_x,\vec{a}_y,\vec{a}_zとします。長さの単位はBohrです。ここから逆格子ベクトル(reciprocal lattice vector)\vec{b}_x,\vec{b}_y,\vec{b}_zを次の様に定義します。

\vec{b}_x = 2\pi \frac{ \vec{a}_y \times \vec{a}_z}{V}

ただし、V=\vec{a}_x \cdot ({ \vec{a}_y \times \vec{a}_z})はセル体積(シミュレーションボックスの体積)です。\vec{b}_y,\vec{b}_zに関しても同様に定義されます。

ここで、格子ベクトルや逆格子ベクトルが必ずしも直交していない場合があるので、直交成分のみを次のように抽出します。

l_x = \vec{b}_x \cdot \frac{\vec{b}_y \times \vec{b}_z}{\left|\vec{b}_y \times \vec{b}_z \right|}

これを用いて実空間グリッド数n_x,n_y,n_zは

n_x \ge \frac{k_l}{|l_x|}

を満たす整数として算出されます。また、実際上はグリッド数n_xが2,3,5の公倍数になるように決めているようです。

また、より簡単な見積もりとして、格子ベクトルが全て直交している場合を考えます。すると、

|l_x|=|\vec{b}_x|=\frac{2\pi}{ |\vec{a}_x|}

ですので、

n_x \ge \frac{ k_l |\vec{a}_x| }{2\pi}

となります。

グリッド数への変換 (OpenMXの場合)

OpenMXの場合にはscf.energycutoffによってカットオフエネルギーE_c[Ry]を指定します。上述のxTAPPではカットオフ波数k_cの2倍となるk_lを用いたように、OpenMXでも

k_l = 2 \sqrt{E_c}

とします。あとはxTAPPの場合と同様に上述の計算から実空間のグリッド数を算出します。

OpenMXにおいてもグリッド数n_xが2,3,5,7の公倍数となるように設定されます。ただし、OpenMXとxTAPPではグリッド数の整数化に伴う計算手順の違いから、同じカットオフエネルギー(波数)を設定しても、実際に利用されるグリッド数は多少異なります。

OpenMXでは、整数化されたグリッド数を求めた後に、そこから逆算して実際のカットオフエネルギーを算出し、標準出力に出力してくれていますので、計算の確認に使ってください。

実空間グリッドはどこで使うのか

平面波基底や原子基底を使っているはずなのに、実空間になぜグリッドが必要なのでしょうか。OpenMXやxTAPPでは電子密度を一度実空間グリッド状のデータとして求めています。そこからフーリエ変換を駆使してHartreeポテンシャルや相関交換ポテンシャルを求めています。おそらく多くのDFTコードが電子密度に関しては実空間グリッド状のデータとして計算を行っています。例えば平面波の場合には全て波数空間で計算を行うことも不可能ではないですが、電子密度の算出およびそこからのポテンシャル算出の際に非常に膨大な計算量が必要になります。フーリエ変換を用いて一度実空間に焼き直すことで、計算量を大幅に縮小できるのです。そもそも、ここでの説明ではグリッド数を「実空間の」という枕詞を付けましたが、波数空間においても波数に関するグリッドが必要で、離散化されています。つまり、平面波を基底として「何個」使うかという数が、フーリエ変換を通して実空間のグリッド数に対応します。非常にざっくりとした説明ですので、詳細は各種文献等を参考にしてください。