A paradigm system for strong correlation and charge transfer competition

In chemistry and condensed matter physics the solution of simple paradigm systems, such as the hydrogen atom and the uniform electron gas, plays a critical role in understanding electron behaviors and developing electronic structure methods. The H$_2$ molecule is a paradigm system for strong correlation with a spin-singlet ground state that localizes the two electrons onto opposite protons at dissociation. We extend H$_2$ to a new paradigm system by using fractional nuclear charges to break the left-right nuclear symmetry, thereby enabling the competition between strong correlation and charge transfer that drives the exotic properties of many materials. This modification lays a foundation for improving practical electronic structure theories and provides an extendable playground for analyzing how the competition appears and evolves.

an equal pull on the two electrons. This left-right symmetry can be broken by replacing one proton with a helium nucleus, thus enabling charge transfer, driving both electrons to localize on the heavier He nucleus, and suppressing the strong correlation. This full replacement of strong correlation by charge transfer is unfortunate, since the competition between them drives exotic properties in many materials. It would be desirable to have a 2-electron paradigm system that can facilitate this competition, but the above analysis indicates that no such system can exist in nature given that all nuclei have integer charges.
Here, we fill this role by extending H 2 to a new paradigm system in which the competition between strong correlation and charge transfer can be tuned continuously.
To engineer such a system we create an asymmetric nuclear potential by replacing the hydrogen nuclei with fictitious fractional nuclear charges, Z A and Z B . This enables charge transfer as Z A /Z B moves away from 1 without completely suppressing strong correlation.
Much like the uniform electron gas that fertilized enormous concepts in condensed matter physics, it is of no importance that the proposed system is fictional, since the physics it captures remain deeply relevant to real correlated materials. To make the two-electron system charge neutral we choose Z A + Z B = 2 and, without loss of generality, we require Z A ≤ Z B for simplicity. We term this as "H FNC 2 ", where FNC stands for "fractional nuclear charge". Related systems were investigated by Cohen and Mori-Sánchez in the context of the DFT derivative discontinuity and delocalization error in Ref 12. For the H FNC 2 paradigm system to be of any use it requires an exact solution. Hartree-Fock (HF) theory is exact for single electron systems such as H and H + 2 (up to the chosen basis set), but it is insufficient for systems of multiple electrons. It is difficult in general to obtain exact solutions for multi-electron systems, normally requiring exponentially scaling methods such as full configuration-interaction (FCI) or QMC. Fortunately, the two electron systems of interest here are small enough such that coupled-cluster at the singles-doubles (CCSD) level is equivalent to FCI, considering all possible excitations. We note however, that CCSD is not generally reliable for strongly correlated systems with more than two electrons. It is precisely this easy availability of exact solutions that make paradigm systems valuable assets for accessing the underlying physics of complex problems. Depending on the ratio Z A /Z B , the ground state is either the |Ψ SO or |Ψ DO , or the configurations can become degenerate. This is summarized in Figure 1  The appearance of degeneracy between |Ψ SO and |Ψ DO has profound implications. It is well known that neither the 1D Hubbard model nor its material realization in the infinite 1D hydrogen chain shows charge transfer physics [13][14][15]. At short inter-atomic distances the hydrogen chain is weakly correlated and metallic, while at larger inter-atomic distances it undergoes a phase transition to a strongly correlated insulating phase [13][14][15]. The hydrogen chain is therefore a prototypical system illustrating the Mott-Hubbard metal insulator transition. Now, consider a hydrogen chain that has a sufficiently large inter-atomic distance such that the electron density overlap between atomic sites is negligible. Analogous to H FNC 2 , we can allow the nuclear charge for a pair of hydrogen atoms in the chain to be fractional under the constraint the whole system remains charge neutral. Then, following the previous analysis, the fractional nuclear charges can be tuned to make the pair close to the SO and DO degeneracy discussed above. Around this degeneracy a small perturbation, e.g., an electric field that enhances the potential at the more positive nucleus of the pair, can easily drive the electron from the less positive to the more positive nuclear site. This charge transfer capability under small perturbation emerging from the insulating hydrogen chain highlights the rich physics brought by the fractional nuclear charge. If more fractionally charged pairs are present then more degenerate states can be generated by tuning the fractional nuclear charges, potentially leading to exotic properties including superconductivity, even within this 1D model [16].
The exact solutions established above can highlight important deficiencies in common approximate electronic structure techniques. Here we focus on DFT, which has become a mainstay of computational materials studies. In principle, DFT is exact for the ground state energy and electron density through an efficient mapping of the interacting-electron problem onto an auxiliary non-interacting electron system described by a single determinant. In practice however, DFT methods must approximate the exchange-correlation energy functional that carries the many-electron effects. Paradigm systems have played critical roles in the development of exchange-correlation approximations [3,[17][18][19], with each greatly enhancing the functional's predictive power when smoothly incorporated with other constraints. This role is played by the uniform electron gas for LSDA [2][3][4] and the hydrogen atom for the strongly-constrained and appropriately-normed (SCAN) density functional and its r 2 SCAN revision [18,19]. Hartree-Fock (HF) theory is also included here for comparison as it uses a single determinant to directly approximate the correlated wave function and is a base upon which many more sophisticated methods are built.
The exact ground state wave function of the H  [20]. Recently, the symmetry broken solutions revealing the strong correlation have been interpreted as "freezing" a fluctuation in the exact correlated ground state wave function [8]. Given this interpretation and the improved energies for strongly correlated systems we adopt the spin symmetry breaking strategy for DFT approximations to be discussed below.
We have selected four non-empirical DFT exchange-correlation functionals as examples from different levels of the Perdew-Schmidt hierarchy [21]. The Perdew-Burke-Ernzerhof (PBE) functional [17] is a standard at the generalized gradient approximation (GGA) level and is the simplest semi-local functional featured, taking only the spin density and its gradient as inputs. The meta-GGA level (the most sophisticated semi-local level) is represented by our recent r 2 SCAN functional [18,19], which includes the non-negative kinetic energy density as an additional ingredient that can be used to satisfy more exact constraints. Beyond the semi-local functionals we take the PBE0 functional [17,22], which replaces 25% of the PBE exchange with 25% of the non-local exact exchange of HF. Finally, we include the r 2 SCAN functional with the Perdew-Zunger self-interaction correction (PZ-SIC+r 2 SCAN c ) [23] in which the self-interaction error is removed on an orbital-by-orbital basis, equivalent to   Given the good performance of PZ-SIC+r 2 SCAN c at infinite separation, this highlights the challenging problem of delivering accuracy for both regions dictated by self-interaction errors and multi-center non-local strong correlation. We therefore expect that H FNC 2 can be a powerful tool for developing the non-local density functionals that have been the focus of much recent DFT development [28][29][30][31][32][33][34].
It is well accepted that DFT with sophisticated exchange-correlation approximations have better accuracy than HF, and that accuracy generally improved when climbing up the the Perdew-Schmidt hierarchy, e.g., from PBE, to r 2 SCAN, and to PBE0. This is consistent with the observation in Figure 2 that general performance is improved from HF, to PBE, to r 2 SCAN, and to PBE0, shown by smaller error scales and overall smaller regions of error.
Similarly, PZ-SIC has been shown as an effective correction to DFT approximations for correlated materials due to the removal of self-interaction errors [35,36]. Correcting DFT with PZ-SIC deteriorates accuracy for normal materials however, an effect which has been called "the paradox for PZ-SIC" [37]. This agrees with the increased error found around the   site correlation localizing electrons into d bands [38]. Density functional methods have typically struggled with such materials, suffering from the self-interaction error that leads to a spurious charge delocalization between the metal and oxygen ions [39].  Table I with the range of errors observed in Figure 2, we see the large region of delocalization error for PBE are reflected in underestimated magnetic moments and qualitatively incorrect band gaps. The region of delocalization error is smaller for r 2 SCAN and correspondingly the material predictions are improved, with all materials correctly insulating though significant underestimation of band gaps remains. The partial self-interaction error correction from the presented can be easily extended to more complex multi-orbital systems, offering a clear and practical sandbox for one of the largest problems remaining in the physical sciences.
wrote the manuscript. J.S. provided computational resources. All authors contributed to editing the manuscript.

DATA AVAILABILITY STATEMENT
Data for Figures 1 and 2 is available from the authors by request.

Fractional Nuclear Charges
Fractional nuclear charges are implemented under the Born-Oppenheimer approximation by assigning desired Z i ∈ R + to each nucleus and evaluating the nuclear-electron attraction and nuclear repulsion integrals in the standard way. This modification is trivial for most existing electronic structure codes and is available in the standard Turbomole release used for this work [41,42]. Total energies calculated using Turbomole V7.4 as the sum of two independent atomic calculations with fractional nuclear charges. The d-aug-cc-pV5Z hydrogen basis functions [26] were used for all atomic calculations. Fractional electron occupation was determined numerically by adjusting the occupation fraction on each atomic fragment (fixed such that the total system contains two electrons) to minimize self-consistent total energy. Total energies calculated for fractional nuclei at finite separation compared to CCSD references, all calculated using Turbomole V7.4. The d-aug-cc-pVQZ hydrogen basis functions [26] were used for all calculations, no basis set superposition error (BSSE) corrections were applied. Coulson-Fisher points were evaluated at regular steps in Z A /Z B by numerically searching for the minimum bond length R (±0.05 Bohr) where the spin restricted and spin unrestricted total energies differed by > 10 −4 eV. No point was recorded if no separation was found < 9 Bohr.

Table I
All materials are in the G-type AFM phase. Calculations use the pseudopotential projector-augmented wave method [43] as implemented in the Vienna ab initio simulation package (VASP) [44,45]. A high-energy cutoff of 500 eV was used to truncate the plane-wave basis set. * jfurness@tulane.edu † jsun@tulane.edu