tlm_adjoint.block_system

Mixed space linear algebra utilities.

Given a linear problem with a potentially singular matrix \(A\)

\[A u = b,\]

a LinearSolver instead solves the linear problem

\[\left[ (I - M U (U^* M U)^{-1} U^*) A (I - V (V^* C V)^{-1} V^* C) + M U S V^* C \right] u = (I - M U (U^* M U)^{-1} U^*) b.\]

Here

  • \(U\) is a full rank matrix whose columns span the left nullspace for a modified system matrix \(\tilde{A}\).

  • \(V\) is a full rank matrix with the same number of columns as \(U\), whose columns span the nullspace for \(\tilde{A}\).

  • \(V^* C V\) and \(S\) are invertible matrices.

  • \(M\) is a Hermitian positive definite matrix.

Here the left nullspace for a matrix is defined to be the nullspace for its Hermitian transpose, and the modified system matrix \(\tilde{A}\) is defined

\[\tilde{A} = (I - M U (U^* M U)^{-1} U^*) A (I - V (V^* C V)^{-1} V^* C).\]

This has two primary use cases:

  1. Where a matrix \(A\) and right-hand-side \(b\) are constructed via finite element assembly on superspaces of the test space and trial space. The typical example is in the application of homogeneous essential Dirichlet boundary conditions.

  2. Where the matrix \(A\) is singular and \(b\) is orthogonal to the left nullspace of \(A\). Typically one would then choose \(U\) and \(V\) so that their columns respectively span the left nullspace and nullspace of \(A\), and the LinearSolver then seeks a solution to the original problem subject to the linear constraints \(V^* C u = 0\).

Spaces are defined via backend spaces or TypedSpace objects, and Sequence objects containing backend spaces, TypedSpace objects, or similar Sequence objects. Similarly variables are defined via backend variables, or Sequence objects containing backend variables, or similar Sequence objects. This defines a basic tree structure which is useful e.g. when defining block matrices in terms of sub-block matrices.

Elements of the tree are accessed in a consistent order using a depth first search. Hence e.g.

((u_0, u_1), u_2)

and

(u_0, u_1, u_2)

where u_0, u_1, and u_2 are backend variables, are both valid representations of a mixed space solution.

Module Contents

class tlm_adjoint.block_system.TypedSpace(space, *, space_type=None)

A space with an associated space type.

Parameters:
spacespace

The backend space.

space_typestr

The space type.

Attributes:
commcommunicator

The communicator associated with the space.

spacespace

The backend space.

space_typestr

The space type.

new()

Create a new variable in the space.

Returns:
variable

The new variable.

class tlm_adjoint.block_system.MixedSpace(spaces)

Used to map between different versions of a mixed space.

This class defines two representations for the space:

  1. As a ‘split space’: A tree defining the mixed space. Stored using TypedSpace and tuple objects, each corresponding to a node in the tree. TypedSpace objects correspond to leaf nodes, and tuple objects to other nodes in the tree.

  2. As a ‘flattened space’: A Sequence containing leaf nodes of the split space with an ordering determined using a depth first search.

Provides methods to allow data to be copied to and from a compatible petsc4py.PETSc.Vec. This allows, for example, the construction (with Firedrake):

u_0 = Function(space_0, name='u_0')
u_1 = Function(space_1, name='u_1')
u_2 = Function(space_2, name='u_2')

mixed_space = MixedSpace(((space_0, space_1), space_2))

and then data can be copied to a compatible petsc4py.PETSc.Vec via

mixed_space.to_petsc(u_petsc, ((u_0, u_1), u_2))

and from a compatible petsc4py.PETSc.Vec via

mixed_space.from_petsc(u_petsc, ((u_0, u_1), u_2))
Parameters:

spaces – Defines the split space.

property split_space

The split space representation.

property flattened_space

The flattened space representation.

tuple_sub(u)
Parameters:

u – An Iterable.

Returns:

A tuple storing elements in u using the tree structure of the split space.

new()
Returns:

A new element in the split space.

class tlm_adjoint.block_system.Nullspace

Represents a nullspace and left nullspace for a square matrix.

abstract apply_nullspace_transformation_lhs_right(x)

Apply the nullspace transformation associated with a matrix action on \(x\),

\[x \rightarrow (I - V (V^* C V)^{-1} V^* C) x.\]
Parameters:

x – Defines \(x\).

abstract apply_nullspace_transformation_lhs_left(y)

Apply the left nullspace transformation associated with a matrix action,

\[y \rightarrow (I - M U (U^* M U)^{-1} U^*) y.\]
Parameters:

y – Defines \(y\).

abstract constraint_correct_lhs(x, y)

Add the linear constraint term to \(y\),

\[y \rightarrow y + M U S V^* C x.\]
Parameters:

  • x – Defines \(x\).

  • y – Defines \(y\).

abstract pc_constraint_correct_soln(u, b)

Add the preconditioner linear constraint term to \(u\),

\[u \rightarrow u + V \tilde{S}^{-1} U^* b,\]

with

\[\tilde{S}^{-1} = \left( V^* C V \right)^{-1} S^{-1} \left( U^* M U \right)^{-1}.\]
Parameters:

  • u – Defines \(u\).

  • b – Defines \(b\).

correct_soln(x)

Correct the linear system solution so that it is orthogonal to space spanned by the columns of \(V\).

Parameters:

x – The linear system solution, to be corrected.

pre_mult_correct_lhs(x)

Apply the pre-left-multiplication nullspace transformation.

Parameters:

x – Defines the vector on which the matrix action is computed.

post_mult_correct_lhs(x, y)

Apply the post-left-multiplication nullspace transformation, and add the linear constraint term.

Parameters:

  • x – Defines the vector on which the matrix action is computed, and used to add the linear constraint term. If None is supplied then the linear constraint term is not added.

  • y – Defines the result of the matrix action on x.

correct_rhs(b)

Correct the linear system right-hand-side so that it is orthogonal to the space spanned by the columns of \(U\).

Parameters:

b – The linear system right-hand-side, to be corrected.

pc_pre_mult_correct(b)

Apply the pre-preconditioner-application nullspace transformation.

Parameters:

b – Defines the vector on which the preconditioner action is computed.

pc_post_mult_correct(u, b)

Apply the post-preconditioner-application left nullspace transformation, and add the linear constraint term.

Parameters:

  • u – Defines the result of the preconditioner action on b.

  • b – Defines the vector on which the preconditioner action is computed, and used to add the linear constraint term. If None is supplied then the linear constraint term is not added.

class tlm_adjoint.block_system.NoneNullspace

An empty nullspace and left nullspace.

apply_nullspace_transformation_lhs_right(x)

Apply the nullspace transformation associated with a matrix action on \(x\),

\[x \rightarrow (I - V (V^* C V)^{-1} V^* C) x.\]
Parameters:

x – Defines \(x\).

apply_nullspace_transformation_lhs_left(y)

Apply the left nullspace transformation associated with a matrix action,

\[y \rightarrow (I - M U (U^* M U)^{-1} U^*) y.\]
Parameters:

y – Defines \(y\).

constraint_correct_lhs(x, y)

Add the linear constraint term to \(y\),

\[y \rightarrow y + M U S V^* C x.\]
Parameters:

  • x – Defines \(x\).

  • y – Defines \(y\).

pc_constraint_correct_soln(u, b)

Add the preconditioner linear constraint term to \(u\),

\[u \rightarrow u + V \tilde{S}^{-1} U^* b,\]

with

\[\tilde{S}^{-1} = \left( V^* C V \right)^{-1} S^{-1} \left( U^* M U \right)^{-1}.\]
Parameters:

  • u – Defines \(u\).

  • b – Defines \(b\).

class tlm_adjoint.block_system.BlockNullspace(nullspaces)

Nullspaces for a square BlockMatrix.

Parameters:

nullspaces – A Nullspace or a Sequence of Nullspace objects defining the nullspace. None indicates a NoneNullspace.

apply_nullspace_transformation_lhs_right(x)

Apply the nullspace transformation associated with a matrix action on \(x\),

\[x \rightarrow (I - V (V^* C V)^{-1} V^* C) x.\]
Parameters:

x – Defines \(x\).

apply_nullspace_transformation_lhs_left(y)

Apply the left nullspace transformation associated with a matrix action,

\[y \rightarrow (I - M U (U^* M U)^{-1} U^*) y.\]
Parameters:

y – Defines \(y\).

constraint_correct_lhs(x, y)

Add the linear constraint term to \(y\),

\[y \rightarrow y + M U S V^* C x.\]
Parameters:

  • x – Defines \(x\).

  • y – Defines \(y\).

pc_constraint_correct_soln(u, b)

Add the preconditioner linear constraint term to \(u\),

\[u \rightarrow u + V \tilde{S}^{-1} U^* b,\]

with

\[\tilde{S}^{-1} = \left( V^* C V \right)^{-1} S^{-1} \left( U^* M U \right)^{-1}.\]
Parameters:

  • u – Defines \(u\).

  • b – Defines \(b\).

class tlm_adjoint.block_system.Matrix(arg_space, action_space)

Represents a matrix defining a mapping \(A\) mapping \(V \rightarrow W\).

Parameters:

  • arg_space – Defines the space \(V\).

  • action_space – Defines the space \(W\).

property arg_space

The space defining \(V\).

property action_space

The space defining \(W\).

mult(x, y)

Compute \(y = A x\).

Parameters:

  • x – Defines \(x\). Should not be modified.

  • y – Defines \(y\).

mult_add(x, y)

Add \(A x\) to \(y\).

Parameters:

  • x – Defines \(x\). Should not be modified.

  • y – Defines \(y\).

class tlm_adjoint.block_system.MatrixFreeMatrix(arg_space, action_space, mult)

Represents a matrix defining a mapping \(A\) mapping \(V \rightarrow W\).

Parameters:

  • arg_space – Defines the space \(V\).

  • action_space – Defines the space \(W\).

mult(x, y)

Compute \(y = A x\).

Parameters:

  • x – Defines \(x\). Should not be modified.

  • y – Defines \(y\).

class tlm_adjoint.block_system.BlockMatrix(arg_space, action_space, blocks=None)

A matrix defining a mapping \(A\) mapping \(V \rightarrow W\), where \(V\) and \(W\) are defined by mixed spaces.

Parameters:

  • arg_space – Defines the space \(V\).

  • action_space – Defines the space \(W\).

  • block – A Mapping defining the blocks of the matrix. Items are ((i, j), block) where the block in the i th and j th column is defined by block. Each block is a tlm_adjoint.block_system.Matrix, or None to indicate a zero block.

mult_add(x, y)

Add \(A x\) to \(y\).

Parameters:

  • x – Defines \(x\). Should not be modified.

  • y – Defines \(y\).

class tlm_adjoint.block_system.LinearSolver(A, *, nullspace=None, solver_parameters=None, pc_fn=None, comm=None)

Solver for a linear system

\[A u = b,\]

using PETSc.

Parameters:

  • A – A tlm_adjoint.block_system.Matrix defining \(A\).

  • nullspace – A Nullspace or a Sequence of Nullspace objects defining the nullspace and left nullspace of \(A\). None indicates a NoneNullspace.

  • solver_parameters – A Mapping defining Krylov solver parameters.

  • pc_fn –

    Defines the application of a preconditioner. A callable

    def pc_fn(u, b):

    The preconditioner is applied to b, and the result stored in u. Defaults to an identity.

  • comm – Communicator.

property ksp

The petsc4py.PETSc.KSP used to solve the linear problem.

solve(u, b, *, correct_initial_guess=True, correct_solution=True)

Solve the linear system.

Parameters:

  • u – Defines the solution \(u\).

  • b – Defines the right-hand-side \(b\).

  • correct_initial_guess – Whether to apply a nullspace correction to the initial guess.

  • correct_solution – Whether to apply a nullspace correction to the solution.

class tlm_adjoint.block_system.Eigensolver(A, B=None, *, B_inv=None, nullspace=None, solver_parameters=None, comm=None)

Solver for an eigenproblem

\[A v = \lambda B v\]

using SLEPc.

Parameters:
property eps

The slepc4py.SLEPc.EPS used to solve the eigenproblem.

is_hermitian_and_positive()
Returns:

Whether the eigenproblem is Hermitian with positive semi-definite \(B\).

solve()

Solve the eigenproblem.

eigenvalues()

Return converged eigenvalues.

Returns:

A numpy.ndarray containing eigenvalues.

eigenpairs()

Return converged eigenpairs.

Returns:

A tuple (Lam, V), where

  • Lam is a numpy.ndarray containing eigenvalues.

  • V is a Sequence of length two tuple objects, (v_r, v_i), defining corresponding eigenvectors. Each eigenvectors is defined by v_r + 1.0j v_i. v_i is None, to indicate a value of zero, for Hermitian eigenvalue problems or with a complex PETSc build.

B_orthonormality_test()

Test \(B\) orthonormality of the eigenvectors for a Hermitian eigenvalue problem.

Returns:

\(\left| V^* B V - I \right|_\infty\) where \(V\) is the matrix whose columns are the eigenvectors.

class tlm_adjoint.block_system.MatrixFunctionSolver(A, *, nullspace=None, solver_parameters=None, comm=None)

Matrix function action evaluation

\[v = f ( A ) u\]

using SLEPc.

Parameters:
property mfn

The slepc4py.SLEPc.MFN used to compute the matrix function action.

solve(u, v)

Compute the matrix function action.

Parameters:

  • u – Defines the argument \(u\).

  • v – Defines the result \(v\).