DMRG runs
In this notebook, we use the SimpleDMRGEngine
class from tenpy_toycodes/d_dmrg.py
to run DMRG and find MPS ground states.
[ ]:
# standard imports and cosmetics
import numpy as np
import matplotlib.pyplot as plt
np.set_printoptions(precision=5, suppress=True, linewidth=100, threshold=50)
plt.rcParams['figure.dpi'] = 150
In previous notebooks, we learned how to initialize SimpleMPS
…
[ ]:
from tenpy_toycodes.a_mps import SimpleMPS, init_FM_MPS, init_Neel_MPS
[ ]:
L = 12
[ ]:
psi_FM = init_FM_MPS(L=L, d=2, bc='finite')
print(psi_FM)
SigmaZ = np.diag([1., -1.])
print(psi_FM.site_expectation_value(SigmaZ))
[ ]:
def init_PM_MPS(L, bc='finite'):
"""Return a paramagnetic MPS (= product state with all spins pointing in +x direction)"""
d = 2
B = np.zeros([1, d, 1], dtype=float)
B[0, 0, 0] = 1./np.sqrt(2)
B[0, 1, 0] = 1./np.sqrt(2)
S = np.ones([1], dtype=float)
Bs = [B.copy() for i in range(L)]
Ss = [S.copy() for i in range(L)]
return SimpleMPS(Bs, Ss, bc=bc)
… and how to initialize the TFIModel
.
[ ]:
from tenpy_toycodes.b_model import TFIModel
g = 1.2
model = TFIModel(L=L, J=1., g=g, bc='finite')
print("<H_bonds> = ", psi_FM.bond_expectation_value(model.H_bonds))
print("energy:", np.sum(psi_FM.bond_expectation_value(model.H_bonds)))
# (make sure the model and state have the same length and boundary conditions!)
For small enough system size \(L \lesssim 16\), you can compare the energies to exact diagonalization in the full Hilbert space (which is exponentially expensive!):
[ ]:
from tenpy_toycodes.tfi_exact import finite_gs_energy
if L <= 16:
energy_exact = finite_gs_energy(L=L, J=1., g=g)
[ ]:
The DMRG algorithm
The file tenpy_toycodes/d_dmrg.py
implements the DMRG algorithm. It can be called like this:
[ ]:
from tenpy_toycodes.d_dmrg import SimpleDMRGEngine, SimpleHeff2
[ ]:
chi_max = 15
psi = init_FM_MPS(model.L, model.d, model.bc)
eng = SimpleDMRGEngine(psi, model, chi_max=chi_max, eps=1.e-10)
for i in range(10):
E_dmrg = eng.sweep()
E = np.sum(psi.bond_expectation_value(model.H_bonds))
print("sweep {i:2d}: E = {E:.13f}".format(i=i + 1, E=E))
print("final bond dimensions: ", psi.get_chi())
mag_x = np.mean(psi.site_expectation_value(model.sigmax))
mag_z = np.mean(psi.site_expectation_value(model.sigmaz))
print("magnetization in X = {mag_x:.5f}".format(mag_x=mag_x))
print("magnetization in Z = {mag_z:.5f}".format(mag_z=mag_z))
if model.L <= 16:
E_exact = finite_gs_energy(L=model.L, J=model.J, g=model.g)
print("err in energy = {err:.3e}".format(err=E - E_exact))
Exercise: read d_dmrg.py
Read the code of tenpy_toycodes/d_dmrg.py
and try to undertstand the general structure of how it works.
[ ]:
Exercise: measure correlation functions
Just looking at expectation values of local operators is not enough. Also measure correlation functions \(\langle Z_{L/4} Z_{3L/4} \rangle\).
[ ]:
[ ]:
Exercise: DMRG runs
Try running DMRG for various different parameters:
Change the bond dimension
chi
and truncation threasholdeps
.Change the system size
L
Change the model parameter
g
(at fixed \(J=1\)) to both the ferromagnetic phase \(g<J\) and the paramagnetic phase \(g>1\).Change the initial state.
[ ]:
Exercise: Phase diagram
To map out the phase diagram, it can be convenient to define a function that just runs DMRG for a given model. Fill in the below template. Use it obtain and plot the energy, correlations and magnetizations for different \(L=16, 32, 64\) against \(g\).
[ ]:
def run_DMRG(model, chi_max=50):
print(f"runnning DMRG for L={model.L:d}, g={model.g:.2f}, bc={model.bc}, chi_max={chi_max:d}")
raise NotImplementedError("TODO: this is an Exercise!")
return psi
[ ]:
results_all = {}
# TIP: comment this out after the first run to avoid overriding your data!
[ ]:
L = 8
mag_X = []
mag_Z = []
E = []
corr = []
max_chi = []
gs = [0.1, 0.5, 0.8, 0.9, 1.0, 1.1, 1.2, 1.5]
for g in gs:
model = TFIModel(L=L, J=1., g=g, bc='finite')
psi = run_DMRG(model)
mag_X.append(np.mean(psi.site_expectation_value(model.sigmax)))
mag_Z.append(np.mean(psi.site_expectation_value(model.sigmaz)))
E.append(np.sum(psi.bond_expectation_value(model.H_bonds)))
corr.append(psi.correlation_function(model.sigmaz, model.L//4, model.sigmaz, model.L * 3 // 4))
max_chi.append(max(psi.get_chi()))
results_all[L] = {
'g': gs,
'mag_X': mag_X,
'mag_Z': mag_Z,
'E': E,
'corr': corr,
'max_chi': max_chi
}
[ ]:
key = 'corr'
plt.figure()
for L in results_all:
res_L = results_all[L]
plt.plot(res_L['g'], res_L[key], marker='o', label="L={L:d}".format(L=L))
plt.xlabel('g')
plt.ylabel(key)
plt.legend(loc='best')
[ ]:
[ ]:
[ ]:
[ ]:
Infinite DMRG
The given DMRG code also works with bc='infinite'
boundary conditions of the model and state. The given SimpleDMRG
code also allows to run infinite DMRG, simply by replacing the bc='finite'
for both the model and the MPS.
Look at the implementation of
d_dmrg.py
(anda_mps.py
) to see where the differences are.
Again, we can compare to analytic calculations possible for the Transverse Field Ising model.
[ ]:
from tenpy_toycodes.tfi_exact import infinite_gs_energy
The L
parameter now just indices the number of tensors insite the unit cell of the infinite MPS. It has to be at least 2
, since we optimize 2 tensors at once in our DMRG code. Note that we now use the mean
to calculate densities of observables instead of extensive quantities:
[ ]:
model = TFIModel(L=2, J=1., g=0.8, bc='infinite') # just change bc='infinite' here
chi_max = 10
psi = init_FM_MPS(model.L, model.d, model.bc)
eng = SimpleDMRGEngine(psi, model, chi_max=chi_max, eps=1.e-7)
for i in range(10):
E_dmrg = eng.sweep()
E = np.mean(psi.bond_expectation_value(model.H_bonds))
#if i % 10 == 9:
print("sweep {i:2d}: E/L = {E:.13f}".format(i=i + 1, E=E))
print("final bond dimensions: ", psi.get_chi())
mag_x = np.mean(psi.site_expectation_value(model.sigmax))
mag_z = np.mean(psi.site_expectation_value(model.sigmaz))
print("magnetization density in X = {mag_x:.5f}".format(mag_x=mag_x))
print("magnetization density in Z = {mag_z:.5f}".format(mag_z=mag_z))
E_exact = infinite_gs_energy(model.J, model.g)
print("err in energy = {err:.3e}".format(err=E - E_exact))
Exercise: Infinite DMRG
Try running the infinite DMRG code. How many sweeps do you now need to perform to converge now? Does it depend on g
and chi
?
[ ]:
[ ]:
[ ]:
From the exercise, you should see that you need significantly more sweeps when the correlation length gets larger, in particular at the critical point! At the critical point, the correlation length (of infinite MPS) scales as \(\xi \propto \chi^\kappa\), so you again need more sweeps to converge at larger bond dimensions.
[ ]:
[ ]:
Advanced exercises for the experts (and those who want to become them ;-) )
Obtain the ground state of the transverse field ising model at the critical point with DMRG for large
L
. Try to plot the corrlation function as a function ofj-i
. What form does it have? Is an MPS a good ansatz for that?
[ ]:
Compare running DMRG and imaginary time evolution with TEBD from
tenpy_toycodes/c_tebd.py
for various parameters ofL
,J
,g
, andbc
. Which one is faster? Do they always agree?
[ ]:
[ ]:
[ ]:
[ ]:
[ ]:
[ ]: