Fieldsplit with Static Condensation in inner blocks

[volpatto]

Hi, Firedrakers.

I have a question, and maybe it is a basic question. I have a code solving a double permeability Darcy's problem. I will skip the details since it is not required to understand my question (but feel free to contact me if you have interest). This problem has 4-fields: 2 pressures and 2 velocities, each one associated with a scale (macro or micro) in the porous medium.

The thing is: I use a hybrid method to solve this problem, introducing a Lagrange multiplier in each scale (similar to HDG). But, so far, I was unable to solve (to be honest, I forever post-pone it) it with a single "condensed system". So, in the past, with the help of Thomas Gibson, we devised a strategy using the fieldsplit facility. We split the system according to the scales (macro and micro, and their associated unknowns), and then apply static condensation in each bock separately. Everything is put inside a Krylov method (it can be GMRES, for instance), which iterates until the convergence is achieved.

Here is the code that configures this setup:

solver_parameters = {
    "snes_type": "ksponly",
    "pmat_type": "matfree",
    "ksp_type": "gmres",
    "ksp_view": None,
    "ksp_monitor_true_residual": None,
    "snes_monitor": None,
    "ksp_rtol": 1.0e-12,
    "ksp_atol": 1.0e-12,
    "pc_type": "fieldsplit",
    "pc_fieldsplit_0_fields": "0,1,2",  # macro
    "pc_fieldsplit_1_fields": "3,4,5",  # micro
    "fieldsplit_0": {
        "pmat_type": "matfree",
        "ksp_type": "preonly",
        "pc_type": "python",
        "pc_python_type": "firedrake.SCPC",
        "pc_sc_eliminate_fields": "0, 1",
        "condensed_field": {
            "ksp_type": "preonly",
            "pc_type": "lu",
            "pc_factor_mat_solver_type": "mumps",
        },
    },
    "fieldsplit_1": {
        "pmat_type": "matfree",
        "ksp_type": "preonly",
        "pc_type": "python",
        "pc_python_type": "firedrake.SCPC",
        "pc_sc_eliminate_fields": "0, 1",
        "condensed_field": {
            "ksp_type": "preonly",
            "pc_type": "lu",
            "pc_factor_mat_solver_type": "mumps",
        },
    },
}

However, I am not sure how the iteration occurs. We solve both condensed systems, retrieve the pressures and velocities block-wisely, push them back into the "monolithic" system (the whole system), and what? How is everything used in the KSP? Just "Picard iterations" would be enough? Or the PC with fieldsplit is just solved once and then the KSP method performs on the monolithic system? I am asking because I have to provide more information in a submission (well, shame on me).

Thanks in advance. Cheers!
Diego.

[wence-]

Can you show the whole of the ksp_view output? I think by default -pc_type fieldsplit gives you a multiplicative PC (ie block forward gauss-seidel), so let's continue with that.

You run GMRES on the outside, split into macro and micro fields, these are both considered the same, so let's just consider the macro fields.

Each outer iteration approximates the inverse with (see https://petsc.org/release/docs/manual/ksp/#sec-block-matrices):

$$ P^{-1} := \begin{bmatrix} I & 0 \\ 0 & A_{11}^{-1} \end{bmatrix} \begin{bmatrix} I & 0 \\ -A_{10} & I \end{bmatrix} \begin{bmatrix} A_{00}^{-1} & 0 \\ 0 & I \end{bmatrix} $$

Where the block inverses are applied using SCPC. That is itself a fieldsplit-analogous preconditioner that performs a (user-guided) element-wise schur complement elimination (in this case, eliminating the first two fields onto the third one) and then assembles that third field into a global system that is inverted with mumps.

So, writing the macro matrix blockwise as:

$$ \begin{bmatrix} C & K \\ K^T & D \end{bmatrix} $$

SCPC, assembles (elementwise using slate) $S := D - K^T C^{-1} K$ (I might have my $K$ and $K^T$ the wrong way round) into an AIJ matrix, and then that is inverted using sparse LU (mumps), followed by back-substitution to reconstruct the fields corresponding to $C$.

Does that help?