firedrake.slate.static_condensation package

Submodules

firedrake.slate.static_condensation.hybridization module

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

Bases: firedrake.slate.static_condensation.sc_base.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.

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.

firedrake.slate.static_condensation.la_utils module

class firedrake.slate.static_condensation.la_utils.LAContext

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)

property field_idx

Alias for field number 2

property lhs

Alias for field number 0

property rhs

Alias for field number 1

firedrake.slate.static_condensation.la_utils.backward_solve(A, b, x, 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.

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

firedrake.slate.static_condensation.la_utils.condense_and_forward_eliminate(A, b, elim_fields)[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.

firedrake.slate.static_condensation.sc_base module

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

Bases: firedrake.preconditioners.base.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.sc_base.timed_region()

Log.Event(type cls, name, klass=None)

firedrake.slate.static_condensation.scpc module

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

Bases: firedrake.slate.static_condensation.sc_base.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)[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.

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)[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.

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