Define QN on a sublattice

Question

Hi Miles,

I'm working on a model with conserved quantum number on a sublattice. To be specific, say we have the following Hamiltonian:
$$
H = S{i}^z S{i+1}^z + S_{2i+1}^x
$$
So the total Sz is conserved only for even sites. Is it possible to define a QN for this kind of case? Thanks a lot!

Best,
Chengshu

Answer 1

Hi Chengshu,
Interesting question. My guess is that it should be possible. The main thing that the IQTensor system requires is that some QN be conserved, but I think (or hope) that it should be flexible enough to accomodate a case like this.

Here's how I would approach it: make a new type of SiteSet, similar to SpinHalf (see itensor/mps/sites for examples; definitely use the very latest version of the master branch for this, post 2.0.11 because I've simplified the SiteSet system recently). For the even sites, define the IQIndex's in the same way as the SpinHalf class. For the odd sites, define the IQIndex for each odd site as follows:

IQIndex(nameint("site ",j), Index("0",2,Site), QN("Sz=",0))

In the above, I've merged the usual two IQIndex blocks for a S=1/2 spin (up and down) into a single block of size 2 (indicated by the single Index of size 2) and given this block a QN of Sz=0.

Hopefully this will work, but you may run into an issue somewhere. Please give it a try and let me know if something seems off with your results.

Miles

Comments

Thanks for your fast reply! I tried the following code in the head file

    SpinHalfSubSite(int n, Args const& args = Args::global())
    {
        if (n % 2 == 0)
            s = IQIndex{nameint("S=1/2 ",n),
               Index(nameint("Up ",n),1,Site),QN("Sz=",+1),
               Index(nameint("Dn ",n),1,Site),QN("Sz=",-1)};
        else
            s = IQIndex{nameint("site ",n), Index("0",2,Site), QN("Sz=",0)};
    }

and the Hamiltonian

    for (j = 1; j < N; ++j)
    {
        ampo += 1.,"Sz",j,"Sz",j+1;
        if (j % 2 == 1)
            ampo += g,"Sx",j;
    }
    ampo += 1.,"Sz",N,"Sz",1;

with the initial state

    for(int i = 1; i <= N; ++i)
    {
        if (i%2 == 1)
            state.set(i,"Up");
        else
            state.set(i,"Dn");
    }

It turns out that while it works with MPO & MPS, it does't work for IQMPS and IQMPO. When running it throws the error message "doTask not defined for task CalcDiv and storage type QMixed<Real>". What do you think is happening?

---

Hm, yes this is one of the tricky parts I anticipated. It may not be hard to fix though.

In the code defining the operators (the "op" function) look for where "Sx" is defined, and just comment out the line

    Op = mixedIQTensor(s,sP);

then re-run the code and see if that fixes it.

---

A better version of the code would check which site you are requesting "Sx" on and not let you get that operator on the even numbered sites. ("Sx" is not a well-defined IQTensor for the case of a usual S=1/2 site, which is what the mixed IQTensor business is about.)

---

Yes that fixes the problem. Thank you very much!

---

Great. Does it seem to work with DMRG and give reasonable energies too?

---

Yes. All three approaches (original head file, new head file+MPO+MPS, new head file+IQMPO+IQMPS) give the same result for DMRG energy.

---

Now I have another problem... It seems that after DMRG I can't regauge the MPS by (say) psi.position(1).

Test Hamiltonian:

    k = 0;

    ampo += -2.*t1,"Sz",k+1;
    ampo += -2.*t1,"Sz",k+2;
    ampo += -2.*t1,"Sz",k+3;
    ampo += -4.*t1,"Sx",k+1,"Sx",k+2;
    ampo += -4.*t1,"Sx",k+2,"Sx",k+3;

    ampo += -2.*h,"Sx",k+4;
    ampo += -2.*h,"Sx",k+5;
    ampo += -2.*h,"Sx",k+6;
    ampo += -2.*h,"Sx",k+7;
    ampo += -2.*h,"Sx",k+8;
    ampo += -2.*h,"Sx",k+9;

initial state:

    auto state = InitState(sites);
    for (j = 1; j <= N; ++j)
    {
        if (j % 2 == 1)
            state.set(j,"Up");
        else
            state.set(j,"Dn");
    }

with head file:

    SpinHalfSite(int n, Args const& args = Args::global())
    {
        if ((n % 9 >= 1) && (n % 9 <= 3))
            s = IQIndex{nameint("S=1/2 ",n),
               Index(nameint("Up ",n),1,Site),QN({1,2}),
               Index(nameint("Dn ",n),1,Site),QN({-1,2})};
        else
            s = IQIndex{nameint("site ",n), Index("0",2,Site),QN({0,2})};
    }

I think the problem might have something to do with the site setup, but can't find a way out. What do you think? Thanks!

---

What do you mean by you can't gauge it? Do you mean the code crashes?

---

Yes. When doing psi.position(1) I get the error message "Segmentation fault (core dumped)"

---

That's unfortunate. Can you please compile and re-run the code in debug mode? That way you may get a more helpful error message. Probably I will need to look at the code eventually, because this bug could be related to a deficiency currently in part of the IQTensor system (it's not a bug per se, but just a case that it doesn't handle well).

---

In the debug mode even the DMRG is not working, which is really surprising to me. The output now reads

Using approx/svd conversion of AutoMPO->IQMPO
nindex() = 1
i = 1
From line 191, file iqindex.cc

IQIndex::index arg out of range

IQIndex::index arg out of range
Aborted (core dumped)

---

So to make things simpler I tried the following "minimum" Hamiltonian
    
    ampo += -4.,"Sp",1,"Sm",2;
    ampo += -4.,"Sm",1,"Sp",2;

And it doesn't work with the following site-setting head file

    s = IQIndex{nameint("S=1/2 ",n),
               Index(nameint("Up ",n),1,Site),QN({1,2}),
               Index(nameint("Dn ",n),1,Site),QN({-1,2})};

I think the problem might be in this parity-conserving QN definition. What do you think about it? Thanks!

---

Problem solved! I should have defined the QN as {1,2} and {0,2}. Sorry about the mess!

---

Ok great! I meant to write back sooner but this is very good news, as I was initially worried there might be a bug in our IQTensor system. Let me know if you run into any more issues.