firedrake.slate.static_condensation package

Submodules

firedrake.slate.static_condensation.hybridization module

class firedrake.slate.static_condensation.hybridization.HybridizationPC[source]

Bases: SCBase

Create a PC context suitable for PETSc.

Matrix free preconditioners should inherit from this class and implement:

backward_substitution(pc, y)[source]

Perform the backwards recovery of eliminated fields.

Parameters:
  • pc – a Preconditioner instance.

  • y – a PETSc vector for placing the resulting fields.

forward_elimination(pc, x)[source]

Perform the forward elimination of fields and provide the reduced right-hand side for the condensed system.

Parameters:
  • pc – a Preconditioner instance.

  • x – a PETSc vector containing the incoming right-hand side.

getSchurComplementBuilder()[source]
initialize(pc)[source]

Set up the problem context. Take the original mixed problem and reformulate the problem as a hybridized mixed system.

A KSP is created for the Lagrange multiplier system.

needs_python_pmat = True

A Slate-based python preconditioner that solves a mixed H(div)-conforming problem using hybridization. Currently, this preconditioner supports the hybridization of the RT and BDM mixed methods of arbitrary degree.

The forward eliminations and backwards reconstructions are performed element-local using the Slate language.

sc_solve(pc)[source]

Solve the condensed linear system for the condensed field.

Parameters:

pc – a Preconditioner instance.

update(pc)[source]

Update by assembling into the operator. No need to reconstruct symbolic objects.

view(pc, viewer=None)[source]

Viewer calls for the various configurable objects in this PC.

class firedrake.slate.static_condensation.hybridization.SchurComplementBuilder(prefix, Atilde, K, KT, pc, vidx, pidx, non_zero_saddle_mat=None)[source]

Bases: object

A Slate-based Schur complement expression builder. The expression is used in the trace system solve and parts of it in the reconstruction calls of the other two variables of the hybridised system. How the Schur complement if constructed, and in particular how the local inverse of the mixed matrix is built, is controlled with PETSc options. All corresponding PETSc options start with hybridization_localsolve. The following option sets are valid together with the usual set of hybridisation options:

{'localsolve': {'ksp_type': 'preonly',
                'pc_type': 'fieldsplit',
                'pc_fieldsplit_type': 'schur'}}

A Schur complement is requested for the mixed matrix inverse which appears inside the Schur complement of the trace system solve. The Schur complements are then nested. For details see defition of build_schur(). No fieldsplit options are set so all local inverses are calculated explicitly.

'localsolve': {'ksp_type': 'preonly',
               'pc_type': 'fieldsplit',
               'pc_fieldsplit_type': 'schur',
               'fieldsplit_1': {'ksp_type': 'default',
                                'pc_type': 'python',
                                'pc_python_type': __name__ + '.DGLaplacian'}}

The inverse of the Schur complement inside the Schur decomposition of the mixed matrix inverse is approximated by a default solver (LU in the matrix-explicit case) which is preconditioned by a user-defined operator, e.g. a DG Laplacian, see build_inner_S_inv(). So \(P_S * S * x = P_S * b\).

'localsolve': {'ksp_type': 'preonly',
                'pc_type': 'fieldsplit',
                'pc_fieldsplit_type': 'schur',
                'fieldsplit_1': {'ksp_type': 'default',
                                'pc_type': 'python',
                                'pc_python_type': __name__ + '.DGLaplacian',
                                'aux_ksp_type': 'preonly'}
                                'aux_pc_type': 'jacobi'}}}}

The inverse of the Schur complement inside the Schur decomposition of the mixed matrix inverse is approximated by a default solver (LU in the matrix-explicit case) which is preconditioned by a user-defined operator, e.g. a DG Laplacian. The inverse of the preconditioning matrix is approximated through the inverse of only the diagonal of the provided operator, see build_Sapprox_inv(). So \(diag(P_S).inv * S * x = diag(P_S).inv * b\).

'localsolve': {'ksp_type': 'preonly',
               'pc_type': 'fieldsplit',
               'pc_fieldsplit_type': 'schur',
               'fieldsplit_0': {'ksp_type': 'default',
                                'pc_type': 'jacobi'}

The inverse of the \(A_{00}\) block of the mixed matrix is approximated by a default solver (LU in the matrix-explicit case) which is preconditioned by the diagonal matrix of \(A_{00}, see :meth:`build_A00_inv\). So \(diag(A_{00}).inv * A_{00} * x = diag(A_{00}).inv * b\).

'localsolve': {'ksp_type': 'preonly',
               'pc_type': 'fieldsplit',
               'pc_fieldsplit_type': 'None',
               'fieldsplit_0':  ...
               'fieldsplit_1':  ...

All the options for fieldsplit_ are still valid if 'pc_fieldsplit_type': 'None'. In this case the mixed matrix inverse which appears inside the Schur complement of the trace system solve is calculated explicitly, but the local inverses of \(A_{00}\) and the Schur complement in the reconstructions calls are still treated according to the options in fieldsplit_.

build_A00_inv()[source]

Calculates the inverse of \(A_{00}\), the (0,0)-block of the mixed matrix Atilde. The inverse is potentially approximated through a solve which is potentially preconditioned with jacobi.

build_Sapprox_inv()[source]

Calculates the inverse of preconditioner to the Schur complement, which can be either the schur complement approximation provided by the user or jacobi. The inverse is potentially approximated through a solve which is potentially preconditioned with jacobi.

build_inner_S()[source]

Build the inner Schur complement.

build_inner_S_inv()[source]

Calculates the inverse of the schur complement. The inverse is potentially approximated through a solve which is potentially preconditioned with the preconditioner P.

build_schur(rhs, non_zero_saddle_rhs=None)[source]

The Schur complement in the operators of the trace solve contains the inverse on a mixed system. Users may want this inverse to be treated with another Schur complement.

Let the mixed matrix Atilde be called A here. Then, if a nested schur complement is requested, the inverse of Atilde is rewritten with help of a a Schur decomposition as follows.

A.inv=[[I, -A00.inv * A01] * [[A00.inv, 0    ] * [[I,            0]
      [0,  I             ]]  [0,        S.inv]]  [-A10* A00.inv, I]]
      --------------------   -----------------   ------------------
             block1                block2              block3
with the (inner) schur complement S = A11 - A10 * A00.inv * A01
inv(A, P, prec, preonly=False)[source]

Calculates the inverse of an operator A. The inverse is potentially approximated through a solve which is potentially preconditioned with the preconditioner P if prec is True. The inverse of A may be just approximated with the inverse of P if prec and replace.

retrieve_user_S_approx(pc, usercode)[source]

Retrieve a user-defined :class:firedrake.preconditioners.AuxiliaryOperator from the PETSc Options, which is an approximation to the Schur complement and its inverse is used to precondition the local solve in the reconstruction calls (e.g.).

firedrake.slate.static_condensation.la_utils module

class firedrake.slate.static_condensation.la_utils.LAContext(lhs, rhs, field_idx)

Bases: tuple

Context information for systems of equations after applying algebraic transformation via Slate-supported operations. This object provides the symbolic expressions for the transformed linear system of equations.

Parameters:
  • lhs – The resulting expression for the transformed left-hand side matrix.

  • rhs – The resulting expression for the transformed right-hand side vector.

  • field_idx – An integer or iterable of integers (if the system is mixed) denoting which field(s) the resulting solution is defined on.

Create new instance of LAContext(lhs, rhs, field_idx)

field_idx

Alias for field number 2

lhs

Alias for field number 0

rhs

Alias for field number 1

class firedrake.slate.static_condensation.la_utils.SchurComplementBuilder(prefix, Atilde, K, KT, pc, vidx, pidx, non_zero_saddle_mat=None)[source]

Bases: object

A Slate-based Schur complement expression builder. The expression is used in the trace system solve and parts of it in the reconstruction calls of the other two variables of the hybridised system. How the Schur complement if constructed, and in particular how the local inverse of the mixed matrix is built, is controlled with PETSc options. All corresponding PETSc options start with hybridization_localsolve. The following option sets are valid together with the usual set of hybridisation options:

{'localsolve': {'ksp_type': 'preonly',
                'pc_type': 'fieldsplit',
                'pc_fieldsplit_type': 'schur'}}

A Schur complement is requested for the mixed matrix inverse which appears inside the Schur complement of the trace system solve. The Schur complements are then nested. For details see defition of build_schur(). No fieldsplit options are set so all local inverses are calculated explicitly.

'localsolve': {'ksp_type': 'preonly',
               'pc_type': 'fieldsplit',
               'pc_fieldsplit_type': 'schur',
               'fieldsplit_1': {'ksp_type': 'default',
                                'pc_type': 'python',
                                'pc_python_type': __name__ + '.DGLaplacian'}}

The inverse of the Schur complement inside the Schur decomposition of the mixed matrix inverse is approximated by a default solver (LU in the matrix-explicit case) which is preconditioned by a user-defined operator, e.g. a DG Laplacian, see build_inner_S_inv(). So \(P_S * S * x = P_S * b\).

'localsolve': {'ksp_type': 'preonly',
                'pc_type': 'fieldsplit',
                'pc_fieldsplit_type': 'schur',
                'fieldsplit_1': {'ksp_type': 'default',
                                'pc_type': 'python',
                                'pc_python_type': __name__ + '.DGLaplacian',
                                'aux_ksp_type': 'preonly'}
                                'aux_pc_type': 'jacobi'}}}}

The inverse of the Schur complement inside the Schur decomposition of the mixed matrix inverse is approximated by a default solver (LU in the matrix-explicit case) which is preconditioned by a user-defined operator, e.g. a DG Laplacian. The inverse of the preconditioning matrix is approximated through the inverse of only the diagonal of the provided operator, see build_Sapprox_inv(). So \(diag(P_S).inv * S * x = diag(P_S).inv * b\).

'localsolve': {'ksp_type': 'preonly',
               'pc_type': 'fieldsplit',
               'pc_fieldsplit_type': 'schur',
               'fieldsplit_0': {'ksp_type': 'default',
                                'pc_type': 'jacobi'}

The inverse of the \(A_{00}\) block of the mixed matrix is approximated by a default solver (LU in the matrix-explicit case) which is preconditioned by the diagonal matrix of \(A_{00}, see :meth:`build_A00_inv\). So \(diag(A_{00}).inv * A_{00} * x = diag(A_{00}).inv * b\).

'localsolve': {'ksp_type': 'preonly',
               'pc_type': 'fieldsplit',
               'pc_fieldsplit_type': 'None',
               'fieldsplit_0':  ...
               'fieldsplit_1':  ...

All the options for fieldsplit_ are still valid if 'pc_fieldsplit_type': 'None'. In this case the mixed matrix inverse which appears inside the Schur complement of the trace system solve is calculated explicitly, but the local inverses of \(A_{00}\) and the Schur complement in the reconstructions calls are still treated according to the options in fieldsplit_.

build_A00_inv()[source]

Calculates the inverse of \(A_{00}\), the (0,0)-block of the mixed matrix Atilde. The inverse is potentially approximated through a solve which is potentially preconditioned with jacobi.

build_Sapprox_inv()[source]

Calculates the inverse of preconditioner to the Schur complement, which can be either the schur complement approximation provided by the user or jacobi. The inverse is potentially approximated through a solve which is potentially preconditioned with jacobi.

build_inner_S()[source]

Build the inner Schur complement.

build_inner_S_inv()[source]

Calculates the inverse of the schur complement. The inverse is potentially approximated through a solve which is potentially preconditioned with the preconditioner P.

build_schur(rhs, non_zero_saddle_rhs=None)[source]

The Schur complement in the operators of the trace solve contains the inverse on a mixed system. Users may want this inverse to be treated with another Schur complement.

Let the mixed matrix Atilde be called A here. Then, if a nested schur complement is requested, the inverse of Atilde is rewritten with help of a a Schur decomposition as follows.

A.inv=[[I, -A00.inv * A01] * [[A00.inv, 0    ] * [[I,            0]
      [0,  I             ]]  [0,        S.inv]]  [-A10* A00.inv, I]]
      --------------------   -----------------   ------------------
             block1                block2              block3
with the (inner) schur complement S = A11 - A10 * A00.inv * A01
inv(A, P, prec, preonly=False)[source]

Calculates the inverse of an operator A. The inverse is potentially approximated through a solve which is potentially preconditioned with the preconditioner P if prec is True. The inverse of A may be just approximated with the inverse of P if prec and replace.

retrieve_user_S_approx(pc, usercode)[source]

Retrieve a user-defined :class:firedrake.preconditioners.AuxiliaryOperator from the PETSc Options, which is an approximation to the Schur complement and its inverse is used to precondition the local solve in the reconstruction calls (e.g.).

firedrake.slate.static_condensation.la_utils.backward_solve(A, b, x, schur_builder, reconstruct_fields)[source]

Returns a sequence of linear algebra contexts containing Slate expressions for backwards substitution.

Parameters:
  • A – a slate.Tensor corresponding to the mixed UFL operator.

  • b – a firedrake.Function corresponding to the right-hand side.

  • x – a firedrake.Function corresponding to the solution.

  • schur_builder – a SchurComplementBuilder

  • reconstruct_fields – a tuple of indices denoting which fields to reconstruct.

Returns:

a list of LAContext for the reconstruction

firedrake.slate.static_condensation.la_utils.condense_and_forward_eliminate(A, b, elim_fields, prefix, pc)[source]

Returns Slate expressions for the operator and right-hand side vector after eliminating specified unknowns.

Parameters:
  • A – a slate.Tensor corresponding to the mixed UFL operator.

  • b – a firedrake.Function corresponding to the right-hand side.

  • elim_fields – a tuple of indices denoting which fields to eliminate.

  • prefix – an option prefix for the condensed field.

  • pc – a Preconditioner instance.

Returns:

a tuple of LAContext and SchurComplementBuilder

firedrake.slate.static_condensation.sc_base module

class firedrake.slate.static_condensation.sc_base.SCBase[source]

Bases: PCBase

A general-purpose base class for static condensation interfaces.

Create a PC context suitable for PETSc.

Matrix free preconditioners should inherit from this class and implement:

apply(pc, x, y)[source]

Applies the static condensation preconditioner.

Parameters:
  • pc – a Preconditioner instance.

  • x – A PETSc vector containing the incoming right-hand side.

  • y – A PETSc vector for the result.

applyTranspose(pc, x, y)[source]

Apply the transpose of the preconditioner.

abstract backward_substitution(pc, y)[source]

Perform the backwards recovery of eliminated fields.

Parameters:
  • pc – a Preconditioner instance.

  • y – a PETSc vector for placing the resulting fields.

abstract forward_elimination(pc, x)[source]

Perform the forward elimination of fields and provide the reduced right-hand side for the condensed system.

Parameters:
  • pc – a Preconditioner instance.

  • x – a PETSc vector containing the incoming right-hand side.

abstract sc_solve(pc)[source]

Solve the condensed linear system for the condensed field.

Parameters:

pc – a Preconditioner instance.

firedrake.slate.static_condensation.scpc module

class firedrake.slate.static_condensation.scpc.SCPC[source]

Bases: SCBase

Create a PC context suitable for PETSc.

Matrix free preconditioners should inherit from this class and implement:

backward_substitution(pc, y)[source]

Perform the backwards recovery of eliminated fields.

Parameters:
  • pc – a Preconditioner instance.

  • y – a PETSc vector for placing the resulting fields.

condensed_system(A, rhs, elim_fields, prefix, pc)[source]

Forms the condensed linear system by eliminating specified unknowns.

Parameters:
  • A – A Slate Tensor containing the mixed bilinear form.

  • rhs – A firedrake function for the right-hand side.

  • elim_fields – An iterable of field indices to eliminate.

  • prefix – an option prefix for the condensed field.

  • pc – a Preconditioner instance.

forward_elimination(pc, x)[source]

Perform the forward elimination of fields and provide the reduced right-hand side for the condensed system.

Parameters:
  • pc – a Preconditioner instance.

  • x – a PETSc vector containing the incoming right-hand side.

initialize(pc)[source]

Set up the problem context. This takes the incoming three-field system and constructs the static condensation operators using Slate expressions.

A KSP is created for the reduced system. The eliminated variables are recovered via back-substitution.

local_solver_calls(A, rhs, x, elim_fields, schur_builder)[source]

Provides solver callbacks for inverting local operators and reconstructing eliminated fields.

Parameters:
  • A – A Slate Tensor containing the mixed bilinear form.

  • rhs – A firedrake function for the right-hand side.

  • x – A firedrake function for the solution.

  • elim_fields – An iterable of eliminated field indices to recover.

  • schur_builder – a SchurComplementBuilder.

needs_python_pmat = True

A Slate-based python preconditioner implementation of static condensation for problems with up to three fields.

sc_solve(pc)[source]

Solve the condensed linear system for the condensed field.

Parameters:

pc – a Preconditioner instance.

update(pc)[source]

Update by assembling into the KSP operator. No need to reconstruct symbolic objects.

view(pc, viewer=None)[source]

Viewer calls for the various configurable objects in this PC.

Module contents