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カウンター
カテゴリ
ユーザータグ

LTspiceAkaiKKRmachikaneyamaScilabKKRPSoCOPアンプCPA強磁性PICモンテカルロ解析常微分方程式odeトランジスタecalj状態密度DOSインターフェース定電流スイッチング回路PDS5022半導体シェルスクリプト乱数レベルシフトHP6632A温度解析分散関係I2Cトランジスタ技術R6452A可変抵抗ブレッドボードセミナーバンドギャップ数値積分確率論反強磁性偏微分方程式バンド構造絶縁熱設計非線形方程式ソルバフォトカプラシュミットトリガLEDLM358カオスISO-I2C三端子レギュレータGW近似A/Dコンバータカレントミラーアナログスイッチ数値微分マフィンティン半径TL431発振回路サーボPC817CUSB直流動作点解析74HC4053補間FFTBSch開発環境パラメトリック解析2ちゃんねるチョッパアンプ量子力学bzqlty電子負荷イジング模型LDA標準ロジックアセンブラ基本並進ベクトルブラべ格子単振り子熱伝導位相図TLP621キュリー温度繰り返し状態方程式MaximaVESTAスイッチト・キャパシタ相対論FETランダムウォークスピン軌道相互作用SMP六方最密充填構造抵抗不規則合金ewidthスレーターポーリング曲線GGAラプラス方程式cygwingfortranQSGW失敗談コバルト条件分岐TLP521テスタLM555Writer509TLP552格子比熱マントルデータロガー自動計測詰め回路ガイガー管ダイヤモンドQNAPMCUFXA-7020ZR過渡解析三角波UPSNE555固有値問題熱力学ブラウン運動フェルミ面awk起電力第一原理計算OpenMPfsolveubuntu最大値xcrysden最小値最適化仮想結晶近似VCA差し込みグラフスーパーセル井戸型ポテンシャル平均場近似シュレディンガー方程式FSMフラクタルOPA2277固定スピンモーメント2SC1815全エネルギー合金multiplotgnuplotc/aTeX結晶磁気異方性interp1ウィグナーザイツ胞初期値マンデルブロ集合疎行列面心立方構造fcc不純物問題非線型方程式ソルバフィルタL10構造PGA半金属二相共存SICZnOウルツ鉱構造BaO重積分クーロン散乱磁気モーメント電荷密度三次元CIF岩塩構造CapSenseノコギリ波デバイ模型ハーフメタル正規分布フォノンquantumESPRESSOルチル構造スワップ領域リジッドバンド模型edelt縮退キーボード軸ラベルグラフの分割凡例トラックボールPC不規則局所モーメント片対数グラフトランス両対数グラフCK1026MAS830L直流解析Excel円周率パラメータ・モデルヒストグラム日本語最小二乗法等価回路モデルGimp線種シンボルTS-110TS-112PIC16F785LMC662化学反応文字列specx.f入出力ifortマテリアルデザインヒストグラム確率論Realforce等高線ジバニャン方程式P-10Ubuntuナイキスト線図Crank-Nicolson法陰解法熱拡散方程式HiLAPWAACircuit連立一次方程式負帰還安定性境界条件EAGLEMBE関数フィッティング

最新コメント
リンク

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