Bond Breaking and Bond Formation: How Electron Correlation is Captured in Many-Body Perturbation Theory and Density-Functional Theory

For the paradigmatic case of H2-dissociation we compare state-of-the-art many-body perturbation theory (MBPT) in the GW approximation and density-functional theory (DFT) in the exact-exchange plus random-phase approximation for the correlation energy (EX+cRPA). For an unbiased comparison and to prevent spurious starting point effects both approaches are iterated to full self-consistency (i.e. sc-RPA and sc-GW). The exchange-correlation diagrams in both approaches are topologically identical, but in sc-RPA they are evaluated with non-interacting and in sc-GW with interacting Green functions. This has a profound consequence for the dissociation region, where sc-RPA is superior to sc-GW. We argue that for a given diagrammatic expansion, the DFT framework outperforms the many-body framework when it comes to bond-breaking. We attribute this to the difference in the correlation energy rather than the treatment of the kinetic energy.

For the paradigmatic case of H2-dissociation we compare state-of-the-art many-body perturbation theory (MBPT) in the GW approximation and density-functional theory (DFT) in the exactexchange plus random-phase approximation for the correlation energy (EX+cRPA). For an unbiased comparison and to prevent spurious starting point effects both approaches are iterated to full selfconsistency (i.e. sc-RPA and sc-GW ). The exchange-correlation diagrams in both approaches are topologically identical, but in sc-RPA they are evaluated with non-interacting and in sc-GW with interacting Green functions. This has a profound consequence for the dissociation region, where sc-RPA is superior to sc-GW . We argue that for a given diagrammatic expansion, sc-RPA outperforms sc-GW when it comes to bond-breaking. We attribute this to the difference in the correlation energy rather than the treatment of the kinetic energy.
First-principles electronic-structure calculations have become indispensable in many fields of science, because they yield atomistic insight and are complementary to purely experimental studies. Since the full many-body problem of interacting electrons and nuclei is intractable for all but the simplest systems, different strategies for approximate approaches have been developed over the years. The most prominent are densityfunctional theory (DFT) [1-3], many-body perturbation theory (MBPT) [4][5][6], coupled-cluster theory [7] and quantum Monte Carlo methods [8]. Each approach has its strengths and weaknesses in terms of accuracy, applicability, and computational efficiency and no consensus has been reached regarding the optimal approach for current and future challenges in electronic-structure theory. In this work we address the difference between DFT and MBPT for the total energy and ask the questions: Given a fixed set of diagrams for the electron-electron interaction, will the DFT and the MBPT framework give the same result? And if not, which one is better? To answer these questions we consider the paradigmatic case of H 2 dissociation. Other diatomic molecules are presented in the Supplemental Material.
In the past, DFT and MBPT have been compared directly in the exchange-only case [11]. In MBPT this corresponds to the Hartree-Fock approach, whereas in DFT a multiplicative Kohn-Sham (KS) potential is constructed by means of the optimized effective potential approach (OEP) [12]. As we will demonstrate in this Letter, the comparison between DFT and MBPT can be extended to encompass correlation using exactexchange plus correlation in the random-phase approximation to DFT (EX+cRPA), [13][14][15][16] referred to as RPA in the following, and the GW approach to MBPT [6,17]. The exchange-correlation diagrams in both approaches are topologically identical (see Fig. 1), but in RPA they are evaluated with a non-interacting KS and in GW with an interacting Green function. To illustrate the impact of these differences we consider the bondbreaking/formation regimes in the binding curves of H 2 . To avoid starting point effects both approaches are iterated to self-consistency, which we denote as sc-RPA and sc-GW . The extension of this study to higher order correlation diagrams is, in principle, possible and will be pursued in future work. Here, we focus on sc-RPA and sc-GW as they provide the simplest (and currently only computationally tractable) way for more complex systems to compare density-functional and many-body theory.
Let us start with the ground-state total-energy expression for an interacting electron system obtained with the  Hartree-Fock (HF) and exact-exchange OEP (OEPx) are identical for H2 and are included for comparison. All calculations were performed using a Gaussian cc-pVQZ [19] basis set.
adiabatic-connection (AC) technique (see e.g. Ref. [16]): Here v(r, r ′ ) is the Coulomb interaction, Tr [AB] denotes drdr ′ A(r, r ′ )B(r ′ , r) and E 0 = T s + E H + E ext . T s is the kinetic energy of the KS independent-particle system, E H the Hartree and E ext the external energy. Along the AC path (i.e. at each value of λ), the electron density n(r) is assumed to be fixed at its physical value, and χ λ , G λ , and Σ λ are the polarizability, the single-particle (timeordered) Green function, and the self-energy, respectively. In the following we adopt the notation G s ≡ G λ=0 for the non-interacting KS Green function and G ≡ G λ=1 for the fully interacting one.
The RPA for the total energy can be most conveniently introduced in Eq. 1 through the approximation Within this approximation, the integrand in Eq. 1 is assumed to depend on λ only through the scaled Coulomb interaction λv. Alternatively the RPA total energy can also be obtained through Eq. 2 by introducing the GW approximation for the proper self-energy [16,20]: where W λ [G λ ] = λv(1 + iλvG λ G λ ) −1 and η is a positive infinitesimal. The RPA total energy is retrieved by omitting the λ-dependence of G λ , i.e. replacing G λ by the KS non-interacting Green function G s = G λ=0 and Combining Eqs. 1 and 4 allows us to express the RPA total energy functional as: We now come to the differences in the evaluation of the total energy in the context of KS-DFT and MBPT. In MBPT, the Green function G λ represents an interacting electron system for λ = 0, and has to satisfy the Dyson equation: (6) with v λ ext being the external potential of the λ-dependent system (chosen to keep the density fixed), and v xc the exchange-correlation potential of the KS non-interacting particle reference system. Making use of Eqs. 3 and 6, the λ-integration in Eq. 2 can be carried out and one arrives at the following expression for the total energy Details for the derivation of Eq. 7 can be found in the supplemental material [21]. In Eq. 7, the functional Φ[G] is defined as [9,10] where Σ (n) is the sum of all self-energy diagrams that contain n explicit Coulomb interaction lines. We note that since Σ = δΦ/δG , an approximation for Φ directly translates into a corresponding approximation for Σ. The diagrammatic representation of Φ in the GW approximation is illustrated in Fig. 1.
In the KS framework, the sc-RPA total energy is obtained by requiring G s in Eq. 5 to satisfy the Dyson equa- is determined by the optimized effective potential equation (also known as the linearized Sham-Schlüter equation) [22,23]: Alternatively, the sc-RPA energy can be obtained by minimizing E RPA in Eq. 5 with respect to the non-interacting input KS Green functions G s . Regarding the energy expression of Eq. 7 as a functional of G yields the well-known Klein functional [24]. This functional is stationary (i.e., δE[G]/δG = 0) at the self-consistent G of the Dyson equation [24]. It has further been shown [23,25] that evaluating the Klein functional (using the GW approximation for Φ) with the KS reference Green function G s one obtains the RPA total energy in Eq. 5 (see also the supplemental material [21]). This offers a second way to look at the difference between sc-RPA and sc-GW : the sc-RPA energy corresponds to a mininum of the Klein functional within a variational subspace of non-interacing KS Green functions, whereas the sc-GW total energy corresponds to a stationary point of the Klein functional in a larger variational space including both noninteracting and interacting Green functions. However, we emphasize that this stationary point is not necessarily a minimum [10,24]. In practical calculations, the sc-GW total energy is actually above the sc-RPA energy as we will show in this Letter.
For a quantitative comparison between sc-GW and sc-RPA, we choose the Galitskii-Migdal (GM) formula [26] for the computation of the sc-GW total energy. At selfconsistency, the GM formula is coincides with the Klein functional (Eq. 7) -as for instance discussed in Refs. [25,27,28]. The GM formula can be expressed as [29]: where all terms on the right hand side of Eq. 10 are regarded as functionals of the Green function G. T is the kinetic energy of the interacting system and E GW c the so-called GW correlation energy defined as Here Σ GW c is the correlation part of the GW self-energy. The evaluation of Eq. 10 with a HF (PBE) Green function is referred to as G 0 W 0 @HF (G 0 W 0 @PBE) total energy.
An inspection of Eq. 4 and 11 reveals that the difference in the sc-GW correlation energy and the sc-RPA correlation energy is twofold. First, the sc-GW expression is evaluated with an interacting Green function as opposed to a Kohn-Sham one in sc-RPA. Second, the kinetic correlation energy -i.e., the difference between the full kinetic energy and that of the non-interacting KS system -are included in Eq. 4 through the coupling constant integration, whereas in sc-GW the correlation term is purely Coulombic. To facilitate a term-by-term comparison between sc-GW and sc-RPA total energies, we separate E RPA c into the Coulomb correlation energy U RPA c and the kinetic correlation energy T RPA c : and The kinetic energy in sc-RPA is then given by: With this reorganization of terms, the kinetic energy in sc-GW can be directly compared to T RPA , and similarly E GW c to U RPA c . Now the only factor responsible for the difference in these different pairs of terms arises from the difference in the input Green functions used to evaluate them.
We determined the RPA correlation potential following the direct minimization scheme of Yang et al. [30]. The resulting orbitals and eigenvalues were used to evaluate the sc-RPA total energy from Eq. 5. We refer to a previous publication for details of the sc-RPA implementation [31]. The sc-GW method -based on the iterative solution of Eqs. 3 and 6 at λ = 1 -has been implemented in the all-electron localized basis code FHI-aims [32], as explained in more detail in Refs. [29,33]. The sc-GW total energy was then obtained from Eq. 10.
We now turn to an assessment of sc-RPA and sc-GW for the potential energy curve of H 2 . In the Supplemental Material, we also show data for other covalently bonded dimers, such as LiH and Li 2 [21]. In the following, we explicitly refer to non-self-consistent calculations by appending the suffix @input to label the Green function used as input. Figure 2 reports the total energy of H 2 for different flavors of GW and RPA. For comparison we reproduce the full configuration interaction (CI) potential energy curve of H 2 [18], that provides an exact reference for this system. We also report the total energy of H 2 evaluated from a beyond-GW /RPA approach that incorporates second-order screened exchange (SOSEX) and renormalized single-excitations in the self-energy [34], referred to in the following as renormalized second order perturbation theory (rPT2) [16,35]. As reported previously [16,[36][37][38], non-self-consistent RPA overestimates the total energy of H 2 at the equilibrium bond length. Around the equilibrium distance, the RPA total energy based on exact exchange (OEPx) and sc-RPA are almost identical and overestimate the total energy by approximately 0.8 eV, compared to full-CI. At intermediate bond distances and in the dissociation region we see a lowering of the sc-RPA energy compared to RPA@OEPx. The spurious "bump" [16,36], present in all RPA calculations for H 2 and other covalently bonded molecules, is reduced in sc-RPA but is still present. The total energy stays below the full-CI energy throughout, indicating a general overestimation of the bonding and dissociation regions.
In agreement with Stan et al. [39], sc-GW provides an accurate total energy for H 2 close to equilibrium. For the Galitskii-Migdal framework, self-consistency is crucial as G 0 W 0 @HF and G 0 W 0 @PBE largely overestimate the total energy. In contrast, the Klein functional evaluated with the HF Green function (RPA@HF) yields results similar to sc-GW . sc-RPA and sc-GW thus provide a qualitatively similar description of the energetics of the covalent bond of H 2 , which results in a slight overestimation of the total energy (see Table I). However, sc-GW is in better agreement with full-CI. Most interestingly, the sc-GW energy is higher than the sc-RPA one. This is in contrast to the exchange only case, in which the HF total energy is always lower than (or equal for a two electron system) the OEPx energy [12]. This is expected, as HF is variational and the local potential in OEPx provides an additional constraint that increases the energy. Conversely, the total energy in sc-GW has to be higher than in sc-RPA, because the variational procedure yields a maximum at the self-consistent Green function [10,24].
In the dissociation region, sc-RPA and sc-GW deviate markedly. For sc-RPA the dissociation energy is below the full-CI energy and is in rather good agreement with the reference. sc-GW , on the other hand, fails dramatically in the dissociation limit and with -24.5 eV underestimates the total energy considerably. On the plus side, sc-GW dissociates monotonically and therefore does not show the unphysical "bump" present in all RPA-based approaches. Again, both non-self-consistent G 0 W 0 @HF and G 0 W 0 @PBE energies give better agreement with the reference curve than sc-GW .
One could surmise that this qualitatively different behavior originates from the different treatment of the kinetic energies that we discussed earlier. Figure 3, however, shows that this is not the case. At equilibrium the kinetic energy in sc-GW differs only slightly from the sc-RPA kinetic energy defined in Eq. 13. This indicates that in the bonding regime the AC framework correctly reproduces the kinetic energy of an interacting system. At larger bond distances, the kinetic energy differs increasingly in the two approaches. However, this effect is averaged out by an opposing change in the external energy that arises from an increasing deviation in the electron densities. The same is observed for the Hartree and the exchange energy, although the absolute magnitude of the effect is smaller. The total energy difference between sc-GW and sc-RPA can be finally ascribed to the Coulomb correlation energy , as the left panel of Fig. 3 demonstrates. Close to equilibrium ∆E c is of the order of 1 eV, but increases to approximately 4 eV at larger bond lengths. This illustrates that it matters decisively whether the correlation energy is evaluated with the interacting sc-GW or the non-interacting sc-RPA Green function.
Why the difference is so pronounced at dissociation is still an open question. A potential explanation can be found in the inverse dependence of the RPA Coulomb correlation energy on the gap between the highest occupied and the lowest unoccupied molecular orbital (HOMO and LUMO, respectively). This is exemplified by the right panel of Fig. 4, which shows the inverse of U c as a function of the gap for a simplified two level system. The large value of U c obtained from sc-RPA for H 2 at dissociation can therefore be traced back to the small HOMO-LUMO gap (left panel of Fig. 4) of the RPA Green function, as also illustrated in the Supplemental Material of Ref. 38. In contrast, due to the spatial non-locality of the self-energy, the HOMO-LUMO gap of the HF and GW Green functions is much larger at every given bond distance. This leads in turn to a smaller Coulomb cor-  relation energy for sc-GW and HF-based perturbative methods.
In conclusion, we have compared MBPT in the GW approximation to DFT in the RPA. We found that the density functional description is superior at dissociation, yielding a total energy in qualitative agreement with the exact energy along the entire dissociation curve. These results illustrate how MBPT and DFT based approaches deal with multi-reference ground-states. We demonstrated that in a DFT-based framework the closure of the (KS) HOMO-LUMO gap is in part responsible for the improved description at dissociation, i.e., static correlation is better accounted for in sc-RPA, than in sc-GW . The same effect in Green function theory has to be achieved by the right (potentially infinite) set of diagrams. We conclude that static and local approximations of exchange-correlation potentials -as opposed to nonlocal, frequency dependent self-energy approximationsare more effective in describing the dissociation regime of covalently bonded molecules.