diff --git a/base/deprecated.jl b/base/deprecated.jl index 1e4853fd6e05a..669df26b57e81 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -856,7 +856,7 @@ end for f in (:sin, :sinh, :sind, :asin, :asinh, :asind, :tan, :tanh, :tand, :atan, :atanh, :atand, :sinpi, :cosc, :ceil, :floor, :trunc, :round, :real, :imag, - :log1p, :expm1, :abs, :abs2, :conj, + :log1p, :expm1, :abs, :abs2, :log, :log2, :log10, :exp, :exp2, :exp10, :sinc, :cospi, :cos, :cosh, :cosd, :acos, :acosd, :cot, :coth, :cotd, :acot, :acotd, diff --git a/base/linalg/factorization.jl b/base/linalg/factorization.jl index fcc53b2108b28..045d88212a755 100644 --- a/base/linalg/factorization.jl +++ b/base/linalg/factorization.jl @@ -39,8 +39,7 @@ function (\){T<:BlasReal}(F::Factorization{T}, B::VecOrMat{Complex{T}}) end for (f1, f2) in ((:\, :A_ldiv_B!), - (:Ac_ldiv_B, :Ac_ldiv_B!), - (:At_ldiv_B, :At_ldiv_B!)) + (:Ac_ldiv_B, :Ac_ldiv_B!)) @eval begin function $f1(F::Factorization, B::AbstractVecOrMat) TFB = typeof(one(eltype(F)) / one(eltype(B))) @@ -57,6 +56,10 @@ for f in (:A_ldiv_B!, :Ac_ldiv_B!, :At_ldiv_B!) $f(A, copy!(Y, B)) end +# fallback methods for transposed solves +At_ldiv_B{T<:Real}(F::Factorization{T}, B::AbstractVecOrMat) = Ac_ldiv_B(F, B) +At_ldiv_B(F::Factorization, B) = conj.(Ac_ldiv_B(F, conj.(B))) + """ A_ldiv_B!([Y,] A, B) -> Y diff --git a/base/sparse/cholmod.jl b/base/sparse/cholmod.jl index 1f127579aa3c6..c35a2653a1494 100644 --- a/base/sparse/cholmod.jl +++ b/base/sparse/cholmod.jl @@ -1504,6 +1504,19 @@ Ac_ldiv_B(L::Factor, B::VecOrMat) = convert(Matrix, solve(CHOLMOD_A, L, Dense(B) Ac_ldiv_B(L::Factor, B::Sparse) = spsolve(CHOLMOD_A, L, B) Ac_ldiv_B(L::Factor, B::SparseVecOrMat) = Ac_ldiv_B(L, Sparse(B)) +for f in (:\, :Ac_ldiv_B) + @eval function ($f)(A::Union{Symmetric{Float64,SparseMatrixCSC{Float64,SuiteSparse_long}}, + Hermitian{Float64,SparseMatrixCSC{Float64,SuiteSparse_long}}, + Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},SuiteSparse_long}}}, B::StridedVecOrMat) + try + return ($f)(cholfact(A), B) + catch e + isa(e, LinAlg.PosDefException) || rethrow(e) + return ($f)(ldltfact(A) , B) + end + end +end + ## Other convenience methods function diag{Tv}(F::Factor{Tv}) f = unsafe_load(get(F.p)) diff --git a/base/sparse/linalg.jl b/base/sparse/linalg.jl index 8dd9b869621dd..18570640c655d 100644 --- a/base/sparse/linalg.jl +++ b/base/sparse/linalg.jl @@ -864,12 +864,7 @@ function (\)(A::SparseMatrixCSC, B::AbstractVecOrMat) return UpperTriangular(A) \ B end if ishermitian(A) - try - return cholfact(Hermitian(A)) \ B - catch e - isa(e, PosDefException) || rethrow(e) - return ldltfact(Hermitian(A)) \ B - end + Hermitian(A) \ B end return lufact(A) \ B else diff --git a/test/linalg/lu.jl b/test/linalg/lu.jl index 1c0ea293fbf41..6e5bbbffa7ab1 100644 --- a/test/linalg/lu.jl +++ b/test/linalg/lu.jl @@ -44,7 +44,7 @@ for eltya in (Float32, Float64, Complex64, Complex128, BigFloat, Int) εb = eps(abs(float(one(eltyb)))) ε = max(εa,εb) -debug && println("(Automatic) Square LU decomposition") +debug && println("(Automatic) Square LU decomposition. eltya: $eltya, eltyb: $eltyb") κ = cond(a,1) lua = factorize(a) @test_throws KeyError lua[:Z] diff --git a/test/sparse/cholmod.jl b/test/sparse/cholmod.jl index d86f2937a5a2e..64527dbfd2241 100644 --- a/test/sparse/cholmod.jl +++ b/test/sparse/cholmod.jl @@ -381,6 +381,7 @@ for elty in (Float64, Complex{Float64}) @test_throws DimensionMismatch F\CHOLMOD.Sparse(sparse(ones(elty, 4))) @test F'\ones(elty, 5) ≈ Array(A1pd)'\ones(5) @test F'\sparse(ones(elty, 5)) ≈ Array(A1pd)'\ones(5) + @test F.'\ones(elty, 5) ≈ conj(A1pd)'\ones(elty, 5) @test logdet(F) ≈ logdet(Array(A1pd)) @test det(F) == exp(logdet(F)) let # to test supernodal, we must use a larger matrix @@ -670,3 +671,19 @@ let m = 400, n = 500 @test s.is_super == 0 @test F\b ≈ ones(m + n) end + +# Test that \ and '\ and .'\ work for Symmetric and Hermitian. This is just +# a dispatch exercise so it doesn't matter that the complex matrix has +# zero imaginary parts +let Apre = sprandn(10, 10, 0.2) - I + for A in (Symmetric(Apre), Hermitian(Apre), + Symmetric(Apre + 10I), Hermitian(Apre + 10I), + Hermitian(complex(Apre)), Hermitian(complex(Apre) + 10I)) + + x = ones(10) + b = A*x + @test x ≈ A\b + @test A.'\b ≈ A'\b + end +end +