From 0fcc986322dea32e799080dd9b96958b30c2face Mon Sep 17 00:00:00 2001 From: JordiManyer Date: Thu, 28 Dec 2023 18:48:30 +1100 Subject: [PATCH] Add nonlinear numerical setup update --- src/BlockSolvers/BlockTriangularSolvers.jl | 12 ++++++++++-- src/LinearSolvers/Krylov/CGSolvers.jl | 7 +++++++ src/LinearSolvers/Krylov/FGMRESSolvers.jl | 10 ++++++++++ src/LinearSolvers/Krylov/GMRESSolvers.jl | 12 ++++++++++++ src/LinearSolvers/Krylov/MINRESSolvers.jl | 11 +++++++++++ src/NonlinearSolvers/NewtonRaphsonSolver.jl | 2 +- src/NonlinearSolvers/NonlinearSolvers.jl | 1 + src/SolverInterfaces/GridapExtras.jl | 7 +------ 8 files changed, 53 insertions(+), 9 deletions(-) diff --git a/src/BlockSolvers/BlockTriangularSolvers.jl b/src/BlockSolvers/BlockTriangularSolvers.jl index a2afb5fd..0fa9b622 100644 --- a/src/BlockSolvers/BlockTriangularSolvers.jl +++ b/src/BlockSolvers/BlockTriangularSolvers.jl @@ -91,7 +91,11 @@ function Gridap.Algebra.numerical_setup!(ns::BlockTriangularSolverNS,mat::Abstra solver = ns.solver mat_blocks = blocks(mat) block_caches = map(update_block_cache!,ns.block_caches,solver.blocks,mat_blocks) - map(numerical_setup!,ns.block_ns,diag(block_caches)) + map(diag(solver.blocks),ns.block_ns,diag(block_caches)) do bi, nsi, ci + if is_nonlinear(bi) + numerical_setup!(nsi,ci) + end + end return ns end @@ -102,7 +106,11 @@ function Gridap.Algebra.numerical_setup!(ns::BlockTriangularSolverNS,mat::Abstra block_caches = map(CartesianIndices(solver.blocks)) do I update_block_cache!(ns.block_caches[I],mat_blocks[I],vec_blocks[I[2]]) end - map(numerical_setup!,ns.block_ns,diag(block_caches),vec_blocks) + map(diag(solver.blocks),ns.block_ns,diag(block_caches),vec_blocks) do bi, nsi, ci, xi + if is_nonlinear(bi) + numerical_setup!(nsi,ci,xi) + end + end return ns end diff --git a/src/LinearSolvers/Krylov/CGSolvers.jl b/src/LinearSolvers/Krylov/CGSolvers.jl index 12a69c3a..081da65d 100644 --- a/src/LinearSolvers/Krylov/CGSolvers.jl +++ b/src/LinearSolvers/Krylov/CGSolvers.jl @@ -46,6 +46,13 @@ end function Gridap.Algebra.numerical_setup!(ns::CGNumericalSetup, A::AbstractMatrix) numerical_setup!(ns.Pl_ns,A) ns.A = A + return ns +end + +function Gridap.Algebra.numerical_setup!(ns::CGNumericalSetup, A::AbstractMatrix, x::AbstractVector) + numerical_setup!(ns.Pl_ns,A,x) + ns.A = A + return ns end function Gridap.Algebra.solve!(x::AbstractVector,ns::CGNumericalSetup,b::AbstractVector) diff --git a/src/LinearSolvers/Krylov/FGMRESSolvers.jl b/src/LinearSolvers/Krylov/FGMRESSolvers.jl index 31adf264..c5de93a1 100644 --- a/src/LinearSolvers/Krylov/FGMRESSolvers.jl +++ b/src/LinearSolvers/Krylov/FGMRESSolvers.jl @@ -62,6 +62,16 @@ function Gridap.Algebra.numerical_setup!(ns::FGMRESNumericalSetup, A::AbstractMa numerical_setup!(ns.Pl_ns,A) end ns.A = A + return ns +end + +function Gridap.Algebra.numerical_setup!(ns::FGMRESNumericalSetup, A::AbstractMatrix, x::AbstractVector) + numerical_setup!(ns.Pr_ns,A,x) + if !isa(ns.Pl_ns,Nothing) + numerical_setup!(ns.Pl_ns,A,x) + end + ns.A = A + return ns end function Gridap.Algebra.solve!(x::AbstractVector,ns::FGMRESNumericalSetup,b::AbstractVector) diff --git a/src/LinearSolvers/Krylov/GMRESSolvers.jl b/src/LinearSolvers/Krylov/GMRESSolvers.jl index b2e5de7b..14d603e9 100644 --- a/src/LinearSolvers/Krylov/GMRESSolvers.jl +++ b/src/LinearSolvers/Krylov/GMRESSolvers.jl @@ -63,6 +63,18 @@ function Gridap.Algebra.numerical_setup!(ns::GMRESNumericalSetup, A::AbstractMat numerical_setup!(ns.Pl_ns,A) end ns.A = A + return ns +end + +function Gridap.Algebra.numerical_setup!(ns::GMRESNumericalSetup, A::AbstractMatrix, x::AbstractVector) + if !isa(ns.Pr_ns,Nothing) + numerical_setup!(ns.Pr_ns,A,x) + end + if !isa(ns.Pl_ns,Nothing) + numerical_setup!(ns.Pl_ns,A,x) + end + ns.A = A + return ns end function Gridap.Algebra.solve!(x::AbstractVector,ns::GMRESNumericalSetup,b::AbstractVector) diff --git a/src/LinearSolvers/Krylov/MINRESSolvers.jl b/src/LinearSolvers/Krylov/MINRESSolvers.jl index fcf4504b..8feb9fb1 100644 --- a/src/LinearSolvers/Krylov/MINRESSolvers.jl +++ b/src/LinearSolvers/Krylov/MINRESSolvers.jl @@ -62,6 +62,17 @@ function Gridap.Algebra.numerical_setup!(ns::MINRESNumericalSetup, A::AbstractMa ns.A = A end +function Gridap.Algebra.numerical_setup!(ns::MINRESNumericalSetup, A::AbstractMatrix, x::AbstractVector) + if !isa(ns.Pr_ns,Nothing) + numerical_setup!(ns.Pr_ns,A,x) + end + if !isa(ns.Pl_ns,Nothing) + numerical_setup!(ns.Pl_ns,A,x) + end + ns.A = A + return ns +end + function Gridap.Algebra.solve!(x::AbstractVector,ns::MINRESNumericalSetup,b::AbstractVector) solver, A, Pl, Pr, caches = ns.solver, ns.A, ns.Pl_ns, ns.Pr_ns, ns.caches V, W, zr, zl, H, g, c, s = caches diff --git a/src/NonlinearSolvers/NewtonRaphsonSolver.jl b/src/NonlinearSolvers/NewtonRaphsonSolver.jl index 08e7bdac..a2c17471 100644 --- a/src/NonlinearSolvers/NewtonRaphsonSolver.jl +++ b/src/NonlinearSolvers/NewtonRaphsonSolver.jl @@ -60,7 +60,7 @@ function _solve_nr!(x,A,b,dx,ns,nls,op) if !done # Update jacobian and solver jacobian!(A, op, x) - numerical_setup!(ns,A) + numerical_setup!(ns,A,x) end end diff --git a/src/NonlinearSolvers/NonlinearSolvers.jl b/src/NonlinearSolvers/NonlinearSolvers.jl index 0ee37902..235769a7 100644 --- a/src/NonlinearSolvers/NonlinearSolvers.jl +++ b/src/NonlinearSolvers/NonlinearSolvers.jl @@ -10,6 +10,7 @@ module NonlinearSolvers using PartitionedArrays using GridapDistributed + using GridapSolvers.SolverInterfaces using GridapSolvers.MultilevelTools using GridapSolvers.SolverInterfaces diff --git a/src/SolverInterfaces/GridapExtras.jl b/src/SolverInterfaces/GridapExtras.jl index a6d5d24f..88a549ee 100644 --- a/src/SolverInterfaces/GridapExtras.jl +++ b/src/SolverInterfaces/GridapExtras.jl @@ -1,11 +1,6 @@ # LinearSolvers that depend on the non-linear solution -""" + function Gridap.Algebra.numerical_setup!(ns::Gridap.Algebra.LinearSolver,A::AbstractMatrix,x::AbstractVector) numerical_setup!(ns,A) end - -function allocate_solver_caches(ns::Gridap.Algebra.LinearSolver,args...;kwargs...) - @abstractmethod -end -""" \ No newline at end of file