Implementing local gauge symmetries

Question

Hi everybody,

I'm trying to set up the lattice Schwinger model which is a 1+1D theory with a U(1) gauge symmetry. I have defined a mixed site set consisting of a spin half site (the matter field) and a boson site (the gauge field) with a cutoff. The Hamiltonian is as given in Eq. 2 of this paper: https://arxiv.org/pdf/1312.6654.pdf.

In addition to a global U(1) symmetry (which I know how to implement natively in iTensor), the physical Hilbert space satisfies a Gauss' law that commutes with the Hamiltonian, is site-dependent, and has the following form:

G_i = E_(i+1) - E_i + (-1)^i n_i = 0

where n_i is the matter number operator at site i and E_i is the gauge electric field (which is diagonal in the site basis). Is there a way to implement this sort of local, site-dependent symmetry to further constrain the block structure of the tensors? I haven't been able to find a way, and the only resource I have found so far is an earlier post discussing how to implement separate symmetries on two sublattices of the lattice.

While the Schwinger model is special in that you can eliminate the gauge degree of freedom, being able to implement gauge symmetries would provide a large speedup when we move to more complicated systems, so it would be great to know whether or not this is possible.

Any pointers would be greatly appreciated!

Thanks.

Comments

Hi, do you think this article is relevant for what you're asking? Their method looks very nice, but may be too different of a setting or system from the one you're studying; I'm not an expert enough to be sure:

https://journals.aps.org/prx/abstract/10.1103/PhysRevX.4.041024

---

I am also working this Hamiltonian, so I would also be very interested in figuring out how to do this. Out of curiosity, I see that the tj.jl code in the Julia version of the code has this segment:

if conserve_sz && conserve_nf
  em = QN(("Nf",0,-1),("Sz", 0)) => 1
  up = QN(("Nf",1,-1),("Sz",+1)) => 1
  dn = QN(("Nf",1,-1),("Sz",-1)) => 1
  return [Index(em,up,dn;tags="Site,tJ,n=$n") for n=1:N]

If I am understanding this correctly, this code shows the presence of multiple conserved quantities at once (please correct me if I am wrong). Does this mean that, for the Schwinger model or other cases with a local gauge symmetry, one could effectively define a conserved quantity on each site, and then stitch them together for force all of them to be conserved at once? It seems to be that this would be a nightmare to write up for a general number of spatial sites, but I am wondering whether this would at least work in theory, if one had the patience to write out this code for a specific number of spatial sites.

---

Hi Sujay, I had missed this comment - I think maybe I didn't get a notification email from the forum?

It's a good question what you're proposing. In principle, the idea of putting a different QN type on each site would work *except* that the way our system currently works, these all get accumulated together each time tensors are added or contracted. So after a while, it could be the case that more than 4 different QN values accumulate together, which goes beyond the hard limit we currently have set.

One approach could be to lift this limit by using dynamically sizeable containers inside of ITensors.jl to hold the different QN values. This could certainly work, however I believe the code that handles the arithmetic on the QN objects might get rather slow so we'd have to do some testing and profiling to be sure. It might work well, but also might get rather involved.

Answer 1

Here is an idea that I had (written using the Julia code). I still highly doubt it will work, for some reasons I will explain, but I have some hope that something like it may work. In particular, I think understanding why it doesn't work would still be helpful. Anyway, here is the code:

ferm_full = QN("Sz",-1) => 1
ferm_vac = QN("Sz",0) => 1
antiferm_full = QN("Sz",+1) => 1
antiferm_vac = QN("Sz",0) => 1
even_link_p1 = QN("Sz",+2) => 1 # p1 = "positive 1"
even_link_0 = QN("Sz",0) => 1
even_link_n1 = QN("Sz",-2) => 1 # n1 = "negative 1"
odd_link_p1 = QN("Sz",-2) => 1
odd_link_0 = QN("Sz",0) => 1
odd_link_n1 = QN("Sz",+2) => 1

function Kogut_Susskind_sites(num_spatial_sites)
  sites = Index[]
  for i = 1:num_spatial_sites
    push!(sites, Index(ferm_full,ferm_vac; tags="Site,S=1/2,n=$(4*i-3)"))
    push!(sites, Index(odd_link_p1,odd_link_0,odd_link_n1; tags="Site,S=1,n=$(4*i-2)"))
    push!(sites, Index(antiferm_full,antiferm_vac; tags="Site,S=1/2,n=$(4*i-1)"))
    push!(sites, Index(even_link_p1,even_link_0,even_link_n1; tags="Site,S=1,n=$(4*i)"))
  end
  return sites
end

function Schwinger_Hamiltonian(sites, num_spatial_sites, x, μ)
  N = num_spatial_sites * 2
  ampo_H = AutoMPO()
  for j=0:N-1
    add!(ampo_H, 1, "Sz2", 2*j+2)
    add!(ampo_H, (μ/2)*(-1)^j * 2, "Sz", 2*j+1)
    if j != N-1
      add!(ampo_H, x * 2*2/sqrt(2), "S+", 2*j+1, "S-", 2*j+2, "S-", 2*j+3)
      add!(ampo_H, x * 2*2/sqrt(2), "S+", 2*j+3, "S+", 2*j+2, "S-", 2*j+1)
    end
  end
  add!(ampo_H, x * 2*2/sqrt(2), "S+", 2*N-1, "S-", 2*N, "S-", 1)
  add!(ampo_H, x * 2*2/sqrt(2), "S+", 1, "S+", 2*N, "S-", 2*N-1)
  H = MPO(ampo_H, sites)
  return H
end

function Schwinger_DMRG(num_spatial_sites)
  sites = Kogut_Susskind_sites(num_spatial_sites)
  H = Schwinger_Hamiltonian(sites, num_spatial_sites, x, μ)
  sweeps = Sweeps(10)
  maxdim!(sweeps, 10,20,100,100,200)
  cutoff!(sweeps, 1E-10)

  # construct initial state as completely empty
  init_state = ["Up" for n = 1:4*num_spatial_sites] # just to initialize
  for n = 1:4*num_spatial_sites
    if n % 4 == 1
      init_state[n] = "Dn" # empty fermion sites
    elseif n % 4 == 3
      init_state[n] = "Up" # empty antifermion states
    else # n % 4 == 2 || n % 4 == 0
      init_state[n] = "Up" # flux = 0 at all links
    end
  end
  psi0 = randomMPS(sites, init_state, 100)
  energy, psi = dmrg(H, psi0, sweeps)
  println("Ground state energy = $energy")
  return energy
end

In case you're wondering why I haven't already tested it, basically it's because I ran into some other more fundamental issue with the randomMPS function in Julia. This issue is mentioned on the Support Q&A (http://itensor.org/support/2186/questions-random-initial-conservation-randommps-linkdim) and on the GitHub page (https://github.com/ITensor/ITensors.jl/issues/424), so I think it will be fixed soon.

The idea is that, by assigning these particular quantum numbers, we have now made it so that each term in the Hamiltonian indeed has zero QN flux. So I think this means that, when the DMRG process varies the state, it will only end up doing things that respect the gauge symmetry? But maybe I'm misunderstanding.

The first problem I see with this is that I have still only defined one QN, so I'm pretty sure this QN conservation is not actually capturing the symmetry at each site. I believe there would be ways for DMRG to change the state in a way that respects this QN but doesn't respect the local constraints, in some way that I don't really understand. If this is the case, then I wonder whether setting up this type of quantum number for every site could be the way to go.

The second problem I see with this is that it is unclear to me how DMRG will work properly. Thanks to the gauge symmetries, the smallest "legal" operation is creating a fermion-antifermion pair on adjacent sites and changing the flux between the two particles by 1 (I'm assuming the fermion is negative, so +1 if the antifermion comes first and -1 if the fermion comes first), or the opposite (destroying a fermion-antifermion pair and making the corresponding flux change for the site between them). This means that one needs to change three adjacent sites (two particle sites and the link between them) to make a legal change to the state, and DMRG only isolates two sites at a time when performing the local optimization. This seems to me like something that would lead to trouble.

Anyway, I'm sorry to disappoint if you were expecting a definitive answer, but I thought this was extensive enough to put as an answer to at least get some thoughts going. I would love to hear people's thoughts on this.