Computing Inverse of a stack of matrices

[bhaveshshrimali]

@stevengj

Use an array of SMatrix (from StaticArrays.jl), which is about 25× faster ( update: 75× faster if I fix my type declaration) than your compute_inv on my machine:

Thanks @stevengj for a superbly informative reply. I wasn’t aware of StaticArrays.jl. This is indeed very fast.

@stevengj

It’s also worthwhile reconsidering what you are doing here. If you are coming from Python/Numpy, then it may have been ingrained that you should break every computation into a sequence of “vectorized” operations like inverting every 3×3 matrix at once.

Absolutely !!

@stevengj

Presumably, you don’t just want the inverses for their own sake, but are planning on using them somehow, for example to solve a huge set of 3×3 equation systems, and potentially do then do further processing. For each 3×3 matrix, you might want to do all processing associated with that matrix (e.g. solve the system of equations, plus any pre/post-processing) before moving on to the next 3×3 matrix.

So this basically ends up in a Finite Element calculation in elasticity where I have a deformation gradient tensor at each quadrature point inside my domain, hence the shape (3,3,1000,4). And since the calculations from hereon depend on the inverse of these (3,3) matrices, I was looking to increase the inversion speed.

@stevengj

(But still use StaticArrays for working with 3×3 systems — that will always be a win.)

Absolutely. Thanks a lot for such an enlightening response.