Re: [eigen] Why not Matrix<T, const int rows, const int cols> ??

[Benoît Jacob <jacob@xxxxxxxxxxxxxxx>]

My reply was obfuscated by an unwanted paste near the end, making it 
incomprehensible, so I resend, please ignore previous reply.

Le vendredi 05 janvier 2007 13:59, Andre Krause a écrit :
> Why is EIGEN restricted to square matrices?

Ah, very good question! At the beginning, Eigen had rectangular matrices,
 with template arguments Width and Height as you suggest. Then it created
complications with the multiplication operators. Eigen does not only provide
fixed-size matrices, it also provides dynamic-size matrices (class MatrixX)
and both inherit the same MatrixBase using a 'curiously recurring template
pattern'. The problem is that multiplication of rectangular matrices works
differently for fixed-size matrices than it does for dynamic-size matrices.
Indeed, for fixed-size matrices this operator takes an additional template
parameter, int other_width. The type of the 'other' matrix is not the same as
the type of *this. On the other hand, with dynamic-size matrices, there is no
need for an additional template parameter and all matrices are of the same
type, MatrixX<T>.

The only solution is to take the matrix-matrix product out of MatrixBase and
implement it separately in Matrix and in MatrixX. Good, that's what I did.
But there are other operator* and multiply() methods in MatrixBase, for
matrix-vector product. They were now hidden by the operator* and multiply()
implemented in the child class. So I had to remove them also from MatrixBase
and implement them also in Matrix and MatrixX separately.

At this point I started wondering whether all that was really useful. It was
really starting to make Eigen less elegant, fatter, and partly defeated the
point of having this curiously recurring template pattern. The goal of Eigen
being to cover the needs of applications, I started to ask people around me
what matrix width and height they planned to actually use in their apps. All
answers mentionned only square matrices. Same for me. Only one person
mentionned rectangular matrices, that was for linear regression, but since
then I have implemented linear regression in Eigen using only square matrices
(by unfolding the rectangular matrix product one needs to compute at some
point. By the way that allowed for a +50% speed optimization).

So I removed all the rectangular stuff, and was amazed at how much it
simplified the source code, the API, the documentation, and the tests (the
test app already weighs 2 MB, think how much it'd be with all the
width/height combinations).

Do you have a real need for rectangular matrices?

> by the way - is there any code that shows how to implement
> pseudo-inverse calculation using eigen?

Hmm, the best way to do that is by singular values decomposition, which Eigen
doesn't do. However, Eigen does LU decomposition with complete pivoting. I
don't know, maybe that could help here. It decomposes your matrix as a product

PLUQ

where P and Q are permutation matrices, hence unitary, hence it is enough to
compute the pseudoinverse of LU. As L and U are triangular, I guess their
pseudoinverses are easily computed, and it remains to deduce from that the
pseudoinverse of LU. I don't know if that's feasible.

Hope that helps,
Benoit

-------------------------------------------------------

Attachment: pgpbnIcqF28Y0.pgp
Description: PGP signature