Phonons & Phonopy: “Pro Tips”

Phonons & Phonopy: “Pro Tips” [Version 2.0] J. M. Skelton Bath CompChem 14th April 2015 Phonons and La*ce Dynamics Crystallography is generally concerned with the sta3c proper3es of crystals, describing features such as the average posi3ons of atoms and the symmetry of a crystal. Solid state physics takes a similar line as far as elementary electronic proper3es are concerned. We know, however, that atoms actually move around inside the crystal structure, since it is these mo3ons that give the concept of temperature […]. The sta3c laBce model, which is only concerned with the average posi3ons of atoms and neglects their mo3ons, can explain a large number of material features […]. There are, however, a number of proper3es that cannot be explained by a sta3c model… MarEn Dove, “IntroducEon to LaNce Dynamics” Bath CompChem, April 2015 | Slide 2 What Can Phonopy Do? “Anharmonicity” Proper:es phonopy
phonopy-gruneisen
phonopy-qha
phono3py
Bath CompChem, April 2015 | Slide 3 Overview •  “PracEcal” Theory •  phonopy
§ 
§ 
§ 
§ 
Workflow/setup CalculaEng forces: opEons and “gotchas” Post processing Input/output files; using phonopy with other force-­‐constant calculators •  phonopy-qha
§  Workflow and post-­‐processing §  Example applicaEons: which DFT funcEonal? •  phono3py
§  Workflow, setup and post processing §  Examples •  Some current developments •  Summary Bath CompChem, April 2015 | Slide 4 “Prac:cal” Theory Force-­‐constant matrix: From finite differences: Dynamical matrix: A_er diagonalisaEon: Bath CompChem, April 2015 | Slide 5 “Prac:cal” Theory Bath CompChem, April 2015 | Slide 6 phonopy: Workflow Input Structure Create Displacements phonopy -d [--dim=“1 1 1”]
Calculate Forces Extract Forces Post-­‐Process phonopy -f vasprun-{001..XXX}.xml
phonopy --fc vasprun.xml
phonopy [-t] [-p] [-s] Settings.conf
Bath CompChem, April 2015 | Slide 7 phonopy: Setup phonopy -d
[--dim=“1 1 1”]
[--tolerance=1e-5]
[-c POSCAR]
Set up the calculaEons on an XxYxZ supercell PosiEon tolerance for symmetry detecEon Path to POSCAR-­‐format structure (“cell”), if not “POSCAR” disp.yaml
POSCAR-001
POSCAR-002
POSCAR-003
...
SPOSCAR
InformaEon about the structure, supercell and displaced atoms POSCAR files containing single atomic displacements “Unperturbed” supercell for DFPT phonon calculaEon Bath CompChem, April 2015 | Slide 8 phonopy: Calcula:ng Forces Sample finite-­‐displacement INCAR: ADDGRID = .TRUE.
EDIFF = 1E-8
ENCUT = 500-800 eV
LREAL = .FALSE.
PREC = High | Accurate
Sample VASP force-­‐constants INCAR: ADDGRID = .TRUE.
EDIFF = 1E-8
ENCUT = 500-800 eV
IBRION = 5|6|7|8
LREAL = .FALSE.
NSW = 1
PREC = High | Accurate
•  Accurate forces are essenEal -­‐> crank the standard seNngs right up •  LREAL = .FALSE. is essenEal, unless you manually adjust ROPT
•  ADDGRID = .TRUE. doesn’t seem to be essenEal, but doesn’t cost much either •  For finite-­‐difference/DFPT phonon calculaEons in VASP, set NSW = 1
Bath CompChem, April 2015 | Slide 9 Bath CompChem, April 2015 | Slide 10 ADDGRID = .TRUE.
LREAL = .FALSE.
ADDGRID = .TRUE.
LREAL = Auto
ADDGRID = .FALSE. ADDGRID = .FALSE.
LREAL = .FALSE.
LREAL = Auto
phonopy: Calcula:ng Forces (PbTe) phonopy: Calcula:ng Forces (ZnS) Bath CompChem, April 2015 | Slide 11 phonopy: Calcula:ng Forces (ZnS) •  Symmetrising the force constants during the post processing with phonopy corrects the PBE dispersion •  SymmetrisaEon enforces the acousEc-­‐sum rule in the TPSS dispersion, but does not remove the artefacts away from Γ •  Symmetrising within phonopy gives the same result as calcula3ng the force constants with the VASP internal DFPT/finite-­‐
displacement rou3nes -­‐ as with most codes, VASP enforces symmetry by default Bath CompChem, April 2015 | Slide 12 phonopy: Calcula:ng Forces (ZnS) Bath CompChem, April 2015 | Slide 13 phonopy: Calcula:ng Forces Bath CompChem, April 2015 | Slide 14 phonopy: Calcula:ng Forces Bath CompChem, April 2015 | Slide 15 phonopy: A Few “Top Tips” If using VASP FD/DFPT, set NWRITE = 3 in the INCAR file, and you can run this bash script on the OUTCAR to obtain a simulated IR spectrum “for free”: hrp://homepage.univie.ac.at/david.karhanek/downloads.html#Entry02 If using DFPT with an LDA/GGA funcEonal, set LEPSILON = .TRUE. in the INCAR file to obtain the staEc dielectric constant, in parEcular the ionic-­‐
relaxaEon part, for a small added cost When using FD/DFPT, VASP tries to change the k-­‐point set internally, which requires NPAR = #Cores to be set in the INCAR file; seNng ISYM = -1 avoids this, and although the number of displacements which need to be evaluated may increase, the performance gained by using band parallelism can quite easily offset this for low-­‐symmetry systems (!) Bath CompChem, April 2015 | Slide 16 phonopy: Post Processing phonopy -p -s Settings.conf
“Plot” SeNngs file “Save” PbTe Sample phonon DOS se*ngs file: DIM = 4 4 4
MP = 48 48 48
GAMMA_CENTER = .TRUE.
Bath CompChem, April 2015 | Slide 17 phonopy: Post Processing phonopy -p -s Settings.conf
“Plot” SeNngs file “Save” PbTe Sample phonon DOS se*ngs file: DIM = 4 4 4
MP = 48 48 48
GAMMA_CENTER = .TRUE.
EIGENVECTORS = .TRUE.
PDOS = 1, 2
Bath CompChem, April 2015 | Slide 18 phonopy: Post Processing SeNngs file phonopy -p -s -t Settings.conf
“Plot” “Save” “[Calculate] thermal properEes” PbTe Sample phonon DOS se*ngs file: DIM = 4 4 4
MP = 48 48 48
GAMMA_CENTER = .TRUE.
Bath CompChem, April 2015 | Slide 19 phonopy: Post Processing phonopy -p -s Settings.conf
“Plot” SeNngs file “Save” PbTe Sample phonon band structure se*ngs file: DIM = 4 4 4
BAND = 0.0 0.0 0.0 0.5 0.25 0.75
0.5 0.0 0.5 0.0 0.0 0.0
0.5 0.5 0.5
BAND_POINTS = 101
BAND_LABELS = \Gamma W X \Gamma L
[EIGENVECTORS = .TRUE.]
Bath CompChem, April 2015 | Slide 20 phonopy: Post Processing phonopy -p -s Settings.conf
“Plot” SeNngs file “Save” PbTe Sample phonon band structure se*ngs file: DIM = 4 4 4
BAND = 0.0 0.0 0.0 0.5 0.25 0.75
0.5 0.0 0.5 0.0 0.0 0.0
0.5 0.5 0.5
BAND_POINTS = 101
BAND_LABELS = \Gamma W X \Gamma L
BAND_CONNECTION = .TRUE.
Bath CompChem, April 2015 | Slide 21 phonopy: Non-­‐Analy:cal Correc:ons •  To apply a non-­‐analyEcal correcEon (LO/TO spliNng) to the phonon frequencies, phonopy needs the Born effecEve charges and electronic-­‐polarisaEon contribuEon to the macroscopic dielectric constant •  In VASP, for LDA/GGA funcEonals these can be computed using DFPT; for others, they need to be computed from the response to an electric field INCAR for Born charges using DFPT: EDIFF = 1E-8
ENCUT = 500-800 eV
LEPSILON = .TRUE.
LREAL = .FALSE.
NSW = 0
PREC = High | Accurate
Bath CompChem, April 2015 | Slide 22 INCAR for Born charges using LCALCEPS: EDIFF = 1E-8
ENCUT = 500-800 eV
LCALCEPS = .TRUE.
LREAL = .FALSE.
NSW = 0
PREC = High | Accurate
[EFIELD_PEAD = Ex Ey Ez]
phonopy: Non-­‐Analy:cal Correc:ons outcar-born > BORN
Sample BORN file: <Conversion Factor>
εxx εxy εxz εyx εyy
Zxx Zxy Zxz Zyx Zyy
Zxx Zxy Zxz Zyx Zyy
Dielectric tensor εyz
Zyz
Zyz
εzx
Zzx
Zzx
εzy
Zzy
Zzy
εzz
Zzz
Zzz
Born charge tensors for unique atoms •  CorrecEons are enabled by seNng NAC = .TRUE. in the configuraEon file, or passing --nac as a command-­‐line argument •  When this opEon is used, phonopy expects to find a BORN file in the working directory Bath CompChem, April 2015 | Slide 23 phonopy: Non-­‐Analy:cal Correc:ons NAC = .FALSE.
PbTe Bath CompChem, April 2015 | Slide 24 NAC = .TRUE.
phonopy: Force-­‐Constant Symmetrisa:on FC_SYMMETRY = 0
FC_SYMMETRY = 1
ZnS •  Force-­‐constant symmetrisaEon is enabled by seNng FC_SYMMETRY = > 0 in the configuraEon file Bath CompChem, April 2015 | Slide 25 phonopy: Output Files Sample mesh.yaml file: Sample band.yaml file: mesh: [ mx, my, mz ]
nqpoint: 32000
natom: 8
phonon:
- q-position: [ qx, qy, qz ]
weight: w1
band:
- # 1
frequency: ω1
...
nqpoint: 808
npath: 8
natom: 8
phonon:
- q-position: [ qx, qy, qz ]
distance: d1
band:
- # 1
frequency: ω1
...
•  If EIGENVECTORS = .TRUE. is set in the configuraEon file, the mode eigenvectors will also appear in these files •  With BAND_CONNECTION = .TRUE., the frequencies for each band in band.yaml are ordered so that they connect across the band structure Bath CompChem, April 2015 | Slide 26 phonopy: Output Files Sample total_dos.dat file: # Sigma = 0.053821
-0.5372... 0.0000...
-0.5103... 0.0000...
-0.4834... 0.0000...
-0.4564... 0.0000...
-0.4295... 0.0000...
-0.4026... 0.0000...
-0.3757... 0.0000...
-0.3488... 0.0000...
-0.3219... 0.0000...
...
•  The “parEal_dos.dat” file generated with EIGENVECTORS = .TRUE. contains one column for each atom in the primiEve cell Bath CompChem, April 2015 | Slide 27 Sample thermal_proper:es.yaml file: unit:
temperature: K
...
natom:
8
zero_point_energy: 18.9108676
high_T_entropy:
847.3220815
thermal_properties:
- temperature:
0.0000000
- free_energy:
18.9108676
- entropy:
0.0000000
- heat_capacity: 0.0000000
- energy:
18.9108676
...
phonopy: Force-­‐Constant Calculators
Sample FORCE_SETS file: 128
2
1
d1x
F1x
F2x
...
2
d2x
F1x
F2x
...
d1y
F1y
F2y
d1z
F1z
F2z
d2y
F1y
F2y
d2z
F1z
F2z
Bath CompChem, April 2015 | Slide 28 Sample FORCE_CONSTANTS file: 128
1 1
Φxx
Φyx
Φzx
1 2
Φxx
Φyx
Φzx
...
Φxy Φxz
Φyy Φyz
Φzy Φzz
Φxy Φxz
Φyy Φyz
Φzy Φzz
phonopy: Force-­‐Constant Calculators
IRMOF-­‐10 *Tinker calculaEons by J. K. Bristow and D. Tiana Bath CompChem, April 2015 | Slide 29 phonopy-qha: Workflow EoS Curve E/V curve -­‐> “e-­‐v.dat” Harmonic Phonopy phonopy -d --dim=“...”
Thermal ProperEes phonopy -f vasprun-{001..XXX}.xml
phonopy -t -p -s Settings.conf
Post-­‐Process phonopy-qha
e-v.dat
thermal_properties-{001..XXX}.yaml
--tmax=980
Needs to be two points smaller [--tstep=10]
than the maximum temperature [--pressure=<p/GPa>]
in the YAML files [-p] [-s]
Bath CompChem, April 2015 | Slide 30 phonopy-qha: Output *J. M. Skelton et al., Phys. Rev. B 89, 205203 (2014) Bath CompChem, April 2015 | Slide 31 phonopy-qha: Output bulk_modulus-­‐temperature.dat <-­‐ B is temperature dependent(!) Cp-­‐temperature.dat Cp-­‐temperature_polyfit.dat Cv-­‐volume.dat <-­‐ CV at each volume, at each temperature dsdv-­‐temperature.dat entropy-­‐volume.dat <-­‐ SV at each volume, at each temperature gibbs-­‐temperature.dat gruneisen-­‐temperature.dat <-­‐ Average Gruneisen parameter (?) helmholtz-­‐volume.dat <-­‐ A at each volume, at each temperature thermal_expansion.dat volume_expansion.dat volume-­‐temperature.dat Bath CompChem, April 2015 | Slide 32 phonopy-qha: Examples *J. M. Skelton et al., in preparaEon Bath CompChem, April 2015 | Slide 33 phonopy-qha: Examples *J. M. Skelton et al., in preparaEon Bath CompChem, April 2015 | Slide 34 phonopy-qha: Examples *J. M. Skelton et al., in preparaEon Bath CompChem, April 2015 | Slide 35 phonopy-qha: Examples ~4 K 300-­‐310 K ~105 K 450-­‐475 K ~150 K 550-­‐600 K *J. M. Skelton et al., Phys. Rev. B 89, 205203 (2014) Bath CompChem, April 2015 | Slide 36 phonopy-qha: Examples SnS: 10 K -­‐ 350 K PbTe: 0 K -­‐ 600 K *PbTe: J. M. Skelton et al., Phys. Rev. B 89, 205203 (2014) Bath CompChem, April 2015 | Slide 37 phono3py: Workflow Input Structure Create Displacements phono3py -d --dim=“2 2 2”
Calculate Forces Extract Forces phono3py --cf3 vasprun-{001..XXX}.xml
phono3py --dim=“2 2 2”
Post-­‐Process phono3py --tmax=1000 --tstep=10 --fc2 --fc3
--dim=“2 2 2” -v --mesh=“24 24 24”
--br --sigma=“0.1”
Bath CompChem, April 2015 | Slide 38 phono3py: Setup phono3py -d --dim=“...”
[--cutoff_pair=4]
[--dim_fc2=“4 4 4”]
phono3py --cf3 vasprun-{001..XXX}.xml
Extract forces from VASP output phono3py --dim=“...”
[--dim_fc2...]
The --cutoff_pair tag uses the same numbering for the displaced POSCAR files as the full calculaEon; this means the cutoff can be increased, and the extra displacements added, systemaEcally, to converge w.r.t. the interacEon range Bath CompChem, April 2015 | Slide 39 phono3py: Post Processing Read in pre-­‐calculated force constants phono3py --tmax=1000 --tstep=10 --fc2 --fc3
--dim=“...” -v --mesh=“24 24 24”
--br --sigma=“0.1” [--nac]
[--dim2=“...”/--dim_fc2=“...”]
# of interacEons per q-­‐point becomes larger with mesh size; cannot easily “max out” as for DOS calculaEons, but needs to be converged •  The post-­‐processing (mainly the phonon-­‐lifeEme calculaEons) takes a very long Eme for large supercells/large or low-­‐symmetry structures •  It is possible to run the calculaEon on (ranges of) q-­‐points separately, and then combine them a_erwards •  Various post-­‐processing tags can be applied, e.g. to incorporate isotope effects Bath CompChem, April 2015 | Slide 40 phono3py: Examples *J. M. Skelton et al., APL Materials 3, 041102 (2015) Bath CompChem, April 2015 | Slide 41 phono3py: Examples Simulated IR/Raman spectra with… • 10 K linewidths • 300 K linewidths … and an instrument broadening of 3.5 cm-­‐1, compared to expt. *J. M. Skelton et al., APL Materials 3, 041102 (2015) Bath CompChem, April 2015 | Slide 42 Some Current Developments *E. L. da Silva et al., Phys. Rev. B 91, 144107 (2015) Bath CompChem, April 2015 | Slide 43 A Few Closing Remarks Bath CompChem, April 2015 | Slide 44