MPO construction from two site custom op

Hi,
I want to build an MPO using a custom defined ITensor operator. The MPO can be constructed using
MPO(::sites, ::list of op) (reference to definition : https://itensor.github.io/ITensors.jl/dev/MPSandMPO.html).
However, the caveat is that the list of operators are two-body. For example, for spin half systems, this two body operator could be Z \otimes Z or the one shown below.

using ITensors
using LinearAlgebra

sites = siteinds("Qubit", 6)

function ITensors.op(::OpName"TwoSiteOp", ::SiteType"Qubit"; θ1, θ2)
    
    M = cos(θ1)*Diagonal(ones(4)) - im*sin(θ2)*[1 0 0 0; 0 -1 0 0; 0 0 -1 0; 0 0 0 1]

    return M
end

Now, the question is how can one construct an MPO for the defined ‘siteinds’, sites?
The op list will be made of various operators with different \theta 1 and \theta 2.

Any suggestion for this will be very helpful.

See Miles’ answer from here:

Essentially the OpSum system can’t do that at the moment, however if you can write it as a string (O_i \otimes O_j \dots ) or a sum of strings as in the post then that can be used. And for generic 2 site matrices you would have to decompose it by hand into that format

Hi Ryan,
Thank you very much for the reply. To clarify, the two-site operator I mentioned above, e^{-i \theta_1 \mathbb Z_{1} \otimes \mathbb Z_{2}}, can be decomposed as a product of single-site tensors as
\begin{pmatrix} \cos(\theta_1) \mathbb{I}_{1} &-i\sin(\theta_1) \mathbb Z_{1} \end{pmatrix} \otimes \begin{pmatrix} \mathbb{I}_{2} \\ \mathbb Z_{2} \end{pmatrix}.
Is there a way to construct an MPO by using such tensors as input?

If you have the operator O = \cos(\theta_1) I_1 \otimes I_2 -i\sin(\theta_1) Z_1 \otimes Z_2 then the corresponding OpSum call would be

os = OpSum()
os += cos(θ1), "I",1,"I",2      # the second I is not needed really
os += -1im*sin(θ2),"Z",1,"Z",2
O = MPO(os,sites)

Hi Ryan,
thank you so much! There might still be a miscommunication though. Just to be sure: I need to build an MPO for the operator
e^{-i\theta\mathbb Z_{1} \otimes \mathbb Z_{2}} \otimes e^{-i\theta\mathbb Z_{3} \otimes \mathbb Z_{4}} \otimes ....... \otimes e^{-i\theta\mathbb Z_{N-1} \otimes \mathbb Z_{N}} .
This is a product rather than sum of two-site operators (unlike the link with Miles’ answer).
Taking your suggestion to use OpSum, I tried the following code :

using ITensors

N = 8
sites = siteinds("Qubit", N)
U = MPO(sites)

θ = 0.5
os = OpSum()
os += cos(θ), "I", 1, "I", 2
os += -1im * sin(θ), "Z", 1, "Z", 2

for j in 1:2:N-1
    U[j:j+1] = MPO(os, sites[j:j+1])
end

However, this does not seem to produce the correct results. Could you help me correct it?

You’ll need to decompose your operator into a list of single site operators
O = \sum_n O_1^n \otimes O_2^n \otimes \dots, which for each n you would write the os+= line as above. You’d then call MPO just once to form the MPO.

Then one trick you can do, if you want say all "Z"s like A*Z_1 \otimes Z_2 \otimes \dots \otimes Z_N where A is a constant, then you can write

A = 0.1 # etc
term = [("Z",i) for i=1:N]
add!(os, A, (term...)...)

Thank you very much.

1 Like