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

only results left to clean up

Browse Source
This commit is contained in:
Pierre-Francois Loos 2022-09-21 14:17:43 +02:00
parent 007a26ab5e
commit 2b76d91354
1 changed files with 48 additions and 50 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

96
g.tex
View File

@ -232,17 +232,16 @@ For example, in the important case of the energy, one can rely on the following
\label{eq:E0} \label{eq:E0}
E_0 = \lim_{N \to \infty } \frac{\mel{\PsiT}{H T^N}{\Psi_0}}{\mel{\PsiT}{T^N}{\Psi_0}}, E_0 = \lim_{N \to \infty } \frac{\mel{\PsiT}{H T^N}{\Psi_0}}{\mel{\PsiT}{T^N}{\Psi_0}},
\ee \ee
where $|\PsiT\rangle$ is some trial vector (some approximation of the true ground-state wave function) on which $T^N \ket{\Psi_0}$ is projected out. where $\ket{\PsiT}$ is a trial wave function (some approximation of the true ground-state wave function) on which $T^N \ket{\Psi_0}$ is projected out.
To proceed further we introduce the time-dependent Green's matrix $G^{(N)}$ defined as To proceed further we introduce the time-dependent Green's matrix $G^{(N)}$ defined as
\be \be
G^{(N)}_{ij}=\mel{j}{T^N}{i}. G^{(N)}_{ij}=\mel{j}{T^N}{i}.
\ee \ee
where $\ket{i}$ and $\ket{j}$ are basis vectors. where $\ket{i}$ and $\ket{j}$ are basis vectors.
The denomination ``time-dependent Green's matrix'' is used here since $G$ may be viewed as a short-time approximation of the (time-imaginary) evolution operator, The denomination ``time-dependent Green's matrix'' is used here since $G$ may be viewed as a short-time approximation of the (time-imaginary) evolution operator $e^{-N\tau H}$, which is usually referred to as the imaginary-time dependent Green's function.
$e^{-N\tau H}$ which is usually referred to as the imaginary-time dependent Green's function.
\titou{Introducing the set of $N-1$ intermediate states, $\{ \ket{i_k} \}_{1 \le k \le N-1}$, in the $N$th product of $T$,} $G^{(N)}$ can be written in the following expanded form Introducing the set of $N-1$ intermediate states, $\{ \ket{i_k} \}_{1 \le k \le N-1}$, in $T^N$, $G^{(N)}$ can be written in the following expanded form
\be \be
\label{eq:cn} \label{eq:cn}
G^{(N)}_{i_0 i_N} = \sum_{i_1} \sum_{i_2} \cdots \sum_{i_{N-1}} \prod_{k=0}^{N-1} T_{i_{k} i_{k+1}}, G^{(N)}_{i_0 i_N} = \sum_{i_1} \sum_{i_2} \cdots \sum_{i_{N-1}} \prod_{k=0}^{N-1} T_{i_{k} i_{k+1}},
@ -250,8 +249,8 @@ $e^{-N\tau H}$ which is usually referred to as the imaginary-time dependent Gree
where $T_{ij} =\mel{i}{T}{j}$. where $T_{ij} =\mel{i}{T}{j}$.
Here, each index $i_k$ runs over all basis vectors. Here, each index $i_k$ runs over all basis vectors.
In quantum physics, Eq.~\eqref{eq:cn} is referred to as the path-integral representation of the Green's matrix (function). In quantum physics, Eq.~\eqref{eq:cn} is referred to as the path-integral representation of the Green's matrix (or function).
The series of states $\ket{i_0}, \ldots,\ket{i_N}$ is interpreted as a ``path'' in the Hilbert space starting at vector $\ket{i_0}$ and ending at vector $\ket{i_N}$ where $k$ plays the role of a time index. The series of states $\ket{i_0}, \ldots,\ket{i_N}$ is interpreted as a ``path'' in the Hilbert space starting at vector $\ket{i_0}$ and ending at vector $\ket{i_N}$, where $k$ plays the role of a time index.
Each path is associated with a weight $\prod_{k=0}^{N-1} T_{i_{k} i_{k+1}}$ and the path integral expression of $G$ can be recast in the more suggestive form as follows: Each path is associated with a weight $\prod_{k=0}^{N-1} T_{i_{k} i_{k+1}}$ and the path integral expression of $G$ can be recast in the more suggestive form as follows:
\be \be
\label{eq:G} \label{eq:G}
@ -261,7 +260,7 @@ Each path is associated with a weight $\prod_{k=0}^{N-1} T_{i_{k} i_{k+1}}$ and
This expression allows a simple and vivid interpretation of the solution. This expression allows a simple and vivid interpretation of the solution.
In the limit $N \to \infty$, the $i_N$th component of the ground-state wave function (obtained as $\lim_{N \to \infty} G^{(N)}_{i_0 i_N})$ is the weighted sum over all possible paths arriving at vector $\ket{i_N}$. In the limit $N \to \infty$, the $i_N$th component of the ground-state wave function (obtained as $\lim_{N \to \infty} G^{(N)}_{i_0 i_N})$ is the weighted sum over all possible paths arriving at vector $\ket{i_N}$.
This result is independent of the initial vector $\ket{i_0}$, apart from some irrelevant global phase factor. This result is independent of the initial vector $\ket{i_0}$, apart from some irrelevant global phase factor.
When the size of the linear space is small the explicit calculation of the full sums involving $M^N$ terms (where $M$ is the size of the Hilbert space) can be performed. When the size of the linear space is small, the explicit calculation of the full sums involving $M^N$ terms (where $M$ is the size of the Hilbert space) can be performed.
In such a case, we are in the realm of what one would call the ``deterministic'' power methods, such as the Lancz\`os or Davidson approaches. In such a case, we are in the realm of what one would call the ``deterministic'' power methods, such as the Lancz\`os or Davidson approaches.
If not, probabilistic techniques for generating only the paths contributing significantly to the sums are to be preferred. If not, probabilistic techniques for generating only the paths contributing significantly to the sums are to be preferred.
This is the central theme of the present work. This is the central theme of the present work.
@ -271,26 +270,26 @@ This is the central theme of the present work.
\label{sec:proba} \label{sec:proba}
%=======================================% %=======================================%
In order to derive a probabilistic expression for the Green's matrix we introduce a so-called guiding vector, $\ket{\PsiG}$, having strictly positive components, \ie, $\PsiG_i > 0$, in order to apply a similarity transformation to the operators $G^{(N)}$ and $T$ as follows: To derive a probabilistic expression for the Green's matrix, we introduce a guiding wave function, $\ket{\PsiG}$, having strictly positive components, \ie, $\PsiG_i > 0$, in order to perform a similarity transformation to the operators $G^{(N)}$ and $T$, as follows:
\begin{align} \begin{align}
\label{eq:defT} \label{eq:defT}
\bar{T}_{ij} & = \frac{\PsiG_j}{\PsiG_i} T_{ij}, \bar{T}_{ij} & = \frac{\PsiG_j}{\PsiG_i} T_{ij},
& &
\bar{G}^{(N)}_{ij}& = \frac{\PsiG_j}{\PsiG_i} G^{(N)}_{ij}. \bar{G}^{(N)}_{ij}& = \frac{\PsiG_j}{\PsiG_i} G^{(N)}_{ij}.
\end{align} \end{align}
Note that, under this similarity transformation, the path integral expression relating $G^{(N)}$ and $T$ [see Eq.~\eqref{eq:G}] remains unchanged for $\bar{G}^{(N)}$ and $\bar{T}$. Note that, thanks to the properties of similarity transformations, the path integral expression relating $G^{(N)}$ and $T$ [see Eq.~\eqref{eq:G}] remains unchanged for $\bar{G}^{(N)}$ and $\bar{T}$.
Next, the matrix elements of $\bar{T}$ are expressed as those of a stochastic matrix multiplied by some residual weight, namely, Next, the matrix elements of $\bar{T}$ are expressed as those of a stochastic matrix multiplied by some residual weight, namely,
\be \be
\label{eq:defTij} \label{eq:defTij}
\bar{T}_{ij}= p_{i \to j} w_{ij}. \bar{T}_{ij}= p_{i \to j} w_{ij}.
\ee \ee
Here, we recall that a stochastic matrix is defined as a matrix with positive entries and obeying Here, we recall that a stochastic matrix is defined as a matrix with positive entries that obeys
\be \be
\label{eq:sumup} \label{eq:sumup}
\sum_j p_{i \to j}=1. \sum_j p_{i \to j}=1.
\ee \ee
To build the transition probability density the following operator is introduced To build the transition probability density, the following operator is introduced
%As known, there is a natural way of associating a stochastic matrix to a matrix having a positive ground-state vector (here, a positive vector is defined here as %As known, there is a natural way of associating a stochastic matrix to a matrix having a positive ground-state vector (here, a positive vector is defined here as
%a vector with all components positive). %a vector with all components positive).
\be \be
@ -330,22 +329,22 @@ We are now in the position to define the stochastic matrix as
\tau \frac{\PsiG_{j}}{\PsiG_{i}} \abs{H_{ij}} \ge 0, & \qif* i\neq j. \tau \frac{\PsiG_{j}}{\PsiG_{i}} \abs{H_{ij}} \ge 0, & \qif* i\neq j.
\end{cases} \end{cases}
\ee \ee
As readily seen in Eq.~\eqref{eq:pij}, the off-diagonal terms of the stochastic matrix are positive, while the diagonal ones can be made positive if $\tau$ is chosen sufficiently small via the condition As readily seen in Eq.~\eqref{eq:pij}, the off-diagonal terms of the stochastic matrix are positive, while the diagonal terms can be made positive if $\tau$ is chosen sufficiently small via the condition
\be \be
\label{eq:cond} \label{eq:cond}
\tau \leq \frac{1}{\max_i\abs{H^+_{ii}-(\EL^+)_{i}}}. \tau \leq \frac{1}{\max_i\abs{H^+_{ii}-(\EL^+)_{i}}}.
\ee \ee
The sum-over-states condition [see Eq.~\eqref{eq:sumup}] The sum-over-states condition [see Eq.~\eqref{eq:sumup}],
\be \be
\sum_j p_{i \to j}= \frac{\mel{i}{T^+}{\PsiG}}{\PsiG_{i}} = 1. \sum_j p_{i \to j}= \frac{\mel{i}{T^+}{\PsiG}}{\PsiG_{i}} = 1,
\ee \ee
follows from the fact that $\ket{\PsiG}$ is eigenvector of $T^+$ [see Eq.~\eqref{eq:relT+}]. follows from the fact that $\ket{\PsiG}$ is eigenvector of $T^+$, as evidenced by Eq.~\eqref{eq:relT+}.
This ensures that $p_{i \to j}$ is indeed a stochastic matrix. This ensures that $p_{i \to j}$ is indeed a stochastic matrix.
At first sight, the condition defining the maximum value of $\tau$ allowed, Eq.~\eqref{eq:cond}, may appear rather tight since, for very large matrices, it may impose an extremely small value of the time step. At first sight, the condition defining the maximum value of $\tau$ [see Eq.~\eqref{eq:cond}] may appear rather tight since, for very large matrices, it may impose an extremely small value of the time step.
However, in practice, during the simulation only a (tiny) fraction of the linear space is sampled, and the maximum absolute value of $H^+_{ii}-(\EL^+)_{i}$ for the sampled states turns out to be not too large. However, in practice, during the simulation only a (tiny) fraction of the linear space is sampled, and the maximum absolute value of $H^+_{ii}-(\EL^+)_{i}$ for the sampled states turns out to be not too large.
Hence, reasonable values of $\tau$ can be selected without violating the positivity of the transition probability matrix. Hence, reasonable values of $\tau$ can be selected without violating the positivity of the transition probability matrix.
\titou{Note that one can eschew this condition via a simple generalization of the transition probability matrix:} Note that one can eschew this condition via a simple generalization of the transition probability matrix:
\be \be
p_{i \to j} p_{i \to j}
= \frac{ \frac{\PsiG_{j}}{\PsiG_{i}} \abs{\mel{i}{T^+}{j}} } = \frac{ \frac{\PsiG_{j}}{\PsiG_{i}} \abs{\mel{i}{T^+}{j}} }
@ -393,7 +392,7 @@ In this framework, the energy defined in Eq.~\eqref{eq:E0} is given by
\frac{ \expval{ \prod_{k=0}^{N-1} w_{i_{k}i_{k+1}} {(H\PsiT)}_{i_N}} } \frac{ \expval{ \prod_{k=0}^{N-1} w_{i_{k}i_{k+1}} {(H\PsiT)}_{i_N}} }
{ \expval{ \prod_{k=0}^{N-1} w_{i_{k}i_{k+1}} {\PsiT}_{i_N} }}. { \expval{ \prod_{k=0}^{N-1} w_{i_{k}i_{k+1}} {\PsiT}_{i_N} }}.
\ee \ee
Note that, instead of using a single walker, it is possible to introduce a population of independent walkers and to calculate the averages over this population. Note that, instead of using a single walker, it is common to introduce a population of independent walkers and to calculate the averages over this population.
In addition, thanks to the ergodic property of the stochastic matrix (see, for example, Ref.~\onlinecite{Caffarel_1988}), a unique path of infinite length from which sub-paths of length $N$ can be extracted may also be used. In addition, thanks to the ergodic property of the stochastic matrix (see, for example, Ref.~\onlinecite{Caffarel_1988}), a unique path of infinite length from which sub-paths of length $N$ can be extracted may also be used.
We shall not here insist on these practical details that are discussed, for example, in Refs.~\onlinecite{Foulkes_2001,Kolorenc_2011}. We shall not here insist on these practical details that are discussed, for example, in Refs.~\onlinecite{Foulkes_2001,Kolorenc_2011}.
@ -465,7 +464,7 @@ Details of the implementation of this effective dynamics can be in found in Refs
Let us now extend the results of Sec.~\ref{sec:single_domains} to a general domain. Let us now extend the results of Sec.~\ref{sec:single_domains} to a general domain.
To do so, we associate to each state $\ket{i}$ a set of states, called the domain of $\ket{i}$ denoted $\cD_i$, consisting of the state $\ket{i}$ plus a certain number of states. To do so, we associate to each state $\ket{i}$ a set of states, called the domain of $\ket{i}$ denoted $\cD_i$, consisting of the state $\ket{i}$ plus a certain number of states.
No particular constraints on the type of domains are imposed. No particular constraints on the type of domains are imposed.
For example, domains associated with different states can be identical, and they may or may not have common states. For example, domains associated with different states may be identical, and they may or may not have common states.
The only important condition is that the set of all domains ensures the ergodicity property of the effective stochastic dynamics, that is, starting from any state, there is a non-zero probability to reach any other state in a finite number of steps. The only important condition is that the set of all domains ensures the ergodicity property of the effective stochastic dynamics, that is, starting from any state, there is a non-zero probability to reach any other state in a finite number of steps.
In practice, it is not difficult to impose such a condition. In practice, it is not difficult to impose such a condition.
@ -525,7 +524,7 @@ corresponding to the last move connecting the inside and outside regions of the
\ee \ee
Physically, $F$ may be seen as a flux operator through the boundary of ${\cal D}_{I}$. Physically, $F$ may be seen as a flux operator through the boundary of ${\cal D}_{I}$.
\titou{Now}, the probability of being trapped $n$ times within ${\cal D}_{I}$ is given by \titou{Now}, the probability of being trapped $n$ times in ${\cal D}_{I}$ is given by
\be \be
\label{eq:PiN} \label{eq:PiN}
P_{I}(n) = \frac{1}{\PsiG_{I}} \mel{ I }{ \qty(T^+_{I})^{n-1} F^+_{I} }{ \PsiG }. P_{I}(n) = \frac{1}{\PsiG_{I}} \mel{ I }{ \qty(T^+_{I})^{n-1} F^+_{I} }{ \PsiG }.
@ -578,7 +577,7 @@ It follows that
\label{eq:Gt} \label{eq:Gt}
G^{(N)}_{I_0 I_N}= G^{(N),{\cal D}}_{I_0 I_N} + G^{(N)}_{I_0 I_N}= G^{(N),{\cal D}}_{I_0 I_N} +
\sum_{p=1}^{N} \sum_{p=1}^{N}
\sum_{|I_1\rangle \notin {\cal D}_{I_0}, \ldots , |I_p\rangle \notin {\cal D}_{I_{p-1}} } \sum_{\ket{I_1} \notin \cD_{I_0}, \ldots , \ket{I_p} \notin \cD_{I_{p-1}} }
\sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1} \sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1}
\delta_{\sum_{k=0}^p n_k,N+1} \delta_{\sum_{k=0}^p n_k,N+1}
\\ \\
@ -602,7 +601,7 @@ To express the fundamental equation of $G$ under the form of a probabilistic ave
\sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1} \sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1}
\\ \\
\times \times
\delta_{\sum_k n_k,N+1} \qty{ \prod_{k=0}^{p-1} \qty[ \frac{\PsiG_{I_{k+1}}}{\PsiG_{I_k}} \mel{ I_k }{ \qty(T_{I_k})^{n_k-1} F_{I_k} }{ I_{k+1} } ] } \delta_{\sum_{k=0}^p n_k,N+1} \qty{ \prod_{k=0}^{p-1} \qty[ \frac{\PsiG_{I_{k+1}}}{\PsiG_{I_k}} \mel{ I_k }{ \qty(T_{I_k})^{n_k-1} F_{I_k} }{ I_{k+1} } ] }
\bar{G}^{(n_p-1),\cD}_{I_p I_N}. \bar{G}^{(n_p-1),\cD}_{I_p I_N}.
\end{multline} \end{multline}
Introducing the weights Introducing the weights
@ -622,7 +621,7 @@ where, in this context, the average of a given function $F$ is defined as
\expval{F} \expval{F}
= \sum_{\ket{I_1} \notin \cD_{I_0}, \ldots , \ket{I_p} \notin \cD_{I_{p-1}}} = \sum_{\ket{I_1} \notin \cD_{I_0}, \ldots , \ket{I_p} \notin \cD_{I_{p-1}}}
\sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1} \sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1}
\delta_{\sum_k n_k,N+1} \delta_{\sum_{k=0}^p n_k,N+1}
\\ \\
\times \times
\prod_{k=0}^{N-1}\cP_{I_k \to I_{k+1}}(n_k-1) F(I_0,n_0;\ldots;I_N,n_N). \prod_{k=0}^{N-1}\cP_{I_k \to I_{k+1}}(n_k-1) F(I_0,n_0;\ldots;I_N,n_N).
@ -658,14 +657,15 @@ Let us define the energy-dependent Green's matrix
G^E_{ij}= \tau \sum_{N=0}^{\infty} \mel{ i }{ T^N }{ j} = \mel{i}{ \qty( H-E \Id )^{-1} }{j}, G^E_{ij}= \tau \sum_{N=0}^{\infty} \mel{ i }{ T^N }{ j} = \mel{i}{ \qty( H-E \Id )^{-1} }{j},
\ee \ee
which does not depends on the time step. which does not depends on the time step.
Note that, \titou{in a continuous space}, this quantity is essentially the \titou{Laplace transform of the time-dependent Green's function}. Note that, in a continuous space, this quantity is essentially the \titou{Laplace transform of the time-dependent Green's function}.
We shall use the same denomination in the following. We shall use the same denomination in the following.
The remarkable property is that, thanks to the summation over $N$ up to the infinity the constrained multiple sums appearing in Eq.~\eqref{eq:Gt} can be factorized in terms of a product of unconstrained single sums as follows The remarkable property is that, thanks to the summation over $N$ up to infinity, the constrained multiple sums appearing in Eq.~\eqref{eq:Gt} can be factorized in terms of a product of unconstrained single sums, as follows
\be \begin{multline}
\sum_{N=1}^\infty \sum_{p=1}^N \sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1} \delta_{n_0+...+n_p,N+1} \titou{??} \sum_{N=1}^\infty \sum_{p=1}^N \sum_{n_0 \ge 1} \cdots \sum_{n_p \ge 1} \delta_{\sum_{k=0}^p n_k,N+1} F(n_0,\ldots,n_N)
= \sum_{p=1}^{\infty} \sum_{n_0=1}^{\infty} \cdots \sum_{n_p=1}^{\infty} \titou{??}. \\
\ee = \sum_{p=1}^{\infty} \sum_{n_0=1}^{\infty} \cdots \sum_{n_p=1}^{\infty} F(n_0,\ldots,n_N).
It is then a trivial matter to integrate out exactly the $n_k$ variables, leading to \end{multline}
It is then trivial to integrate out exactly the $n_k$ variables, leading to
\begin{multline} \begin{multline}
\label{eq:eqfond} \label{eq:eqfond}
\mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ I_N } \mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ I_N }
@ -677,7 +677,7 @@ It is then a trivial matter to integrate out exactly the $n_k$ variables, leadin
\times \times
\mel{ I_p }{ P_p \qty( H-E \Id)^{-1} P_p }{ I_N } \mel{ I_p }{ P_p \qty( H-E \Id)^{-1} P_p }{ I_N }
\end{multline} \end{multline}
As an illustration, Appendix \ref{app:A} reports the exact derivation of this formula in the case of a two-state system. As a didactical example, Appendix \ref{app:A} reports the exact derivation of this formula in the case of a two-state system.
%----------------------------% %----------------------------%
\subsubsection{Dyson equation} \subsubsection{Dyson equation}
@ -694,7 +694,7 @@ where $H_0$ is some arbitrary reference Hamiltonian, we have the Dyson equation
\ee \ee
with $G^E_{0,ij} = \mel{i}{ \qty( H_0-E \Id )^{-1} }{j}$. with $G^E_{0,ij} = \mel{i}{ \qty( H_0-E \Id )^{-1} }{j}$.
Let us choose $H_0$ such that $\mel{ i }{ H_0 }{ j } = \mel{ i }{ P_i H P_i }{ j }$ for all $i$ and $j$. Let us choose $H_0$ such that $\mel{ i }{ H_0 }{ j } = \mel{ i }{ P_i H P_i }{ j }$ for all $i$ and $j$.
Then, The Dyson equation \eqref{eq:GE} becomes Then, the Dyson equation \eqref{eq:GE} becomes
\begin{multline} \begin{multline}
\mel{ i }{ \qty(H-E \Id)^{-1} }{ j } \mel{ i }{ \qty(H-E \Id)^{-1} }{ j }
= \mel{ i }{ P_i \qty(H-E \Id)^{-1} P_i }{ j } = \mel{ i }{ P_i \qty(H-E \Id)^{-1} P_i }{ j }
@ -741,24 +741,25 @@ Finally, the probabilistic expression reads
%----------------------------% %----------------------------%
\subsubsection{Energy estimator} \subsubsection{Energy estimator}
%----------------------------% %----------------------------%
To calculate the energy we introduce the following energy estimator To calculate the energy, we introduce the following estimator
\be \be
\cE(E) = \frac{ \mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ H\PsiT } } {\mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ \PsiT } }, \cE(E) = \frac{ \mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ H\PsiT } } {\mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ \PsiT } },
\label{calE} \label{calE}
\ee \ee
and search for the solution $E=E_0$ of $\cE(E)= E$. and search for the solution of the non-linear equation $\cE(E)= E$.
Using the spectral decomposition of $H$, we have Using the spectral decomposition of $H$, we have
\be \be
\label{eq:calE} \label{eq:calE}
\cE(E) = \frac{ \sum_i \frac{E_i c_i}{E_i-E}}{\sum_i \frac{c_i}{E_i-E}} \cE(E) = \frac{ \sum_i \frac{E_i c_i}{E_i-E}}{\sum_i \frac{c_i}{E_i-E}},
\ee \ee
with $c_i = \braket{ I_0 }{ \Phi_i } \braket{ \Phi_i}{ \PsiT }$, \titou{where $\Phi_i$ are eigenstates of $H$}. where $c_i = \braket{ I_0 }{ \Phi_i } \braket{ \Phi_i}{ \PsiT }$ and $\Phi_i$ are eigenstates of $H$, \ie, $H \ket{\Phi_i} = E_i \ket{\Phi_i}$.
It is easy to check that, in the vicinity of $E = E_0$, $\cE(E)$ is a linear function of $E - E_0$.
Therefore, in practice, we compute the value of $\cE(E)$ for several values of $E$, and fit these data using a linear or quadratic function of $E$ in order to obtain, via extrapolation, the exact value $E_0$. It is easy to check that, in the vicinity of the exact energy $E_0$, $\cE(E)$ is a linear function of $E - E_0$.
Therefore, in practice, we compute the value of $\cE(E)$ for several values of $E$, and fit these data using a linear or quadratic function of $E$ in order to obtain, via extrapolation, an estimate of $E_0$.
In order to have a precise extrapolation of the energy, it is best to compute $\cE(E)$ for values of $E$ as close as possible to the exact energy. In order to have a precise extrapolation of the energy, it is best to compute $\cE(E)$ for values of $E$ as close as possible to the exact energy.
However, as $E \to E_0$, both the numerators and denominators of Eq.~\eqref{eq:calE} diverge. However, as $E \to E_0$, both the numerators and denominators of Eq.~\eqref{eq:calE} diverge.
This is reflected by the fact that one needs to compute more and more $p$-components with an important increase of statistical fluctuations. \titou{This is reflected by the fact that one needs to compute more and more $p$-components with an important increase of statistical fluctuations.}
Thus, from a practical point of view, a trade-off has to be found between the \titou{quality} of the extrapolation and the amount of statistical fluctuations. Thus, from a practical point of view, a trade-off has to be found between the quality of the extrapolation and the amount of statistical fluctuations.
%Let us describe the 3 functions used here for the fit. %Let us describe the 3 functions used here for the fit.
@ -800,7 +801,7 @@ Let us consider the one-dimensional Hubbard Hamiltonian for a chain of $N$ sites
H= -t \sum_{\expval{ i j } \sigma} \hat{a}^+_{i\sigma} \hat{a}_{j\sigma} 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}, + U \sum_i \hat{n}_{i\uparrow} \hat{n}_{i\downarrow},
\ee \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. where $\expval{ i j }$ 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$. 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 In the site representation, a general vector of the Hilbert space can be written as
\be \be
@ -910,7 +911,7 @@ Let us begin with a small chain of 4 sites with $U=12$.
From now on, we shall take $t=1$. 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 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...$. 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$. In what follows $\ket{I_0}$ will be systematically chosen as one of the two N\'eel states, {\it e.g.} $\ket{I_0} = \ket{\uparrow,\downarrow, \uparrow,\ldots}$.
Figure \ref{fig1} shows the convergence of $H_p$ as a function of $p$ for different values of the reference energy $E$. Figure \ref{fig1} shows the convergence of $H_p$ as a function of $p$ for different values of the reference energy $E$.
We consider the simplest case where a single-state domain is associated to each state. We consider the simplest case where a single-state domain is associated to each state.
@ -949,7 +950,7 @@ Table \ref{tab1} illustrates the dependence of the Monte Carlo results upon the
The reference energy is $E=-1$. The reference energy is $E=-1$.
The first column indicates the various domains consisting of the union of some elementary domains as explained above. The first column indicates the various domains consisting of the union of some elementary domains as explained above.
The first line of the table gives the results when using a minimal single-state domain for all states, and the last one for the maximal domain containing the full linear space. The first line of the table gives the results when using a minimal single-state domain for all states, and the last one for the maximal domain containing the full linear space.
The size of the various domains is given in column 2, the average trapping time for the state $|I_0\rangle$ in the third column, and an estimate of the speed of convergence of the $p$-expansion for the energy in the fourth column. The size of the various domains is given in column 2, the average trapping time for the state $\ket{I_0}$ in the third column, and an estimate of the speed of convergence of the $p$-expansion for the energy in the fourth column.
To quantify the rate of convergence, we report the quantity, $p_{conv}$, defined as the smallest value of $p$ for which the energy is converged with six decimal places. To quantify the rate of convergence, we report the quantity, $p_{conv}$, defined as the smallest value of $p$ for which the energy is converged with six decimal places.
The smaller $p_{conv}$, the better the convergence is. The smaller $p_{conv}$, the better the convergence is.
Although this is a rough estimate, it is sufficient here for our purpose. Although this is a rough estimate, it is sufficient here for our purpose.
@ -1160,15 +1161,12 @@ This project has received funding from the European Union's Horizon 2020 --- Res
\label{app:A} \label{app:A}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
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 For the simplest case of a system containing only two states, $\ket{1}$ and $\ket{2}$, the fundamental equation given in Eq.~\eqref{eq:eqfond} simplifies to
\be \be
\begin{split}
\cI \cI
& = \mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ \Psi } = \mel{ I_0 }{ \qty(H-E \Id)^{-1} }{ \Psi }
\\ = \mel{ I_0 }{ P_0 \qty(H-E \Id)^{-1} P_0 }{ \Psi }
& = \mel{ I_0 }{ P_0 \qty(H-E \Id)^{-1} P_0 }{ \Psi }
+ \sum_{p=1}^{\infty} \cI_p, + \sum_{p=1}^{\infty} \cI_p,
\end{split}
\ee \ee
with with
\begin{multline} \begin{multline}