The time-evolving block decimation (TEBD) algorithm is a numerical scheme used to simulate one-dimensional quantum many-body systems, characterized by at most nearest-neighbour interactions. It is dubbed Time-evolving Block Decimation because it dynamically identifies the relevant low-dimensional Hilbert subspaces of an exponentially larger original Hilbert space. The algorithm, based on the Matrix Product States formalism, is highly efficient when the amount of entanglement in the system is limited, a requirement fulfilled by a large class of quantum many-body systems in one dimension.
Contents
- Introduction
- Introducing the decomposition of State
- The Schmidt decomposition
- Building the decomposition of state
- Rationale
- The update of the decomposition
- One qubit gates acting on qubit k
- Two qubit gates acting on qubits k k1
- The computational cost
- The numerical simulation
- The Suzuki Trotter expansion
- Simulation of the time evolution
- Adaptive Schmidt dimension
- References
Introduction
There is nowadays a considerable interest in the field of quantum theory for computational methods well-suited to the physics of many-body systems. Considering the inherent difficulties of simulating general quantum many-body systems, the exponential increase in parameters with the size of the system, and correspondingly, the high computational costs, one solution would be to look for numerical methods that deal with special cases, where one can profit from the physics of the system. The raw approach, by directly dealing with all the parameters used to fully characterize a quantum many-body system is seriously impeded by the lavishly exponential buildup with the system size of the amount of variables needed for simulation, which leads, in the best cases, to unreasonably long computational times and extended use of memory. To get around this problem a number of various methods have been developed and put into practice in the course of time, one of the most successful ones being the quantum Monte Carlo method (QMC). Also the density matrix renormalization group (DMRG) method, next to QMC, is a very reliable method, with an expanding community of users and an increasing number of applications to physical systems.
When the first quantum computer will be plugged in and functioning, the perspectives for the field of computational physics will look rather promising, but until that day one has to restrict oneself to the mundane tools offered by classical computers. While experimental physicists are putting a lot of effort in trying to build the first quantum computer, theoretical physicists are searching, in the field of quantum information theory (QIT), for genuine quantum algorithms, appropriate for problems that would perform badly when trying to be solved on a classical computer, but pretty fast and successful on a quantum one. The search for such algorithms is still going, the best-known (and almost the only ones found) being the Shor's algorithm, for factoring large numbers, and Grover's search algorithm.
In the field of QIT one has to identify the primary resources necessary for genuine quantum computation. Such a resource may be responsible for the speedup gain in quantum versus classical, identifying them means also identifying systems that can be simulated in a reasonably efficient manner on a classical computer. Such a resource is quantum entanglement; hence, it is possible to establish a distinct lower bound for the entanglement needed for quantum computational speedups.
Guifré Vidal, then at the Institute for Quantum Information, Caltech, has recently proposed a scheme useful for simulating a certain category of quantum systems. He asserts that "any quantum computation with pure states can be efficiently simulated with a classical computer provided the amount of entanglement involved is sufficiently restricted" . This happens to be the case with generic Hamiltonians displaying local interactions, as for example, Hubbard-like Hamiltonians. The method exhibits a low-degree polynomial behavior in the increase of computational time with respect to the amount of entanglement present in the system. The algorithm is based on a scheme that exploits the fact that in these one-dimensional systems the eigenvalues of the reduced density matrix on a bipartite split of the system are exponentially decaying, thus allowing us to work in a re-sized space spanned by the eigenvectors corresponding to the eigenvalues we selected.
One can also estimate the amount of computational resources required for the simulation of a quantum system on a classical computer, knowing how the entanglement contained in the system scales with the size of the system. The classically (and quantum, as well) feasible simulations are those that involve systems only a trifle entangled—the strongly entangled ones being, on the other hand, good candidates only for genuine quantum computations.
The numerical method is efficient in simulating real-time dynamics or calculations of ground states using imaginary-time evolution or isentropic interpolations between a target Hamiltonian and a Hamiltonian with an already-known ground state. The computational time scales linearly with the system size, hence many-particles systems in 1D can be investigated.
A useful feature of the TEBD algorithm is that it can be reliably employed for time evolution simulations of time-dependent Hamiltonians, describing systems that can be realized with cold atoms in optical lattices, or in systems far from equilibrium in quantum transport. From this point of view, TEBD had a certain ascendance over DMRG, a very powerful technique, but until recently not very well suited for simulating time-evolutions. With the Matrix Product States formalism being at the mathematical heart of DMRG, the TEBD scheme was adopted by the DMRG community, thus giving birth to the time dependent DMRG [6], t-DMRG for short.
Around the same time, other groups have developed similar approaches in which quantum information plays a predominant role as, for example, in DMRG implementations for periodic boundary conditions [7], and for studying mixed-state dynamics in one-dimensional quantum lattice systems,. Those last approaches actually provide a formalism that is more general than the original TEBD approach, as it also allows to deal with evolutions with matrix product operators; this enables the simulation of nontrivial non-infinitesimal evolutions as opposed to the TEBD case, and is a crucial ingredient to deal with higher-dimensional analogues of matrix product states.
Introducing the decomposition of State
Let us consider a chain of N qubits, described by the function
where M is the on-site dimension.
The trick of TEBD is to re-write the coefficients
This may seem rather esoteric in the beginning but let's try take a look at how this decomposition is obtained and what it is good for.
For this we have to make use of another piece of mathematics called the Schmidt decomposition of a state, or singular value decomposition.
The Schmidt decomposition
Let us consider the state of a bipartite system
where
For example, the two-qubit state:
has the following SD:
with
On the other hand, the state:
is a product state:
Building the decomposition of state
At this point we probably know enough to try to see how we explicitly build the decomposition (let's call it D).
Consider the bipartite splitting
The process can be decomposed in three steps, iterated for each bond (and, correspondingly, SD) in the chain:
Step 1: express the
The vectors
Step 2: write each vector
Step 3: make the substitutions and obtain:
Repeating the steps 1 to 3, one can construct the whole decomposition of state D. The last
The Schmidt eigenvalues, are given explicitly in D:
The Schmidt eigenvectors are simply:
and
Rationale
Now, looking at D, instead of
Therefore, it is possible to take into account only some of the Schmidt coefficients (namely the largest ones), dropping the others and consequently normalizing again the state:
where
Let's get away from this abstract picture and refresh ourselves with a concrete example, to emphasize the advantage of making this decomposition. Consider for instance the case of 50 fermions in a ferromagnetic chain, for the sake of simplicity. A dimension of 12, let's say, for the
Even if the Schmidt eigenvalues don't have this exponential decay, but they show an algebraic decrease, we can still use D to describe our state
The update of the decomposition
One can proceed now to investigate the behaviour of the decomposition D when acted upon with one-qubit gates (OQG) and two-qubit gates (TQG) acting on neighbouring qubits. Instead of updating all the
One-qubit gates acting on qubit k
The OQGs are affecting only the qubit they are acting upon, the update of the state
Two-qubit gates acting on qubits k, k+1
The changes required to update the
Following Vidal's original approach,
The subspace J is spanned by the eigenvectors of the reduced density matrix
In a similar way, the subspace K is spanned by the eigenvectors of the reduced density matrix:
The subspaces
Using the same reasoning as for the OQG, the applying the TQG V to qubits k, k + 1 one needs only to update
We can write
where
To find out the new decomposition, the new
The square roots of its eigenvalues are the new
From the left-hand eigenvectors,
after expressing them in the basis
The computational cost
The dimension of the largest tensors in D is of the order
The numerical simulation
The numerical simulation is targeting (possibly time-dependent) Hamiltonians of a system of
It is useful to decompose
Any two-body terms commute:
The Suzuki-Trotter expansion
The Suzuki-Trotter expansion of the first order (ST1) represents a general way of writing exponential operators:
or, equivalently
The correction term vanishes in the limit
For simulations of quantum dynamics it is useful to use operators that are unitary, conserving the norm (unlike power series expansions), and there's where the Trotter-Suzuki expansion comes in. In problems of quantum dynamics the unitarity of the operators in the ST expansion proves quite practical, since the error tends to concentrate in the overall phase, thus allowing us to faithfully compute expectation values and conserved quantities. Because the ST conserves the phase-space volume, it is also called a symplectic integrator.
The trick of the ST2 is to write the unitary operators
where
Simulation of the time-evolution
The operators
since any two operators
The time-evolution can be made according to
For each "time-step"
Our goal is to make the time evolution of a state
It is rather troublesome, if at all possible, to construct the decomposition
i) construct the decomposition
ii) relate
iii) make an imaginary-time evolution towards the ground state of the Hamiltonian
or, alternatively, simulate an isentropic evolution using a time-dependent Hamiltonian, which interpolates between the Hamiltonian
iv)finally, make the time-evolution of the state
"Adaptive" Schmidt dimension
One thing that can save a lot of computational time without loss of accuracy is to use a different Schmidt dimension for each bond instead of a fixed one for all bonds, keeping only the necessary amount of relevant coefficients, as usual. For example, taking the first bond, in the case of qubits, the Schmidt dimension is just two. Hence, at the first bond, instead of futilely diagonalizing, let us say, 10 by 10 or 20 by 20 matrices, we can just restrict ourselves to ordinary 2 by 2 ones, thus making the algorithm generally faster. What we can do instead is set a threshold for the eigenvalues of the SD, keeping only those that are above the threshold.
TEBD also offers the possibility of straightforward parallelization due to the factorization of the exponential time-evolution operator using the Suzuki-Trotter expansion. A parallel-TEBD has the same mathematics as its non-parallelized counterpart, the only difference is in the numerical implementation.