PHononは、フォノンに関する計算を行うパッケージである。 フォノンの分散や状態密度などを計算できる。

フォノンの計算では、原子が平衡位置、すなわち力の釣り合いの位置にないとおかしな結果が出てしまうので注意する。 事前に構造緩和計算が必要である。

Siを例に、フォノンのもっとも簡単な計算を実行してみる。 事前にpw.xによる電子状態計算が必要である。

si.scf.in

&control
   calculation='scf',
   prefix='si'
   pseudo_dir = './',
   outdir='./'
/
&system
   ibrav=2, celldm(1)=10.2, nat=2, ntyp=1,
   ecutwfc=18.0
/
&electrons
/
ATOMIC_SPECIES
 Si  28.086  Si.pz-vbc.UPF
ATOMIC_POSITIONS (alat)
 Si 0.00 0.00 0.00
 Si 0.25 0.25 0.25
K_POINTS automatic
 8 8 8 1 1 1

$ pw.x < si.scf.in > si.scf.out

電子状態の計算が終われば、フォノンの振動数を計算できるようになる。 まず、k=(0,0,0)のΓ点におけるフォノンを計算してみる。 以下のファイルを用意する。

si.phG.in

phonons of Si at Gamma
 &inputph
  tr2_ph=1.0d-14,
  prefix='si',
  epsil=.true.,
  amass(1)=28.08,
  outdir='./',
  fildyn='si.dynG',
 /
0.0 0.0 0.0

第1行はコメント行とみなされ無視される。 また、最後の行は計算する波数（デカルト座標で、単位は格子定数aを使って2π/a）を表す。 各引数の意味は次のとおりである。 詳細はマニュアルを参照のこと。

実行は次のようにする。

$ ph.x < si.phG.in > si.phG.out

結果はsi.dynGで見ることができる。

si.dynG

（中略）
 **************************************************************************
     freq (    1) =      -0.170404 [THz] =      -5.684069 [cm-1]
 (  0.274931  0.000000 -0.504256  0.000000 -0.412479  0.000000 ) 
 (  0.274931  0.000000 -0.504256  0.000000 -0.412479  0.000000 ) 
     freq (    2) =      -0.170404 [THz] =      -5.684069 [cm-1]
 ( -0.161796  0.000000  0.380828  0.000000 -0.573404  0.000000 ) 
 ( -0.161796  0.000000  0.380828  0.000000 -0.573404  0.000000 ) 
     freq (    3) =      -0.170404 [THz] =      -5.684069 [cm-1]
 (  0.631059  0.000000  0.317327  0.000000  0.032690  0.000000 ) 
 (  0.631059  0.000000  0.317327  0.000000  0.032690  0.000000 ) 
     freq (    4) =      15.296893 [THz] =     510.249417 [cm-1]
 ( -0.463984  0.000000  0.270574  0.000000  0.459900  0.000000 ) 
 (  0.463984  0.000000 -0.270574  0.000000 -0.459900  0.000000 ) 
     freq (    5) =      15.296893 [THz] =     510.249417 [cm-1]
 (  0.527053  0.000000  0.327503  0.000000  0.339053  0.000000 ) 
 ( -0.527053  0.000000 -0.327503  0.000000 -0.339053  0.000000 ) 
     freq (    6) =      15.296893 [THz] =     510.249417 [cm-1]
 ( -0.083269  0.000000  0.565271  0.000000 -0.416575  0.000000 ) 
 (  0.083269  0.000000 -0.565271  0.000000  0.416575  0.000000 ) 
 **************************************************************************

• 6つの振動モードがあって、それぞれのΓ点における振動数がわかる。
• 音響モード（-0.170404 THz, あるいは 0 Thz）と光学モード（15.296893 THz）がそれぞれ三重に縮退している。
• 原子の変位ベクトルは、例えばfreq 4であれば1番目のSiが(-0.463984, 0.270574, 0.459900)方向、2番目のSiが(0.463984, -0.270574, -0.459900)である（デカルト座標）。

音響モードの振動数を0にしたければ構造緩和を行うとよい。

今度はX点(2π/a,0,0)のフォノンを調べる。 手順はさきほどと同じであるが、誘電率が計算できないことに注意する。

si.phX.in

phonons of Si at X
 &inputph
  tr2_ph=1.0d-14,
  prefix='si',
  amass(1)=28.08,
  outdir='./',
  fildyn='si.dynX',
 /
1.0 0.0 0.0

$ ph.x < si.phX.in > si.phX.out

結果はsi.dynXに出力される。

（中略）
 **************************************************************************
     freq (    1) =       4.207901 [THz] =     140.360482 [cm-1]
 ( -0.000000  0.000000 -0.023972  0.000000  0.706700  0.000000 ) 
 ( -0.000000  0.000000  0.706700  0.000000 -0.023972  0.000000 ) 
     freq (    2) =       4.207901 [THz] =     140.360482 [cm-1]
 ( -0.000000  0.000000 -0.706700  0.000000 -0.023972  0.000000 ) 
 ( -0.000000  0.000000 -0.023972  0.000000 -0.706700  0.000000 ) 
     freq (    3) =      12.236850 [THz] =     408.177391 [cm-1]
 ( -0.999882  0.000000  0.000000  0.000000  0.000000  0.000000 ) 
 ( -0.015351  0.000000 -0.000000  0.000000  0.000000  0.000000 ) 
     freq (    4) =      12.236850 [THz] =     408.177391 [cm-1]
 ( -0.015351  0.000000  0.000000  0.000000  0.000000  0.000000 ) 
 (  0.999882  0.000000 -0.000000  0.000000  0.000000  0.000000 ) 
     freq (    5) =      13.746712 [THz] =     458.540957 [cm-1]
 (  0.000000  0.000000 -0.000000  0.000000  0.707107  0.000000 ) 
 ( -0.000000  0.000000 -0.707107  0.000000  0.000000  0.000000 ) 
     freq (    6) =      13.746712 [THz] =     458.540957 [cm-1]
 (  0.000000  0.000000  0.707107  0.000000  0.000000  0.000000 ) 
 (  0.000000  0.000000  0.000000  0.000000 -0.707107  0.000000 ) 
 **************************************************************************

2重に縮退したものが3つあることがわかる。