- 1次元
- 2次元
- 3次元
規格化されたGaussianは次の形で与えられる。
規格化されている(積分値が1になっている)ことは、次のようにして確かめられる。
を用いた表記は、統計学的な処理を行う際に便利である(と想像される)が、直感的にどれくらいの幅を持つのかが分かりにくい。
そこで、半値半幅を用いた表記を考えることにする。
の定義は、次で与えられる。
したがって、の変換は、
数値的に見積もると、となる。半値半幅はよりも気持ちちょっと大きい程度である。
元のGaussianの式をで表すと、
(おそらくLorentz関数でもそうだと思うが、)Gaussian型のピーク形状の場合、半値全幅(文献によっては、を半値全幅の表記に用いる場合もあるため、注意が必要)を分解能として扱うため、Gaussian broadening等で理論値に実験誤差を入れ込む場合には、半値幅の用いた表記の方が便利である。
通常の結晶構造であれば、データベースか何かから引っ張って来れるが、それをベースにした超周期構造を作ろうと思うと、対称性(空間群)が変わるため、手でやるには歯が立たない。
そこで、入力した原子座標に対して対称性を見つけてくれる spglib を使って、反強磁性でのNiOの対称性を決定することを試みる。
spglibはC言語、FORTRAN, pythonに対応しているが、いずれにおいても「実行ファイル(今風に言えばアプリケーション)」ではなく、(~libとあるように)ライブラリーとして提供されている。そのためGUIのようなものに結晶構造を入れるようなものではないため、多少プログラミングに慣れている必要がある。
Spglib — spglib 1.12.0 documentation
個人的には、pythonが簡便でオススメである。なので、ここではpythonでの例を紹介する。
Spglib for Python — spglib 1.12.0 documentation
pythonバージョンのインストールは、pipですぐ入る。
pip install spglib
これを使って、まずは通常のNiOの結晶構造で試してみる。
import numpy as np import spglib # lattice constant a = 4.19 # structural information ( [lattice constants], [coordination], [atomic numbers] ) # Caution: it should be a tuple nio = \ ( [ [ a, 0, 0 ], [ 0, a, 0 ], [ 0, 0, a ] ], [ [ 0, 0, 0 ], # Ni [ 0, 0.5, 0.5 ], # Ni [ 0.5, 0, 0.5 ], # Ni [ 0.5, 0.5, 0 ], # Ni [ 0.5, 0, 0 ], # O [ 0, 0.5, 0 ], # O [ 0, 0, 0.5 ], # O [ 0.5, 0.5, 0.5 ] # O ], [28,] * 4 + [8,] * 4 ) # just space group print( spglib.get_spacegroup( nio ) ) #Fm-3m (225) # primitive cell print( spglib.find_primitive( nio ) ) #(array([[0. , 2.095, 2.095], # [2.095, 0. , 2.095], # [2.095, 2.095, 0. ]]), array([[0. , 0. , 0. ], # [0.5, 0.5, 0.5]]), array([28, 8], dtype=int32))
のsuper cell でもちゃんと元の構造の対称性を返してくれる。
# some coordinates are just added 1 to discribe them in the other unit cell # therefore, lattice constans and coordinates should be multipled and divided by 2, respectively nio_sc \ = [ [ ( a, 0, 0 ), ( 0, a, 0 ), ( 0, 0, a ) ], [ [ 0, 0, 0 ], # Ni 0, 0, 0 [ 0, 0.5, 0.5 ], [ 0.5, 0, 0.5 ], [ 0.5, 0.5, 0 ], [ 1, 0, 0 ], # Ni 1, 0, 0 [ 1, 0.5, 0.5 ], [ 1.5, 0, 0.5 ], [ 1.5, 0.5, 0 ], [ 0, 1, 0 ], # Ni 0, 1, 0 [ 0, 1.5, 0.5 ], [ 0.5, 1, 0.5 ], [ 0.5, 1.5, 0 ], [ 0, 0, 1 ], # Ni 0, 0, 1 [ 0, 0.5, 1.5 ], [ 0.5, 0, 1.5 ], [ 0.5, 0.5, 1 ], [ 1, 1, 0 ], # Ni 1, 1, 0 [ 1, 1.5, 0.5 ], [ 1.5, 1, 0.5 ], [ 1.5, 1.5, 0 ], [ 1, 0, 1 ], # Ni 1, 0, 1 [ 1, 0.5, 1.5 ], [ 1.5, 0, 1.5 ], [ 1.5, 0.5, 1 ], [ 0, 1, 1 ], # Ni 0, 1, 1 [ 0, 1.5, 1.5 ], [ 0.5, 1, 1.5 ], [ 0.5, 1.5, 1 ], [ 1, 1, 1 ], # Ni 1, 1, 1 [ 1, 1.5, 1.5 ], [ 1.5, 1, 1.5 ], [ 1.5, 1.5, 1 ], [ 0.5, 0, 0 ], # O 0, 0, 0 [ 0, 0.5, 0 ], [ 0, 0, 0.5 ], [ 0.5, 0.5, 0.5 ], [ 1.5, 0, 0 ], # O 1, 0, 0 [ 1, 0.5, 0 ], [ 1, 0, 0.5 ], [ 1.5, 0.5, 0.5 ], [ 0.5, 1, 0 ], # O 0, 1, 0 [ 0, 1.5, 0 ], [ 0, 1, 0.5 ], [ 0.5, 1.5, 0.5 ], [ 0.5, 0, 1 ], # O 0, 0, 1 [ 0, 0.5, 1 ], [ 0, 0, 1.5 ], [ 0.5, 0.5, 1.5 ], [ 1.5, 1, 0 ], # O 1, 1, 0 [ 1, 1.5, 0 ], [ 1, 1, 0.5 ], [ 1.5, 1.5, 0.5 ], [ 1.5, 0, 1 ], # O 1, 0, 1 [ 1, 0.5, 1 ], [ 1, 0, 1.5 ], [ 1.5, 0.5, 1.5 ], [ 0.5, 1, 1 ], # O 0, 1, 1 [ 0, 1.5, 1 ], [ 0, 1, 1.5 ], [ 0.5, 1.5, 1.5 ], [ 1.5, 1, 1 ], # O 1, 1, 1 [ 1, 1.5, 1 ], [ 1, 1, 1.5 ], [ 1.5, 1.5, 1.5 ] ], [28,] * 32 + [8,] * 32 ] # check whether there is dupulication for i, p1 in enumerate( nio_sc[1] ): for p2 in nio_sc[1][ i+1 : ]: if p1 == p2: print( p1 ) # lattice constants should be twice in 2*2*2 super cell nio_sc[0] = [ list( a_ ) for a_ in np.array( nio_sc[0] ) * 2 ] # coordinate should be divided by 2 nio_sc[1] = [ list( p_ ) for p_ in np.array( nio_sc[1] ) / 2 ] nio_sc = tuple( nio_sc ) print( spglib.get_spacegroup( nio_sc ) ) print( spglib.find_primitive( nio_sc ) )
座標の重複があると、'None'が返って来るため、チェックする必要がある。目で見て一生懸命探すのには限界がある。
最後に反強磁性だが、格子定数と座標はsuper cellのままでよく、原子番号で反強磁性秩序の情報を入れる。仮に磁化軸方向と平行な原子をNi(28)として、反平行なものを仮にPd(46)とする。酸素はそのままで良い。
NiOの反強磁性磁化軸は[111]方向なので、それを考慮して原子番号を対応させていく。
nio_af \ = [ [ ( a, 0, 0 ), ( 0, a, 0 ), ( 0, 0, a ) ], [ [ 0, 0, 0 ], # Ni 0, 0, 0 [ 0, 0.5, 0.5 ], [ 0.5, 0, 0.5 ], [ 0.5, 0.5, 0 ], [ 1, 0, 0 ], # Ni 1, 0, 0 [ 1, 0.5, 0.5 ], [ 1.5, 0, 0.5 ], [ 1.5, 0.5, 0 ], [ 0, 1, 0 ], # Ni 0, 1, 0 [ 0, 1.5, 0.5 ], [ 0.5, 1, 0.5 ], [ 0.5, 1.5, 0 ], [ 0, 0, 1 ], # Ni 0, 0, 1 [ 0, 0.5, 1.5 ], [ 0.5, 0, 1.5 ], [ 0.5, 0.5, 1 ], [ 1, 1, 0 ], # Ni 1, 1, 0 [ 1, 1.5, 0.5 ], [ 1.5, 1, 0.5 ], [ 1.5, 1.5, 0 ], [ 1, 0, 1 ], # Ni 1, 0, 1 [ 1, 0.5, 1.5 ], [ 1.5, 0, 1.5 ], [ 1.5, 0.5, 1 ], [ 0, 1, 1 ], # Ni 0, 1, 1 [ 0, 1.5, 1.5 ], [ 0.5, 1, 1.5 ], [ 0.5, 1.5, 1 ], [ 1, 1, 1 ], # Ni 1, 1, 1 [ 1, 1.5, 1.5 ], [ 1.5, 1, 1.5 ], [ 1.5, 1.5, 1 ], [ 0.5, 0, 0 ], # O 0, 0, 0 [ 0, 0.5, 0 ], [ 0, 0, 0.5 ], [ 0.5, 0.5, 0.5 ], [ 1.5, 0, 0 ], # O 1, 0, 0 [ 1, 0.5, 0 ], [ 1, 0, 0.5 ], [ 1.5, 0.5, 0.5 ], [ 0.5, 1, 0 ], # O 0, 1, 0 [ 0, 1.5, 0 ], [ 0, 1, 0.5 ], [ 0.5, 1.5, 0.5 ], [ 0.5, 0, 1 ], # O 0, 0, 1 [ 0, 0.5, 1 ], [ 0, 0, 1.5 ], [ 0.5, 0.5, 1.5 ], [ 1.5, 1, 0 ], # O 1, 1, 0 [ 1, 1.5, 0 ], [ 1, 1, 0.5 ], [ 1.5, 1.5, 0.5 ], [ 1.5, 0, 1 ], # O 1, 0, 1 [ 1, 0.5, 1 ], [ 1, 0, 1.5 ], [ 1.5, 0.5, 1.5 ], [ 0.5, 1, 1 ], # O 0, 1, 1 [ 0, 1.5, 1 ], [ 0, 1, 1.5 ], [ 0.5, 1.5, 1.5 ], [ 1.5, 1, 1 ], # O 1, 1, 1 [ 1, 1.5, 1 ], [ 1, 1, 1.5 ], [ 1.5, 1.5, 1.5 ] ], [ 28, 46, 46, 46 ] + [ 46, 28, 28, 28 ] * 3 + [ 28, 46, 46, 46 ] * 3 + [ 46, 28, 28, 28 ] + [ 8, ] * 32 ] # check whether there is dupulication for i, p1 in enumerate( nio_af[1] ): for p2 in nio_af[1][ i+1 : ]: if p1 == p2: print( p1 ) # lattice constants should be twice in 2*2*2 super cell nio_af[0] = [ list( a_ ) for a_ in np.array( nio_af[0] ) * 2 ] # coordinate should be divided by 2 nio_af[1] = [ list( p_ ) for p_ in np.array( nio_af[1] ) / 2 ] nio_af = tuple( nio_af ) print( spglib.get_spacegroup( nio_af ) ) #R-3m (166) print( spglib.find_primitive( nio_af ) ) #(array([[ 1.48138871, 0.85528017, 4.83819526], # [-1.48138871, 0.85528017, 4.83819526], # [ 0. , -1.71056034, 4.83819526]]), array([[-3.70074342e-17, 0.00000000e+00, 0.00000000e+00], # [ 5.00000000e-01, 5.00000000e-01, 5.00000000e-01], # [ 7.50000000e-01, 7.50000000e-01, 7.50000000e-01], # [ 2.50000000e-01, 2.50000000e-01, 2.50000000e-01]]), array([28, 46, 8, 8], dtype=int32))
対称性がFm-3mからR-3mに落ちているのが分かる。
ここで、結晶軸がベクトル表記のままだと結晶構造作成プログラムでは不便な場合があるので、格子定数の形式に変換すると、
prim_nio_af = spglib.find_primitive( nio_af ) axes = prim_nio_af[0] abc = [ np.linalg.norm( axis ) for axis in axes ] print( abc ) #[5.131681011130759, 5.131681011130759, 5.131681011130759] # alpha beta gamma abg =[] abg.append( np.degrees( np.arccos( np.dot( axes[1], axes[2] ) / np.linalg.norm(axes[1]) / np.linalg.norm(axes[2]) ) ) ) abg.append( np.degrees( np.arccos( np.dot( axes[2], axes[0] ) / np.linalg.norm(axes[2]) / np.linalg.norm(axes[0]) ) ) ) abg.append( np.degrees( np.arccos( np.dot( axes[0], axes[1] ) / np.linalg.norm(axes[0]) / np.linalg.norm(axes[1]) ) ) ) print( abg ) #[33.55730976192072, 33.55730976192072, 33.55730976192072]
これで、第一原理計算プログラムで反強磁性が計算できるようになった。
何も考えずにプログラムを作っていたら、配列の次元が7個を超えるとコンパイルできないらしい。
program test_array7 real( kind( 0d0 ) ) :: r( 1, 2, 3, 4, 5, 6, 7, 8 ) end program test_array7
これを gfortran でコンパイルすると
test_array7.f90:3:46: real( kind( 0d0 ) ) :: r( 1, 2, 3, 4, 5, 6, 7, 8 ) 1 Error: Array specification at (1) has more than 7 dimensions
添字が11個ある変数の計算しようとしていたので、これはちょっと困った。。。
構造体とかを使って、何とか逃げるしかないか。。。
ちなみに構造体なら7個以上の要素があっても大丈夫であった。
program test_type7 implicit none type :: test real( kind( 0d0 ) ) :: r1 real( kind( 0d0 ) ) :: r2 real( kind( 0d0 ) ) :: r3 real( kind( 0d0 ) ) :: r4 real( kind( 0d0 ) ) :: r5 real( kind( 0d0 ) ) :: r6 real( kind( 0d0 ) ) :: r7 real( kind( 0d0 ) ) :: r8 end type type( test ) :: t end program test_type7
自分で毎回忘れるので、まとめておく。
Fourier変換の符号や係数は、変換と逆変換とで対になっていれば良いため(語弊があるが)不定性があり、分野によって異なる(気がしている)。
一応、自分ルールとして、
と考える。
量子力学的な観点でここでは考えていく。
運動量およびエネルギーの固有関数は、それぞれ次のような平面波で表される。
したがって、任意のやの関数をこれらの固有関数の重なりで表現すると思うと、
という位相を採用するのが適当であろう。
注意するべきは、これらは「逆変換」の式である。何となく「逆」とついていると後から決まるような感じがするが、ここではそうではない。
これによって、「順」方向側の変換の位相が定まる。
「基準となる関数から、知りたい平面波の係数を抜き出す」というイメージになるだろうか。
残るは係数であるが、
という関係であるため、一意には決まらない。
ここで、Fourier変換のモチベーションに戻ると、基準であるの関数はわかっているから、知りたいのはである。
そのため、, という式をパッと使いたくなるし、頻出することが予想される。それゆえの「順」方向である。
したがって、順方向側の係数を簡単にして、逆方向側にメンドクサイのを押し付ける方が便利(なように感じるの)である。
よって、
最終形は、
あるセル内の散乱によるセルT行列を と置くと、全てのサイトによる散乱Green関数 は次のように書ける。
ここで、 は自由散乱状態におけるGreen関数である。また、セルに属するサイトをという風に表した。
セルT行列 はサイトT行列 を用いて次のように書ける。
セルT行列に関する詳しい話は前回まとめた。
koideforest.hatenadiary.com
ここで、系が並進対称性を持ち、ユニットセルが周期的に繰り返しているとすると、セルT行列はセルに依存しなくなる(, )。
Fourier級数展開によって、セルの添字を結晶運動量に変換する。
はサイト間の差にしか依存しないため、
したがって、
よって
逆行列計算に計算コストが掛かるため、ユニットセル内の原子だけに抑えられるのはかなりのメリットだと思われる。
T行列は元々サイトポテンシャルの展開で以下の様に定義されている。
しかし、「違うサイトへの伝搬」と「同じサイト内での多重散乱」が混在しており、実際の計算には非常に不便である。
そこで、同じサイトのポテンシャル散乱を繰り込んだサイトT行列を定義すると、サイト内の多重散乱とサイト間の伝搬を分けることができる。これがいわゆる多重散乱展開である。
サイトT行列間の伝搬の際には、かならず異なるサイトへ行くため、散乱のイメージが付きやすい。
注意として、あるサイトから他サイトに伝搬した後、もう一度同じサイトに戻って来ることは可能である。
この考えを更に進めて、サイトを各集合に分類し、それをセルと呼ぶことにすると、セル内の多重散乱とセル間の散乱に分けることができる。セル内の多重散乱を繰り込んだものをセルT行列と定義すると、
ここで、セルに含まれてるサイトは、和を取ると全サイトの和に等しくなるように定義している。
注意として、セル間の散乱といえども、伝搬はあくまでサイト間である。「サイト同士が必ず異なるセルに属するような散乱」を「セル間の散乱」とここでは呼んでいる。
まとめると、
以下、もう少し細かく説明する。
サイトの和から、あるサイトだけを抜き出してを作ることを考える。
したがって、三次までで整理すると、
級数の中でがやがてになっていくことが何となく感じられたと思う。
ポイントは、から移るときには必ず異なるサイトへ行くようになっていることである。
他のサイトに関しても同様にを作るように抜き出していくと、最終的にサイトT行列展開が得られる。
サイトT行列展開の時と同様に、サイトの和からあるセルに含まれるサイト集合を抜き出してを作ることを考える。
異なるセルに属していれば、サイトが被ることは無い。
したがって、三次までで整理すると、
級数の中でがやがてになっていくことが何となく感じられたと思う。
セルのみを抜き出したため、の和はとなる場合があることに注意。そのために同じサイトへの伝搬が無い様に和が複雑になっている。
他のセルに関しても同様にを作るように抜き出していくと、最終的にセルT行列展開が得られる。