Phonons from density-functional-perturbation theory: Difference between revisions

From VASP Wiki
VisualWikitext
No edit summary
No edit summary
 
(12 intermediate revisions by 3 users not shown)
Line 1: Line 1:
'''Density-functional-perturbation theory''' provides a way to compute the second-order linear response to ionic displacement, strain, and electric fields. The equations are derived as follows.
The phonon calculations using [[Phonons:_Theory#Density_functional_perturbation_theory|density-functional-perturbation theory (DFPT)]] are carried out by setting [[IBRION#ibrion78|'''IBRION'''=7 or 8]] in the {{FILE|INCAR}} file.
{{NB|mind|Only zone-center (Γ-point) frequencies are calculated.}}


In density-functional theory, we solve the Kohn-Sham (KS) equations
In general, the DFPT routines in VASP are somewhat rudimentary and only support displacements commensurate with the supercell, i.e., so-called q=0 phonons.
Therefore, this approach offers few advantages over computing [[Phonons_from_finite_differences|phonons from finite differences]].
In particular, the DFPT routines are limited to LDA and GGA functionals, and they do not determine the elastic tensors, since the perturbation with respect to the strain tensor is not implemented.
The only advantage of the linear response routines is that they eliminate the need to choose the magnitude of the finite displacement {{TAG|POTIM}}.
Therefore, it might be helpful to first calculate phonon frequencies within DFPT and then switch to the [[Phonons from finite differences|finite differences approach]] in order to determine the largest displacement that will produce results compatible with the linear response routines.


:<math>
A few technical comments: VASP solves the [[Phonons:_Theory#IonSternheimer|linear Sternheimer equation]] to determine the linear response of the orbitals.
H(\mathbf{k}) | \psi_{n\mathbf{k}} \rangle=
Hence, unoccupied orbitals are not required.
e_{n\mathbf{k}}S(\mathbf{k}) | \psi_{n\mathbf{k}} \rangle,
Internally, the VASP routines for linear response rely on finite differences in two places:
</math>
#  The first place is the determination of the second derivative of the [[:Category:Exchange-correlation_functionals|exchange-correlation functional]]: Since most functionals do not support an algebraic determination of second derivatives, VASP always resorts to finite differences to determine the second-order change of the exchange correlation-potential and the PAW one-center terms for each atomic displacement.  
where ... (MTH: please define all quantities.)
# Second, after VASP has determined the first-order change of the orbitals, it computes all [[Phonons:_Theory#ForceConstantsDFPT|second derivatives using finite displacements]].  


Taking the derivative with respect to the ionic positions <math>R_i^a</math>, we obtain the Sternheimer equations
To this end, VASP displaces the selected atom in the selected directions [[Phonons:_Theory#FiniteDiffWF|adds the calculated linear response to the orbitals]], and finally determines the differences in the forces and the stress tensor for positive and negative displacements.
It can be shown that this yields precisely the [[Phonons:_Theory#InternalStrainDFPT|second-order force constants]] and the [[Phonons:_Theory#InternalStrainDFPT|internal strain tensor]], respectively.


:<math>
== Input ==
\left[ H(\mathbf{k}) - e_{n\mathbf{k}}S(\mathbf{k}) \right]
| \partial_{u_i^a}\psi_{n\mathbf{k}} \rangle
=
-\partial_{u_i^a} \left[ H(\mathbf{k}) - e_{n\mathbf{k}}S(\mathbf{k})\right]
| \psi_{n\mathbf{k}} \rangle
</math>


Once the derivative of the KS orbitals is computed from the Sternheimer equations, we can write
To use DFPT, set the tag [[IBRION#ibrion78|'''IBRION'''=7 or 8]] in the {{TAG|INCAR}} file.
:<math>
There are two options for using the DFPT routines to compute the second-order force-constants
| \psi^{u^a_i}_\lambda \rangle =  
* {{TAG|IBRION}}=7, all the atoms are displaced in all three Cartesian directions,
| \psi \rangle +
* {{TAG|IBRION}}=8, uses symmetry to reduce the number of displacements.
\lambda | \partial_{u^a_i}\psi \rangle.
</math>


The second-order response to ionic displacement, i.e., the force constants or Hessian matrix, are then computed using
If {{TAG|LEPSILON}}=.TRUE. or {{TAG|LCALCEPS}}=.TRUE., additional dielectric properties are computed.  
:<math>
\Phi^{ab}_{ij}=
\frac{\partial^2E}{\partial u^a_i \partial u^b_j}=
-\frac{\partial F^a_i}{\partial u^b_j}
\approx
-\frac{
\left(
  \mathbf{F}[\{\psi^{u^b_j}_{\lambda/2}\}]-
  \mathbf{F}[\{\psi^{u^b_j}_{-\lambda/2}\}]
\right)^a_i}{\lambda},
</math>
where <math>\mathbf{F}</math> yields the [[:Category:Forces|forces]] for a given set of KS orbitals.


MTH: Here, it would be good to explicitly write the eigenvalue equation that is solved to obtain phonon frequencies.
== Output ==


<!--
The second derivates of the total energy with respect to ionic displacements (interatomic force constants) are computed,
MTH: The following is not concerning the calculation of phonons:
the [[Phonons:_Theory#DynamicalMatrix|dynamical matrix]] is constructed, diagonalized and the phonon modes and frequencies of the system are reported after the following lines:


The internal strain tensor is computed using
  Eigenvectors and eigenvalues of the dynamical matrix
:<math>
  ----------------------------------------------------
\Xi^a_{il}=\frac{\partial^2 E}{\partial \eta^a_i \partial u^b_j}=
\frac{\partial \sigma_l}{\partial u^a_i}
\approx
\frac{
    \left(
        \sigma[\{\psi^{u^a_i}_{\lambda/2}\}]-
        \sigma[\{\psi^{u^a_i}_{-\lambda/2}\}]
    \right)_l
}{\lambda}.
</math>


At the end of the calculation if {{TAG|LEPSILON}}=.TRUE., the Born effective charges are computed using Eq. (42) of Ref. {{cite|gonze:prb:1997}}.
The mixed second derivative with respect to the strain and the ionic displacement (internal strain tensor) are evaluated and reported.
Although the contributions from the ionic relaxations to the elastic tensor are calculated, the ion-clamped elastic tensor (rigid ion) is not determined because the perturbation with respect to the strain tensor is not implemented.


:<math>
Furthermore, the [[Phonons:_Theory#BornChargeDFPT|Born effective charges]] are determined analytically by contracting the linear response of the orbitals over the "polarization" vector Eq. (30) in Ref. {{cite|gajdos:prb:2006}}.
Z^{a*}_{ij} =
These should agree well with the Born effective charges that were previously determined when the linear response with respect to external fields {{TAG|LEPSILON}}=.TRUE. was calculated (there are two different routes to calculate mixed derivatives).
2 \frac{\Omega_0}{(2\pi)^3}
The final summary output towards the end of the {{TAG|OUTCAR}} file writes the Born effective charges determined from the linear response with respect to external fields.
\int_\text{BZ}
\sum_n
\langle
  \partial_{u^a_i}\psi_{n\mathbf{k}} | \nabla_{\mathbf{k}} \tilde{u}_{n\mathbf{k}}
\rangle d\mathbf{k}
</math>


where <math>a</math> is the atom index, <math>i</math> the direction of the displacement of atom and <math>j</math> the polarization direction.
It is possible to [[Computing_the_phonon_dispersion |obtain the phonon dispersion at different '''q''' points]] by computing the force constants on a sufficiently large supercell and Fourier interpolating the dynamical matrices in the primitive cell.
The results should be equivalent to the ones obtained using {{TAG|LCALCEPS}} and {{TAG|LEPSILON}}.


MTH: This should be part of a more extensive discussion in a how-to article. The beginning of the article seemed more like a theory article, so I propose to move it in that category.
It is also possible to use phonopy{{cite|phonopy}} to use the results of a density-functional-perturbation theory calculation done with VASP.{{cite|phonopy_dfpt}}
{{NB|mind|{{TAG|IBRION}}{{=}}7 or 8 are supported by VASP.5.1 and later versions.}}


When {{TAG|IBRION}}=7 or 8 VASP solves the Sternheimer equation above with an ionic displacement perturbation.
== Related tags and sections ==
If {{TAG|IBRION}}=7 no symmetry is used and the displacement of all the ions is computed.
{{TAG|IBRION}},
When {{TAG|IBRION}}=8 only the irreducible displacements are computed with the physical quantities being reconstructed by symmetry.
{{TAG|LEPSILON}},
-->
[[Phonons: Theory]],
[[Phonons from finite differences]]


== References ==
== References ==
<references/>
<references/>


[[Category:Phonons]][[Category:Theory]]
[[Category:Phonons]][[Category:Howto]]

Latest revision as of 09:41, 23 October 2023

The phonon calculations using density-functional-perturbation theory (DFPT) are carried out by setting IBRION=7 or 8 in the INCAR file.

Mind: Only zone-center (Γ-point) frequencies are calculated.

In general, the DFPT routines in VASP are somewhat rudimentary and only support displacements commensurate with the supercell, i.e., so-called q=0 phonons. Therefore, this approach offers few advantages over computing phonons from finite differences. In particular, the DFPT routines are limited to LDA and GGA functionals, and they do not determine the elastic tensors, since the perturbation with respect to the strain tensor is not implemented. The only advantage of the linear response routines is that they eliminate the need to choose the magnitude of the finite displacement POTIM. Therefore, it might be helpful to first calculate phonon frequencies within DFPT and then switch to the finite differences approach in order to determine the largest displacement that will produce results compatible with the linear response routines.

A few technical comments: VASP solves the linear Sternheimer equation to determine the linear response of the orbitals. Hence, unoccupied orbitals are not required. Internally, the VASP routines for linear response rely on finite differences in two places:

  1. The first place is the determination of the second derivative of the exchange-correlation functional: Since most functionals do not support an algebraic determination of second derivatives, VASP always resorts to finite differences to determine the second-order change of the exchange correlation-potential and the PAW one-center terms for each atomic displacement.
  2. Second, after VASP has determined the first-order change of the orbitals, it computes all second derivatives using finite displacements.

To this end, VASP displaces the selected atom in the selected directions adds the calculated linear response to the orbitals, and finally determines the differences in the forces and the stress tensor for positive and negative displacements. It can be shown that this yields precisely the second-order force constants and the internal strain tensor, respectively.

Input

To use DFPT, set the tag IBRION=7 or 8 in the INCAR file. There are two options for using the DFPT routines to compute the second-order force-constants

  • IBRION=7, all the atoms are displaced in all three Cartesian directions,
  • IBRION=8, uses symmetry to reduce the number of displacements.

If LEPSILON=.TRUE. or LCALCEPS=.TRUE., additional dielectric properties are computed.

Output

The second derivates of the total energy with respect to ionic displacements (interatomic force constants) are computed, the dynamical matrix is constructed, diagonalized and the phonon modes and frequencies of the system are reported after the following lines: 

 Eigenvectors and eigenvalues of the dynamical matrix
 ----------------------------------------------------

The mixed second derivative with respect to the strain and the ionic displacement (internal strain tensor) are evaluated and reported. Although the contributions from the ionic relaxations to the elastic tensor are calculated, the ion-clamped elastic tensor (rigid ion) is not determined because the perturbation with respect to the strain tensor is not implemented.

Furthermore, the Born effective charges are determined analytically by contracting the linear response of the orbitals over the "polarization" vector Eq. (30) in Ref. [1]. These should agree well with the Born effective charges that were previously determined when the linear response with respect to external fields LEPSILON=.TRUE. was calculated (there are two different routes to calculate mixed derivatives). The final summary output towards the end of the OUTCAR file writes the Born effective charges determined from the linear response with respect to external fields.

It is possible to obtain the phonon dispersion at different q points by computing the force constants on a sufficiently large supercell and Fourier interpolating the dynamical matrices in the primitive cell.

It is also possible to use phonopy[2] to use the results of a density-functional-perturbation theory calculation done with VASP.[3]

Mind: IBRION=7 or 8 are supported by VASP.5.1 and later versions.

Related tags and sections

IBRION, LEPSILON, Phonons: Theory, Phonons from finite differences

References

  1. ↑ M. Gajdoš, K. Hummer, G. Kresse, J. Furthmüller, and F. Bechstedt, Phys. Rev. B 73, 045112 (2006).
  2. ↑ http://phonopy.github.io/phonopy/index.html (2022).
  3. ↑ http://phonopy.github.io/phonopy/vasp-dfpt.html (2022).
Retrieved from "https://www.vasp.at/wiki/index.php?title=Phonons_from_density-functional-perturbation_theory&oldid=22573"