Update cpr.jl

sintefmath · Feb 21, 2023 · 92983d8 · 92983d8
1 parent feefc52
commit 92983d8
Showing 1 changed file with 18 additions and 12 deletions.
diff --git a/src/cpr.jl b/src/cpr.jl
@@ -188,25 +188,31 @@ function apply!(x, cpr::CPRPreconditioner, r, arg...)
         # Apply preconditioner to pressure part
         @tic "p apply" begin
             p_rtol = cpr.p_rtol
-            apply!(Δp, cpr.pressure_precond, r_p)
-            if !isnothing(p_rtol)
-                A_p = cpr.A_p
-                if isnothing(cpr.psolver)
-                    cpr.psolver = FgmresSolver(A_p, r_p)
-                end
-                psolve = cpr.psolver
-                warm_start!(psolve, Δp)
-                M = Jutul.PrecondWrapper(linear_operator(cpr.pressure_precond))
-                fgmres!(psolve, A_p, r_p, M = M, rtol = p_rtol, itmax = 20)
-                @. Δp = psolve.x
-            end
+            p_precond = cpr.pressure_precond
+            cpr_p_apply!(Δp, cpr, p_precond, r_p, p_rtol)
         end
         @tic "r update" correct_residual_for_dp!(y, x, Δp, bz, cpr.buf, cpr.A_ps)
         @tic "s apply" apply!(x, cpr.system_precond, y)
         @tic "Δp" increment_pressure!(x, Δp, bz)
     end
 end
 
+
+function cpr_p_apply!(Δp, cpr, p_precond, r_p, p_rtol)
+    apply!(Δp, p_precond, r_p)
+    if !isnothing(p_rtol)
+        A_p = cpr.A_p
+        if isnothing(cpr.psolver)
+            cpr.psolver = FgmresSolver(A_p, r_p)
+        end
+        psolve = cpr.psolver
+        warm_start!(psolve, Δp)
+        M = Jutul.PrecondWrapper(linear_operator(p_precond))
+        fgmres!(psolve, A_p, r_p, M = M, rtol = p_rtol, itmax = 20)
+        @. Δp = psolve.x
+    end
+end
+
 reservoir_residual(lsys) = lsys.r
 reservoir_jacobian(lsys) = lsys.jac