Skip to content

imaginary_time_evolve

Computes the imaginary-time evolution,

\[\vert \psi(\tau) \rangle = e^{-(H - \delta) \tau} \vert \psi_0\rangle,\]

of a State \(\vert \psi_0 \rangle\) and a Hermitian operator \(H\) using an iterative algorithm. \(\delta\) here denotes a real number which can be chosen as the ground state energy \(\delta=E_0\) of \(H\).

Sources
imaginary_time_evolve.hpp
imaginary_time_evolve.cpp
imaginary_time_evolve.jl


Definition

The method is provided in two variants:

  1. Returning a new state while the input state remains untouched. This variant is safe to use and simple to code.

    State imaginary_time_evolve(OpSum const &H, State psi0, double time,
                                double precision = 1e-12, double shift = 0.);
    
    imaginary_time_evolve(ops::OpSum, psi0::State, time::Float64; precision::Float64 = 1e-12, 
                          shift::Float64=0.0)::State
    
  2. An inplace variant imaginary_time_evolve_inplace, where the input state is overwritten and contains the time evolved state upon exit. This version is more memory efficient than imaginary_time_evolve.

    void imaginary_time_evolve_inplace(OpSum const &H, State &psi0, double time,
                                       double precision = 1e-12, shift = 0.);
    
    imaginary_time_evolve_inplace(ops::OpSum, psi0::State, time::Float64; 
                                  precision::Float64 = 1e-12, shift::Float64=0.0)
    

Parameters

Name Description Default
H OpSum defining the hermitian operator \(H\) for time evolution
psi0 initial State \(\vert \psi_0 \rangle\) of the time evolution
time time \(\tau\) until which the state is evolved
precision accuracy of the computed time evolved state \(\vert \psi(\tau) \rangle\) 1e-12
shift the offset \(\delta\) when computing \(\vert \psi(t) \rangle = e^{-(H - \delta) \tau} \vert \psi_0\rangle\) 0.0

The routine calls the subroutine evolve_lanczos implementing a Lanczos algorithm to perform the evolution. This routine can also be called explicitly if more control is desired. Please also confer to the page evolve_lanczos for further details on the specifics of the algorithm. The parameter shift can be used to turn all eigenvalues of the matrix \(H - \delta \;\textrm{Id}\) positive whenever \(\delta < E_0\), where \(E_0\) denotes the ground state energy of \(H\).


Usage Example

int N = 8;
int nup = N / 2;
auto block = Spinhalf(N, nup);

// Define the nearest-neighbor Heisenberg model
auto ops = OpSum();
for (int i=0; i<N; ++i) {
  ops += Op("SdotS", {i, (i+1) % N});
}

// Compute ground state energy
double e0 = eigval0(ops, block);

auto psi0 = product_state(block, {"Up", "Dn", "Up", "Dn", "Up", "Dn", "Up", "Dn"});
double time = 1.0;
double precision = 1e-12;
auto psi = imaginary_time_evolve(ops, psi0, time, precision, e0);
imaginary_time_evolve_inplace(ops, psi0, time, precision, e0);
XDIAG_SHOW(isapprox(psi0, psi));
let
    N = 8
    nup = N รท 2
    block = Spinhalf(N, nup)

    # Define the nearest-neighbor Heisenberg model
    ops = OpSum()
    for i in 1:N
        ops += Op("SdotS", [i, mod1(i+1, N)])
    end

    # Compute ground state energy
    e0 = eigval0(ops, block)

    psi0 = product_state(block, ["Up", "Dn", "Up", "Dn", "Up", "Dn", "Up", "Dn"])
    time = 1.0
    psi = imaginary_time_evolve(ops, psi0, time,
                                precision=1e-12, shift=e0)
    imaginary_time_evolve_inplace(ops, psi0, time, precision=1e-12, shift=e0)
    @show isapprox(psi0, psi)
end