saving work

Manuscript/SRGGW.tex

@@ -186,7 +186,7 @@ where $\eta$ is a positive infinitesimal and the screened two-electron integrals
 \end{equation}
 where $\bX$ and $\bY$ are the components of the eigenvectors of the particle-hole direct RPA problem defined as 
 \begin{equation}
-  \mqty( \bA & \bB \\ -\bA & -\bB ) \mqty( \bX \\ \bY ) =  \mqty( \bX \\ \bY ) \boldsymbol{\Omega},
+  \mqty( \bA & \bB \\ -\bA & -\bB ) \mqty( \bX \\ \bY ) =  \mqty( \bX \\ \bY ) \mqty( \boldsymbol{\Omega} & \bO \\ \bO & -\boldsymbol{\Omega} ),
 \end{equation}
 with
 \begin{subequations}
@@ -196,7 +196,7 @@ with
   B_{ia,jb} & = \eri{ij}{ab}.
 \end{align}
 \end{subequations}
-The diagonal matrix $\boldsymbol{\Omega}$ contains the eigenvalues and its elements $\Omega_\nu$ appear in Eq.~\eqref{eq:GW_selfenergy}.
+The diagonal matrix $\boldsymbol{\Omega}$ contains the positive eigenvalues and its elements $\Omega_\nu$ appear in Eq.~\eqref{eq:GW_selfenergy}.
 \titou{The TDA case is discussed in Appendix \ref{sec:nonTDA}.}
 Throughout the manuscript, the indices $p,q,r,s$ are general orbitals while $i,j,k,l$ and $a,b,c,d$ refers to occupied and virtual orbitals, respectively.
 The indices $\mu$ and $\nu$ are composite indices, \eg $\nu=(ia)$, referring to neutral excitations.
@@ -253,51 +253,52 @@ Various other regularisers are possible and in particular one of us has shown th
 But it would be more rigorous, and more instructive, to obtain this regulariser from first principles by applying the SRG formalism to many-body perturbation theory.
 This is the aim of the rest of this work.
 
+\ANT{The matrix element expressions should be changed}
 Applying the SRG to $GW$ could gradually remove the coupling between the quasi-particle and the satellites resulting in a renormalized quasi-particle.
 However, to do so one needs to identify the coupling terms in Eq.~\eqref{eq:quasipart_eq}, which is not straightforward.
 The way around this problem is to transform Eq.~\eqref{eq:quasipart_eq} to its upfolded version and the coupling terms will elegantly appear in the process.
-The upfolded $GW$ quasi-particle equation is
+The upfolded $GW$ quasi-particle equation is \cite{Bintrim_2021,Tolle_2023}
 \begin{equation}
   \label{eq:GWlin}
   \begin{pmatrix}
-    \bF				&	\bV^{\text{2h1p}}		&	\bV^{\text{2p1h}}	\\
+    \bF				&	\bW^{\text{2h1p}}		&	\bW^{\text{2p1h}}	\\
-    \T{(\bV^{\text{2h1p}})}	&	\bC^{\text{2h1p}}		&	\bO					\\
+    \T{(\bW^{\text{2h1p}})}	&	\bC^{\text{2h1p}}		&	\bO					\\
-    \T{(\bV^{\text{2p1h}})}	&	\bO				&	\bC^{\text{2p1h}}	\\
+    \T{(\bW^{\text{2p1h}})}	&	\bO				&	\bC^{\text{2p1h}}	\\
   \end{pmatrix}
   \begin{pmatrix}
-    \bX	\\
+    \bZ^{\text{1h/1p}}	\\
-    \bY^{\text{2h1p}}	\\
+    \bZ^{\text{2h1p}}	\\
-    \bY^{\text{2p1h}}	\\
+    \bZ^{\text{2p1h}}	\\
   \end{pmatrix}
   =
   \begin{pmatrix}
-    \bX	\\
+    \bZ^{\text{1h/1p}}	\\
-    \bY^{\text{2h1p}}	\\
+    \bZ^{\text{2h1p}}	\\
-    \bY^{\text{2p1h}}	\\
+    \bZ^{\text{2p1h}}	\\
   \end{pmatrix}
   \boldsymbol{\epsilon},
 \end{equation}
 where $\boldsymbol{\epsilon}$ is a diagonal matrix collecting the quasi-particle and satellite energies, the 2h1p and 2p1h matrix elements are
 \begin{subequations}
   \begin{align}
-    C^\text{2h1p}_{ija,klc} & = \qty[ \qty( \epsilon_{i} + \epsilon_{j} - \epsilon_{a} ) \delta_{jl} \delta_{ac} - \eri{jc}{al} ] \delta_{ik} ,
+    C^\text{2h1p}_{i\nu,j\mu} & =  \left(\epsilon_i - \Omega_\nu\right)\delta_{ij}\delta_{\nu\mu},
     \\
-    C^\text{2p1h}_{iab,kcd} & = \qty[ \qty( \epsilon_{a} + \epsilon_{b} - \epsilon_{i} ) \delta_{ik} \delta_{ac} + \eri{ak}{ic} ] \delta_{bd},
+    C^\text{2p1h}_{a\nu,b\mu} & = \left(\epsilon_a + \Omega_\nu\right)\delta_{ab}\delta_{\nu\mu},
   \end{align}
 \end{subequations}
-and the corresponding coupling blocks read
+and the corresponding coupling blocks read [see Eq.~(\ref{eq:GW_sERI})]
 \begin{align}
-  V^\text{2h1p}_{p,klc} & = \eri{pc}{kl},
+  W^\text{2h1p}_{p,i\nu} & = W_{p,i\nu},
   &
-    V^\text{2p1h}_{p,kcd} & = \eri{pk}{dc}.
+    W^\text{2p1h}_{p,a\nu} & = W_{p,a\nu}.
 \end{align}
 The usual $GW$ non-linear equation can be obtained by applying L\"odwin partitioning technique \cite{Lowdin_1963} to Eq.~\eqref{eq:GWlin} which gives the following expression for the self-energy \cite{Bintrim_2021} 
 \begin{equation}
 \begin{split}
   \bSig(\omega) 
-  & = \bV^{\hhp} \qty(\omega \mathbb{1} - \bC^{\hhp})^{-1}  \T{(\bV^{\hhp})}
+  & = \bW^{\hhp} \qty(\omega \bI - \bC^{\hhp})^{-1}  \T{(\bW^{\hhp})}
   \\
-  & + \bV^{\pph} \qty(\omega \mathbb{1} - \bC^{\pph})^{-1}  \T{(\bV^{\pph})}, 
+  & + \bW^{\pph} \qty(\omega \bI - \bC^{\pph})^{-1}  \T{(\bW^{\pph})}, 
 \end{split}
 \end{equation}
 which can be further developed to give exactly Eq.~(\ref{eq:GW_selfenergy}).
@@ -339,7 +340,7 @@ where $\boldsymbol{\eta}(s)$, the flow generator, is defined as
 \boldsymbol{\eta}(s) = \dv{\bU(s)}{s}	\bU^\dag(s)	= - \boldsymbol{\eta}^\dag(s).
 \end{equation}
 
-To solve \titou{this equation} at a lower cost than the one associated with the diagonalization of the initial Hamiltonian, one must introduce an approximate form for $\boldsymbol{\eta}(s)$.
+To solve the flow equation at a lower cost than the one associated with the diagonalization of the initial Hamiltonian, one must introduce an approximate form for $\boldsymbol{\eta}(s)$.
 In this work, we consider Wegner's canonical generator \cite{Wegner_1994}
 \begin{equation}
   \boldsymbol{\eta}^\text{W}(s) = \comm{\bH^\text{d}(s)}{\bH(s)} = \comm{\bH^\text{d}(s)}{\bH^\text{od}(s)},
@@ -387,9 +388,9 @@ As hinted at the end of Sec.~\ref{sec:gw}, the diagonal and off-diagonal parts a
   & \\
     \bH^\text{od}(s)  &= 
   \begin{pmatrix}
-    \bO			                         	&     \bV^{\text{2h1p}} &     \bV^{\text{2p1h}} \\
+    \bO			                         	&     \bW^{\text{2h1p}} &     \bW^{\text{2p1h}} \\
-    (\bV^{\text{2h1p}})^{\mathrm{T}}	&	\bO	   	        &	\bO                    	\\
+    (\bW^{\text{2h1p}})^{\mathrm{T}}	&	\bO	   	        &	\bO                    	\\
-    (\bV^{\text{2p1h}})^{\mathrm{T}}	&	\bO		        &	\bO	                \\
+    (\bW^{\text{2p1h}})^{\mathrm{T}}	&	\bO		        &	\bO	                \\
   \end{pmatrix}
 \end{align}
 where we have omitted the $s$ dependence of the matrix elements for the sake of brevity.
@@ -408,8 +409,8 @@ Then, the aim is to solve order by order the flow equation [see Eq.~\eqref{eq:fl
                  \bO & \bO 
                \end{pmatrix} &
                                \bHod{1}(0) &=   \begin{pmatrix}
-                                                  \bO & \bV{}{} \\
+                                                  \bO & \bW{}{} \\
-                                                  \bV{}{\dagger} & \bO \notag
+                                                  \bW^{\dagger} & \bO \notag
                                                 \end{pmatrix}   \notag 
 \end{align}
 where we have defined the matrices $\bC$ and $\bV$ that collects the 2h1p and 2p1h channels for the sake of conciseness.
@@ -417,7 +418,7 @@ Once the analytical low-order perturbative expansions are known they can be inse
 In particular, in this manuscript, the focus will be on the second-order renormalized quasi-particle equation.
 
 %///////////////////////////%
-\subsection{Zeroth-order matrix elements}
+\subsection{Zeroth-order matrix elements \ANT{This subsection is false, I'll rewrite it soon}}
 %///////////////////////////%
 
 There is only one zeroth order term in the right-hand side of the flow equation
@@ -429,32 +430,31 @@ and performing the block matrix products gives the following system of equations
   \begin{align}
     \dv{\bF^{(0)}}{s}  &= \bO \\
     \dv{\bC^{(0)}}{s}  &= \bO \\
-    \dv{\bV^{(0),\dagger}}{s} &= 2 \bC^{(0)}\bV^{(0),\dagger}\bF^{(0)} - \bV^{(0),\dagger}(\bF^{(0)})^2 - (\bC^{(0)})^2\bV^{(0),\dagger} \\
+    \dv{\bW^{(0),\dagger}}{s} &= 2 \bC^{(0)}\bW^{(0),\dagger}\bF^{(0)} - \bW^{(0),\dagger}(\bF^{(0)})^2 - (\bC^{(0)})^2\bW^{(0),\dagger} \\
-    \dv{\bV^{(0)}}{s} &= 2 \bF^{(0)}\bV^{(0)}\bC^{(0)} - (\bF^{(0)})^2\bV^{(0)} - \bV^{(0)}(\bC^{(0)})^2 
+    \dv{\bW^{(0)}}{s} &= 2 \bF^{(0)}\bW^{(0)}\bC^{(0)} - (\bF^{(0)})^2\bW^{(0)} - \bW^{(0)}(\bC^{(0)})^2 
   \end{align}
 \end{subequations}
-where the $s$ dependence of $\bV^{(0)}$ and $\bV^{(0),\dagger}$  has been dropped in the last two equations.
+where the $s$ dependence of $\bW^{(0)}$ and $\bW^{(0),\dagger}$  has been dropped in the last two equations.
 $\bF^{(0)}$ and $\bC^{(0)}$ do not depend on $s$ as a consequence of the first two equations.
-The last equation can be solved by introducing $\bU$ the matrix that diagonalizes $\bC^{(0)} = \bU \bD^{(0)} \bU^{-1}$ such that the differential equation for $\bV^{(0)}$ becomes
+% The last equation can be solved by introducing $\bU$ the matrix that diagonalizes $\bC^{(0)} = \bU \bD^{(0)} \bU^{-1}$ such that the differential equation for $\bV^{(0)}$ becomes
-\begin{equation}
+% \begin{equation}
-  \label{eq:eqdiffW0}
+%   \label{eq:eqdiffW0}
-  \dv{\bW^{(0)}}{s} = 2 \bF^{(0)}\bW^{(0)} \bD^{(0)} - (\bF^{(0)})^2\bW^{(0)} - \bW^{(0)} (\bD^{(0)})^2
+%   \dv{\bW^{(0)}}{s} = 2 \bF^{(0)}\bW^{(0)} \bD^{(0)} - (\bF^{(0)})^2\bW^{(0)} - \bW^{(0)} (\bD^{(0)})^2
-\end{equation}
+% \end{equation}
-where $\bW^{(0)}(s)= \bV^{(0)}(s) \bU$.
+% where $\bW^{(0)}(s)= \bV^{(0)}(s) \bU$.
-The matrix elements of $\bU$ and $\bD^{(0)}$ are
+% The matrix elements of $\bU$ and $\bD^{(0)}$ are
-\begin{align}
+% \begin{align}
-  U_{p\nu,q\mu} &= \delta_{pq} \bX_{\nu\mu} \\
+%   U_{p\nu,q\mu} &= \delta_{pq} \bX_{\nu\mu} \\
-  D_{p\nu,q\mu}^{(0)} &= \left(\epsilon_p + \sgn(\epsilon_p-\epsilon_\text{F})\Omega_\nu\right)\delta_{pq}\delta_{\nu\mu}
+%   D_{p\nu,q\mu}^{(0)} &= \left(\epsilon_p + \sgn(\epsilon_p-\epsilon_\text{F})\Omega_\nu\right)\delta_{pq}\delta_{\nu\mu}
-\end{align}
+% \end{align}
-where $\epsilon_\text{F}$ is the Fermi level.
+% where $\epsilon_\text{F}$ is the Fermi level.
-Note that the matrix $\bU$ is also used in the downfolding process of Eq.~\eqref{eq:GWlin}. \cite{Bintrim_2021}
+% Note that the matrix $\bU$ is also used in the downfolding process of Eq.~\eqref{eq:GWlin}. \cite{Bintrim_2021}
 
-Thanks to the diagonal structure of $\bF^{(0)}$ and $\bD^{(0)}$, Eq.~\eqref{eq:eqdiffW0} can be easily solved and give 
+Thanks to the diagonal structure of $\bF^{(0)}$ and $\bC^{(0)}$, the last equation can be easily solved and give 
 \begin{equation}
   W_{p,q\nu}^{(0)}(s) = W_{p,q\nu}^{(0)}(0) e^{- (F_{pp}^{(0)} - D_{q\nu,q\nu}^{(0)})^2 s}
 \end{equation}
-Due to the initial conditions  $\bV^{(0)}(0) = \bO$, we have $\bW^{(0)}(s)=\bO$ and therefore $\bV^{(0)}(s)=\bO=\bV^{(0)}(0) $.
+The initial condition $\bW^{(0)}(0) = \bO$ implies $\bW^{(0)}(s)=\bO$ and therefore the zeroth order Hamiltonian is
-Therefore, the zeroth order Hamiltonian is
 \begin{equation}
   \bH^{(0)}(s) = \bH^{(0)}(0),
 \end{equation}
@@ -468,7 +468,7 @@ Knowing that $\bHod{0}(s)=\bO$, the first order flow equation is
 \begin{equation}
   \dv{\bH^{(1)}}{s} = \comm{\comm{\bHd{0}}{\bHod{1}}}{\bHd{0}}
 \end{equation}
-which gives the same system of equations as in the previous subsection except that $\bV^{(0)}$ and $\bV^{(0),\dagger}$ should be replaced by $\bV^{(1)}$ and $\bV^{(1),\dagger}$.
+which gives the same system of equations as in the previous subsection except that $\bW^{(0)}$ and $\bW^{(0),\dagger}$ should be replaced by $\bW^{(1)}$ and $\bW^{(1),\dagger}$.
 
 Once again the two first equations are easily solved
 \begin{align}
@@ -496,21 +496,21 @@ with
 \begin{align}
   \widetilde{\bF}(s) &= \bF^{(0)}+\bF^{(2)}(s)\\
   \label{eq:srg_sigma}
-  \widetilde{\bSig}(\omega; s) &= \bV^{(1)}(s) \left(\omega \mathbb{1} - \bC^{(0)}\right)^{-1}  (\bV^{(1)}(s))^{\dagger}
+  \widetilde{\bSig}(\omega; s) &= \bV^{(1)}(s) \left(\omega \bI - \bC^{(0)}\right)^{-1}  (\bV^{(1)}(s))^{\dagger}
 \end{align}
 
 As can be readily seen above, $\bF^{(2)}$ is the only second-order block of the effective Hamiltonian contributing to the second-order SRG quasi-particle equation.
 Collecting every second-order terms in the flow equation and performing the block matrix products results in the following differential equation for $\bF^{(2)}$
-\begin{equation}
+\begin{multline}
   \label{eq:diffeqF2}
-  \dv{\bF^{(2)}}{s}  = \bF^{(0)}\bV^{(1)}\bV^{(1),\dagger} + \bV^{(1)}\bV^{(1),\dagger}\bF^{(0)} - 2 \bV^{(1)}\bC^{(0)}\bV^{(1),\dagger} .
+  \dv{\bF^{(2)}}{s}  = \bF^{(0)}\bW^{(1)}\bW^{(1),\dagger} + \bW^{(1)}\bW^{(1),\dagger}\bF^{(0)} \\ - 2 \bW^{(1)}\bC^{(0)}\bW^{(1),\dagger} .
-\end{equation}
+\end{multline}
 This can be solved by simple integration along with the initial condition $\bF^{(2)}=\bO$ to give
 \begin{multline}
-  F_{pq}^{(2)}(s) =  \sum_{r\mu} \frac{\Delta_{pr\mu}+ \Delta_{qr\mu}}{\Delta_{pr\mu}^2 + \Delta_{qr\mu}^2} W_{p,r\mu} W^{\dagger}_{r\mu,q} \times \\
+  F_{pq}^{(2)}(s) =  \sum_{r\nu} \frac{\Delta_{pr\nu}+ \Delta_{qr\nu}}{\Delta_{pr\nu}^2 + \Delta_{qr\nu}^2} W_{p,r\nu} W^{\dagger}_{r\nu,q} \times \\
-  \left(1 - e^{-(\Delta_{pr\mu}^2 + \Delta_{qr\mu}^2) s}\right),
+  \left(1 - e^{-(\Delta_{pr\nu}^2 + \Delta_{qr\nu}^2) s}\right),
 \end{multline}
-with $\Delta_{prv} = \epsilon_p - \epsilon_r - \sgn(\epsilon_r-\epsilon_F)\Omega_v$.
+with $\Delta_{pr\nu} = \epsilon_p - \epsilon_r - \sgn(\epsilon_r-\epsilon_F)\Omega_\nu$.
 
 At $s=0$, this second-order correction is null while for $s\to\infty$ it tends towards the following static limit
 \begin{equation}
@@ -524,7 +524,7 @@ y
 Interestingly, the static limit, \ie $s\to\infty$ limit, of Eq.~\eqref{eq:GW_renorm} defines an alternative qs$GW$ approximation to the one defined by Eq.~\eqref{eq:sym_qsgw} which matrix elements read as
 \begin{equation}
   \label{eq:sym_qsGW}
-    \Sigma_{pq}^{\text{qs}}(\eta) = \frac{1}{2} \sum_{r\mu}  \left( \frac{\Delta_{pr\mu}}{\Delta_{pr\mu}^2 + \eta^2} +\frac{\Delta_{qr\mu}}{\Delta_{qr\mu}^2 + \eta^2}  \right) W_{p,r\mu} W^{\dagger}_{r\mu,q}.
+    \Sigma_{pq}^{\text{qs}}(\eta) = \frac{1}{2} \sum_{r\nu}  \left( \frac{\Delta_{pr\nu}}{\Delta_{pr\nu}^2 + \eta^2} +\frac{\Delta_{qr\nu}}{\Delta_{qr\nu}^2 + \eta^2}  \right) W_{p,r\nu} W^{\dagger}_{r\nu,q}.
 \end{equation}
 Yet, both approximation are closely related as they share the same diagonal terms when $\eta=0$.
 Also, note that the SRG static approximation is naturally Hermitian as opposed to the symmetrized case where it is enforced by symmetrization.
@@ -536,7 +536,7 @@ Therefore, we will define the SRG-qs$GW$ static effective Hamiltonian as
 
 \begin{multline}
   \label{eq:SRG_qsGW}
-  \Sigma_{pq}^{\text{SRG}}(s) = \sum_{r\mu} \frac{\Delta_{pr\mu}+ \Delta_{qr\mu}}{\Delta_{pr\mu}^2 + \Delta_{qr\mu}^2} W_{p,r\mu} W_{q,r\mu} \times \\ \left(1 - e^{-(\Delta_{pr\mu}^2 + \Delta_{qr\mu}^2) s}\right)  
+  \Sigma_{pq}^{\text{SRG}}(s) = \sum_{r\nu} \frac{\Delta_{pr\nu}+ \Delta_{qr\nu}}{\Delta_{pr\nu}^2 + \Delta_{qr\nu}^2} W_{p,r\nu} W_{q,r\nu} \times \\ \left(1 - e^{-(\Delta_{pr\nu}^2 + \Delta_{qr\nu}^2) s}\right)  
 \end{multline}
 which depends on one regularising parameter $s$ analogously to $\eta$ in the usual case.
 The fact that the $s\to\infty$ static limit does not always converge when used in a qs$GW$ calculation could have been predicted because in this limit even the intruder states have been included in $\tilde{\bF}$.
@@ -612,7 +612,7 @@ This behavior as a function of $s$ can be \ant{approximately} streamlined by app
 Through second order, the principal IP is
 \begin{equation}
   \label{eq:2nd_order_IP}
-  I_k(2) = - \epsilon_k - \sum_{i\nu} \frac{W_{k,i\nu}^2}{\epsilon_k - \epsilon_i + \Omega_\nu} - \sum_{a\nu} \frac{W_{k,a\nu}^2}{\epsilon_k - \epsilon_a - \Omega_\nu} + \order{3}
+  I_k = - \epsilon_k - \sum_{i\nu} \frac{W_{k,i\nu}^2}{\epsilon_k - \epsilon_i + \Omega_\nu} - \sum_{a\nu} \frac{W_{k,a\nu}^2}{\epsilon_k - \epsilon_a - \Omega_\nu} + \order{3}
 \end{equation}
 where $k$ is the index of the highest molecular MO (HOMO).
 The first term is the zeroth order IP and the two following terms come from the 2h1p and 2p1h coupling, respectively.
@@ -695,7 +695,7 @@ The matrix elements of the $GW$ self-energy within the TDA are the same as in Eq
   \label{eq:GWnonTDA_sERI}
   W_{p,q\nu} = \sum_{ia}\eri{pi}{qa}\qty( \bX_{\nu})_{ia},
 \end{equation}
-where $\bX$ ais the eigenvector matrix of the TDA particle-hole dRPA problem defined as
+where $\bX$ is the eigenvector matrix of the TDA particle-hole dRPA problem defined as
 \begin{equation}
   \label{eq:full_dRPA}
   \bA \bX = \bX \boldsymbol{\Omega}
@@ -704,7 +704,7 @@ with
 \begin{align}
   A_{ia,jb} &= (\epsilon_i - \epsilon_a) \delta_{ij}\delta_{ab} + \eri{ib}{aj} \\
 \end{align}
-$\boldsymbol{\Omega}$ is the diagonal matrix of eigenvalues. Note that $\boldsymbol{\Omega}$ has the same size as without the TDA because in RPA we consider only the positive excitations of the full dRPA problem.
+$\boldsymbol{\Omega}$ is the diagonal matrix of eigenvalues.
 
 Defining an unfold version of this equation that does not require a diagonalization of the RPA problem before unfolding is a tricky task (see supplementary material of Ref.~\onlinecite{Bintrim_2021}).
 However, because we will eventually downfold again the upfolded matrix, we can use the following matrix \cite{Tolle_2022}