To increase accuracy of DMRG code

DMRG and Numerical Methods
1

Hi ,
I am using DMRG code for XXZ spin- 1/2 AFM spin chain in presence of longitudinal field. In AFM phase I expect degenerate ground state energy(E0 and E1). Whereas I am getting fluctuations in E1-E0 as a function of hz at few values below second critical magnetic field(hc2). Also as I am using randomMPS I get different answers for E0/N(when I write upto 10 decimal places, answer is same upto 6 decimal places ) at each run. How can I control the DMRG parameters?

Length of the chain

 int N = 32;

Hamiltonian used

 for(int j = 1; j < N; ++j)
           {
            ampo += 1,"Sx",j,"Sx",j+1;
            ampo += 1,"Sy",j, "Sy",j+1;
            ampo += 3, "Sz", j ,"Sz",j+1;
             }
            ampo += 1,"Sx",N,"Sx",1;
            ampo += 1,"Sy",N, "Sy",1;
            ampo += 3, "Sz", N ,"Sz",1;

 for(int j=1;j<=N;++j)
               {
                 ampo += hz, "Sz",j;
               }

DMRG parameters fixed

 auto sweeps=Sweeps(6);
 sweeps.maxdim()=10,20,30,40,50,80;
 sweeps.cutoff()=1E-6,1E-7,1E-8,1E-9;
 sweeps.niter()= 2 ;
 sweeps.noise()=1E-8,1E-10,0.0;

Ground state and corrosponding energy using dmrg

auto [E0,psi0] = dmrg(H,randomMPS(sites,2),sweeps,{"Quiet=",true});
2

Hi, so if I understand correctly you’re saying the energy is consistent each time up to about 6 decimal places? If so, that’s quite good and is a sign that things are converging well.

If you then want precision beyond that, there are two main things to do:

  1. increase the number of sweeps. So you could do 10 or even 20 sweeps, say.
  2. lower the cutoff, from 1E-9 to lower such as 1E-10 or even 1E-12
3

Using these two things it works and I get precision upto 10 decimal places for 20 sweeps and 1E-20 cutoff. Thank you so much.

4

Glad to hear it!

5

Hi,
I did same thing s for L=16, and anisotropy=2, and calculated local quantity <S^z_i> at each site at hz=0. Since anisotropy is in z-direction I expect some non-zero local quantity at each site but answer I got is,

#bond      hz                local_z         real_local_z                imag_local_z          
1                  0       (-6.0697280535e-12,0)   -6.0697280535e-12                   0 
2                   0      (5.809099729e-12,0)     5.809099729e-12                   0 
3                   0      (-5.8841074374e-12,0)   -5.8841074374e-12                   0
4                   0     (5.9981750469e-12,0)    5.9981750469e-12                   0
5                   0    (-6.0050848842e-12,0)   -6.0050848842e-12                   0
6                   0     (6.0311209152e-12,0)    6.0311209152e-12                   0
7                   0    (-6.0352209341e-12,0)   -6.0352209341e-12                   0
8                   0     (6.0492166831e-12,0)    6.0492166831e-12                   0
9                   0    (-6.0547677982e-12,0)   -6.0547677982e-12                   0
10                   0      (6.062005932e-12,0)     6.062005932e-12                   0
11                   0    (-6.0593585271e-12,0)   -6.0593585271e-12                   0
12                   0     (6.0610986716e-12,0)    6.0610986716e-12                   0
13                   0    (-6.0466466903e-12,0)   -6.0466466903e-12                   0
14                   0     (6.0527607232e-12,0)    6.0527607232e-12                   0
15                   0    (-6.0349919506e-12,0)   -6.0349919506e-12                   0
16                   0       (6.12790374e-12,0)      6.12790374e-12                   0
6

It could be for many reasons. You might have measured S^z incorrectly; your results might not be converged; the system might not develop order; there could be a bug in how you inputted your Hamiltonian. Those are some of the possibilities. Please check over each of those, and compare to exact results if available or check on small system sizes. Try adjusting various things like number of sweeps to see if the results change.

If you spend some more days on it and can’t figure it out, please post again here with more details about what you tried.

7

I tried different ways to modify answer,

  1. Set initial Neel state
    2.increase number of sweeps upto 28
    3.increase maxdim upto 1000
  2. lowered cutoff upto 1E-40
    5.checked for isotropic XXX model - all expectation values are getting of the order of 1E-12.
    6.took OBC

still answer is not getting what I expect.

7.tried for different length

#L               <S^z_1>
16               1.28e-11
32                3.83e-5               
48                0.094
64                0.35
128              0.43
  1. Decrease sweeps upto 1
#swp           <Sz_1>
22                    1.28e-11
5                      5.42e-05
4                       0.0027 
3                        0.039
2                         0.026
1                          0.38
8

Isn’t the true ground state of this model on a finite chain a superposition of Neel-like states? That would explain your observations, i.e. for a short chain with a large number of sweeps you manage to resolve the exponentially small energy difference between the true symmetric ground state and the symmetry broken states.

9

I agree with Vladimir. I did some experiments myself and here’s what I think is going on:

  1. in the infinite L limit, it’s known that for the XXZ model the Neel order only develops for \Delta > 1 where the Z term is \Delta S^z_j S^z_{j+1}.
  2. however, there is no guarantee for a finite-sized system whether this ordering will happen; it is only a limiting statement for large enough system size L
  3. therefore your system could be on the other side of the phase transition (or not) and so the only way to tell is to study various system sizes and extrapolate to the thermodynamic limit (alternatively you can try even larger \Delta)

When I tried \Delta=1.5 I found that for L=100 I could see some staggered \langle S^z_j\rangle order developing.

Some comments on your approach:

  1. with DMRG it is much, much better to use open boundary conditions when possible. So I would strongly recommend against periodic boundary conditions which could slow down your calculations and make the results less accurate

  2. it’s better to write the Hamiltonian as \frac{1}{2} S^+_j S^-_{j+1} + \frac{1}{2} S^-_j S^+_{j+1} + \Delta S^z_j S^z_{j+1} since the “Sy” operator involves complex numbers so will make ITensor use complex data and will slow down your calculation considerably. Also it won’t allow you to use quantum numbers to speed up your calculation if you use the “Sx” and “Sy” operators.

10

This isn’t exactly what I had in mind. Let me go into more detail:

  1. Whenever a Hamiltonian has some symmetry, its eigenstates can always be chosen to also be eigenstates of the symmetry operator. Hence, its non-degenerate eigenstates must be eigenstates of the symmetry, i.e. exhibit no symmetry breaking.

  2. The XXZ model has a Z2 symmetry. For infinite \Delta it becomes the classical Ising model, whose ground states are the two symmetry-related Neel states. Since they’re degenerate, one of them is as good a ground state as their symmetric or anti-symmetric superpositions. However, at any finite \Delta the S^+_jS^-_{j+1}​ terms introduce mixing between the Neel states, opening a gap and leaving a Z2 symmetric state with \langle S^z_i \rangle =0 as the unique ground state (and a Z2 anti-symmetric state as the first excited state) for any finite L and \Delta. This can be easily checked using exact diagonalization for L=16. The gap decays exponentially as system size is increased since mixing the Neel states requires an L order term in perturbation theory, thus allowing for spontaneous symmetry breaking in the thermodynamic limit.

  3. Specifically for L=16 and \Delta=2, the gap between the ground state and the first excited state is about 0.07, or 4e-3 per site. The energy difference between the ground state and the best symmetry broken state would be of that order of magnitude, i.e. well within the resolution of a high quality DMRG run. When the system size is increased or the run quality decreased, DMRG finds a symmetry broken state instead of the exact ground state.