0 votes
asked by (200 points)


I am trying to calculate the density of up-spin particles in an hexagonal chain with 6 sites considering only first-neighbor tight binding interactions.

When I set the chain to have one electron-up only, the density is equally distributed between all 6 sites. However, when more electrons are added, the density is no longer uniform, regardless of the DMRG parameters I choose (maxm, niter, cutoff, noise and number of sweeps). This is very weird, as I don't see any reason the density shouldn't be the same in all sites.

Also, the results obtained can change depending on the position I choose to initialize the electrons. For instance, if I choose


the density of "up" states in each site is different than in the case


Here is an example of my complete code for the case with two electrons. Am I doing something wrong? Isn't the expected behavior that electrons should be evenly distributed between all sites?

#include "itensor/all.h"
#include <iostream>
#include <fstream>

using namespace itensor;

int main()

    int N = 6;
    double t = 1.0;
    int nsweeps = 20;
    int maxm = 200;
    double cutoff = 1E-15;
    int niter = 2;
    double noise = 1E-5;

    auto sites = Hubbard(N);

    auto ampo = AutoMPO(sites);

   //Tight-binding: first neighbors in the linear chain (sites 1-2 and 5-6)
    for(int j = 1; j <= N-1; j+=4)
       ampo += -t,"Cdagup",j,"Cup",j+1;
       ampo += -t,"Cdagup",j+1,"Cup",j;
       ampo += -t,"Cdagdn",j,"Cdn",j+1;
       ampo += -t,"Cdagdn",j+1,"Cdn",j;

   //Tight-binding: second neighbors in the linear chain (sites 1-3, 2-4, 3-5, 4-6)
    for(int j=1; j <= N-2; j+=1)
       ampo += -t,"Cdagup",j,"Cup",j+2;
       ampo += -t,"Cdagup",j+2,"Cup",j;
       ampo += -t,"Cdagdn",j,"Cdn",j+2;
       ampo += -t,"Cdagdn",j+2,"Cdn",j;

    //Use IQMPO or MPO to generate the hamiltonian
    auto H = IQMPO(ampo);

    //Define initial wavefunction MPS for IQMPS
    auto state1 = InitState(sites);


    //Create IQMPS
    auto psi = IQMPS(state1);

    //Define DMRG sweeps
    auto sweeps = Sweeps(nsweeps);
    sweeps.maxm() = maxm/20,maxm/10,maxm/2,maxm/2,maxm;
    sweeps.cutoff() = cutoff;
    sweeps.niter() = niter;
    sweeps.noise() = noise;

    auto energy = dmrg(psi,H,sweeps,{"Quiet",true});

    //Calculating density in each site

    std::ofstream file1;
    file1.open ("nup.txt");

    for(int ii = 1; ii <= N; ++ii)
      auto Nup = sites.op("Nup",ii);
      ITensor ket = psi.A(ii);
      ITensor bra = dag(prime(ket,Site));
      auto nup = (bra*Nup*ket).real();
      file1 << nup << std::endl;


    return 0;
    } //End of main    

P.S.: Since I am working with an hexagon, I needed to number each site and displace them in a linear chain in order to conduct the calculations, which ended up creating second neighbor interactions in the linear chain. That is why the code includes first and second neighbors interactions for the linear chain. To make things clearer, here is how I chose to number the sites in the hexagon:

          /       \
         1         5      
         |         |
         2         6
          \       /

1 Answer

0 votes
answered by (70.1k points)

Hi rmagaldi,
Thanks for the question. So I think I may have figured out what's going wrong. The good news is that it does not appear to be the DMRG optimization itself, but the definition of the Hamiltonian.

What I did to test this is to observe that your system is just a periodic ring, so one can make the Hamiltonian in a simpler way:

for(int j = 1; j <= N-1; j+=1)
    ampo += -t,"Cdagup",j,"Cup",j+1;
    ampo += -t,"Cdagup",j+1,"Cup",j;
    ampo += -t,"Cdagdn",j,"Cdn",j+1;
    ampo += -t,"Cdagdn",j+1,"Cdn",j;
ampo += +t,"Cdagup",1,"Cup",N;
ampo += +t,"Cdagup",N,"Cup",1;
ampo += +t,"Cdagdn",1,"Cdn",N;
ampo += +t,"Cdagdn",N,"Cdn",1;

This is just a nearest-neighbor hopping with an additional "long bond" connecting site 1 to site N. Note that the coupling has the opposite sign for the long bond. Choosing this sign gives a uniform density, regardless of the initial state. Choosing a -t coupling gives the behavior you reported.

So probably if one works through this carefully, the reason the +t is the correct coupling is because really one should think of the long bond as going from N to 1 with a -t, which is equivalent to going from 1 to N with a +t, due to the fact that C and Cdagger anti commute.

Hope that helps!


commented by (200 points)
Hello Miles,

Thank you for your answer.

I understand your solution and it does work perfectly for the case of 2 electrons in a 6-site chain. However, I found two limitations to this approach:

1) It only works for an even number of electrons (2, 4 or 6). If I choose an odd number of electrons, the density once again becomes non-uniform. This happens even for only one electron, which is very weird because it is the simplest case and even my original approach was able to get an uniform density in this situation.

2) I intend to work with larger hexagonal chains, and there seems to be no clear best way to number the sites to make the Hamiltonian "simpler", as you did for the 6-site case. I tried out some possible numbering layouts for a 10 site chain, and the density ended up not being uniform in any of my tests.

Do you have any ideas or suggestions to solve these issues?
Welcome to ITensor Support Q&A, where you can ask questions and receive answers from other members of the community.

Formatting Tips:
  • To format code, indent by four spaces
  • To format inline LaTeX, surround it by @@ on both sides
  • To format LaTeX on its own line, surround it by $$ above and below
  • For LaTeX, it may be necessary to backslash-escape underscore characters to obtain proper formatting. So for example writing \sum\_i to represent a sum over i.
If you cannot register due to firewall issues (e.g. you cannot see the capcha box) please email Miles Stoudenmire to ask for an account.

To report ITensor bugs, please use the issue tracker.