diff --git a/base/deprecated.jl b/base/deprecated.jl index 50f8b30b86268..5d03feff16a09 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -146,12 +146,6 @@ end @deprecate expr(hd, a...) Expr(hd, a...) @deprecate expr(hd, a::Array{Any,1}) Expr(hd, a...) -@deprecate cholfact chol -@deprecate cholpfact cholp -@deprecate lufact lu -@deprecate qrfact qr -@deprecate qrpfact qrp - @deprecate logb exponent @deprecate ref_shape index_shape @deprecate assign_shape_check setindex_shape_check diff --git a/base/exports.jl b/base/exports.jl index 163102b592357..622920171fbac 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -550,7 +550,11 @@ export # linear algebra chol, + cholfact, + cholfact!, cholp, + cholpfact, + cholpfact!, cond, cross, ctranspose, diff --git a/base/linalg/dense.jl b/base/linalg/dense.jl index a5ab2bb56a45b..6cfc015a0e19b 100644 --- a/base/linalg/dense.jl +++ b/base/linalg/dense.jl @@ -415,7 +415,7 @@ function det(A::Matrix) end det(x::Number) = x -logdet(A::Matrix) = 2.0 * sum(log(diag(chol(A)[:U]))) +logdet(A::Matrix) = 2.0 * sum(log(diag(cholfact(A)[:U]))) function inv(A::StridedMatrix) if istriu(A) return inv(Triangular(A, 'U')) end diff --git a/base/linalg/factorization.jl b/base/linalg/factorization.jl index bcf7217e5851f..b462991ca4c41 100644 --- a/base/linalg/factorization.jl +++ b/base/linalg/factorization.jl @@ -13,10 +13,14 @@ type CholeskyDense{T<:BlasFloat} <: Factorization{T} end CholeskyDense{T<:BlasFloat}(A::Matrix{T}, uplo::Char) = CholeskyDense{T}(A, uplo) -chol(A::Matrix, uplo::Symbol) = CholeskyDense(copy(A), string(uplo)[1]) -chol(A::Matrix) = chol(A, :U) -chol{T<:Integer}(A::Matrix{T}, args...) = chol(float64(A), args...) -chol(x::Number) = imag(x) == 0 && real(x) > 0 ? sqrt(x) : error("Argument not positive-definite") +cholfact!(A::Matrix, uplo::Symbol) = CholeskyDense(A, string(uplo)[1]) +cholfact(A::Matrix, uplo::Symbol) = cholfact!(copy(A), uplo) +cholfact!(A::Matrix) = cholfact!(A, :U) +cholfact(A::Matrix) = cholfact(A, :U) +cholfact{T<:Integer}(A::Matrix{T}, args...) = cholfact(float(A), args...) +cholfact(x::Number) = imag(x) == 0 && real(x) > 0 ? sqrt(x) : error("Argument not positive-definite") + +chol(A) = cholfact(A, :U)[:U] size(C::CholeskyDense) = size(C.UL) size(C::CholeskyDense,d::Integer) = size(C.UL,d) @@ -59,10 +63,13 @@ function CholeskyPivotedDense{T<:BlasFloat}(A::Matrix{T}, uplo::Char, tol::Real) CholeskyPivotedDense{T}(uplo == 'U' ? triu!(A) : tril!(A), uplo, piv, rank, tol, info) end -cholp(A::Matrix, uplo::Symbol, tol::Real) = CholeskyPivotedDense(copy(A), string(uplo)[1], tol) -cholp(A::Matrix, tol::Real) = cholp(A, :U, tol) -cholp(A::Matrix) = cholp(A, -1.) -cholp{T<:Int}(A::Matrix{T}, args...) = cholp(float64(A), args...) +cholpfact!(A::Matrix, uplo::Symbol, tol::Real) = CholeskyPivotedDense(A, string(uplo)[1], tol) +cholpfact(A::Matrix, uplo::Symbol, tol::Real) = cholpfact!(copy(A), uplo, tol) +cholpfact!(A::Matrix, tol::Real) = cholpfact!(A, :U, tol) +cholpfact(A::Matrix, tol::Real) = cholpfact(A, :U, tol) +cholpfact!(A::Matrix) = cholpfact!(A, -1.) +cholpfact(A::Matrix) = cholpfact(A, -1.) +cholpfact{T<:Int}(A::Matrix{T}, args...) = cholpfact(float(A), args...) size(C::CholeskyPivotedDense) = size(C.UL) size(C::CholeskyPivotedDense,d::Integer) = size(C.UL,d) @@ -127,9 +134,16 @@ function LUDense{T<:BlasFloat}(A::Matrix{T}) LUDense{T}(LU, ipiv, info) end -lu(A::Matrix) = LUDense(copy(A)) -lu{T<:Integer}(A::Matrix{T}) = lu(float(A)) -lu(x::Number) = (one(x), x, [1]) +lufact!(A::Matrix) = LUDense(A) +lufact(A::Matrix) = lufact!(copy(A)) +lufact!{T<:Integer}(A::Matrix{T}) = lufact!(float(A)) +lufact{T<:Integer}(A::Matrix{T}) = lufact(float(A)) +lufact(x::Number) = (one(x), x, [1]) + +function lu(A::Matrix) + F = lufact(A) + return (F[:L], F[:U], F[:P]) +end size(A::LUDense) = size(A.LU) size(A::LUDense,n) = size(A.LU,n) @@ -182,9 +196,15 @@ type QRDense{S} <: Factorization{S} end QRDense(A::Matrix) = QRDense(LAPACK.geqrt3!(A)...) -qr(A::Matrix) = QRDense(copy(A)) -qr{T<:Integer}(A::Matrix{T}) = qr(float(A)) -qr(x::Number) = (one(x), x) +qrfact!(A::Matrix) = QRDense(A) +qrfact(A::Matrix) = qrfact!(copy(A)) +qrfact{T<:Integer}(A::Matrix{T}) = qrfact(float(A)) +qrfact(x::Number) = (one(x), x) + +function qr(A::Matrix) + F = qrfact(A) + return (F[:Q], F[:R]) +end size(A::QRDense, args::Integer...) = size(A.vs, args...) @@ -254,7 +274,9 @@ type QRPivotedDense{T} <: Factorization{T} end end QRPivotedDense{T<:BlasFloat}(A::Matrix{T}) = QRPivotedDense{T}(LAPACK.geqp3!(A)...) -qrp(A::Matrix) = QRPivotedDense(copy(A)) + +qrpfact!(A::Matrix) = QRPivotedDense(A) +qrpfact(A::Matrix) = qrpfact!(copy(A)) # QRDenseQ(A::QRPivotedDense) = QRDenseQ(A.hh, A.tau) size(A::QRPivotedDense, args::Integer...) = size(A.hh, args...) diff --git a/base/linalg/linalg.jl b/base/linalg/linalg.jl index b8da036136514..335ee84ff81d4 100644 --- a/base/linalg/linalg.jl +++ b/base/linalg/linalg.jl @@ -26,7 +26,11 @@ export # Functions chol, + cholfact, + cholfact!, cholp, + cholpfact, + cholpfact!, cond, copy!, cross, diff --git a/test/linalg.jl b/test/linalg.jl index 985018d61a9e4..be6fbefd18d1b 100644 --- a/test/linalg.jl +++ b/test/linalg.jl @@ -8,7 +8,7 @@ for elty in (Float32, Float64, Complex64, Complex128) apd = a'*a # symmetric positive-definite b = convert(Vector{elty}, b) - capd = chol(apd) # upper Cholesky factor + capd = cholfact(apd) # upper Cholesky factor r = capd[:U] @test_approx_eq r'*r apd @test_approx_eq b apd * (capd\b) @@ -16,10 +16,10 @@ for elty in (Float32, Float64, Complex64, Complex128) @test_approx_eq a*(capd\(a'*b)) b # least squares soln for square a @test_approx_eq det(capd) det(apd) - l = chol(apd, :L)[:L] # lower Cholesky factor + l = cholfact(apd, :L)[:L] # lower Cholesky factor @test_approx_eq l*l' apd - cpapd = cholp(apd) # pivoted Choleksy decomposition + cpapd = cholpfact(apd) # pivoted Choleksy decomposition @test rank(cpapd) == n @test all(diff(diag(real(cpapd.UL))).<=0.) # diagonal should be non-increasing @test_approx_eq b apd * (cpapd\b) @@ -32,21 +32,21 @@ for elty in (Float32, Float64, Complex64, Complex128) @test_approx_eq inv(bc2) * apd eye(elty, n) @test_approx_eq apd * (bc2\b) b - lua = lu(a) # LU decomposition + lua = lufact(a) # LU decomposition l,u,p = lua[:L], lua[:U], lua[:p] @test_approx_eq l*u a[p,:] @test_approx_eq l[invperm(p),:]*u a @test_approx_eq a * inv(lua) eye(elty, n) @test_approx_eq a*(lua\b) b - qra = qr(a) # QR decomposition + qra = qrfact(a) # QR decomposition q,r = qra[:Q], qra[:R] @test_approx_eq q'*full(q, false) eye(elty, n) @test_approx_eq q*full(q, false)' eye(elty, n) @test_approx_eq q*r a @test_approx_eq a*(qra\b) b - qrpa = qrp(a) # pivoted QR decomposition + qrpa = qrpfact(a) # pivoted QR decomposition q,r,p = qrpa[:Q], qrpa[:R], qrpa[:p] @test_approx_eq q'*full(q, false) eye(elty, n) @test_approx_eq q*full(q, false)' eye(elty, n)