この節では炭化チタンTiCを例に電子状態計算を行います。
* [[https://ja.wikipedia.org/wiki/%E7%82%AD%E5%8C%96%E3%83%81%E3%82%BF%E3%83%B3|炭化チタン - Wikipedia]]
第一原理計算の基本的な流れは
- 結晶構造データの用意
- 計算方法の設定
- 実際の電子状態計算（self-consistent field計算、SCF計算）
- 各種物理量の計算
となります。
ここでは3番目までの方法を解説します。
====== SCF計算 ======
===== 結晶構造データの作成 =====
WIEN2kでは結晶構造はcase.structというファイルに保存されます。
caseはプロジェクト名であり自由に名付けることができますが、今回はTiC.structというファイルを作りましょう。
作業ディレクトリTiCを作成し、その中で作業します。
作業ディレクトリ名は必ずプロジェクト名と同じにします。
ターミナルを立ち上げて、次のように入力します。
$ mkdir TiC
$ cd TiC
structファイルの作り方は、プログラムmakestruct_lapwを使って自分で結晶構造を入力する方法と、結晶構造の統一フォーマットであるCIFファイルを変換する方法の2通りがあります。
既知の物質を調べる場合、結晶構造は論文か物質科学系のデータベースサイトを参照すると良いでしょう。
今回は以下のサイトに記載されているTiCの情報を使います。
* [[http://www.crystallography.net/cod/5910091.html|Crystallography Open Database: Information card for entry 5910091]]
==== makestruct_lapwを使う方法 ====
基本的には質問に答えていくだけで完了します。
SPHERE RadIIやRMTの設定はとりあえずデフォルトのままにします。
$ makestruct_lapw
************************************************
* *
********** Terminal struct maker **********
********** (C) 2012 by Morteza Jamal **********
* *
************************************************
TITLE :TiC
This Program accepts a SPACE GROUP (symbol or number) or a
LATTICE TYPE (P, F, B, H, R, CXY, CXZ, CYZ ).
But, with LATTICE TYPE YOU HAVE TO put in all the atomic positions by hand.
Would you like to enter Spacegroup or Lattice (S/L)(def=S)? S
SPACE GROUP: (type ENTER or give first LETTER for a list)
give SPACE GROUP as SYMBOL or NUMBER: 225
Info: space group is : 225 F Fm-3m -F4;2;3
Units of lattice parameters (Bohr/Angstrom) (b/A) (def=ANG):A
Lattice PARAMETERS as a b c (3 numbers):4.29 4.29 4.29
ANGLES BETWEEN lattice vectors, as alpha beta gamma (def=90.0 90.0 90.0):90.0 90.0 90.0
NUMBER INEQUEVALENT ATOMS :2
ATOM 1 (ELEMENT): Ti
POSITION OF ATOM Ti as X,Y,Z (def=0 0 0) :0 0 0
ATOM 2 (ELEMENT): C
POSITION OF ATOM C as X,Y,Z (def=0 0 0) :0.5 0.5 0.5
Now, 'datastruct' file is ready.
We Run 'Tmaker' for making WIEN2k struct file.
Fm-3m
8.10692830565198 8.10692830565198 8.10692830565198
90.0000000000000 90.0000000000000 90.0000000000000
ATOM NAME:Ti
Ti/1 0.00000000 0.00000000 0.00000000
ATOM NAME:C
C/1 0.50000000 0.50000000 0.50000000
'init.struct' file is ready
SETTING UP SPHERE RadII:
SPECIFY possible REDUCTION of SPHERE RADII in % (def=0)
0
specify nn-bondlength factor: (usually=2) [and optionally dlimit, dstmax (about
1.d-5, 20)]
DSTMAX: 20.0000000000000
iix,iiy,iiz 5 5 5 40.5346400000000
40.5346400000000 40.5346400000000
ATOM 1 Ti ATOM 2 C
RMT( 1)=2.00000 AND RMT( 2)=2.00000
SUMS TO 4.00000 LT. NN-DIST= 4.05346
ATOM 2 C ATOM 1 Ti
RMT( 2)=2.00000 AND RMT( 1)=2.00000
SUMS TO 4.00000 LT. NN-DIST= 4.05346
NN ENDS
0.0u 0.0s 0:00.02 0.0% 0+0k 2240+32io 9pf+0w
atom Z RMT-max RMT
1 22.0 2.22 2.22
2 6.0 1.81 1.81
file init.struct_setrmt generated
rerun setrmt ?(y,N) (def=N):
N
The file init.struct has been created
for modifications of your input you can also edit file datastruct and run
Tmaker / setrmt init -r X individually
$
init.structというファイルを (プロジェクト名).struct に変更しておきます。
$ mv init.struct TiC.struct
==== CIFファイルからstructファイルを生成する方法 ====
実際にWIEN2kを使った第一原理計算を行うときは、対話形式でstructファイルを作成するよりも、結晶構造データの統一フォーマットであるCIFファイルから生成することが多くなると思います。
上記のCODのサイトからCIFファイルをダウンロードし、ファイル名をTiC.cifに変更して作業ディレクトリに保存します。
用意ができたら、プログラムcif2structを使ってstructファイルに変換します。
$ cif2struct TiC.cif
これで入力ファイルTiC.structができます。
ただし、このままだとマフィンチン半径$R_{\rm MT}$が正しく設定されていないので計算します。
拡張子の.structは不要です。
$ setrmt_lapw TiC
マフィンチン半径が正しく設定された入力ファイルTiC.struct_setrmtというファイルができるので、さきほどの入力ファイルに上書きします。
$ mv -f TiC.struct_setrmt TiC.struct
===== 計算の設定の初期化 =====
SCF計算に入る前に、計算の設定の初期化をおこないます。
ここではk点数とRKmaxを設定しましょう。
k点数やRKmaxを増やせば計算の精度はあがりますが、計算時間もかかるようになります。
最終的にほしい物理量がこれらのユーザー側が設定した量に対して収束しているか確認する必要があります。
なお、RKmaxは波動関数のカットオフを決めます。原子の種類によって必要な大きさが異なるので、選び方は http://susi.theochem.tuwien.ac.at/reg_user/faq/rkmax.html を参照してください。
ここではk点数を1000に、RKmaxを7.5にしてみましょう。
$ init_lapw -b -numk 1000 -rkmax 7.5
* ここでフラグ -b はバッチモードで、これをつけることで指定のもの以外はすべてまとめてデフォルトの値に設定されます。
* 最後に init_lapw finished ok と出たらうまく行っています。このメッセージが出てこなかった場合、どこかで問題が起こっているので修正します。
オプションの意味は次のとおりです
^オプション^初期値^説明^
|-rkmax|7.0|RKmax|
|-numk|1000|k点数|
他の設定はフラグ-hにより参照できます。
$ init_lapw -h
===== SCF計算の実行 =====
これで準備が整いましたので計算を実行します。
run_lapw を使用します。
$ run_lapw -cc 0.0001 -ec 0.00001 -i 100
^オプション^初期値^説明^
|-cc|なし|電荷の収束。単位はe。指定しなければ電荷の収束はチェックされないが、見ておいたほうが無難|
|-ec|0.0001|エネルギーの収束。単位はRy。デフォルトだと少し荒い|
|-i|40|SCF計算を行うときの、繰り返し回数の最大値。40回だと足りないこともあるので少し多めに。|
他のオプションはフラグ-hで確認できます。
$ run_lapw -h
実行後、いろいろなファイルができますが、計算の経過や重要な情報は TiC.scf で確認できます。
$ grep ":ENE" TiC.scf # 全エネルギー
（略）
:ENE : ********** TOTAL ENERGY IN Ry = -1783.96056300
$ grep ":FER" TiC.scf # フェルミ準位。単位は Ry です。
（略）
:FER : F E R M I - ENERGY(TETRAH.M.)= 0.7492912035
$ grep ":GAP" TiC.scf # バンドギャップの大きさ。金属の場合は0になります。
（略）
:GAP (global) : 0.0 Ry = 0.0 eV (metal)
これで電子状態が得られましたので、これをもとに状態密度やバンド分散などを計算することになります。
===== SCF計算の再実行 =====
SCF計算を再度実行したい場合は、ブロイデン法の履歴を削除してから再実行します。
（電子密度の分布は前回の計算のものが引き継がれます）
$ rm *.broyd*
$ run_lapw -cc 0.0001 -ec 0.00001 -i 100
フラグ-NIをつけて再実行することも可能ですが、計算条件を変更する場合は使えません（エラーが出ます）。
$ run_lapw -cc 0.0001 -ec 0.00001 -i 100 -NI
====== k点数の変更 ======
k点数が十分であったかを見るために、k点数を変更したいとします。
再び init_lapw を走らせると全ての設定や計算結果が初期化されますので、k点数だけ変更するために x_lapw kgen を使います。
$ x_lapw kgen
NUMBER OF K-POINTS IN WHOLE CELL: (0 allows to specify 3 divisions of G)
2000
length of reciprocal lattice vectors: 1.342 1.342 1.342 12.599 12.599 12.599
72 k-points generated, ndiv= 12 12 12
KGEN ENDS
0.0u 0.0s 0:05.99 0.6% 0+0k 0+304io 0pf+0w
計算の再実行は上で述べたようブロイデン法の履歴を削除してから行います。
最初は少ないk点数（10くらい）から初めて、徐々に増やしていきましょう。
k点数に対して全エネルギーやフェルミ準位（半導体の場合はエネルギーギャップの大きさ）をプロットしてみて、収束しているかどうか必ず確認しましょう。
====== 計算結果の保存 ======
計算の入出力ファイルを残しておきたい場合は save_lapw を使用します。
例えばk点を増やして計算する際に、k点数を1000で計算したときの入力ファイルや結果を k1000 という名前のディレクトリに保存しておきたい場合は次のようにします。
$ save_lapw -d k1000
また、この保存したファイルを復活させる場合は次のようにします。-fは同じ名前のファイルがあった時に、強制的に上書きするオプションです。
$ restore_lapw -d k1000 -f
====== ジョブの一時停止 ======
ジョブを一旦止めたい時には、作業ディレクトリに.stopという名前の空のファイルを置けばよいようになっています。
$ touch .stop
この空のファイルは、ジョブが止まった後に自動的に削除されます。
====== プログラム名について ======
プログラム名の最後に _lapw をつけていましたが、これを省略しても良いようになっています。
例えば x_lapw は単に x だけで実行可能です。
====== エラーや警告への対処 ======
===== atom 1 has a large sphere , consider setting HDLOs and/or larger LVNS =====
init_lapw を実行したとき、lstartの段階でこのようなメッセージが出たときは、オプション -LVNS X （Xは4から10までの整数）をつけて、エラーが消えるまでXを増やしてください。
$ init_lapw -b -numk 1000 -rkmax 7.5 -lvns 5