Reload MPO

[cfengno1]

Dear Huanchen,

I am trying to save Hamiltonian MPO and reload it in the successive DMRG calculations. In the initial run, I save the MPO using the following code:
``
mpo = driver.get_mpo(b.finalize(adjust_order=True, fermionic_ops="cdCD"),

algo_type=MPOAlgorithmTypes.FastBipartite, iprint=1)

mpo.save_data('Ham.bin')
``

and then in the successive run, I reload the MPO by:
``
from block2.sany import *
mpo=MPO(0)

mpo.load_data('Ham.bin')

energy = driver.dmrg(mpo, mps, n_sweeps=sweep_num, bond_dims=[M] * sweep_num,
noises=[0], thrds=[1e-10] * 8, dav_max_iter=30, iprint=2, forward=forward, sweep_start=sweep_start)
``

But I get the following error:
python: /superscratch/fchen/packages/block2/src/core/expr.hpp:988: std::shared_ptr<block2::OpExpr<S> > block2::load_expr(std::istream&) [with S = SAnyT<>; FL = double; std::istream = std::basic_istream<char>]: Assertion `false' failed.

Can you help?

Best,
Feng

[hczhai]

In order to reproduce the problem and save time guessing the context, please provide the complete Python script for saving and loading MPO, respectively.

[cfengno1]

Here they are.

I first run the following code to produce MPO and initial MPS:

from pyblock2.driver.core import DMRGDriver, SymmetryTypes, MPOAlgorithmTypes
from block2.sany import *
import numpy as np
LX, LY = 4, 6
L = LX * LY
J1 = 1
t1 = 3
J2 = 0.06
t2 = t1 * (J2 ** 0.5)
theta = np.pi/12
Nh = 2
NE = L - Nh
TWO_S = 0
NTHREADS = 32
Mem = 100
driver = DMRGDriver(scratch="./tmp", symm_type=SymmetryTypes.SAny | SymmetryTypes.CPX , n_threads=NTHREADS, stack_mem=Mem << 30)

driver.set_symmetry_groups("U1Fermi", "U1") # electron number, doublon number, total spin, total spin
Q = driver.bw.SX

site_basis, site_ops = [], []

for k in range(L):
    basis = [(Q(0, 0), 1), (Q(1, 1), 1), (Q(1, -1), 1)] # [0 up down]
    ops = {
        "": np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1]]),        # identity
        "c": np.array([[0, 0, 0], [1, 0, 0], [0, 0, 0]]),       # alapha+
        "d": np.array([[0, 1, 0], [0, 0, 0], [0, 0, 0]]),  # alpha
        "C": np.array([[0, 0, 0], [0, 0, 0], [1, 0, 0]]),       # beta+
        "D": np.array([[0, 0, 1], [0, 0, 0], [0, 0, 0]]),       # beta
        "z": np.array([[0, 0, 0], [0, 0.5, 0], [0, 0, -0.5]]),       # Sz
        "p": np.array([[0, 0, 0], [0, 0, 1], [0, 0, 0]]),       # S+
        "m": np.array([[0, 0, 0], [0, 0, 0], [0, 1, 0]]),       # S-
        "n": np.array([[0, 0, 0], [0, 1, 0], [0, 0, 1]]),  # electron number
    }
    site_basis.append(basis)
    site_ops.append(ops)

driver.initialize_system(n_sites=L, vacuum=Q(0, 0), target=Q(NE, TWO_S), hamil_init=False)
driver.ghamil = driver.get_custom_hamiltonian(site_basis, site_ops)
b = driver.expr_builder()
la = np.loadtxt("Lattice.dat",delimiter=",",dtype=int)
for i in range(0,la.shape[0]):
    bd = la[i][0]
    s1 = la[i][1] - 1
    s2 = la[i][2] - 1
    if bd == 1:
        th = theta
        t = t1
        J = J1
    else:
        th = 0
        t = t2
        J = J2

    b.add_term("pm", [s1, s2], J/2 * np.exp(-2j*th))
    b.add_term("pm", [s2, s1], J/2 * np.exp(2j*th))
    b.add_term("zz", [s1, s2], J)
    b.add_term("nn", [s1, s2], -J/4)
    b.add_term("cd", [s1, s2], -t*np.exp(-1j*th))
    b.add_term("cd", [s2, s1], -t*np.exp(1j*th))
    b.add_term("CD", [s1, s2], -t*np.exp(1j*th))
    b.add_term("CD", [s2, s1], -t*np.exp(-1j*th))

mpo = driver.get_mpo(b.finalize(adjust_order=True, fermionic_ops="cdCD"), algo_type=MPOAlgorithmTypes.FastBipartite, iprint=1)
mpo.save_data('Ham.bin')
loc_h = np.round(np.linspace(0,L-1,Nh)).astype(int)
Occu = np.ones(L)
Occu[loc_h] = 0

mps = driver.get_random_mps(tag="KET", bond_dim=1000, nroots=1, occs=Occu)

Then I reload MPO and MPS to do DMRG:

from pyblock2.driver.core import DMRGDriver, SymmetryTypes, MPOAlgorithmTypes
from block2.sany import *
import numpy as np

LX, LY = 4, 6
L = LX * LY
J1 = 1 
t1 = 3 
J2 = 0.06
t2 = t1 * (J2 ** 0.5)
theta = np.pi/12
Nh = 2 
NE = L - Nh
TWO_S = 0 
NTHREADS = 40
Nmem = 120 

driver = DMRGDriver(scratch="./tmp", symm_type=SymmetryTypes.SAny | SymmetryTypes.CPX , n_threads=NTHREADS, stack_mem=Nmem << 30, restart_dir="./MPS", clean_scratch=False)

driver.set_symmetry_groups("U1Fermi", "U1") # electron number, doublon number, total spin, total spin
Q = driver.bw.SX

site_basis, site_ops = [], []

for k in range(L):
    basis = [(Q(0, 0), 1), (Q(1, 1), 1), (Q(1, -1), 1)] # [0 up down]
    ops = { 
        "": np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1]]),        # identity
        "c": np.array([[0, 0, 0], [1, 0, 0], [0, 0, 0]]),       # alapha+
        "d": np.array([[0, 1, 0], [0, 0, 0], [0, 0, 0]]),  # alpha
        "C": np.array([[0, 0, 0], [0, 0, 0], [1, 0, 0]]),       # beta+
        "D": np.array([[0, 0, 1], [0, 0, 0], [0, 0, 0]]),       # beta
        "z": np.array([[0, 0, 0], [0, 0.5, 0], [0, 0, -0.5]]),       # Sz
        "p": np.array([[0, 0, 0], [0, 0, 1], [0, 0, 0]]),       # S+
        "m": np.array([[0, 0, 0], [0, 0, 0], [0, 1, 0]]),       # S-
        "n": np.array([[0, 0, 0], [0, 1, 0], [0, 0, 1]]),  # electron number
    }   
    site_basis.append(basis)
    site_ops.append(ops)

driver.initialize_system(n_sites=L, vacuum=Q(0, 0), target=Q(NE, TWO_S), hamil_init=False)
driver.ghamil = driver.get_custom_hamiltonian(site_basis, site_ops)
mpo=MPO(0)
mpo.load_data('Ham.bin')

mps = driver.load_mps(tag="KET", nroots=1)

sweep_start = 0
sweep_num = 12

M = 8000
energy = driver.dmrg(mpo, mps, n_sweeps=sweep_num, bond_dims=[M] * sweep_num,
        noises=[0], thrds=[1e-10] * 8, dav_max_iter=30, iprint=2, forward=forward, sweep_start=sweep_start)
print("DMRG energy = %20.15f (per site = %10.6f)" % (energy, energy / L))

This gives the error:
python: /superscratch/fchen/packages/block2/src/core/expr.hpp:988: std::shared_ptr<block2::OpExpr<S> > block2::load_expr(std::istream&) [with S = SAnyT<>; FL = double; std::istream = std::basic_istream<char>]: Assertion `false' failed.

[hczhai]

Thanks for providing the details.

The problem is that the symmetry mode used in from block2.sany import *; mpo = MPO(0) (real custom symmetry) does not match the actual symmetry mode in your script symm_type=SymmetryTypes.SAny | SymmetryTypes.CPX (complex custom symmetry).

To solve the problem, you can remove from block2.sany import * and change mpo = MPO(0) to mpo = driver.bw.bs.MPO(0).