第一原理計算によるフォノン計算についてまとめる。
計算に使用するソフトウェアはQuantum Espressoである。

フォノン計算方法

フォノン計算の流れ

フォノン計算は次のように行う。[1]
1. scf計算 (pw.x)
2. 粗いq点メッシュ上のdynamical matrixを計算(ph.x)
3. 2.で得られた結果を実空間にフーリエ変換し、force constantを計算(q2r.x)
4. 3.で得られた結果を任意の波数にフーリエ変換(matdyn.x)
3., 4.の過程ではacoustic sum ruleを用いることで、計算誤差によるimaginary phonon modeが出ないようにすることができる。

Siでの計算例は[1A]に、Acoustic Sum Ruleに関しては[1B]にある。

scf計算

フォノン計算する際に、誘電率も計算したい場合、occupationsでsmearingを設定しないこと。
半導体であってもsmearingに設定していると、金属と判断されてしまう。

フォノン計算

フォノン計算の入力ファイルを用意する。
バンド分散を計算する場合は、ldisp=.true.とし、q点のメッシュを指定する。

誘電率を計算したい場合、epsil=.true.とするが、Γ点でしか計算できないため、ldisp=.flase.とし、Γ点(0 0 0)を指定する。また、上述したように金属では計算できないため注意が必要である。

実空間へのフーリエ変換

ph.xで計算したdynamical matrixをq2r.xを使って実空間にフーリエ変換する。

波数空間へのフーリエ変換

任意の波数のフォノンの情報をmatdyn.xで計算する。

フォノンの状態密度(DOS: Density Of States)を計算する場合は、dos=.true.とする。

フォノンのバンド分散を計算する際には、dos=.false.とし、q_in_band_form=.true.とする。

計算するq点の情報を入力する必要があるが、逆格子ベクトルで指定しても良いし、Γ=gG、X, Aのように記号で指定することもできる。(個人的には記号で指定するほうが間違いがなくて良い。)

誘電率の計算

GaAsのフォノン計算の例に誘電率についての記載がある。[2]
ph.x のinputファイルにおいてepsil = .true.とすればoutputファイルに、Dielectric constant in cartesian axisが出力される。
※GaAsのような分極をもつ材料の場合、LOフォノンとTOフォノンがΓ点で分裂する。
[2]ではそのことに言及されていないが、正しいフォノン分散を得るためには必要不可欠である。

さて、誘電率は格子系誘電率と電子系誘電率に分けられる。
周波数0の誘電率のみで良ければ、電子系誘電率はepsilon.xから、格子系誘電率はph.xから求めることができる。[3]

Siのepsilon.xの例は[4]にある。Si, GaAsの例は[5]にある。