placeholders

Author: vpuri3

Project.toml

@@ -9,8 +9,10 @@ IterativeSolvers = "42fd0dbc-a981-5370-80f2-aaf504508153"
 Krylov = "ba0b0d4f-ebba-5204-a429-3ac8c609bfb7"
 KrylovKit = "0b1a1467-8014-51b9-945f-bf0ae24f4b77"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+Pardiso = "46dd5b70-b6fb-5a00-ae2d-e8fea33afaf2"
 RecursiveFactorization = "f2c3362d-daeb-58d1-803e-2bc74f2840b4"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
+Requires = "ae029012-a4dd-5104-9daa-d747884805df"
 SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"

src/LinearSolve.jl

@@ -10,11 +10,13 @@ using SparseArrays
 using SciMLBase: AbstractDiffEqOperator, AbstractLinearAlgorithm
 using Setfield
 using UnPack
+using Requires
 
 # wrap
 import Krylov
-import KrylovKit
+import KrylovKit # TODO
 import IterativeSolvers
+import Pardiso
 
 using Reexport
 @reexport using SciMLBase
@@ -38,14 +40,16 @@ function __init__()
   else
     IS_OPENBLAS[] = occursin("openblas", BLAS.get_config().loaded_libs[1].libname)
   end
+
+  @require Pardiso="46dd5b70-b6fb-5a00-ae2d-e8fea33afaf2" PardisoJL, PardisoJL_Default
 end
 
 export LUFactorization, SVDFactorization, QRFactorization, GenericFactorization,
        RFLUFactorizaation
-export KrylovJL, KrylovJL_CG, KrylovJL_GMRES, KrylovJL_BICGSTAB,
-       KrylovJL_MINRES,
+export KrylovJL, KrylovJL_CG, KrylovJL_GMRES, KrylovJL_BICGSTAB, KrylovJL_MINRES,
        IterativeSolversJL, IterativeSolversJL_CG, IterativeSolversJL_GMRES,
-       IterativeSolversJL_BICGSTAB, IterativeSolversJL_MINRES
+       IterativeSolversJL_BICGSTAB, IterativeSolversJL_MINRES,
+       PardisoJL, PardisoJL_Default
 export DefaultLinSolve
 
 end

src/wrappers.jl

@@ -297,3 +297,63 @@ function SciMLBase.solve(cache::LinearCache, alg::IterativeSolversJL; kwargs...)
 
     return cache.u
 end
+
+## Paradiso
+
+struct PardisoJL{F,S,P,A,K} <: SciMLLinearSolveAlgorithm
+    alg_selector::F
+    schur::S
+    parallel::P
+    args::A
+    kwargs::K
+end
+
+function PardisoJL(args...; alg_selector=:default_alg,
+                          schur=:smart,
+                          parallel=:fancy,
+                          kwargs...)
+
+    return PardisoJL(alg_selector, schur, parallel,
+                   args, kwargs)
+end
+
+PardisoJL_Default(args...; kwargs...) = PardisoJL(args...;
+                                             alg_selector=:this_one,
+                                             kwargs...)
+
+function init_cacheval(alg::PardisoJL, cache::LinearCache)
+
+    solver =
+    if alg.alg_selector == :default_alg
+        MKLPardisoSolver()
+    elseif alg_selector == :this_one
+        PardisoSolver()
+    else
+        PardisoSolver()
+    end
+
+    alg.verbose && set_msglvl!(solver, Pardiso.MESSAGE_LEVEL_ON)
+
+    return solver
+end
+
+function SciMLBase.solve(cache::LinearCache, alg::PardisoJL; kwargs...)
+    @unpack cacheval, A, b, u = cache
+
+    if cache.isfresh
+        solver = init_cacheval(alg, cache)
+        cache = set_cacheval(cache, solver)
+    end
+
+    abstol  = cache.abstol
+    reltol  = cache.reltol
+    verbose = cache.verbose
+
+    kwargs = (abstol=abstol, reltol=reltol,
+              alg.kwargs...)
+
+    Pardiso.solve!(cacheval, u, A, b)
+
+    return cache.u
+end
+