DMRG of the Heisenberg model

Within the study of the quantum many-body problem in physics, the DMRG analysis of the Heisenberg model is an important theoretical example applying techniques of the density matrix renormalization group (DMRG) to the Heisenberg model of a chain of spins. This article presents the "infinite" DMRG algorithm for the $S=1$ antiferromagnetic Heisenberg chain, but the recipe can be applied for every translationally invariant one-dimensional lattice.

DMRG is a renormalization-group technique because it offers an efficient truncation of the Hilbert space of one-dimensional quantum systems.

The algorithm

The Starting Point

To simulate an infinite chain, starting with four sites. The first is the block site, the last the universe-block site and the remaining are the added sites, the right one is added to the universe-block site and the other to the block site.

The Hilbert space for the single site is ${\mathfrak {H}}$ with the base $\{|S,S_{z}\rangle \}\equiv \{|1,1\rangle ,|1,0\rangle ,|1,-1\rangle \}$ . With this base the spin operators are $S_{x}$ , $S_{y}$ and $S_{z}$ for the single site. For every block, the two blocks and the two sites, there is its own Hilbert space ${\mathfrak {H}}_{b}$ , its base $\{|w_{i}\rangle \}$ ( $i:1\dots \dim({\mathfrak {H}}_{b})$ )and its own operators $O_{b}:{\mathfrak {H}}_{b}\rightarrow {\mathfrak {H}}_{b}$ :

block: ${\mathfrak {H}}_{B}$ , $\{|u_{i}\rangle \}$ , $H_{B}$ , $S_{x_{B}}$ , $S_{y_{B}}$ , $S_{z_{B}}$
left-site: ${\mathfrak {H}}_{l}$ , $\{|t_{i}\rangle \}$ , $S_{x_{l}}$ , $S_{y_{l}}$ , $S_{z_{l}}$
right-site: ${\mathfrak {H}}_{r}$ , $\{|s_{i}\rangle \}$ , $S_{x_{r}}$ , $S_{y_{r}}$ , $S_{z_{r}}$
universe: ${\mathfrak {H}}_{U}$ , $\{|r_{i}\rangle \}$ , $H_{U}$ , $S_{x_{U}}$ , $S_{y_{U}}$ , $S_{z_{U}}$

At the starting point all four Hilbert spaces are equivalent to ${\mathfrak {H}}$ , all spin operators are equivalent to $S_{x}$ , $S_{y}$ and $S_{z}$ and $H_{B}=H_{U}=0$ . This is always (at every iterations) true only for left and right sites.

Step 1: Form the Hamiltonian matrix for the Superblock

The ingredients are the four block operators and the four universe-block operators, which at the first iteration are $3\times 3$ matrices, the three left-site spin operators and the three right-site spin operators, which are always $3\times 3$ matrices. The Hamiltonian matrix of the superblock (the chain), which at the first iteration has only four sites, is formed by these operators. In the Heisenberg antiferromagnetic S=1 model the Hamiltonian is:

$\mathbf {H} _{SB}=-J\sum _{\langle i,j\rangle }\mathbf {S} _{x_{i}}\mathbf {S} _{x_{j}}+\mathbf {S} _{y_{i}}\mathbf {S} _{y_{j}}+\mathbf {S} _{z_{i}}\mathbf {S} _{z_{j}}$

These operators live in the superblock state space: ${\mathfrak {H}}_{SB}={\mathfrak {H}}_{B}\otimes {\mathfrak {H}}_{l}\otimes {\mathfrak {H}}_{r}\otimes {\mathfrak {H}}_{U}$ , the base is $\{|f\rangle =|u\rangle \otimes |t\rangle \otimes |s\rangle \otimes |r\rangle \}$ . For example: (convention):

$|1000\dots 0\rangle \equiv |f_{1}\rangle =|u_{1},t_{1},s_{1},r_{1}\rangle \equiv |100,100,100,100\rangle$

$|0100\dots 0\rangle \equiv |f_{2}\rangle =|u_{1},t_{1},s_{1},r_{2}\rangle \equiv |100,100,100,010\rangle$

The Hamiltonian in the DMRG form is (we set $J=-1$ ):

$\mathbf {H} _{SB}=\mathbf {H} _{B}+\mathbf {H} _{U}+\sum _{\langle i,j\rangle }\mathbf {S} _{x_{i}}\mathbf {S} _{x_{j}}+\mathbf {S} _{y_{i}}\mathbf {S} _{y_{j}}+\mathbf {S} _{z_{i}}\mathbf {S} _{z_{j}}$

The operators are $(d*3*3*d)\times (d*3*3*d)$ matrices, $d=\dim({\mathfrak {H}}_{B})\equiv \dim({\mathfrak {H}}_{U})$ , for example:

$\langle f|\mathbf {H} _{B}|f'\rangle \equiv \langle u,t,s,r|H_{B}\otimes \mathbb {I} \otimes \mathbb {I} \otimes \mathbb {I} |u',t',s',r'\rangle$

$\mathbf {S} _{x_{B}}\mathbf {S} _{x_{l}}=S_{x_{B}}\mathbb {I} \otimes \mathbb {I} S_{x_{l}}\otimes \mathbb {I} \mathbb {I} \otimes \mathbb {I} \mathbb {I} =S_{x_{B}}\otimes S_{x_{l}}\otimes \mathbb {I} \otimes \mathbb {I}$

Step 2: Diagonalize the superblock Hamiltonian

At this point you must choose the eigenstate of the Hamiltonian for which some observables is calculated, this is the target state . At the beginning you can choose the ground state and use some advanced algorithm to find it, one of these is described in:

The Iterative Calculation of a Few of the Lowest Eigenvalues and Corresponding Eigenvectors of Large Real-Symmetric Matrices, Ernest R. Davidson; Journal of Computational Physics 17, 87-94 (1975)

This step is the most time-consuming part of the algorithm.

If $|\Psi \rangle =\sum \Psi _{i,j,k,w}|u_{i},t_{j},s_{k},r_{w}\rangle$ is the target state, expectation value of various operators can be measured at this point using $|\Psi \rangle$ .

Step 3: Reduce density matrix

Form the reduced density matrix $\rho$ for the first two block system, the block and the left-site. By definition it is the $(d*3)\times (d*3)$ matrix: $\rho _{i,j;i',j'}\equiv \sum _{k,w}\Psi _{i,j,k,w}\Psi _{i',j',k,w}$ Diagonalize $\rho$ and form the $m\times (d*3)$ matrix $T$ , which rows are the $m$ eigenvectors associated with the $m$ largest eigenvalues $e_{\alpha }$ of $\rho$ . So $T$ is formed by the most significant eigenstates of the reduced density matrix. You choose $m$ looking to the parameter $P_{m}\equiv \sum _{\alpha =1}^{m}e_{\alpha }$ : $1-P_{m}\cong 0$ .

Step 4: New block and universe-block operators

Form the $(d*3)\times (d*3)$ matrix representation of operators for the system composite of the block and left-site, and for the system composite of right-site and universe-block, for example:

$H_{B-l}=H_{B}\otimes \mathbb {I} +S_{x_{B}}\otimes S_{x_{l}}+S_{y_{B}}\otimes S_{y_{l}}+S_{z_{B}}\otimes S_{z_{l}}$

$S_{x_{B-l}}=\mathbb {I} \otimes S_{x_{l}}$

$H_{r-U}=\mathbb {I} \otimes H_{U}+S_{x_{r}}\otimes S_{x_{U}}+S_{y_{r}}\otimes S_{y_{U}}+S_{z_{r}}\otimes S_{z_{U}}$

$S_{x_{r-U}}=S_{x_{r}}\otimes \mathbb {I}$

Now, form the $m\times m$ matrix representations of the new block and universe-block operators, form a new block by changing basis with the transformation $T$ , for example:

{\begin{matrix}&H_{B}=TH_{B-l}T^{\dagger }&S_{x_{B}}=TS_{x_{B-l}}T^{\dagger }\end{matrix}}

At this point the iteration is ended and the algorithm goes back to step 1. The algorithm stops successfully when the observable converges to some value.