two sites operator and apply function with MPO

ITensor Julia Questions
, ,
1

Hello everyone,
I’m a beginner with ITensor and I have a question about the usage of apply function in two-site operator case.
My code is following.

using ITensorMPS
using ITensors

s = siteinds("Qubit",2) 
qubit_state_origin_1 = MPS(s,"0")
rho_id = outer(qubit_state_origin_1',qubit_state_origin_1) #create density matrix
rho_id[1][4] = 1 # create identity at each site
rho_id[2][4] = 1 
rho_mult = apply(rho_id,rho_id; alg="zipup") #contract by apply function
rho_ele_mult = rho_id' .* rho_id #contract by boardcasting
sz_test = op("Sz",s[1])*op("Id",s[2])
state_evolution = apply(sz_test,rho_id)

First, I create a two-site MPO called rho_id with each site is 2 by 2 identity matrix. And what I want to do is applying the rho_id to another rho_id to find whether each site is still identity matrix after contraction. The result expected to be [1 0; 0 1] for each site. However, the result is

@show rho_mult[1]
rho_mult[1] = ITensor ord=3
Dim 1: (dim=2|id=766|"Qubit,Site,n=1")'
Dim 2: (dim=2|id=766|"Qubit,Site,n=1")
Dim 3: (dim=1|id=958|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}
 2×2×1
[:, :, 1] =
 -1.4142135623730947   0.0
  0.0                 -1.414213562373095
ITensor ord=3 (dim=2|id=766|"Qubit,Site,n=1")' (dim=2|id=766|"Qubit,Site,n=1") (dim=1|id=958|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}

I have no idea why there is a additional minus sign before the strange number.
And if I use boardcast “.*”, the result is reasonable with the cost of adding a dimension.

 @show rho_ele_mult[2]
rho_ele_mult[2] = ITensor ord=4
Dim 1: (dim=2|id=722|"Qubit,Site,n=2")''
Dim 2: (dim=1|id=865|"Link,l=1")'
Dim 3: (dim=2|id=722|"Qubit,Site,n=2")
Dim 4: (dim=1|id=865|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}
 2×1×2×1
[:, :, 1, 1] =
 1.0
 0.0

[:, :, 2, 1] =
 0.0
 1.0
ITensor ord=4 (dim=2|id=722|"Qubit,Site,n=2")'' (dim=1|id=865|"Link,l=1")' (dim=2|id=722|"Qubit,Site,n=2") (dim=1|id=865|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}

The reason I play with this kind of thing is that I’d like to create a global operator like sz_test defined in the code. And applying it to MPO state to evolute the qubit state. The result is also strange.

@show state_evolution[1]
state_evolution[1] = ITensor ord=3
Dim 1: (dim=2|id=766|"Qubit,Site,n=1")'
Dim 2: (dim=2|id=766|"Qubit,Site,n=1")
Dim 3: (dim=4|id=8|"Link,n=1")
NDTensors.Dense{Float64, Vector{Float64}}
 2×2×4
[:, :, 1] =
 -0.7071067811865472  0.0
  0.0                 0.7071067811865475

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

[:, :, 3] =
 0.0  1.0
 0.0  0.0

[:, :, 4] =
 0.7071067811865475  0.0
 0.0                 0.7071067811865475
ITensor ord=3 (dim=2|id=766|"Qubit,Site,n=1")' (dim=2|id=766|"Qubit,Site,n=1") (dim=4|id=8|"Link,n=1")
NDTensors.Dense{Float64, Vector{Float64}}

I’d sincerely appreciate it if someone could help me to solve this problem or share the method of applying global operator to MPO state.

2

Maybe the error comes from MPO structure, since what I was doing in my code is applying a MPO with link dimension to itself, rather than applying MPO without link dimension, though I still don’t get the difference.
For the implementation, I was thinking there are two kinds of methods to implement this global operator. Like the following graphs showing.

Methods
Methods1310×1079 56.6 KB

The first method’s code is something like this

sz_test = op("Sz",s[1])*op("Id",s[2])
state_evolution = apply(sz_test,rho_id)

And the code of second method is

gates = ITensor[]
for j in 1:2
    s1 = s[j]
    hj =
        op("Sz", s1)
    push!(gates, hj)
end
state = apply(gates,rho_id)

I don’t know whether the methods I’m think is correct or not, since the final output are strange.
Output of 1st method

@show state_evolution
state_evolution = MPO
[1] ((dim=2|id=988|"Qubit,Site,n=1")', (dim=2|id=988|"Qubit,Site,n=1"), (dim=4|id=166|"Link,n=1"))
[2] ((dim=4|id=166|"Link,n=1"), (dim=2|id=599|"Qubit,Site,n=2")', (dim=2|id=599|"Qubit,Site,n=2"))

@show state_evolution[1]
state_evolution[1] = ITensor ord=3
Dim 1: (dim=2|id=988|"Qubit,Site,n=1")'
Dim 2: (dim=2|id=988|"Qubit,Site,n=1")
Dim 3: (dim=4|id=166|"Link,n=1")
NDTensors.Dense{Float64, Vector{Float64}}
 2×2×4
[:, :, 1] =
 -0.7071067811865472  0.0
  0.0                 0.7071067811865475

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

[:, :, 3] =
 0.0  1.0
 0.0  0.0

[:, :, 4] =
 0.7071067811865475  0.0
 0.0                 0.7071067811865475
ITensor ord=3 (dim=2|id=988|"Qubit,Site,n=1")' (dim=2|id=988|"Qubit,Site,n=1") (dim=4|id=166|"Link,n=1")
NDTensors.Dense{Float64, Vector{Float64}}

@show state_evolution[2]
state_evolution[2] = ITensor ord=3
Dim 1: (dim=4|id=166|"Link,n=1")
Dim 2: (dim=2|id=599|"Qubit,Site,n=2")'
Dim 3: (dim=2|id=599|"Qubit,Site,n=2")
NDTensors.Dense{Float64, Vector{Float64}}
 4×2×2
[:, :, 1] =
 -0.7071067811865476  0.0
  0.0                 0.0
  0.0                 0.0
  0.0                 0.0

[:, :, 2] =
 0.0  -0.7071067811865475
 0.0   0.0
 0.0   0.0
 0.0   0.0
ITensor ord=3 (dim=4|id=166|"Link,n=1") (dim=2|id=599|"Qubit,Site,n=2")' (dim=2|id=599|"Qubit,Site,n=2")
NDTensors.Dense{Float64, Vector{Float64}}

Output of 2nd method

@show state
state = MPO
[1] ((dim=2|id=988|"Qubit,Site,n=1")', (dim=2|id=988|"Qubit,Site,n=1"), (dim=1|id=384|"Link,l=1"))
[2] ((dim=2|id=599|"Qubit,Site,n=2")', (dim=2|id=599|"Qubit,Site,n=2"), (dim=1|id=384|"Link,l=1"))

@show state[1]
state[1] = ITensor ord=3
Dim 1: (dim=2|id=988|"Qubit,Site,n=1")'
Dim 2: (dim=2|id=988|"Qubit,Site,n=1")
Dim 3: (dim=1|id=384|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}
 2×2×1
[:, :, 1] =
 -0.7071067811865472  0.0
  0.0                 0.7071067811865475
ITensor ord=3 (dim=2|id=988|"Qubit,Site,n=1")' (dim=2|id=988|"Qubit,Site,n=1") (dim=1|id=384|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}

@show state[2]
state[2] = ITensor ord=3
Dim 1: (dim=2|id=599|"Qubit,Site,n=2")'
Dim 2: (dim=2|id=599|"Qubit,Site,n=2")
Dim 3: (dim=1|id=384|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}
 2×2×1
[:, :, 1] =
 -0.3535533905932738  -0.0
  0.0                  0.3535533905932738
ITensor ord=3 (dim=2|id=599|"Qubit,Site,n=2")' (dim=2|id=599|"Qubit,Site,n=2") (dim=1|id=384|"Link,l=1")
NDTensors.Dense{Float64, Vector{Float64}}
3

If I understand what you’re doing, I think you’re seeing a gauge degree of freedom from the MPO. You unfortunately can’t always inspect each tensor to know if its equal due to this. However we can try something else:

julia> norm(rho_mult - rho_id)
4.965068306494546e-16

julia> rho_mult ≈ rho_id
true

And you can see that I^2 = I.

Then for the second part, the only thing I can say is that usually one evolves a state by applying the unitary on both sides. What you’ve done is take |\psi\rangle\langle\psi| \to U|\psi\rangle\langle\psi|. I would also try to focus on observables, or at least things not affected by gauge degrees of freedom

4

hi ryan,

Thanks a lot for your reply. It’s really helpful and your suggestion is very useful, directly to the point.
Additionally, is there any difference between using a MPO with 1-dim link and a MPO without linkdimension as a evolution in TEBD method? (maybe the front MPO is just svd form of later)

two_kind_MPO
two_kind_MPO1839×823 54.2 KB

And futhermore, are the two methods really different in practice (or for the perspective of ITensor)?

5

The with links is potentially a compression of the second (without links). They can always potentially equal each other but the second (without links) is exponentially large always, and the first (with links) can be smaller. This is the idea of “rank”.

When contracting all the indices (prod) you can turn any ‘with links’ into a single tensor without links, however due to the exponential scaling one has to be careful about the size, as most are intractable.

6

Thank you very much!

7

This topic was automatically closed 10 days after the last reply. New replies are no longer allowed.