Writing the following Hamiltonian as MPOs for DMRG calculation

ITensor Julia Questions
, , ,
1

Dear Admin,

I’m trying to find the ground state for the following Hamiltonian via DMRG (for a spin-half chain of 2N sites):

H = \sum_{k=1}^{N} c_k (Z_{2k-1} - Z_{2k})^2 + \frac{J}{2N} ( H_{odd/Left} - H_{even/Right} )^2

where H_{odd/Left} are Hamiltonians defined on odd sites and H_{even/Right} are Hamiltonians defined on even sites. The (H_{odd/Left} - H_{even/right})^2 will be highly nonlocal and for now that’s fine for me. For simplicity I fix all c_k coupling to be the same

My ITensors Julia code is in the following form:

using ITensors, ITensorMPS

# Some previous code block turning ampol, ampor,
# OpSum() objects, into Hamiltonians HLeft & HRight

HLeft = MPO(ampol, sites);
HRight = MPO(ampor, sites);

# Here we construct the MPO of (H_L - H_R)^2

H_diff = HLeft - HRight;
H_sq_mpo = apply(H_diff, H_diff);
H_sq_mpo = truncate(H_sq_mpo; maxdim=maxdim);
H_sq_mpo = (BigJ/(2*N))*H_sq_mpo;

# Here we construct sum of the MPO of (Z_L - Z_R)^2

c_k_sqrt = sqrt(100);

OL = OpSum();
OR = OpSum();
add!(OL, c_k_sqrt, "Z", 1);
add!(OR, c_k_sqrt, "Z", 2);
OL = MPO(OL, sites);
OR = MPO(OR, sites);
O_diff = OL - OR;
O_sq = apply(O_diff, O_diff);
O_sq = truncate(O_sq; maxdim=maxdim);
Obscoupling = O_sq;

for j in 2:N
    temp1 = OpSum()
    temp2 = OpSum()
    add!(temp1, c_k_sqrt, "Z", 2j-1)
    add!(temp2, c_k_sqrt, "Z", 2j)
    temp1 = MPO(temp1, sites)
    temp2 = MPO(temp2, sites)
    temp = temp1 - temp2
    temp_sq = apply(temp, temp)
    temp_sq = truncate(temp_sq; maxdim=maxdim)
    Obscoupling = Obscoupling + temp_sq
end;

# Then do DMRG
H_list = [Obscoupling, H_sq_mpo]

# After initializing random MPS, perform DMRG
energy, gs = dmrg(H_list, psi0;nsweeps=nsweeps, maxdim=maxdims,cutoff=cutoffs)

Did I implement both terms correctly in my DMRG calculation? The calculation proceeds bug-free but I’m not sure if what I have written is the term I want.

Best,

Brian.

2

Hi Brian,
I would suggest trying the following: convert your MPO’s to full tensors and add them by doing Hfull = prod(Obscoupling) + prod(H_sq_mpo). Of course you’ll only be able to do this for very small systems.

Then use either a power method or a full diagonalization to get the ground state and check against exact results. Could that help you to be sure you inputted H correctly?