Data assignments between arrays of arrays

New to Julia
, ,
1

As an old Matlab user I am still sometimes struggling to understand how Julia assigns data in “arrays of arrays”:
For example, I want to copy a page of data from a 3D data array. The result is a Matrix, as expected:

A = zeros(3,3,2);
b = a[:,:,1]
julia> b
3×3 Matrix{Float64}:
** 0.0 0.0 0.0**
** 0.0 0.0 0.0**
** 0.0 0.0 0.0**

Similarly, I can replace a page of A:

A = zeros(3,3,2);
b = rand(3,3)
A[:,:,end] .= b
julia> A
3×3×2 Array{Float64, 3}:
[:, :, 1] =
** 0.0 0.0 0.0**
** 0.0 0.0 0.0**
** 0.0 0.0 0.0**

[:, :, 2] =
** 0.169528 0.718634 0.5193**
** 0.452352 0.807573 0.204051**
** 0.226577 0.223511 0.633584**

However, it doesn’t seem to work for “Arrays in Arrays” as I have naively thought. There, the extracted part b is still identical(!) to A

A = fill(rand(2,2,2),2,2)
b = A[:,:][:,:,1]

julia> A == b
true

So how can I copy a subset from A?
Analogously, it doesn’t seem to be possible to replace a subset of A by other data, like:

A = fill(rand(2,2,2),2,2)
julia> A[:,:][:,:,1] .= rand(2,2)
ERROR: MethodError: Cannot convert an object of type Float64 to an object of type Array{Float64, 3}

I know that this is a beginners question, but it really bothers me for some time now and all my web searches were unsucessful. I would appreciate any hints. Thanks!
Alex

2

In this case you query for the full array using A[:,:] and then ask for A[:,:,1], which works, as the index allows higher dimensions when using 1 as index.

julia> A[:,:] == A
true

julia> A[:,:,1] == A
true

julia> A[:,:,1,1] == A
true

What you want is the [:,:,1] for every index of A. This you could do using getindex, but this returns a Matrix per index of A.

julia> getindex.(A,:,:,1)
2×2 Matrix{Matrix{Float64}}:
 [0.310872 0.881689; 0.779912 0.138566]  [0.310872 0.881689; 0.779912 0.138566]
 [0.310872 0.881689; 0.779912 0.138566]  [0.310872 0.881689; 0.779912 0.138566]

A better way might be to use a 5 dimensional data structure.

julia> A = zeros(2,2,2,2,2);

julia> subA = rand(2,2,2)
2×2×2 Array{Float64, 3}:
[:, :, 1] =
 0.975893  0.584884
 0.463343  0.217259

[:, :, 2] =
 0.948323  0.414374
 0.381338  0.156916

julia> for i in 1:size(A,1), j in 1:size(A,2)
           A[i,j,:,:,:] .= subA
       end

julia> A[:,:,:,:,1]
2×2×2×2 Array{Float64, 4}:
[:, :, 1, 1] =
 0.975893  0.975893
 0.975893  0.975893

[:, :, 2, 1] =
 0.463343  0.463343
 0.463343  0.463343

[:, :, 1, 2] =
 0.584884  0.584884
 0.584884  0.584884

[:, :, 2, 2] =
 0.217259  0.217259
 0.217259  0.217259

To make your posts on discourse easier to read you can surround your code blocks by ``` .

4 Likes
3

There is also the package TensorCast.jl that might help with more complicated transforms.

Others to look at are Tullio.jl, TensorOperations.jl, OMEinsum.jl, Einsum.jl

4

As @feanor12 said, A[:, :] doesn’t do what you think, it just creates a copy of A. And since it is a copy of A, you would normally not see changes to A[:, :] in A, because you were changing the copy, not the original. Interestingly, in this case you do see the changes, because A and copy(A) are both pointing to the same array internally.

And that brings me to a different point. This:

A = fill(rand(2,2,2),2,2)

probably does not do what you expect. It creates a single array (rand(2,2,2)) and then creates a 2x2 array filled with references to the very same array. So this happens:

1.7.2> A = fill(rand(2,2,2),2,2);

1.7.2> A[1,1][1,1,1] = 100; # only make one change

1.7.2> A[1,1][:,:,1]  # expect one element to be 100.0
2×2 Matrix{Float64}:
 100.0       0.336786
   0.142365  0.36947

1.7.2> A[2,1][:,:,1]  # did not expect any element = 100.0, but there it is!
2×2 Matrix{Float64}:
 100.0       0.336786
   0.142365  0.36947

This is because, as explained, A is just a bunch of pointers all pointing to the same array.

Perhaps this is what you wanted:

A = [rand(2,2,2) for i in 1:2, j in 1:2]  # creates 4 different random arrays.

Now, if you want to set the first slice of each array inside A to the same values, you can try this:

1.7.2> setindex!.(A, Ref(fill(100.0, 2,2)), :, :, 1)
2×2 Matrix{Array{Float64, 3}}:
 [100.0 100.0; 100.0 100.0;;; 0.305401 0.883966; 0.541377 0.487306]  …  [100.0 100.0; 100.0 100.0;;; 0.218273 0.0704365; 0.527901 0.199329]
 [100.0 100.0; 100.0 100.0;;; 0.267275 0.0161309; 0.800307 0.14134]     [100.0 100.0; 100.0 100.0;;; 0.51063 0.903457; 0.907086 0.70358]

1.7.2> A[1,1]
2×2×2 Array{Float64, 3}:
[:, :, 1] =
 100.0  100.0
 100.0  100.0

[:, :, 2] =
 0.305401  0.883966
 0.541377  0.487306

1.7.2> A[2,1]
2×2×2 Array{Float64, 3}:
[:, :, 1] =
 100.0  100.0
 100.0  100.0

[:, :, 2] =
 0.267275  0.0161309
 0.800307  0.14134

Is that what you wanted?

6 Likes