This page was generated from 01_dmrg.ipynb.

A first finite DMRG Example

Like examples/d_dmrg.py, this notebook shows the basic interface for DMRG. It initalized the transverse field Ising model \(H = J XX + g Z\) at the critical point \(J=g=1\), and a finite MPS in the all-up state \(|\uparrow\cdots \uparrow\rangle\). It then runs DMRG to find the ground state. Finally, we look at the profile of the entanglement-cuts.

[1]:
import numpy as np
import scipy
import matplotlib.pyplot as plt
[2]:
import tenpy
import tenpy.linalg.np_conserved as npc
from tenpy.algorithms import dmrg
from tenpy.networks.mps import MPS
from tenpy.models.tf_ising import TFIChain
[3]:
L = 100
[4]:
model_params = {
    'J': 1. , 'g': 1.,  # critical
    'L': L,
    'bc_MPS': 'finite',
}

M = TFIChain(model_params)
Reading 'bc_MPS'='finite' for config TFIChain
Reading 'L'=100 for config TFIChain
Reading 'J'=1.0 for config TFIChain
Reading 'g'=1.0 for config TFIChain
[5]:
psi = MPS.from_lat_product_state(M.lat, [['up']])
[6]:
dmrg_params = {
    'mixer': None,  # setting this to True helps to escape local minima
    'max_E_err': 1.e-10,
    'trunc_params': {
        'chi_max': 100,
        'svd_min': 1.e-10,
    },
    'verbose': True,
    'combine': True
}
eng = dmrg.TwoSiteDMRGEngine(psi, M, dmrg_params)
E, psi = eng.run() # the main work; modifies psi in place
Reading 'combine'=True for config TwoSiteDMRGEngine
Subconfig 'trunc_params'=Config(<3 options>, 'trunc_params') for config TwoSiteDMRGEngine
Reading 'max_E_err'=1e-10 for config TwoSiteDMRGEngine
Reading 'mixer'=None for config TwoSiteDMRGEngine
================================================================================
sweep 1, age = 100
Energy = -126.9290280127265333, S = 0.3766760098039985, norm_err = 1.2e-01
Current memory usage 97.4 MB, time elapsed: 2.3 s
Delta E = nan, Delta S = nan (per sweep)
max_trunc_err = 0.0000e+00, max_E_trunc = 1.8474e-13
MPS bond dimensions: [2, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 2]
================================================================================
sweep 2, age = 100
Energy = -126.9618018107068309, S = 0.5108558254154341, norm_err = 7.0e-03
Current memory usage 99.5 MB, time elapsed: 5.6 s
Delta E = -3.2774e-02, Delta S = 1.3418e-01 (per sweep)
max_trunc_err = 0.0000e+00, max_E_trunc = 1.9895e-13
MPS bond dimensions: [2, 4, 8, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 16, 8, 4, 2]
================================================================================
sweep 3, age = 100
Energy = -126.9618767383882130, S = 0.5277643082687588, norm_err = 1.1e-05
Current memory usage 107.0 MB, time elapsed: 10.1 s
Delta E = -7.4928e-05, Delta S = 1.6908e-02 (per sweep)
max_trunc_err = 4.7789e-20, max_E_trunc = 2.8422e-13
MPS bond dimensions: [2, 4, 8, 16, 27, 35, 39, 43, 46, 49, 49, 52, 55, 60, 62, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 64, 61, 59, 60, 55, 51, 47, 44, 42, 42, 38, 35, 32, 23, 16, 8, 4, 2]
================================================================================
sweep 4, age = 100
Energy = -126.9618767396799086, S = 0.5277994024302479, norm_err = 5.8e-10
Current memory usage 115.4 MB, time elapsed: 16.9 s
Delta E = -1.2917e-09, Delta S = 3.5094e-05 (per sweep)
max_trunc_err = 1.3388e-18, max_E_trunc = 2.8422e-13
MPS bond dimensions: [2, 4, 8, 16, 22, 30, 34, 36, 41, 45, 47, 48, 52, 54, 57, 61, 61, 63, 64, 66, 67, 70, 71, 72, 74, 76, 76, 77, 79, 79, 80, 82, 83, 83, 83, 84, 84, 84, 84, 86, 86, 89, 88, 89, 89, 88, 88, 89, 88, 88, 88, 87, 86, 87, 87, 87, 87, 86, 86, 85, 85, 84, 84, 84, 83, 83, 83, 83, 81, 80, 79, 77, 77, 77, 75, 73, 71, 70, 67, 66, 64, 62, 61, 61, 57, 54, 52, 48, 47, 45, 41, 36, 34, 30, 22, 16, 8, 4, 2]
================================================================================
sweep 5, age = 100
Energy = -126.9618767396792691, S = 0.5277994025033066, norm_err = 2.1e-13
Current memory usage 115.4 MB, time elapsed: 20.4 s
Delta E = 6.3949e-13, Delta S = 7.3059e-11 (per sweep)
max_trunc_err = 6.5986e-20, max_E_trunc = 3.1264e-13
MPS bond dimensions: [2, 4, 8, 16, 22, 30, 34, 36, 41, 45, 47, 48, 52, 54, 57, 61, 61, 63, 64, 66, 67, 71, 71, 74, 77, 77, 78, 79, 81, 81, 81, 83, 83, 84, 85, 85, 85, 85, 86, 88, 89, 90, 90, 90, 90, 90, 90, 91, 91, 91, 90, 90, 90, 90, 90, 90, 90, 90, 89, 88, 86, 86, 85, 85, 85, 84, 83, 83, 81, 81, 81, 79, 78, 78, 77, 73, 71, 71, 67, 66, 64, 63, 61, 61, 57, 54, 52, 48, 47, 45, 41, 36, 34, 30, 22, 16, 8, 4, 2]
================================================================================
DMRG finished after 5 sweeps.
total size = 100, maximum chi = 91
================================================================================

Expectation Values

[7]:
# the ground state energy was directly returned by dmrg.run()
print("ground state energy = ", E)

# there are other ways to extract the energy from psi:
E1 = M.H_MPO.expectation_value(psi)  # based on the MPO
E2 = np.sum(M.bond_energies(psi))  # based on bond terms of H, works only for a NearestNeighborModel
assert abs(E-E1) < 1.e-10 and abs(E-E2) < 1.e-10
ground state energy =  -126.96187673967927
[8]:
# onsite expectation values

X = psi.expectation_value("Sigmax")
Z = psi.expectation_value("Sigmaz")
x = np.arange(psi.L)
plt.figure()
plt.plot(x, Z, label="Z")
plt.plot(x, X, label="X")   # note: it's clear that this is zero due to charge conservation!
plt.xlabel("site")
plt.ylabel("onsite expectation value")
plt.legend()
plt.show()
../_images/notebooks_01_dmrg_9_0.png
[9]:
# correlation functions

i0 = psi.L // 4  # for fixed `i`
j = np.arange(i0 + 1, psi.L)
XX = psi.term_correlation_function_right([("Sigmax", 0)], [("Sigmax", 0)], i_L=i0, j_R=j)
XX_disc = XX - X[i0] * X[j]
ZZ = psi.term_correlation_function_right([("Sigmaz", 0)], [("Sigmaz", 0)], i_L=i0, j_R=j)
ZZ_disc = ZZ - Z[i0] * Z[j]

dx = j - i0
plt.figure()
plt.plot(dx, XX_disc, label="X X")
plt.plot(dx, ZZ_disc, label="Z Z")
plt.xlabel(r"distance $|i-j|$")
plt.ylabel(r"disconnected correlations $\langle A_i A_j\rangle - \langle A_i \rangle\langle A_j\rangle$")
plt.legend()
plt.loglog()
plt.show()
../_images/notebooks_01_dmrg_10_0.png

We find power-law decaying correlations, as expected for a critical model. For a gapped model, we would expect exponentially decaying correlations.


We now look at the entanglement entropy. The transverse-field Ising model is critical at \(g=J\). Conformal field theory, Calabrese,Cardy 2004, predicts an entanglement entropy profile of

\[ \begin{align}\begin{aligned}S(l, L) = \frac{c}{6} \log\left(\frac{2L}{\pi a} \sin\left(\frac{\pi l}{L}\right)\right) + \textrm{const}\\where :math:`c=0.5` is the central charge, :math:`a` is the lattice spacing (we set :math:`a=1`), :math:`L` is the total size of the system and we consider subsystems of size :math:`l` and :math:`L-l` as left/right. Note that this yields the familiar :math:`S(L/2, L) = \frac{c}{6} \log(L) + \textrm{const}` for the half-chain entropy as a function of system size.\end{aligned}\end{align} \]
[12]:
S = psi.entanglement_entropy()

bonds = np.arange(0.5, psi.L-1)
plt.plot(bonds, S, 'o', label="S")

# preform fit to extract the central charge
central_charge, const, res = tenpy.tools.fit.central_charge_from_S_profile(psi)
fit = tenpy.tools.fit.entropy_profile_from_CFT(bonds + 0.5, psi.L, central_charge, const)
print(f"extraced central charge {central_charge:.5f} with residuum {res:.2e}")
print("(Expect central charge = 0.5 for the transverse field Ising model.)")
plt.plot(bonds, fit, label=f"fit with $c={central_charge:.3f}$")
plt.xlabel("bond")
plt.ylabel("entanglement entropy $S$")
plt.legend()
plt.show()
extraced central charge 0.50788 with residuum 1.95e-09
(Expect central charge = 0.5 for the transverse field Ising model.)
../_images/notebooks_01_dmrg_12_1.png
[ ]: