DMRG2 + MPO_ham_ising(..., cyclic=True) gives ~0 energy and almost-zero state norm

[TianHao-Yang-tim]

What happened?

Summary

When using quimb.tensor.DMRG2 with a periodic Ising Hamiltonian constructed via MPO_ham_ising(..., cyclic=True), the DMRG2 energy is extremely close to zero, while the same Hamiltonian converted to a dense matrix and diagonalized gives a reasonable negative ground state energy. The returned DMRG state also has (almost) zero norm.

This suggests that, at least in this combination (Ising + PBC + DMRG2), the resulting MPS/state and/or energy are not consistent with <psi|H|psi>.

---

Minimal reproducible example

import quimb as qu
import quimb.tensor as qtn
import numpy as np

L = 6
J = -1.0
h = 0.5

# Target Pauli model:
#   H = J sum σ^z_i σ^z_{i+1} + h sum σ^x_i, with PBC
#
# In terms of spin operators S^α with S^z = σ^z / 2, S^x = σ^x / 2,
# this corresponds to:
#   H = 4 J sum S^z_i S^z_{i+1} - 2 h sum S^x_i
#
# So we use:
#   j  = 4 * J
#   bx = -2 * h

H = qtn.MPO_ham_ising(
    L,
    j=4 * J,
    bx=-2 * h,
    S=1/2,
    cyclic=True,
)

# 1) Dense Hamiltonian ground state energy (works as expected)
H_dense = H.to_dense()
evals = np.linalg.eigvalsh(H_dense)
E_exact = float(evals[0])
print("Exact GS from dense diag =", E_exact)

# 2) DMRG2 on the same MPO
dmrg = qtn.DMRG2(
    H,
    bond_dims=[10, 20, 100, 100, 200],
    cutoffs=1e-10,
    which='SA',
)

dmrg.solve(tol=1e-6, verbosity=0)
print("E_quimb =", float(dmrg.energy))

gs0 = dmrg.state
print("state norm <psi|psi> =", float((gs0.H @ gs0)))
print("state.norm()=", gs0.norm())

### What did you expect to happen?

Exact GS from dense diag = -6.384694563603669
E_quimb = -1.6783801434103012e-09
state norm <psi|psi> = 2.947477296255865e-10
state.norm()=1.7168218592084225e-05

### Minimal Complete Verifiable Example

```Python

Relevant log output

Anything else we need to know?

I found that using MPS+DMRG can not solve PBC TFIM even for very small system. Hope can fix it. Ty!

Environment

quimb 1.11.2 python 3.12.12

[jcmgray]

Hi @TianHao-Yang-tim, yes a quirk of periodic DMRG is that it works less well for small system sizes (since the loop is too short). Can you see if the problem remains as you push the system size up? For very small system sizes you can just use the periodic hamtiltonian embedded into a non periodic MPO.

[jcmgray]

You can also use a autograd optimizer directly like so, this is a much more stable algorithm than periodic DMRG:

mps = qtn.MPS_rand_state(L, bond_dim=20, phys_dim=2, cyclic=True)

def norm_fn(mps):
    return mps / mps.expec(mps.H, compress=False)**0.5

def energy(mps, H):
    return mps.expec(H, mps.H, compress=False)

tnopt = qtn.TNOptimizer(
    mps,
    loss_fn=energy,
    norm_fn=norm_fn,
    loss_constants={'H': H},
)

tnopt.optimize(100)
# -6.384694562864 [best: -6.384694562864] :  40%|████      | 40/100 [00:01<00:02, 24.33it/s]

[TianHao-Yang-tim]

You can also use a autograd optimizer directly like so, this is a much more stable algorithm than periodic DMRG:

mps = qtn.MPS_rand_state(L, bond_dim=20, phys_dim=2, cyclic=True)

def norm_fn(mps): return mps / mps.expec(mps.H, compress=False)**0.5

def energy(mps, H): return mps.expec(H, mps.H, compress=False)

tnopt = qtn.TNOptimizer( mps, loss_fn=energy, norm_fn=norm_fn, loss_constants={'H': H}, )

tnopt.optimize(100)

-6.384694562864 [best: -6.384694562864] : 40%|████ | 40/100 [00:01<00:02, 24.33it/s]

Thanks for your explanation!! I find above method works well for small size but don't work when L increase to 1e2. So I wonder if you recommend me to use TNOptimizer solving small system while DMRG solving large system? I'd appreciate for your reply! (And DMRG test takes a bit more time so I'm not sure if it works right now. I would reply once I get the results! Still thanks a lot!)

[TianHao-Yang-tim]

Hi @jcmgray, I’ve followed your advice and obtained a much more stable result for the 50-qubit chain. Thank you so much!

I still have a couple of questions:

L = 50 ... psi_opt = tnopt.get_tn_opt() psi_opt.expec(H, psi_opt, compress=False)

energy = -47.84711788082012

I extract the final state from TNOptimizer, but I find that the energy computed from psi_opt does not coincide with the value reported by tnopt during optimization:

[best: -49.638034820557] : 13%|█▎ | 67/500 [20:33<2:12:54, 18.42s/it]

Is this discrepancy between the final MPS energy and the best loss value a common phenomenon, or am I using TNOptimizer incorrectly?

2. When I increase the system size to 50 or 100 sites, the DMRG method still fails to converge and the energy stays close to 0. At the same time, the runtime becomes extremely long. Do you have any suggestions for stabilizing DMRG in this regime (e.g. particular options or settings I should be using)?

[jcmgray]

Regarding TNOptimizer, it should match... I can't replicate that behavior here https://colab.research.google.com/drive/1DgZrGizDACL1C2vDSx_wMxvr01PnrzLi?usp=sharing.

Perhaps the optimizer is defaulting to using jax which defaults to single precision and there is some conditioning problem? You could check the norm of psi_opt, or as above specfiy autograd to compute gradients to be sure.

There is a PBC DMRG example notebook here with the various things to explore: https://quimb.readthedocs.io/en/latest/examples/ex_dmrg_periodic.html. There actually might be a bug affecting performance of this currently, you could install the dev versions of quimb and cotengra to check.

[TianHao-Yang-tim]

Thank you so much! It works after specifying 'autograd' with TNOptimizer.

I have tried PBC DMRG on Heissenberg model as notebook and get great result before testing Ising. Now I use dev versions and find it works well for large size.