AkaiKKRで凝集エネルギー(未完)

凝集エネルギーは(ゼロ点振動のエネルギーなどを無視すれば)、凝集状態の全エネルギーと孤立原子の全エネルギーの差で表すことができます。仮に、孤立原子のエネルギーを、結晶の格子定数を大きくしていった極限として見積もることが可能なら、通常の第一原理計算パッケージでも凝集エネルギーを計算できることになります。

今回はAkaiKKR(machikaneyama)を用いてダイヤモンド構造のシリコンの格子定数を大きくしていって、全エネルギーの変化を調べる事により、シリコンの凝集エネルギーを見積もることを試みました。しかしながら、格子定数を大きくするに従って収束させるのが難しくなり、凝集エネルギーとしてリーズナブルな値を得ることはできませんでした。

CohesiveSi.png

Fig.1: ダイヤモンド型構造シリコンの全エネルギーの格子定数依存性。



凝集エネルギー


全エネルギーって何だよ?では、第一原理計算パッケージの出力する全エネルギーと固体物理の熱力学で言うところの凝集エネルギーの関係について書きました。
話の結論としては、全エネルギーは、エネルギーの原点に物理的な意味がないので、凝集エネルギーとは直接比較することができないが、よく似た概念であるということでした。

そこで、凝集エネルギーの原点である孤立原子の全エネルギーを、第一原理計算パッケージで計算しておけば、凝集状態の全エネルギーとの差から、凝集エネルギーを計算できるのではないでしょうか?

CohesiveOsaka2002.png

固体中の原子を引き離していったら…(PDF)で計算されているシリコンの格子定数と全エネルギー


このような試みはosaka2002のテクニカルノートにも見つけることができます(参考: 固体中の原子を引き離していったら…(PDF))。そこで今回は、AkaiKKR(machikaneyama)を用いてダイヤモンド構造のシリコンの格子定数を大きくしていって、全エネルギーの変化を調べました。

計算方法


これまで通り、入力ファイルのテンプレートからシェルスクリプトを利用して連続的に入力ファイルを作成していきます。交換・相関ポテンシャルにはmjwを使いました。これはMoruzzi Janak and Williamsに原子のエネルギーと凝集エネルギーの値が載っている事を意識してのことでしたが、この中にシリコンの値はありませんでした。Osaka2002のテクニカルノートにはスピン分極に関する言及もあるので、スピン分極ありの計算(mag)としました。相対論の効果は考慮せず nrl を指定しています。

c------------------------------------------------------------
go data/CohesiveSi_ABOHR
c------------------------------------------------------------
c brvtyp a c/a b/a alpha beta gamma
fcc ABOHR , , , , , ,
c------------------------------------------------------------
c edelt ewidth reltyp sdftyp magtyp record
0.001 1.0 nrl mjw mag 2nd
c------------------------------------------------------------
c outtyp bzqlty maxitr pmix
update 4 1000 0.01
c------------------------------------------------------------
c ntyp
2
c------------------------------------------------------------
c type ncmp rmt field mxl anclr conc
Si 1 1 0.0 2
14 100
Vc 1 1 0.0 2
0 100
c------------------------------------------------------------
c natm
4
c------------------------------------------------------------
c atmicx(in the unit of a) atmtyp
0 0 0 Si
1/4 1/4 1/4 Si
1/2 1/2 1/2 Vc
3/4 3/4 3/4 Vc
c------------------------------------------------------------


#!/bin/csh -f

## *** プロジェクト名 ***
set PROJECT="CohesiveSi"
## *** 格子定数 ***
set ABOHR_LIST=( 10.26 11.97 13.68 15.39 17.10 18.81 20.52 22.23 23.94 )

## *** 繰り返し計算 ****
foreach ABOHR ( ${ABOHR_LIST} )
## ポテンシャルが存在しなければコピー
if ( ! -e data/${PROJECT}_${ABOHR} ) then
if ( -e data/${PROJECT} ) then
cp data/${PROJECT} data/${PROJECT}_${ABOHR}
endif
endif
## 入力ファイルの作成
sed 's/'ABOHR'/'${ABOHR}'/g' template/${PROJECT}_template.in > in/${PROJECT}_${ABOHR}.in
## 第一原理計算
specx < in/${PROJECT}_${ABOHR}.in > out/${PROJECT}_${ABOHR}.out
## 次のポテンシャルをコピー
cp data/${PROJECT}_${ABOHR} data/${PROJECT}
end


格子定数を大きくするほどバンドの分散が小さくなるので、ewidthも格子定数を大きくするのに従って小さくするほうが収束しやすくなるのだと思いますが、すべての計算で同じ値としました。

結果


Fig.2-3に格子定数 a0 = 10.26 Bohr のときのバンド構造と a = 15.39 Bohr (a = 1.5 a0, a0/a = 0.67)のときのバンド構造を示します。格子定数が大きくなることによってバンドの分散が小さくなっていることが分かります。

CohesiveSi1026.png
Fig.2: ダイヤモンド構造シリコンのバンド構造(a0 = 10.26 Bohr)

CohesiveSi1539.png
Fig.3: ダイヤモンド構造シリコンのバンド構造(a = 15.39 Bohr)


Fig.1に a0 = 10.26 Bohr としたときの格子定数の変化に対する全エネルギーをプロットしたものを示します。格子定数を大きくするほど全エネルギーが大きくなりましたが a = 20.52 Bohr (a0/a = 0.5) よりも大きい辺りから収束させるのが難しくなっていきました。これはbzqlty=4としましたが、bzqlty=9としても同様です。

計算できなくなる直前の値が、ほぼ孤立原子と同じエネルギーになっていると考えると凝集エネルギーは Ecoh = 0.47 Ry (6.37 eV) 程度となります。実験値は Ecoh = 4.63 eV との事なので(スピン分極を考慮した計算であるにもかかわらず)過大評価をしています。

議論(反省)


実を言うと凝集エネルギーに関する議論は、AkaiKKRの掲示板に投稿があります(How can I calculate cohesive energy? (in Japanese))。その内容によると、孤立原子のエネルギーは計算できるようにする予定が、少なくとも当時は、あったようですが、現在のバージョンにも実装はされていないようです。

前述の通り、格子定数を大きくするとバンドの分散が小さくなるのでewidthを小さくするほうがよさそうな気がします。
その一方で、同じ基準でエネルギーを比較するためにはewidthを固定する必要があるようにも思います。

格子定数 a = 10.26 Bohr では ewidth=1.0 Ry 程度が必要です。通常のgo計算では ewidth の範囲を source/specx.f の mse0 で指定されたエネルギー分解能で計算するので、この値を大きくすることにより、大きな ewidth でも高分解能で起算できるかもしれません。

関連エントリ




参考URL




付録


このエントリで使用したファイルを添付します。ファイル名末尾の".txt"を削除して、"_"を"."に変更すれば使えるはずです。(参考:ねがてぃぶろぐの付録)


参考文献/使用機器




フィードバック



にほんブログ村 その他趣味ブログ 電子工作へ

 ↑ 電子工作ブログランキング参加中です。1クリックお願いします。


コメント・トラックバックも歓迎です。 ↓      


 ↓ この記事が面白かった方は「拍手」をお願いします。


tag: AkaiKKR machikaneyama KKR 

comment

Secret

FC2カウンター
カテゴリ
ユーザータグ

LTspiceAkaiKKRmachikaneyamaScilabKKRPSoC強磁性CPAPICOPアンプecalj常微分方程式モンテカルロ解析状態密度トランジスタodeDOSインターフェース定電流スイッチング回路PDS5022半導体シェルスクリプト分散関係レベルシフト乱数HP6632AR6452A可変抵抗トランジスタ技術温度解析ブレッドボードI2C反強磁性確率論数値積分セミナーバンドギャップバンド構造偏微分方程式非線形方程式ソルバ熱設計絶縁ISO-I2Cカオス三端子レギュレータLM358GW近似マフィンティン半径A/DコンバータフォトカプラシュミットトリガLEDPC817C発振回路数値微分直流動作点解析サーボカレントミラーTL431アナログスイッチUSB74HC4053bzqltyVESTA補間電子負荷アセンブライジング模型BSch量子力学単振り子2ちゃんねるチョッパアンプLDA開発環境基本並進ベクトルFFT標準ロジックブラべ格子パラメトリック解析抵抗SMPMaxima失敗談ラプラス方程式繰り返し位相図スイッチト・キャパシタ熱伝導状態方程式キュリー温度gfortranコバルトTLP621不規則合金Quantum_ESPRESSO六方最密充填構造ランダムウォーク相対論ewidthスピン軌道相互作用FETQSGWVCAcygwinスレーターポーリング曲線GGA仮想結晶近似PWscfシュレディンガー方程式LM555ハーフメタル固有値問題NE555最小値ガイガー管QNAPUPS自動計測ダイヤモンドマントルTLP552格子比熱最適化MCU井戸型ポテンシャル最大値xcrysdenCIF条件分岐詰め回路フェルミ面差し込みグラフスーパーセルfsolveブラウン運動awk過渡解析起電力三角波第一原理計算FXA-7020ZRWriter509Ubuntuテスタ熱力学データロガーTLP521OpenMPubuntu平均場近似MAS830LトランスCK1026PIC16F785PGA2SC1815EAGLEノコギリ波負帰還安定性ナイキスト線図MBEOPA2277P-10フィルタCapSenseAACircuitLMC662文字列固定スピンモーメントFSMTeX結晶磁気異方性全エネルギーc/a合金multiplotgnuplot非線型方程式ソルバL10構造正規分布等高線ジバニャン方程式初期値interp1fcc面心立方構造ウィグナーザイツ胞半金属デバイ模型電荷密度重積分SIC二相共存磁気モーメント不純物問題PWgui擬ポテンシャルゼーベック係数ZnOウルツ鉱構造edeltquantumESPRESSOフォノンリジッドバンド模型スワップ領域BaO岩塩構造ルチル構造ヒストグラム確率論マテリアルデザインフラクタルマンデルブロ集合キーボードRealforceクーロン散乱三次元疎行列縮退化学反応関数フィッティング最小二乗法Excel直流解析PCTS-110TS-112日本語パラメータ・モデル等価回路モデルcif2cell入出力陰解法熱拡散方程式HiLAPW両対数グラフCrank-Nicolson法連立一次方程式specx.fifort境界条件片対数グラフグラフの分割円周率ヒストグラム不規則局所モーメントGimpシンボル軸ラベル凡例線種トラックボール

最新コメント
リンク

にほんブログ村 その他趣味ブログ 電子工作へ