From 501e6dd9e9f80afe88fc93898b015a53fc8768d6 Mon Sep 17 00:00:00 2001 From: "Steven G. Johnson" Date: Wed, 15 Apr 2015 21:13:52 -0400 Subject: [PATCH] =?UTF-8?q?support=20user-specified=20permutation=20in=20C?= =?UTF-8?q?HOLMOD=20factorizations,=20change=20=CE=B2=20optional=20argumen?= =?UTF-8?q?t=20to=20a=20shift=20keyword=20arg,=20some=20code=20consolidati?= =?UTF-8?q?on;=20change=20&=20to=20new=20Ref=20syntax=20in=20cholmod?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- NEWS.md | 6 ++ base/deprecated.jl | 3 + base/sparse/cholmod.jl | 144 +++++++++++++++----------------------- doc/stdlib/linalg.rst | 14 +++- test/sparsedir/cholmod.jl | 26 +++---- 5 files changed, 90 insertions(+), 103 deletions(-) diff --git a/NEWS.md b/NEWS.md index 83a67235fb396..aa5237b51e370 100644 --- a/NEWS.md +++ b/NEWS.md @@ -132,6 +132,10 @@ Library improvements * Added generic Cholesky factorization, and the Cholesky factorization is now parametrized on the matrix type ([#7236]). + * Sparse `cholfact` and `ldltfact` functions now accept a `perm` keyword + for user-provided permutations and a `shift` keyword to factorize + a shifted matrix ([#10844]). + * Add `svds` for sparse truncated SVD. ([#9425]) * Symmetric and Hermitian immutables are now parametrized on matrix type ([#7992]). @@ -1334,3 +1338,5 @@ Too numerous to mention. [#10543]: https://github.com/JuliaLang/julia/issues/10543 [#10659]: https://github.com/JuliaLang/julia/issues/10659 [#10709]: https://github.com/JuliaLang/julia/issues/10709 +[#10747]: https://github.com/JuliaLang/julia/issues/10747 +[#10844]: https://github.com/JuliaLang/julia/issues/10844 diff --git a/base/deprecated.jl b/base/deprecated.jl index d64330fe820e5..5987fc7f92af1 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -439,6 +439,9 @@ end @deprecate dlsym_e Libdl.dlsym_e @deprecate find_library Libdl.find_library +@deprecate cholfact(A::AbstractMatrix, β::Number) cholfact(A, shift=β) +@deprecate ldltfact(A::AbstractMatrix, β::Number) ldltfact(A, shift=β) + # 0.4 discontinued functions @noinline function subtypetree(x::DataType, level=-1) diff --git a/base/sparse/cholmod.jl b/base/sparse/cholmod.jl index 6fcf3f37ccf50..74a939318753a 100644 --- a/base/sparse/cholmod.jl +++ b/base/sparse/cholmod.jl @@ -33,6 +33,8 @@ const common_final_ll = (1:4) + cholmod_com_offsets[7] const common_print = (1:4) + cholmod_com_offsets[13] const common_itype = (1:4) + cholmod_com_offsets[18] const common_dtype = (1:4) + cholmod_com_offsets[19] +const common_nmethods = (1:4) + cholmod_com_offsets[15] +const common_postorder = (1:4) + cholmod_com_offsets[17] ## macro to generate the name of the C function according to the integer type macro cholmod_name(nm,typ) string("cholmod_", eval(typ) == SuiteSparse_long ? "l_" : "", nm) end @@ -283,7 +285,7 @@ function allocate_dense(nrow::Integer, ncol::Integer, d::Integer, ::Type{Complex d end -free_dense!{T}(p::Ptr{C_Dense{T}}) = ccall((:cholmod_l_free_dense, :libcholmod), Cint, (Ptr{Ptr{C_Dense{T}}}, Ptr{Void}), &p, common(Cint)) +free_dense!{T}(p::Ptr{C_Dense{T}}) = ccall((:cholmod_l_free_dense, :libcholmod), Cint, (Ref{Ptr{C_Dense{T}}}, Ptr{Void}), p, common(Cint)) function zeros{T<:VTypes}(m::Integer, n::Integer, ::Type{T}) d = Dense(ccall((:cholmod_l_zeros, :libcholmod), Ptr{C_Dense{T}}, @@ -568,7 +570,7 @@ for Ti in IndexTypes A end - function sdmult!{Tv<:VTypes}(A::Sparse{Tv,$Ti}, transpose::Bool, α::Tv, β::Tv, X::Dense{Tv}, Y::Dense{Tv}) + function sdmult!{Tv<:VTypes}(A::Sparse{Tv,$Ti}, transpose::Bool, α::Number, β::Number, X::Dense{Tv}, Y::Dense{Tv}) m, n = size(A) nc = transpose ? m : n nr = transpose ? n : m @@ -576,10 +578,10 @@ for Ti in IndexTypes throw(DimensionMismatch("incompatible dimensions, $nc and $(size(X,1))")) end @isok ccall((@cholmod_name("sdmult", $Ti),:libcholmod), Cint, - (Ptr{C_Sparse{Tv,$Ti}}, Cint, Ptr{Cdouble}, Ptr{Cdouble}, - Ptr{C_Dense{Tv}}, Ptr{C_Dense{Tv}}, Ptr{UInt8}), - A.p, transpose, &α, &β, - X.p, Y.p, common($Ti)) + (Ptr{C_Sparse{Tv,$Ti}}, Cint, + Ref{Complex128}, Ref{Complex128}, + Ptr{C_Dense{Tv}}, Ptr{C_Dense{Tv}}, Ptr{UInt8}), + A.p, transpose, α, β, X.p, Y.p, common($Ti)) Y end @@ -605,7 +607,7 @@ for Ti in IndexTypes end # cholmod_cholesky.h - # For analyze, factorize and factorize_p, the Common argument must be + # For analyze, analyze_p, and factorize_p!, the Common argument must be # supplied in order to control if the factorization is LLt or LDLt function analyze{Tv<:VTypes}(A::Sparse{Tv,$Ti}, cmmn::Vector{UInt8}) f = Factor(ccall((@cholmod_name("analyze", $Ti),:libcholmod), Ptr{C_Factor{Tv,$Ti}}, @@ -614,18 +616,22 @@ for Ti in IndexTypes finalizer(f, free!) f end - function factorize!{Tv<:VTypes}(A::Sparse{Tv,$Ti}, F::Factor{Tv,$Ti}, cmmn::Vector{UInt8}) - @isok ccall((@cholmod_name("factorize", $Ti),:libcholmod), Cint, - (Ptr{C_Sparse{Tv,$Ti}}, Ptr{C_Factor{Tv,$Ti}}, Ptr{UInt8}), - A.p, F.p, cmmn) - F + function analyze_p{Tv<:VTypes}(A::Sparse{Tv,$Ti}, perm::Vector{$Ti}, + cmmn::Vector{UInt8}) + length(perm) != size(A,1) && throw(BoundsError()) + f = Factor(ccall((@cholmod_name("analyze_p", $Ti),:libcholmod), Ptr{C_Factor{Tv,$Ti}}, + (Ptr{C_Sparse{Tv,$Ti}}, Ptr{$Ti}, Ptr{$Ti}, Csize_t, Ptr{UInt8}), + A.p, perm, C_NULL, 0, cmmn)) + finalizer(f, free!) + f end - function factorize_p!{Tv<:VTypes}(A::Sparse{Tv,$Ti}, β::Tv, fset::Vector{$Ti}, F::Factor{Tv,$Ti}, cmmn::Vector{UInt8}) + function factorize_p!{Tv<:VTypes}(A::Sparse{Tv,$Ti}, β::Real, F::Factor{Tv,$Ti}, cmmn::Vector{UInt8}) + # note that β is passed as a complex number (double beta[2]), + # but the CHOLMOD manual says that only beta[0] (real part) is used @isok ccall((@cholmod_name("factorize_p", $Ti),:libcholmod), Cint, - (Ptr{C_Sparse{Tv,$Ti}}, Ptr{Cdouble}, Ptr{$Ti}, Csize_t, + (Ptr{C_Sparse{Tv,$Ti}}, Ref{Complex128}, Ptr{$Ti}, Csize_t, Ptr{C_Factor{Tv,$Ti}}, Ptr{UInt8}), - A.p, &β, fset, length(fset), - F.p, cmmn) + A.p, β, C_NULL, 0, F.p, cmmn) F end @@ -984,42 +990,45 @@ Ac_mul_B(A::Sparse, B::VecOrMat) = Ac_mul_B(A, Dense(B)) ## Factorization methods -function cholfact(A::Sparse) + +function fact_{Tv<:VTypes,Ti<:ITypes,Ti2<:Integer}(A::Sparse{Tv,Ti}, cm::Array{UInt8}; + shift::Real=0.0, + perm::AbstractVector{Ti2}=Int[], + postorder::Bool=true, + userperm_only::Bool=true) sA = unsafe_load(A.p) sA.stype == 0 && throw(ArgumentError("sparse matrix is not symmetric/Hermitian")) - cm = common(indtype(A)) + if !postorder + cm[common_postorder] = reinterpret(UInt8, [zero(Cint)]) + end - # Hack! makes it a llt - cm[common_final_ll] = reinterpret(UInt8, [one(Cint)]) - F = analyze(A, cm) - factorize!(A, F, cm) - s = unsafe_load(F.p) - s.minor < size(A, 1) && throw(Base.LinAlg.PosDefException(s.minor)) + if isempty(perm) + F = analyze(A, cm) + else # user permutation provided + if userperm_only # use perm even if it is worse than AMD + cm[common_nmethods] = reinterpret(UInt8, [one(Cint)]) + end + F = analyze_p(A, Ti[p-1 for p in perm], cm) + end + + factorize_p!(A, shift, F, cm) return F end -function cholfact{Tv<:VTypes,Ti<:ITypes}(A::Sparse{Tv,Ti}, β::Tv) - sA = unsafe_load(A.p) - sA.stype == 0 && throw(ArgumentError("sparse matrix is not symmetric/Hermitian")) - - cm = common(Ti) +function cholfact(A::Sparse; kws...) + cm = common(indtype(A)) # Hack! makes it a llt cm[common_final_ll] = reinterpret(UInt8, [one(Cint)]) - F = analyze(A, cm) - factorize_p!(A, β, Ti[0:sA.ncol - 1;], F, cm) - + F = fact_(A, cm; kws...) s = unsafe_load(F.p) s.minor < size(A, 1) && throw(Base.LinAlg.PosDefException(s.minor)) return F end -function ldltfact(A::Sparse) - sA = unsafe_load(A.p) - sA.stype == 0 && throw(ArgumentError("sparse matrix is not symmetric/Hermitian")) - +function ldltfact(A::Sparse; kws...) cm = common(indtype(A)) # Hack! makes it a ldlt @@ -1028,72 +1037,29 @@ function ldltfact(A::Sparse) # Hack! really make sure it's a ldlt by avoiding supernodal factorisation cm[common_supernodal] = reinterpret(UInt8, [zero(Cint)]) - F = analyze(A, cm) - factorize!(A, F, cm) - - # Check if decomposition failed - s = unsafe_load(F.p) - s.minor < size(A, 1) && throw(Base.LinAlg.ArgumentError("matrix has one or more zero pivots")) - - return F -end - -function ldltfact{Tv<:VTypes,Ti<:ITypes}(A::Sparse{Tv,Ti}, β::Tv) - sA = unsafe_load(A.p) - sA.stype == 0 && throw(ArgumentError("sparse matrix is not symmetric/Hermitian")) - - cm = common(Ti) - - # Hack! makes it a ldlt - cm[common_final_ll] = reinterpret(UInt8, [zero(Cint)]) - - # Hack! really make sure it's a ldlt by avoiding supernodal factorisation - cm[common_supernodal] = reinterpret(UInt8, [zero(Cint)]) - - F = analyze(A, cm) - factorize_p!(A, β, Ti[0:sA.ncol - 1;], F, cm) - - # Check if decomposition failed + F = fact_(A, cm; kws...) s = unsafe_load(F.p) s.minor < size(A, 1) && throw(Base.LinAlg.ArgumentError("matrix has one or more zero pivots")) - return F end -cholfact{T<:VTypes}(A::Sparse{T}, β::Number) = cholfact(A, convert(T, β)) -cholfact(A::SparseMatrixCSC) = cholfact(Sparse(A)) -cholfact(A::SparseMatrixCSC, β::Number) = cholfact(Sparse(A), β) -cholfact{Ti}(A::Symmetric{Float64,SparseMatrixCSC{Float64,Ti}}) = cholfact(Sparse(A)) -cholfact{Ti}(A::Symmetric{Float64,SparseMatrixCSC{Float64,Ti}}, β::Number) = cholfact(Sparse(A), β) -cholfact{Ti}(A::Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},Ti}}) = cholfact(Sparse(A)) -cholfact{Ti}(A::Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},Ti}}, β::Number) = cholfact(Sparse(A), β) - -ldltfact{T<:VTypes}(A::Sparse{T}, β::Number) = ldltfact(A, convert(T, β)) -ldltfact(A::SparseMatrixCSC) = ldltfact(Sparse(A)) -ldltfact(A::SparseMatrixCSC, β::Number) = ldltfact(Sparse(A), β) -ldltfact{Ti}(A::Symmetric{Float64,SparseMatrixCSC{Float64,Ti}}) = ldltfact(Sparse(A)) -ldltfact{Ti}(A::Symmetric{Float64,SparseMatrixCSC{Float64,Ti}}, β::Number) = ldltfact(Sparse(A), β) -ldltfact{Ti}(A::Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},Ti}}) = ldltfact(Sparse(A)) -ldltfact{Ti}(A::Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},Ti}}, β::Number) = ldltfact(Sparse(A), β) - -function update!{Tv<:VTypes,Ti<:ITypes}(F::Factor{Tv,Ti}, A::Sparse{Tv,Ti}) - cm = common(Ti) - s = unsafe_load(F.p) - if Bool(s.is_ll) - cm[common_final_ll] = reinterpret(UInt8, [one(Cint)]) # Hack! makes it a llt +for f in (:cholfact, :ldltfact) + @eval begin + $f(A::SparseMatrixCSC; kws...) = $f(Sparse(A); kws...) + $f{Ti}(A::Symmetric{Float64,SparseMatrixCSC{Float64,Ti}}; kws...) = $f(Sparse(A); kws...) + $f{Ti}(A::Hermitian{Complex{Float64},SparseMatrixCSC{Complex{Float64},Ti}}; kws...) = $f(Sparse(A); kws...) end - factorize!(A, F, cm) end -function update!{Tv<:VTypes,Ti<:ITypes}(F::Factor{Tv,Ti}, A::Sparse{Tv,Ti}, β::Tv) + +function update!{Tv<:VTypes,Ti<:ITypes}(F::Factor{Tv,Ti}, A::Sparse{Tv,Ti}; shift::Real=0.0) cm = common(Ti) s = unsafe_load(F.p) if Bool(s.is_ll) cm[common_final_ll] = reinterpret(UInt8, [one(Cint)]) # Hack! makes it a llt end - factorize_p!(A, β, Ti[0:size(F, 1) - 1;], F, cm) + factorize_p!(A, shift, F, cm) end -update!{T<:VTypes}(F::Factor{T}, A::SparseMatrixCSC{T}) = update!(F, Sparse(A)) -update!{T<:VTypes}(F::Factor{T}, A::SparseMatrixCSC{T}, β::Number) = update!(F, Sparse(A), convert(T, β)) +update!{T<:VTypes}(F::Factor{T}, A::SparseMatrixCSC{T}; kws...) = update!(F, Sparse(A); kws...) ## Solvers diff --git a/doc/stdlib/linalg.rst b/doc/stdlib/linalg.rst index 91ecc2b351a31..ded3af51ab127 100644 --- a/doc/stdlib/linalg.rst +++ b/doc/stdlib/linalg.rst @@ -95,10 +95,15 @@ Linear algebra functions in Julia are largely implemented by calling functions f Compute the Cholesky factorization of a dense symmetric positive (semi)definite matrix ``A`` and return either a ``Cholesky`` if ``pivot==Val{false}`` or ``CholeskyPivoted`` if ``pivot==Val{true}``. ``LU`` may be ``:L`` for using the lower part or ``:U`` for the upper part. The default is to use ``:U``. The triangular matrix can be obtained from the factorization ``F`` with: ``F[:L]`` and ``F[:U]``. The following functions are available for ``Cholesky`` objects: ``size``, ``\``, ``inv``, ``det``. For ``CholeskyPivoted`` there is also defined a ``rank``. If ``pivot==Val{false}`` a ``PosDefException`` exception is thrown in case the matrix is not positive definite. The argument ``tol`` determines the tolerance for determining the rank. For negative values, the tolerance is the machine precision. -.. function:: cholfact(A) -> CHOLMOD.Factor +.. function:: cholfact(A; shift=0, perm=Int[]) -> CHOLMOD.Factor Compute the Cholesky factorization of a sparse positive definite matrix ``A``. A fill-reducing permutation is used. The main application of this type is to solve systems of equations with ``\``, but also the methods ``diag``, ``det``, ``logdet`` are defined. The function calls the C library CHOLMOD and many other functions from the library are wrapped but not exported. + Setting optional ``shift`` keyword argument computes the factorization + of ``A+shift*I`` instead of ``A``. If the ``perm`` argument is nonempty, + it should be a permutation of `1:size(A,1)` giving the ordering to use + (instead of CHOLMOD's default AMD ordering). + .. function:: cholfact!(A [,LU=:U [,pivot=Val{false}]][;tol=-1.0]) -> Cholesky ``cholfact!`` is the same as :func:`cholfact`, but saves space by overwriting the input ``A``, instead of creating a copy. ``cholfact!`` can also reuse the symbolic factorization from a different matrix ``F`` with the same structure when used as: ``cholfact!(F::CholmodFactor, A)``. @@ -107,10 +112,15 @@ Linear algebra functions in Julia are largely implemented by calling functions f Compute a factorization of a positive definite matrix ``A`` such that ``A=L*Diagonal(d)*L'`` where ``L`` is a unit lower triangular matrix and ``d`` is a vector with non-negative elements. -.. function:: ldltfact(A) -> CHOLMOD.Factor +.. function:: ldltfact(A; shift=0, perm=Int[]) -> CHOLMOD.Factor Compute the LDLt factorization of a sparse symmetric or Hermitian matrix ``A``. A fill-reducing permutation is used. The main application of this type is to solve systems of equations with ``\``, but also the methods ``diag``, ``det``, ``logdet`` are defined. The function calls the C library CHOLMOD and many other functions from the library are wrapped but not exported. + Setting optional ``shift`` keyword argument computes the factorization + of ``A+shift*I`` instead of ``A``. If the ``perm`` argument is nonempty, + it should be a permutation of `1:size(A,1)` giving the ordering to use + (instead of CHOLMOD's default AMD ordering). + .. function:: qr(A [,pivot=Val{false}][;thin=true]) -> Q, R, [p] Compute the (pivoted) QR factorization of ``A`` such that either ``A = Q*R`` or ``A[:,p] = Q*R``. Also see ``qrfact``. The default is to compute a thin factorization. Note that ``R`` is not extended with zeros when the full ``Q`` is requested. diff --git a/test/sparsedir/cholmod.jl b/test/sparsedir/cholmod.jl index 82915f68e02bf..0635ef888b954 100644 --- a/test/sparsedir/cholmod.jl +++ b/test/sparsedir/cholmod.jl @@ -98,6 +98,7 @@ chma = ldltfact(A) # LDL' form @test unsafe_load(chma.p).is_ll == 0 # check that it is in fact an LDLt x = chma\B @test_approx_eq x ones(size(x)) +@test nnz(ldltfact(A, perm=1:size(A,1))) > nnz(chma) chma = cholfact(A) # LL' form @test CHOLMOD.isvalid(chma) @@ -105,6 +106,7 @@ chma = cholfact(A) # LL' form x = chma\B @test_approx_eq x ones(size(x)) @test nnz(chma) == 489 +@test nnz(cholfact(A, perm=1:size(A,1))) > nnz(chma) #lp_afiro example afiro = CHOLMOD.Sparse(27, 51, @@ -360,13 +362,13 @@ for elty in (Float64, Complex{Float64}) # Factor @test_throws ArgumentError cholfact(A1) @test_throws Base.LinAlg.PosDefException cholfact(A1 + A1' - 2eigmax(full(A1 + A1'))I) - @test_throws Base.LinAlg.PosDefException cholfact(A1 + A1', -2eigmax(full(A1 + A1'))) + @test_throws Base.LinAlg.PosDefException cholfact(A1 + A1', shift=-2eigmax(full(A1 + A1'))) @test_throws ArgumentError ldltfact(A1 + A1' - 2real(A1[1,1])I) - @test_throws ArgumentError ldltfact(A1 + A1', -2real(A1[1,1])) + @test_throws ArgumentError ldltfact(A1 + A1', shift=-2real(A1[1,1])) @test_throws ArgumentError cholfact(A1) - @test_throws ArgumentError cholfact(A1, 1.0) + @test_throws ArgumentError cholfact(A1, shift=1.0) @test_throws ArgumentError ldltfact(A1) - @test_throws ArgumentError ldltfact(A1, 1.0) + @test_throws ArgumentError ldltfact(A1, shift=1.0) F = cholfact(A1pd) tmp = IOBuffer() show(tmp, F) @@ -391,16 +393,16 @@ for elty in (Float64, Complex{Float64}) if elty <: Real @test CHOLMOD.issym(Sparse(A1pd, 0)) @test CHOLMOD.Sparse(cholfact(Symmetric(A1pd, :L))) == CHOLMOD.Sparse(cholfact(A1pd)) - @test CHOLMOD.Sparse(cholfact(Symmetric(A1pd, :L), 2)) == CHOLMOD.Sparse(cholfact(A1pd, 2)) + @test CHOLMOD.Sparse(cholfact(Symmetric(A1pd, :L), shift=2)) == CHOLMOD.Sparse(cholfact(A1pd, shift=2)) @test CHOLMOD.Sparse(ldltfact(Symmetric(A1pd, :L))) == CHOLMOD.Sparse(ldltfact(A1pd)) - @test CHOLMOD.Sparse(ldltfact(Symmetric(A1pd, :L), 2)) == CHOLMOD.Sparse(ldltfact(A1pd, 2)) + @test CHOLMOD.Sparse(ldltfact(Symmetric(A1pd, :L), shift=2)) == CHOLMOD.Sparse(ldltfact(A1pd, shift=2)) else @test !CHOLMOD.issym(Sparse(A1pd, 0)) @test CHOLMOD.ishermitian(Sparse(A1pd, 0)) @test CHOLMOD.Sparse(cholfact(Hermitian(A1pd, :L))) == CHOLMOD.Sparse(cholfact(A1pd)) - @test CHOLMOD.Sparse(cholfact(Hermitian(A1pd, :L), 2)) == CHOLMOD.Sparse(cholfact(A1pd, 2)) + @test CHOLMOD.Sparse(cholfact(Hermitian(A1pd, :L), shift=2)) == CHOLMOD.Sparse(cholfact(A1pd, shift=2)) @test CHOLMOD.Sparse(ldltfact(Hermitian(A1pd, :L))) == CHOLMOD.Sparse(ldltfact(A1pd)) - @test CHOLMOD.Sparse(ldltfact(Hermitian(A1pd, :L), 2)) == CHOLMOD.Sparse(ldltfact(A1pd, 2)) + @test CHOLMOD.Sparse(ldltfact(Hermitian(A1pd, :L), shift=2)) == CHOLMOD.Sparse(ldltfact(A1pd, shift=2)) end @@ -414,17 +416,17 @@ for elty in (Float64, Complex{Float64}) @test size(F, 3) == 1 @test_throws ArgumentError size(F, 0) - F = cholfact(A1pdSparse, 2) + F = cholfact(A1pdSparse, shift=2) @test isa(CHOLMOD.Sparse(F), CHOLMOD.Sparse{elty, CHOLMOD.SuiteSparse_long}) - @test_approx_eq CHOLMOD.Sparse(CHOLMOD.update!(copy(F), A1pd, 2.0)) CHOLMOD.Sparse(F) # surprisingly, this can cause small ulp size changes so we cannot test exact equality + @test_approx_eq CHOLMOD.Sparse(CHOLMOD.update!(copy(F), A1pd, shift=2.0)) CHOLMOD.Sparse(F) # surprisingly, this can cause small ulp size changes so we cannot test exact equality F = ldltfact(A1pd) @test isa(CHOLMOD.Sparse(F), CHOLMOD.Sparse{elty, CHOLMOD.SuiteSparse_long}) @test_approx_eq CHOLMOD.Sparse(CHOLMOD.update!(copy(F), A1pd)) CHOLMOD.Sparse(F) # surprisingly, this can cause small ulp size changes so we cannot test exact equality - F = ldltfact(A1pdSparse, 2) + F = ldltfact(A1pdSparse, shift=2) @test isa(CHOLMOD.Sparse(F), CHOLMOD.Sparse{elty, CHOLMOD.SuiteSparse_long}) - @test_approx_eq CHOLMOD.Sparse(CHOLMOD.update!(copy(F), A1pd, 2.0)) CHOLMOD.Sparse(F) # surprisingly, this can cause small ulp size changes so we cannot test exact equality + @test_approx_eq CHOLMOD.Sparse(CHOLMOD.update!(copy(F), A1pd, shift=2.0)) CHOLMOD.Sparse(F) # surprisingly, this can cause small ulp size changes so we cannot test exact equality @test isa(CHOLMOD.factor_to_sparse!(F), CHOLMOD.Sparse) @test_throws CHOLMOD.CHOLMODException CHOLMOD.factor_to_sparse!(F)