Re: [eigen] Divide & conquer SVD algorithm

Thank you for the quick and developped answers. We discussed our options and decided to start with functionality C and continue with B.

Le 27/05/2014 11:51, Gael Guennebaud a écrit :

You're very welcome to contribute to Eigen.

On Tue, May 27, 2014 at 11:23 AM, Jitse Niesen <jitseniesen@xxxxxxxxx> wrote:

On Tue, 27 May 2014, Tristan Sahel wrote:

Then, we are planning to continue the work of our fellow students of the Ensimag from last year, who were trying to improve the SVD by developping a divide & conquer algorithm. Their work was very promising and their algorithm was slightly faster than the Jocobi SVD algorithm. However, they did not have time to finish the implementation of their divide & conquer: they were missing one step, the last one, which they explained step by step in their todoBDcsvd.txt document. If we manage to finish the job, the calculation time should become o(n^2).

I hate to scupper your plan, but this is already done. At the moment we have a complete but basic implementation of the divide and conquer algorithm in the experimental section. It does need some improvements (there seems to be at least one bug and several inefficiencies).

Indeed, the D&C SVD implementation is now sufficiently advanced, and a student project to polish it would not make much sense.

We're abandoning the idea of another SVD algorithm then, and instead we'll be focusing on your suggestions.

I do not know exactly what is required for your project, but here are some ideas along similar lines. What I would consider the standard algorithm for the SVD based on Golub and Kahan using bi-diagonalization is not implemented in Eigen (we did have an implementation in the past but we removed it as we were not happy with it), so you could do that. You could also try implementing a divide-and-conquer algorithm for computing eigenvalues. Another idea is to implement an method for the polar decomposition.

I'm not sure that a third SVD implementation would be really useful. I've more useful ideas in mind:

A - Improve the column-pivoting QR factorization class: In my opinion, column-pivoting QR is the method of choice for solving least-square and rank-deficient problems. It is rank revealing, in theory it allows for solving for the minimal norm solution, and it is significantly faster than SVD. However, the current implementation is somewhat naive and two important improvements are needed to make is really useful:

A1 - Regarding performance, we need to leverage matrix-matrix operations (as in the non-pivoting variant). This basically means rewriting the algorithm to work on sub-blocks of the matrix.

A2 - For rank-deficient problem, we need to add the possibility to perform a full orthogonalization: A*P = Q*T*Z where Z is unitary, and T is of the following form:

| T1 0 |

| 0 0 |

with T1 upper triangular and full-rank.

A3 - bonus: add functions to extract the kernel and image of the matrix from the column-pivoting factorization.

We may try and implement this functionality if we have spare time at the end of our project.

B - Implement eigenvalue decompositions for sparse matrices. A must to have for a library called Eigen ;), and really useful!

This will become the main part of our project. We want to use a block Lanczos algorithm (we have barely researched the subject at this point but it seems to be a good solution to tackle the sparse matrices problem).

C - Implement new eigenvalue decompositions for dense matrices allowing to compute only a subset of the eigenvalues/vectors. Useful for large eigenvalue problems. This requires completely different algorithms than the one we currently have which can only compute all eigenvalues.

We intend to do this first, as a way to get used to the library. We want to use the power iteration method to find the eigenpair with the largest eigenvalue.

Thank you for your help,
Imène, Bediss and Tristan.