Requests for technical support from the VASP group should be posted in the VASP-forum.

Low-scaling GW: The space-time formalism

Available as of VASP.6 are low-scaling algorithms for ACFDT/RPA.[1] This page describes the formalism of the corresponding low-scaling GW approach.[2] A theoretical description of the ACFDT/RPA total energies is found here. A brief summary regarding GW theory is given below, while a practical guide can be found here.

Theory

The GW implementations in VASP described in the papers of Shishkin et al.[3][4] avoid storage of the Green's function ${\displaystyle G}$ as well as Fourier transformations between time and frequency domain entirely. That is, all calculations are performed solely on the real frequency axis using Kramers-Kronig transformations for convolutions in the equation of ${\displaystyle \chi }$ and ${\displaystyle \Sigma }$ in reciprocal space and results in a relatively high computational cost that scales with ${\displaystyle N^{4}}$ (number of electrons).

The scaling with system size can, however, be reduced to ${\displaystyle N^{3}}$ by performing a so-called Wick-rotation to imaginary time ${\displaystyle t\to i\tau }$.[5]

Following the low scaling ACFDT/RPA algorithms the space-time implementation determines first, the non-interacting Green's function on the imaginary time axis in real space

${\displaystyle G({\bf {r}},{\bf {r}}',i\tau )=-\sum _{n{\bf {k}}}\phi _{n{\bf {k}}}^{*(0)}({\bf {r}})\phi _{n{\bf {k}}}^{(0)}({\bf {r}}')e^{-(\epsilon _{n{\bf {k}}}-\mu )\tau }\left[\Theta (\tau )(1-f_{n{\bf {k}}})-\Theta (-\tau )f_{n{\bf {k}}}\right]}$

Here ${\displaystyle \Theta }$ is the step function and ${\displaystyle f_{n{\bf {k}}}}$ the occupation number of the state ${\displaystyle \phi _{n{\bf {k}}}^{(0)}}$. Because the Green's function is non-oscillatory on the imaginary time axis it can be represented on a coarse grid ${\displaystyle \tau _{m}}$, where the number of time points can be selected in VASP via the NOMEGA tag. Usually 12 to 16 points are sufficient for insulators and small band gap systems.[6]

Subsequently, the irreducible polarizability is calculated from a contraction of two imaginary time Green's functions

${\displaystyle \chi ({\bf {r}},{\bf {r}}',i\tau _{m})=-G({\bf {r}},{\bf {r}}',i\tau _{m})G({\bf {r}}',{\bf {r}},-i\tau _{m})}$

Afterwards, the same compressed Fourier transformation as for the low scaling ACFDT/RPA algorithms is employed to obtain the irreducible polarizability in reciprocal space on the imaginary frequency axis ${\displaystyle \chi ({\bf {r}},{\bf {r}}',i\tau _{m})\to \chi _{{\bf {G}}{\bf {G}}'}({\bf {q}},i\omega _{n})}$.[6][2]

The next step is the computation of the screened potential

${\displaystyle W_{{\bf {G}}{\bf {G}}'}({\bf {q}},i\omega _{m})=\left[\delta _{{\bf {G}}{\bf {G}}'}-\chi _{{\bf {G}}{\bf {G}}'}({\bf {q}},i\omega _{m})V_{{\bf {G}}{\bf {G}}'}({\bf {q}})\right]^{-1}V_{{\bf {G}}{\bf {G}}'}({\bf {q}})}$

followed by the inverse Fourier transform ${\displaystyle W_{{\bf {G}}{\bf {G}}'}({\bf {q}},i\omega _{n})\to \chi ({\bf {r}},{\bf {r}}',i\tau _{m})}$ and the calculation of the self-energy

${\displaystyle \Sigma ({\bf {r}},{\bf {r}}',i\tau _{m})=-G({\bf {r}},{\bf {r}}',i\tau _{m})W({\bf {r}}',{\bf {r}},i\tau _{m})}$

From here, several routes are possible including all approximations mentioned above, that is the single-shot, EVG0 and QPEVG0 approximation. All approximations have one point in common.

In contrast to the real-frequency implementation, the low-scaling GW algorithms require an analytical continuation of the self-energy from the imaginary frequency axis to the real axis. In general, this is an ill-defined problem and usually prone to errors, since the self-energy is known on a finite set of points. VASP determines internally a Padé approximation of the self-energy ${\displaystyle \Sigma (z)}$ from the calculated set of NOMEGA points ${\displaystyle \Sigma (i\omega _{n})}$ and solves the non-linear eigenvalue problem

${\displaystyle \left[T+V_{ext}+V_{h}+\Sigma (z)\right]\left|\phi _{n{\bf {k}}}\right\rangle =z\left|\phi _{n{\bf {k}}}\right\rangle }$

on the real frequency axis ${\displaystyle z=\omega }$.

Because, preceding Fourier transformations have been carried out with exponentially suppressed errors, the analytical continuation ${\displaystyle \Sigma (z)}$ of the self-energy can be determined with high accuracy. The analytical continuation typically yields energies that differ less than 20 meV from quasi-particle energies obtained from the real-frequency calculation.[2]

In addition, the space-time formulation allows to solve the full Dyson equation for ${\displaystyle G({\bf {r,r'}},i\tau )}$ with decent computational cost.[7] This approach is known as the self-consistent GW approach (scGW) and is available as of VASP6.

Finite temperature formalism

The zero-temperature formalism of many-body perturbation theory breaks down for metals (systems with zero energy band-gap) as pointed out by Kohn and Luttinger.[8] This conundrum is lifted by considering diagrammatic perturbation theory at finite temperature ${\displaystyle T>0}$, which may be understood by an analytical continuation of the real-time ${\displaystyle t}$ to the imaginary time axis ${\displaystyle -i\tau }$. Matsubara has shown that this Wick-rotation in time ${\displaystyle t\to -i\tau }$ reveals an intriguing connection to the inverse temperature ${\displaystyle \beta =1/T}$ of the system.[9] More precisely, Matsubara has shown that all terms in perturbation theory at finite temperature can be expressed as integrals of imaginary time quantities (such as the polarizability ${\displaystyle \chi (-i\tau )}$) over the fundamental interval ${\displaystyle -\beta \leq \tau \leq \beta }$.

As a consequence, one decomposes imaginary time quantities into a Fourier series with period ${\displaystyle \beta }$ that determines the spacing of the Fouier modes. For instance the imaginary polarizability can be written as

${\displaystyle \chi (-i\tau )={\frac {1}{\beta }}\sum _{m=-\infty }^{\infty }{\tilde {\chi }}(i\nu _{m})e^{-i\nu _{m}\tau },\quad \nu _{m}={\frac {2m}{\beta }}\pi }$

and the corresponding random phase approximation of the correlation energy at finite temperature becomes a series over (in this case bosonic) Matsubara frequencies

${\displaystyle \Omega _{c}^{\rm {RPA}}={\frac {1}{2}}{\frac {1}{\beta }}\sum _{m=-\infty }^{\infty }{\rm {Tr}}\left\lbrace \ln \left[1-{\tilde {\chi }}(i\nu _{m})V\right]-{\tilde {\chi }}(i\nu _{m})V\right\rbrace ,\quad \nu _{m}={\frac {2m}{\beta }}\pi }$

The Matsubara formalism has the advantage that all contributions to the Green's function and the polarizability are mathematically well-defined, including contributions from states close to the chemical potential ${\displaystyle \epsilon _{n{\bf {k}}}\approx \mu }$, such that Matsubara series converge also for metallic systems.

Although formally convenient, Matsubara series converge poorly with the number of considered terms in practice. The VASP code, therefore, uses a compressed representation of the Fourier modes by employing the Minimax Isometry method.[10] This approach converges exponentially with the number of considered frequency points.