Skip to content

Mean-field calculations: tricks and hints

Tricks and hints


Vxc0 is the G=0 component of the exchange-correlation potential V_{xc}. Vxc0 is added to the eigenvalues in the wavefunction files produced by PARATEC and ESPRESSO. Vxc0 is also included in the VXC and vxc.dat files produced by PARATEC and ESPRESSO (the VXC file contains Vxc(G) and the vxc.dat file contains matrix elements of V_{xc}). The two Vxc0 terms cancel out when calculating the quasiparticle corrections in the Sigma code.

Vacuum level

To correct the DFT eigenvalues for the vacuum level take the average of the electrostatic potential on the faces of the unit cell in the non-periodic directions and subtract it from the DFT eigenvalues. The electrostatic potential is defined as V_{elec} = V_{bare} + V_{hartree}, while the total potential is V_{tot} = V_{bare} + V_{hartree} + V_{xc}, hence V_{tot} contains Vxc0 and V_{elec} does not. The average of V_{elec} is fed into the Sigma code using keywords avgpot and avgpot_outer. The potentials can be generated with PARATEC and ESPRESSO. In PARATEC, use elecplot for V_{elec} or potplot for V_{tot}, then convert from a3dr to cube using the Visual/ script. In ESPRESSO, use iflag=3, output_format=6, and plot_num=11 for V_{elec} or plot_num=1 for V_{tot}. The averaging is done with the Visual/ script.

Unit cell size

If you truncate the Coulomb interaction, make sure that the size of the unit cell in non-periodic directions is at least two times larger than the size of the charge distribution. This is needed to avoid spurious interactions between periodic replicas but at the same time not to alter interactions within the same unit cell. Run Visual/surface.x to plot an isosurface that contains 99% of the charge density (see Visual/README for instructions on how to do this). The code will print the size of the box that contains the isosurface to stdout. Multiply the box dimensions in non-periodic directions by two to get the minimum size of the unit cell.

Inversion symmetry

When using the real flavor of the code, make sure the inversion symmetry has no associated fractional translation (if it does, shift the coordinate origin). Otherwise WFN, VXC, RHO (and VSC in SAPO) have non-vanishing imaginary parts which are simply dropped in the real flavor of the code. This won't do any good to your calculation.


There are some notorious cases where typical DFT calculations might not provide a good starting point for one-shot perturbation-theory calculations. Examples include:

  • Germanium crystal, which is often predicted to be metallic at DFT. This issue can be remedied by either using another mean-field starting point, or performing some sort of self-consistent iteration, for instance, based on the static COHSEX approximation.

  • Molecules such as silane (SiH4), in which semilocal DFT often yields a LUMO orbital that bound, where in reality it is unbound. These systems can also be remedied by using a different starting mean-field point, or performing self-consistent GW calculations. For molecules, another common approach is known as best G, best W, where one picks one mean-field starting point to write the Green's function G in \Sigma=iGW, such as Hartree-Fock, and another one to compute the polarizability \chi^0 used to construct W, such as LDA.

  • Some strongly correlated systems, especially those with partially filled f orbitals. Systems such as transition-metal oxides are often tackled with a Hubbard-type of correction scheme, such as DFT+U.