Re: [eigen] Band Lanczos Algorith

On Tue, Dec 23, 2014 at 9:55 PM, Micah Chambers <micahc.vt@xxxxxxxxx> wrote:

Hey everyone. As part of my PhD work I needed an implementation of the Band Lanczos Algorithm, so I implemented it from the "Templates for the Solution of Algebraic Eigenvalue Problems" book (link) using eigen, but its not really a full blown template class. I figured it might be worth expanding and adding to the Unsupported set of libraries. I have a couple of questions about contributing, since it will still be a bit of work to make as generic as the rest of the decomposition classes.

1) Does this already exist somewhere?

AFAIK, nothing like this yet.

2) How does Eigen deal with Hermitian/non-Hermitian naming, since usually a very similar algorithm exists for both. Should I make a separate class or just add a flag for whether the matrix is hermitian?

So far we used to have different classes because the underlying algorithms where really different. Implementation-wise, if the code is the same, you can factorize it through a common internal free function, and if the API is also common, then a common base class will do it. We use the word "SelfAdjoint" instead of Hermitian though.

3) In these sorts of problems its sometimes possible that the Matrix being decomposed can't actually be formulated, so a function needs to passed that takes a vector and returns a vector. Are there any classes that already do this so I can make sure my API looks the same?

What does this function do? A matrix vector product? In general, we request that the abstract matrix class follows the part of the Eigen API required by the algorithm, like implementing rows(), cols(), and operator*(const Eigen::MatrixBase<Derived> &) for iterative solvers as BiCGSTAB.

4) Which branch should I check out to add stuff to?

The default one.

5) Can I just put the tests in unsupported/test?

That's the right place ! You might look at the test/bicgstab.cpp and test/conjugate_gradient.cpp files for some examples which makes use of a common set of testing functions for sparse solvers (in test/sparse_solver.h)

Feel free to ask more questions as the internal mechanisms of Eigen are not yet documented.

cheers,

Gael

Thanks

-Micah

--
Micah Chambers, M.S.
Ahmanson-Lovelace Brain Mapping Center
University of California Los Angeles
635 Charles E. Young Drive SW, Suite 225
Los Angeles, CA 90095-7334
310-206-2101 (Office)
310-206-5518 (Fax)
703-307-7692 (Cell)
email: micahcc@xxxxxxxx