This page was generated from 12_dmrg_mixer.ipynb (download).

Why you need the Mixer in DMRG

This notebooks demonstrates why you need the “Mixer” in DMRG.

The problem: getting stuck in local minima

As you probably know, DMRG performs a local optimization of the MPS on one (Single-Site DMRG) or two (Two-Site DMRG) sites of the MPS, sweeping left and right through the system. During such a local update, say of sites i, i+1, the MPS outside of the one or two sites to be optimized is kept fix. The effective Hamiltonian is given in the product basis of Schmidt states left of i, local basis on sites i and i+1, and the Schmidt states right of i+1. The truncation to a maximal bond dimension of the MPS hence corresponds to a projection of the Hamiltonian into the (kept) Schmidt bases on the left and right. In particular when you have a long-range hopping (e.g. if you consider 2D cylinder models mapped to a 1D MPS), this can project out some relevant terms of the Hamiltonian, such that DMRG never gets to consider them in the local optimization: it’s “stuck in a local minima”.

The solution: enable the mixer to turn on a subspace expansion

The “Mixer” performs what is called a “subspace expansion” in the literature (see also the documention of the mixer classes). In other words, it artifically perturbs the SVD used to get Schmidt states, trying to include the relevant terms of the Hamiltonian into the basis of Schmidt states that you keep. The strength of the perturbation is controlled by the amplitude of the mixer.

The usual strategy in TeNPy is to enable the mixer initially and then switch it off after a couple of sweeps to perform final sweeps of convergence once we are in the optimal subspace of Schmidt states. If you already have a very good initial guess of the ground state, choose the initial amplitude fairly small (but still larger than the smallest Schmidt values you keep during truncation…)

The SingleSiteDMRG by construction stays in the subspace of a given bond-dimension MPS when you don’t use a mixer. Hence you should always use a mixer for SingleSiteDMRG, unless you start with an intial MPS with both the correct final bond dimension and distribution of charge blocks you want to have.

While the TwoSiteDMRG formally avoids this problem and is able to both grow the bond dimension and redistribute charges as necessary (at least if your hamiltonian has nearest-neighbor hoppings), it can still get stuck in a local minima as well. This is demonstrate below by the extreme example of a single long-range coupling over the whole chain when you consider periodic boundary conditions of a finite system.

Even for TwoSiteDMRG it is hence recommended to keep the mixer enabled for the initial sweeps.

[1]:
import numpy as np
import scipy
import matplotlib.pyplot as plt

import yaml
import time

np.set_printoptions(precision=5, suppress=True, linewidth=100)
plt.rcParams['figure.dpi'] = 150
[2]:
import tenpy
import tenpy.linalg.np_conserved as npc
from tenpy.algorithms import dmrg
from tenpy.algorithms.dmrg import TwoSiteDMRGEngine, SingleSiteDMRGEngine
from tenpy.networks.mps import MPS
from tenpy.models.spins import SpinModel
from tenpy.models.spins_nnn import SpinChainNNN2

from tenpy.algorithms.dmrg import SubspaceExpansion, DensityMatrixMixer

tenpy.tools.misc.setup_logging(to_stdout="WARNING")  # disable lengthy status reports

Finite DMRG with periodic boundary conditions

Tenpy supports DMRG using an (open-boundary) MPS on a model with periodic boundary condtions. Even for a nearest-neighbor model, for example here a spin-1 Heisenberg chain with L=100, the periodic boundary conditions get mapped to a single long-range coupling between site 0 and site 99.

[ ]:

[3]:
L = 100
model_params = {
    'L': L, 'bc_MPS': 'finite',
    'bc_x': 'periodic',
    'S': 1.,
}
M = SpinModel(model_params)
print(M.H_MPO.chi)
E_exact = -140.14840390392   # taken from https://arxiv.org/abs/cond-mat/0508709, Fig. 3

# NOTE: running this notebook with L=100 takes quite a bit of time.
# and you won't see progress due to disabled log.
# you can also experiment with smaller systems, but you need a different reference energy.
[2, 5, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 2]
[4]:
dmrg_params = {
    'mixer': None,
    'mixer_params': {
        'amplitude': 1.e-6,
        'decay': 1.,
        'disable_after': 18,
    },
    'max_sweeps': 20,
    'min_sweeps': 20,
    'chi_list': {
        0: 50,
        4: 100,
        8: 200,
       12: 400,
       16: 600,
    },
}
[5]:
def run_dmrg(DMRGEngine, mixer):
    psi1 = MPS.from_product_state(M.lat.mps_sites(), ['0.0'] * L, M.lat.bc_MPS)
    options = dmrg_params.copy()
    options['mixer'] = mixer
    t0 = time.time()
    eng = DMRGEngine(psi1, M, options)
    print('running ', DMRGEngine, 'with mixer ', mixer)
    eng.run()
    print("took {t:.1f} seconds".format(t=time.time()-t0))
    return eng
[6]:
engs = {}
for mixer in [None, DensityMatrixMixer, SubspaceExpansion]:
    eng = run_dmrg(TwoSiteDMRGEngine, mixer)
    key = "None" if mixer is None else mixer.__name__
    engs["TwoSiteDMRG, " + key] = eng
# takes long to run!
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  None
took 460.0 seconds
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.DensityMatrixMixer'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

took 497.0 seconds
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.SubspaceExpansion'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

took 564.9 seconds
[7]:
for mixer in [SubspaceExpansion, DensityMatrixMixer]:
    eng = run_dmrg(SingleSiteDMRGEngine, mixer)
    key = mixer.__name__
    engs["SingleSiteDMRG, " + key] = eng
# takes long to run!
running  <class 'tenpy.algorithms.dmrg.SingleSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.SubspaceExpansion'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:1852: UserWarning: H is zero in the given block, nothing to diagonalize.We just return the initial state again.
  warnings.warn("H is zero in the given block, nothing to diagonalize."

took 354.4 seconds
running  <class 'tenpy.algorithms.dmrg.SingleSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.DensityMatrixMixer'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:1852: UserWarning: H is zero in the given block, nothing to diagonalize.We just return the initial state again.
  warnings.warn("H is zero in the given block, nothing to diagonalize."

took 368.5 seconds
[8]:
plt.figure(dpi=150)
ax = plt.gca()
for key, eng in engs.items():
    eng.plot_sweep_stats(ax, 'sweep', 'E', y_exact=E_exact, label=key)
for s, chi in eng.options['chi_list'].items():
    ax.axvline(s+0.5, linestyle='--', color='k')
    ax.text(s+1, 0.03, f"$\\chi={chi:d}$")
ax.set_ylim(1.e-9, 1.e-1)
plt.legend(loc='lower left')
[8]:
<matplotlib.legend.Legend at 0x7fd9de0b5a00>
../_images/notebooks_12_dmrg_mixer_10_1.png
[ ]:

Reordering the MPS

The 'default' geometry of the chain used above just has one very-long range coupling:

  F -- E
 /      \             .--------------------.                .--------------------.
A        D    =>     /                      \      =>      /                      \
 \      /           A -- B -- C -- D -- E -- F            0 -- 1 -- 2 -- 3 -- 4 -- 5
  B -- C

 physical couplings           stretched                 chosen MPS orering

In this particular case, another even better solution is to re-order the MPS to use the order='folded' of the Chain. In this case, we get at most next-nearest neighbor interactions in the MPS language, which DMRG can handle much better than the one very-long-range coupling:

  F -- E
 /      \               F -- E -- D                1 -- 3 -- 5
A        D        =>    |         |      =>        |         |
 \      /               A -- B -- C                0 -- 2 -- 4
  B -- C

 physical couplings     stretched            chosen MPS orering

The example below shows that that TwoSiteDMRG does not get stuck even with disabled mixer. Nevertheless, I’d recommended to use the mixer in this case as well, since most couplings go beyond the local two sites you optimize at once.

[9]:
L = 100
model_params = {
    'L': L, 'bc_MPS': 'finite',
    'bc_x': 'periodic',
    'order': 'folded',
    'S': 1.,
}
M = SpinModel(model_params)
print(M.lat.order[:, 0])  # the order of MPS sites relative to original lattice indices
[ 0 99  1 98  2 97  3 96  4 95  5 94  6 93  7 92  8 91  9 90 10 89 11 88 12 87 13 86 14 85 15 84 16
 83 17 82 18 81 19 80 20 79 21 78 22 77 23 76 24 75 25 74 26 73 27 72 28 71 29 70 30 69 31 68 32 67
 33 66 34 65 35 64 36 63 37 62 38 61 39 60 40 59 41 58 42 57 43 56 44 55 45 54 46 53 47 52 48 51 49
 50]
[10]:
engs_reordered = {}
for mixer in [None, SubspaceExpansion]:
    eng = run_dmrg(TwoSiteDMRGEngine, mixer)
    key = "None" if mixer is None else mixer.__name__
    engs_reordered["TwoSiteDMRG, " + key] = eng
# takes long to run!
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  None
took 412.2 seconds
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.SubspaceExpansion'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

took 537.7 seconds
[11]:
plt.figure(dpi=150)
ax = plt.gca()
for key, eng in engs.items():
    if not key.startswith('Two'):
        continue
    eng.plot_sweep_stats(ax, 'sweep', 'E', y_exact=E_exact, label='default' + key, linestyle='--')
for key, eng in engs_reordered.items():
    eng.plot_sweep_stats(ax, 'sweep', 'E', y_exact=E_exact, label='folded ' + key, marker='s')
for s, chi in eng.options['chi_list'].items():
    ax.axvline(s+0.5, linestyle='--', color='k')
    ax.text(s+1, 0.03, f"$\\chi={chi:d}$")
ax.set_ylim(1.e-9, 1.e-1)
plt.legend(loc='lower left')
[11]:
<matplotlib.legend.Legend at 0x7fd9ddfd29d0>
../_images/notebooks_12_dmrg_mixer_15_1.png

Infinite DMRG for next-nearest neighbors

[12]:
L = 4
J, Jp = 0., 0.5  # purely next-nearest neighbors
model_params = {
    'Jx': J, 'Jy': J, 'Jz': 1.,
    'Jxp': Jp, 'Jyp': Jp, 'Jzp': Jp,
    'L': L, 'bc_MPS': 'infinite',
    'conserve': 'Sz',
    'S': 1.
}
M = SpinChainNNN2(model_params)
E_exact_inf = -0.7610925710077936 # from running 2-site DMRG with mixer at chi=2000
[13]:
dmrg_params = {
    'mixer': mixer,
    'mixer_params': {
        'amplitude': 1.e-6,
        'decay': 1.,
        'disable_after': 90,
    },
    'update_env': 0,
    'N_sweeps_check': 5,
    'max_sweeps': 100,
    'min_sweeps': 100,
    'chi_list': {
        0: 50,
        20: 100,
        40: 200,
        60: 400,
        80: 600,
    },
}
[ ]:

[14]:
engs_inf = {}
for mixer in [None, DensityMatrixMixer, SubspaceExpansion]:
    eng = run_dmrg(TwoSiteDMRGEngine, mixer)
    key = "None" if mixer is None else mixer.__name__
    engs_inf["TwoSiteDMRG, " + key] = eng
# takes long to run!
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  None
took 5.4 seconds
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.DensityMatrixMixer'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

took 115.4 seconds
running  <class 'tenpy.algorithms.dmrg.TwoSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.SubspaceExpansion'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

took 153.6 seconds
[15]:
for mixer in [SubspaceExpansion, DensityMatrixMixer]:
    eng = run_dmrg(SingleSiteDMRGEngine, mixer)
    key = mixer.__name__
    engs_inf["SingleSiteDMRG, " + key] = eng
# takes long to run!
running  <class 'tenpy.algorithms.dmrg.SingleSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.SubspaceExpansion'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

took 69.0 seconds
running  <class 'tenpy.algorithms.dmrg.SingleSiteDMRGEngine'> with mixer  <class 'tenpy.algorithms.dmrg.DensityMatrixMixer'>
WARNING : /home/johannes/postdoc/tenpy/tenpy/algorithms/dmrg.py:173: UserWarning: Mixer with decay=1. doesn't decay
  warnings.warn("Mixer with decay=1. doesn't decay")

WARNING : final DMRG state not in canonical form up to norm_tol=1.00e-05: norm_err=1.44e-04
WARNING : norm_err=1.44e-04 still too high after environment_sweeps, call psi.canonical_form()
took 101.4 seconds
[16]:
plt.figure(dpi=150)
ax = plt.gca()
for key, eng in engs_inf.items():
    eng.plot_sweep_stats(ax, 'sweep', 'E', y_exact=E_exact_inf, label=key)
for s, chi in eng.options['chi_list'].items():
    ax.axvline(s, linestyle='--', color='k')
    ax.text(s+1, 0.03, f"$\\chi={chi:d}$")
ax.set_ylim(1.e-9, 1.e-1)
plt.legend(loc='lower left')
[16]:
<matplotlib.legend.Legend at 0x7fd9dd8278e0>
../_images/notebooks_12_dmrg_mixer_22_1.png

Note: the fluctuations for the DensityMatrixMixer with SingleSIteDMRG are on the order of the mixer amplitude…

[ ]: