* アルミニウム結晶 [#wa27a0e1] このチュートリアルでは金属系に関するk点およびスメアリング関数を用いてフェルミ面を取り扱う際のぼかし幅について収束性のテストに関する説明を行います。 バンド構造および状態密度の計算方法も説明します。 ** 収束性のテスト [#kd001903] 金属系では一般的に半導体・絶縁体に比べて多くのk点が必要になります。 さらにSTATEではフェルミ面を取り扱うために様々な関数系によるスメアリング法、あるいはテトラへドロン法の使用が可能です。 *** k点に関する収束性の調査 [#lda8499b] ここではMethfessel-Paxton法とテトラへドロン法を用いて全エネルギーのk点に関する収束を調べます。 スメアリング法を使用するためには負のぼかし幅 (WIDTH)を指定します。ぼかし幅の単位はRyです。 正のWIDTHを使用すると多項式のスメアリング関数が使用され、電子のエントロピーの項は考慮されません。 テトラへドロン法を使用する場合にはWIDTHを-10よりも小さい任意の数値にします。 - スメアリング関数(Methfessel-Paxton法)を使用する場合の入力ファイル 0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY) 4.00 8.00 1 1 1 : GMAX GMAXP NTYP NATM NATM2 221 2 : NUM_SPACE_GROUP TYPE 7.5967 7.5967 7.5967 90.0 90.0 90.0 : A B C ALPHA BETA GAMMA 06 06 06 1 1 1 : N1 N2 N3 M1 M2 M3 0 0 : NCORD, NINV 0.00 0.00 0.00 1 0 1 : CPS(1,1:3) IWEI IMDTYP ITYP 13 0.50 26.98 6 1 0.2 : IATOMN ALFA AMION ILOC IVAN 0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC 0 1 : IPRE IPRI 30 30 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP 6 1 : WAY_MIX MIX_WHAT 0 20 0.60 : ITER_START KBXMIX MIX_ALPHA 0.20 0.30 0.20 0.20 0.20 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM 300.00 4 1 0.50D-09 : DTIO IMDALG IEXPL EDELTA -0.0020 0.50D+03 0 : WIDTH FORCCR ISTRESS ggapbe 1 : XCTYPE KSPIN 2.00 : DESTM 101 : NBZTYPE 4 4 4 : NKX NKY NKZ (DUMMY) 4 4 4 : NKX2 NKY2 NKZ2 (DUMMY) 6 : KEG 1 : NEXTST 0 : (DUMMY) 2 : IMSD 0 : EVALUATE_EKO_DIFF 0 : NPDOSAO 0 0.000 : SM_N DOPPING (DUMMY) - テトラへドロン法を用いる場合の入力ファイル 0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY) 4.00 8.00 1 1 1 : GMAX GMAXP NTYP NATM NATM2 221 2 : NUM_SPACE_GROUP TYPE 7.5967 7.5967 7.5967 90.0 90.0 90.0 : A B C ALPHA BETA GAMMA 06 06 06 1 1 1 : N1 N2 N3 M1 M2 M3 0 0 : NCORD, NINV 0.00 0.00 0.00 1 0 1 : CPS(1,1:3) IWEI IMDTYP ITYP 13 0.50 26.98 6 1 0.2 : IATOMN ALFA AMION ILOC IVAN 0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC 0 1 : IPRE IPRI 30 30 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP 6 1 : WAY_MIX MIX_WHAT 0 20 0.60 : ITER_START KBXMIX MIX_ALPHA 0.20 0.30 0.20 0.20 0.20 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM 300.00 4 1 0.50D-09 : DTIO IMDALG IEXPL EDELTA -10.02 0 0.50D+03 0 : WIDTH FORCCR ISTRESS ggapbe 1 : XCTYPE KSPIN 2.00 : DESTM 101 : NBZTYPE 4 4 4 : NKX NKY NKZ (DUMMY) 4 4 4 : NKX2 NKY2 NKZ2 (DUMMY) 6 : KEG 1 : NEXTST 0 : (DUMMY) 2 : IMSD 0 : EVALUATE_EKO_DIFF 0 : NPDOSAO 0 0.000 : SM_N DOPPING (DUMMY) 上記のファイルの違いはWIDTHの値だけです。~ k点メッシュを変えながら計算を実行すると、k点の関数として全エネルギーを以下のように計算できます。 #ref(http://www-cp.prec.eng.osaka-u.ac.jp/puki_state/graph/etot_nk_Al_fcc.png,center,nolink,30%) *** スメアリング幅に関する収束性のテスト [#pbe549d0] 金属系の全エネルギー(電子のエントロピーを含む自由エネルギー)はスメアリング関数のぼかし幅に敏感です。 特にガウス関数とフェルミ-ディラック関数の場合は顕著です。 以下ではde Gironcoliによって最初に示された全エネルギーのぼかし幅 依存性をアルミニウムについて確認します[Phys. Rev. B51, 6773(R) (1995)]。Marzariによる[[記事>http://theossrv1.epfl.ch/Main/ElectronicTemperature]]も大変参考になります。 デフォルトのMethfessel-Paxton (一次のHermite-Gaussian関数) 以外のスメアリング関数を使用する場合は&OTHERS ... &ENDというブロックを入力ファイルに追加し、使用したいスメアリング関数名を指定します。 フェルミ-ディラック分布関数の場合 &OTHERS FERMI_DIRAC &END ガウシアン関数の場合 &OTHERS GAUSSIAN &END MarzariとVanderbiltのcold smearingを使用する場合 &OTHERS COLD &END とします。 全エネルギーをぼかし幅の関数として計算すると以下が得られます。 #ref(http://www-cp.prec.eng.osaka-u.ac.jp/puki_state/graph/etot_sigma_Al_fcc.png,center,nolink,30%) 先に述べたようにフェルミ-ディラック関数とガウス関数を使用した場合、全エネルギーのぼかし幅依存性は著しく、Hermite-Gaussianあるいはcold smearingを使用した場合、全エネルギーはぼかし幅にほぼ依存性しないことが分かります。 一般に我々の計算は0Kで行われており、ぼかし幅0の極限での全エネルギーが我々の求めたい値と考えられます。 そのような目的のもとでHermite-Gaussianあるいはcold smearingが適切な選択肢と考えられます。 また一般にぼかし幅が大きければk点サンプリング数は少なくて済みます(k点数に関する収束性のテストは常に必要)。 ぼかし幅とk点数のバランスをうまくとり、効率的に計算を進めることができるパラメーターを選ぶ必要があります。 * バンド計算 [#c7d42ba2] ここではバンド構造の計算方法を説明します。 バンド計算を行うためには、まず最初に電子密度をSCF計算により求める必要があります。 以下はテスト用に比較的少ないk点数を用いたSCF計算用の入力ファイルの例を示します。 - 入力ファイル (nfinp_scf) 0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY) 4.00 8.00 1 1 1 : GMAX GMAXP NTYP NATM NATM2 221 2 : NUM_SPACE_GROUP TYPE 7.5967 7.5967 7.5967 90.0 90.0 90.0 : A B C ALPHA BETA GAMMA 06 06 06 1 1 1 : N1 N2 N3 M1 M2 M3 0 0 : NCORD, NINV 0.00 0.00 0.00 1 0 1 : CPS(1,1:3) IWEI IMDTYP ITYP 13 0.50 26.98 6 1 0.2 : IATOMN ALFA AMION ILOC IVAN 0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC 0 1 : IPRE IPRI 30 30 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP 6 1 : WAY_MIX MIX_WHAT 0 20 0.60 : ITER_START KBXMIX MIX_ALPHA 0.20 0.30 0.20 0.20 0.20 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM 300.00 4 1 0.50D-09 : DTIO IMDALG IEXPL EDELTA -0.0020 0.50D+03 0 : WIDTH FORCCR ISTRESS ggapbe 1 : XCTYPE KSPIN 2.00 : DESTM 101 : NBZTYPE 4 4 4 : NKX NKY NKZ (DUMMY) 4 4 4 : NKX2 NKY2 NKZ2 (DUMMY) 6 : KEG 1 : NEXTST 0 : (DUMMY) 2 : IMSD 0 : EVALUATE_EKO_DIFF 0 : NPDOSAO 0 0.000 : SM_N DOPPING (DUMMY) 次にSCF計算により得られた電荷密度を用いて、(電荷密度を固定した)非SCF (non-SCF) 計算を実行します。 バンド計算においてはブリルアンゾーン内の高対称点をつなぐ線上のk点を選びKohn-Sham方程式を解きます。 以下はnon-SCFバンド計算用の入力ファイルの例です。 - 入力ファイル (nfinp_band) 0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY) 4.00 8.00 1 1 1 : GMAX GMAXP NTYP NATM NATM2 221 2 : NUM_SPACE_GROUP TYPE 7.5967 7.5967 7.5967 90.0 90.0 90.0 : A B C ALPHA BETA GAMMA 06 06 06 1 1 1 : N1 N2 N3 M1 M2 M3 0 0 : NCORD, NINV 0.00 0.00 0.00 1 0 1 : CPS(1,1:3) IWEI IMDTYP ITYP 13 0.50 26.98 6 1 0.2 : IATOMN ALFA AMION ILOC IVAN 22 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC 0 1 : IPRE IPRI 30 30 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP 6 1 : WAY_MIX MIX_WHAT 0 20 0.60 : ITER_START KBXMIX MIX_ALPHA 0.20 0.30 0.20 0.20 0.20 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM 300.00 4 1 0.50D-09 : DTIO IMDALG IEXPL EDELTA -0.0020 0.50D+03 0 : WIDTH FORCCR ISTRESS ggapbe 1 : XCTYPE KSPIN 2.00 : DESTM 101 : NBZTYPE 4 4 4 : NKX NKY NKZ (DUMMY) 4 4 4 : NKX2 NKY2 NKZ2 (DUMMY) 6 : KEG 1 : NEXTST 0 : (DUMMY) 2 : IMSD 0 : EVALUATE_EKO_DIFF 0 : NPDOSAO 0 0.000 : SM_N DOPPING (DUMMY) &KPOINTS_BAND NKSEG 4 KMESH 40 20 20 20 KPOINTS 0.000 0.000 0.000 0.000 0.500 0.500 0.250 0.500 0.750 0.500 0.500 0.500 0.000 0.000 0.000 &END バンド計算を実行する時は以下のようにICONDを22に指定します。 22 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC ICONDに加えて以下のように&KPOINTS_BAND ... &ENDというブロックを追加し、バンドを計算するブリルアンゾーン内の線を定義する高対称点を指定します。 &KPOINTS_BAND NKSEG 4 KMESH 40 20 20 20 KPOINTS 0.000 0.000 0.000 0.000 0.500 0.500 0.250 0.500 0.750 0.500 0.500 0.500 0.000 0.000 0.000 &END ここで NKSEG はk点のセグメント(区分、線)の数を指定します。 各セグメントの刻み幅は KMESH でセグメントの数だけ指定します。 高対称点は KPOINTS というキーワードに続いて指定します。 指定するk点の数はNKSEG+1となります。~ この入力ファイルを用いてnon-SCFバンド計算を実行すると、出力ファイルに加えてk点と固有値が"energy.data"というファイルに出力されます。このデータは直接プロットすることができません。 バンドをプロットするために"energy2band"というプログラムが用意されています。~ "energy2band"は以下のように実行します。 > energy2band $ energy2band 計算に用いられたバンドの数、プロットするバンドの数、バンド計算の際に使用されたk点の数、エネルギーの基準(Hartree単位、一般にはFermi levelを指定)が聞かれるのでそれらを正しく入力すると"band.data"というファイルが生成されます。 これはgnuplotやxmgraceでプロット可能です。 以下はgnuplotスクリプトの例です。 set terminal postscript eps 'Helvetica' 22 set output 'band.eps' emin=-12.5 emax=12.5 xmin=0.0 xmax=2.54176932 G=0.0 X=0.82709404 W=1.24064106 L=1.82548487 G2=2.54176932 offset=0.75 set xrange [xmin:xmax] set yrange [emin:emax] unset key unset xtics unset xlabel set ylabel 'E-E_F (eV)' set xzeroaxis set arrow from X,emin to X,emax nohead set arrow from W,emin to W,emax nohead set arrow from L,emin to L,emax nohead set label '{/Symbol G}' at G,emin-offset center set label 'X' at X,emin-offset center set label 'W' at W,emin-offset center set label 'L' at L,emin-offset center set label '{/Symbol G}' at G2,emin-offset center plot 'band.data' using 1:2 with lines lw 3 これを"band.gp"というファイルに保存して $ gnuplot band.gp を実行すると"band.eps"というファイルが得られます。 #ref(http://www-cp.prec.eng.osaka-u.ac.jp/puki_state/graph/band_Al_fcc.png,,center,nolink,30%) * 練習問題 [#lf33eb74] - 格子定数のk点およびスメアリング幅依存性を調べてみましょう。またスメアリング法とテトラへドロン法を用いて得られた格子定数を比較してみましょう。