Skip to content

Commit

Permalink
Deprecate eig(...) in favor of eigfact(...) and factorization destruc…
Browse files Browse the repository at this point in the history
…turing.
  • Loading branch information
Sacha0 committed May 13, 2018
1 parent 76e77d8 commit 78344d1
Show file tree
Hide file tree
Showing 14 changed files with 108 additions and 103 deletions.
21 changes: 10 additions & 11 deletions stdlib/LinearAlgebra/docs/src/index.md
Original file line number Diff line number Diff line change
Expand Up @@ -198,16 +198,16 @@ Legend:

### Matrix factorizations

| Matrix type | LAPACK | [`eig`](@ref) | [`eigvals`](@ref) | [`eigvecs`](@ref) | [`svd`](@ref) | [`svdvals`](@ref) |
|:------------------------- |:------ |:------------- |:----------------- |:----------------- |:------------- |:----------------- |
| [`Symmetric`](@ref) | SY | | ARI | | | |
| [`Hermitian`](@ref) | HE | | ARI | | | |
| [`UpperTriangular`](@ref) | TR | A | A | A | | |
| [`LowerTriangular`](@ref) | TR | A | A | A | | |
| [`SymTridiagonal`](@ref) | ST | A | ARI | AV | | |
| [`Tridiagonal`](@ref) | GT | | | | | |
| [`Bidiagonal`](@ref) | BD | | | | A | A |
| [`Diagonal`](@ref) | DI | | A | | | |
| Matrix type | LAPACK | [`eigfact`](@ref) | [`eigvals`](@ref) | [`eigvecs`](@ref) | [`svd`](@ref) | [`svdvals`](@ref) |
|:------------------------- |:------ |:----------------- |:----------------- |:----------------- |:------------- |:----------------- |
| [`Symmetric`](@ref) | SY | | ARI | | | |
| [`Hermitian`](@ref) | HE | | ARI | | | |
| [`UpperTriangular`](@ref) | TR | A | A | A | | |
| [`LowerTriangular`](@ref) | TR | A | A | A | | |
| [`SymTridiagonal`](@ref) | ST | A | ARI | AV | | |
| [`Tridiagonal`](@ref) | GT | | | | | |
| [`Bidiagonal`](@ref) | BD | | | | A | A |
| [`Diagonal`](@ref) | DI | | A | | | |

Legend:

Expand Down Expand Up @@ -334,7 +334,6 @@ LinearAlgebra.lqfact
LinearAlgebra.lq
LinearAlgebra.bkfact
LinearAlgebra.bkfact!
LinearAlgebra.eig
LinearAlgebra.eigvals
LinearAlgebra.eigvals!
LinearAlgebra.eigmax
Expand Down
1 change: 0 additions & 1 deletion stdlib/LinearAlgebra/src/LinearAlgebra.jl
Original file line number Diff line number Diff line change
Expand Up @@ -80,7 +80,6 @@ export
diagind,
diagm,
dot,
eig,
eigfact,
eigfact!,
eigmax,
Expand Down
40 changes: 40 additions & 0 deletions stdlib/LinearAlgebra/src/deprecated.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1282,3 +1282,43 @@ function lu(x::Number)
"`LU` object (`lup = lufact(x)`)."), :lu)
return first.((lufact(x)...,))
end

# deprecate eig(...) in favor of eigfact and factorization destructuring
export eig
function eig(A::Union{Number, StridedMatrix}; permute::Bool=true, scale::Bool=true)
depwarn(string("`eig(A[, permute, scale])` has been deprecated in favor of ",
"`eigfact(A[, permute, scale])`. Whereas `eig(A[, permute, scale])` ",
"returns a tuple of arrays, `eigfact(A[, permute, scale])` returns ",
"an `Eigen` object. So for a direct replacement, use ",
"`(eigfact(A[, permute, scale])...,)`. But going forward, consider ",
"using the direct result of `eigfact(A[, permute, scale])` instead, ",
"either destructured into its components ",
"(`vals, vecs = eigfact(A[, permute, scale])`) ",
"or as an `Eigen` object (`eigf = eigfact(A[, permute, scale])`)."), :eig)
return (eigfact(A, permute=permute, scale=scale)...,)
end
function eig(A::AbstractMatrix, args...)
depwarn(string("`eig(A, args...)` has been deprecated in favor of ",
"`eigfact(A, args...)`. Whereas `eig(A, args....)` ",
"returns a tuple of arrays, `eigfact(A, args...)` returns ",
"an `Eigen` object. So for a direct replacement, use ",
"`(eigfact(A, args...)...,)`. But going forward, consider ",
"using the direct result of `eigfact(A, args...)` instead, ",
"either destructured into its components ",
"(`vals, vecs = eigfact(A, args...)`) ",
"or as an `Eigen` object (`eigf = eigfact(A, args...)`)."), :eig)
return (eigfact(A, args...)...,)
end
eig(A::AbstractMatrix, B::AbstractMatrix) = _geneig(A, B)
eig(A::Number, B::Number) = _geneig(A, B)
function _geneig(A, B)
depwarn(string("`eig(A::AbstractMatrix, B::AbstractMatrix)` and ",
"`eig(A::Number, B::Number)` have been deprecated in favor of ",
"`eigfact(A, B)`. Whereas the former each return a tuple of arrays, ",
"the latter returns a `GeneralizedEigen` object. So for a direct ",
"replacement, use `(eigfact(A, B)...,)`. But going forward, consider ",
"using the direct result of `eigfact(A, B)` instead, either ",
"destructured into its components (`vals, vecs = eigfact(A, B)`), ",
"or as a `GeneralizedEigen` object (`eigf = eigfact(A, B)`)."), :eig)
return (eigfact(A, B)...,)
end
97 changes: 33 additions & 64 deletions stdlib/LinearAlgebra/src/eigen.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,13 @@ end
GeneralizedEigen(values::AbstractVector{V}, vectors::AbstractMatrix{T}) where {T,V} =
GeneralizedEigen{T,V,typeof(vectors),typeof(values)}(values, vectors)

# iteration for destructuring into factors
Base.start(::Union{Eigen,GeneralizedEigen}) = Val(:values)
Base.next(F::Union{Eigen,GeneralizedEigen}, ::Val{:values}) = (F.values, Val(:vectors))
Base.next(F::Union{Eigen,GeneralizedEigen}, ::Val{:vectors}) = (F.vectors, Val(:done))
Base.done(F::Union{Eigen,GeneralizedEigen}, ::Val{:done}) = true
Base.done(F::Union{Eigen,GeneralizedEigen}, ::Any) = false

isposdef(A::Union{Eigen,GeneralizedEigen}) = isreal(A.values) && all(x -> x > 0, A.values)

"""
Expand Down Expand Up @@ -94,6 +101,20 @@ julia> F.values
18.0
julia> F.vectors
3×3 Array{Float64,2}:
1.0 0.0 0.0
0.0 1.0 0.0
0.0 0.0 1.0
julia> vals, vecs = F; # destructuring via iteration
julia> vals
3-element Array{Float64,1}:
1.0
3.0
18.0
julia> vecs
3×3 Array{Float64,2}:
1.0 0.0 0.0
0.0 1.0 0.0
Expand All @@ -107,38 +128,6 @@ function eigfact(A::StridedMatrix{T}; permute::Bool=true, scale::Bool=true) wher
end
eigfact(x::Number) = Eigen([x], fill(one(x), 1, 1))

function eig(A::Union{Number, StridedMatrix}; permute::Bool=true, scale::Bool=true)
F = eigfact(A, permute=permute, scale=scale)
F.values, F.vectors
end

"""
eig(A::Union{SymTridiagonal, Hermitian, Symmetric}, irange::UnitRange) -> D, V
eig(A::Union{SymTridiagonal, Hermitian, Symmetric}, vl::Real, vu::Real) -> D, V
eig(A, permute::Bool=true, scale::Bool=true) -> D, V
Computes eigenvalues (`D`) and eigenvectors (`V`) of `A`.
See [`eigfact`](@ref) for details on the
`irange`, `vl`, and `vu` arguments
(for [`SymTridiagonal`](@ref), [`Hermitian`](@ref), and
[`Symmetric`](@ref) matrices)
and the `permute` and `scale` keyword arguments.
The eigenvectors are returned columnwise.
# Examples
```jldoctest
julia> eig([1.0 0.0 0.0; 0.0 3.0 0.0; 0.0 0.0 18.0])
([1.0, 3.0, 18.0], [1.0 0.0 0.0; 0.0 1.0 0.0; 0.0 0.0 1.0])
```
`eig` is a wrapper around [`eigfact`](@ref), extracting all parts of the
factorization to a tuple; where possible, using [`eigfact`](@ref) is recommended.
"""
function eig(A::AbstractMatrix, args...)
F = eigfact(A, args...)
F.values, F.vectors
end

"""
eigvecs(A; permute::Bool=true, scale::Bool=true) -> Matrix
Expand Down Expand Up @@ -378,6 +367,18 @@ julia> F.values
0.0 - 1.0im
julia> F.vectors
2×2 Array{Complex{Float64},2}:
0.0-1.0im 0.0+1.0im
-1.0-0.0im -1.0+0.0im
julia> vals, vecs = F; # destructuring via iterationw
julia> vals
2-element Array{Complex{Float64},1}:
0.0 + 1.0im
0.0 - 1.0im
julia> vecs
2×2 Array{Complex{Float64},2}:
0.0-1.0im 0.0+1.0im
-1.0-0.0im -1.0+0.0im
Expand All @@ -390,38 +391,6 @@ end

eigfact(A::Number, B::Number) = eigfact(fill(A,1,1), fill(B,1,1))

"""
eig(A, B) -> D, V
Computes generalized eigenvalues (`D`) and vectors (`V`) of `A` with respect to `B`.
`eig` is a wrapper around [`eigfact`](@ref), extracting all parts of the
factorization to a tuple; where possible, using [`eigfact`](@ref) is recommended.
# Examples
```jldoctest
julia> A = [1 0; 0 -1]
2×2 Array{Int64,2}:
1 0
0 -1
julia> B = [0 1; 1 0]
2×2 Array{Int64,2}:
0 1
1 0
julia> eig(A, B)
(Complex{Float64}[0.0+1.0im, 0.0-1.0im], Complex{Float64}[0.0-1.0im 0.0+1.0im; -1.0-0.0im -1.0+0.0im])
```
"""
function eig(A::AbstractMatrix, B::AbstractMatrix)
F = eigfact(A,B)
F.values, F.vectors
end
function eig(A::Number, B::Number)
F = eigfact(A,B)
F.values, F.vectors
end

"""
eigvals!(A, B) -> values
Expand Down
2 changes: 1 addition & 1 deletion stdlib/LinearAlgebra/src/triangular.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2156,7 +2156,7 @@ function log(A0::UpperTriangular{T}) where T<:BlasFloat
R[i,i+1] = i / sqrt((2 * i)^2 - 1)
R[i+1,i] = R[i,i+1]
end
x,V = eig(R)
x,V = eigfact(R)
w = Vector{Float64}(undef, m)
for i = 1:m
x[i] = (x[i] + 1) / 2
Expand Down
4 changes: 2 additions & 2 deletions stdlib/LinearAlgebra/test/bidiag.jl
Original file line number Diff line number Diff line change
Expand Up @@ -231,8 +231,8 @@ srand(1)

@testset "Eigensystems" begin
if relty <: AbstractFloat
d1, v1 = eig(T)
d2, v2 = eig(map(elty<:Complex ? ComplexF64 : Float64,Tfull))
d1, v1 = eigfact(T)
d2, v2 = eigfact(map(elty<:Complex ? ComplexF64 : Float64,Tfull))
@test (uplo == :U ? d1 : reverse(d1)) d2
if elty <: Real
test_approx_eq_modphase(v1, uplo == :U ? v2 : v2[:,n:-1:1])
Expand Down
14 changes: 7 additions & 7 deletions stdlib/LinearAlgebra/test/eigen.jl
Original file line number Diff line number Diff line change
Expand Up @@ -28,12 +28,12 @@ aimg = randn(n,n)/2

α = rand(eltya)
β = rand(eltya)
eab = eig(α,β)
eab = (eigfact(α,β)...,)
@test eab[1] == eigvals(fill(α,1,1),fill(β,1,1))
@test eab[2] == eigvecs(fill(α,1,1),fill(β,1,1))

@testset "non-symmetric eigen decomposition" begin
d, v = eig(a)
d, v = eigfact(a)
for i in 1:size(a,2)
@test a*v[:,i] d[i]*v[:,i]
end
Expand Down Expand Up @@ -70,7 +70,7 @@ aimg = randn(n,n)/2
@test eigvecs(f) === f.vectors
@test_throws ErrorException f.Z

d,v = eig(asym_sg, a_sg'a_sg)
d,v = eigfact(asym_sg, a_sg'a_sg)
@test d == f.values
@test v == f.vectors
end
Expand All @@ -89,7 +89,7 @@ aimg = randn(n,n)/2
@test eigvecs(a1_nsg, a2_nsg) == f.vectors
@test_throws ErrorException f.Z

d,v = eig(a1_nsg, a2_nsg)
d,v = eigfact(a1_nsg, a2_nsg)
@test d == f.values
@test v == f.vectors
end
Expand All @@ -98,11 +98,11 @@ end

@testset "eigenvalue computations with NaNs" begin
for eltya in (NaN16, NaN32, NaN)
@test_throws(ArgumentError, eig(fill(eltya, 1, 1)))
@test_throws(ArgumentError, eig(fill(eltya, 2, 2)))
@test_throws(ArgumentError, eigfact(fill(eltya, 1, 1)))
@test_throws(ArgumentError, eigfact(fill(eltya, 2, 2)))
test_matrix = rand(typeof(eltya),3,3)
test_matrix[2,2] = eltya
@test_throws(ArgumentError, eig(test_matrix))
@test_throws(ArgumentError, eigfact(test_matrix))
end
end

Expand Down
2 changes: 1 addition & 1 deletion stdlib/LinearAlgebra/test/lu.jl
Original file line number Diff line number Diff line change
Expand Up @@ -51,7 +51,7 @@ dimg = randn(n)/2
-eps(real(one(eltya)))/4 eps(real(one(eltya)))/2 -1.0 0;
-0.5 -0.5 0.1 1.0])
F = eigfact(A, permute=false, scale=false)
eig(A, permute=false, scale=false)
(eigfact(A, permute=false, scale=false)...,)
@test F.vectors*Diagonal(F.values)/F.vectors A
F = eigfact(A)
# @test norm(F.vectors*Diagonal(F.values)/F.vectors - A) > 0.01
Expand Down
2 changes: 1 addition & 1 deletion stdlib/LinearAlgebra/test/schur.jl
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@ aimg = randn(n,n)/2
view(apd, 1:n, 1:n)))
ε = εa = eps(abs(float(one(eltya))))

d,v = eig(a)
d,v = eigfact(a)
f = schurfact(a)
@test f.vectors*f.Schur*f.vectors' a
@test sort(real(f.values)) sort(real(d))
Expand Down
14 changes: 7 additions & 7 deletions stdlib/LinearAlgebra/test/symmetric.jl
Original file line number Diff line number Diff line change
Expand Up @@ -218,29 +218,29 @@ end

@testset "symmetric eigendecomposition" begin
if eltya <: Real # the eigenvalues are only real and ordered for Hermitian matrices
d, v = eig(asym)
d, v = eigfact(asym)
@test asym*v[:,1] d[1]*v[:,1]
@test v*Diagonal(d)*transpose(v) asym
@test isequal(eigvals(asym[1]), eigvals(asym[1:1,1:1]))
@test abs.(eigfact(Symmetric(asym), 1:2).vectors'v[:,1:2]) Matrix(I, 2, 2)
eig(Symmetric(asym), 1:2) # same result, but checks that method works
(eigfact(Symmetric(asym), 1:2)...,) # same result, but checks that method works
@test abs.(eigfact(Symmetric(asym), d[1] - 1, (d[2] + d[3])/2).vectors'v[:,1:2]) Matrix(I, 2, 2)
eig(Symmetric(asym), d[1] - 1, (d[2] + d[3])/2) # same result, but checks that method works
(eigfact(Symmetric(asym), d[1] - 1, (d[2] + d[3])/2)...,) # same result, but checks that method works
@test eigvals(Symmetric(asym), 1:2) d[1:2]
@test eigvals(Symmetric(asym), d[1] - 1, (d[2] + d[3])/2) d[1:2]
# eigfact doesn't support Symmetric{Complex}
@test Matrix(eigfact(asym)) asym
@test eigvecs(Symmetric(asym)) eigvecs(asym)
end

d, v = eig(aherm)
d, v = eigfact(aherm)
@test aherm*v[:,1] d[1]*v[:,1]
@test v*Diagonal(d)*v' aherm
@test isequal(eigvals(aherm[1]), eigvals(aherm[1:1,1:1]))
@test abs.(eigfact(Hermitian(aherm), 1:2).vectors'v[:,1:2]) Matrix(I, 2, 2)
eig(Hermitian(aherm), 1:2) # same result, but checks that method works
(eigfact(Hermitian(aherm), 1:2)...,) # same result, but checks that method works
@test abs.(eigfact(Hermitian(aherm), d[1] - 1, (d[2] + d[3])/2).vectors'v[:,1:2]) Matrix(I, 2, 2)
eig(Hermitian(aherm), d[1] - 1, (d[2] + d[3])/2) # same result, but checks that method works
(eigfact(Hermitian(aherm), d[1] - 1, (d[2] + d[3])/2)...,) # same result, but checks that method works
@test eigvals(Hermitian(aherm), 1:2) d[1:2]
@test eigvals(Hermitian(aherm), d[1] - 1, (d[2] + d[3])/2) d[1:2]
@test Matrix(eigfact(aherm)) aherm
Expand Down Expand Up @@ -365,7 +365,7 @@ end
end

@testset "Issues #8057 and #8058. f=$f, A=$A" for f in
(eigfact, eigvals, eig),
(eigfact, eigvals, eigfact),
A in (Symmetric([0 1; 1 0]), Hermitian([0 im; -im 0]))
@test_throws ArgumentError f(A, 3, 2)
@test_throws ArgumentError f(A, 1:4)
Expand Down
2 changes: 1 addition & 1 deletion stdlib/LinearAlgebra/test/triangular.jl
Original file line number Diff line number Diff line change
Expand Up @@ -249,7 +249,7 @@ for elty1 in (Float32, Float64, BigFloat, ComplexF32, ComplexF64, Complex{BigFlo

# eigenproblems
if !(elty1 in (BigFloat, Complex{BigFloat})) # Not handled yet
vals, vecs = eig(A1)
vals, vecs = eigfact(A1)
if (t1 == UpperTriangular || t1 == LowerTriangular) && elty1 != Int # Cannot really handle degenerate eigen space and Int matrices will probably have repeated eigenvalues.
@test vecs*diagm(0 => vals)/vecs A1 atol=sqrt(eps(float(real(one(vals[1])))))*(norm(A1,Inf)*n)^2
end
Expand Down
10 changes: 5 additions & 5 deletions stdlib/LinearAlgebra/test/tridiag.jl
Original file line number Diff line number Diff line change
Expand Up @@ -243,7 +243,7 @@ end
w, iblock, isplit = LAPACK.stebz!('V', 'B', -infinity, infinity, 0, 0, zero, b, a)
evecs = LAPACK.stein!(b, a, w)

(e, v) = eig(SymTridiagonal(b, a))
(e, v) = eigfact(SymTridiagonal(b, a))
@test e w
test_approx_eq_vecs(v, evecs)
end
Expand All @@ -266,10 +266,10 @@ end
end
@testset "eigenvalues/eigenvectors of symmetric tridiagonal" begin
if elty === Float32 || elty === Float64
DT, VT = @inferred eig(A)
@inferred eig(A, 2:4)
@inferred eig(A, 1.0, 2.0)
D, Vecs = eig(fA)
DT, VT = @inferred eigfact(A)
@inferred eigfact(A, 2:4)
@inferred eigfact(A, 1.0, 2.0)
D, Vecs = eigfact(fA)
@test DT D
@test abs.(VT'Vecs) Matrix(elty(1)I, n, n)
test_approx_eq_modphase(eigvecs(A), eigvecs(fA))
Expand Down
1 change: 0 additions & 1 deletion stdlib/SparseArrays/src/linalg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1009,7 +1009,6 @@ function factorize(A::LinearAlgebra.RealHermSymComplexHerm{Float64,<:SparseMatri
end

chol(A::SparseMatrixCSC) = error("Use cholfact() instead of chol() for sparse matrices.")
eig(A::SparseMatrixCSC) = error("Use IterativeEigensolvers.eigs() instead of eig() for sparse matrices.")

function Base.cov(X::SparseMatrixCSC; dims::Int=1, corrected::Bool=true)
vardim = dims
Expand Down
1 change: 0 additions & 1 deletion stdlib/SparseArrays/test/sparse.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1780,7 +1780,6 @@ end
C, b = A[:, 1:4], fill(1., size(A, 1))
@test !Base.USE_GPL_LIBS || factorize(C)\b Array(C)\b
@test_throws ErrorException chol(A)
@test_throws ErrorException eig(A)
@test_throws ErrorException inv(A)
end

Expand Down

0 comments on commit 78344d1

Please sign in to comment.