Error in using MPS loaded from mps.h5 file in julia

Question

Hello Miles,

I have been trying to calculate correlations from the ground state MPS saved as mps.h5 file. Here is a minimal example for free fermions:
Step 1: Obtaining and saving the ground state MPS:

using ITensors, HDF5, JLD, DelimitedFiles  
let
    N = 4
    numb_sweeps = 5
    sweeps = Sweeps(numb_sweeps)
    maxdim!(sweeps,10,20,100)
    cutoff!(sweeps,1e-10)

    sites = siteinds("Fermion",N)
    states = ["Occ","Occ","Emp","Emp"]
    psi0 = productMPS(sites,states)
    ampo = AutoMPO()
    for j = 1:N-1
            ampo += "C",j,"Cdag",j+1
            ampo += "C",j+1,"Cdag",j
    end
    ampo += "C",N,"Cdag",1
    ampo += "C",1,"Cdag",N
    H = MPO(ampo,sites)
    energy,psi = dmrg(H,psi0,sweeps)

    f = h5open("mps.h5","w")
    write(f,"mps",psi)
    close(f)
    return
end

Step 2: use the save mps.h5 file to calculate correlation

using ITensors, HDF5, JLD, DelimitedFiles
let
N = 4
sites = siteinds("Fermion",N)
ampo = AutoMPO()
ampo += "C",1,"Cdag",2  
corr = MPO(ampo,sites)

f = h5open("mps.h5","r")
psi=read(f,"mps",MPS)
close(f)
println("correlation = ",inner(psi,corr,psi))
return
end

Everything is pretty standard, but I got the error: "ERROR: LoadError: MPO A and MPS B must share site indices. On site 1, A has site indices IndexSet{2} (dim=2|id=347|"Fermion,Site,n=1")' (dim=2|id=347|"Fermion,Site,n=1") while B has site indices IndexSet{1} (dim=2|id=806|"Fermion,Site,n=1")."

It seems that when I define the correlation operator corr in the 2nd code, the site indices are regenerated and no longer match with that of the MPS obtained in the 1st code. If this is indeed the case, then is there a way around this issue?

Thank you for your time,
-Mason

Answer 1

Hi Mason,
Yes, your diagnosis is correct about what is causing the issue. There are two ways around this:

1) you can write the sites array to your HDF5 file too, and read it back in instead of making a new sites array. Then when you make your MPO using AutoMPO, it will use your old sites array and things will match.

2) we have a function called replace_siteinds(M::MPS,sites) which given an MPS M and an array of site indices sites, will replace all of the site indices of M (it figures these out for you) by the new ones in the provided array. So in your case you would just call psi = replace_siteinds(psi,sites) before using psi in your code. There is an in-place version too called replace_siteinds!.

Hope that helps!

Miles

Comments

Dear miles,
I ran into a similar issue and would like to store the SiteSet into hdf5, but it seems h5_write doesn't support SiteSet object in C++ v3.