diff --git a/base/deprecated.jl b/base/deprecated.jl index 58046dafaa64c5..828b3ce5c75388 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -1491,6 +1491,18 @@ export conv, conv2, deconv, filt, filt!, xcorr @deprecate cov(X::AbstractVector, Y::AbstractVector, corrected::Bool) cov(X, Y, corrected=corrected) @deprecate cov(X::AbstractVecOrMat, Y::AbstractVecOrMat, vardim::Int, corrected::Bool) cov(X, Y, vardim, corrected=corrected) +# bkfact +function bkfact(A::StridedMatrix, uplo::Symbol, symmetric::Bool = issymmetric(A), rook::Bool = false) + depwarn("bkfact with uplo and symmetric arguments deprecated. Please use bkfact($(symmetric ? "Symmetric(" : "Hermitian(")A, :$uplo))", + :bkfact) + return bkfact(symmetric ? Symmetric(A, uplo) : Hermitian(A, uplo), rook) +end +function bkfact!(A::StridedMatrix, uplo::Symbol, symmetric::Bool = issymmetric(A), rook::Bool = false) + depwarn("bkfact! with uplo and symmetric arguments deprecated. Please use bkfact!($(symmetric ? "Symmetric(" : "Hermitian(")A, :$uplo))", + :bkfact!) + return bkfact!(symmetric ? Symmetric(A, uplo) : Hermitian(A, uplo), rook) +end + # END 0.7 deprecations # BEGIN 1.0 deprecations diff --git a/base/linalg/bunchkaufman.jl b/base/linalg/bunchkaufman.jl index 3b674df501d44a..248be82756835b 100644 --- a/base/linalg/bunchkaufman.jl +++ b/base/linalg/bunchkaufman.jl @@ -22,36 +22,38 @@ BunchKaufman(A::AbstractMatrix{T}, ipiv::Vector{BlasInt}, uplo::Char, symmetric: `bkfact!` is the same as [`bkfact`](@ref), but saves space by overwriting the input `A`, instead of creating a copy. """ -function bkfact!(A::StridedMatrix{<:BlasReal}, uplo::Symbol = :U, - symmetric::Bool = issymmetric(A), rook::Bool = false) - - if !symmetric - throw(ArgumentError("Bunch-Kaufman decomposition is only valid for symmetric matrices")) +function bkfact!(A::HermOrSym{T,S} where {T<:BlasReal,S<:StridedMatrix{T}}, rook::Bool = false) + if rook + LD, ipiv, info = LAPACK.sytrf_rook!(A.uplo, A.data) + else + LD, ipiv, info = LAPACK.sytrf!(A.uplo, A.data) end + BunchKaufman(LD, ipiv, A.uplo, true, rook, info) +end +function bkfact!(A::Symmetric{T,S} where {T<:BlasComplex,S<:StridedMatrix{T}}, rook::Bool = false) if rook - LD, ipiv, info = LAPACK.sytrf_rook!(char_uplo(uplo), A) + LD, ipiv, info = LAPACK.sytrf_rook!(A.uplo, A.data) else - LD, ipiv, info = LAPACK.sytrf!(char_uplo(uplo), A) + LD, ipiv, info = LAPACK.sytrf!(A.uplo, A.data) end - BunchKaufman(LD, ipiv, char_uplo(uplo), symmetric, rook, info) + BunchKaufman(LD, ipiv, A.uplo, true, rook, info) end -function bkfact!(A::StridedMatrix{<:BlasComplex}, uplo::Symbol=:U, - symmetric::Bool=issymmetric(A), rook::Bool=false) - +function bkfact!(A::Hermitian{T,S} where {T<:BlasComplex,S<:StridedMatrix{T}}, rook::Bool = false) if rook - if symmetric - LD, ipiv, info = LAPACK.sytrf_rook!(char_uplo(uplo), A) - else - LD, ipiv, info = LAPACK.hetrf_rook!(char_uplo(uplo), A) - end + LD, ipiv, info = LAPACK.hetrf_rook!(A.uplo, A.data) else - if symmetric - LD, ipiv, info = LAPACK.sytrf!(char_uplo(uplo), A) - else - LD, ipiv, info = LAPACK.hetrf!(char_uplo(uplo), A) - end + LD, ipiv, info = LAPACK.hetrf!(A.uplo, A.data) + end + BunchKaufman(LD, ipiv, A.uplo, false, rook, info) +end +function bkfact!(A::StridedMatrix{<:BlasFloat}, rook::Bool = false) + if ishermitian(A) + return bkfact!(Hermitian(A), rook) + elseif issymmetric(A) + return bkfact!(Symmetric(A), rook) + else + throw(ArgumentError("Bunch-Kaufman decomposition is only valid for symmetric or Hermitian matrices")) end - BunchKaufman(LD, ipiv, char_uplo(uplo), symmetric, rook, info) end """ @@ -69,13 +71,8 @@ The following functions are available for [^Bunch1977]: J R Bunch and L Kaufman, Some stable methods for calculating inertia and solving symmetric linear systems, Mathematics of Computation 31:137 (1977), 163-179. [url](http://www.ams.org/journals/mcom/1977-31-137/S0025-5718-1977-0428694-0/). """ -bkfact(A::StridedMatrix{<:BlasFloat}, uplo::Symbol=:U, symmetric::Bool=issymmetric(A), - rook::Bool=false) = - bkfact!(copy(A), uplo, symmetric, rook) -bkfact(A::StridedMatrix{T}, uplo::Symbol=:U, symmetric::Bool=issymmetric(A), - rook::Bool=false) where {T} = - bkfact!(convert(Matrix{promote_type(Float32, typeof(sqrt(one(T))))}, A), - uplo, symmetric, rook) +bkfact(A::AbstractMatrix{T}, rook::Bool=false) where {T} = + bkfact!(copy_oftype(A, typeof(sqrt(one(T)))), rook) convert(::Type{BunchKaufman{T}}, B::BunchKaufman{T}) where {T} = B convert(::Type{BunchKaufman{T}}, B::BunchKaufman) where {T} = @@ -90,9 +87,9 @@ ishermitian(B::BunchKaufman) = !B.symmetric issuccess(B::BunchKaufman) = B.info == 0 -function Base.show(io::IO, F::BunchKaufman) - println(io, summary(F)) - print(io, "successful: $(issuccess(F))") +function Base.show(io::IO, mime::MIME{Symbol("text/plain")}, B::BunchKaufman) + println(io, summary(B)) + print(io, "successful: $(issuccess(B))") end function inv(B::BunchKaufman{<:BlasReal}) diff --git a/base/linalg/symmetric.jl b/base/linalg/symmetric.jl index 7a669199ab61c8..13101861467e04 100644 --- a/base/linalg/symmetric.jl +++ b/base/linalg/symmetric.jl @@ -330,7 +330,6 @@ for T in (:Symmetric, :Hermitian), op in (:+, :-, :*, :/) @eval ($op)(A::$T, x::$S) = ($T)(($op)(A.data, x), Symbol(A.uplo)) end -bkfact(A::HermOrSym) = bkfact(A.data, Symbol(A.uplo), issymmetric(A)) factorize(A::HermOrSym) = bkfact(A) det(A::RealHermSymComplexHerm) = real(det(bkfact(A))) diff --git a/test/linalg/bunchkaufman.jl b/test/linalg/bunchkaufman.jl index 42289b4a7732f5..46e77b3886c4b9 100644 --- a/test/linalg/bunchkaufman.jl +++ b/test/linalg/bunchkaufman.jl @@ -23,16 +23,17 @@ bimg = randn(n,2)/2 a = eltya == Int ? rand(1:7, n, n) : convert(Matrix{eltya}, eltya <: Complex ? complex.(areal, aimg) : areal) a2 = eltya == Int ? rand(1:7, n, n) : convert(Matrix{eltya}, eltya <: Complex ? complex.(a2real, a2img) : a2real) @testset for atype in ("Array", "SubArray") - asym = a'+a # symmetric indefinite - apd = a'*a # symmetric positive-definite + asym = a.'+ a # symmetric indefinite + aher = a' + a # Hermitian indefinite + apd = a' * a # Positive-definite if atype == "Array" - a = a + a = a a2 = a2 else - a = view(a, 1:n, 1:n) - a2 = view(a2, 1:n, 1:n) - asym = view(asym, 1:n, 1:n) - apd = view(apd, 1:n, 1:n) + a = view(a , 1:n, 1:n) + a2 = view(a2 , 1:n, 1:n) + aher = view(aher, 1:n, 1:n) + apd = view(apd , 1:n, 1:n) end ε = εa = eps(abs(float(one(eltya)))) @@ -48,8 +49,12 @@ bimg = randn(n,2)/2 εb = eps(abs(float(one(eltyb)))) ε = max(εa,εb) - @testset "(Automatic) Bunch-Kaufman factor of indefinite matrix" begin - bc1 = factorize(asym) + # check that factorize gives a Bunch-Kaufman + @test isa(factorize(asym), LinAlg.BunchKaufman) + @test isa(factorize(aher), LinAlg.BunchKaufman) + + @testset "$uplo Bunch-Kaufman factor of indefinite matrix" for uplo in (:L, :U) + bc1 = bkfact(Hermitian(aher, uplo)) @test LinAlg.issuccess(bc1) @test logabsdet(bc1)[1] ≈ log(abs(det(bc1))) if eltya <: Real @@ -57,27 +62,27 @@ bimg = randn(n,2)/2 else @test logabsdet(bc1)[2] ≈ sign(det(bc1)) end - @test inv(bc1)*asym ≈ eye(n) - @test asym*(bc1\b) ≈ b atol=1000ε + @test inv(bc1)*aher ≈ eye(n) + @test aher*(bc1\b) ≈ b atol=1000ε @testset for rook in (false, true) - @test inv(bkfact(a.'+a, :U, true, rook))*(a.'+a) ≈ eye(n) + @test inv(bkfact(Symmetric(a.' + a, uplo), rook))*(a.' + a) ≈ eye(n) @test size(bc1) == size(bc1.LD) - @test size(bc1,1) == size(bc1.LD,1) - @test size(bc1,2) == size(bc1.LD,2) + @test size(bc1, 1) == size(bc1.LD, 1) + @test size(bc1, 2) == size(bc1.LD, 2) if eltya <: BlasReal @test_throws ArgumentError bkfact(a) end end end - @testset "Bunch-Kaufman factors of a pos-def matrix" begin - @testset for rook in (false, true) - bc2 = bkfact(apd, :U, issymmetric(apd), rook) + @testset "$uplo Bunch-Kaufman factors of a pos-def matrix" for uplo in (:U, :L) + @testset "rook pivoting: $rook" for rook in (false, true) + bc2 = bkfact(Hermitian(apd, uplo), rook) @test LinAlg.issuccess(bc2) @test logdet(bc2) ≈ log(det(bc2)) @test logabsdet(bc2)[1] ≈ log(abs(det(bc2))) @test logabsdet(bc2)[2] == sign(det(bc2)) - @test inv(bc2)*apd ≈ eye(n) + @test inv(bc2)*apd ≈ eye(eltyb, n) @test apd*(bc2\b) ≈ b atol=150000ε @test ishermitian(bc2) == !issymmetric(bc2) end @@ -104,11 +109,13 @@ end end @testset for rook in (false, true) - F = bkfact(As, :U, issymmetric(As), rook) - @test !LinAlg.issuccess(F) - @test det(F) == 0 - @test_throws LinAlg.SingularException inv(F) - @test_throws LinAlg.SingularException F \ ones(size(As, 1)) + @testset for uplo in (:L, :U) + F = bkfact(issymmetric(As) ? Symmetric(As, uplo) : Hermitian(As, uplo), rook) + @test !LinAlg.issuccess(F) + @test det(F) == 0 + @test_throws LinAlg.SingularException inv(F) + @test_throws LinAlg.SingularException F \ ones(size(As, 1)) + end end end end