Jutho / KrylovKit.jl

Krylov methods for linear problems, eigenvalues, singular values and matrix functions
Other
284 stars 37 forks source link

Geneigsolve failing for matrices with larger elements #91

Open JoeyT1994 opened 1 month ago

JoeyT1994 commented 1 month ago

The following code:

using Random
using LinearAlgebra: isposdef, ishermitian
using KrylovKit: geneigsolve

Random.seed!(1234)
n = 8
ω = 100
B = randn((n,n))
B = ω *(B * B')

A = randn((n,n))
A = ω *(A + A')

@assert ishermitian(A)
@assert isposdef(B)

x0 = randn((n,))
vals, vecs, info = geneigsolve((A,B), x0, 1, :SR; isposdef = true, ishermitian = true)

@show first(vals)

fails for ω > 50 due to error: ERROR: LoadError: PosDefException: matrix is not positive definite; Cholesky factorization failed.

It seems that for matrices with sufficiently large elements the posdef flag gets raised even though the matrices involved are of the correct form (A is symmetric and B is posdef)? I assume this is related to some sort of numerical precision/ tolerance in the check that is not relative to the norm of the matrix itself.

Is there any way to make the algorithm more tolerant with respect to that check? Any help would be greatly appreciated.

ryanlevy commented 1 month ago

I just ran into this problem as well. I suspect the issue is only when krylovdim==n (the default is 30), then something is wrong with the HA/HB matrices. For example

vals, vecs, info = geneigsolve((A,B), x0, 1, :SR; isposdef = true, ishermitian = true, krylovdim=7)

for your example works, just anything under 8.

JoeyT1994 commented 1 month ago

Yeah I noticed that 8 seems to be the threshold for a range of cases which is rather peculiar! I am keeping the krylovdim lower for now