MultiIndices.jl

Introduction

The idea is provide means for creating tensor product spaces that can be efficiently indexed into; a 3d problem, maybe described in Cartesian coordinates or spherical coordinates, etc,, all of which maybe discretized or represented using basis functions, etc. Depending on the choice of coordinate system, but also on the storage pattern of the different degrees of freedom, the sparsity structure of various operators will be different. MultiIndex.jl’s purpose is to provide an efficient means of setting up these tensor product spaces and to use “logical” indices, which are mapped to linear indices in the underlying vector in the product space.

Usage pattern

using MultiIndices

Cartesian grids

A simple Cartesian 1d grid can be setup like this:

𝕽 = Cartesian(1000)

MultiIndices.Cartesian{1,(1,)}((1000,))

while a Cartesian 2d grid can be setup like this:

𝖃𝖄 = CartesianXY(100,500)

MultiIndices.Cartesian{2,(1, 2)}((100, 500))

CartesianXY indicates that X is the “fast dimension”, i.e. the dimension whose elements are found next to each other in vector, while subsequent elements in Y are strided.

We see then that the second element along the X dimension has linear index 2:

𝖃𝖄[2,1]

2

while the second element along the Y dimension has the linear index 101:

𝖃𝖄[1,2]

101

Tensor product spaces

As a more complicated example, the argon atom consists of 18 electron, each of which has a position along the radial direction r, an angular momentum ℓ, a magnetic quantum number m, and a spin. Its tensor product space can be created the following way:

1. We first define the orbital angular momentum basis, which can be built up from a direct sum of subspaces which are invariant under rotation. The first subspace, with ℓ=0 (called s in spectroscopic notation) has one element, m=0; the second subspace, ℓ=1 (p) has 2ℓ+1=3 elements, m=-1,0,1; and so on. If we want to include all ℓ∈[0,10], we can form the space thus:

   𝕷 = RotationSpace(0)
   for ℓ = 1:10
       𝕷 = 𝕷 ⊕ RotationSpace(ℓ)
   end
   𝕷
       

   MultiIndices.SumSpace{MultiIndices.RotationSpace}(MultiIndices.RotationSpace[Rot{ℓ = s}, Rot{ℓ = p}, Rot{ℓ = d}, Rot{ℓ = f}, Rot{ℓ = g}, Rot{ℓ = h}, Rot{ℓ = i}, Rot{ℓ = k}, Rot{ℓ = l}, Rot{ℓ = m}, Rot{ℓ = n}])
       

   We can query the size of this space:

   size(𝕷)
       

   121
       

   i.e., there are 121 partial waves. We can also find the linear index of a partial wave, given its values of ℓ and m (no bounds checking, so make sure you obey m ≤ |ℓ|):

   𝕷[ℓm(1,1)...]
       

   3 (wrong at the moment)
       

   Spectroscopic notation is also supported:

   𝕷[ℓm('d',1)...]
       

   4 (wrong at the moment)
       

2. We then define the electron spin momentum space and the combined angular momentum space:

   𝔖 = Spin{1//2}()
   𝕵 = 𝕷 ⊗ 𝔖
       

   MultiIndices.ProductSpace{2,(1, 2)}(MultiIndices.SumSpace{MultiIndices.RotationSpace}(MultiIndices.RotationSpace[Rot{ℓ = s}, Rot{ℓ = p}, Rot{ℓ = d}, Rot{ℓ = f}, Rot{ℓ = g}, Rot{ℓ = h}, Rot{ℓ = i}, Rot{ℓ = k}, Rot{ℓ = l}, Rot{ℓ = m}, Rot{ℓ = n}]), MultiIndices.Spin{1//2}())
       

3. Next, we define a 1d coordinate system for the radial dimension and the total product space for one electron

   𝕽 = Cartesian(100)
   𝖂 = 𝕽 ⊗ 𝕵
       

   MultiIndices.ProductSpace{3,(1, 2, 3)}(MultiIndices.Cartesian{1,(1,)}((100,)), MultiIndices.ProductSpace{2,(1, 2)}(MultiIndices.SumSpace{MultiIndices.RotationSpace}(MultiIndices.RotationSpace[Rot{ℓ = s}, Rot{ℓ = p}, Rot{ℓ = d}, Rot{ℓ = f}, Rot{ℓ = g}, Rot{ℓ = h}, Rot{ℓ = i}, Rot{ℓ = k}, Rot{ℓ = l}, Rot{ℓ = m}, Rot{ℓ = n}]), MultiIndices.Spin{1//2}()))
       

   the size of which is

   size(𝖂)
       

   (100, 121, 2)
       

   as expected.

4. With the product space for one electron ready, we can easily replicate it 18 times:

   𝕱 = 𝖂 ⊗ 18
       

   MultiIndices.CopySpace{4,(1, 2, 3, 4)}(MultiIndices.ProductSpace{3,(1, 2, 3)}(MultiIndices.Cartesian{1,(1,)}((100,)), MultiIndices.ProductSpace{2,(1, 2)}(MultiIndices.SumSpace{MultiIndices.RotationSpace}(MultiIndices.RotationSpace[Rot{ℓ = s}, Rot{ℓ = p}, Rot{ℓ = d}, Rot{ℓ = f}, Rot{ℓ = g}, Rot{ℓ = h}, Rot{ℓ = i}, Rot{ℓ = k}, Rot{ℓ = l}, Rot{ℓ = m}, Rot{ℓ = n}]), MultiIndices.Spin{1//2}())), 18)
       

   size(𝕱)
       

   (100, 121, 2, 18)
       

   We now try to find the linear index of the second element along each dimension:

   [𝕱[2,1,1,1] 𝕱[1,2,1,1] 𝕱[1,1,2,1] 𝕱[1,1,1,2]]
       

   2

   101

   12101

   24201

   which is correct.

Wish list

• Interval indexing, e.g. 𝕱[:,1,2,3]
• Efficient iteration, via Base.Cartesian