# Operations on elements of a common index of an array of matrices

**URL:** <https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715>\
**Category:** Performance\
**Tags:** performance, arrays, matrices\
**Created:** [October 3, 2020, 10:29pm UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715 "2020-10-03T22:29:43Z")\
**Posts on this page:** 8\
**Page:** 1

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [October 3, 2020, 10:29pm UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/1 "2020-10-03T22:29:43Z")

</div>

I have an array of matrices. Is it possible to update groups of common indices in place? The code below doesn’t work but I think conveys what I am trying to do. Yes it would be easier to do this with a 3D matrix, but for other reasons it would be very convenient to keep the array of matrices:

```julia
using FFTW

A = [rand(ComplexF64,(100,100)) for _ in 1:1000]

fft!(getindex.(A,1)) #apply fft in-place to the 1000 element grouping consisting of the 1st element of each matrix

```

---

<div class="post-metadata">

**Author:** ![tomerarnon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tomerarnon/32/3170_2.png) [@tomerarnon](https://discourse.julialang.org/u/tomerarnon)\
**Post date:** [October 3, 2020, 10:42pm UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/2 "2020-10-03T22:42:39Z")

</div>

No, I don’t think you can do this with an array of matrices straight out of the box. You can try with `view.(A, 1)`, but considering this is a 1000-element array of 1x1 views with uncorrelated memory address, etc., even if it works (i.e if `fft!` is written generically enough for this sort of input) I expect it will be pretty slow.

---

<div class="post-metadata">

**Author:** ![Henrique\_Becker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/henrique_becker/32/15443_2.png) [@Henrique\_Becker](https://discourse.julialang.org/u/Henrique_Becker)\
**Post date:** [October 3, 2020, 11:22pm UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/3 "2020-10-03T23:22:16Z")

</div>

I think this is a clear case of [_copying data is not always bad_](https://docs.julialang.org/en/v1/manual/performance-tips/#Copying-data-is-not-always-bad).

What I would recommend is:

```julia
using FFTW

A = [rand(ComplexF64,(100,100)) for _ in 1:1000]

setindex.((A,), 1, fft(getindex.(A,1))

```

The broadcast in the `setindex` may, at least, help avoid allocating intermediary data structs.

---

<div class="post-metadata">

**Author:** ![jleman](https://avatars.discourse-cdn.com/v4/letter/j/df788c/32.png) [@jleman](https://discourse.julialang.org/u/jleman)\
**Post date:** [October 4, 2020, 3:49am UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/4 "2020-10-04T03:49:30Z")

</div>

Your suggestion was a good compromise. In the end it was still a significant improvement in my overall speed to start with an array of matrices for initial operations, copy each slice into a 3D matrix for the fft function, then copy back to an array of matrices for subsequent operations.

---

<div class="post-metadata">

**Author:** ![tomerarnon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tomerarnon/32/3170_2.png) [@tomerarnon](https://discourse.julialang.org/u/tomerarnon)\
**Post date:** [October 4, 2020, 10:23am UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/5 "2020-10-04T10:23:40Z")

</div>

I don’t think the broadcasting could prevent allocating in this case, since the `getindex.(...)` is inherently creating a new array (although it’s nice and compact to do it this way!). It seems about the same (maybe a bit slower even) as allocating the intermediate:

```julia
julia> @btime setindex!.($A, fft!(getindex.($A,1)), 1);
  94.378 μs (28 allocations: 25.94 KiB)

julia> @btime let A = $A
           ffta1 = fft!(getindex.(A, 1))
           @inbounds for i in eachindex(A)
               A[i][1] = ffta1[i]
           end
       end
  82.882 μs (27 allocations: 18.00 KiB)

```

---

<div class="post-metadata">

**Author:** ![Henrique\_Becker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/henrique_becker/32/15443_2.png) [@Henrique\_Becker](https://discourse.julialang.org/u/Henrique_Becker)\
**Post date:** [October 4, 2020, 1:06pm UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/6 "2020-10-04T13:06:00Z")

</div>

It is kinda unfair you use `@inbounds` in one and not in the other. For me:

```julia
julia> @btime let A = $A
           ffta1 = fft!(getindex.(A, 1))
           @inbounds for i in eachindex(A)
               A[i][1] = ffta1[i]
           end
       end
  49.373 μs (44 allocations: 18.38 KiB)

julia> @btime @inbounds setindex!.($A, fft!(getindex.($A,1)), 1);
  48.683 μs (45 allocations: 26.31 KiB)

```

~~Also, considering the number of allocations in both your benchmark and mine, it is clear the broadcast remove a single allocation, that is actually, 30% of the memory used… The difference in time is small, but my suggestion has a slightly better time in my machine, uses less memory~~ and is more compact…

I now perceive my first code was wrong, this is the reason @jleman found it to be considerably faster. Explain in my next post.

---

<div class="post-metadata">

**Author:** ![tomerarnon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tomerarnon/32/3170_2.png) [@tomerarnon](https://discourse.julialang.org/u/tomerarnon)\
**Post date:** [October 4, 2020, 2:41pm UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/7 "2020-10-04T14:41:11Z")

</div>

I’m seeing the same relationship annotating or not annotating both… `@inbounds` is actually a microsecond slower on average even… I’m running on version 1.5.1.

```julia
julia> begin 
       @btime let A = $A
           ffta1 = fft!(getindex.(A, 1))
           @inbounds for i in eachindex(A)
               A[i][1] = ffta1[i]
           end
       end

       @btime @inbounds setindex!.($A, fft!(getindex.($A,1)), 1);
       end;
  82.740 μs (27 allocations: 18.00 KiB)
  93.370 μs (28 allocations: 25.94 KiB)

julia> begin 
       @btime let A = $A
           ffta1 = fft!(getindex.(A, 1))
           for i in eachindex(A)
               A[i][1] = ffta1[i]
           end
       end

       @btime setindex!.($A, fft!(getindex.($A,1)), 1);
       end;
  81.483 μs (27 allocations: 18.00 KiB)
  92.707 μs (28 allocations: 25.94 KiB)

```

> [@Henrique\_Becker](#):
>
> Also, considering the number of allocations in both your benchmark and mine, it is clear the broadcast remove a single allocation, that is actually, 30% of the memory used

Am I misunderstanding something? From the looks of it, the `setindex!.` method is the one allocating more memory.

---

<div class="post-metadata">

**Author:** ![Henrique\_Becker](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/henrique_becker/32/15443_2.png) [@Henrique\_Becker](https://discourse.julialang.org/u/Henrique_Becker)\
**Post date:** [October 4, 2020, 3:05pm UTC](https://discourse.julialang.org/t/operations-on-elements-of-a-common-index-of-an-array-of-matrices/47715/8 "2020-10-04T15:05:06Z")

</div>

In my first post I wrote:

```julia
setindex!.((A,), 1, fft(getindex.(A,1))

```

what is wrong. This will write the value 1 in the position given by the fft, what makes no sense:

```julia
A[fft(getindex.(A,1))] .= 1 

```

@jleman has probably corrected the order of the parameters to:

```julia
setindex!.((A,), fft(getindex.(A,1), 1)

```

But this is:

```julia
A[1] .= fft(getindex.(A,1)

```

What explains some gain in performance, it is overwriting contiguous positions in the first matrix with the the `fft` values, instead of overwriting the first value in each of the matrices. (I am not sure, however, how it did not have a type problem, as the `fft` return is a vector, not a matrix.)

The correct version is:

```julia
@btime @inbounds setindex!.($A, fft!(getindex.($A,1)), 1);

```

That somehow is slightly faster than looping to me, even if it makes one extra allocation. (I will need to investigate further, it does not make sense to me why it does make an extra allocation nor why it is faster/equivalent despite that.)

I am sorry for the noise.