Broadcast over a specific dimension of array

Hello, Julians!

I have a 3 dimensional complex F64 array. I want to run the inv() over the inner two dimensions and broadcast along the outer dimension. When I try simple broadcasting inv.(A) I get different results than if I do a for-loop. How do I specify the dimension I want to broadcast over?

The easiest solution would be to have a Vector of Matrix objects. That said, it is very rare for inv to be the best method for solving problems. What’s the usecase?

2 Likes

Broadcast is completely point-wise. That is, it applies to every element in the array, not slices thereof. Thus if you want to use broadcast across slices like this, you need to use some sort of array-of-arrays structure.

2 Likes

The matrix represents the mapping of inputs to outputs of an electrical circuit. The outer dimension represents how that matrix changes at different frequency points. So to flip the mapping from outputs to inputs I just need to invert the matrix, but I want to do it at each frequency so I need to invert the matrix for each matrix along the array. It sounds like a 1D array of 2D arrays would be the best approach.

It looks like mapslices might be the method I was looking for. But, perhaps there is a better way to do what I want than iterating over the outer dimension?

1 Like

It sounds like a 1D array of 2D arrays would be the best approach.

map or broadcast with eachslice can yield a 1D array of 2D arrays:

julia> x = rand(3,3,4);

julia> inv.(eachslice(x, dims=3))
4-element Array{Array{Float64,2},1}:
 [-1.6590312226570239 2.0793019606694854 0.5372029679273337; -0.41725584494784745 -0.16227671804840568 1.4356241190356016; 3.1581583049406534 -1.8836845469856938 -1.3709070403858319]
 [2.2504272336124984 -1.3387872100956475 -1.3175125593090484; -0.6955200736976385 2.3904260659293395 -0.15231503786728723; -1.6588907341843613 -5.504678105886818 6.354210163417449]
 [-6.522265014454815 11.007596635562658 -0.8997442890195426; 5.918653872368978 -6.235209223176827 -0.20055712370467133; 3.1672517888989637 -7.134447496548201 2.0670158780550607]
 [1.257897188987869 0.5534734488273791 -3.7933556798265524; -1.3639997740291203 1.2156395483916482 3.2530995922802717; -1.1021239862109267 -3.2404755836203387 21.985361866828296]

julia> map(inv, eachslice(x, dims=3))
4-element Array{Array{Float64,2},1}:
 [-1.6590312226570239 2.0793019606694854 0.5372029679273337; -0.41725584494784745 -0.16227671804840568 1.4356241190356016; 3.1581583049406534 -1.8836845469856938 -1.3709070403858319]
 [2.2504272336124984 -1.3387872100956475 -1.3175125593090484; -0.6955200736976385 2.3904260659293395 -0.15231503786728723; -1.6588907341843613 -5.504678105886818 6.354210163417449]
 [-6.522265014454815 11.007596635562658 -0.8997442890195426; 5.918653872368978 -6.235209223176827 -0.20055712370467133; 3.1672517888989637 -7.134447496548201 2.0670158780550607]
 [1.257897188987869 0.5534734488273791 -3.7933556798265524; -1.3639997740291203 1.2156395483916482 3.2530995922802717; -1.1021239862109267 -3.2404755836203387 21.985361866828296]

julia> mapslices(inv, x, dims=(1,2))
3×3×4 Array{Float64,3}:
[:, :, 1] =
 -1.65903    2.0793     0.537203
 -0.417256  -0.162277   1.43562
  3.15816   -1.88368   -1.37091

[:, :, 2] =
  2.25043  -1.33879  -1.31751
 -0.69552   2.39043  -0.152315
 -1.65889  -5.50468   6.35421

[:, :, 3] =
 -6.52227  11.0076   -0.899744
  5.91865  -6.23521  -0.200557
  3.16725  -7.13445   2.06702

[:, :, 4] =
  1.2579    0.553473  -3.79336
 -1.364     1.21564    3.2531
 -1.10212  -3.24048   21.9854
3 Likes

Thanks, everyone, for the help.

1 Like

I found this thread because I had the same question, essentialy. I was not satisfied with the answers because I wanted

  • a view
  • high performance

So I created a typestable, zero allocation (so it is fast) method for selectdim. Like selectdim it returns a view. Unlike Julia’s version, its second argument must be a Val.

@generated function Base.selectdim(a,::Val{d},i) where{d}
    precols = ()
    pstcols = ()
    for i = 1:d-1
        precols = (precols...,:)
    end
    for i = 1:ndims(a)-d
        pstcols = (pstcols...,:)
    end
    return quote
        return view(a,$(precols...),i,$(pstcols...))
    end
end

a = randn(2,4,5)
g = selectdim(a,Val(3),2)  # view to a[:,:,2]