Could you file an issue and also the matrix? That way, this won't get lost.
-viral On Saturday, July 11, 2015 at 1:04:47 AM UTC+5:30, Stef Kynaston wrote: > > Dear all, > > I am experiencing some problems with eigs, and wondered if anyone had any > suggestions for a workaround. > > I need to find the smallest positive eigenvalue, lambda, of the > generalised eigenvalue problem > A*v = lambda*B*v > Here, A and B are large, sparse symmetric stiffness matrices from a > (fairly involved) finite element model. I have exactly the same > implementation in Matlab, and have checked that the matrices A and B I > obtain in Matlab match up with the Julia ones. In Matlab, doing > > [v, lambda] = eigs(A, B, 1, sigma), > > > where 1 specifies the number of eigenvalues, and sigma specifies the shift > (I use 1e-12 here, to avoid the zero eigenvalues) returns the "correct" > solution. > > Now, using Julia's eigs, as > > lambda, v = eigs(A, B, nev=1, sigma=1e-12, ritzvec=true) > > > gives me the error "ERROR: Base.LinAlg.PosDefException(0) in eigs at > linalg/arnoldi.jl:38". > > The issue seems to be that B is positive semi-definite, and is numerically > producing some very small negative eigenvalues. Matlab seems able to work > around this and obtains sensible results, but Julia is not playing ball > with this. > > Does anyone have an idea for a way around this? Also interesting is that > there appears to be some randomness in the Julia implementation. > Occasionally eigs doesn't fail and returns a very incorrect answer. Even > more occasionally, it returns the correct answer. > > Many thanks in advance for any insight here, > Stef >