set ITensor elements Array-wise

Question

HI!
I was wondering if it is possible to set ITensor elements Array-wise instead of element-wise, something like:

using ITensors
i = Index(2)
j = Index(2)
k = Index(2)
l = Index(2)

A = randomITensor(i,j,k)
A[i=>1, j=>2, :, :] = [0 0; 0 0]

The use case is that I would like to "diagonalize" a Tensor with (2+n) indices, in the sense that all entries of the Tensor should become zero when the values of the first two indices (which have the same dimensions) are not equal.

I am using Julia v1.5.2 and ITensors v0.1.24.

Thanks for any reply!

Answer 1

Hi,

Good question. Unfortunately we don't have very good support right now for slicing ITensors like that. However, in this case you can make use of the NDTensors.jl library, which is a "basic tensor" library that is used by ITensors for basic tensor operations. Here is an example for how to use it:

using ITensors
using ITensors.NDTensors

i = Index(2)
j = Index(2)
k = Index(2)
l = Index(2)

A = randomITensor(i, j, k, l)
tensor(A)[1, 2, :, :] = [0 0; 0 0]

Note that when you index into a Tensor object (which is what you get when you call tensor(A)), it will use the current Index ordering of A, so you will have to be careful about the ordering of the indices. You can print the current Index ordering with @show inds(A).

The ITensor functionality you suggest would be pretty easy to implement, I'm working on it now and I'll add a comment here when it is implemented.

-Matt

Comments

Hi Matt

Thanks for the answer, very insightful. I didn't know about this way of indexing, but it does exactly what I was looking for. Is there a way to display the result of such indexing (@show tensor(A)[1, 2, :, :] does throw a MethodError)?
Is there a general rule regarding index ordering of an ITensor (by id e.g.)?

And thanks for right away working on this feature - that's incredible!

Christian

---

Hi Christian,

If I understand your first question, I think what you are asking is what:

    tensor(A)[1, 2, :, :]

actually returns. In fact, this operation isn't defined right now (just as an oversight on our side), but this similar slicing operation is defined:

    julia> tensor(A)[1, 2, 1:end, 1:end]
    Dim 1: 1
    Dim 2: 1
    Dim 3: 2
    Dim 4: 2
    Dense{Float64,Base.ReshapedArray{Float64,1,SubArray{Float64,4,Array{Float64,4},NTuple{4,UnitRange{Int64}},false},Tuple{Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64},Base.MultiplicativeInverses.SignedMultiplicativeInverse{Int64}}}}
     1×1×2×2
    [:, :, 1, 1] =
     0.0
    
    [:, :, 2, 1] =
     0.0
    
    [:, :, 1, 2] =
     0.0
    
    [:, :, 2, 2] =
     0.0

This looks complicated, but basically it is a view into the specified slice of the Tensor. However, for the operation:

    tensor(A)[1, 2, :, :] = [0 0; 0 0]

it doesn't actually require creating the tensor `tensor(A)[1, 2, :, :]`, since Julia directly turns it into the function call `setindex!(tensor(A), [0 0; 0 0], 1, 2, :, :)` (this is a slightly technical detail about how Julia parses those types of expressions).

As for the Index ordering of an ITensor, at first it is ordered based on how it was defined. Afterward, in general it is ordered in such a way that data permutations are minimized (for example after a contraction). Note that you can permute the ITensor to change the ordering of the ITensor, using the function:

    permute(A, l, i, k, j)

Note that this performs a permutation of the ITensor data as well, so the data matches the new index ordering.

-Matt

---

I've added support for slicing operations here: https://github.com/ITensor/ITensors.jl/pull/535

It will be included in the next version of ITensors.jl (either today or tomorrow).

---

Awesome! Thanks a lot!