How to translate from a state/density matrix formalism to matrix product state representation?

Question

From what I understand, MPS is just a simpler way to write out a state, compared to the density matrix. But how do I get those $A_i$ matrices? From all the examples I read, people just somehow "have" the matrices in their pocket already. Is there a way to generate them from some random state, say $|\Psi\rangle = \alpha|\uparrow\rangle^{\otimes N}+\beta|\downarrow\rangle^{\otimes N}$, or some random density matrix? Feel free to set $N$ to small numbers if it becomes too complicated, I just really want to see an example of how this works in general, not just a GHZ state! Thanks a lot!

Crosspost from: https://quantumcomputing.stackexchange.com/questions/17590/how-to-find-the-a-i-in-the-matrix-product-state-representation?noredirect=1#comment25914_17590

Answer 1

The review does a good job of answering the question, but just wanted to present the material using less forbidding notation and to make the connection to density matrices, entanglement entropy, and bond dimension (another quantity you might come across) more clear.

Suppose you are interested in calculating the MPS for the state $|\Psi\rangle = \sum_{\{s_1,...,s_N\}} \psi(s_1,...,s_N)|s_1,...,s_N\rangle$. For example, this might be a spin-1/2 system on N lattice sites. Because each of the qubits can take on one of $d$ values (in the spin example above, $d=2$), the total number of unique $\{s_1,...,s_N\}$ configurations is $d^N$. Thus, there are $d^N$ terms in this sum, and you would need to store $d^N$ variables--the coefficients $\psi(s_1,...,s_N)$--in order to reconstruct this state from basis elements.

Cool! Now let's begin forming the MPS. We start by cutting up the full state space into two subsystems, which we will call $A$ and $B$. We can do this anywhere on the chain; for now, let's be perfectly general and snip the chain apart on the $n^\text{th}$ bond.

Another thing we're going to do is come up with a better way of expressing the coefficients. Because there are $d^N$ terms in the sum over basis states, we can think of $\psi(s_1,s_2,...,s_N)$ as a vector with $d^N$ components. All we're going to do is reshape this vector into a matrix, $\psi \in \mathbb{R}^{d^n \times d^{N-n}}$.

After all this headache, we've just ended up with an equivalent way of writing the original state:

$$\begin{align}|\Psi\rangle_{AB} &= \sum_{\{s_1,...,s_n\},\{s_{n+1},...,s_N\}}\psi_{(s_1,...,s_n),(s_{n+1},...,s_N)}|s_1,...,s_n\rangle_A \otimes |s_{n+1},...,s_N\rangle_B\\ &\equiv \sum_{a,b}\psi_{a,b}|a\rangle_A |b\rangle_B \end{align}$$

In the last step I've let $a$ and $b$ correspond to a particular configuration of states in the $A$ subsystem and $B$ subsystem, respectively. This is just for ease of notation; $\psi$ is still a $d^n \times d^{N-n}$ matrix, and $\psi_{a,b}$ is still just some scalar coefficient.

Whenever we have a big matrix, a common thing to try is decomposing it. Generally speaking, a decomposition might allow us to do two things:

1. Observe patterns in the "data" (the coefficient matrix)
2. Reduce the rank of the state space (in this case, the rank of the Hilbert space in which $|\Psi\rangle$ lives). I'm not going to focus too much on this point, but suffice to say that it comes in very useful for many-body physics simulations, where the dimension of the Hilbert space corresponding to the full system is simply too vast to permit manipulation and storage.

Anyhow. One of the matrix decompositions we can try is called the Schmidt decomposition. I won't go into too much detail because it's covered in the review, but basically it allows us to rewrite the coefficient matrix in a useful way:

$$\begin{align} |\Psi\rangle_{AB} &= \sum_{a,b}\psi_{a,b}|a\rangle_A \otimes|b\rangle_B\\ &= \sum_{a,b} [U\sigma V^T]_{a,b}|a\rangle_A \otimes|b\rangle_B\\ &= \sum_{a,b}\sum_{k=1}^D U_{i,k}\sigma_k V_{k,j}^T|a\rangle_A \otimes|b\rangle_B\\ &= \sum_{k=1}^D \sigma_k U|k\rangle_A \otimes V|k\rangle_B\\ &\equiv \sum_{k=1}^D \sqrt{\lambda_k} |\alpha_k\rangle_A \otimes |\beta_k\rangle_B \end{align}$$ In the last line, I've just redefined some things.

Whew! That was a lot of work. Now for the good stuff.

1. To answer your original question: The algorithm I just described is the basic building block for constructing an MPS. If we want to work from right to left, for example, we could choose system $A$ to be sites $\{1,...,N-1\}$, and assign system $B$ to site $N$, as shown below:

I won't get into the messy details, but the basic ideas is that at each step you absorb $\sigma_k V_{k,j}^T$ into a single matrix, and then perform SVD again. The right-unitary matrix, $U$, is the $n^{th}$ matrix in the final product form of the MPS that you are familiar with.

2. To make the connection between density matrices and matrix product states: As is turns out, we can very easily compute the reduced density matrices $\rho_A$ and $\rho_B$ from the SVD we performed earlier. Namely,

$$\begin{align} \rho_A = \sum_{k=1}^D \lambda_k |\alpha_k\rangle\langle\alpha_k|\\ \rho_B = \sum_{k=1}^D \lambda_k |\beta_k\rangle\langle\beta_k| \end{align}$$

This is actually really cool! It tells us that the entanglement entropy (AKA Von Neumann entropy), which is only a function of the density matrix eigenvalues, is shared across the two subsystems. The entanglement entropy is (as the name suggests) a measurement of the entanglement between systems $A$ and $B$. The following statements are equivalent:

• The rank of the SVD, $D$, is equal to $1$
• There is only one nonzero eigenvalue $\lambda$ (and it's going to be equal to $1$)
• $\rho_A$ and $\rho_B$ represent pure states
• Entanglement entropy is zero
• The mutual information of systems $A$ and $B$ is zero
• Systems $A$ and $B$ are not entangled

Otherwise, take the converse of all the above statements, and you get a situation where systems $A$ and $B$ are entangled.

3. The rank of the Schmidt decomposition, $D$, is commonly known as the bond dimension, and corresponds to $\min(\dim A, \dim B)$. Although a lot of sources will refer to the bond dimension of a system, it's important to keep in mind that $D$ actually depends on which "bond" you're considering. In the 1D chain we've been considering, the bond dimension between systems $A = \{1,2,...,N-1\}$ and $B = \{N\}$ is just $d$; the bond dimension between systems $A = \{1,2,...,N/2\}$ and $B = \{N/2+1,...,N\}$ is $d^{N/2}$. This is still really big! Which brings us to my last point

4. A lot of the time, we will fix a maximum value for the bond dimension; say, $D=16$. Obviously, we lose a lot information, because the SVD is no longer exact. But this is an easy bullet to bite if it means being able to simulate complex many-body systems numerically.

Answer 2

Any state can be written as a matrix product state. There are systematic procedures to construct such a description, based on sequential SVDs, see e.g. Section 4.1.3 of this review.

On the other hand, this description is usually of interest if the resulting MPS description has much less parameters than the $2^N$ parameters needed to describe a general quantum state. This is, e.g., the case for the example you give.