Iterative projection into ground state In this paper, it says that one can "determine the ground state wave function by applying the projection operator $\exp(-\tau H)$ to an arbitrary initial state $|\Psi\rangle$," and that in the limit of $\tau \to \infty$, we have that $\exp(-\tau H)|\Psi\rangle$ converges to the true ground state of the system. It then goes on to state that that this projection cannot be done in a single step since terms in the Hamiltonian $H$ do not commute with each other. My question is two-fold:

*

*How do we know that $\exp(-\tau H)|\Psi\rangle$ converges to the true ground state of the system in the limit $\tau \to \infty$? I don't recall this coming up in any of my QM courses, but I might be missing something very self-evident here.

*Why does a Hamiltonian with non-commuting terms mean that we cannot do this projection in a single step? I understand how non-commuting terms cause difficulty when we want to split $\tau$ into many small pieces and "build up" to the full $\tau$ by applying the projection operator iteratively. However, the way this sentence is written in the paper makes me think that there is some a priori reason for why $\exp(-\tau H)|\Psi\rangle$ does not give the true ground state, after a single application, when $H$ contains non-commuting terms.

Thanks for the help!
 A: Meng Cheng pretty much answered the question already, but it might still be helpful to elaborate on the meaning of "too hard to calculate $e^{-\tau H}$ for a finite $\tau$".
What hardness refers to is the computational complexity of the problem/proposed numerical algorithm. That is to say, given a problem of size $N$ (e.g. here in the paper referenced $N$ is the number of spins wished to simulate), how much memory and how much time do you need.
Why this is important is because we typically want to solve problems in which $N$ is big: obviously, computing $e^{-\tau H}$ for say a 1024x1024 Hamiltonian $H$ describing 10 spins is not a problem -- just type exp(H) in MATLAB or Mathematica. But ask MATLAB to exponentiate $H$ for 100 spins and your computer will crash.
Since exponentiation of a Hamiltonian typically results in a dense matrix, this has bad memory scaling: exponential in $N$. For scaling in time, see:
https://mathoverflow.net/questions/239073/what-is-the-time-complexity-of-the-matrix-exponential
Therefore, naively trying to compute $e^{-\tau H}$ is folly.
The key insights that the paper uses are:

*

*Exponential over finite time can be written as product of many exponentials over small times $e^{-\tau H} = e^{-d \tau H}e^{-d \tau H}\cdots e^{-d \tau H}$ (this is not an approximation).


*For small time exponentials $e^{-d \tau H}$, it can be approximated (Trotter decomposition) as in Eq. 4 of the paper: $e^{-d\tau H} = e^{-d\tau H_z} e^{-d\tau H_y} e^{-d\tau H_x} + O(d\tau^2)$


*Naively you might think, $H_z, H_y, H_x$ are still big matrices, so exponentiating them is still hard, like in the general case. But what's nice about each $H_z, H_y, H_x$ is that they are each made of sums of local terms which commute. So the full matrix $e^{-d\tau H_z}$ is itself simply made of products of local exponentiations $e^{-d\tau H_{z,local}}$ where $H_{z,local}$ only acts on 4 sites (assuming a nearest-neighbor model) and is thus easy to do.
Furthermore, note that we intend to apply $e^{-\tau H}$ to a wavefunction. So we don't actually have to form the full $e^{-d \tau H_z}$ etc. which would be humongous. Instead we just apply $e^{-d\tau H_{z,local}}$ to each local part of the wavefunction. This is very fast and memory efficient!
A: 

*

*How do we know that $\exp(−\tau H)|Ψ\rangle$ converges to the true ground state
of the system in the limit $\tau \to \infty$?


You can write the Hamiltonian in its diagonalized form
$$H=\sum_{n=0}^\infty E_n|E_n\rangle\langle E_n| \tag{1}$$
where the $E_n$ are the eigenvalues and $|E_n\rangle$
the eigenvectors of $H$.
Especially, $E_0$ is the lowest eigenvalue, and $|E_0\rangle$
is the ground state.
From this you can get
$$e^{-\tau H}=\sum_{n=0}^\infty e^{-\tau E_n}|E_n\rangle\langle E_n| \tag{2}$$
You can decompose any vector $|\Psi\rangle$
as a linear combination of the eigenvectors $|E_n\rangle$
$$|\Psi\rangle=\sum_{m=0}^\infty c_m|E_m\rangle \tag{3}$$
Now apply the operator (2) to the vector (3):
$$\begin{align}
   e^{-\tau H}|\Psi\rangle
=& \sum_{n=0}^\infty e^{-\tau E_n}|E_n\rangle\langle E_n|
   \sum_{m=0}^\infty c_m|E_m\rangle \\
=& \sum_{n=0}^\infty e^{-\tau E_n}|E_n\rangle
   \sum_{m=0}^\infty c_m\langle E_n|E_m\rangle \\
=& \sum_{n=0}^\infty e^{-\tau E_n}|E_n\rangle
   \sum_{m=0}^\infty c_m\delta_{nm} \\
=& \sum_{n=0}^\infty e^{-\tau E_n} c_n |E_n\rangle \\
=&\ e^{-\tau E_0} \sum_{n=0}^\infty e^{-\tau (E_n-E_0)} c_n |E_n\rangle \\
=&\ e^{-\tau E_0} \left(c_0|E_0\rangle
    + \sum_{n=1}^\infty \underbrace{e^{-\tau (E_n-E_0)}}_{\to 0} c_n |E_n\rangle\right) \\
\to&\ e^{-\tau E_0} c_0|E_0\rangle
\end{align}$$
So for $\tau \to \infty$ you finally got a vector proportional to the ground state $|E_0\rangle$.



*Why does a Hamiltonian with non-commuting terms mean
that we cannot do this projection in a single step?


I don't understand either about which non-commuting terms the authors
are talking in this statement.
Certainly $H$ commutes with $H$, so this not what they mean.
A: For your question 2, I think the authors were saying nothing more profound or deeper than "it is too hard to calculate $e^{-\tau H}$ directly for a finite $\tau$" because of the non-commutativity between different terms in $H$. This is why they split it into many small pieces, and for small $\tau$ one can apply the approximation in equation (4) of the paper.
A: Since I do not have access to the paper I cannot say much about 2. Regarding 1, the statement as it stands is trivially false, even if it is true when, say, roughly speaking.
What it is mathematically true is that, for every non-vanishing vector $\psi$,
$$\lim_{t\to +\infty} \frac{e^{-tH} \psi}{||e^{-tH}\psi||} = \psi_0$$
where  $\psi_0$ is a normalized eigenvector of $H$ with eigenvalue given by  $$\min \sigma(H) \cap \sigma_\psi(H)\:,$$ where $\sigma(H)$ is the spectrum of $H$ and $\sigma_\psi(H)$ is the minimal set of eigenvalues whose eigenvectors span $\psi$. (I am assuming to deal with Hermitian matrices and $H$ is such a matrix so everything is finite dimensional.)
In particular it is easy to pick out $\psi$ such that the limit above does not converge to the ground state of $H$ (assuming that the  eigenspace of minimal energy is one-dimensional).
If, up to phases, $\phi$ represents such ground state, it is sufficient to chose $\psi \perp \phi$.
The proof of my assertions is a straightforward re-adaptation of Thomas Fritsch's answer.
