概要

GaAsを例に、ガンマ点におけるフォノンの縦光学(LO)フォノンと横光学フォノン(TO)の分裂を扱います
極性半導体ではLO-TO分裂が起こるので、この計算は必須です。

計算手順

準備

GaAsのCIFファイルは以下のサイトのものを使い、格子定数はvc-relax計算で得られたものに変更しました

https://materials.springer.com/isp/crystallographic/docs/sd_0311662
vc-relax計算で得られた格子定数は 5.663 Å = 10.702 Bohr

擬ポテンシャルはPS Libraryからウルトラソフト型を利用しました。

擬ポテンシャル … https://www.quantum-espresso.org/pseudopotentials/ps-library/
- As.pbesol-n-rrkjus_psl.1.0.0.UPF
- Ga.pbesol-dnl-rrkjus_psl.1.0.0.UPF

第一原理計算にはQEを使用します。入力ファイルは以下のように作成します。

$ cif2cell -p pwscf -f GaAs.cif -o GaAs.scf.in

ボルン有効電荷の計算

DFPTを用いてボルン有効電荷を計算します。 SCF計算で基底状態を求めたあと、ph.xを使って計算します

GaAs.scf.in

&control
  prefix      = 'GaAs',
  calculation = 'scf',
  pseudo_dir  = './pseudo/',
  outdir      = './tmp/',
/

&SYSTEM
  ibrav = 0
  nat   = 2
  ntyp  = 2
  ecutwfc =  40.0
  ecutrho = 240.0
/

&electrons
  conv_thr = 1.0d-10
/

CELL_PARAMETERS {angstrom}
   2.8314   2.8314   0.0000
   2.8314   0.0000   2.8314
   0.0000   2.8314   2.8314

ATOMIC_SPECIES
  As   74.92100  As.pbesol-n-rrkjus_psl.1.0.0.UPF
  Ga   69.72300  Ga.pbesol-dnl-rrkjus_psl.1.0.0.UPF

ATOMIC_POSITIONS {crystal}
As            0.2500000000        0.2500000000        0.2500000000
Ga            0.0000000000        0.0000000000        0.0000000000

K_POINTS automatic
  8 8 8 1 1 1

GaAs.phG.in

phonons of GaAs
&inputph
  tr2_ph = 1.0d-14
  prefix = 'GaAs'
  outdir = './tmp/'
  fildyn = 'GaAs.dynG'
  epsil	 = .true.
  zeu    = .true.
  trans  = .false.
/
0.0  0.0  0.0

epsilは誘電率を、zeuは原子に加わる力からボルン有効電荷を計算するように指定します。このほか、分極からボルン有効電荷を計算するzueもありますが、どちらを使ってもほぼ同じ値になります。

また、フォノンの振動数の計算は必要ないので、transは.false.にしておきます。

$ pw.x < GaAs.scf.in > GaAs.scf.out
$ ph.x < GaAs.phG.in > GaAs.phG.out

誘電率およびボルン有効電荷は次のように出力されます。

GaAs.phG.out

（略）
          Dielectric constant in cartesian axis 

          (      14.192331204       0.000000000      -0.000000000 )
          (       0.000000000      14.192331204      -0.000000000 )
          (      -0.000000000      -0.000000000      14.192331204 )

          Effective charges (d Force / dE) in cartesian axis

           atom      1   As 
      Ex  (       -2.13215       -0.00000        0.00000 )
      Ey  (       -0.00000       -2.13215       -0.00000 )
      Ez  (       -0.00000       -0.00000       -2.13215 )
           atom      2   Ga 
      Ex  (        2.14310        0.00000       -0.00000 )
      Ey  (       -0.00000        2.14310        0.00000 )
      Ez  (        0.00000        0.00000        2.14310 )
（略）

Dynamical matrixの計算

QE用の超格子の入力ファイルを作成します

$ cif2cell --supercell=[3,3,3] -p pwscf -f GaAs.cif -o GaAs_333.pw.in

超格子のインプットファイルに足りない部分を書き加えます。 cif2cellのデフォルトで格子定数はA（オングストローム単位）で指定されますが、celldm(1)（ボーア半径単位）に直しておきます。また、原子に加わる力の大きさを出力するために tprnfor=.true. を忘れないようにします。

GaAs_333.harmonic.in

&control
  prefix      = 'GaAs',
  calculation = 'scf',
  pseudo_dir  = './pseudo/',
  outdir      = './tmp/',
  tprnfor=.true.
/

&SYSTEM
  ibrav = 0
  celldm(1) = 10.702
  nat = 54
  ntyp = 2
/

&electrons
  conv_thr = 1.0d-8
/

CELL_PARAMETERS {alat}
  1.500000000000000   1.500000000000000   0.000000000000000 
  1.500000000000000   0.000000000000000   1.500000000000000 
  0.000000000000000   1.500000000000000   1.500000000000000 

ATOMIC_SPECIES
  Ga   69.72300  Ga.pbesol-dnl-rrkjus_psl.1.0.0.UPF
  As   74.92100  As.pbesol-n-rrkjus_psl.1.0.0.UPF

K_POINTS automatic
  4 4 4 1 1 1

ATOMIC_POSITIONS {crystal}
As   0.083333333333333   0.083333333333333   0.083333333333333 
As   0.416666666666667   0.083333333333333   0.083333333333333 
（略）
Ga   0.666666666666667   0.333333333333333   0.333333333333333

また、ALAMODE用の原子間力定数の計算の入力ファイルを作成します

GaAs_333.pw.in

&general
  PREFIX = GaAs_333
  MODE = suggest
  NAT = 54; NKD = 2
  KD = As Ga
/

&interaction
  NORDER = 1
/

&cell
  10.702
  1.500000000000000   1.500000000000000   0.000000000000000 
  1.500000000000000   0.000000000000000   1.500000000000000 
  0.000000000000000   1.500000000000000   1.500000000000000 
/

&cutoff 
  *-* None
/

&position
1   0.083333333333333   0.083333333333333   0.083333333333333 
1   0.416666666666667   0.083333333333333   0.083333333333333 
（略）
2   0.666666666666667   0.333333333333333   0.333333333333333 
/

できたらALMを実行します。

$ alm GaAs_333.harmonic.in > GaAs_333.harmonic.log

GaAs_333.pattern_HARMONIC というファイルができるので、ALAMODEに付属のPythonスクリプトを使ってQEの入力ファイルを作成します。

$ python -m displace --QE GaAs_333.pw.in --mag 0.01 --prefix disp -pf GaAs_333.pattern_HARMONIC

今回の場合はQEの入力ファイルが2つできます。出てきた入力ファイルをすべてQEで実行します

$ pw.x < disp1.pw.in > disp1.pw.out
$ pw.x < disp2.pw.in > disp2.pw.out

計算は重いのでOpenMPやMPIによる並列計算をおすすめします

これをALAMODEに付属のPythonスクリプトで整理します。

$ python -m extract --QE GaAs_333.pw.in *.pw.out > DFSET_harmonic

ALMを使って原子間力定数を計算します。 GaAs_333.harminic.inをコピーして GaAs.harmonic_opt.in を生成し、次の点を変更します

&generalフィールドで MODE=optimize に変更
&optimizeフィールドを作成して、中に DFSET = DFSET_harmonic と記述

Si_333.harmonic_opt.in

&general
  PREFIX = GaAs_333
  MODE = optimize
  NAT = 54; NKD = 2
  KD = As Ga
/

&interaction
  NORDER = 1
/

&cell
  10.702
  1.500000000000000   1.500000000000000   0.000000000000000 
  1.500000000000000   0.000000000000000   1.500000000000000 
  0.000000000000000   1.500000000000000   1.500000000000000 
/

&cutoff 
  *-* None
/

&optimize
  DFSET = DFSET_harmonic
/

&position
1   0.083333333333333   0.083333333333333   0.083333333333333 
1   0.416666666666667   0.083333333333333   0.083333333333333 
（略）
2   0.666666666666667   0.333333333333333   0.333333333333333  
/

ALMを実行します

$ alm GaAs_333.harmonic_opt.in > GaAs.harmonic_opt.log

GaAs_333.fcs と GaAs_333.xml が出力されたら成功です。

フォノンの分散

DFPT計算で得た GaAs.phG.out から誘電率とボルン有効電荷を記録したファイルを作ります

GaAs.born

14.200395680      -0.000000000      -0.000000000
-0.000000000      14.200395680      -0.000000000
-0.000000000      -0.000000000      14.200395680
-2.12591       -0.00000        0.00000
-0.00000       -2.12591       -0.00000
 0.00000        0.00000       -2.12591
 2.14222        0.00000        0.00000
 0.00000        2.14222        0.00000
 0.00000        0.00000        2.14222

誘電率、各原子のボルン有効電荷の順に書きます

各原子の有効ボルン電荷の順番は、以下のanphonを実行することで得られる GaAs.phband.log 内のものに合わせてください

GaAs.phband.log

（略）
  Atomic positions in the primitive cell (fractional):
   1:   2.500000e-01   2.500000e-01   2.500000e-01  As
   2:   0.000000e+00   0.000000e+00   0.000000e+00  Ga
（略）

ボルン有効電荷なし

GaAs.phband.in

&general
  PREFIX = GaAs_333
  MODE   = phonons
  FCSXML = GaAs_333.xml

  NKD = 2
  KD = As Ga
/

&cell
  10.702
  0.500000000000000   0.500000000000000   0.000000000000000 
  0.500000000000000   0.000000000000000   0.500000000000000 
  0.000000000000000   0.500000000000000   0.500000000000000 
/

&kpoint
  1
  G 1.00 1.00 1.00  X 0.50 0.50 1.00 51
  X 0.50 0.50 1.00  G 0.00 0.00 0.00 51
  G 0.00 0.00 0.00  L 0.50 0.50 0.50 51
  L 0.50 0.50 0.50  X 0.50 0.50 0.00 51
  X 0.50 0.50 0.00  W 0.75 0.50 0.25 51
  W 0.75 0.50 0.25  L 0.50 0.50 0.50 51
/

ボルン有効電荷あり

GaAs.phband.in

&general
  PREFIX = GaAs_333
  MODE   = phonons
  FCSXML = GaAs_333.xml
  NONANALYTIC = 3
  BORNINFO = GaAs.born

  NKD = 2
  KD = As Ga
/

&cell
  10.702
  0.500000000000000   0.500000000000000   0.000000000000000 
  0.500000000000000   0.000000000000000   0.500000000000000 
  0.000000000000000   0.500000000000000   0.500000000000000 
/

&kpoint
  1
  G 1.00 1.00 1.00  X 0.50 0.50 1.00 51
  X 0.50 0.50 1.00  G 0.00 0.00 0.00 51
  G 0.00 0.00 0.00  L 0.50 0.50 0.50 51
  L 0.50 0.50 0.50  X 0.50 0.50 0.00 51
  X 0.50 0.50 0.00  W 0.75 0.50 0.25 51
  W 0.75 0.50 0.25  L 0.50 0.50 0.50 51
/

$ anphon GaAs.phband.in > GaAs.phband.log

赤い実線がボルン有効電荷あり、青い破線がなしです。ガンマ点でTOフォノンとLOフォノンが分裂するのがわかります。

Table of Contents

概要

計算手順

準備

ボルン有効電荷の計算

Dynamical matrixの計算

フォノンの分散