Heisenberg Model
[1]:
!pip install block2==0.5.2rc13 -qq --progress-bar off --extra-index-url=https://block-hczhai.github.io/block2-preview/pypi/
Introduction
In this tutorial we explain how to solve the Heisenberg model using the python interface of block2.
First, we have to define the “site-spin” of the model. The parameter heis_twos represents two times the spin in each site. Namely, for \(S = 1/2\) Heisenberg model, heis_twos = 1, for \(S = 1\) Heisenberg model, heis_twos = 2, etc. Note that arbitrary non-negative half-integer and integer \(S\) can be supported in block2.
Second, we can solve the model using the SU2 symmetry (SU2 mode) or U1 symmetry (SGB mode). The SU2 symmetry can be more efficient (for large S) and can generate states with well-defined total spin symmetry, but requires some additional rearrangement of the Hamiltonian terms.
The SGB Mode
The Hamiltonian is
We can solve this model using the following code:
[2]:
import numpy as np
from pyblock2.driver.core import DMRGDriver, SymmetryTypes
L = 100
heis_twos = 1
driver = DMRGDriver(scratch="./tmp", symm_type=SymmetryTypes.SGB, n_threads=4)
driver.initialize_system(n_sites=L, heis_twos=heis_twos, heis_twosz=0)
b = driver.expr_builder()
for i in range(L - 1):
b.add_term("PM", [i, i + 1], 0.5)
b.add_term("MP", [i, i + 1], 0.5)
b.add_term("ZZ", [i, i + 1], 1.0)
heis_mpo = driver.get_mpo(b.finalize(), iprint=0)
def run_dmrg(driver, mpo):
ket = driver.get_random_mps(tag="KET", bond_dim=250, nroots=1)
bond_dims = [250] * 4 + [500] * 4
noises = [1e-5] * 4 + [1e-6] * 2 + [0]
thrds = [1e-6] * 4 + [1e-8] * 4
return driver.dmrg(
mpo,
ket,
n_sweeps=8,
bond_dims=bond_dims,
noises=noises,
thrds=thrds,
cutoff=1E-24,
iprint=1,
)
energies = run_dmrg(driver, heis_mpo)
print('DMRG energy = %20.15f' % energies)
Sweep = 0 | Direction = forward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 32.844 | E = -44.1201487994 | DW = 1.12e-09
Sweep = 1 | Direction = backward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 38.582 | E = -44.1277383773 | DE = -7.59e-03 | DW = 1.21e-13
Sweep = 2 | Direction = forward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 41.882 | E = -44.1277383773 | DE = -4.01e-12 | DW = 7.28e-14
Sweep = 3 | Direction = backward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 44.918 | E = -44.1277383773 | DE = -1.34e-12 | DW = 2.79e-14
Sweep = 4 | Direction = forward | Bond dimension = 500 | Noise = 1.00e-06 | Dav threshold = 1.00e-08
Time elapsed = 60.141 | E = -44.1277398826 | DE = -1.51e-06 | DW = 4.18e-17
Sweep = 5 | Direction = backward | Bond dimension = 500 | Noise = 1.00e-06 | Dav threshold = 1.00e-08
Time elapsed = 79.821 | E = -44.1277398826 | DE = -1.92e-13 | DW = 3.84e-17
Sweep = 6 | Direction = forward | Bond dimension = 500 | Noise = 0.00e+00 | Dav threshold = 1.00e-08
Time elapsed = 94.075 | E = -44.1277398826 | DE = -1.35e-13 | DW = 1.61e-20
DMRG energy = -44.127739882610513
The SU2 Mode
To solve the model in the SU2 mode, we define the following (triplet) spin tensor operator:
Then we have
Note that in the above calculation, we have used the following CG factors:
[3]:
from block2 import SU2CG
print('<Jp=1,Jzp=+1; Jq=1,Jzq=-1|J=0,Jz=0> = ', SU2CG().cg(2, 2, 0, 2, -2, 0))
print('<Jp=1,Jzp=-1; Jq=1,Jzq=+1|J=0,Jz=0> = ', SU2CG().cg(2, 2, 0, -2, 2, 0))
print('<Jp=1,Jzp= 0; Jq=1,Jzq= 0|J=0,Jz=0> = ', SU2CG().cg(2, 2, 0, 0, 0, 0))
print(1 / 3 ** 0.5)
<Jp=1,Jzp=+1; Jq=1,Jzq=-1|J=0,Jz=0> = 0.5773502691896257
<Jp=1,Jzp=-1; Jq=1,Jzq=+1|J=0,Jz=0> = 0.5773502691896257
<Jp=1,Jzp= 0; Jq=1,Jzq= 0|J=0,Jz=0> = -0.5773502691896257
0.5773502691896258
So the Hamiltonian in SU2 notation is
We can solve this model using the following code:
[5]:
import numpy as np
from pyblock2.driver.core import DMRGDriver, SymmetryTypes
L = 100
heis_twos = 1
driver = DMRGDriver(scratch="./tmp", symm_type=SymmetryTypes.SU2, n_threads=4)
driver.initialize_system(n_sites=L, heis_twos=heis_twos, spin=0)
b = driver.expr_builder()
for i in range(L - 1):
b.add_term("(T+T)0", [i, i + 1], -np.sqrt(3) / 2)
heis_mpo = driver.get_mpo(b.finalize(adjust_order=False), iprint=0)
energies = run_dmrg(driver, heis_mpo)
print('DMRG energy = %20.15f' % energies)
Sweep = 0 | Direction = forward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 48.587 | E = -44.1275743858 | DW = 7.23e-12
Sweep = 1 | Direction = backward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 53.973 | E = -44.1277393752 | DE = -1.65e-04 | DW = 6.90e-14
Sweep = 2 | Direction = forward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 59.107 | E = -44.1277393752 | DE = -3.84e-13 | DW = 3.59e-16
Sweep = 3 | Direction = backward | Bond dimension = 250 | Noise = 1.00e-05 | Dav threshold = 1.00e-06
Time elapsed = 63.144 | E = -44.1277393752 | DE = -1.28e-13 | DW = 8.81e-17
Sweep = 4 | Direction = forward | Bond dimension = 500 | Noise = 1.00e-06 | Dav threshold = 1.00e-08
Time elapsed = 75.269 | E = -44.1277398850 | DE = -5.10e-07 | DW = 3.06e-22
Sweep = 5 | Direction = backward | Bond dimension = 500 | Noise = 1.00e-06 | Dav threshold = 1.00e-08
Time elapsed = 94.287 | E = -44.1277398850 | DE = -1.99e-13 | DW = 6.38e-23
Sweep = 6 | Direction = forward | Bond dimension = 500 | Noise = 0.00e+00 | Dav threshold = 1.00e-08
Time elapsed = 103.865 | E = -44.1277398850 | DE = 0.00e+00 | DW = 1.42e-23
DMRG energy = -44.127739885005937