Merge pull request #255 from j-fu/sparspak_defaults

More on sparspak

Author: ChrisRackauckas

README.md

@@ -19,6 +19,7 @@ maximum efficiency. Specifically, LinearSolve.jl includes:
 - KLU for faster sparse LU factorization on unstructured matrices
 - UMFPACK for faster sparse LU factorization on matrices with some repeated structure
 - MKLPardiso wrappers for handling many sparse matrices faster than SuiteSparse (KLU, UMFPACK) methods
+- Sparspak.jl for sparse LU factorization in pure Julia for generic number types and for non-GPL distributions
 - GPU-offloading for large dense matrices
 - Wrappers to all of the Krylov implementations (Krylov.jl, IterativeSolvers.jl, KrylovKit.jl) for easy
   testing of all of them. LinearSolve.jl handles the API differences, especially with the preconditioner

src/LinearSolve.jl

@@ -75,13 +75,19 @@ end
         A = sprand(4, 4, 0.3) + I
         b = rand(4)
         prob = LinearProblem(A, b)
-        sol = solve(prob)
         sol = solve(prob, KLUFactorization())
         sol = solve(prob, UMFPACKFactorization())
-        sol = solve(prob, SparspakFactorization())
     end
 end
 
+SnoopPrecompile.@precompile_all_calls begin
+    A = sprand(4, 4, 0.3) + I
+    b = rand(4)
+    prob = LinearProblem(A, b)
+    sol = solve(prob)
+    sol = solve(prob, SparspakFactorization())
+end
+
 export LUFactorization, SVDFactorization, QRFactorization, GenericFactorization,
        GenericLUFactorization, SimpleLUFactorization, RFLUFactorization,
        UMFPACKFactorization, KLUFactorization, FastLUFactorization, FastQRFactorization,

src/default.jl

@@ -43,18 +43,20 @@ function defaultalg(A::Diagonal, b, ::OperatorAssumptions{Nothing})
     DiagonalFactorization()
 end
 
+function defaultalg(A::SparseMatrixCSC{Tv, Ti}, b,
+                    ::OperatorAssumptions{true}) where {Tv, Ti}
+    SparspakFactorization()
+end
+
 @static if INCLUDE_SPARSE
-    function defaultalg(A::SparseMatrixCSC, b, ::OperatorAssumptions{true})
+    function defaultalg(A::SparseMatrixCSC{<:Union{Float64,ComplexF64}, Ti}, b,
+                        ::OperatorAssumptions{true}) where {Ti}
         if length(b) <= 10_000
             KLUFactorization()
         else
             UMFPACKFactorization()
         end
     end
-else
-    function defaultalg(A::SparseMatrixCSC, b, ::OperatorAssumptions{true})
-        SparspakFactorization()
-    end
 end
 
 function defaultalg(A::GPUArraysCore.AbstractGPUArray, b, ::OperatorAssumptions{true})

src/factorization.jl

@@ -511,7 +511,7 @@ function init_cacheval(::SparspakFactorization, A, b, u, Pl, Pr, maxiters::Int,
                        reltol,
                        verbose::Bool, assumptions::OperatorAssumptions)
     A = convert(AbstractMatrix, A)
-    return sparspaklu(A, factorize=false)
+    return sparspaklu(A, factorize = false)
 end
 
 function SciMLBase.solve(cache::LinearCache, alg::SparspakFactorization; kwargs...)