• View
• Edit
• Rename
• History
• Print
• Backlinks

TDVP

In this tutorial, we will introduce the tools used to perform time-evolution simulations using the time-dependent variational principle (TDVP) algorithm for matrix product states. The main advantage of using TDVP instead of the time-evolving block decimation (TEBD) algorithm is that TDVP work with Hamiltonian with any arbitrary interaction lengths, while TEBD can only normally handle nearest-neighbour interactions

For our example calculation, we will look at a chain of spin-1/2s, where the initial state is a domain wall of spins pointing up and down in {$z$} direction, {$$|\Psi(t=0)\rangle = |\cdots\uparrow\uparrow\uparrow\uparrow\uparrow\downarrow\downarrow\downarrow\downarrow\downarrow\cdots\rangle.$$} We then time evolve this state under the XX Hamiltonian, which is a version of the anisotropic Heisenberg model where the spin coupling in the {$z$} dimension is zero, and the couplings in the {$x$} and {$y$} dimensions have equal magnitude, {$$H_\text{XX} = \sum_j (S^x_j S^x_{j+1} + S^y_j S^y_{j+1}) = \frac{1}{2} \sum_j (S^-_j S^+_{j+1} + S^+_j S^-_{j+1}).$$} From the second form of this Hamiltonian, we can see that the effect of each term is to swap neighbouring up and down spins (this model is mathematically equivalent to a free fermion model on a lattice, where up and down spins map to occupied and unoccupied sites, respectively). In fact, the dynamics of this magnetization of this state under this Hamiltonian is known exactly [1], placing the domain wall between sites {$-1$} and {$0$}, we have {$$S^z(t, x) = -\frac{1}{2} \sum_{n=-x}^x J^2_n(t), \qquad x \geq 0,$$} {$$S^z(t, x) = -S^z(t, -x-1), \qquad x < 0,$$} where {$J_n(t)$} is a Bessel function of the first kind.

For this simulation, we use the spin chain lattice file with {$U(1)$} symmetry (the Hamiltonian and initial state break the full {$SU(2)$} symmetry, so we have to use a smaller symmetry group), and we can generate the initial state for a finite system of 20 sites using mp-construct:

Now we are ready to perform time evolution on our state using mp-tdvp:

here, the -H (--Hamiltonian) flag sets the Hamiltonian to use for the time evolution, for the XX Hamiltonian we use lattice:H_J1t (other possible operators can be seen using the mp-lattice-info command), -w (--wavefunction) specifies the initial state, -o (--output) species the prefix for output files, which will have filenames of the form <prefix>.t<time>, -t (--timestep) specifies the size of the timestep and -n (--num-timesteps) the number of timesteps.

In the terminal output, we can see that the maximum bond dimension in the chain (shown as the MaxStates number in the output) increases during the simulation, and the program slows down as MaxStates increases. This is because the program uses bond expansion, which is necessary for TDVP to give accurate results. This is since during time evolution, the state will becomes more entangled over time, and a larger bond dimension is needed to represent the state to the same accuracy.

To tune the speed and accuracy of the evolution, we can use the -d (--eigen-cutoff) flag to specify the eigenvalue cutoff for truncating the state. Using a smaller cutoff will For example, if we use a cutoff of {$10^{-10}$}:

The bond dimension will be much smaller (roughly half when I tested it), the program will run faster, and if we plot the results, they will be visually identical (unless we plot the errors on a logarithmic scale: see the exercises at the bottom of the page).

After this command finishes running (or we cancel it with Ctrl+C), we will want to calculate the results of the time-evolution simulation: for simplicity, first, we’ll only examine the {$z$}-component of the spin to the left of the domain wall. To calculate expectation values from the full time evolution simulation, we could use a bash loop, such as

But if we have access to multiple CPU cores, then we can speed this up by calculating each expectation value in parallel using GNU Parallel:

where the -k (--keep-order) flag makes sure the output is printed in order. To plot these results, we can copy the output from the terminal, or just redirect them into an output file, which can be done by appending > data.dat right at the end of the line containing the for loop or parallel command, then we ise gnuplot (or your favourite plotting program):

As we can see, the numerical results lie right on top of the exact results.

(On a slightly technical note, there is a neat property of TDVP in that it conserves the energy of a time-independent Hamiltonian and any constants of motion which can be written as the sum of single-site operators (i.e. which commute with the single-site projector). However, this property only holds exactly when we do not truncate the state: with truncation, we will see some slight deviation from the original value (e.g. in the energy E in the program output), but this will deviation will generally be around the order of magnitude of the truncation threshold, which should be fairly small in the first place.)

To calculate the expectation values for the full chain, we can use a double loop, which can be done easily using GNU Parallel

We can then plot this in gnuplot using

If we continued to evolve for further times, the wavefront will reflect off the boundary, which will interfere with the results. Hence, if we wanted to look at the evolution for further times, we would need to redo the time evolution with a larger system size. A better alternative is to use infinite boundary condition (IBC) wavefunctions, where we work with an infinite-size system, but only perturb a finite-sized region, which can be dynamically expanded during the time evolution simulation as the wavefront spreads outwards.

Infinite TDVP

We can also consider a similar quench with a translation-invariant initial state {$$|\Psi(t=0)\rangle = |\cdots\uparrow\downarrow\uparrow\downarrow\uparrow\downarrow\uparrow\downarrow\uparrow\downarrow\cdots\rangle.$$} In this case, the expectation value of {$S^z$} is given by [2] {$$S^z(t, x) = -\frac{1}{2} J_0(2t).$$} We can generate this state as an iMPS with a two-site unit cell using mp-construct

We can then evolve this state using mp-itdvp and calculate the expectation values:

We can plot the results in gnuplot with

Something to note here is that the growth of the entanglement is much greater in the translation-invariant state than the domain wall state. This is since the domain wall is the only part of the latter state which is nonstationary in time evolution, while the state is close to an eigenstate far away from the domain wall wavefront. On the other hand, the translation-invariant state essentially has a domain wall between every neighbouring site, and so the state is nonstationary at every point at the lattice right from the beginning, and the entanglement increases even more as the domain wall wavefronts move through the chain. As a result, the domain wall state has a linear growth of bond dimension in time evolution, while the translation-invariant state has an exponential growth of bond dimension.

Exercises

1. We can plot the absolute error of the expectation value of {$S_z$} in the infinite simulation in gnuplot using

References