tlm_adjoint.hessian_system

Hessian matrix utilities.

Module Contents

class tlm_adjoint.hessian_system.HessianLinearSolver(H, M, *args, **kwargs)

Solver for linear systems involving a Hessian matrix.

Solves the linear system

\[H u = b\]

for \(u\) using PETSc.

Parameters:
HHessian

Defines \(H\).

Mvariable or Sequence[variable, …]

Defines the control and its value.

args, kwargs

Passed to the tlm_adjoint.block_system.LinearSolver constructor.

solve(u, b, *args, **kwargs)

Solve a linear system involving a Hessian matrix.

Solves the linear system

\[H u = b\]

for \(u\).

Parameters:
uvariable or Sequence[variable, …]

The solution \(u\). Also defines the initial guess.

bvariable or Sequence[variable, …]

The conjugate of the right-hand-side \(b\).

args, kwargs

Passed to the tlm_adjoint.block_system.LinearSolver.solve() method.

class tlm_adjoint.hessian_system.HessianEigensolver(H, M, B_action, B_inv_action, *args, **kwargs)

Solver for generalized eigenproblems involving a Hessian matrix.

Solves the eigenproblem

\[H v = \lambda B v\]

using SLEPc.

Parameters:
HHessian

Defines \(H\).

Mvariable or Sequence[variable, …]

Defines the control and its value.

B_actioncallable

Accepts one or more variables as arguments, defining the direction, and returns a variable or a Sequence of variables defining the conjugate of the action of \(B\) on this direction. Arguments should not be modified.

B_inv_actioncallable

Accepts one or more variables as arguments, defining the direction, and returns a variable or a Sequence of variables defining the action of \(B^{-1}\) on the conjugate of this direction. Arguments should not be modified.

args, kwargs

Passed to the Eigensolver constructor.

spectral_approximation_solve(b)

\((H + B)^{-1}\) action approximation.

Computes an approximation for

\[(H + B)^{-1} b\]

where \(H\) and \(B\) define the eigenproblem solved by this HessianEigensolver. The approximation is constructed using a partial eigenspectrum – see HessianEigensolver.spectral_pc_fn().

Parameters:
bvariable or Sequence[variable, …]

The conjugate of the right-hand-side \(b\).

Returns:
variable or tuple[variable, …]

The approximation for the action on \(b\).

spectral_pc_fn()

Construct a partial eigenspectrum preconditioner.

Constructs a matrix preconditioner using a partial eigenspectrum. Specifically for a matrix

\[C = H + B,\]

where \(H\) and \(B\) define the eigenproblem solved by this HessianEigensolver, the approximation is defined via

\[C^{-1} \approx B^{-1} - V \Lambda \left( I + \Lambda \right)^{-1} V^*\]

where

\[H V = B V \Lambda,\]

and where \(\Lambda\) is a diagonal matrix and \(V\) has \(B\)-orthonormal columns, \(V^* B V = I\).

This low rank update approximation for the Hessian matrix inverse is described in

  • Tobin Isaac, Noemi Petra, Georg Stadler, and Omar Ghattas, ‘Scalable and efficient algorithms for the propagation of uncertainty from data through inference to prediction for large-scale problems, with application to flow of the Antarctic ice sheet’, Journal of Computational Physics, 296, pp. 348–368, 2015, doi: 10.1016/j.jcp.2015.04.047

See in particular their equation (20).

Returns:
callable

Suitable for use as the pc_fn argument to HessianLinearSolver.solve().