Extracting QN information from an index

ITensor Julia Questions
1

Suppose I have an index with quantum numbers, e.g.

i
(dim=9|id=689|"Link,u") <Out>
 1: QN("Sz",0) => 4
 2: QN("Sz",-2) => 4
 3: QN("Sz",-4) => 1

How do I extract information about which index values are in which QN blocks? For example, I can see i=>6 is in the QN(“Sz”,-2) block, but I can’t find a function that tells me this, say qn(i,6) == QN(“Sz”,-2). Alternatively, say I wanted the range of index values that correspond to QN(“Sz”,-2), here 5:9. Is there something like qnrange(i,QN(“Sz”,-2)) == 5:9 ?

A specific example is examining the Schmidt decomposition for a link: one might like to print out or work with the Schmidt values organized by QN. Also, one might like to extract some of the Schmidt vectors by QN.

2

Hi Steve,
It’s a good question. We should have more of an intro in the documentation to the QN system in terms of the Julia implementation. Right now, though, Matt is reworking a lot of the internals and specifically a range-based block sparse system which could trickle up into the internals of QN Index eventually too.

But the good news is that most of the functions you asked about are there, except for your qnrange suggestion. For that one, I’ve provided an implementation below.

See below for a complete example showing qn, blockdim (for obtaining dimensions of each block), and qnrange.

using ITensors
using Test

function qnrange(i::Index, q::QN)
  j = findfirst(b->(b.first==q),space(i))
  isnothing(j) && error("Index $i\n does not have block with QN = $q")
  start = sum(k->blockdim(i,k),1:(j-1); init=1)
  return start:(start+blockdim(i,j)-1)
end

@testset "QN Index" begin

  # Index with three "QN blocks"
  i = Index(QN("Sz",0) => 4, QN("Sz",-2) => 4, QN("Sz",-4) => 1)

  @test nblocks(i) == 3
  @test dim(i) == 9

  # Get individual QNs of each QN block
  @test qn(i,1) == QN("Sz",0)
  @test qn(i,2) == QN("Sz",-2)
  @test qn(i,3) == QN("Sz",-4)

  # Get dimensions of each QN block
  @test blockdim(i,1) == 4
  @test blockdim(i,2) == 4
  @test blockdim(i,3) == 1

  # Get the underlying array of block info
  @test space(i) == [QN("Sz",0) => 4, QN("Sz",-2) => 4, QN("Sz",-4) => 1]

  return
end
1 Like
3

To get Schmidt values after running a DMRG calculation, here’s an example I just put together:

using ITensors

function qnrange(i::Index, q::QN)
  j = findfirst(b->(b.first==q),space(i))
  isnothing(j) && error("Index $i\n does not have block with QN = $q")
  start = sum(k->blockdim(i,k),1:(j-1); init=1)
  return start:(start+blockdim(i,j)-1)
end

let
  N = 10
  s = siteinds("S=1", N; conserve_qns=true)

  hterms = OpSum()
  for j=1:N-1
    hterms += 1/2,"S+",j,"S-",j+1
    hterms += 1/2,"S-",j,"S+",j+1
    hterms += "Sz",j,"Sz",j+1
  end
  H = MPO(hterms,s)

  nsweeps = 8
  cutoff = 1E-6
  maxdim = [20,40,80,160,200]

  psi0 = randomMPS(s, [isodd(j) ? "Up" : "Dn" for j=1:N]; linkdims=4)

  energy, psi = dmrg(H,psi0; nsweeps, cutoff, maxdim)

  c = N÷2
  orthogonalize!(psi,c)
  l = commonindex(psi[c],psi[c-1])

  U,S,V = svd(psi[c],(l,s[c]))
  u = commonindex(S,U)

  q = QN("Sz",0) # could also get this from getting qn(u,1), qn(u,2), ...
  println("Schmidt values for q=$q")
  for j=qnrange(u,q)
    @show S[j,j]
  end

  q = QN("Sz",2) # could also get this from getting qn(u,1), qn(u,2), ...
  println("Schmidt values for q=$q")
  for j=qnrange(u,q)
    @show S[j,j]
  end
2 Likes
4

Hi Miles,
Excellent, thanks for the very quick response.
The other type of qn(i,j) would find the quantum number associated with a particular index value j, with 1 <= j <= dim(i). Here is a version of that, called qnind:

function qnind(u,i)
    @assert 0 < i <= dim(u)
    for iq=1:nblocks(u)
        bd = blockdim(u,iq)
        i <= bd && return qn(u,iq)
        i -= bd
    end
    error("Shouldnt get here")
end

Here are two ways to print out all the Schmidt values, the first being what you referred to above:

for iq=1:nblocks(u)
    q = qn(u,iq)
    println("Schmidt values for q=$q")
    for j=qnrange(u,q)
        printsp(j,S[j,j])
    end
end

for j=1:dim(u)
    printsp(j,qnind(u,j),S[j,j])
end

Here is sample output for these to versions:

Schmidt values for q=QN("Sz",2)
1 2.6469291817725617e-5
Schmidt values for q=QN("Sz",0)
2 0.009569622735109015
3 0.0013009691142406793
Schmidt values for q=QN("Sz",-2)
4 0.9999533634622407
1 QN("Sz",2) 2.6469291817725617e-5
2 QN("Sz",0) 0.009569622735109015
3 QN("Sz",0) 0.0013009691142406793
4 QN("Sz",-2) 0.9999533634622407
2 Likes
5

I forgot to mention that my code snippet uses my utility function printsp to automatically insert spaces between the variables, and flush stdout. It is defined with the lines:

printlnf(a...) = (println(a...); flush(stdout))
printsp(a) = printlnf(a)
printsp(a,b...) = (print(a," "); printsp(b...))
6

I like printsp, I’ll have to use that!

I see what you were asking now in your first post about i=>6 being one of the Index values versus a block number. We don’t have that function provided unfortunately, so it’s helpful that you posted an implementation. Matt and I will have to come back to this thread as he reworks the QN and blocksparse system. (Also Olivier Gauthe is working on QNs specifically right now in the context of non-Abelian symmetries.) Though I think much of that work is pretty orthogonal to adding functions like these.