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近似 

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 半導体 バンドギャップ バンド構造 分散関係 

AkaiKKRでバンド構造(分散関係)

最新版(August 26, 2015)のAkaiKKR(machikaneyama)でバンド構造が簡単に描ける機能が追加されました。そのまま外部向けの文章(プレゼンテーションや論文など)に使うにはそっけないですが、計算の確認には便利です。

002_20150912192936f79.png

003_2015091219293582c.png
Fig.1-2: 面心立方構造(fcc)のNi90Fe10合金のバンド構造(upスピンとdownスピン)


今回は、この機能の使い方について書きます。


バンド構造計算の新機能


AkaiKKR(machikaneyama)では、バンド構造をプロットする代わりに、ブロッホスペクトル関数をプロットします。これまでもこの方法でバンド構造を描くことが可能でしたが、ブロッホスペクトル関数を計算するk点のパスを「全て」手動で設定しなければならなかったり、描画のためのgnuplotのスクリプトを自分で書かなければならなかったりと、不便がありました。(参考:AkaiKKRで不規則NiMn合金の分散関係AkaiKKRでSrTiO3ペロフスカイト)

August 26, 2015版のcpa2002v009cでは、これらの手間が随分と軽減される機能が新たに実装されました。この機能に関しては付属のReleaseNoteV009cに以下のように書かれています。(タイポと考えられる点は勝手に修正してあります。)

Another new feature provided by this version is the "spc" program that plots a dispersion relation (or broadened dispersion/Bloch spectrum function for alloy cases) if gnuplot is installed. If "spc" is specified instead of "go", then the program create data/file_up.spc and data/file_dn.spc, where data/file is the file specified as a data file in the input. Using the command

spc data/file_up.spc

will plot the dispersion on a gnuplot window. In the input data some additonal data specifying the symmetry points are needed in addition to the number of k-points that will be plotted.
For example, input data like

spc data/nife
fcc 6.65 , , , , , ,
0.001 1.0 nrl mjw mag 2nd
update 4 50 0.03
1
NiFe 2 1 0.0 2 26 40
28 60
1
0 0 0 NiFe
c------- below needed only for "spc" cases
c---number of k-point used for the plot
300
c--- W-point
1 0.5 0
c--- L-point
0.5 0.5 0.5
c--- Gamma-point
0 0 0
c--- X-point
1 0 0
c--- W-point
1 0.5 0
c--- K-point
0.75 0.75 0

will create data/nife_up.spc and data/nife_dn.spc.
Then

spc data/nife_up.spc
spc data/nife_dn.spc

plot the broadend dispersion curves for spin up and down states of Ni60Fe40 disorderd alloy along the W-L-Gamma-X-W-K symmetry lines using 300 k-points. To create the "spc" command, just type

make spc

This is a similar procedure to obtain "gdp" command that plot the DOS curve on a gnuplot window using "make gdp".


かなりいい加減な意訳をすると、以下のような感じでしょうか。

このバージョンのもう一つの新機能は"spc"というプログラムです。これはgnuplotがインストールされていれば、分散関係(合金の場合はにじんだ分散/ブロッホスペクトル関数)をプロットします。"go"の代わりに"spc"を指定することにより、specxプログラムはdata/file_up.spcとdata/file_dn.spcを作成します。ここでdata/fileは入力ファイルで指定されたデータファイルです。以下のコマンドを利用してgnuplotのウインドウ上に分散関係をプロットすることができます。

spc data/file_up.spc

入力ファイルには、通常の計算に加えていくつかのデータが必要になります。プロットする対称性の高いk点とその分割数です。
例えば、入力ファイルは以下のようになります。

spc data/nife
fcc 6.65 , , , , , ,
0.001 1.0 nrl mjw mag 2nd
update 4 50 0.03
1
NiFe 2 1 0.0 2 26 40
28 60
1
0 0 0 NiFe
c------- below needed only for "spc" cases
c---number of k-point used for the plot
300
c--- W-point
1 0.5 0
c--- L-point
0.5 0.5 0.5
c--- Gamma-point
0 0 0
c--- X-point
1 0 0
c--- W-point
1 0.5 0
c--- K-point
0.75 0.75 0

これはdata/nife_up.spcとdata/nife_dn.spcを作成します。

次に以下のコマンドを実行します。
spc data/nife_up.spc
spc data/nife_dn.spc

するとNi60Fe40合金のアップスピンとダウンスピンのにじんだ分散関係をW-L-Γ-X-W-Kの経路を300個のk点でプロットします。

表示のための"spc"実行ファイルを作成するには、以下のように端末に入力します。

make spc

これは状態密度曲線をgnuplotで表示する実行ファイル"gpd"を作成するために実行するmake gpdの手続きと同様です。


もはや上記の通りなのですがNi90Fe10合金のバンド構造を計算しています。

表示用プログラムspcのコンパイル


第一原理計算を実行する前に、結果をプロットするための実行ファイルをコンパイルします。コンパイラはifortでもgfortranでもコンパイル可能です。specxを作成したときと同様にAkaiKKRのインストールディレクトリに移動してmakeを利用します。

make spc

これで実行ファイルspcが作成されるはずです。

ポテンシャルファイルの作成


状態密度の計算と同様に、あらかじめgo計算を収束せてポテンシャルファイルを作成しておく必要があります。

./specx < in/nife


この際の入力ファイルは、AkaiKKRのインストールディレクトリにあるin/nifeを利用しました。

c--------------------NiFe------------------------------------
go data/nife
c------------------------------------------------------------
c brvtyp a c/a b/a alpha beta gamma
fcc 6.65 , , , , , ,
c------------------------------------------------------------
c edelt ewidth reltyp sdftyp magtyp record
0.001 1.0 nrl mjw mag 2nd
c------------------------------------------------------------
c outtyp bzqlty maxitr pmix
update 4 50 0.03
c------------------------------------------------------------
c ntyp
1
c------------------------------------------------------------
c type ncmp rmt field mxl anclr conc
NiFe 2 1 0.0 2
28 90
26 10
c------------------------------------------------------------
c natm
1
c------------------------------------------------------------
c atmicx atmtyp
0 0 0 NiFe
c------------------------------------------------------------


maxitr=50だとerr=-6まで達するのに2回ぐらい実行する必要があると思います。(参考:AkaiKKRでテスト計算)

バンド構造の計算


次にバンド構造を計算するための入力ファイルを作成します。
以下の内容のファイルをin/nife-spcとして保存しました。

c--------------------NiFe------------------------------------
spc data/nife
c------------------------------------------------------------
c brvtyp a c/a b/a alpha beta gamma
fcc 6.65 , , , , , ,
c------------------------------------------------------------
c edelt ewidth reltyp sdftyp magtyp record
0.001 1.0 nrl mjw mag 2nd
c------------------------------------------------------------
c outtyp bzqlty maxitr pmix
update 4 50 0.03
c------------------------------------------------------------
c ntyp
1
c------------------------------------------------------------
c type ncmp rmt field mxl anclr conc
NiFe 2 1 0.0 2
28 90
26 10
c------------------------------------------------------------
c natm
1
c------------------------------------------------------------
c atmicx atmtyp
0 0 0 NiFe
c------------------------------------------------------------
c------- below needed only for "spc" cases
c---number of k-point used for the plot
300
c--- W-point
1 0.5 0
c--- L-point
0.5 0.5 0.5
c--- Gamma-point
0 0 0
c--- X-point
1 0 0
c--- W-point
1 0.5 0
c--- K-point
0.75 0.75 0


go計算用の入力ファイルからの変更点はgoをspcにしたことと、ファイルの末尾にk点の分割数と計算するパスの対称性の高い点の座標を追加したことです。

これをspecxへ渡して.spcファイルを作成します。

./specx < in/nife-spc

するとdata/nife_up.spcとdata/nife_dn.spcが作成されます。
これらをそれぞれ実行ファイルspcへ渡せば、バンド構造が表示されます。

./spc data/nife_up.spc


補足


AkaiKKR(machikaneyama)のバージョンは、現在cpa2002v009cですが、この同じ名前のバージョンの中で頻繁に更新がされています。そんなわけで、バージョン名の代わりに公開された日付で、バージョンの違いを表しています。今回の記事ではAugust 26, 2015のバージョンを使って計算しています。

ブロッホスペクトル関数の計算における、入力ファイルと出力結果のフォーマットは、バージョンごとによく変更されるようです。今回の更新の前もAkaiKKRで不規則NiMn合金の分散関係のエントリを書いていた時期とAkaiKKRでSrTiO3ペロフスカイトのエントリを書いていた時期では出力のフォーマットが異なります。

今後、ねがてぃぶろぐでどのバージョンを使うかは、悩ましい問題ではありますが、おそらく今回(August 26, 2015)のバージョンよりも前の古いフォーマットのものを使い続けることになると思います。その理由は、AkaiKKRで銅と銅亜鉛合金のフェルミ面のような計算をする際に、任意のk点でブロッホスペクトル関数を指定できると便利だからです。

関連エントリ




参考URL




参考文献/使用機器




フィードバック



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

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


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


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

tag: AkaiKKR machikaneyama KKR CPA 分散関係 バンド構造 不規則合金 

AkaiKKRでニッケル・鉄・コバルト

AkaiKKRで不規則NiFe合金の強磁性AkaiKKRで不規則NiMn合金の分散関係では、AkaiKKR(Machikaneyama)を用いて面心立法構造(face-centered cubic; fcc)を持つ金属の電子構造を計算し、AkaiKKRでダイヤモンド型構造半導体ではfcc構造を基本とする半導体の計算を行いました。

今回は、fcc構造に加えて基本的な結晶構造である体心立方構造(body centered cubic; bcc)と六方最密充填構造(hexagonal closed package; hcp)の構造をもつ金属として強磁性のニッケル・鉄・コバルトの計算を行いました。

NiDOS_20130210233144.png FeDOS.png CoDOS.png



specx.fの設定


今回の計算では、specx.fを下記の設定にしてmakeしました。
     & (natmmx=4, ncmpmx=4, msizmx=198, mxlmx=3, nk1x=2200, nk3x=2688,
& msex=201, ngmx=15, nrpmx=650, ngpmx=650, npmx=350, msr=400)


fcc構造ニッケルの入力ファイル


まずは、fcc構造を持つニッケルでこれまでの復習をします。入力ファイルにおいて、それぞれの結晶構造の個性を表現しているのは、fcc, bcc, hcpといった構造の指定に加えて、原子位置の指定ブロッホスペクトル関数を計算するk点の座標です。

原子位置の指定に関して言えば、例えば、同じfcc構造を指定したとしても、単位格子の中で実空間の(0 0 0)にだけ原子をおくとfcc構造となり、実空間の(0 0 0)と(0.25 0.25 0.25)におくとfcc構造を1/4ずつずらして作られたダイヤモンド構造になります。(AkaiKKRでダイヤモンド型構造半導体参照)

以下に示すのが、VESTAを用いて作図したfcc構造ニッケルの原子配置です。cifファイルは結晶構造ギャラリーのものを利用させていただきました。AkaiKKRの入力ファイルでは、実空間の(0 0 0)にだけ原子を置きます。


fccNi.png
Fig.1: fcc構造ニッケルの原子配置


次に、ブロッホスペクトル関数を計算するk点の座標です。この波数空間の座標は、第一Brillouin zone(B.Z.)の特徴的なk点をつなぐように選ばれます。上手な選び方は一通りではないのですが、今回はtable.1の選び方にしました。table.1の中の波数空間の座標の単位は2π/aです。

FCC.png
Fig.2: fcc構造の第一B.Z.(Wikipediaより転載)


kxkykz
W1/201
L1/21/21/2
Γ000
X001
W1/201
K3/403/4
table.1: fcc構造の計算するk点のパス、座標の単位は2π/a


これらを踏まえて作成した入力ファイルがNi_in.txtです。ブロッホスペクトル関数を計算するk点は、table.1のパスに沿った第一B.Z.の特徴的なk点の間をそれぞれ100等分したものを使いました。

fcc構造ニッケルの全状態密度


以下にニッケルの全状態密度(total DOS)の計算結果を示します。
ニッケルの全状態密度は、アウトプットファイルのtotal DOSのところに書き込まれます。今回扱う金属はいずれも強磁性状態の計算なので、上向きスピンと下向きスピンの両方の2つが見つかるはずです。

この部分を取り出して、下向きスピンの状態密度を-1倍したものがNiDOS_dat.txtです。更にこれをgnuplotでプロットしたものがFig.3です。plotにつかったgnuplot用のスクリプトはNiDOS_plt.txtです。


NiDOS_20130210233144.png
Fig.3: fcc構造ニッケルの状態密度


fcc構造ニッケルのバンド構造


次にニッケルのブロッホスペクトル関数の計算結果から、バンド構造(E-k曲線)を描くことにします。
AkaiKKRで不規則NiMn合金の分散関係で書いたとおりgnuplotでカラーマップを使った2次元プロット(pm3d map)の方法を用いてE-k曲線を描きます。そのためには、生成されたspcファイルの行頭に波数空間のパスの道のりの数値を補う必要があります。

W0
L0.70710678
Γ1.57313218
X2.57313218
W3.07313218
K3.42668558
table.2: W点からの逆格子空間の道のり


AkaiKKRで不規則NiMn合金の分散関係のときは、Excelでちまちまと計算したと書きましたが、今回はxyzzyのlispで簡単なマクロを書きました。xyzzyでAkaiKKR-modeに従って導入し、spcファイルを開いた状態で

M-x AkaiKKR-spc2ek


とやれば処理してくれるはずです。かなりいい加減に書いたマクロなので動作は無保証です。また、ファイルは破壊的に編集されるのでバックアップは自分でとってください。あと、処理速度は遅いです。

このようにして出来たファイルをgnuplotでグラフにしたのがFig.4です。グラフのプロットに利用したgnuplotのスクリプトはNiSPC_plt.txtです。


003_20130120185338.png
004_20130120185337.pngFig.4: fcc構造ニッケルのバンド構造


bcc構造鉄の入力ファイル


鉄は常温常圧でbcc構造をとります。Fig.5に示したのがbcc構造の鉄の原子配置です。VESTAを用いて結晶構造ギャラリーのcifファイルを描画しています。fcc構造ニッケルと同様に実空間の(0 0 0)にだけ原子を置きます。


bccFe.png
Fig.5: bcc構造鉄の原子配置


Fig.6に示したのがbcc構造の第一B.Z.です。ブロッホスペクトル関数を計算するk点のパスは、table.3のようにとりました。table.3の中の波数空間の座標の単位は2π/aです。

BCC.png
Fig.6: bcc構造の第一B.Z.(Wikipediaより転載)


kxkykz
Γ000
H001
N1/201/2
P1/21/21/2
Γ000
N1/201/2
table.3: bcc構造の第一B.Z.、座標の単位は2π/a


これらを踏まえて作成した入力ファイルがFe_in.txtです。ブロッホスペクトル関数を計算するk点は、table.3のパスに沿った第一B.Z.の特徴的なk点の間をそれぞれ100等分したものを使いました。

bcc構造鉄の全状態密度


以下に鉄の全状態密度(total DOS)の計算結果を示します。
ニッケルのときと同様な処理をしたデータファイルがFeDOS_dat.txtで、gnuplotのスクリプトファイルがFeDOS_plt.txtです。


FeDOS.png
Fig.7: bcc構造鉄の状態密度


bcc構造鉄のバンド構造


続いて、鉄のバンド構造に関してです。これも基本的にニッケルと同じで、第一B.Z.の形が異なるので逆格子空間での道のりも当然ながら違います(table.4)。

Γ0
H1
N1.70710678
P2.20710678
Γ3.07313218
N3.78023897
table.4: Γ点からの逆格子空間の道のり


gnuplotでE-k曲線を描くための前処理ですが、道のりを計算するだけなのでAkaiKKR-spc2ekは、bcc構造に対してもそのまま使えます。グラフをプロットするためのgnuplotスクリプトはFeSPC_plt.txtです。


FeUp.png
FeDown.pngFig.8: bcc構造鉄のバンド構造


hcp構造コバルトの入力ファイル


最後にhcp構造を持つコバルトについての計算を行います。
まず、Fig.9にコバルトの原子配置をVESTAを用いて描画したものを示します。結晶構造ギャラリーにはコバルトのcifファイルが存在しないので、同じhcp構造をもつマグネシウムのcifファイルの原子の名前と格子定数の部分を編集したものを利用しました。


hcpCo.png
Fig.9: hcp構造コバルトの原子配置


hcp構造の入力ファイルの考え方は、fccやbccのものと比較して多少複雑です。その理由は、下記の2点が主なものです。
  • 基本単位格子の中に2つの原子を含む
  • 格子定数が a ≠ c

fccやbccでは、c/a = 1は明らかなので省略してきました。hcpであっても、省略は可能で、その場合は c/a = 1.633となります。
この値は、理想的なhcp構造のc/a = 2√2/√3ですが、現実のhcp構造のc/aは必ずしも1.633とはなりません。
差し当たり、今回は省略してc/a=1.633を使うことにします。

基本単位格子の中に2この原子を含むということは、nmatm=2ということです。これはAkaiKKRでダイヤモンド型構造半導体で扱ったダイヤモンド構造の結晶も同様でした。このときは、実空間での原子の座標atmicxを直交座標で指定しました。hcp構造では、基本ベクトルで指定するほうが簡単だと思います。計算機マテリアルデザイン入門 (大阪大学新世紀レクチャー)には、以下のように書かれています。

原子の位置. a を単位とする. 直交座標で指定するときは0.5, 0.5, 0.5 などとかき, 基本ベクトルで指定するときは0.5a, 0.5b, 0.5cなどと書く.


hcp構造の計算では、実空間の(0.33333a 0.66667b 0.25c)と(0.66667a 0.33333b 0.75c)の2点に原子をおくことにします。原子位置の指定は結局、以下のようになります。

   0.33333a  0.66667b  0.25c    Co
0.66667a 0.33333b 0.75c Co


Fig.10に示したのがhcp構造の第一B.Z.で、table.5がブロッホスペクトル関数を計算するk点のパスです。
fcc構造やbcc構造での波数空間の座標は 2π/a を単位としていました。hcp構造の場合、多くの教科書やwebサイトではkx,kyは 2π/a を単位としていますが、kzだけは 2π/c を単位としています。今回のtable.5でもそのようになっていますし、入力ファイルでもそのようにします。

HEX.png
Fig.10: hcp構造の第一B.Z.(Wikipediaより転載)


kxkykz
Γ000
M01/√30
K1/31/√30
Γ000
A001/2
L01/√31/2
H1/31/√31/2
A001/2
table.5: hcp構造の第一B.Z. 単位は1-2列目は 2π/a , 3列目は 2π/c


これらを踏まえて作成した入力ファイルがCo_in.txtです。ブロッホスペクトル関数を計算するk点は、table.5のパスに沿った第一B.Z.の特徴的なk点の間をそれぞれ100等分したものを使いました。

hcp構造コバルトの全状態密度


以下にコバルトの全状態密度(total DOS)の計算結果を示します。
ニッケルのときと同様な処理をしたデータファイルがCoDOS_dat.txtで、gnuplotのスクリプトファイルがCoDOS_plt.txtです。アウトプットファイルに出力されるtotal DOSの単位は states/Ry/unit-cell です(だと思います)。bcc構造やfcc構造の金属と比較するためには、states/Ry/atom に統一するほうがいいと思います。hcpは原子を二つ配置したのでDOSをプロットする際に2で割ります。この処理はCoDOS_plt.txtのなかで行っています。


CoDOS.png
Fig.11: hcp構造コバルトの状態密度


hcp構造コバルトのバンド構造


コバルトのバンド構造に関しても、hcp構造特有の注意点があります。
入力ファイルを作る際に、k空間においてkzだけ単位が 2π/c であったことを思い出してください。逆格子空間での道のりを考えるときにこのことを考慮しなければいけません。

計算によって出力されたspcファイルには # の後にk空間の座標が書かれています。
実を言うと、この座標の単位はkx, ky, kzの全てに関して 2π/a です。

重要なことなのでもう一度書きます。
hcp構造では、入力ファイルのk空間の座標はkzだけ 2π/c の単位を持ちますが、出力ファイルのk空間の座標の単位は全て 2π/a です。

一見すると統一感がないようにも見えますが、この仕様はとても便利です。なぜなら、入力の時にはc/aの値を気にする事無く教科書どおりにk点の座標を与えることが出来、出力ファイルを処理するときもc/aの値を気にする事無く道のりの長さを求めることが出来るからです。したがって、fcc構造やbcc構造と同様にAkaiKKR-spc2ekでバンド図を描くための前処理が出来ます。

Γ0
M0.577350269
K0.910683603
Γ1.577350269
A1.883536487
L2.460886756
H2.79422009
A3.460886756
table.6: Γ点からの逆格子空間の道のり


table.6がc/a=1.633(理想的なhcp構造)のときの逆格子空間の特徴的なk点までの道のりです。以下に示したバンド図作成に利用したgnuplotスクリプトはCoSPC_plt.txtです。


CoUp.png
CoDown.pngFig.12: hcp構造コバルトのバンド構造


小口多美夫著 遷移金属のバンド理論に強磁性bcc構造Feと強磁性hcp構造Coの全状態密度の計算例が載っています。

関連エントリ




参考URL




付録


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


参考文献/使用機器




フィードバック



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

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


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


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

tag: AkaiKKR machikaneyama KKR 強磁性 状態密度 DOS バンド構造 分散関係 

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

LTspiceAkaiKKRScilabmachikaneyamaKKRPSoCOPアンプPICCPA強磁性モンテカルロ解析常微分方程式トランジスタode状態密度インターフェースDOSPDS5022ecaljスイッチング回路定電流半導体シェルスクリプトレベルシフト乱数HP6632A温度解析ブレッドボードI2CR6452A分散関係トランジスタ技術可変抵抗確率論数値積分反強磁性セミナー非線形方程式ソルバ絶縁バンドギャップ熱設計偏微分方程式バンド構造GW近似カオス三端子レギュレータLEDフォトカプラシュミットトリガISO-I2CA/DコンバータLM358USBカレントミラーTL431マフィンティン半径PC817C数値微分アナログスイッチ発振回路サーボ直流動作点解析74HC40532ちゃんねる標準ロジックチョッパアンプLDAアセンブラFFTbzqltyイジング模型ブラべ格子開発環境補間量子力学電子負荷BSchパラメトリック解析単振り子基本並進ベクトル熱伝導繰り返しGGAMaximaTLP621ewidthSMP相対論抵抗位相図ランダムウォークスピン軌道相互作用六方最密充填構造不規則合金FETコバルト失敗談QSGWcygwinスレーターポーリング曲線スイッチト・キャパシタラプラス方程式gfortranキュリー温度状態方程式条件分岐格子比熱TLP552LM555TLP521三角波NE555過渡解析FXA-7020ZRWriter509テスタ詰め回路MCUマントルダイヤモンドQNAPデータロガーガイガー管自動計測UPS井戸型ポテンシャルawk第一原理計算仮想結晶近似ブラウン運動差し込みグラフ平均場近似fsolve起電力熱力学OpenMPスーパーセル固有値問題最適化最小値VCAシュレディンガー方程式VESTAubuntu最大値面心立方構造PGAOPA2277L10構造非線型方程式ソルバ2SC1815fccフェルミ面等高線ジバニャン方程式ヒストグラム確率論マテリアルデザイン正規分布結晶磁気異方性interp1フィルタ初期値ウィグナーザイツ胞c/aルチル構造岩塩構造スワップ領域リジッドバンド模型edeltBaOウルツ鉱構造重積分SIC二相共存ZnOquantumESPRESSOCapSensegnuplotmultiplot全エネルギー固定スピンモーメントFSM合金ノコギリ波フォノンデバイ模型ハーフメタル半金属TeXifortTS-110不規則局所モーメントTS-112等価回路モデルパラメータ・モデルヒストグラムExcel円周率GimpトラックボールPC直流解析入出力文字列マンデルブロ集合キーボードフラクタル化学反応三次元Realforce縮退日本語最小二乗法関数フィッティング疎行列シンボル線種ナイキスト線図陰解法負帰還安定性熱拡散方程式EAGLECrank-Nicolson法連立一次方程式P-10クーロン散乱Ubuntu境界条件MBEHiLAPW軸ラベルトランスCK1026MAS830L凡例PIC16F785LMC662AACircuit両対数グラフ片対数グラフグラフの分割specx.f

最新コメント
リンク

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