AkaiKKRのspc計算の入出力フォーマット その2

AkaiKKRのspc計算の入出力フォーマットでは、AkaiKKR(machikaneyama)の26 Augsut 2015のバージョンを利用して、それよりも以前の 22 May 2015 のバージョンと同じspcの入出力フォーマットに戻す方法を書きました。入力フォーマットの変更は簡単だったのですが、出力フォーマットの指定方法の変更は、少し面倒くさい感じでした。最新版の 22 July 2016 バージョンのAkaiKKRでは、この出力フォーマットの指定方法も簡単になりました。

具体的には source/spmain.firdfmt の値を変更することによって入力のフォーマットを、また、 iwrtfmt の値を変更することによって出力のフォーマットを変更できます。22 May 2015 のバージョンと同じspcの入出力フォーマットに戻すには、両方とも 1 を指定します。デフォルトの 3 を指定しておくほうがバンド構造のテスト計算には便利かもしれません(参考:AkaiKKRでバンド構造(分散関係))。しかしながら、ねがてぃぶろぐでは、今後とも基本的に22 May 2015 のバージョンと同じ入出力フォーマットで計算をすることにします。

関連エントリ




参考URL




フィードバック



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

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


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


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

tag: AkaiKKR machikaneyama 分散関係 

AkaiKKRとecaljでCuGaTe2 その1

第一原理計算パッケージには、それぞれ特徴があり、計算したい物質によって適切に使い分ける必要に迫られることがあります。AkaiKKR(machikaneyama)は不規則系に適しており、ecaljは半導体のバンドギャップを求めるのに適しています。

例えば、不規則を含む半導体の計算をAkaiKKRで行いたいと考えたとき、不規則を含まない端成分の計算をecaljの結果と比較しておくことは有用です。今回はCuGaTe2を対象として、AkaiKKRで状態密度の計算をおこないました。

CuGaTe2DOS.png
Fig.1: CuGaTe2の状態密度



AkaiKKRとecaljの長所


AkaiKKR(machikaneyama)は、コヒーレントポテンシャル近似(CPA)を導入することによって、合金などの不規則性を扱うことが可能であるという特徴があります。
またecaljはGW近似を用いて、半導体のバンドギャップの見積もりを局所密度近似(LDA)から改善できる長所があります。

他にもさまざまな第一原理計算パッケージが、それぞれ特有の長所を持っています。このため、しばしば複数のコードでの計算結果を比較するということが起こります。

今回と次回では、AkaiKKRの掲示板に投稿された CuGaTe2 のバンドギャップをこれら二つのコードで計算し、バンドギャップと状態密度の比較を行います。今回はAkaiKKRでの計算です。

計算手法


入力ファイルはCannot reproduce the bandgap of CuGaTe2に投稿されているものとほとんど同じですが、少しだけ変更してあります。一つ目の変更点は、スピン軌道相互作用を(計算が重いので)はずした事。二つ目はewidthを小さくしたことです。

c--------------------CuGaTe2---------------------------------
go data/cugate2
c------------------------------------------------------------
c brvtyp a c/a b/a alpha beta gamma
bct 11.5388 1.992 1 90 90 90
c------------------------------------------------------------
c edelt ewidth reltyp sdftyp magtyp record
0.001 0.7 sra pbe nmag 2nd
c------------------------------------------------------------
c outtyp bzqlty maxitr pmix
update 4 500 0.015
c------------------------------------------------------------
c ntyp
5
c------------------------------------------------------------
c type ncmp rmt field mxl anclr conc
Cu 1 0 0.0 2 29 100
Ga 1 0 0.0 2 31 100
Te 1 0 0.0 2 52 100
Es1 1 0 0.0 0 0 100
Es2 1 0 0.0 0 0 100
c------------------------------------------------------------
c natm
16
c------------------------------------------------------------
c atmicx atmtyp
0.23703x 1/4y 1/8z Te
0.76297x 3/4y 1/8z Te
3/4x 0.23703y 3/8z Te
1/4x 0.76297y 3/8z Te
1/2x 1/2y 0.0z Ga
1/2x 0.0y 1/4z Ga
0.0x 0.0y 0.0z Cu
0.0x 1/2y 1/4z Cu
c
0.75x 1/4y 1/8z Es1
0.25x 3/4y 1/8z Es1
3/4x 0.75y 3/8z Es1
1/4x 0.25y 3/8z Es1
c
0.0x 0.0y 0.25z Es2
1/2x 1/2y 0.25z Es2
0.0x 1/2y 0.0z Es2
1/2x 0.0y 0.0z Es2
c------------------------------------------------------------


結果


Fig.1に状態密度を示します。
AkaiKKRでの状態密度やバンド構造(ブロッホスペクトル関数)のエネルギー分解能は source/specx.f の msex で指定することが可能で、デフォルトでは msex=201 となっています。したがって、状態密度を計算するために ewidth = 0.8 Ry とした場合の分解能は 4 mRy 程度になります。その結果、状態密度の図だけを見ると、バンドギャップが存在するか否かが微妙です。

AkaiKKRでバンドギャップの測り方では、バンドギャップを決める場合、状態密度から値を読むよりも、バンド構造から見るほうが良さそうであると書きました。CuGaTe2は、伝導帯の上端(CBM)と価電子帯の下端(VBM)が共にΓ点に存在する直接遷移型の半導体であるとの事なので、その付近のバンド構造をプロットしたのがFig.2です。

CuGaTe2band.png
Fig.2: Γ点周辺のCuGaTe2のバンド構造


GaAsの場合と異なり、CBMにフェルミ準位(というか計算上のエネルギー基準点)が張り付いてしまっていますが、電子の数を足し上げるときの数値計算上の誤差と思うので、いまは気にしないことにします。

ローレンツ関数へのフィッティングは、あまりきれいにいかなかったので、目視で読むと、バンドギャップの大きさはおよそ 30 mRy 程度でしょうか。換算すると 0.4 eV 程度となるので、Cannot reproduce the bandgap of CuGaTe2に書かれている通り 1 eV 程度存在するはずのバンドギャップから見ると過小評価です。

AkaiKKRに限らず密度汎関数理論(DFT)に局所密度近似(LDA)や一般化勾配近似(GGA)を組み合わせた第一原理計算パッケージは、バンドギャップを過小評価してしまう問題が広く知られています。
ecaljで利用できるGW近似は、この問題に対する回答のひとつです。AkaiKKRとecaljでCuGaTe2 その2では、ecaljを用いてCuGaTe2の状態密度とバンドギャップを計算します。

関連エントリ




参考URL




付録


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


参考文献/使用機器




フィードバック



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

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


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


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

tag: AkaiKKR machikaneyama KKR CPA ecalj 半導体 バンドギャップ バンド構造 分散関係 GW近似 

AkaiKKRのspc計算の入出力フォーマット

AkaiKKR(machikaneyama)は26 August 2015からブロッホスペクトル関数の入出力形式が変更されました。
簡単にバンド図が描けるようになり便利なのですが、一方で、これまでのほうが便利だったケースもありえます。
今回は、これを22 May 2015のバージョンと同じフォーマットに戻す方法について書きます。

入力フォーマットは source/spmain.f でサブルーチン readk を呼び出す際の引数 iread を1に指定します。
出力フォーマットは source/spmain.f でサブルーチン wrtspc を呼び出している行をコメントアウトしてしまい、代わりに22 May 2015で出力をしている部分のソースコードをコピペします。


AkaiKKRのブロッホスペクトル関数


AkaiKKR(machikaneyama)は通常のバンド図を出力する代わりにブロッホスペクトル関数を出力します。cpa2002v009cの22 May 2015までは、計算するk点をすべて手動で入力しなければなりませんでした。そのためAkaiKKRでSrTiO3ペロフスカイトなどではgnuplotでバンド構造をプロットするための補助的なScilabスクリプトを使っています。

これに対して26 Augst 2015では、他の一般的なDFTパッケージで見られるように、特徴的なk点の座標とその分割数を指定することで、バンド図を描くことができるようになりました(参考: AkaiKKRでバンド構造(分散関係))。

便利になった反面、特徴的なk点をつなぐパスに沿った計算以外を行う場合、例えばAkaiKKRで銅と銅亜鉛合金のフェルミ面のような計算はやりにくくなりました。また、それ以外にもこれまで書いたスクリプトを再利用するのにも不便です。
そんな折、26 August 2015のソースコードを編集して、以前の入出力形式に戻す方法を教えていただいたので、今回は、この方法についてまとめます。

入力ファイルの形式


spc計算のためのk点指定方法は source/spmain.f の382行目(ぐらい)にある call readk(vkp,coa,boa,nk1,nk3,nk3x,nkmx,kcrt,25,kblst,korder,3) の最後の引数の3を変更することで変えることができます。引数は1, 2, 3の三種類から選ぶことができ、その内容はsource/readk.fのコメントに書かれています。
c-----------------------------------------------------------------------
c This program read-in the k-points according either format (a),
c (b) or (c).
c
c (a) simply input all the k-point to be calculated.
c
c (b) input k-point corresponding to the Gamma, K, M, etc., with
c the number of devision of the k-mesh between that point and
c the previous point. For example, input data such as
c 0.0 0.0 0.0 0
c 1.0 0.0 0.0 50
c
c means that the line connecting (0, 0, 0) and (1, 0, 0)
c is divided into 50 pieces, and on each grid points
c including both endpoints the Bloch spectrum function
c will be calculated.
c
c (c) input k-point corresponding to the Gamma, K, M, etc., with
c the total number of k-point to be used as the first data.
c For example, input data such as
c 101
c 0.0 0.0 0.0
c 1.0 0.0 0.0
c 1.0 0.5 0.0
c
c means that the line connecting (0, 0, 0), (1, 0, 0), and
c (1, 0.5, 0) is divided into more or less equidistant 100
c pieces, and on each grid points including both endpoints
c the Bloch spectrum function will be calculated.
c
c In all formats, data will be read-in untill any data that
c do not follow the above format be met. Which format should be
c adopted is controled by a parameter iread. If iread=1,
c format (a), iread=2 for (b), and iread=3 for (c)
c All the data shoulbe be given in the unit of 2pi*(1/a, 1/b, 1/c),
c i.e. actual k vector is 2pi*(kx/a, ky/b, kz/c), where a, b, c
c are the lattice constant along x, y, and z direction.
c
c Coded by H. Akai, 4 June 2015, Tokyo
c-----------------------------------------------------------------------

今回は、すべてのk点を手動で指定する(a)に戻すので、最後の引数である iread の3を1に変更します。
c     call readk(vkp,coa,boa,nk1,nk3,nk3x,nkmx,kcrt,25,kblst,korder,3)
call readk(vkp,coa,boa,nk1,nk3,nk3x,nkmx,kcrt,25,kblst,korder,1)


出力ファイルの形式


出力形式の指定は source/wrtspc.f の中で行われますが source/readk.f ほどには、整備されていないようです。このsource/readk.fが追加されたのは26 August 2015からで22 May 2015の段階では、まだsource/wrtspc.fは作成されておらずファイルの出力も source/spmain.f の中で行われています。26 August 2015では wrtspcのサブルーチンは source/spmain.f の735行目(ぐらい)の call wrtspc(spctrl,is,e,mse,kcrt,kblst,nk3,ef) で呼び出されています。そこでこの呼び出しの部分をコメントアウトしてしまい、代わりに22 May 2015の source/spmain.f の else if(ids .eq. 4) then (748行目ぐらい)と endif (784行目ぐらい)の間をコピー&ペーストします。

こうすることによって26 August 2015を使いながら、ブロッホスペクトル関数の入出力を22 May 2015以前の形式にすることができます。

関連エントリ




参考URL




参考文献/使用機器




フィードバック



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

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


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


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

tag: AkaiKKR machikaneyama 分散関係 

ecaljでシリコンのバンド構造(GW近似)

LDA計算では、半導体のバンドギャップを過小評価するという欠点が広く認知されており、GW近似はこの問題点を改善することが期待されます。今回はecaljでシリコンのバンド構造を計算し、GW近似によりバンドギャップの大きさが改善されていることを確認しました。


001_20150917001717f58.png
Fig.1: シリコンのバンド構造。緑はGW近似による計算、赤はLDAによる計算。GW近似によりバンドギャップの大きさが改善されていることが分かる。



GW近似計算に必要なファイル


ecaljでシリコンのバンド構造(LDA計算)ではecaljを用いてダイヤモンド構造のシリコンのバンド図を局所密度近似(LDA)の範囲で描画しました。ecaljは、更にGW近似を用いた計算も可能です。

この際にecaljは、LDA計算に利用した結晶構造ファイル ctrls.si とバンド図のためのk点指定ファイル syml.si を無編集でそのまま使うことができます。

STRUC   ALAT=10.26
PLAT=0.0 1/2 1/2
1/2 0.0 1/2
1/2 1/2 0.0
SITE ATOM=Si POS=0.0 0.0 0.0
ATOM=Si POS=1/4 1/4 1/4

# num  from            to                name
41 0.5 0.5 0.5 0.0 0.0 0.0 L Gamma
41 0.0 0.0 0.0 1.0 0.0 0.0 Gamma X
21 1.0 0.0 0.0 1.0 0.5 0.0 X W
41 1.0 0.5 0.0 0.0 0.0 0.0 W Gamma
0


これらのファイルを使ってecaljの実行手順(GW近似)に従ってGW近似計算を行いました。

結果


以下のFig.2に緑のラインで描かれたものが、QSGW計算によって得られたシリコンのバンド構造です。対比のためにecaljでシリコンのバンド構造(LDA計算)で計算したLDAの結果をFig.3に示します。


002_20150917001716541.png
Fig.2: GW近似を用いたダイヤモンド構造のシリコンのバンド図

001_20150916211725117.png
Fig.3: LDAを用いたダイヤモンド構造のシリコンのバンド図


更にこれらを同時にプロットしたものが冒頭のFig.1です。

密度汎関数理論(DFT)と局所密度近似(LDA)を組み合わせた第一原理計算パッケージには、半導体のバンドギャップを過小評価するなどの問題点がある事が広く知られています。
この問題を克服するための試みもまた、広く行われており「LDAを超える試み」のキャッチフレーズで色々な方法論が提案されています。
GW近似は、これらの方法の中でもっとも有名なもののひとつで、実際に色々な半導体のバンドギャップの計算結果が、実験により得られている値に対して、LDAによる計算よりもはるかに近くなることが知られています。

実際、今回行ったecaljの計算ではGW近似によるバンドギャップはLDA計算よりも大きく、実験値である1.11eVに近い値となっていることが確認できます。

関連エントリ




参考URL




付録


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


参考文献/使用機器




フィードバック



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

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


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


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

tag: ecalj GW近似 半導体 バンドギャップ バンド構造 分散関係  

ecaljでシリコンのバンド構造(LDA計算)

これまでecaljのインストール(Ubuntu + gfortran)ではecaljをインストールし、ecaljの実行手順(LDA計算)で計算の具体的な手順を確認しました。
今回は、最初の計算例としてダイヤモンド構造を持ったシリコンのバンド構造の計算を局所密度近似(LDA)の範囲で行いました。


001_20150916211725117.png
Fig.1: LDA計算によるシリコンのバンド構造



ecaljの計算に必要なファイル


ecaljは最低限、結晶構造を記述する ctrls.si だけユーザーが作成すれば第一原理計算が可能です。バンド構造を描画するところまで含めても、ユーザーが追加で作成しなければならないファイルは計算するk点のパスを指定する syml.si だけです。
それ以外は ctrlgenM1.py などのスクリプトがテンプレートを半自動的に作成してくれるので、ユーザー側ではそれをそのままコピーするか、多少編集する程度です。

今回は、早速、ダイヤモンド構造のシリコンのバンド構造を計算してみます。

ダイヤモンド構造の結晶構造ファイル


ecaljの結晶構造ファイル ctrls.si は以下の3つの情報だけで作成することができます。
  • 格子定数 a (Bohr)
  • 基本並進ベクトル
  • 基底の原子位置


AkaiKKRでダイヤモンド型構造半導体で書いた通り、ダイヤモンド構造は面心立方格子のそれぞれの格子点に(0 0 0)と(1/4 1/4 1/4)の2原子の基底をもつ結晶構造です。

AkaiKKR(machikaneyama)では格子を指定する際に、ブラべ格子を気ワード指定する方法と基本並進ベクトルを使う方法の2通りがありました。(参考: AkaiKKRのブラベ格子, AkaiKKRの基本並進ベクトル その1, その2)
ecaljでは、基本並進ベクトルを使います。

面心立方格子の基本並進ベクトルは以下のように表すことができます。

\begin{equation}
\begin{pmatrix}
0 & 1/2 & 1/2 \\
1/2 & 0 & 1/2 \\
1/2 & 1/2 & 0
\end{pmatrix}
\end{equation}

ecaljでもAkaiKKRと同様に、格子定数の単位はBohr(原子単位系)で与えます。今回はAkaiKKRでダイヤモンド型構造半導体のときと同様に a = 10.26 Bohr としました。
ecaljのパッケージに付属しているシリコンの結晶構造ファイルは、もう少しいろいろと書いてありますが、今回のエントリでは最も簡素な入力ファイルの例として、以下のものを用意しました。

STRUC   ALAT=10.26
PLAT=0.0 1/2 1/2
1/2 0.0 1/2
1/2 1/2 0.0
SITE ATOM=Si POS=0.0 0.0 0.0
ATOM=Si POS=1/4 1/4 1/4


以降はecaljの実行手順(LDA計算)の2.以降を実行することでLDA計算を行うことができます。

バンド構造の描画


バンド構造の計算も状態密度の計算と同様にポスト処理ということになります。
バンド構造の計算には、あらかじめk点のパスを指定する syml.si を用意しておく必要があります。
内容は、左から順に分割数、k点の開始座標、終了座標、」これらの座標に付ける名前です。
今回は、以下のようなファイルを作成しました。

# num  from            to                name
41 0.5 0.5 0.5 0.0 0.0 0.0 L Gamma
41 0.0 0.0 0.0 1.0 0.0 0.0 Gamma X
21 1.0 0.0 0.0 1.0 0.5 0.0 X W
41 1.0 0.5 0.0 0.0 0.0 0.0 W Gamma
0


なお syml.si の方では、分数を使わない方が良さげです。

gnuplot用ファイルの編集


バンド構造の計算が終わるとgnuplotのウインドウが立ち上がってバンド図が描画されます。この描画スクリプトはbandplot.isp1.gltというなまえで保存されているので、後から以下のようなコマンドで簡単に再描画できます。

gnuplot -persist bandplot.isp1.glt


また、ファイルの中身はただのgnuplotのスクリプトなので、自分で編集することも簡単です。
冒頭のFig.1も自分で編集して、多少見栄えをよくしてあります。

関連エントリ




参考URL




付録


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


参考文献/使用機器




フィードバック



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

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


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


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

tag: ecalj 半導体 バンドギャップ バンド構造 分散関係 

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

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

最新コメント
リンク

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