The 0.5 behavior is the same as Fortran, and is technically correct. ARPACK cannot create more (orthogonal) Krylov vectors than the rank of the matrix A, so your example needs the additional keyword ncv=2. This doesn't seem to be adequately documented in arpack-ng either. The 0.4 behavior looks unreliable.
On Saturday, August 6, 2016 at 2:03:21 AM UTC-4, Madeleine Udell wrote: > > Actually: I'm not sure we should chalk the error up to ARPACK. Julia 0.4 > is able to produce a (correct, I think) answer to > > eigs(A, C) > > but 0.5 gives an ARPACK error. I don't suppose ARPACK changed between > Julia versions...!? > > On Sat, Aug 6, 2016 at 1:54 AM, Madeleine Udell <[email protected] > <javascript:>> wrote: > >> Andreas, thanks for the investigation. I'll use 0.5 for now, and hope the >> real problems I encounter are within the capabilities of ARPACK. >> >> It's embarrassing to be bested by Matlab... >> >> On Fri, Aug 5, 2016 at 9:23 PM, Andreas Noack <[email protected] >> <javascript:>> wrote: >> >>> I've looked a bit into this. I believe there is a bug in the Julia >>> wrappers on 0.4. The good news is that the bug appears to be fixed on 0.5. >>> The bad news is the ARPACK seems to have a hard time with the problem. I get >>> >>> julia> eigs(A,C,nev = 1, which = :LR)[1] >>> ERROR: Base.LinAlg.ARPACKException("unspecified ARPACK error: -9999") >>> in aupd_wrapper(::Type{T}, >>> ::Base.LinAlg.#matvecA!#67{Array{Float64,2}}, ::Base.LinAlg.##64#71{Array >>> {Float64,2}}, ::Base.LinAlg.##65#72, ::Int64, ::Bool, ::Bool, ::String, >>> ::Int64, ::Int64, ::String, : >>> :Float64, ::Int64, ::Int64, ::Array{Float64,1}) at ./linalg/arpack.jl:53 >>> in #_eigs#60(::Int64, ::Int64, ::Symbol, ::Float64, ::Int64, ::Void, >>> ::Array{Float64,1}, ::Bool, ::B >>> ase.LinAlg.#_eigs, ::Array{Float64,2}, ::Array{Float64,2}) at >>> ./linalg/arnoldi.jl:271 >>> in (::Base.LinAlg.#kw##_eigs)(::Array{Any,1}, ::Base.LinAlg.#_eigs, >>> ::Array{Float64,2}, ::Array{Floa >>> t64,2}) at ./<missing>:0 >>> in #eigs#54(::Array{Any,1}, ::Function, ::Array{Float64,2}, >>> ::Array{Float64,2}) at ./linalg/arnoldi. >>> jl:80 >>> in (::Base.LinAlg.#kw##eigs)(::Array{Any,1}, ::Base.LinAlg.#eigs, >>> ::Array{Float64,2}, ::Array{Float6 >>> 4,2}) at ./<missing>:0 >>> >>> and since SciPy ends up with the same conclusion I conclude that the >>> issue is ARPACK. Matlab is doing something else because they are able to >>> handle this problem. >>> >>> Given that 0.5 is almost release, I'll not spend more time on the issue >>> on 0.4. Thought, if anybody is able to figure out what is going on, please >>> let us know. >>> >>> On Friday, August 5, 2016 at 8:47:26 AM UTC-4, Madeleine Udell wrote: >>>> >>>> Setting `which=:LR, nev=1` does not return the generalized eigenvalue >>>> with the largest real parts, and does not give a warning or error: >>>> >>>> n = 10 >>>> C = eye(n) >>>> A = zeros(n,n) >>>> A[1] = 100 >>>> A[end] = -100 >>>> @assert eigs(A, C, nev=1, which=:LR)[1][1] == maximum(eigs(A, C)[1]) >>>> >>>> Am I expected to set nev greater than the number of eigenvalues I truly >>>> desire, based on my intuition as a numerical analyst? Or has eigs broken >>>> its implicit guarantee? >>>> >>>> >> >> >> -- >> Madeleine Udell >> Assistant Professor, Operations Research and Information Engineering >> Cornell University >> https://people.orie.cornell.edu/mru8/ >> (415) 729-4115 >> > > > > -- > Madeleine Udell > Assistant Professor, Operations Research and Information Engineering > Cornell University > https://people.orie.cornell.edu/mru8/ > (415) 729-4115 >
