problem with the last test in test/arnoldi.jl in the 0.6.2 version

[aitzkora]

I run multiple times the following test

# Symmetric generalized with singular B
let n = 10
    k = 3
    A = randn(n,n); A = A'A
    B = randn(n,k);  B = B*B'
    @test sort(eigs(A, B, nev = k, sigma = 1.0)[1]) ≈ sort(eigvals(A, B)[1:k])
end

and I get an error, to reproduce the error I stored the matrices B and A in files and rexecute the test test_arnoldi

using Base.Test
A = readdlm("A.txt")
B = readdlm("B.txt")
m,n  = size(A)
@assert m == n
k = 3
@assert p == n
@assert n == 10
@test sort(eigs(A, B, nev = k, sigma = 1.0)[1]) ≈ sort(eigvals(A, B)[1:k])

A.txt B.txt

[RalphAS]

You should select entries 1:k from eigvals after sorting.

There will still be problems, since IIUC eigvals(A,B) for symmetric arguments calls LAPACK.sygvd!, which (according to LAPACK documentation) assumes that B is positive definite.

Let's ask @andreasnoack whether this is a bug or something that should be documented.

[andreasnoack]

For once, this isn't an issue with our ARPACK wrappers. As @RalphAS said, you'd have to sort before you extract the values. See

julia> eigvals(A, B)
10-element Array{Float64,1}:
 Inf
   0.905373
   0.0877561
   0.047194
 Inf
 Inf
 Inf
 Inf
 Inf
 Inf

Even though A and B are symmetric eigvals(A,B) uses the general solver. As @RalphAS also points out, the symmetric solver requires that B is PD so

julia> eigvals(Symmetric(A), Symmetric(B))
ERROR: Base.LinAlg.PosDefException(15)
Stacktrace:
 [1] chkposdef at ./linalg/lapack.jl:46 [inlined]
 [2] sygvd!(::Int64, ::Char, ::Char, ::Array{Float64,2}, ::Array{Float64,2}) at ./linalg/lapack.jl:4910
 [3] eigvals!(::Symmetric{Float64,Array{Float64,2}}, ::Symmetric{Float64,Array{Float64,2}}) at ./linalg/symmetric.jl:448
 [4] eigvals(::Symmetric{Float64,Array{Float64,2}}, ::Symmetric{Float64,Array{Float64,2}}) at ./linalg/eigen.jl:409

[aitzkora]

thanks for the reply, the aim of my question was to modify this test, because it generate errors when a friend try to make package for guix. It seems it was solved in the developpement versions fixing the seed in the file stdlib/IterativeEigensolvers/test/runtests.jl

@testset "Symmetric generalized with singular B" begin
        srand(127)
        n = 10
        k = 3
        A = randn(n,n); A = A'A
        B = randn(n,k); B = B*B'
        @test sort(eigs(A, B, nev = k, sigma = 1.0)[1]) ≈ sort(eigvals(A, B)[1:k])
    end