Hello there, I’m trying to figure out how to do matrix operations for 2D slices of a 3D array whose shape is (3,3,N) where N is a large number (e.g. N=10^6). I’d like to write a vectorized code of C=A+B’ using loop fusion as explained in More Dots: Syntactic Loop Fusion in Julia for the first 1024 …

Since you’re working with such small slices, the cost of creating all those views (and doing the associated bookkeeping on transposes) is hard to amortize. Adding @inbounds speeds up the expanded loop even further, but if you can use StaticArrays, I think they give you the syntax you’re looking for: …

Thanks for the instant reply and the cool demonstration of StaticArrays! I got similar results with your code. However, I’m a bit confused about its applicability. I thought StaticArrays are stack allocated and so they have limited memory quota. Typing ulimit -a in my MacBook suggests that my stack …

Good question. StaticArrays can be stored on the stack but since they’re isbits , they can also be stored immediately inline within an array. This means that the array [@SVector[1,2],@SVector[3,4]] is stored in memory simply as the sequence of the four Ints smashed together. This is in contrast to…

Thanks for the explanation. So in this case the entire s1 is located on heap which has no 8MB limit, right? The remaining issue is that I usually obtain s1 from hdf5 files where data has the format of a 3D array (e.g. s1 = h5read("data_3darray.h5", "group")). Therefore, I need to convert this 3x3x10…

This can be done without copying http://juliaarrays.github.io/StaticArrays.jl/latest/pages/api.html#Arrays-of-static-arrays-1

Worked like a charm! For reference, the syntax is s1 = reinterpret(SMatrix{3,3,Float64,9}, reshape(v1,9,1000000)). It doesn’t seem to be compatible with the mutable static matrix though (it fails if I replace SMatrix with MMatrix in the above line), which is a pity… but I guess for an array of matri…

I asked a question a while ago, which has some similar themes to your question. You might find some of the answers I got to be useful. [image] Multiply many-matrices by many-vectors Performance All, I have a critical portion of code that essentially involves the multipl…

Thanks for the reference. It seems that my best option is to initialize a large array of MMatrix m = zeros(MMatrix{3,3,Float64,9},N) and then assign each element with a for-loop.

I would have opened a new thread for this. Could you ensure that your code is runnable? Which package is h5write from, for example? As generic advice, have you checked the functions for type stability (e.g. @code_warntype read1())? I kind of expect h5read to be type-unstable.

Just opened a new topic! Thanks!

Note that these will (afaik) not be stored inline. Is use the SMatrices. Note that while an SMatrix is immutable, so is a Float64, but you can still change the contents of an Array holding (either of) them.

Efficient vectorized matrix operations for 2D matrix slices of a 3D array?

New to Julia

GregVernon May 1, 2019, 5:47pm 8

I asked a question a while ago, which has some similar themes to your question. You might find some of the answers I got to be useful.

1 Like

Topic		Replies	Views
Optimizing iteration over slices of multiple matrices Performance question , array	14	1986	December 30, 2018
Matrix-by-(slice of)vector multiplication with limited allocation New to Julia question	6	857	September 11, 2020
Convert a 3D Array to 2D array of vectors General Usage	7	1353	August 4, 2021
Recent developments in StaticArrays.jl land Package Announcements performance , array , linearalgebra , arrays	12	3797	March 11, 2021
Performance of Matrix{T} vs Vector{Vector{T}} New to Julia question , performance	14	768	June 8, 2023

Efficient vectorized matrix operations for 2D matrix slices of a 3D array?

Related topics