mat: Redesign Eigen

[btracey]

Right now, Eigen stores the eigenvalue matrix as a *Dense. Conceptually, this is wrong, because in the general case, it's really a complex matrix. The code instead uses extra rows and columns to represent complex values, which means that the physical size of the matrix is larger than its conceptual size. It would be better to just have the values just be complex. Acyclic imports mean we can't have mat128 depend on mat64 and vis-versa, but we can have mat64 import blas128. Eigen can return the blas structs, which can then be converted into mat128 matrices. This will take some work to interface with lapack routines, but not doing this work means we just push the problem onto the user (making them do the conversion from 2x2 float64 blocks to complex values). The eigenvalues by blas128.Banded

We should also reconsider the function signatures. First of all, it seems easier to me to have it be
Eigen(a *Dense, epsilon) (eigenvalues blas128.Banded, eigenvectors *Dense)
I'm not sure what the EigenFactors struct helps with
Secondly, we should also have
EigenSym(a *Symmetric, epsilon) (eigenvalues *Diagonal, eigenvectors *Dense)
as the eigenvalues of a symmetric matrix are real.
The asymmetry between the two is unfortunate (Diagonal vs. Banded), however there is no Diagonal matrix in BLAS. We could add a definition both into blas64 and blas128, or we could plausibly have a non-blas package containing the definition of RawDiagonal (name tbd) for both complex and real diagonal matrices. We could also leave as-is, and solve the problem in a different manner (there are a few choices).

[kortschak]

Other choices? I don't like the idea of giving a result vallue that cannot be used.

[btracey]

I'm not sure what you're referring to by "result value that cannot be used". Which can't be used? Arguably the current return from Eigen can't be used either, as to construct the actual eigenvalues you have to look through the *Dense matrix, see which blocks are 1x1 and which are 2x2, and construct the cmplx128 yourself.

[kortschak]

Moving from one non-usable to another doesn't seem like a win.

Can we return two *mat64.Dense, one real and the other imaginary, then provide a helper in cma128 (when it happens) takes two *mat64.Dense and returns a *cmat128.Dense?

[btracey]

How about we just return []cmplx128 and []float64 for the eigenvalues.
There are many cases where you'd rather have them as a list than a matrix
(stability, spectral analysis), and it solves the import issue. If the
user needs it as a diagonal matrix, it's easy to martial it into one
On Jan 13, 2015 12:55 PM, "Dan Kortschak" notifications@github.com wrote:

Moving from one non-usable to another doesn't seem like a win.

Can we return two *mat64.Dense, one real and the other imaginary, then
provide a helper in cma128 (when it happens) takes two *mat64.Dense and
returns a *cmat128.Dense?

—
Reply to this email directly or view it on GitHub
https://github.com/gonum/matrix/issues/79#issuecomment-69818112.

[kortschak]

SGTM

[jonlawlor]

I had a thought earlier - the mat64 interfaces don't allow complex128s to be returned, but that doesn't prevent you from having a []complex128 backing vector. It would be possible to create a matrix type with a backing of complex, that can still satisfy the matrix interfaces, except that when things like "At" are called, it would return either the real component or the magnitude instead of the complex values. If you add a conversion function to cmat128.Dense (or whatever) then you could get at the imaginary values directly. Alternatively, instead of just having an At, you could also have an AtComplex.

[btracey]

Partially fixed by gonum/matrix#309 . Full fix needs gonum/matrix#308

[kortschak]

So what is needed here?