MPO matrix representations

ITensor Julia Questions
, ,
1

Hi,
is there a way to print the matrices of an MPO/MPS? What I mean by that is, an MPO is a formed by multiplying matrices on each site. However, ITensors, stores these matrices as 4 indexed objects in case of MPOs and 3 indexed objects in case of MPS. Now, what I want is to print these 4 indexed and 3 indexed objects as matrices. In other words, I want to print the constituent matrices of MPO/MPS.

2

Hello,

You can easily print the constituent matrices using Julia array access notation. For example

julia> siteinds("S=1/2", 10)
julia> psi = randomMPS(s)
julia> @show psi[5]
psi[5] = ITensor ord=3
    Dim 1: (dim=1|id=416|"Link,l=4")
    Dim 2: (dim=2|id=192|"S=1/2,Site,n=5")
    Dim 3: (dim=1|id=967|"Link,l=5")
    Dense{Float64, Vector{Float64}}
     1Ă—2Ă—1
    [:, :, 1] =
     -0.8970638949108982  0.4419008581653681
    ITensor ord=3 (dim=1|id=416|"Link,l=4") (dim=2|id=192|"S=1/2,Site,n=5") (dim=1|id=967|"Link,l=5")
    Dense{Float64, Vector{Float64}}
3

Hi Bhargava,
Karl’s answer is right, but you may want to see the output organized so that the printing always shows the “matrices” in the sense of how a “matrix product state” is a collection of d matrices on each site, one for every setting of the “physical index” of dimension d. And similar for an MPO thought of as a collection of d^2 indices on each site.

So I made a small function that permutes the index ordering of each tensor to put the physical or Site indices last. This way the printout shows each matrix grouped together. Please give it a try:

function print_matrices(A::AbstractMPS)
  s = siteinds(A)
  for j=1:length(s)
    Aj = A[j]
    links = uniqueinds(Aj,s[j])
    Aj = permute(Aj,links...,collect(s[j])...)
    println(Aj)
  end
end

For example, when I use it to print a Hamiltonian MPO (XY model for S=1/2) I get:

Dim 1: (dim=4|id=324|"Link,l=1")
Dim 2: (dim=2|id=640|"Qubit,Site,n=1")'
Dim 3: (dim=2|id=640|"Qubit,Site,n=1")
NDTensors.Dense{Float64, Vector{Float64}}
 4Ă—2Ă—2
[:, :, 1] =
 0.0  0.0
 0.0  0.0
 0.0  0.5
 1.0  0.0

[:, :, 2] =
 0.0  0.0
 0.5  0.0
 0.0  0.0
 0.0  1.0
ITensor ord=4
Dim 1: (dim=4|id=324|"Link,l=1")
Dim 2: (dim=4|id=576|"Link,l=2")
Dim 3: (dim=2|id=762|"Qubit,Site,n=2")'
Dim 4: (dim=2|id=762|"Qubit,Site,n=2")
NDTensors.Dense{Float64, Vector{Float64}}
 4Ă—4Ă—2Ă—2
[:, :, 1, 1] =
 1.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  1.0

[:, :, 2, 1] =
 0.0  0.0  0.0  0.0
 1.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.5  0.0

[:, :, 1, 2] =
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 1.0  0.0  0.0  0.0
 0.0  0.5  0.0  0.0

[:, :, 2, 2] =
 1.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0
 0.0  0.0  0.0  1.0

...
2 Likes
4

Hi Karl and Miles,
Thank you very much for the detailed reply. Your answers were very useful. :slight_smile: