What is the best way to set up block preconditioner with inner solves

1

I would like to use Fenics in python to assemble and solve a block system with a block diagonal preconditioner where the blocks in the preconditioner have their own inner solvers. Let’s take a mixed Poisson equation as an example with 2\times 2 block structure. (For my application I am in fact interested in general n\times n block structure.) After discretization we end up with a block matrix

A = \begin{bmatrix} M & B^t \\ B & 0\\ \end{bmatrix}.

Let’s say that we have a reasonable preconditioner for A given by

P = \begin{bmatrix} M & \\ & S\\ \end{bmatrix}.

To use this preconditioner I have the following code for setup and solve

PETScOptions.set("ksp_type", "minres")  #fgmres with inner solvers
PETScOptions.set("pc_type", "lu")
PETScOptions.set("pc_factor_mat_solver_type", "mumps")
PETScOptions.set("ksp_converged_reason")

solver = PETScKrylovSolver()
solver.set_operators(A, P)
solver.set_from_options()
x = Function(W)
solver.solve(x.vector(), b)

I would like to use different inner solvers to solve the systems M and S in the preconditioner. Typically using preconditioned conjugate gradients. I know that this can be done in PETSc using PCFIELDSPLIT which allows to use different KSPs with their inner PCs as solvers of the diagonal blocks. Does Fenics provide some interface to this code? Is there some recommended way of dealing with block systems?
The options I have found so far are petsc4py and cbc.block but neither has a nice documentation. I’d be grateful for any recommendation on what would be the right tool/approach to use.

2

See, for example, the Stokes Taylor-Hood demo.

2 Likes