$GW$ density matrix for estimation of self-consistent $GW$ total energies in solids

The $GW$ approximation is a well-established method for calculating ionization potentials and electron affinities in solids and molecules. For numerous years, obtaining self-consistent $GW$ total energies in solids has been a challenging objective that is not accomplished yet. However, it was shown recently that the linearized $GW$ density matrix permits a reliable prediction of the self-consistent $GW$ total energy for molecules [F. Bruneval, M. Rodríguez-Mayorga, P. Rinke, and M. Dvorak, J. Chem. Theory Comput. 17, 2126 (2021)] for which self-consistent $GW$ energies are available. Here we implement, test, and benchmark the linearized $GW$ density matrix for several solids. We focus on the total energy, lattice constant, and bulk modulus obtained from the $GW$ density matrix and compare our findings to more traditional results obtained within the random-phase approximation (RPA). We conclude on the improved stability of the total energy obtained from the linearized $GW$ density matrix with respect to the mean-field starting point. We bring compelling clues that the RPA and the $GW$ density matrix total energies are certainly close to the self-consistent $GW$ total energy in solids if we use hybrid functionals with enriched exchange as a starting point.

Figure 1

HF total energies as a function of lattice constant with different codes and pseudopotentials for bulk silicon. The equilibrium lattice constant for each calculation is given in the legend.

Figure 2

Spin-summed natural occupations obtained from PBE, ${\gamma }^{\mathrm{HF}}$@PBE, and ${\gamma }^{GW}$@PBE density matrices in bulk silicon for k point (1/8, 0, 0) in reciprocal lattice vectors. The natural occupations are ordered by descending values, from 1 to 150, which is $N_b$, the dimension of the matrix. The $x$ axis has been cut to show the first and the last values.

Figure 3

Electron count deviation $\mathrm{\Delta }N\left(\mathbf{k}\right)$ as a function of the wave vector $\mathbf{k}$ in the first Brillouin zone for bulk silicon. $\mathbf{k}=k_1{\mathbf{b}}_1+k_2{\mathbf{b}}_2+k_3{\mathbf{b}}_3$ is reported in reduced coordinates (${\mathbf{b}}_i$ are the reciprocal lattice vectors). The plane $k_3=0$ is represented. The isoline $\mathrm{\Delta }N\left(\mathbf{k}\right)=0$ is drawn with a bold black line. Special points $\mathrm{\Gamma }$ and $X$ are marked.

Figure 4

Total energies for the water molecule in the gas phase (left-hand panel) and crystalline silicon (right-hand panel) as a function of the gKS starting point. The water results were extracted from a previous study [21]. The silicon QMC value comes from Ref. [64].

Figure 5

Energy-lattice constant curves for crystalline silicon for RPA functional (left-hand panel) and ${\gamma }^{GW}$ energy functional (central panel). Equilibrium lattices are summarized in the right-hand panel. The experimental lattice constant is given as a reference.

Figure 6

Energy as a function of the spacing between the layers for h-BN with different gKS approximations and different starting points for RPA.