ok with appendix

This commit is contained in:
Pierre-Francois Loos 2022-09-21 11:07:57 +02:00
parent db7e8089c5
commit 9a5df4440e

343
g.tex
View File

@ -53,6 +53,7 @@
\newcommand{\cC}{\mathcal{C}}
\newcommand{\cD}{\mathcal{D}}
\newcommand{\cE}{\mathcal{E}}
\newcommand{\cI}{\mathcal{I}}
\newcommand{\EPT}{E_{\PT}}
\newcommand{\laPT}{\lambda_{\PT}}
@ -674,7 +675,7 @@ It is then a trivial matter to integrate out exactly the $n_k$ variables, leadin
\qty[ \prod_{k=0}^{p-1} \mel{ I_k }{ P_k \qty( H-E \Id )^{-1} P_k (-H)(1-P_k) }{ I_{k+1} } ]
\mel{ I_p }{ P_p \qty( H-E \Id)^{-1} P_p }{ I_N }
\end{multline}
As an illustration, Appendix \ref{app:A} reports the exact derivation of this formula in the case of a \titou{two-level} system.
As an illustration, Appendix \ref{app:A} reports the exact derivation of this formula in the case of a two-state system.
%----------------------------%
\subsubsection{Dyson equation}
@ -793,140 +794,116 @@ Thus, from a practical point of view, a trade-off has to be found between the \t
Let us consider the one-dimensional Hubbard Hamiltonian for a chain of $N$ sites
\be
\hat{H}= -t \sum_{\langle i j\rangle \sigma} \hat{a}^+_{i\sigma} \hat{a}_{j\sigma}
+ U \sum_i \hat{n}_{i\uparrow} \hat{n}_{i\downarrow}
\hat{H}= -t \sum_{\expval{ i j } \sigma} \hat{a}^+_{i\sigma} \hat{a}_{j\sigma}
+ U \sum_i \hat{n}_{i\uparrow} \hat{n}_{i\downarrow},
\ee
where $\langle i j\rangle$ denotes the summation over two neighboring sites,
$\hat{a}_{i\sigma} (\hat{a}_{i\sigma})$ is the fermionic creation (annihilation) operator of
an spin-$\sigma$ electron (with $\sigma$ = $\uparrow$ or $\downarrow$) on site $i$, $\hat{n}_{i\sigma} = \hat{a}^+_{i\sigma} \hat{a}_{i\sigma}$ the
number operator, $t$ the hopping amplitude and $U$ the on-site Coulomb repulsion.
We consider a chain with an even number of sites and open boundary conditions
at half-filling, that is, $N_{\uparrow}=N_{\downarrow}=N/2$.
In the site-representation, a general vector of the Hilbert space will be written as
where $\langle i j\rangle$ denotes the summation over two neighboring sites, $\hat{a}_{i\sigma}$ ($\hat{a}_{i\sigma}$) is the fermionic creation (annihilation) operator of an spin-$\sigma$ electron (with $\sigma$ = $\uparrow$ or $\downarrow$) on site $i$, $\hat{n}_{i\sigma} = \hat{a}^+_{i\sigma} \hat{a}_{i\sigma}$ the number operator, $t$ the hopping amplitude, and $U$ the on-site Coulomb repulsion.
We consider a chain with an even number of sites and open boundary conditions at half-filling, that is, $N_{\uparrow}=N_{\downarrow}=N/2$.
In the site representation, a general vector of the Hilbert space can be written as
\be
|n\rangle = |n_{1 \uparrow},...,n_{N \uparrow},n_{1 \downarrow},...,n_{N \downarrow}\rangle
\ket{n} = \ket{n_{1 \uparrow},\ldots,n_{N \uparrow},n_{1 \downarrow},\ldots,n_{N \downarrow}},
\ee
where $n_{i \sigma}=0,1$ is the number of electrons of spin $\sigma$ at site $i$.
where $n_{i \sigma}=0$ or $1$ is the number of electrons of spin $\sigma$ on site $i$.
For the one-dimensional Hubbard model with open boundary conditions, the components of the ground-state vector have the same sign (say, $c_i \ge 0$).
It is then possible to identify the guiding and trial vectors, that is, $|c^+\rangle=|c_T\rangle$.
As trial wave function we shall employ a generalization of the Gutzwiller wave function\cite{Gutzwiller_1963} written under the form
It is then possible to equate the guiding and trial vectors, \ie, $\ket{c^+} = \ket{c_\text{T}}$.
As trial wave function, we shall employ a generalization of the Gutzwiller wave function \cite{Gutzwiller_1963}
\be
\langle n|c_T\rangle= e^{-\alpha n_D(n)-\beta n_A(n)}
\braket{ n }{ c_\text{T} }= e^{-\alpha n_D(n)-\beta n_A(n)},
\ee
where $n_D(n)$ is the number of doubly occupied sites for the configuration $|n\rangle$ and $n_A(n)$
the number of nearest-neighbor antiparallel pairs defined as
where $n_D(n)$ is the number of doubly occupied sites for the configuration $\ket{n}$ and $n_A(n)$ the number of nearest-neighbor antiparallel pairs defined as
\be
n_A(n)= \sum_{\langle i j\rangle} n_{i\uparrow} n _{j\downarrow} \delta(n_{i\uparrow}-1) \delta(n_{j\downarrow}-1)
n_A(n)= \sum_{\expval{ i j }} n_{i\uparrow} n _{j\downarrow} \delta_{n_{i\uparrow},1} \delta_{n_{j\downarrow},1}.
\ee
The parameters $\alpha$ and $\beta$ are optimized by minimizing the variational energy of the trial wave function, \ie,
\be
E_\text{T}(\alpha,\beta) = \frac{\mel{ c_\text{T} }{H }{c_\text{T}}} {\braket{ c_\text{T} }{ c_\text{T} }}.
\ee
where the kronecker function $\delta(k)$ is equal to 1 when $k=0$ and 0 otherwise.
The parameters $\alpha, \beta$ are optimized by minimizing the variational energy,
$E_v(\alpha,\beta)=\frac{\langle c_T|H|c_T\rangle} {\langle c_T|c_T\rangle}$.\\
{\it Domains}. The efficiency of the method depends on the choice of states forming each domain.
As a general guiding principle, it is advantageous to build domains associated with a large average trapping time in order to integrate out the most important
part of the Green's matrix. Here, as a first illustration of the method, we shall
consider the large-$U$ regime of the Hubbard model where the construction of such domains is rather natural.
At large $U$ and half-filling,
the Hubbard model approaches the Heisenberg limit where only the $2^N$ states with no double occupancy, $n_D(n)=0$, have a significant weight in the wave function.
The contribution of the other states vanishes as $U$ increases with a rate increasing sharply with $n_D(n)$. In addition, for a given
number of double occupations, configurations with large values of $n_A(n)$ are favored due to their high kinetic energy (electrons move
more easily). Therefore, we build domains associated with small $n_D$ and large $n_A$ in a hierarchical way as described below.
For simplicity and decreasing the number of diagonalizations to perform,
we shall consider only one non-trivial domain called here the main domain and denoted as ${\cal D}$.
As discussed above, the efficiency of the method depends on the choice of states forming each domain.
As a general guiding principle, it is advantageous to build domains associated with a large average trapping time in order to integrate out the most important part of the Green's matrix.
Here, as a first illustration of the method, we shall consider the large-$U$ regime of the Hubbard model where the construction of such domains is rather natural.
At large $U$ and half-filling, the Hubbard model approaches the Heisenberg limit where only the $2^N$ states with no double occupancy, $n_D(n)=0$, have a significant weight in the wave function.
The contribution of the other states vanishes as $U$ increases with a rate increasing sharply with $n_D(n)$.
In addition, for a given number of double occupations, configurations with large values of $n_A(n)$ are favored due to their high kinetic energy (electrons move more easily).
Therefore, we build domains associated with small $n_D$ and large $n_A$ in a hierarchical way as described below.
For simplicity and decreasing the number of diagonalizations to perform, we shall consider only one non-trivial domain called here the main domain and denoted as $\cD$.
This domain will be chosen common to all states belonging to it, that is
\be
{\cal D}_i= {\cal D} \;\; {\rm for } \; |i\rangle \in {\cal D}
\cD_i= \cD \qq{for} \ket{i} \in \cD.
\ee
For the other states we choose a single-state domain
For the other states, we choose a single-state domain
\be
{\cal D}_i= \{ |i\rangle \} \;\; {\rm for} \; |i\rangle \notin {\cal D}
\cD_i= \qty{ \ket{i} } \qq{for} \ket{i} \notin \cD.
\ee
To define ${\cal D}$, let us introduce the following set of states
To define $\cD$, let us introduce the following set of states
\be
{\cal S}_{ij} = \{ |n\rangle; n_D(n)=i {\rm \; and\;} n_A(n)=j \}.
\titou{\cS_{ij} = \qty{ \ket{n} : n_D(n)=i \land n_A(n)=j }}.
\ee
${\cal D}$ is defined as the set of states having up to $n_D^{\rm max}$ double occupations and, for each state with a number
of double occupations equal to $m$, a number of nearest-neighbor antiparallel pairs between $n_A^{\rm min}(m)$ and $n_A^{\rm max}(m)$.
Here, $n_A^{\rm max}(m)$ will not be varied and taken to be the maximum possible for a given $m$,
$n_A^{\rm max}(m)= {\rm Max}(N-1-2m,0)$.
$\cD$ is defined as the set of states having up to $n_D^\text{max}$ double occupations and, for each state with a number of double occupations equal to $m$, a number of nearest-neighbor antiparallel pairs between $n_A^\text{min}(m)$ and $n_A^\text{max}(m)$.
Here, $n_A^\text{max}(m)$ will not be varied and taken to be the maximum possible for a given $m$, $n_A^\text{max}(m)= \max(N-1-2m,0)$.
Using these definitions, the main domain is taken as the union of some elementary domains
\be
{\cal D} = \cup_{n_D=0}^{n_D^{\rm max}}{\cal D}(n_D,n_A^{\rm min}(n_D))
\cD = \cup_{n_D=0}^{n_D^\text{max}}\cD(n_D,n_A^{\rm min}(n_D))
\ee
where the elementary domain is defined as
\be
{\cal D}(n_D,n_A^{\rm min}(n_D))=\cup_{ n_A^{\rm min}(n_D) \leq j \leq n_A^{\rm max}(n_D)} {\cal S}_{n_D j}
\cD(n_D,n_A^\text{min}(n_D))=\cup_{ n_A^\text{min}(n_D) \leq j \leq n_A^{\rm max}(n_D)} \cS_{n_D j}
\ee
The two quantities defining the main domain are thus $n_D^{\rm max}$ and
$n_A^{\rm min}(m)$.
The two quantities defining the main domain are thus $n_D^\text{max}$ and $n_A^\text{min}(m)$.
To give an illustrative example, let us consider the 4-site case. There are 6 possible elementary domains
$$
{\cal D}(0,3)= {\cal S}_{03}
$$
$$
{\cal D}(0,2)= {\cal S}_{03} \cup {\cal S}_{02}
$$
$$
{\cal D}(0,1)= {\cal S}_{03} \cup {\cal S}_{02} \cup {\cal S}_{01}
$$
$$
{\cal D}(1,1)= {\cal S}_{11}
$$
$$
{\cal D}(1,0)= {\cal S}_{11} \cup {\cal S}_{10}
$$
$$
{\cal D}(2,0)= {\cal S}_{20}
$$
\begin{align*}
\cD(0,3) & = \cS_{03},
&
\cD(0,2) & = \cS_{03} \cup \cS_{02},
\\
\cD(0,1) & = \cS_{03} \cup \cS_{02} \cup \cS_{01},
&
\cD(1,1) & = \cS_{11},
\\
\cD(1,0) & = \cS_{11} \cup \cS_{10},
&
\cD(2,0) & = \cS_{20},
\end{align*}
where
$$
{\cal S}_{03} = \{\; |\uparrow,\downarrow,\uparrow,\downarrow \rangle, |\downarrow,\uparrow,\downarrow,\uparrow \rangle \}
\;\;({\rm the\; two\; N\acute eel\; states})
$$
$$
{\cal S}_{02} = \{\; |\uparrow, \downarrow, \downarrow, \uparrow \rangle, |\downarrow, \uparrow, \uparrow, \downarrow \rangle \}
$$
$$
{\cal S}_{01} = \{\; |\uparrow, \uparrow, \downarrow, \downarrow \rangle, |\downarrow, \downarrow, \uparrow, \uparrow \rangle \}
$$
$$
{\cal S}_{11} = \{\; |\uparrow \downarrow, \uparrow ,\downarrow, 0 \rangle, |\uparrow \downarrow, 0, \uparrow,\uparrow \rangle + ...
\}
$$
$$
{\cal S}_{10} = \{\; |\uparrow \downarrow, \uparrow, 0, \downarrow \rangle, |\uparrow \downarrow, 0, \uparrow, \downarrow \rangle + ... \}
$$
$$
{\cal S}_{20} = \{\; |\uparrow \downarrow, \uparrow \downarrow, 0 ,0 \rangle + ... \}
$$
For the three last cases, the dots indicate the remaining states obtained by permuting the position of the pairs.\\
\begin{align*}
\cS_{03} & = \qty{ \ket{\uparrow,\downarrow,\uparrow,\downarrow }, \ket{\downarrow,\uparrow,\downarrow,\uparrow } },
\qq{(the two N\'eel states)}
\\
Let us now present our QMC calculations for the Hubbard model. In what follows
we shall restrict ourselves to the case of the Green's Function Monte Carlo approach where trapping times are integrated out exactly.
\cS_{02} & = \qty{ \ket{\uparrow, \downarrow, \downarrow, \uparrow }, \ket{\downarrow, \uparrow, \uparrow, \downarrow } },
\\
\cS_{01} & = \qty{ \ket{\uparrow, \uparrow, \downarrow, \downarrow }, \ket{\downarrow, \downarrow, \uparrow, \uparrow } },
\\
\cS_{11} & = \qty{ \ket{\uparrow \downarrow, \uparrow ,\downarrow, 0 }, \ket{\uparrow \downarrow, 0, \uparrow,\uparrow } + \ldots },
\\
\cS_{10} & = \qty{ \ket{\uparrow \downarrow, \uparrow, 0, \downarrow }, \ket{\uparrow \downarrow, 0, \uparrow, \downarrow } + \ldots },
\\
\cS_{20} & = \qty{ \ket{\uparrow \downarrow, \uparrow \downarrow, 0 ,0 } + \ldots }.
\end{align*}
For the three last cases, the dots indicate the remaining states obtained by permuting the position of the pairs.
Following Eqs.(\ref{final_E},\ref{calE}), the practical formula used for computing the QMC energy is written as
Let us now present our QMC calculations for the Hubbard model.
In what follows we shall restrict ourselves to the case of the Green's function Monte Carlo approach where trapping times are integrated out exactly.
Following Eqs.~\eqref{eq:final_E} and \eqref{calE}, the practical formula used for computing the QMC energy is written as
\be
\cE_{QMC}(E,p_{ex},p_{max})= \frac{H_0 +...+ H_{p_{ex}} + \sum_{p=p_{ex}+1}^{p_{max}} H^{QMC}_p }
{S_0 +...+ S_{p_{ex}} + \sum_{p=p_{ex}+1}^{p_{max}} S^{QMC}_p }
\cE^\text{QMC}(E,p_{ex},p_{max})= \frac{H_0 +...+ H_{p_{ex}} + \sum_{p=p_{ex}+1}^{p_{max}} H^\text{QMC}_p }
{S_0 +...+ S_{p_{ex}} + \sum_{p=p_{ex}+1}^{p_{max}} S^\text{QMC}_p }
\label{calE}
\ee
where $p_{ex}+1$ is the number of pairs, ($H_p$, $S_p$), computed analytically. For $p_{ex} < p \le p_{max}$
the Monte Carlo estimates are written as
where $p_{ex}+1$ is the number of pairs, ($H_p$, $S_p$), computed analytically.
For $p_{ex} < p \le p_{max}$ the Monte Carlo estimates are written as
\be
H^{QMC}_p= \Bigg \langle \Big( \prod_{k=0}^{p-1} W^E_{I_k I_{k+1}} \Big) \langle I_p| P_{I_p} \frac{1}{H-E} P_{I_p}|H\PsiT\rangle \Bigg \rangle
H^\text{QMC}_p= \expval{ \qty( \prod_{k=0}^{p-1} W^E_{I_k I_{k+1}} ) \mel{ I_p }{ P_{I_p} \frac{1}{H-E} P_{I_p} }{ H\PsiT } }
\ee
and
\be
S^{QMC}_p= \Bigg \langle \Big(\prod_{k=0}^{p-1} W^E_{I_k I_{k+1}} \Big) \langle I_p| P_{I_p} \frac{1}{H-E} P_{I_p}|\PsiT\rangle \Bigg \rangle.
S^\text{QMC}_p= \expval{ \qty(\prod_{k=0}^{p-1} W^E_{I_k I_{k+1}} ) \mel{ I_p }{ P_{I_p} \frac{1}{H-E} P_{I_p} }{ \PsiT } }.
\ee
Let us begin with a small chain of 4 sites with $U=12$. From now on, we shall take $t=1$.
The size of the linear space is ${\binom{4}{2}}^2= 36$ and the ground-state energy obtained by exact diagonalization
is $E_0=-0.768068...$. The two variational parameters of the trial vector have been optimized and fixed at the values of
$\alpha=1.292$, and $\beta=0.552$ with a variational energy
of $E_v=-0.495361...$. In what follows
Let us begin with a small chain of 4 sites with $U=12$.
From now on, we shall take $t=1$.
The size of the linear space is ${\binom{4}{2}}^2= 36$ and the ground-state energy obtained by exact diagonalization is $E_0=-0.768068...$. The two variational parameters of the trial vector have been optimized and fixed at the values of $\alpha=1.292$, and $\beta=0.552$ with a variational energy of $E_\text{T}=-0.495361...$. In what follows
$|I_0\rangle$ will be systematically chosen as one of the two N\'eel states, {\it e.g.} $|I_0\rangle =|\uparrow,\downarrow, \uparrow,...\rangle$.
Figure \ref{fig1} shows the convergence of $H_p$ as a function of $p$ for different values of the reference energy $E$.
@ -1049,7 +1026,7 @@ we need to use better trial wave functions, better domains, and also larger simu
laptop. Of course, it will also be particularly interesting to take advantage of the fully parallelizable character of the algorithm to get much lower error bars.
All these aspects will be considered in a forthcoming work.
\begin{table}[h!]
\begin{table}
\caption{$N=4$, $U=12$, and $E=-1$. Dependence of the statistical error on the energy with the number of $p$-components calculated
analytically. Same simulation as for Table \ref{tab1}. Results are presented when a single-state domain
is used for all states and when
@ -1073,7 +1050,7 @@ $8$ & $2.2 \times10^{-5}$ &$ 0.05 \times 10^{-8}$\\
\end{table}
\begin{table}[h!]
\begin{table}
\caption{$N=4$, $U=12$, $\alpha=1.292$, $\beta=0.552$. Main domain = ${\cal D}(0,1) \cup {\cal D}(1,0)$. Simulation with 20 independent blocks and $10^6$ paths.
$p_{ex}=4$. The various fits are done with the five values of $E$}
\label{tab3}
@ -1094,7 +1071,7 @@ $E_0$ exact & -0.768068...\\
\end{ruledtabular}
\end{table}
\begin{table}[h!]
\begin{table}
\caption{$N$=4, Domain ${\cal D}(0,1) \cup {\cal D}(1,0)$}
\label{tab4}
\begin{ruledtabular}
@ -1112,7 +1089,7 @@ $U$ & $\alpha,\beta$ & $E_{var}$ & $E_{ex}$ & $\bar{t}_{I_0}$ \\
\end{ruledtabular}
\end{table}
\begin{table*}[h!]
\begin{table*}
\caption{$U=12$. The fits to extrapolate the QMC energies are done using the two-component function}
\label{tab5}
\begin{ruledtabular}
@ -1143,11 +1120,12 @@ To obtain the effective transition probability (which provides the probability o
This is the main computationally intensive step of the present approach.
The efficiency of the method is directly related to the importance of the average time spent by the stochastic trajectories in each domain.
Being able to define domains with large average trapping times is the key aspect of the method since it may lead to some important reduction of the statistical error as illustrated in our numerical applications.
Therefore, a trade-off has to be found between maximizing the average trapping time and minimizing the cost of computing the local Green's functions.
Being able to define domains with large average trapping times is the key aspect of the method since it may lead to some important reduction of the statistical error, as illustrated in our numerical applications.
Therefore, a trade-off has to be found between maximizing the average trapping time and minimizing the cost of computing the domain Green's functions.
In practice, there is no general rule to construct such domains.
For each system at hand, one needs to determine, on physical grounds, which \titou{spatial} regions are preferentially sampled by the stochastic trajectories and to build domains of minimal size enclosing such regions.
In the first application presented here on the 1D Hubbard model we exploit the physics of the large-$U$ regime that is known to approach the Heisenberg limit where double occupations have small weights.
In the first application presented here on the one-dimensional Hubbard model, we exploit the physics of the large-$U$ regime that is known to approach the Heisenberg limit where double occupations have small weights.
This simple example has been chosen to illustrate the various aspects of the approach.
Our goal is, of course, to tackle much larger systems, like those treated by state-of-the-art methods, such as selected CI, \cite{Huron_1973,Harrison_1991,Giner_2013,Holmes_2016,Schriber_2016,Tubman_2020}, FCIQMC, \cite{Booth_2009,Cleland_2010}, AFQMC, \cite{Zhang_2003} or DMRG. \cite{White_1999,Chan_2011}
@ -1158,11 +1136,8 @@ This is outside the scope of the present study and will be presented in a forthc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\acknowledgements{
P.F.L., A.S., and M.C. acknowledge funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (Grant agreement No.~863481).
This work was supported by the European Centre of
Excellence in Exascale Computing TREX --- Targeting Real Chemical
Accuracy at the Exascale. This project has received funding from the
European Union's Horizon 2020 --- Research and Innovation program ---
under grant agreement no.~952165.
This work was supported by the European Centre of Excellence in Exascale Computing TREX --- Targeting Real Chemical Accuracy at the Exascale.
This project has received funding from the European Union's Horizon 2020 --- Research and Innovation program --- under grant agreement no.~952165.
}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@ -1171,91 +1146,105 @@ under grant agreement no.~952165.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Particular case of the $2\times2$ matrix}
\section{Particular case of a $2\times2$ matrix}
\label{app:A}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
For the simplest case of a two-state system the fundamental equation \eqref{eq:eqfond} reads
For the simplest case of a two-state system, $\ket{1}$ and $\ket{2}$, the fundamental equation given in Eq.~\eqref{eq:eqfond} simplifies to
\be
{\cal I}=
\langle I_0|\frac{1}{H-E}|\Psi\rangle = \langle I_0|P_0\frac{1}{H-E} P_0|\Psi\rangle
+ \sum_{p=1}^{\infty} {\cal I}_p
\begin{split}
\cI
& = \mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ \Psi }
\\
& = \mel{ I_0 }{ P_0 \qty(H-E \Id)^{-1} P_0 }{ \Psi }
+ \sum_{p=1}^{\infty} \cI_p,
\end{split}
\ee
with
\be
{\cal I}_p=
\sum_{I_1 \notin {\cal D}_0, \hdots , I_p \notin {\cal D}_{p-1}}
\ee
\be
\Big[ \prod_{k=0}^{p-1} \langle I_k| P_k \frac{1}{H-E} P_k (-H)(1-P_k)|I_{k+1} \rangle \Big]
\langle I_p| P_p \frac{1} {H-E} P_p|\Psi\rangle
\ee
To treat simultaneously the two possible cases for the final state, $|I_N\rangle =|1\rangle$ or $|2\rangle$,
the equation has been slightly generalized to the case of a general vector for the final state written as
\be
|\Psi\rangle = \Psi_1 |1\rangle + \Psi_2 |2\rangle
\ee
where
$|1\rangle$ and $|2\rangle$ denote the two states. Let us choose a single-state domain for both states, namely ${\cal D}_1=\{|1\rangle \}$ and ${\cal D}_2=\{ |2\rangle \}$. Note that the single exit state for each state is the other state. Thus there are only two possible deterministic "alternating" paths, namely either
$|1\rangle \to |2\rangle \to |1\rangle,...$ or $|2\rangle \to |1\rangle \to |2\rangle,...$
We introduce the following quantities
\begin{align}
A_1 & = \langle 1| P_1 \frac{1}{H-E} P_1 (-H)(1-P_1)|2\rangle \;\;\;
\begin{multline}
\cI_p= \sum_{I_1 \notin \cD_0, \ldots , I_p \notin \cD_{p-1}}
\qty[ \prod_{k=0}^{p-1} \mel{ I_k }{ P_k \qty(H-E \Id)^{-1} P_k (-H)(1-P_k) }{ I_{k+1} } ]
\\
A_2 & = \langle 2| P_2 \frac{1}{H-E} P_2 (-H) (1-P_2)|1\rangle
\label{defA}
\times \mel{ I_p }{ P_p \qty(H-E \Id)^{-1} P_p }{ \Psi }.
\end{multline}
To treat simultaneously the two possibilities for the final state, \ie, $\ket{I_N} = \ket{1}$ or $\ket{2}$, \titou{Eq.~\eqref{}} has been slightly generalized to the case of a general vector for the final state written as
\be
\ket{\Psi} = \Psi_1 \ket{1} + \Psi_2 \ket{2}.
\ee
Let us choose a single-state domain for both states, namely $\cD_1 = \qty{ \ket{1} }$ and $\cD_2 = \qty{ \ket{2} }$.
%\titou{Note that the single exit state for each state is the other state.}
Note that, due to the simplicity of the present two-state model, there are only two possible deterministic ``alternating'' paths, namely, $\ket{1} \to \ket{2} \to \ket{1},\ldots$ and $\ket{2} \to \ket{1} \to \ket{2},\ldots$.
For the sake of convenience, we introduce the following quantities
\begin{align}
\label{eq:defA1}
A_1 & = \mel{ 1 }{ P_1 \qty(H-E \Id)^{-1} P_1 (-H)(1-P_1) }{ 2 },
\\
\label{eq:defA2}
A_2 & = \mel{ 2 }{ P_2 \qty(H-E \Id)^{-1} P_2 (-H) (1-P_2) }{ 1 },
\end{align}
and
\begin{align}
C_1 & = \mel{ 1 }{ P_1 \qty(H-E \Id)^{-1} P_1 }{ \Psi },
\\
C_2 & = \mel{ 2 }{ P_2 \qty(H-E \Id)^{-1} P_2 }{ \Psi }.
\end{align}
Without loss of generality, let us choose, for example, $\ket{I_0} = \ket{1}$.
Then, it is easy to show that
%\begin{align}
% \cI_1 & = A_1 C_2,
% &
% \cI_2 & = A_1 A_2 C_1,
% &
% \cI_3 & = A_1 A_2 A_1 C_2,
% &
% \qq{etc.}
%\end{align}
%that is
\begin{align}
\cI_{2k+1} & = \frac{C_2}{A_2} (A_1 A_2)^{2k+1},
&
\cI_{2k} & = C_1 (A_1 A_2)^{2k},
\end{align}
which yields
\be
C_1= \langle 1| P_1 \frac{1}{H-E} P_1 |\Psi\rangle \;\;\;
C_2= \langle 2| P_2 \frac{1}{H-E} P_2 |\Psi\rangle
\ee
Let us choose, for example, $|I_0\rangle =|1\rangle$, we then have
\be
{\cal I}_1= A_1 C_2 \;\; {\cal I}_2 = A_1 A_2 C_1 \;\; {\cal I}_3 = A_1 A_2 A_1 C_2 \;\; etc.
\ee
that is
\be
{\cal I}_{2k+1} = \frac{C_2}{A_2} (A_1 A_2)^{2k+1}
\sum_{p=1}^{\infty}
\cI_p
= \frac{C_2}{A_2} \sum_{p=1}^{\infty} (A_1 A_2)^p + C_1 \sum_{p=1}^{\infty} (A_1 A_2)^p
= A_1 \frac{C_2 + C_1 A_2}{1-A_1 A_2},
\ee
and
\be
{\cal I}_{2k} = C_1 (A_1 A_2)^{2k}
\mel{ 1 }{ \qty(H-E \Id)^{-1} }{ \Psi }
= C_1 + A_1 \frac{C_2 + C_1 A_2}{1-A_1 A_2}.
\ee
Then
For a $2\times2$ matrix of the form
\be
\sum_{p=1}^{\infty}
{\cal I}_p
= \frac{C_2}{A_2} \big[ \sum_{p=1}^{\infty} (A_1 A_2)^p \big] + {C_1} \big[ \sum_{p=1}^{\infty} (A_1 A_2)^p \big]
H =
\begin{pmatrix}
H_{11} & H_{12} \\
H_{12} & H_{22} \\
\end{pmatrix},
\ee
which gives
it is easy to evaluate the $A_i$'s.
Using Eqs.~\eqref{eq:defA1} and \eqref{eq:defA2}, one gets, for $i = 1$ or $2$,
\begin{align}
A_i & = -\frac{H_{12}}{H_{ii}-E},
&
C_i & = \frac{1}{H_{ii}-E} \Psi_i.
\end{align}
which finally yields
\be
\sum_{p=1}^{\infty}
{\cal I}_p
= A_1 \frac{C_2 + C_1 A_2}{1-A_1 A_2}
\mel{ 1 }{ \qty(H-E \Id)^{-1} }{ \Psi}
% = \qty(H_{11}-E)^{-1} \Psi_1
% - H_{12} \frac{ \qty( \Psi_2 - \frac{ H_{12}\Psi_1}{H_{11}-E} )}
% { (H_{11}-E)(H_{22}-E) - H^2_{12}}
% \\
= \frac{ H_{22}-E}{\Delta} \Psi_1 - \frac{H_{12}}{\Delta} \Psi_2,
\ee
and thus
\be
\langle 1 | \frac{1} {H-E} |\Psi\rangle
=
C_1 + A_1 \frac{C_2 + C_1 A_2}{1-A_1 A_2}
\ee
For a two by two matrix it is easy to evaluate the $A_i$'s, Eq.\ref{defA}, we
have
$$
A_i = -\frac{H_{12}}{H_{ii}-E} \;\;\; \;\;\; C_i= \frac{1}{H_{ii}-E} \Psi_i \;\;\; i=1,2
$$
Then
$$
\langle 1 | \frac{1} {H-E} |\Psi\rangle
= \frac{1}{H_{11}-E} \Psi_1 -H_{12} \frac{ \big( \Psi_2 - \frac{ H_{12}\Psi_1}{H_{11}-E} \big)}
{ (H_{11}-E)(H_{22}-E) - H^2_{12}}
$$
$$
= \frac{ (H_{22}-E) \Psi_1}{\Delta} -\frac{H_{12} \Psi_2}{\Delta}
$$
where $\Delta$ is the determinant of the matrix $H$. It is easy to check that the RHS is equal to the LHS, thus validating the fundamental equation
for this particular case.
where $\Delta$ is the determinant of $H$.
\titou{It is easy to check that the RHS is equal to the LHS, thus validating the fundamental equation for this particular case. ??}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\bibliography{g}