EigSolver.jl
``````using IterativeSolvers, KrylovKit, Parameters, Arpack, LinearAlgebra
# In this file, we regroud a way to provide eigen solvers

abstract type EigenSolver end

# the following function returns the n-th eigenvectors computed by an eigen solver. This function is necessary given the different return types each eigensolver has
getEigenVector(eigsolve::EigenSolver, vecs, n::Int) = vecs[:, n]

####################################################################################################
# Solvers for default \ operator (backslash)
####################################################################################################
"""
The struct `Default` is used to  provide the backslash operator to our Package
"""
struct Default_eig <: EigenSolver end

function (l::Default_eig)(J, nev::Int64)
# I put Array so we can call it on small sparse matrices
F = eigen(Array(J))
I = sortperm(F.values, by = x-> real(x), rev = true)
return F.values[I[1:nev]], F.vectors[:, I[1:nev]]
end

# case of sparse matrices
struct Default_eig_sp <: EigenSolver end

function (l::Default_eig_sp)(J, nev::Int64)
λ, ϕ = Arpack.eigs(J, nev = nev, which = :LR)
return λ, ϕ
end
####################################################################################################
# Solvers for IterativeSolvers
####################################################################################################
@with_kw struct eig_IterativeSolvers{T} <: EigenSolver
tol::T = T(1e-4)        # tolerance for solver
restart::Int64 = 200    # number of restarts
maxiter::Int64 = 100
N = 0                   # dimension of the problem
verbose = true
log = true
end

function (l::eig_IterativeSolvers{T})(J, nev::Int64) where T
# for now, we don't have an eigensolver for non hermitian matrices
@assert 1==0 "Not implemented: IterativeSolvers does not have an eigensolver yet!"
return res[1], length(res)>1, res[2].iters
end
####################################################################################################
# Solvers for KrylovKit
####################################################################################################
@with_kw struct eig_KrylovKit{T} <: EigenSolver
dim::Int64 = KrylovDefaults.krylovdim        # Krylov Dimension
tol::T = T(1e-4)        # tolerance for solver
restart::Int64 = 200    # number of restarts
maxiter::Int64 = KrylovDefaults.maxiter
verbose::Int = 0
which = :LR
end

function (l::eig_KrylovKit{T})(J, nev::Int64) where T
@assert typeof(J) <:  AbstractMatrix
vals, vec, info = KrylovKit.eigsolve(J, nev, l.which;  verbosity = l.verbose, krylovdim = l.dim, maxiter = l.maxiter, tol = l.tol)
return vals, vec, true, info.numops
end

getEigenVector(eigsolve::eig_KrylovKit{T}, vecs, n::Int) where T = vecs[n]

# Matrix-Free version, needs to specify an example of rhs x₀
@with_kw struct eig_MF_KrylovKit{T, vectype} <: EigenSolver
dim::Int64 = KrylovDefaults.krylovdim        # Krylov Dimension
tol::T = T(1e-4)        # tolerance for solver
restart::Int64 = 200    # number of restarts
maxiter::Int64 = KrylovDefaults.maxiter
verbose::Int = 0
which = :LR
x₀::vectype
end

function (l::eig_MF_KrylovKit{T, vectype})(J, nev::Int64) where {T, vectype}
vals, vec, info = KrylovKit.eigsolve(J, l.x₀, nev, l.which;  verbosity = l.verbose, krylovdim = l.dim, maxiter = l.maxiter, tol = l.tol)
return vals, vec, true, info.numops
end

getEigenVector(eigsolve::eig_MF_KrylovKit{T, vectype}, vecs, n::Int) where {T, vectype} = vecs[n]
``````