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 バンド構造 分散関係 

AkaiKKRでダイヤモンド型構造半導体

AkaiKKR(Machikaneyama)を利用してダイヤモンド型構造の半導体としてシリコンとダイヤモンド(炭素)の状態密度やバンド構造を計算しました。ダイヤモンド構造の計算は基本的にfcc構造として行いますが、原子のない位置に明示的に空孔を入れて計算しないと計算の精度が下がることがあるようです。

計算されたバンド構造から、シリコンやダイヤモンドは間接遷移型の半導体であることが読み取れます。計算されたバンドギャップの大きさは、実際に測定されているものよりもかなり小さい値を示しました。一般的に局所密度近似では、バンドギャップの大きさを正しく計算することはとても難しく、過小評価する傾向があるということです。

002_20130202154314.png 003_20130202154313.png 004_20130202154313.png 005_20130202154313.png


AkaiKKRで不規則NiFe合金の強磁性AkaiKKRで不規則NiMn合金の分散関係では、AkaiKKR(Machikaneyama)を利用して面心立方(face-centered cubic:fcc)構造の金属の状態密度やバンド分散を計算しました。

ねがてぃぶろぐは電子工作ブログということになっているので(?)、金属だけでなく半導体の状態密度やバンド分散を計算します。
色々な化合物が半導体として知られていますが、今回は単体元素で半導体となるダイヤモンド型構造のシリコンやダイヤモンドを扱います。

ダイヤモンド型結晶構造


以下に示すのはVESTAを利用して、ダイヤモンドの結晶構造を描画したものです。
結晶構造ファイルの作り方は情報科学演習 2011 4. 発展問題 ~氷の構造などに分かりやすい解説がなされていますが、今回は結晶構造ギャラリーに公開されている4-diamond.cifを利用させていただきました。

001_20130202051148.png
Fig.1: ダイヤモンドの結晶構造


ダイヤモンド型の結晶構造は、一見するとたくさんの原子を含む複雑なもののように見えますが、よくよく眺めてみると、2つのfcc構造の結晶をx,y,zの各方向に1/4ずつずらして重ねたものであることが見て取れます。

したがって、AkaiKKR(Machikaneyama)の入力ファイルは、基本的には以下のような物になります。(が、後述する通りもう一工夫必要です。)

 go  data/si
fcc 10.26 1.0 1.0, , , ,
0.001 1.5 sra vwnasa nmag init
update 4 100 0.02
1
Si 1 0 0.0 2 14 100
2
0.00000 0.00000 0.00000 Si
0.25000 0.25000 0.25000 Si


しかしながら、実際にこの入力で状態密度を描いてみると、おおざぱな状態密度の形自体は何となくあっているようにも見えるのですが、本来ならバンドギャップであるはずの領域にも有限の状態密度があるような図が出来てしまいます。

こういった問題に対する処方箋がMachikaneyama2000 の使用に関するメモを注意深く読むと見つかります。

以上がYMnAl の場合であるが,inmn6mn0as0as の場合もinput file の意味を探ることは容易である.ここでは原子番号としてZ=0 の指定があるが,これは原子空孔を意味する.
ダイヤモンド構造やzincblende 構造の結晶の場合オープン構造であり,計算の精度をあげるために原子の無い位置にも空孔を入れているが,必ずしも必要なわけではない.


古いversionであるcpa2002v008bのinフォルダの中のZnSの入力ファイルがあるので、それを参考に原子空孔を明示的に書き込んだシリコンの入力ファイルを作ります。

 go  data/si
fcc 10.26 1.0 1.0, , , ,
0.001 1.5 sra vwnasa nmag init
update 4 100 0.02
2
Si 1 0 0.0 2 14 100
Vc 1 0 0.0 2 0 100
4
0.00000 0.00000 0.00000 Si
0.25000 0.25000 0.25000 Si
0.75000 0.75000 0.75000 Vc
0.50000 0.50000 0.50000 Vc


ダイヤモンド型構造の逆格子


バンド分散関係を計算するためには、第一Brillouin Zone(第一B.Z.)の形を知らなければなりません。
先ほど、ダイヤモンド構造の結晶構造は2つのfcc構造を1/4ずつずらして重ね合わせたものだと書きました。同様に、ダイヤモンド型構造の第一B.Z.は、fcc構造のものと同じです。このことについては金持徹著 固体電子論に噛み砕いて書かれています。

したがって、AkaiKKRで不規則NiMn合金の分散関係でバンド構造を計算するときに指定したものと全く同じk点を使うことが出来ます。

最終的な入力ファイルは、Si_in.txtDiamond_in.txtとなりました。

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)


状態密度とバンド構造


以下に計算された状態密度とバンド構造を示します。

002_20130202154314.png

Fig.2: シリコンの状態密度

003_20130202154313.pngFig.3: シリコンのバンド構造


シリコンは状態密度(Fig.2)からバンドギャップを持ち、バンド構造(Fig.3)から価電子帯の上端と伝導体の底が異なる波数ベクトルのところに存在するため間接遷移が起こるタイプの半導体であることが分かります。

以下に示すダイヤモンド(炭素)でも同様に、間接遷移型の半導体であることが読み取れます。

004_20130202154313.png

Fig.4: ダイヤモンド(炭素)の状態密度

005_20130202154313.png
Fig.5: ダイヤモンド(炭素)のバンド構造


シリコンとダイヤモンドのバンドギャップの大きさを比べると、ダイヤモンドのほうがはるかに大きいことが読み取れます。Wikipediaのバンドギャップの一覧によると、実際に測定されたバンドギャップの大きさでもダイヤモンドのバンドギャップはシリコンの5倍程度大きいです。

元素実際のバンドギャップ(eV)計算されたバンドギャップ(eV)
シリコン1.11~0.4
ダイヤモンド5.5~4
table.1: 計算と実験によるバンドギャップの大きさ


しかし、その定量的な大きさは、どちらの計算結果も現実の値を過小評価しています。また、その過小評価の程度も、シリコンではかなり大きく、ダイヤモンドでは比較的小さいと、扱う物質しだいで変わってしまうようです。この問題は、第一原理計算の業界にとって共通する重要な課題であるようで、本当に色々なところで議論されています。(例えば第一原理計算 (解析編)Wikipdeiaの該当箇所)

AkaiKKR(Machikaneyama)の掲示板にも下記の通り該当する書き込みがあります。

Can KKR methoud give a better band gap for semiconductor compared with LDA?
Posted on : August 10, 2012 (Fri) 22:18:55
by tfl

Dear all,

I want to calculate the band structure of some semiconductor. Using traditional LDA method, it is found that the band gap are obviously underestimated. Would you please to tell me that KKR method can provide more accurate band gap values for semiconductor?


[Re:01] Can KKR methoud give a better band gap for semiconductor compared with LDA?
Posted on : December 09, 2012 (Sun) 10:13:46
by Koji Kobashi

Since Akai-KKR uses the LDA, the calculated bandgaps are lower than the actual values of semiconductors. The ASA improves it but only slightly.


と言うことなので、今回の計算ではASA(Atomic Sphere Approximation, 原子球近似)を利用しています。

関連エントリ




参考URL




付録


このエントリで使用したAkaiKKRの入力ファイルを添付します。



参考文献/使用機器




フィードバック



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

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


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


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

tag: AkaiKKR machikaneyama KKR 半導体 バンドギャップ 状態密度 DOS バンド構造 分散関係 

AkaiKKRで不規則NiMn合金の分散関係

密度汎関数法の発展で紹介されているニッケル-マンガン合金のブロッホスペクトル関数(バンド構造, エネルギー分散関係, E-k曲線)をAkaiKKR(Machikaneyama)を用いて計算し、不規則性の効果によりにじんだバンド分散関係の図が得られました。




AkaiKKRで不規則NiFe合金の強磁性では、AkaiKKR(Machikaneyama)を用いてfcc構造のニッケルと鉄を端成分とする固溶体合金の電子の状態密度を計算し、鉄の濃度の上昇により75%付近で強磁性が消失することを確かめました。

今回は、状態密度だけでなく電子のバンド構造(エネルギー分散関係,E-k曲線)を描画してみます。計算対象の物質は、密度汎関数法の発展で紹介されているニッケル-マンガン合金です。

入力ファイル


セルフコンシステント計算(go)と状態密度計算(dos)のための入力部分はAkaiKKRで不規則NiFe合金の強磁性のものとほぼ同じで固溶させる原子の番号が26(Fe)から25(Mn)に変更された程度です。前回は、たくさんの組成に対して計算を行いましたが、今回は純粋なニッケルとニッケルに15%マンガンを固溶させたものの2つの組成だけを計算しました。

その代わり電子の分散関係を描画するためのブロッホスペクトル関数の計算(spc)を行います。
この計算には、計算を行うk点をk空間での座標で指定する必要があります。このときの座標は 2π/a で規格化されているので格子定数 a が異なる結晶でも同じ結晶構造なら使いまわすことが出来ます。


001_20130120185338.png

Fig.1: fcc構造の第一Brillouin Zone(Wikipediaより転載)

W1/200
L1/21/21/2
Γ000
X001
W1/201
K3/403/4
table.1: 計算するk点のパス


k点のパスの選び方はどう取るのがよいのかいまいち良く分からないのですが、table.1のk点を通るようにし、それぞれの特徴的なk点の間を100点ずつ計算することとしました。(密度汎関数法の発展のパスのとり方と違ってしまいました・・・すみません。)

最終的な入力ファイルは、NiMn_in.txtとなりました。

specx.fの編集


おそらくCygwinでAkaiKKR(Machikaneyama)のspecx.fの設定では、nk1xかnk3xのあたりが不足してしまうと思います。私が今回の入力ファイルを計算したときのspecx.fの設定は下記の通りです。

     & (natmmx=4, ncmpmx=4, msizmx=198, mxlmx=3, nk1x=2200, nk3x=2688,
& msex=201, ngmx=15, nrpmx=650, ngpmx=650, npmx=350, msr=400)

specx.fを編集したらふたたびmakeします。

gnuplotで分散関係(E-k曲線)の描画


以上のような設定で計算を実行するとdataディレクトリにni.spcやni85mn15.spcといったファイルが出来ます。(ファイル容量が5MBと大きいため、fc2blogにはアップロードできませんでした。)

これらのファイルは、入力ファイルで指定したk点の座標(波数ベクトル:k)ごとに、エネルギー(E)とブロッホスペクトル関数(A(k,E))の値が書かれています。これを、横軸に波数ベクトルk、縦軸にエネルギーEをとり、ブロッホスペクトル関数A(k,E)を色の濃淡で表すことを考えます。

これをどのように実現するのがベターなのか、私には自信がありませんが、gnuplotではカラーマップを使った2次元プロット(pm3d map)で描画できそうです。

入力のデータ形式は3次元データのフォーマットの様にx,y,zの組で与える必要があります。ni.spcをみるとyとzに相当するエネルギーEとブロッホスペクトル関数A(k,E)の値は既に組になって書かれているので、1列目にx軸の数値を補ってやる必要があります。

k点のパスは W → L → Γ → X → W → K なので、この順に道のりの長さをx軸の数値にします。
三次元空間の2点間の距離は

\sqrt{(x_2 - x_1)^{2} + (y_2 - y_1)^{2} + (z_2 - z_1)^{2}}


なのでパスの最初の点であるWからの道のりの長さはそれぞれtable.2のようになります。ここの計算は簡単なスクリプトで自動化できるのだろうとは思うのですが、差し当たりExcelでちまちま計算しました。一度計算すれば(結晶構造が同じでエネルギーメッシュやk点の分解能が同じなら)コピー&ペーストで使いまわせます。(fc2blogの容量制限でアップローで出来ずすみません。)時間が出来たらなんらかのスクリプトを書こうと思います。

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


目盛見出しの「目盛見出しを任意の文字に変更する」の方法を使えばx軸の値を文字列に置き換えられます。ラベルにα,βの様なギリシャ文字を使いたいの方法でギリシャ文字を含むラベルがかけます。リンク先の例ではpostscriptのターミナルでしか使えないようにも見える記述ですが、
gnuplot > set terminal windows enhanced

gnuplot > set terminal png enhanced
も可能です。
table.2の対応を考慮して
set xtics   ("{W}" 0.000000, "{L}" 0.707107, "{/Symbol G}" 1.57313, "{X}" 2.57313, "{W}" 3.07313, "{K}" 3.42669)
としました。

計算結果: 純ニッケル


以降は計算結果です。まずは純ニッケルについて。


002_20130120185338.png

Fig.2: Niの状態密度


Fig.2は、AkaiKKRで不規則NiFe合金の強磁性のときと全く同じ計算で、上半分が上向きスピンの状態密度、下半分が下向きスピンの状態密度で、エネルギーの原点がフェルミ面です。

003_20130120185338.pngFig.3: Niの上向きスピンの分散関係004_20130120185337.pngFig.4: Niの下向きスピンの分散関係


Fig.3-4は、ニッケルのバンド分散関係を上向きスピン、下向きスピン別々にプロットしたものです。どちらも似たような形をしていますが、フェルミ面が曲線の傾きが小さいバンドの密集しているところを横切っている(下向きスピン)かそれよりも高いところを横切っているか(上向きスピン)の違いがわかります。これは、状態密度の図(Fig.2)でフェルミ面が、下向きスピンでは状態密度の高いところにあり、上向きスピンでは低いところにあることと対応しています。

計算結果: Ni0.85Mn0.15


次にニッケルマンガン合金の計算結果です。


005_20130120185337.png

Fig.5: Ni0.85Mn0.15の状態密度


ニッケルの状態密度は、鋭いピークを持った形状をしていましたが、ニッケルマンガン合金は俄然ブロードな形となっています。

006_20130120185337.pngFig.6: Ni0.85Mn0.15の上向きスピンの分散関係007_20130120185347.pngFig.7: Ni0.85Mn0.15の下向きスピンの分散関係


ブロードな状態密度の形状は、分散関係の図でバンドを示す線がにじんでしまっていることと対応付けられます。(Fig.3-4の段階で既に線がにじんでしまっているように見えますが、これは本来シャープな曲線のはずです。にじんでしまっているように見えるのは、おそらく、私の技術的な問題です。)

密度汎関数法の発展には、このことに関して以下のようにあります。

Ni0.85Mn0.15合金のエネルギー分散の図が本を水でぬらしてしまったようにところどころにじんでいるのがわかる。これは印刷の失敗でも汚れでもなくて、まさに不規則性の効果が現れたものである。Niの場合には系が規則的なために結晶運動量hk/2πが良い量子数となって固有エネルギーを決定するのに対して、Ni0.84Mn0.15合金では系がもはや規則的でないために結晶運動量hk/2πについて固有状態が得られずブロッホスペクトル関数が広がってしまうのである。


関連エントリ




参考URL




付録


このエントリで使用したAkaiKKRの入力ファイルを添付します。


参考文献/使用機器





フィードバック



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

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


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


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

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

AkaiKKRで不規則NiFe合金の強磁性

計算機マテリアルデザイン入門 (大阪大学新世紀レクチャー)で紹介されている不規則NiFe合金の状態密度をAkaiKKR(Machikaneyama)を利用して計算しました。

NiFe.gif
Fig.1: NiFe合金の状態密度


その結果、NiFe合金は65%以上のFe濃度から急激に(特に上向きスピンの)状態密度の形状が変化し、75%Fe程度で強磁性を失うことが確認できました。


コヒーレントポテンシャル近似(CPA)


第一原理バンド計算法は、結晶構造の周期性を利用して計算を行っているため、ランダムに原子が置換している不規則合金の計算を苦手としています。比較的よく利用されている解決策は、複数の単位格子を合わせて大きな単位格子(スーパーセル)をつくり、その中の原子を不純物原子に置き換えて計算を行うというものです。しかしながら、この方法では希薄な不純物合金の計算をするためには大きなスーパーセルが必要となるため計算量が増えるという問題点もあります。

コヒーレントポテンシャル近似(coherent potential approximation: CPA)は、こういった問題を解決するために考案された近似法で、TB(tight-banding)法やKKR(Korringa-Kohn-Rostoker)法といったバンド計算法と組み合わせて不規則合金の電子構造を計算することが出来ます。CPAの解説記事としてはコヒーレント・ポテンシャル近似 CPA-1--2-CPAの応用(米沢富美子 固体物理 1971, 1972)などがあるようです。AkaiKKR(Machikaneyama)は、大阪大学の赤井久純先生が公開されているKKR-CPAを用いた第一原理計算パッケージです。

CygwinでAkaiKKR(Machikaneyama)ではWindows上でAkaiKKR(Machikaneyama)を利用する方法を書きました。今回はCPAの応用として良く挙げられているfcc構造のニッケル-鉄合金の状態密度計算を行いました。(参考: 計算機マテリアルデザイン入門 (大阪大学新世紀レクチャー))


入力ファイル


入力ファイルには、複数の計算設定を順に書いておくことによって、一挙に計算をさせることが出来ます。下記にNiとNi95Fe5の計算部分を示します。実際の入力ファイルは5%ずつ鉄の濃度を振ってfcc構造の純鉄まで計算を行いました。(NiFe_in.txt)

c ***********************************************
c pure Nickel
c ***********************************************
c *** Self-consistent calculation ***
go data/ni
fcc 6.67 , , , , , ,
0.001 1.2 nrl mjw mag init
update 4 100 0.02
1
Ni 1 0 0 2 28 100
1
0 0 0 Ni
c *** Density of States (DOS) calculation ***
dos data/ni
fcc 6.67 , , , , , ,
0.001 1.2 nrl mjw mag 2nd
update 13 100 0.02
1
Ni 1 0 0 2 28 100
1
0 0 0 Ni

c ***********************************************
c Nickel + 5% Iron
c ***********************************************
c *** Self-consistent calculation ***
go data/ni95fe5
fcc 6.67 , , , , , ,
0.001 1.2 nrl mjw mag init
update 4 100 0.02
1
NiFe 2 0 0 2 28 95
26 5
1
0 0 0 NiFe
c *** Density of States (DOS) calculation ***
dos data/ni95fe5
fcc 6.67 , , , , , ,
0.001 1.2 nrl mjw mag 2nd
update 13 100 0.02
1
NiFe 2 0 0 2 28 95
26 5
1
0 0 0 NiFe


格子定数はCygwinでAkaiKKR(Machikaneyama)のときと同様に a = 6.67 (bohr)としFeの濃度にかかわらず一定となるようにしました。描画する状態密度をきれいにするためにdos計算のbzqltyは前回よりも大きくしてあります。bzqlty=13は、前回のspecx.fの設定を用いたfcc構造の計算では最大の値です。13よりも大きい数を設定すると、計算の途中でspecx.exeが停止します。したがって、これ以上見栄えの良いDOSを描きたければ、specx.fのnk1mxとnk3mxを大きな値に修正した後、再びmakeをする必要があります。

結果


以下に、ニッケルに0-100%の鉄を固溶させた合金を5%ごとに計算した状態密度を示します。


000Fe.png
Fig.2: Pure Nickel

005Fe.png
Fig.3: Nickel + 5% Iron

010Fe.png
Fig.4: Nickel + 10% Iron

015Fe.png
Fig.5: Nickel + 15% Iron

020Fe.png
Fig.6: Nickel + 20% Iron

025Fe.png
Fig.7: Nickel + 25% Iron

030Fe.png
Fig.8: Nickel + 30% Iron

035Fe.png
Fig.9: Nickel + 35% Iron

040Fe.png
Fig.10: Nickel + 40% Iron

045Fe.png
Fig.11: Nickel + 45% Iron

050Fe.png
Fig.12: Nickel + 50% Iron

055Fe.png
Fig.13: Nickel + 55% Iron

060Fe.png
Fig.14: Nickel + 60% Iron

065Fe.png
Fig.15: Nickel + 65% Iron

070Fe.png
Fig.16: Nickel + 70% Iron

075Fe.png
Fig.17: Nickel + 75% Iron

080Fe.png
Fig.18: Nickel + 80% Iron

085Fe.png
Fig.19: Nickel + 85% Iron

090Fe.png
Fig.20: Nickel + 90% Iron

095Fe.png
Fig.21: Nickel + 95% Iron

100Fe.png
Fig.22: Pure fcc Iron


強磁性の消失


ニッケル-鉄合金の強磁性は、実験から約74%Feで消失することが知られています。
今回計算した結果の図は、上半分が上向きスピン、下半分が下向きスピンの状態密度をあらわしています。強磁性の消失は、75%Fe(Fig.17)で上下対称な図になっていることから読み取れます。

不規則合金の電子構造を理解するにあたってリジッドバンドモデル(固定バンドモデル、Rigid band model)が良く持ち出されるが、コヒーレント・ポテンシャル近似と合金の強磁性(金森順次郎 (1972) 固体物理)によると不規則NiFe合金の強磁性の消失は、一見リジッドバンドモデルで説明が付く様に見えるが、平均電子数の変化によるフェルミ準位の位置の変化が本質では無いとのことです。

関連エントリ




参考URL




付録


このエントリで使用したAkaiKKRの入力ファイルを添付します。


参考文献/使用機器






フィードバック



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

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


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


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

tag: AkaiKKR machikaneyama KKR CPA 状態密度 DOS 強磁性 不規則合金 

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

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

最新コメント
リンク

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