loos/ufGW
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

saving work

Browse Source
This commit is contained in:
Pierre-Francois Loos 2022-02-22 15:05:56 +01:00
parent 82a0c98534
commit 51c015560a
9 changed files with 159 additions and 65 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

BIN
Manuscript/fig1.pdf Normal file
View File

Binary file not shown.

BIN
Manuscript/fig1a.pdf
View File

Binary file not shown.

BIN
Manuscript/fig1b.pdf
View File

Binary file not shown.

BIN
Manuscript/fig2a.pdf Normal file
View File

Binary file not shown.

BIN
Manuscript/fig2b.pdf Normal file
View File

Binary file not shown.

BIN
Manuscript/fig3a.pdf Normal file
View File

Binary file not shown.

BIN
Manuscript/fig3b.pdf Normal file
View File

Binary file not shown.

70
Manuscript/ufGW.bib
View File

@ -1,13 +1,62 @@
%% This BibTeX bibliography file was created using BibDesk.
%% https://bibdesk.sourceforge.io/
%% http://bibdesk.sourceforge.net/
%% Created for Pierre-Francois Loos at 2022-02-21 21:38:32 +0100
%% Created for Pierre-Francois Loos at 2022-02-22 14:38:25 +0100
%% Saved with string encoding Unicode (UTF-8)
@misc{Loos_2022,
archiveprefix = {arXiv},
author = {Pierre-Fran{\c c}ois Loos and Pina Romaniello},
date-added = {2022-02-22 14:34:24 +0100},
date-modified = {2022-02-22 14:34:28 +0100},
eprint = {2202.07936},
primaryclass = {physics.chem-ph},
title = {Static and Dynamic Bethe-Salpeter Equations in the $T$-Matrix Approximation},
year = {2022}}
@article{Backhouse_2020b,
author = {Backhouse, Oliver J. and Booth, George H.},
date-added = {2022-02-22 13:31:29 +0100},
date-modified = {2022-02-22 13:32:27 +0100},
doi = {10.1021/acs.jctc.0c00701},
journal = {J. Chem. Theory Comput.},
number = {10},
pages = {6294-6304},
title = {Efficient Excitations and Spectra within a Perturbative Renormalization Approach},
volume = {16},
year = {2020},
bdsk-url-1 = {https://doi.org/10.1021/acs.jctc.0c00701}}
@article{Backhouse_2020a,
author = {Backhouse, Oliver J. and Nusspickel, Max and Booth, George H.},
date-added = {2022-02-22 13:31:26 +0100},
date-modified = {2022-02-22 13:32:17 +0100},
doi = {10.1021/acs.jctc.9b01182},
journal = {J. Chem. Theory Comput.},
number = {2},
pages = {1090-1104},
title = {Wave Function Perspective and Efficient Truncation of Renormalized Second-Order Perturbation Theory},
volume = {16},
year = {2020},
bdsk-url-1 = {https://doi.org/10.1021/acs.jctc.9b01182}}
@article{Backhouse_2021,
author = {Backhouse, Oliver J. and Santana-Bonilla, Alejandro and Booth, George H.},
date-added = {2022-02-22 13:31:19 +0100},
date-modified = {2022-02-22 13:32:57 +0100},
doi = {10.1021/acs.jpclett.1c02383},
journal = {J. Phys. Chem. Lett.},
number = {31},
pages = {7650-7658},
title = {Scalable and Predictive Spectra of Correlated Molecules with Moment Truncated Iterated Perturbation Theory},
volume = {12},
year = {2021},
bdsk-url-1 = {https://doi.org/10.1021/acs.jpclett.1c02383}}
@article{Shee_2021,
author = {Shee, James and Loipersberger, Matthias and Rettig, Adam and Lee, Joonho and Head-Gordon, Martin},
date-added = {2022-02-21 21:37:57 +0100},
@ -25,9 +74,10 @@
abstract = {Low-order scaling GW implementations for molecules are usually restricted to approximations with diagonal self-energy. Here, we present an all-electron implementation of quasiparticle self-consistent GW for molecular systems. We use an efficient algorithm for the evaluation of the self-energy in imaginary time, from which a static non-local exchange-correlation potential is calculated via analytical continuation. By using a direct inversion of iterative subspace method, fast and stable convergence is achieved for almost all molecules in the GW100 database. Exceptions are systems which are associated with a breakdown of the single quasiparticle picture in the valence region. The implementation is proven to be starting point independent and good agreement of QP energies with other codes is observed. We demonstrate the computational efficiency of the new implementation by calculating the quasiparticle spectrum of a DNA oligomer with 1,220 electrons using a basis of 6,300 atomic orbitals in less than 4 days on a single compute node with 16 cores. We use then our implementation to study the dependence of quasiparticle energies of DNA oligomers consisting of adenine-thymine pairs on the oligomer size. The first ionization potential in vacuum decreases by nearly 1 electron volt and the electron affinity increases by 0.4 eV going from the smallest to the largest considered oligomer. This shows that the DNA environment stabilizes the hole/electron resulting from photoexcitation/photoattachment. Upon inclusion of the aqueous environment via a polarizable continuum model, the differences between the ionization potentials reduce to 130 meV, demonstrating that the solvent effectively compensates for the stabilizing effect of the DNA environment. The electron affinities of the different oligomers are almost identical in the aqueous environment.},
author = {F{\"o}rster, Arno and Visscher, Lucas},
date-added = {2022-02-21 21:12:20 +0100},
date-modified = {2022-02-21 21:12:37 +0100},
date-modified = {2022-02-22 14:37:55 +0100},
doi = {10.3389/fchem.2021.736591},
journal = {Front. Chem.},
pages = {736591},
title = {Low-Order Scaling Quasiparticle Self-Consistent GW for Molecules},
volume = {9},
year = {2021},
@ -296,9 +346,10 @@
abstract = {We use the GW100 benchmark set to systematically judge the quality of several perturbation theories against high-level quantum chemistry methods. First of all, we revisit the reference CCSD(T) ionization potentials for this popular benchmark set and establish a revised set of CCSD(T) results. Then, for all of these 100 molecules, we calculate the HOMO energy within second and third-order perturbation theory (PT2 and PT3), and, GW as post-Hartree-Fock methods. We found GW to be the most accurate of these three approximations for the ionization potential, by far. Going beyond GW by adding more diagrams is a tedious and dangerous activity: We tried to complement GW with second-order exchange (SOX), with second-order screened exchange (SOSEX), with interacting electron-hole pairs (W<sub>TDHF</sub>), and with a GW density-matrix (γ<sup>GW</sup>). Only the γ<sup>GW</sup> result has a positive impact. Finally using an improved hybrid functional for the non-interacting Green's function, considering it as a cheap way to approximate self-consistency, the accuracy of the simplest GW approximation improves even more. We conclude that GW is a miracle: Its subtle balance makes GW both accurate and fast.},
author = {Bruneval, Fabien and Dattani, Nike and van Setten, Michiel J.},
date-added = {2022-01-26 15:11:16 +0100},
date-modified = {2022-01-26 15:11:31 +0100},
date-modified = {2022-02-22 14:38:24 +0100},
doi = {10.3389/fchem.2021.749779},
journal = {Front. Chem.},
pages = {749779},
title = {The GW Miracle in Many-Body Perturbation Theory for the Ionization Potential of Molecules},
volume = {9},
year = {2021},
@ -3993,14 +4044,12 @@
@article{Foerster_2011,
author = {Foerster,D. and Koval,P. and S{\'a}nchez-Portal,D.},
date-added = {2020-05-18 21:40:28 +0200},
date-modified = {2020-05-18 21:40:28 +0200},
date-modified = {2022-02-22 14:37:24 +0100},
doi = {10.1063/1.3624731},
eprint = {https://doi.org/10.1063/1.3624731},
journal = {J. Chem. Phys.},
number = {7},
pages = {074105},
title = {An O(N3) implementation of Hedin's GW approximation for molecules},
url = {https://doi.org/10.1063/1.3624731},
volume = {135},
year = {2011},
bdsk-url-1 = {https://doi.org/10.1063/1.3624731}}
@ -5010,16 +5059,11 @@
@article{Romaniello_2012,
author = {Romaniello, Pina and Bechstedt, Friedhelm and Reining, Lucia},
date-added = {2020-05-18 21:40:28 +0200},
date-modified = {2020-05-18 21:40:28 +0200},
date-modified = {2022-02-22 14:33:08 +0100},
doi = {10.1103/PhysRevB.85.155131},
file = {/Users/loos/Zotero/storage/PLHVA772/Romaniello_2012.pdf},
issn = {1098-0121, 1550-235X},
journal = {Phys. Rev. B},
language = {en},
month = apr,
number = {15},
pages = {155131},
shorttitle = {Beyond the {{G W}} Approximation},
title = {Beyond the {{G W}} Approximation: {{Combining}} Correlation Channels},
volume = {85},
year = {2012},

154
Manuscript/ufGW.tex
View File

@ -65,12 +65,13 @@
% orbital energies
\newcommand{\eps}[2]{\epsilon_{#1}^{#2}}
\newcommand{\reps}[2]{\Tilde{\epsilon}_{#1}^{#2}}
\newcommand{\Om}[2]{\Omega_{#1}^{#2}}
% Matrix elements
\newcommand{\Sig}[2]{\Sigma_{#1}^{#2}}
\newcommand{\SigC}[1]{\Sigma^\text{c}_{#1}}
\newcommand{\rSigC}[1]{\widetilde{\Sigma}^\text{c}_{#1}}
\newcommand{\rSigC}[1]{\Tilde{\Sigma}^\text{c}_{#1}}
\newcommand{\SigX}[1]{\Sigma^\text{x}_{#1}}
\newcommand{\SigXC}[1]{\Sigma^\text{xc}_{#1}}
\newcommand{\MO}[1]{\phi_{#1}}
@ -125,8 +126,8 @@
\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, these issues can be directly related to the intruder state problem.
A simple and efficient regularization procedure is proposed to avoid such issues.
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.
%\bigskip
%\begin{center}
% \boxed{\includegraphics[width=0.5\linewidth]{TOC}}
@ -140,7 +141,7 @@ A simple and efficient regularization procedure is proposed to avoid such issues
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
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}
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:
\begin{equation}
@ -171,7 +172,7 @@ 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 gap, total and correlation energies, as well as vertical 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}
@ -185,8 +186,9 @@ In the present article, via an upfolding process of the non-linear $GW$ equation
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.
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.\cite{Hybertsen_1986,Shishkin_2007,Blase_2011,Faber_2011,Rangel_2016,Gui_2018,Faleev_2004,vanSchilfgaarde_2006,Kotani_2007,Ke_2011,Kaplan_2016}
Moreover, we consider a restricted Hartree-Fock (HF) starting point but it can be straightforwardly extended to a Kohn-Sham starting point.
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.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Downfolding: The non-linear $GW$ problem}
@ -197,7 +199,7 @@ Within the {\GOWO} approximation, in order to obtain the quasiparticle energies
\label{eq:qp_eq}
\eps{p}{\HF} + \SigC{p}(\omega) - \omega = 0
\end{equation}
where $\eps{p}{\HF}$ is the $p$th HF orbital energy and the correlation part of the {\GOWO} self-energy is constituted by a hole (h) and a particle (p) term, as follows
where $\eps{p}{\HF}$ is the $p$th HF orbital energy and the correlation part of the {\GOWO} self-energy is constituted by a hole (h) and a particle (p) term as follows
\begin{equation}
\label{eq:SigC}
\SigC{p}(\omega)
@ -221,6 +223,7 @@ and
\ERI{pq}{ia} = \iint \MO{p}(\br_1) \MO{q}(\br_1) \frac{1}{\abs{\br_1 - \br_2}} \MO{i}(\br_2) \MO{a}(\br_2) d\br_1 \dbr_2
\end{equation}
are two-electron integrals over the HF (spatial) orbitals $\MO{p}(\br)$.
Because one must compute all the RPA eigenvalues and eigenvectors to construct the self-energy \eqref{eq:SigC}, the computational cost is $\order*{\Nocc^3 \Nvir^3} = \order*{\Norb^6}$, where $\Nocc$ and $\Nvir$ are the number of occupied and virtual orbitals, respectively, and $\Norb = \Nocc + \Nvir$ is the total number of orbitals.
As a non-linear equation, Eq.~\eqref{eq:qp_eq} has many solutions $\eps{p,s}{\GW}$ and their corresponding weights are given by the value of the following renormalization factor
\begin{equation}
@ -240,7 +243,7 @@ which physically shows that the mean-field solution of unit weight is ``scattere
\section{Upfolding: the linear $GW$ problem}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
The non-linear quasiparticle equation \eqref{eq:qp_eq} can be \textit{exactly} transformed into a larger linear problem via an upfolding process where the 2h1p and 2p1h sectors
are upfolded from the 1h and 1p sectors. \cite{Bintrim_2021a,Riva_2022}
are upfolded from the 1h and 1p sectors. \cite{Backhouse_2020a,Backhouse_2020b,Bintrim_2021a,Backhouse_2021,Riva_2022}
For each orbital $p$, this yields a linear eigenvalue problem of the form
\begin{equation}
\bH^{(p)} \cdot \bc{}{(p,s)} = \eps{p,s}{\GW} \bc{}{(p,s)}
@ -269,13 +272,9 @@ and the corresponding coupling blocks read
&
V^\text{2p1h}_{p,kcd} & = \sqrt{2} \ERI{pd}{kc}
\end{align}
The size of this eigenvalue problem is $1 + O^2 V + O V^2$ (where $O$ and $V$ are the number of occupied and virtual orbitals, respectively), and this eigenvalue problem has to be solved for each orbital that one wishes to correct.
The size of this eigenvalue problem is $1 + \Nocc^2 \Nvir + \Nocc \Nvir^2 = \order*{\Norb^3}$, and it has to be solved for each orbital that one wishes to correct.
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.
Because the renormalization factor \eqref{eq:Z} corresponds to the projection of the vector $\bc{}{(p,s)}$ onto the reference (or internal) space, the weight of a solution $(p,s)$ is given by the first coefficient of their corresponding eigenvector $\bc{}{(p,s)}$, \ie,
\begin{equation}
\label{eq:Z_proj}
Z_{p,s} = \qty[ c_{1}^{(p,s)} ]^{2}
\end{equation}
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}.
This can be further illustrated by expanding the secular equation associated with Eq.~\eqref{eq:Hp}
@ -290,8 +289,13 @@ and comparing it with Eq.~\eqref{eq:qp_eq} by setting
+ \bV{p}{\text{2p1h}} \cdot \qty(\omega \bI - \bC{}{\text{2p1h}} )^{-1} \cdot \T{\qty(\bV{p}{\text{2p1h}})}
\end{multline}
where $\bI$ is the identity matrix.
Because the renormalization factor \eqref{eq:Z} corresponds to the projection of the vector $\bc{}{(p,s)}$ onto the reference (or internal) space, the weight of a solution $(p,s)$ is given by the first coefficient of their corresponding eigenvector $\bc{}{(p,s)}$, \ie,
\begin{equation}
\label{eq:Z_proj}
Z_{p,s} = \qty[ c_{1}^{(p,s)} ]^{2}
\end{equation}
One can see this downfolding process as the construction of a frequency-dependent effective Hamiltonian where the internal space is composed by a single determinant of the 1h or 1p sector and the external (or outer) space by all the 2h1p and 2p1h configurations. \cite{Dvorak_2019a,Dvorak_2019b,Bintrim_2021a}
One can see this downfolding process as the construction of a frequency-dependent effective Hamiltonian where the internal space is composed by a single Slater determinant of the 1h or 1p sector and the external (or outer) space by all the 2h1p and 2p1h configurations. \cite{Dvorak_2019a,Dvorak_2019b,Bintrim_2021a}
The main mathematical difference between the two approaches is that, by diagonalizing Eq.~\eqref{eq:Hp}, one has directly access to the internal and external components of the eigenvectors associated with each quasiparticle and satellite, and not only their projection in the reference space as shown by Eq.~\eqref{eq:Z_proj}.
The element $\eps{p}{\HF}$ of $\bH^{(p)}$ [see Eq.~\eqref{eq:Hp}] corresponds to the (approximate) relative energy of the $(\Ne\pm1)$-electron reference determinant (compared to the $\Ne$-electron HF determinant) while the eigenvalues of the blocks $\bC{}{\text{2h1p}}$ and $\bC{}{\text{2p1h}}$, which are $\eps{i}{\HF} - \Om{m}{\RPA}$ and $\eps{a}{\HF} + \Om{m}{\RPA}$ respectively, provide an estimate of the relative energy of the 2h1p and 2p1h determinants.
@ -303,43 +307,87 @@ As we shall see below, discontinuities, which are ubiquitous in molecular system
\section{An illustrative example}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
In order to illustrate the appearance and the origin of these multiple solutions, we consider the hydrogen molecule in the 6-31G basis set which corresponds to a two-electron system with four spatial orbitals (one occupied and three virtuals).
This example was already considered in our previous work \cite{Veril_2018} but here we provide further insights on the origin of the appearances of these discontinuities.
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}.
We denote as $\ket*{1\Bar{1}}$ the $\Ne$-electron ground-state 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.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FIGURE 1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{figure*}
\includegraphics[width=0.47\linewidth]{fig1a}
\hspace{0.05\linewidth}
\includegraphics[width=0.47\linewidth]{fig1b}
\begin{figure}
\includegraphics[width=\linewidth]{fig1}
\caption{
\label{fig:H2}
Selection of quasiparticle and satellite energies $\eps{p,s}{\GW}$ (solid lines) and their renormalization factor $Z_{p,s}$ (dashed lines) as functions of the internuclear distance $\RHH$ (in \si{\angstrom}) for the HOMO$+2$ ($p=3$) and HOMO$+3$ ($p=4$) orbitals of \ce{H2} at the {\GOWO}@HF/6-31G level.
Quasiparticle energies $\eps{p}{\GW}$ as functions of the internuclear distance $\RHH$ (in \si{\angstrom}) of \ce{H2} at the {\GOWO}@HF/6-31G level.
}
\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FIGURE 2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{figure*}
\includegraphics[width=0.47\linewidth]{fig2a}
\hspace{0.05\linewidth}
\includegraphics[width=0.47\linewidth]{fig2b}
\caption{
\label{fig:H2_zoom}
Selection of quasiparticle and satellite energies $\eps{p,s}{\GW}$ (solid lines) and their renormalization factor $Z_{p,s}$ (dashed lines) as functions of the internuclear distance $\RHH$ (in \si{\angstrom}) for the LUMO$+1$ ($p=3$) and LUMO$+2$ ($p=4$) orbitals of \ce{H2} at the {\GOWO}@HF/6-31G level.
The quasiparticle solution (which corresponds to the solution with the largest weight) is represented as a thicker line.
}
\end{figure*}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Figure \ref{fig:H2} shows the evolution of the quasiparticle energy, energetically close-by satellites and their corresponding weights at the {\GOWO}@HF level as a function on the internuclear distance $\RHH$.
One can easily diagnose two problematic regions showing obvious discontinuities around $\RHH = \SI{1.2}{\angstrom}$ for the HOMO$+2$ ($p = 3$) and $\RHH = \SI{0.5}{\angstrom}$ for the HOMO$+3$ ($p = 4$).
In order to illustrate the appearance and the origin of these multiple solutions, we consider the hydrogen molecule in the 6-31G basis set which corresponds to a two-electron system with four spatial orbitals (one occupied and three virtuals).
This example was already considered in our previous work \cite{Veril_2018} but here we provide further insights on the origin of the appearances of these discontinuities.
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.
Let us first look more closely at the region around $\RHH = \SI{1.2}{\angstrom}$.
As one can see, an avoided crossing is formed between two solutions of the quasiparticle equation.
Inspection of their corresponding eigenvectors reveals that the $(\Ne+1)$-electron determinants involved are the reference 1p determinant $\ket*{1\Bar{1}3}$ and an excited $(\Ne+1)$-electron of configuration $\ket*{12\Bar{2}}$ which becomes lower in energy than the reference determinant for $\RHH > \SI{1.2}{\angstrom}$.
The algorithm diabatically follows the reference determinant $\ket*{1\Bar{1}3}$, while it should be adiabatically switching to the $\ket*{12\Bar{2}}$ determinant which becomes lower in energy for $\RHH > \SI{1.1}{\angstrom}$.
This is not the $\Ne$-electron situation where one has to check if it is multireference, but for the $\Ne\pm1$-electron situations.
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$).
As thoroughly explained in Ref.~\onlinecite{Veril_2018}, if one relies on the linearization of quasiparticle equation \eqref{eq:qp_eq} to compute the quasiparticle energies, \ie, $\eps{p}{\GW} = \eps{p}{\HF} + Z_{p} \SigC{p}(\eps{p}{\HF})$, these discontinuities are transformed into irregularities as the renormalization factor cancels out the singularities of the self-energy.
As similar scenario is at the play in the region around $\RHH = \SI{0.5}{\angstrom}$ but it now involves 3 electronic configurations: the $\ket*{1\Bar{1}4}$ reference determinant as well as two external determinants of configuration $\ket*{1\Bar{?}?}$ and $\ket*{1\Bar{?}?}$.
These forms two avoided crossings in rapid successions, which involves two discontinuities in the energy surface.
Figure \ref{fig:H2_zoom} shows the evolution of the quasiparticle energy, the energetically close-by satellites and their corresponding weights as functions of $\RHH$.
Let us first look more closely at the region around $\RHH = \SI{1.2}{\angstrom}$ involving the LUMO$+1$ (left panel of Fig.~\ref{fig:H2_zoom}).
As one can see, an avoided crossing is formed between two solutions of the quasiparticle equation ($s = 4$ and $s = 5$).
Inspection of their corresponding eigenvectors reveals that the $(\Ne+1)$-electron determinants principally involved are the reference 1p determinant $\ket*{1\Bar{1}3}$ and an excited $(\Ne+1)$-electron determinant of configuration $\ket*{12\Bar{2}}$ that becomes lower in energy than the reference determinant for $\RHH > \SI{1.2}{\angstrom}$.
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.
As 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{?}?}$.
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.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Introducing regularized $GW$ methods}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
One way to alleviate these issues and to massively improve the convergence properties of self-consistent $GW$ calculations is to resort to a regularization of the self-energy without altering too much the quasiparticle energies.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FIGURE 3
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{figure*}
\includegraphics[width=0.47\linewidth]{fig3a}
\hspace{0.05\linewidth}
\includegraphics[width=0.47\linewidth]{fig3b}
\caption{
\label{fig:H2reg_zoom}
Comparison between non-regularized (solid lines) and regularized (dashed lines) energies as functions of the internuclear distance $\RHH$ (in \si{\angstrom}) for the LUMO$+1$ ($p=3$) and LUMO$+2$ ($p=4$) orbitals of \ce{H2} at the {\GOWO}@HF/6-31G level.
The quasiparticle solution is represented as a thicker line.}
\end{figure*}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FIGURE 4
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{figure}
% \includegraphics[width=\linewidth]{fig4}
\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.
}
\end{figure}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
One way to alleviate the issues discussed above and to massively improve the convergence properties of self-consistent $GW$ calculations is to resort to a regularization of the self-energy without altering too much the quasiparticle energies.
From a general perspective, a regularized $GW$ self-energy reads
\begin{equation}
@ -352,9 +400,14 @@ From a general perspective, a regularized $GW$ self-energy reads
\end{equation}
where various choices for the ``regularizer'' $f_\eta$ are possible.
The main purpose of $f_\eta$ is to ensure that $\rSigC{p}(\omega;\eta)$ remains finite even if one of the denominators goes to zero.
The regularized solutions $\reps{p,s}{\GW}$ are then obtained by solving the following regularized quasiparticle equation:
\begin{equation}
\label{eq:rqp_eq}
\eps{p}{\HF} + \rSigC{p}(\omega;\eta) - \omega = 0
\end{equation}
The most common and well-established way of regularizing $\Sigma$ is via the simple regularizer $f_\eta(\Delta) = (\Delta \pm \eta)^{-1}$, a strategy somehow related to the imaginary shift used in multiconfigurational perturbation theory. \cite{Forsberg_1997}
Other choices are legitimate like the regularizers considered by Head-Gordon and coworkers within second-order M{\o}ller-Plesset theory. \cite{Lee_2018a,Shee_2021}
Other choices are legitimate like the regularizers considered by Head-Gordon and coworkers within orbital-optimized second-order M{\o}ller-Plesset theory. \cite{Lee_2018a,Shee_2021}
Our investigations have shown that the following regularizer
\begin{equation}
@ -367,26 +420,23 @@ Of course, by construction, we have
\lim_{\eta\to0} \rSigC{p}(\omega;\eta) = \SigC{p}(\omega)
\end{equation}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FIGURE 2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{figure*}
% \includegraphics[width=\linewidth]{fig1a}
% \includegraphics[width=\linewidth]{fig1b}
\caption{
\label{fig:H2_reg}
Quasiparticle energies (left), correlation part of the self-energy (center) and renormalization factor (right) as functions of the internuclear distance $\RHH$ (in \si{\angstrom}) for various orbitals of \ce{H2} at the {\GOWO}@HF/6-31G level.
}
\end{figure*}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Figure \ref{fig:H2reg_zoom} compares the non-regularized and regularized quasiparticle energies in the two regions of interest for $\eta = 0.1$, $1$, and $10$.
It clearly shows how the regularization of the $GW$ self-energy diabatically linked the two solutions to get rid of the discontinuities.
However, this diabatization is more or less accurate depending on the value of $\eta$.
For $\eta = 10$, the value is clearly too large which induces 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.
Figure \ref{fig:H2_reg} evidences how the regularization of the $GW$ self-energy diabatically linked the two solutions to get rid of the discontinuities.
To further evidence this, Fig.~\ref{fig:H2reg} reports the difference between regularized and non-regularized quasiparticle energies as functions of $\RHH$.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Concluding remarks}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%The $T$-matrix-based formalism as well as second-order Green's function (or second Born) scheme \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 problems.
In the present article, we have provided mathematical and physical explanations behind the appearance of multiple solutions and discontinuities in various physical quantities computed within the $GW$ approximation.
More precisely, we have evidenced that intruder states are the main cause behind these issues and that this downfall of $GW$ is a key signature of strong correlation.
A simple and efficient regularization procedure inspired by the similarity renormalization group has been proposed to remove these discontinuities without altering too much the quasiparticle energies.
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 problems.
We hope that these new physical insights and technical developments will broaden the applicability of Green's function methods in the molecular electronic structure community and beyond.
%%%%%%%%%%%%%%%%%%%%%%%%
\acknowledgements{