GridapSolvers.LinearSolvers

Krylov solvers

GridapSolvers.LinearSolvers.CGSolverType
struct CGSolver <: LinearSolver
  ...
end

CGSolver(Pl;maxiter=1000,atol=1e-12,rtol=1.e-6,flexible=false,verbose=0,name="CG")

Left-Preconditioned Conjugate Gradient solver.

source
GridapSolvers.LinearSolvers.MINRESSolverType
struct MINRESSolver <: LinearSolver 
  ...
end

MINRESSolver(m;Pl=nothing,maxiter=100,atol=1e-12,rtol=1.e-6,verbose=false,name="MINRES")

MINRES solver, with optional left preconditioners Pl. The preconditioner must be symmetric and positive definite.

source
GridapSolvers.LinearSolvers.GMRESSolverType
struct GMRESSolver <: LinearSolver 
  ...
end

GMRESSolver(m;Pr=nothing,Pl=nothing,restart=false,m_add=1,maxiter=100,atol=1e-12,rtol=1.e-6,verbose=false,name="GMRES")

GMRES solver, with optional right and left preconditioners Pr and Pl.

The solver starts by allocating a basis of size m. Then:

  • If restart=true, the basis size is fixed and restarted every m iterations.
  • If restart=false, the basis size is allowed to increase. When full, the solver allocates m_add new basis vectors.
source
GridapSolvers.LinearSolvers.FGMRESSolverType
struct FGMRESSolver <: LinearSolver 
  ...
end

FGMRESSolver(m,Pr;Pl=nothing,restart=false,m_add=1,maxiter=100,atol=1e-12,rtol=1.e-6,verbose=false,name="FGMRES")

Flexible GMRES solver, with right-preconditioner Pr and optional left-preconditioner Pl.

The solver starts by allocating a basis of size m. Then:

  • If restart=true, the basis size is fixed and restarted every m iterations.
  • If restart=false, the basis size is allowed to increase. When full, the solver allocates m_add new basis vectors at a time.
source

Smoothers

Given a linear system $Ax = b$, a smoother is an operator S that takes an iterative solution $x_k$ and its residual $r_k = b - A x_k$, and modifies them in place

\[ S : (x_k,r_k) \rightarrow (x_{k+1},r_{k+1})\]

such that $|r_{k+1}| < |r_k|$.

GridapSolvers.LinearSolvers.RichardsonSmootherType
struct RichardsonSmoother{A} <: LinearSolver
  M     :: A
  niter :: Int64
  ω     :: Float64
end

Iterative Richardson smoother. Given a solution x and a residual r, performs niter Richardson iterations with damping parameter ω using the linear solver M. A Richardson iteration is given by:

dx = ω * inv(M) * r
x  = x + dx
r  = r - A * dx

Updates both the solution x and the residual r in place.

source

Preconditioners

Given a linear system $Ax = b$, a preconditioner is an operator that takes an iterative residual $r_k$ and returns a correction $dx_k$.

GridapSolvers.LinearSolvers.GMGLinearSolverFunction
GMGLinearSolver(
  mh::ModelHierarchy,
  matrices::AbstractArray{<:AbstractMatrix},
  prolongations,
  restrictions;
  pre_smoothers   = Fill(RichardsonSmoother(JacobiLinearSolver(),10),num_levels(mh)-1),
  post_smoothers  = pre_smoothers,
  coarsest_solver = LUSolver(),
  mode::Symbol    = :preconditioner,
  maxiter = 100, atol = 1.0e-14, rtol = 1.0e-08, verbose = false,
)

Creates an instance of GMGLinearSolverFromMatrices from the underlying model hierarchy, the system matrices at each level and the transfer operators and smoothers at each level except the coarsest.

The solver has two modes of operation, defined by the kwarg mode:

  • :solver: The GMG solver takes a rhs b and returns a solution x.
  • :preconditioner: The GMG solver takes a residual r and returns a correction dx.
source
GMGLinearSolver(
  mh::ModelHierarchy,
  trials::FESpaceHierarchy,
  tests::FESpaceHierarchy,
  biforms::AbstractArray{<:Function},
  interp,
  restrict;
  pre_smoothers   = Fill(RichardsonSmoother(JacobiLinearSolver(),10),num_levels(mh)-1),
  post_smoothers  = pre_smoothers,
  coarsest_solver = Gridap.Algebra.LUSolver(),
  mode::Symbol    = :preconditioner,
  is_nonlinear    = false,
  maxiter = 100, atol = 1.0e-14, rtol = 1.0e-08, verbose = false,
)

Creates an instance of GMGLinearSolverFromMatrices from the underlying model hierarchy, the trial and test FEspace hierarchies, the weakform lhs at each level and the transfer operators and smoothers at each level except the coarsest.

The solver has two modes of operation, defined by the kwarg mode:

  • :solver: The GMG solver takes a rhs b and returns a solution x.
  • :preconditioner: The GMG solver takes a residual r and returns a correction dx.
source

Wrappers

PETSc

Building on top of GridapPETSc.jl, GridapSolvers provides specific solvers for some particularly complex PDEs:

GridapSolvers.LinearSolvers.CachedPETScNSType
struct CachedPETScNS <: NumericalSetup
  ...
end

Wrapper around a PETSc NumericalSetup, providing highly efficiend reusable caches:

When converting julia vectors/PVectors to PETSc vectors, we purposely create aliasing of the vector values. This means we can avoid copying data from one to another before solving, but we need to be careful about it.

This structure takes care of this, and makes sure you do not attempt to solve the system with julia vectors that are not the ones you used to create the solver cache.

source
GridapSolvers.LinearSolvers.CachedPETScNSMethod
function CachedPETScNS(ns::PETScLinearSolverNS,x::AbstractVector,b::AbstractVector)

Create a new instance of CachedPETScNS from its underlying properties. Once this structure is created, you can only solve the system with the same vectors you used to create it.

source

IterativeSolvers.jl

GridapSolvers provides wrappers for some iterative solvers from the package IterativeSolvers.jl: