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

シリコン結晶 (SCF計算) †

ここでは結晶の計算例としてダイアモンド構造のシリコンのSCF計算を行いましょう。
現在のワーキングディレクトリが~/STATE/examplesであるとして、先ずはSiディレクトリに移動しましょう。

$ cd Si

ここでlsを実行し、nfで始まる入力ファイル群（nf*）や.shで終わるジョブスクリプト（*.sh）が含まれていることを確認します。

以下を実行して入力ファイルnfinp_scfを見てみましょう。

$ cat nfinp_scf

#
# Crystalline silicon in the diamond structure
#
WF_OPT    DAV
NTYP      1
NATM      2
TYPE      2
NSPG      227
GMAX      4.00
GMAXP     8.00
KPOINT_MESH    8   8   8
KPOINT_SHIFT   OFF OFF OFF
WIDTH     0.0002
EDELTA    0.5000D-09
NEG    8
CELL   10.30  10.30  10.30  90.00  90.00  90.00
&ATOMIC_SPECIES
 Si 28.0900 pot.Si_pbe1
&END
&ATOMIC_COORDINATES CRYSTAL
      0.000000000000      0.000000000000      0.000000000000    1    1    1
      0.250000000000      0.250000000000      0.250000000000    1    1    1
&END

以下を実行するとシリコン結晶の構造を表示するためのSi2.xsfが生成されます。

$ chkinpf nfinp_scf

XCrySDenやVESTAを使用して自分が計算した結晶の構造となっていることを確認します。 同様にsb100で使用するジョブスクリプトqsub_sb100.shも確認しましょう。

$ cat qsub_sb100.sh

#$ -S /bin/sh
#$ -cwd
#$ -q sb.q
#$ -pe x6 6
#$ -N Si

module load intel/2021.2.0
module load intelmpi/2021.2.0
 
# Disable OPENMP parallelism
 
export OMP_NUM_THREADS=1
 
unset I_MPI_TCP_NETMASK
 
# Set the execuable of the STATE code
 
ln -fs ${HOME}/STATE/src/state-5.6.9/src/STATE .
 
# Set the pseudopotential data
 
ln -fs ../gncpp/pot.Si_pbe1
  
# Set the input/output file
 
INPUT_FILE=nfinp_scf
OUTPUT_FILE=nfout_scf
 
# Run!
 
mpirun ./STATE < ${INPUT_FILE} > ${OUTPUT_FILE}
 

このジョブスクリプトを使ってジョブを投入します。

$ qsub qsub_sb100.sh

計算の状態は

$ qstat -u [username]

を実行することで確認できます。
この例では出力ファイルの名前をnfout_scfとしています。計算が終わったらlessなどのページャーを使ってその内容を確認してみましょう。

$ less nfout_scf

格子ベクトルは以下のように表示されます。

 PRIM. LAT. VECTOR(BOHR) :        0.000000       5.150000       5.150000
 PRIM. LAT. VECTOR(BOHR) :        5.150000       0.000000       5.150000
 PRIM. LAT. VECTOR(BOHR) :        5.150000       5.150000       0.000000

原子位置は以下のように表示されます。

 ********************************* ATOMS *******************************
   ATOM    X(BOHR)    Y(BOHR)    Z(BOHR)     TAUX    TAUY    TAUZ IW  IR
   1  1    0.00000    0.00000    0.00000   0.0000  0.0000  0.0000  1   0
   2  1    2.57500    2.57500    2.57500   0.2500  0.2500  0.2500  1   0
 ***********************************************************************

計算で使用するバンドの数は以下のように表示されます。

 NUMBER OF BANDS CONSIDERED :        8

交換相関汎関数は以下のように表示されます。

EXCHANGE CORRELATION FUNCTIONALS : ggapbe

SCF計算中のエネルギーは以下のように表示されます。

 ***********************************************************************
 *                                                                     *
 *                              START SCF                              *
 *                                                                     *
 ***********************************************************************

  NSCF NADR            ETOTAL          EDEL          CDEL CONV      TCPU
     1    0       -6.05513096   0.60551E+01   0.32033E-02    0      0.21
     2    0       -7.84013758   0.17850E+01   0.50625E-02    0      0.11
     3    1       -7.87280128   0.32664E-01   0.23361E-02    1      0.11
     4    2       -7.87352825   0.72697E-03   0.44450E-03    2      0.11
     5    3       -7.87355444   0.26191E-04   0.14565E-03    2      0.11
     6    4       -7.87355827   0.38304E-05   0.17592E-04    3      0.11
     7    5       -7.87355833   0.59844E-07   0.53913E-05    3      0.11
     8    6       -7.87355833   0.17526E-08   0.14864E-05    3      0.11
     9    7       -7.87355833   0.15877E-09   0.19139E-06    4      0.11
    10    8       -7.87355833   0.12939E-10   0.18081E-07    5      0.11
    11    9       -7.87355833   0.84555E-12   0.23101E-07    6      0.11

計算が収束した場合には以下のように全エネルギーとその成分が表示されます。

                     TOTAL ENERGY AND ITS COMPONENTS
                  TOTAL ENERGY     =          -7.87355833 A.U.
                KINETIC ENERGY     =           3.01922419 A.U.
                HARTREE ENERGY     =           0.55014198 A.U.
                     XC ENERGY     =          -2.40098652 A.U.
                  LOCAL ENERGY     =          -0.84294926 A.U.
               NONLOCAL ENERGY     =           0.16885291 A.U.
                  EWALD ENERGY     =          -8.36784162 A.U.
                     PC ENERGY     =           0.00000000 A.U.
               ENTROPIC ENERGY     =           0.00000000 A.U.

計算が終わると波動関数データzaj.data、ポテンシャル（電荷密度）データpotential.dataが出力されます。 それに加えて状態密度のデータを含むdos.dataが出力されます。これはgnuplotなどで可視化することが可能です。