ITensor Basics

This tutorial introduces Index objects in ITensors, how they are used to label tensor legs, and how ITensors contracts tensors by matching indices.

For the original ITensors tutorial on basic ITensor objects, see the ITensor examples page.

For the tensor-network background behind these objects, see Tensor Networks in Physics.

Setup

We load ITensors for basic tensor operations and LinearAlgebra for a dense singular-value check at the end.

using ITensors
import LinearAlgebra

roundreal(x; digits=6) = round(real(x); digits=digits)

roundreal (generic function with 1 method)

Two indices with the same dimension and tags can still be different indices. This is intentional. ITensors contracts tensors only when the indices are the same objects, not merely when they look similar.

This rule becomes extremely important later, especially in Liouville-space vectorization and process-tensor contractions.

Creating an ITensor

We now create a tensor

\[T_{ijk}\]

with dimensions 2 × 3 × 2.

We fill it with the numbers 1, 2, ..., 12.

A = reshape(1:12, 2, 3, 2)
T = ITensor(A, i, j, k)

println("T has indices:")
println(inds(T))

T has indices:
((dim=2|id=461|"i"), (dim=3|id=719|"j"), (dim=2|id=677|"k"))

Accessing a component looks like the mathematical notation:

\[T_{2,3,1}.\]

println("T[i=>2, j=>3, k=>1] = ", T[i => 2, j => 3, k => 1])

@assert T[i => 2, j => 3, k => 1] == 6

T[i=>2, j=>3, k=>1] = 6.0

Pro tip

The order of the index-value pairs does not matter when accessing an ITensor element. The labels identify the legs. This is in contrast with the conventional methods where we need to worry about the order of the indices.

You can also modify tensor entries.

println("Original T[1,1,1] = ", T[i => 1, j => 1, k => 1])

T[i => 1, j => 1, k => 1] = -100
println("Modified T[1,1,1] = ", T[i => 1, j => 1, k => 1])

@assert T[i => 1, j => 1, k => 1] == -100

Original T[1,1,1] = 1.0
Modified T[1,1,1] = -100.0

Convert the ITensor to an ordinary Julia array by specifying the index order.

T_array = Array(T, i, j, k)

println("Size of Array(T, i, j, k): ", size(T_array))
println("Array element [2,3,1]: ", T_array[2, 3, 1])
println("T[i => 2, j => 3, k => 1]: ", T[i => 2, j => 3, k => 1])

@assert T_array[2, 3, 1] == T[i => 2, j => 3, k => 1]

Size of Array(T, i, j, k): (2, 3, 2)
Array element [2,3,1]: 6.0
T[i => 2, j => 3, k => 1]: 6.0

Why specify the index order?

ITensors are not meant to be treated as arrays with a hidden axis order. When you convert an ITensor into an array, you explicitly say which index should become the first, second, third, ... array axis.

Contractions and the SVD

Contracting tensor indices

Tensor contraction means summing over shared indices.

For example, if

\[A_{ij}\]

and

\[B_{jk}\]

share the index j, then their contraction is

\[C_{ik} = \sum_j A_{ij} B_{jk}.\]

In ITensors, this contraction is simply multiplication.

a = Index(2, "a")
b = Index(3, "b")
c = Index(4, "c")

A = random_itensor(a, b)
B = random_itensor(b, c)

C = A * B

println("inds(A) = ", inds(A))
println("inds(B) = ", inds(B))
println("inds(C) = ", inds(C))

@assert a in collect(inds(C))
@assert c in collect(inds(C))
@assert !(b in collect(inds(C)))

inds(A) = ((dim=2|id=0|"a"), (dim=3|id=308|"b"))
inds(B) = ((dim=3|id=308|"b"), (dim=4|id=651|"c"))
inds(C) = ((dim=2|id=0|"a"), (dim=4|id=651|"c"))

The shared index b disappeared because it was summed over.

Tracing two legs with a delta tensor

We can also contract two legs of the same tensor using a Kronecker delta.

Suppose

\[T_{ijk}\]

has dim(i) = dim(k). Then

\[v_j = \sum_{a=1}^{\dim(i)} T_{aja}\]

is obtained by inserting a delta tensor

\[\delta_{ik}.\]

In code:

δik = delta(i, k)
v = δik * T

println("inds(delta(i,k) * T) = ", inds(v))
println("v as an array over j:")
println(Array(v, j))

inds(delta(i,k) * T) = ((dim=3|id=719|"j"),)
v as an array over j:
[-92.0, 13.0, 17.0]

Let us check the formula manually.

manual_v = zeros(Float64, dim(j))
for jj in 1:dim(j)
    manual_v[jj] = T[i => 1, j => jj, k => 1] + T[i => 2, j => jj, k => 2]
end

@assert Array(v, j) ≈ manual_v

Pro tip

A delta tensor is often the cleanest way to express traces, index identifications, and partial traces in ITensor code.

Tensor SVD

Singular value decomposition is the basic operation behind tensor-network compression.

For a matrix,

\[M = U S V^\dagger.\]

For a tensor, we first decide which indices belong to the left group and which belong to the right group.

Here we view

\[T_{ijk}\]

as a matrix

\[T_{(ij),k}.\]

Then we decompose it as

\[T_{ijk} = \sum_{\alpha} U_{ij,\alpha} S_{\alpha} V_{\alpha,k}.\]

U, S, V, spec = svd(T, (i, j))

T_reconstructed = U * S * V

println("Reconstruction error after SVD:")
println(norm(T - T_reconstructed))

@assert norm(T - T_reconstructed) < 1e-10

Reconstruction error after SVD:
2.862998051573224e-14

The singular values are the singular values of the reshaped matrix T[(i,j), k].

T_matrix = reshape(Array(T, i, j, k), dim(i) * dim(j), dim(k))
svals_dense = LinearAlgebra.svdvals(T_matrix)

println("Singular values of T[(i,j),k]:")
println(roundreal.(svals_dense))

Singular values of T[(i,j),k]:
[100.573979, 23.10573]

Why SVD matters

Tensor SVD is the local move behind MPS construction, MPS compression, MPO compression, TEBD, TDVP, and the memory compression ideas that appear later in process-tensor algorithms.

Next: MPS and MPO Basics applies these ideas to many-body states and operators on spin chains.

ITensor Basics

Setup

ITensor indices and objects

Similar-looking indices are not necessarily the same

Creating an ITensor

Contractions and the SVD

Contracting tensor indices

Tracing two legs with a delta tensor

Tensor SVD