Linearsolve2 (#74)

* mutable linsolve * linearsolve 2.0.1 * update project.toml
WIAS-PDELib · Jun 20, 2023 · c6c9ece · c6c9ece · j-fu · Jun 20, 2023
1 parent 7bb7d58
commit c6c9ece
Show file tree

Hide file tree

Showing 7 changed files with 25 additions and 40 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,7 +1,7 @@
 name = "VoronoiFVM"
 uuid = "82b139dc-5afc-11e9-35da-9b9bdfd336f3"
 authors = ["Jürgen Fuhrmann <[email protected]>", "Dilara Abdel", "Jan Weidner", "Alexander Seiler", "Patricio Farrell", "Matthias Liero"]
-version = "1.7.1"
+version = "1.8.0"
 
 [deps]
 BandedMatrices = "aae01518-5342-5314-be14-df237901396f"
@@ -39,7 +39,7 @@ ExtendableSparse = "1.2"
 ForwardDiff = "0.10.35"
 GridVisualize = "0.5.2,0.6.1,1"
 JLD2 = "0.4.29"
-LinearSolve = "1.37"
+LinearSolve = "2.2"
 PrecompileTools = "1"
 RecursiveArrayTools = "^2.14.2"
 RecursiveFactorization = "0.2.18"

diff --git a/examples/Example205_NonlinearPoisson2D.jl b/examples/Example205_NonlinearPoisson2D.jl
@@ -13,7 +13,7 @@ using ILUZero
 
 function main(; n = 10, Plotter = nothing, verbose = false, unknown_storage = :sparse,
               method_linear = nothing, assembly=:edgewise,
-              precon_linear = A -> LinearSolve.Identity())
+              precon_linear = A -> VoronoiFVM.Identity())
     h = 1.0 / convert(Float64, n)
     X = collect(0.0:h:1.0)
     Y = collect(0.0:h:1.0)

diff --git a/pluto-examples/api-update.jl b/pluto-examples/api-update.jl
@@ -225,7 +225,7 @@ md"""
         sys0;
         inival = 0.1,
         method_linear = KrylovJL_BICGSTAB(),
-        precon_linear = A -> factorize(A, SparspakFactorization()),
+        precon_linear = SparspakFactorization(),
     verbose="nlad"
 	)
 

diff --git a/src/VoronoiFVM.jl b/src/VoronoiFVM.jl
@@ -13,7 +13,6 @@ using BandedMatrices
 # using MultidiagonalMatrices
 
 using LinearSolve
-using LinearSolve: Identity
 using Statistics
 
 using ForwardDiff

diff --git a/src/vfvm_linsolve.jl b/src/vfvm_linsolve.jl
@@ -204,45 +204,27 @@ function Random.rand(
 end
 
 
-mutable struct FactorizationPreconditioner{C}
-    cache::C
-end
-
 """
-   factorize(A,method::LinearSolve.AbstractFactorization)
-
-Calculate an LU factorization of A using one of the methods available in LinearSolve.jl. 
+    Make preconditioner constructors from methods
 """
-function LinearAlgebra.factorize(A, method::LinearSolve.AbstractFactorization)
+
+function (method::LinearSolve.AbstractFactorization)(A)
     pr = LinearProblem(A, zeros(eltype(A), size(A, 1)))
-    p = FactorizationPreconditioner(init(pr, method))
-    p
+    init(pr, method)
 end
 
-function LinearAlgebra.factorize(A, method::LinearSolve.SciMLLinearSolveAlgorithm)
+function (method::LinearSolve.SciMLLinearSolveAlgorithm)(A)
     pr = LinearProblem(SparseMatrixCSC(A), zeros(eltype(A), size(A, 1)))
-    p = FactorizationPreconditioner(init(pr, method))
-    p
-end
-
-
-function (f::LinearSolve.AbstractFactorization)(A)
-    factorize(A, f)
-end
-
-function (f::LinearSolve.SciMLLinearSolveAlgorithm)(A)
-    factorize(SparseMatrixCSC(A), f)
+    init(pr, method)
 end
 
 function (f::ExtendableSparse.AbstractFactorization)(A)
     factorize!(f, A)
 end
 
-
-function LinearAlgebra.ldiv!(u, p::FactorizationPreconditioner, b)
-    p.cache = LinearSolve.set_b(p.cache, b)
-    sol = solve(p.cache)
-    p.cache = sol.cache
+function LinearAlgebra.ldiv!(u, cache::LinearSolve.LinearCache, b)
+    cache.b=b
+    sol = solve!(cache)
     copyto!(u, sol.u)
 end
 
@@ -260,19 +242,16 @@ function _solve_linear!(u, system, nlhistory, control, method_linear, A, b)
             Pl,
         )
     else
-        system.linear_cache = LinearSolve.set_A(system.linear_cache, SparseMatrixCSC(A))
-        system.linear_cache = LinearSolve.set_b(system.linear_cache, b)
+        system.linear_cache.A=SparseMatrixCSC(A)
+        system.linear_cache.b=b
         if control.keepcurrent_linear
-            Pl = control.precon_linear(SparseMatrixCSC(A))
             nlhistory.nlu += 1
-            system.linear_cache =
-                LinearSolve.set_prec(system.linear_cache, Pl, LinearSolve.Identity())
+            system.linear_cache.Pl=control.precon_linear(SparseMatrixCSC(A))
         end
     end
 
     try
-        sol = LinearSolve.solve(system.linear_cache)
-        system.linear_cache = sol.cache
+        sol = LinearSolve.solve!(system.linear_cache)
         u .= sol.u
         nliniter = sol.iters
         nlhistory.nlin = sol.iters

diff --git a/src/vfvm_solvercontrol.jl b/src/vfvm_solvercontrol.jl
@@ -1,5 +1,12 @@
 ################################################
 const FactorizationStrategy= Union{Nothing, Function, Type, ExtendableSparse.AbstractFactorization, LinearSolve.AbstractFactorization,LinearSolve.SciMLLinearSolveAlgorithm}
+
+struct Identity end
+Identity(A)=Identity()
+LinearAlgebra.ldiv!(u,I::Identity,v)=u.=v
+LinearAlgebra.ldiv!(I::Identity,u)=nothing
+
+
 """
     SolverControl
     SolverControl(;kwargs...)

diff --git a/src/vfvm_testfunctions.jl b/src/vfvm_testfunctions.jl
@@ -79,7 +79,7 @@ function testfunction(factory::TestFunctionFactory, bc0, bc1)
         method_linear = UMFPACKFactorization()
     end
 
-    p = LinearProblem(factory.tfsystem.matrix, vec(f))
+    p = LinearProblem(SparseMatrixCSC(factory.tfsystem.matrix), vec(f))
     sol = solve(p, method_linear)
     sol.u
 end