diff --git a/Manuscript/fig4.pdf b/Manuscript/fig4.pdf index 8cadb05..d27de77 100644 Binary files a/Manuscript/fig4.pdf and b/Manuscript/fig4.pdf differ diff --git a/Manuscript/ufGW.tex b/Manuscript/ufGW.tex index ae9a6f1..6bcd44e 100644 --- a/Manuscript/ufGW.tex +++ b/Manuscript/ufGW.tex @@ -127,7 +127,7 @@ \begin{abstract} By recasting the non-linear frequency-dependent $GW$ quasiparticle equation into a linear eigenvalue problem, we explain the appearance of multiple solutions and unphysical discontinuities in various physical quantities computed within the $GW$ approximation. Considering the $GW$ self-energy as an effective Hamiltonian, it is shown that these issues are key signatures of strong correlation in the $(N\pm1)$-electron states and can be directly related to the intruder state problem. -A simple and efficient regularization procedure inspired by the similarity renormalization group is proposed to avoid such issues. +A simple and efficient regularization procedure inspired by the similarity renormalization group is proposed to avoid such issues and speed up convergence of partially self-consistent $GW$ calculations. %\bigskip %\begin{center} % \boxed{\includegraphics[width=0.5\linewidth]{TOC}} @@ -140,7 +140,7 @@ A simple and efficient regularization procedure inspired by the similarity renor \section{Introduction} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -The $GW$ approximation of many-body perturbation theory, \cite{Hedin_1965,Martin_2016} allows to compute accurate charged excitation (\ie, ionization potentials, electron affinities and fundamental gaps) in solids and molecules. \cite{Aryasetiawan_1998,Onida_2002,Reining_2017,Golze_2019} +The $GW$ approximation of many-body perturbation theory \cite{Hedin_1965,Martin_2016} allows to compute accurate charged excitation (\ie, ionization potentials, electron affinities and fundamental gaps) in solids and molecules. \cite{Aryasetiawan_1998,Onida_2002,Reining_2017,Golze_2019} Its popularity in the molecular electronic structure community is rapidly growing \cite{Ke_2011,Bruneval_2012,Bruneval_2013,Bruneval_2015,Blase_2016,Bruneval_2016, Bruneval_2016a,Koval_2014,Hung_2016,Blase_2018,Boulanger_2014,Li_2017,Hung_2016,Hung_2017,vanSetten_2015,vanSetten_2018,vanSetten_2015, Maggio_2017,vanSetten_2018,Richard_2016,Gallandi_2016,Knight_2016,Dolgounitcheva_2016,Bruneval_2015,Krause_2015,Govoni_2018,Caruso_2016} thanks to its relatively low computational cost \cite{Foerster_2011,Liu_2016,Wilhelm_2018,Forster_2021,Duchemin_2021} and somehow surprising accuracy for weakly-correlated systems. \cite{Korbel_2014,vanSetten_2015,Caruso_2016,Hung_2017,vanSetten_2018,Bruneval_2021} The idea behind the $GW$ approximation is to recast the many-body problem into a set of non-linear one-body equations. The introduction of the self-energy $\Sigma$ links the non-interacting Green's function $G_0$ to its fully-interacting version $G$ via the following Dyson equation: @@ -172,21 +172,21 @@ Electron correlation is then explicitly incorporated into one-body quantities vi %\end{subequations} %where $v$ is the bare Coulomb interaction, $\delta(12)$ is Dirac's delta function and $(1)$ is a composite coordinate gathering spin, space and time variables $(\sigma_1,\boldsymbol{r}_1,t_1)$. -In recent studies, \cite{Loos_2018b,Veril_2018,Loos_2020e,Berger_2021,DiSabatino_2021} we discovered that one can observe (unphysical) irregularities and/or discontinuities in the energy surfaces of several key quantities (ionization potential, electron affinity, fundamental gap, total and correlation energies, as well as vertical excitation energies) even in the weakly-correlated regime. +In recent studies, \cite{Loos_2018b,Veril_2018,Loos_2020e,Berger_2021,DiSabatino_2021} we discovered that one can observe (unphysical) irregularities and/or discontinuities in the energy surfaces of several key quantities (ionization potential, electron affinity, fundamental and optical gaps, total and correlation energies, as well as excitation energies) even in the weakly-correlated regime. These issues were discovered in Ref.~\onlinecite{Loos_2018b} while studying a model two-electron system \cite{Seidl_2007,Loos_2009a,Loos_2009c} and they were further investigated in Ref.~\onlinecite{Veril_2018}, where we provided additional evidences and explanations of these undesirable features in real molecular systems. In particular, we showed that each branch of the self-energy $\Sigma$ is associated with a distinct quasiparticle solution, and that each switch between solutions implies a significant discontinuity in the quasiparticle energy due to the transfer of weight between two solutions of the quasiparticle equation. \cite{Veril_2018} Multiple solution issues in $GW$ appears frequently, \cite{vanSetten_2015,Maggio_2017,Duchemin_2020} especially for orbitals that are energetically far from the Fermi level, such as in core ionized states. \cite{Golze_2018,Golze_2020} -In addition to obvious irregularities on potential energy surfaces that hampers the accurate determination of properties such as equilibrium bond lengths and harmonic vibrational frequencies, \cite{Loos_2020e,Berger_2021} one direct consequence of these discontinuities is the difficulty to converge (partially) self-consistent $GW$ calculations as the self-consistent procedure jumps erratically from one solution to the other even if convergence accelerator techniques such as DIIS are employed. \cite{Pulay_1980,Pulay_1982,Veril_2018} +In addition to obvious irregularities in potential energy surfaces that hampers the accurate determination of properties such as equilibrium bond lengths and harmonic vibrational frequencies, \cite{Loos_2020e,Berger_2021} one direct consequence of these discontinuities is the difficulty to converge (partially) self-consistent $GW$ calculations as the self-consistent procedure jumps erratically from one solution to the other even if convergence accelerator techniques such as DIIS are employed. \cite{Pulay_1980,Pulay_1982,Veril_2018} Note in passing that the present issues do not only appear in $GW$ as the $T$-matrix \cite{Romaniello_2012,Zhang_2017,Li_2021b,Loos_2022} and second-order Green's function (or second Born) formalisms \cite{SzaboBook,Casida_1989,Casida_1991,Stefanucci_2013,Ortiz_2013,Phillips_2014,Phillips_2015,Rusakov_2014,Rusakov_2016,Hirata_2015,Hirata_2017} exhibit the same drawbacks. It was shown that these problems can be tamed by using a static Coulomb-hole plus screened-exchange (COHSEX) \cite{Hedin_1965,Hybertsen_1986,Hedin_1999,Bruneval_2006} self-energy \cite{Berger_2021} or by considering a fully self-consistent $GW$ scheme, \cite{Stan_2006,Stan_2009,Rostgaard_2010,Caruso_2012,Caruso_2013,Caruso_2013a,Caruso_2013b,Koval_2014,Wilhelm_2018} where one considers not only the quasiparticle solution but also the satellites at each iteration. \cite{DiSabatino_2021} However, none of these solutions is completely satisfying as a static approximation of the self-energy can induce significant loss in accuracy and fully self-consistent calculations can be quite challenging in terms of implementation and cost. In the present article, via an upfolding process of the non-linear $GW$ equation, \cite{Bintrim_2021a} we provide further physical insights into the origin of these discontinuities by highlighting, in particular, the role of intruder states. -Inspired by regularized electronic structure theories, \cite{Lee_2018a,Evangelista_2014b} these new insights allow us to propose a cheap and efficient regularization scheme in order to avoid these issues. +Inspired by regularized electronic structure theories, \cite{Lee_2018a,Evangelista_2014b} these new insights allow us to propose a cheap and efficient regularization scheme in order to avoid these issues and speed up convergence of partially self-consistent $GW$ calculations. -Here, we consider the one-shot {\GOWO} \cite{Strinati_1980,Hybertsen_1985a,Hybertsen_1986,Godby_1988,Linden_1988,Northrup_1991,Blase_1994,Rohlfing_1995,Shishkin_2007} for the sake of simplicity but the same analysis can be performed in the case of (partially) self-consistent schemes such as ev$GW$, \cite{Hybertsen_1986,Shishkin_2007,Blase_2011,Faber_2011,Rangel_2016} where one updates only the quasiparticle energies or qs$GW$, \cite{Gui_2018,Faleev_2004,vanSchilfgaarde_2006,Kotani_2007,Ke_2011,Kaplan_2016} where both quasiparticle energies and orbitals are updated at each iteration. +Here, we consider the one-shot {\GOWO} \cite{Strinati_1980,Hybertsen_1985a,Hybertsen_1986,Godby_1988,Linden_1988,Northrup_1991,Blase_1994,Rohlfing_1995,Shishkin_2007} for the sake of simplicity but the same analysis can be performed in the case of (partially) self-consistent schemes such as ev$GW$, \cite{Hybertsen_1986,Shishkin_2007,Blase_2011,Faber_2011,Rangel_2016} where one updates only the quasiparticle energies, and qs$GW$, \cite{Gui_2018,Faleev_2004,vanSchilfgaarde_2006,Kotani_2007,Ke_2011,Kaplan_2016} where both quasiparticle energies and orbitals are updated at each iteration. Moreover, we consider a Hartree-Fock (HF) starting point but it can be straightforwardly extended to a Kohn-Sham starting point. Throughout this article, $p$ and $q$ are general (spatial) orbitals, $i$, $j$, $k$, and $l$ denotes occupied orbitals, $a$, $b$, $c$, and $d$ are vacant orbitals, while $m$ labels single excitations $i \to a$. Atomic units are used throughout. @@ -211,7 +211,7 @@ Within the Tamm-Dancoff approximation (that we enforce here for the sake of simp \begin{equation} \ERI{pq}{m} = \sum_{ia} \ERI{pq}{ia} X_{ia,m}^\RPA \end{equation} -where $\Om{m}{\RPA}$ and $\bX{m}{\RPA}$ are respectively the $m$th eigenvalue and eigenvector of the random-phase approximation (RPA) problem in the Tamm-Dancoff approximation, \ie, +where $\Om{m}{\RPA}$ and $\bX{m}{\RPA}$ are respectively the $m$th eigenvalue and eigenvector of the random-phase approximation (RPA) problem, \ie, \begin{equation} \bA{}{\RPA} \cdot \bX{m}{\RPA} = \Om{m}{\RPA} \bX{m}{\RPA} \end{equation} @@ -277,7 +277,7 @@ The size of this eigenvalue problem is $1 + \Nocc^2 \Nvir + \Nocc \Nvir^2 = \ord Thus, this step scales as $\order*{\Norb^9}$ with conventional diagonalization algorithms. Note, however, that the blocks $\bC{}{\text{2h1p}}$ and $\bC{}{\text{2p1h}}$ do not need to be recomputed for each orbital. -It is paramount to understand that diagonalizing $\bH^{(p)}$ [see Eq.~\eqref{eq:Hp}] is completely equivalent to solving the quasiparticle equation \eqref{eq:qp_eq}. +It is crucial to understand that diagonalizing $\bH^{(p)}$ [see Eq.~\eqref{eq:Hp}] is completely equivalent to solving the quasiparticle equation \eqref{eq:qp_eq}. This can be further illustrated by expanding the secular equation associated with Eq.~\eqref{eq:Hp} \begin{equation} \det[ \bH^{(p)} - \omega \bI ] = 0 @@ -340,7 +340,7 @@ This example was already considered in our previous work \cite{Veril_2018} but h The downfolded and upfolded {\GOWO} schemes have been implemented in the electronic structure package QuAcK \cite{QuAcK} which is freely available at \url{https://github.com/pfloos/QuAcK}. These calculations are based on restricted HF eigenvalues and orbitals. We denote as $\ket*{1\Bar{1}}$ the $\Ne$-electron ground-state Slater determinant where the orbital 1 is occupied by one spin-up and one spin-down electron. -Similarly notations will be employed for the $(\Ne\pm1)$-electron configurations. +Similar notations will be employed for the $(\Ne\pm1)$-electron configurations. In Fig.~\ref{fig:H2}, we report the variation of the quasiparticle energies of the four orbitals as functions of the internuclear distance $\RHH$. One can easily diagnose two problematic regions showing obvious discontinuities around $\RHH = \SI{1.2}{\angstrom}$ for the LUMO$+1$ ($p = 3$) and $\RHH = \SI{0.5}{\angstrom}$ for the LUMO$+2$ ($p = 4$). @@ -353,7 +353,7 @@ Inspection of their corresponding eigenvectors reveals that the $(\Ne+1)$-electr By construction, the quasiparticle solution diabatically follows the reference determinant $\ket*{1\Bar{1}3}$ through the avoided crossing (thick lines in Fig.~\ref{fig:H2_zoom}) which is precisely the origin of the energetic discontinuity. A similar scenario is at play in the region around $\RHH = \SI{0.5}{\angstrom}$ for the LUMO$+2$ (right panel of Fig.~\ref{fig:H2_zoom}) but it now involves three solutions ($s = 5$, $s = 6$, and $s = 7$). -The electronic configurations of the Slater determinant involved are the $\ket*{1\Bar{1}4}$ reference determinant as well as two external determinants of configuration $\ket*{1\Bar{?}?}$ and $\ket*{1\Bar{?}?}$. +The electronic configurations of the Slater determinant involved are the $\ket*{1\Bar{1}4}$ reference determinant as well as two external determinants of configuration \titou{$\ket*{1\Bar{?}?}$} and \titou{$\ket*{1\Bar{?}?}$}. These forms two avoided crossings in rapid successions, which create two discontinuities in the energy surface (see Fig.~\ref{fig:H2}). In this region, although the ground-state wave function is well described by the $\Ne$-electron HF determinant, a situation that can be safely labeled as single-reference, one can see that the $(\Ne+1)$-electron state involves three Slater determinants and can then be labeled as a multi-reference (or strongly-correlated) situation with near-degenerate electronic configurations. Therefore, one can conclude that this downfall of $GW$ is a key signature of strong correlation in the $(\Ne\pm1)$-electron states that yields a significant redistribution of weights amongst electronic configurations. @@ -380,7 +380,8 @@ Therefore, one can conclude that this downfall of $GW$ is a key signature of str % FIGURE 4 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \begin{figure} - \includegraphics[width=\linewidth]{fig4} +% \includegraphics[width=\linewidth]{fig4a} + \includegraphics[width=\linewidth]{fig4b} \caption{ \label{fig:H2reg} Difference between regularized and non-regularized quasiparticle energies $\reps{p}{\GW} - \eps{p}{\GW}$ as functions of the internuclear distance $\RHH$ (in \si{\angstrom}) of \ce{H2} at the {\GOWO}@HF/6-31G level. @@ -439,6 +440,7 @@ However, this diabatization is more or less accurate depending on the value of $ For $\eta = 10$, the value is clearly too large inducing a large difference between the two sets of quasiparticle energies (purple curves). For $\eta = 0.1$, we have the opposite scenario where the value is too small and some irregularities remain (green curves). We have found that $\eta = 1$ is a good compromise that does not alter significantly the quasiparticle energies while providing a smooth transition between the two solutions. +\titou{This values could be further refined for specfici applications.} To further evidence this, Fig.~\ref{fig:H2reg} reports the difference between regularized and non-regularized quasiparticle energies as functions of $\RHH$.