firedrake.adjoint package

Submodules

firedrake.adjoint.ensemble_reduced_functional module

class firedrake.adjoint.ensemble_reduced_functional.EnsembleReducedFunctional(J, control, ensemble, scatter_control=True, gather_functional=None)[source]

Bases: ReducedFunctional

Enable solving simultaneously reduced functionals in parallel.

Consider a functional \(J\) and its gradient \(\dfrac{dJ}{dm}\), where \(m\) is the control parameter. Let us assume that \(J\) is the sum of \(N\) functionals \(J_i(m)\), i.e.,

\[J = \sum_{i=1}^{N} J_i(m).\]

The gradient over a summation is a linear operation. Therefore, we can write the gradient \(\dfrac{dJ}{dm}\) as

\[\frac{dJ}{dm} = \sum_{i=1}^{N} \frac{dJ_i}{dm},\]

The EnsembleReducedFunctional allows simultaneous evaluation of \(J_i\) and \(\dfrac{dJ_i}{dm}\). After that, the allreduce Ensemble operation is employed to sum the functionals and their gradients over an ensemble communicator.

If gather_functional is present, then all the values of J are communicated to all ensemble ranks, and passed in a list to gather_functional, which is a reduced functional that expects a list of that size of the relevant types.

Parameters:
  • J (pyadjoint.OverloadedType) – An instance of an OverloadedType, usually pyadjoint.AdjFloat. This should be the functional that we want to reduce.

  • control (pyadjoint.Control or list of pyadjoint.Control) – A single or a list of Control instances, which you want to map to the functional.

  • ensemble (Ensemble) – An instance of the Ensemble. It is used to communicate the functionals and their derivatives between the ensemble members.

  • scatter_control (bool) – Whether scattering a control (or a list of controls) over the ensemble communicator Ensemble.ensemble comm.

  • gather_functional (An instance of the pyadjoint.ReducedFunctional.) – that takes in all of the Js.

Notes

The functionals \(J_i\) and the control must be defined over a common ensemble.comm communicator. To understand more about how ensemble parallelism works, please refer to the Firedrake manual.

__call__(values)[source]

Computes the reduced functional with supplied control value.

Parameters:

values (pyadjoint.OverloadedType) – If you have multiple controls this should be a list of new values for each control in the order you listed the controls to the constructor. If you have a single control it can either be a list or a single object. Each new value should have the same type as the corresponding control.

Returns:

The computed value. Typically of instance of pyadjoint.AdjFloat.

Return type:

pyadjoint.OverloadedType

derivative(adj_input=1.0, options=None)[source]

Compute derivatives of a functional with respect to the control parameters.

Parameters:
  • adj_input (float) – The adjoint input.

  • options (dict) – Additional options for the derivative computation.

Returns:

  • dJdm_total (pyadjoint.OverloadedType)

  • The result of Allreduce operations of dJdm_local into dJdm_total over the`Ensemble.ensemble_comm`.

hessian(m_dot, options=None)[source]

The Hessian is not yet implemented for ensemble reduced functional.

Raises:

NotImplementedError – This method is not yet implemented for ensemble reduced functional.

firedrake.adjoint.ufl_constraints module

class firedrake.adjoint.ufl_constraints.UFLConstraint(form, control)[source]

Bases: Constraint

Easily implement scalar constraints using UFL.

The form must be a 0-form that depends on a Function control.

function(m)[source]

Evaluate c(m), where c(m) == 0 for equality constraints and c(m) >= 0 for inequality constraints.

c(m) must return a numpy array or a dolfin Function or Constant.

hessian_action(m, dm, dp, result)[source]

Computes the Hessian action of c(m) in direction dm and dp.

Stores the result in result.

jacobian(m)[source]

Returns the full Jacobian matrix as a list of vector-like objects representing the gradient of the constraint function with respect to the parameter m.

The objects returned must be of the same type as m’s data.

jacobian_action(m, dm, result)[source]

Computes the Jacobian action of c(m) in direction dm.

Stores the result in result.

jacobian_adjoint_action(m, dp, result)[source]

Computes the Jacobian adjoint action of c(m) in direction dp.

Stores the result in result.

output_workspace()[source]

Return an object like the output of c(m) for calculations.

update_control(m)[source]
class firedrake.adjoint.ufl_constraints.UFLEqualityConstraint(form, control)[source]

Bases: UFLConstraint, EqualityConstraint

class firedrake.adjoint.ufl_constraints.UFLInequalityConstraint(form, control)[source]

Bases: UFLConstraint, InequalityConstraint

Module contents

The public interface to Firedrake’s adjoint.

To start taping, run:

from firedrake.adjoint import *
continue_annotation()