[Feature Request] Eigenvalues of small asymmetric matrices without allocation?

[akriegman]

Recently I found myself trying to find the eigenvalues of 1 million 3x3 rotations, and I was very frustrated by how slow it was going, so I turned to StaticArrays. It looks like StaticArrays provides methods to speed this up for hermitian matrices, but not for asymmetric matrices like rotations, so I was out of luck. Would anyone else be interested in seeing this feature added to StaticArrays?

One obstacle that could be a problem is the return type. Eigenvalues of general real matrices can be either real or complex. The standard library handles this by having eigen, eigvals, and eigvecs return Union types over complex and real, e.g. Union{Vector{Float64}, Vector{ComplexF64}} for eigvals. I think this is a problem in general because of runtime type reflection and unnecessary allocations due to the unknown size of the result. If I'm not mistaken this will especially be a problem if we want to add this feature to StaticArrays since we can't have a result of unknown size without allocating.

So there's a few things we can do:

• Have the eigen functions accept an argument for the scalar type in the return value, i.e. eigvals(ComplexF64, A) would return Vector{ComplexF64}. We could make this argument mandatory, or we could have the type default to complex, or have it default to real and throw a DomainError when the eigenvalues are complex.
• We could implement a no-allocation eigen decomposition algorithm just for complex matrices, and continue to default to the standard library for general real matrices. This way someone can cast their matrix to a complex matrix when they don't expect their eigenvalues to be real and they don't want to pay for allocations and type reflection.
• Have separate methods based on the return type: eigvalsreal, eigvalscomplex, etc.

I think the first option would make the most sense. What do people think? I'm also wondering if for static arrays larger than 3x3 there's a way we can call the LAPACK routines without allocating, which would speed things up for medium sized matrices as well. But anyways, this all sounds like a big project. Would people want these features?

[hyrodium]

Recently I found myself trying to find the eigenvalues of 1 million 3x3 rotations, and I was very frustrated by how slow it was going, so I turned to StaticArrays.

Although it does not support general matrix eigenvalue problems, Rotations.jl supports eigenvalues of 3×3 rotation matrix. Please check JuliaGeometry/Rotations.jl#170.

[fredrikekre]

Have the eigen functions accept an argument for the scalar type in the return value, i.e. eigvals(ComplexF64, A) would return Vector{ComplexF64} [...]

Another option would be to require the user to pass the matrix with the correct type, i.e. passing this argument implicitly via the data.

[stevengj]

Duplicate of #753?