- Materials Square

## [MatSQ Tip] Performing Phonon Calculation

2020-05-07 18:01:55

1. Definition of Phonon

Phonon is defined as quanta of crystal vibrational field and it is considered quasiparticle.

In the solid, various interatomic forces continuously change interact with each other. When the movement has periodicity, it called 'Lattice vibration'. This lattice vibration can describe as the phonon dispersion relation curve, representing the phonon frequency with respect to the wave vector.

The phonon dispersion consists of acoustic phonon mode and optical phonon mode. The acoustic phonon is a coherent movement of vibration (in-phase), which means one atom and neighbour atom have vibrate some direction, and the optical phonon is out-of-phase vibration, vice versa. Figure 1 shows (a) phonon eigenvector and (b) phonon dispersion of graphene. At G point, there are 6 phonon modes existed with 3 acoustic modes (blue solid line) and 3 optical modes (red solid line). Also, you may find the phonon eigenvectors of acoustic mode (Figure 1. (a) upper) and optical mode (Figure 1. (a) lower).

The number of phonon modes can calculate simply. A primitive cell has N of atoms, the degree of freedom is 3N with considering x, y, z directions. Among them, acoustic phonon mode is always 3 and 3N-3 is optical phonon mode. For example, the phonon dispersion of silicon described with 6 of phonon modes, because silicon having 2 atoms in primitive cell.

The acoustic phonon has a relatively low frequency, and move to the same phase in the unit cell. It has one longitudinal acoustic mode (LA), having large phonon group velocity, and two transverse acoustic modes (TA). Therefore, acoustic phonons have a great influence on thermal conductivity characteristic which is important to have a large phonon group velocity and a long phonon relaxation time.

The optical phonon has high frequency than acoustic phonon and has opposite vibration direction of one atom and neighbour atom. Thus, it has small phonon group velocity and broad energy range than acoustic phonon modes. Especially, it interacts with photon. The out-of-plane vibration cause the dipole moment. Therefore, IR active and Raman active properties of the structure occur due to optical phonon.

Among thermal properties of lattice, the heat capacity follows Dulong-Petit's rule, which assumes the heat capacity in solid is irrelevant to the temperature of the material, using classical thermodynamics and equipartition theorem, which satisfying the heat capacity of system converged into 0 in classically. But when measuring the heat capacity, it decreases to 0 in low temperature, and approach to Dulong-Petit value in high temperature.

To understand that, Einstein assumes that all of the vibrations of solid have the same frequency and can describe as a harmonic oscillator. In here, the average quantum number can find using Bose-Einstein distribution. For the Einstein model, it corresponds to the Dulong-Petit value in high temperature. However, it has a large difference in low temperature having not sufficient energy to excite from the ground state.

In order to solve this Einstein model, Debye proposed a Debye model (ω = νk) that linearly approximates the dispersion relationship of acoustic phonon and that this dispersion can have several k. He published the Debye T3 law which fits the experimental value better than the exponential behaviour of the Einstein model.

In conclusion, the Debye model appropriate for heat capacity in low temperature, and the Einstein model appropriate in high temperature.

Like this, you can obtain a lot of information such as thermal, optical, electronic properties, from phonon calculation.

However, the initial setting of phonon calculation is complicated, and it takes a quite long calculation time. So, it is challenging to obtain proper results.

Materials Square developed the phonon module that anyone can easily use, to overcome these difficulties.

In present, you can obtain interatomic force constant (IFC), Dielectric constant ε, Effective charge Z*, Electron-phonon coupling coefficient λ, Phonon dispersion, Phonon DOS in MatSQ phonon module.

In this posting, we'll check the method of how can obtain proper phonon results using MatSQ phonon module.

2. The Step of Phonon Calculation

Phonon calculation must start from 'fully relaxed structure (ground structure)'. It means that you should perform very 'accurate' DFT calculation. Follow the steps below to perform phonon calculation.

• Step 1. Structural relaxation
• Step 2. SCF calculation
• Step 3. Phonon calculation
• Step 4. Post-processing
1) Calculate interatomic force constant (IFC)
2-1) Calculate Electron-phonon coupling coefficient (For metal)
2-2) Calculate Dielectric constant & Effective charge (For non-metal)
3) Calculate band structure (phonon dispersion)
4) Calculate phonon DOS

In this time, Step 1 will progress in Quantum Espresso module, and Step 2~4  in Phonon module.

In Step 1, vc-relax (variable cell relaxation) has been performed in order to find the ground state of the structure. Most of the set-up parameters are ideal with normal relaxation calculation, but must use ACCURATE some parameters. You may concern about the following parameters.

• Calculation type : vc-relax
· Minimize ions and cell (full-optimization)
• Force threshold
· Use smaller value than ordinary calculations
· Phonon is a sensitive calculation
(e.g. It is accurately calculated between un-displaced and 0.01 A displaced atoms of the system)
· Recommended value is lower than 1.0e-5
• Ecut (wfc), Ecut (rho)
· Larger value (than ordinary calculation) is mandatory for very accurate calculation
· Recommended Ecut(rho) = 8-12 x Ecut(wfc)
• Cell factor (optional)
· Recommended cell factor = 2 for low-dimensional
• (Electrons) Convergence threshold
· Stop criteria for electron SCF loop
· Recommended value is lower than 1.0e-9
• K-points sampling grid
· Depending on the cell

In order to calculate phonon properties by using Quantum Espresso, you may perform SCF calculation (Step 1) before the phonon calculation (Step 2~4).
However, Materials Square not only provide one-shot calculation combining SCF calculation and phonon calculation simultaneously, but also post-processing. The procedures as follow,

(1) Add Ph.x module connected with previous calculation (Step 1. Relaxation).
(2) Set-up parameters (Check electron-phonon calculation  or Dielectric calculation if you needed)
(3) Add a tab for post-processing
(4) Run simulation

PWscf part is the Quantum Espresso PWscf input script for SCF calculation to obtain more accurate ground state electronic structure.

Ph.x part is for setting input parameters for phonon calculation, such as q-points and convergence threshold.

Quantum Espresso use Density functional perturbation theory (DFPT) for phonon calculation.

Quantum Espresso implements density functional perturbation theory (DFPT) for the calculation of phonon.
In the DFPT calculation, the response properties, such as phonon frequencies, electron-phonon interactions, have been calculated through adding external perturbations. Then, various representations generated without the equivalent ones considering symmetry operations. After generate irreducible representations, Quantum Espresso performs SCF calculation for each image to obtain the atomic forces, and calculate the force constant by gathering them. The convergence condition of each SCF calculation defined by "Convergence threshold". In the Materials Square, 10-12 Ry of convergence threshold set by default in phonon calculation, because the phonon calculation needs extremely accurate calculation. However, you should perform the phonon calculation with lower convergence threshold (i.e. 10-16 Ry) if the phonon dispersion has an imaginary mode (≈ ω < 0) although the structure is very stable.

Click the 'Update' button of the tab. The Dielectric Constant ε & Effective Charge Z* results will be shown in this tab.

Click + button to add the post process tab to perform post-processing calculation. You can obtain band structure, DOS by using this tab.

Please refer to Weekly tip #2, for the definition of k-point and setting method. For band k-path (wave vector, high symmetry point), please refer to Weekly tip #12.

Set the proper initial parameters and click Start Job! to start the calculation.

The Phonon calculation takes a long calculation time. Please check your credits before start calculation.

3. Phonon Calculation Results

After finish calculation, click the 'Update' button to check the job status.

'This Job is has been finished normally' message will appear if the job is finished successfully.

1. Interatomic force constant (IFC)

👉 Check the result

You can find the interatomic force constant information at the 'job.stdout.q2r.x' file of the data page.

The job.stdout.q2r.x file reproduced the stress tensor value of ph.x.dyn* files. If you change the sampled number of q-point, the number of dynamics files will be changed.

2. Band structure & DOS

👉 Check the result

The phonon band structure and DOS can check at the post process tab.

The wave vector sets according to the initial setting.

4. Other Properties with Phonon Calculation

You can select the 'Properties' option at Ph.x part of phonon module.

But you should set the input script of optimization carefully according to the property will obtain.

1. Dielectric Constant (ε) and Effective Charge (Z*)

You can calculate Dielectric Constant (ε) and Effective Charge (Z*) if you check the option when starting phonon calculation.

The system must be a non-metal system. Because the dielectric constant has no physical meaning for metal.

Therefore, you should change the input script like following when performing optimization.

• Occupations = 'fixed'
· To calculate 'Dielectric Constant (ε) and Effective Charge (Z*)', set this tag as 'fixed' to assume the system is non-metal.

The result will appear in the PWscf & Phonon tab.

👉 Check the result

2. Electron-Phonon Coupling Coefficient (λ)

Check 'Electron-Phonon Coupling Coefficient (λ)' option to calculate the electron-phonon coupling coefficient.

You must modify the input script like following when performing optimization calculation. Because the coefficient can obtain from the metal system only, and the 'smearing' option is for the metal system.

• Occupations = 'smearing'
· To calculate 'Electron-Phonon Coupling Coefficient (λ)', set this tag as 'smearing' to assume the system is metal.

Thus, you must add the post process tab to obtain the electron-phonon coupling coefficient. Be careful about it is cannot be calculated if the post process tab does not exist. The LAMBDA section will be added in post process tab if you set to calculate the coefficient.

The calculation result will be shown in the LAMBDA sub-tab of the post process tab.

👉 Check the result

In this post, we have found out the definition of phonon and check the result of phonon band structure, DOS graph and other phonon properties for silicon primitive cell.

For further detail of theory, please refer to our last phonon webinar of April 8. You can watch the recording on YouTube 'Materials Square' Channel.

Make your phonon calculation with Materials Square!

Materials Square