Feature Request -- support for Lanczos about arnoldimethod.jl HOT 7 CLOSED

Thanks!

I have to figure lanczos out, but if I remember correctly some articles advocate using Arnoldi even in the case of symmetric / Hermitian matrices, because of stability. Is Arnoldi too slow for your use case?

from arnoldimethod.jl.

Thanks for the quick response. Admittedly, I am not an expert in numerical matrix algebra but only a user of ARPACK for diagonalizing quantum many-body Hamiltonians (hermitian and sparse matrices). Intuitively, Lanczos uses the information about the symmetric/hermitian structure of the matrix so that supposedly there could only be some gain. But I might be mistaken. I can try testing (speed/correctness) your package for some of my use cases.

from arnoldimethod.jl.

In my experience the real advantage of Lanczos comes from a reduction of memory usage. It comes at cost of numerical stability as haampie points out. Once one eigenvalue converges cancellation effects will come in and lessen the accuracy of the other eigvals. On the other hand, in most engineering applications full precision is not needed.

from arnoldimethod.jl.

Even just for the memory it would be nice to have.

from arnoldimethod.jl.

I would propose to add support this by making Arnoldi and abstract class and let

struct FullArnoldi{T,TV<:StridedMatrix{T},TH<:StridedMatrix{T}} <: Arnoldi
    V::TV 
    H::TH
...

be the standard behavior, and define a new class which only stores the recent vectors

struct Lanczos{T,TV<:StridedVector{T},TH<:StridedMatrix{T}} <: Arnoldi
    v0::TV 
    v1::TV 
    v2::TV 
    H::TH

The expansion could then easily be handled with dispatch:
Arnoldi:

function iterate_arnoldi!(A, arnoldi::FullArnoldi{T}, range::UnitRange{Int}) where {T}

and Lanczos

function iterate_arnoldi!(A, arnoldi::Lanczos{T}, range::UnitRange{Int}) where {T}

from arnoldimethod.jl.

On second thought, there is a decision to made here. 1) The implicitly restarted Lanczos as described in the eigenvalue templates http://www.netlib.org/utk/people/JackDongarra/etemplates/node118.html does not have the advantage of reduced memory, since one needs to store the whole basis to do the restart. It does reduce the orthogonalization CPU-time. 2) If one only wants eigenvalues (not eigenvectors) and do not use restart, one can get the reduced memory advantage. 3) Another solution is to use double sweeps, which would double the computation cost, but have the reduced memory. Support for 1), 2) or 3)? Or all of them?

from arnoldimethod.jl.