reuse cache for multiple solves

Author: vpuri3

src/common.jl

@@ -11,7 +11,7 @@ struct LinearCache{TA,Tb,Tu,Tp,Talg,Tc,Tr,Tl}
 #   k::Tk # iteration count
 end
 
-function set_A(cache, A)
+function set_A(cache, A) # and ! to function name
     @set! cache.A = A
     @set! cache.isfresh = true
 end
@@ -20,6 +20,10 @@ function set_b(cache, b)
     @set! cache.b = b
 end
 
+function set_u(cache, u)
+    @set! cache.u = u
+end
+
 function set_p(cache, p)
     @set! cache.p = p
     # @set! cache.isfresh = true
@@ -57,10 +61,8 @@ function SciMLBase.init(
     Pr = LinearAlgebra.I
     Pl = LinearAlgebra.I
 
-#   @show (A, b, u0, p) |> typeof
-
-    A = alias_A ? A : copy(A)
-    b = alias_b ? b : copy(b)
+    A = alias_A ? A : deepcopy(A)
+    b = alias_b ? b : deepcopy(b)
 
     cache = LinearCache{
         typeof(A),
@@ -95,3 +97,21 @@ function (alg::SciMLLinearSolveAlgorithm)(x,A,b,args...;u0=nothing,kwargs...)
     x = solve(prob,alg,args...;kwargs...)
     return x
 end
+
+# how to initialize cahce?
+
+# use the same cache to solve multiple linear problems
+function (cache::LinearCache)(x,A,b,args...;u0=nothing,kwargs...)
+    set_A(cache, A)
+    set_b(cache, b)
+
+    if u0 == nothing
+        x = zero(x)
+    else
+        x = u0
+    end
+    set_u(cache, x)
+
+    x = solve(cache)
+    return x
+end

test/runtests.jl

@@ -24,20 +24,23 @@ using Test
 
     # Factorization
     for alg in (:LUFactorization, :QRFactorization, :SVDFactorization,
-                :KrylovJL)
+                :KrylovJL,
+#               :KrylovKitJL,
+               )
         @eval begin
             @test $A * solve($prob, $alg();) ≈ $b
-            $alg()($x,$A,$b)
+            $alg()($x, $A, $b)
+            @test $A * $x ≈ $b
+
+            cache = SciMLBase.init($prob, $alg())
+            cache($x, $A, $b)
             @test $A * $x ≈ $b
         end
     end
 
     # test on some ODEProblem
-    using OrdinaryDiffEq
-    using DiffEqProblemLibrary.ODEProblemLibrary
-    ODEProblemLibrary.importodeproblems()
-
-    #   add this problem to DiffEqProblemLibrary
+#   using OrdinaryDiffEq
+#   #   add this problem to DiffEqProblemLibrary
 #   kx = 1
 #   kt = 1
 #   ut(x,t) = sin(kx*pi*x)*cos(kt*pi*t)
@@ -50,8 +53,12 @@ using Test
 #   func = ODEFunction(dudt!)
 #   prob = ODEProblem(func,u0,tspn)
 
-    prob = ODEProblemLibrary.prob_ode_linear
-    sol  = solve(prob, Rodas5(linsolve=SVDFactorization()); saveat=0.1)
-    @show sol.retcode
+#   using OrdinaryDiffEq
+#   using DiffEqProblemLibrary.ODEProblemLibrary
+#   ODEProblemLibrary.importodeproblems()
+#   prob = ODEProblemLibrary.prob_ode_linear
+#   @show prob
+#   sol  = solve(prob, Rodas5(linsolve=KrylovJL()); saveat=0.1)
+#   @show sol.retcode
 
 end