MPO flux not equal when it is symmetric but contains non-zero flux

I would like to act a pair of electron annihilation operator on a state and calculate some properties based on it. Therefore, I create an MPO using the OpSum() and apply it to the MPS I want. However, I find that the obtained MPO does not obey the symmetry.

using ITensors,ITensorMPS
N = 2;
sites = siteinds("Fermion", N; conserve_qns=true);
os = OpSum();
os += "C", 1, "C", 2;
H = MPO(os, sites);
checkflux(H) # will raise error

and these cause some problem. I find that in the post

the OP mt similar problem and received the reply

Unfortunately it is a known issue that the OpSum to MPO constructor doesn’t work when the MPO has nonzero flux (in principle it could, but we haven’t gotten around to implementing it).

As the pair of electrons I am going to remove is rather complicated, it would be hard and computationally inefficient to write down the string of op directly. Meanwhile, I find that that OP received error when he run

cdagmpo=MPO(ampo,sites)

but I do not met this error when I constructed the MPO. So I guess the developers have made changes to the related code, and I would like to ask if there is a possibility that this will be implemented in the near future.
Besides, if I want to fix it myself, are there any suggestions?

The library I wrote ITensorMPOConstruction should be able to handle this just fine, but disclaimer I haven’t tested this behavior very much so let me know how it goes.

using ITensors,ITensorMPS, ITensorMPOConstruction
N = 2;
sites = siteinds("Fermion", N; conserve_qns=true);
os = OpSum();
os += "C", 1, "C", 2;
H = MPO_new(os, sites);
checkflux(H)
2 Likes

ThanksIt works for me.

2 Likes