goalie package¶
Submodules¶
goalie.adjoint module¶
Drivers for solving adjoint problems on sequences of meshes.
- class AdjointMeshSeq(time_partition, initial_meshes, **kwargs)[source]¶
Bases:
MeshSeq
An extension of
MeshSeq
to account for solving adjoint problems on a sequence of meshes.For time-dependent quantities of interest, the solver should access and modify
J
, which holds the QoI value.- Parameters:
time_partition (
TimePartition
) – a partition of the temporal domaininitial_meshes (
list
orMeshGeometry
) – a list of meshes corresponding to the subinterval of the time partition, or a single mesh to use for all subintervalsget_function_spaces – a function as described in
get_function_spaces()
get_initial_condition – a function as described in
get_initial_condition()
get_form – a function as described in
get_form()
get_solver – a function as described in
get_solver()
get_qoi – a function as described in
get_qoi()
- check_qoi_convergence()[source]¶
Check for convergence of the fixed point iteration due to the relative difference in QoI value being smaller than the specified tolerance.
- Returns:
True
if QoI convergence is detected, elseFalse
- Return type:
- get_checkpoints(solver_kwargs=None, run_final_subinterval=False)[source]¶
Solve forward on the sequence of meshes, extracting checkpoints corresponding to the starting fields on each subinterval.
The QoI is also evaluated.
- Parameters:
- Returns:
checkpoints for each subinterval
- Return type:
- get_qoi(subinterval)[source]¶
Get the function for evaluating the QoI, which has either zero or one arguments, corresponding to either an end time or time integrated quantity of interest, respectively. If the QoI has an argument then it is for the current time.
Signature for the function to be returned:
` :arg t: the current time (for time-integrated QoIs) :type t: :class:`float` :return: the QoI as a 0-form :rtype: :class:`ufl.form.Form` `
- Parameters:
solution_map (
dict
withstr
keys and values andfiredrake.function.Function
values) – a dictionary whose keys are the solution field names and whose values are the corresponding solutionssubinterval (
int
) – the subinterval index
- Returns:
the function for obtaining the QoI
- Return type:
see docstring above
- get_solve_blocks(field, subinterval, has_adj_sol=True)[source]¶
Get all blocks of the tape corresponding to solve steps for prognostic solution field on a given subinterval.
- property initial_condition¶
Get the initial conditions associated with the first subinterval.
- Returns:
a dictionary whose keys are field names and whose values are the corresponding initial conditions applied on the first subinterval
- Return type:
AttrDict
withstr
keys andfiredrake.function.Function
values
- solve_adjoint(solver_kwargs=None, adj_solver_kwargs=None, get_adj_values=False, test_checkpoint_qoi=False)[source]¶
Solve an adjoint problem on a sequence of subintervals.
As well as the quantity of interest value, solution fields are computed - see
AdjointSolutionData
for more information.- Parameters:
solver_kwargs (
dict
withstr
keys and values which may take various types) – parameters for the forward solver, as well as any parameters for the QoI, which should be included as a sub-dictionary with key ‘qoi_kwargs’adj_solver_kwargs (
dict
withstr
keys and values which may take various types) – parameters for the adjoint solverget_adj_values (
bool
) – ifTrue
, adjoint actions are also returned at exported timestepstest_checkpoint_qoi – solve over the final subinterval when checkpointing so that the QoI value can be checked across runs
- Returns:
the solution data of the forward and adjoint solves
- Return type:
goalie.error_estimation module¶
Tools to automate goal-oriented error estimation.
- get_dwr_indicator(F, adjoint_error, test_space=None)[source]¶
Given a 1-form and an approximation of the error in the adjoint solution, compute a dual weighted residual (DWR) error indicator.
Note that each term of a 1-form contains only one
firedrake.ufl_expr.TestFunction
. The 1-form most commonly corresponds to the variational form of a PDE. If the PDE is linear, it should be written as in the nonlinear case (i.e., with the solution field in place of anyfiredrake.ufl_expr.TrialFunction
s.- Parameters:
F (
ufl.form.Form
) – the formadjoint_error (
firedrake.function.Function
ordict
withstr
keys andfiredrake.function.Function
values) – a dictionary whose keys are field names and whose values are the approximations to the corresponding components of the adjoint error, or a single such componenttest_space (
firedrake.functionspaceimpl.WithGeometry
) – a dictionary whose keys are field names and whose values are the test spaces for the corresponding fields, or a single such test space (orNone
to determine the test space(s) automatically)
- Returns:
the DWR indicator
- Return type:
goalie.function_data module¶
Nested dictionaries of solution data Function
s.
- class AdjointSolutionData(*args, **kwargs)[source]¶
Bases:
FunctionData
Class representing solution data for general adjoint problems.
For a given exported timestep, the field types are:
'forward'
: the forward solution after taking the timestep;'forward_old'
: the forward solution before taking the timestep (provided the problem is not steady-state)'adjoint'
: the adjoint solution after taking the timestep;'adjoint_next'
: the adjoint solution before taking the timestep backwards (provided the problem is not steady-state).
- Parameters:
time_partition – the
TimePartition
used to discretise the problem in timefunction_spaces – the dictionary of
FunctionSpace
s used to discretise the problem in space
- class ForwardSolutionData(*args, **kwargs)[source]¶
Bases:
FunctionData
Class representing solution data for general forward problems.
For a given exported timestep, the field types are:
'forward'
: the forward solution after taking the timestep;'forward_old'
: the forward solution before taking the timestep (provided the problem is not steady-state).
- Parameters:
time_partition – the
TimePartition
used to discretise the problem in timefunction_spaces – the dictionary of
FunctionSpace
s used to discretise the problem in space
- class IndicatorData(time_partition, meshes)[source]¶
Bases:
FunctionData
Class representing error indicator data.
Note that this class has a single dictionary with the field name as the key, rather than a doubly-nested dictionary.
- Parameters:
time_partition – the
TimePartition
used to discretise the problem in timemeshes – the list of meshes used to discretise the problem in space
goalie.go_mesh_seq module¶
Drivers for goal-oriented error estimation on sequences of meshes.
- class GoalOrientedMeshSeq(*args, **kwargs)[source]¶
Bases:
AdjointMeshSeq
An extension of
AdjointMeshSeq
to account for goal-oriented problems.- Parameters:
time_partition (
TimePartition
) – a partition of the temporal domaininitial_meshes (
list
orMeshGeometry
) – a list of meshes corresponding to the subinterval of the time partition, or a single mesh to use for all subintervalsget_function_spaces – a function as described in
get_function_spaces()
get_initial_condition – a function as described in
get_initial_condition()
get_form – a function as described in
get_form()
get_solver – a function as described in
get_solver()
get_qoi – a function as described in
get_qoi()
- check_estimator_convergence()[source]¶
Check for convergence of the fixed point iteration due to the relative difference in error estimator value being smaller than the specified tolerance.
- Returns:
True
if estimator convergence is detected, elseFalse
- Return type:
- error_estimate(absolute_value=False)[source]¶
Deduce the error estimator value associated with error indicator fields defined over the mesh sequence.
- fixed_point_iteration(adaptor, parameters=None, update_params=None, enrichment_kwargs=None, adaptor_kwargs=None, solver_kwargs=None, indicator_fn=<cyfunction get_dwr_indicator>)[source]¶
Apply goal-oriented mesh adaptation using a fixed point iteration loop approach.
- Parameters:
adaptor – function for adapting the mesh sequence. Its arguments are the mesh sequence and the solution and indicator data objects. It should return
True
if the convergence criteria checks are to be skipped for this iteration. Otherwise, it should returnFalse
.parameters (
GoalOrientedAdaptParameters
) – parameters to apply to the mesh adaptation processupdate_params – function for updating
params
at each iteration. Its arguments are the parameter class and the fixed point iterationenrichment_kwargs (
dict
withstr
keys and values which may take various types) – keyword arguments to pass to the global enrichment methodsolver_kwargs (
dict
withstr
keys and values which may take various types) – parameters to pass to the solveradaptor_kwargs (
dict
withstr
keys and values which may take various types) – parameters to pass to the adaptorindicator_fn – function which maps the form, adjoint error and enriched space(s) as arguments to the error indicator
firedrake.function.Function
- Returns:
solution and indicator data objects
- Rtype1:
:class:`~.AdjointSolutionData
- Rtype2:
:class:`~.IndicatorData
- get_enriched_mesh_seq(enrichment_method='p', num_enrichments=1)[source]¶
Construct a sequence of globally enriched spaces.
The following global enrichment methods are supported: * h-refinement (
enrichment_method='h'
) - refine each mesh elementuniformly in each direction;
p-refinement (
enrichment_method='p'
) - increase the function space polynomial order by one globally.
- indicate_errors(enrichment_kwargs=None, solver_kwargs=None, indicator_fn=<cyfunction get_dwr_indicator>)[source]¶
Compute goal-oriented error indicators for each subinterval based on solving the adjoint problem in a globally enriched space.
- Parameters:
enrichment_kwargs (
dict
withstr
keys and values which may take various types) – keyword arguments to pass to the global enrichment method - seeget_enriched_mesh_seq()
for the supported enrichment methods and optionssolver_kwargs (
dict
withstr
keys and values which may take various types) – parameters for the forward solver, as well as any parameters for the QoI, which should be included as a sub-dictionary with key ‘qoi_kwargs’indicator_fn – function which maps the form, adjoint error and enriched space(s) as arguments to the error indicator
firedrake.function.Function
- Returns:
solution and indicator data objects
- Rtype1:
:class:`~.AdjointSolutionData
- Rtype2:
:class:`~.IndicatorData
- property indicators¶
- Returns:
the error indicator data object
- Return type:
goalie.log module¶
Loggers for Goalie.
Code mostly copied from the Thetis project.
- critical(msg, *args, **kwargs)[source]¶
Log ‘msg % args’ with severity ‘CRITICAL’.
To pass exception information, use the keyword argument exc_info with a true value, e.g.
logger.critical(“Houston, we have a %s”, “major disaster”, exc_info=1)
- debug(msg, *args, **kwargs)[source]¶
Log ‘msg % args’ with severity ‘DEBUG’.
To pass exception information, use the keyword argument exc_info with a true value, e.g.
logger.debug(“Houston, we have a %s”, “thorny problem”, exc_info=1)
- error(msg, *args, **kwargs)[source]¶
Log ‘msg % args’ with severity ‘ERROR’.
To pass exception information, use the keyword argument exc_info with a true value, e.g.
logger.error(“Houston, we have a %s”, “major problem”, exc_info=1)
- info(msg, *args, **kwargs)[source]¶
Log ‘msg % args’ with severity ‘INFO’.
To pass exception information, use the keyword argument exc_info with a true value, e.g.
logger.info(“Houston, we have a %s”, “interesting problem”, exc_info=1)
- pyrint(msg, *args, **kwargs)¶
Log ‘msg % args’ with severity ‘INFO’.
To pass exception information, use the keyword argument exc_info with a true value, e.g.
logger.info(“Houston, we have a %s”, “interesting problem”, exc_info=1)
goalie.math module¶
goalie.mesh_seq module¶
Sequences of meshes corresponding to a TimePartition
.
- class MeshSeq(time_partition, initial_meshes, **kwargs)[source]¶
Bases:
object
A sequence of meshes for solving a PDE associated with a particular
TimePartition
of the temporal domain.- Parameters:
time_partition (
TimePartition
) – a partition of the temporal domaininitial_meshes (
list
orMeshGeometry
) – a list of meshes corresponding to the subinterval of the time partition, or a single mesh to use for all subintervalsget_function_spaces – a function as described in
get_function_spaces()
get_initial_condition – a function as described in
get_initial_condition()
get_form – a function as described in
get_form()
get_solver – a function as described in
get_solver()
transfer_method (
str
) – the method to use for transferring fields between meshes. Options are “project” (default) and “interpolate”. Seeanimate.interpolation.transfer()
for detailstransfer_kwargs (
dict
withstr
keys and values which may take various types) – kwargs to pass to the chosen transfer method
- check_element_count_convergence()[source]¶
Check for convergence of the fixed point iteration due to the relative difference in element count being smaller than the specified tolerance.
- fixed_point_iteration(adaptor, parameters=None, update_params=None, solver_kwargs=None, adaptor_kwargs=None)[source]¶
Apply mesh adaptation using a fixed point iteration loop approach.
- Parameters:
adaptor – function for adapting the mesh sequence. Its arguments are the mesh sequence and the solution data object. It should return
True
if the convergence criteria checks are to be skipped for this iteration. Otherwise, it should returnFalse
.parameters (
AdaptParameters
) – parameters to apply to the mesh adaptation processupdate_params – function for updating
params
at each iteration. Its arguments are the parameter class and the fixed point iterationsolver_kwargs (
dict
withstr
keys and values which may take various types) – parameters to pass to the solveradaptor_kwargs (
dict
withstr
keys and values which may take various types) – parameters to pass to the adaptor
- Returns:
solution data object
- Return type:
- property form¶
See
get_form()
.
- property function_spaces¶
Get the function spaces associated with the mesh sequence.
- Returns:
a dictionary whose keys are field names and whose values are the corresponding function spaces
- Return type:
AttrDict
withstr
keys andfiredrake.functionspaceimpl.FunctionSpace
values
- get_checkpoints(run_final_subinterval=False, solver_kwargs=None)[source]¶
Get checkpoints corresponding to the starting fields on each subinterval.
- Parameters:
- Returns:
checkpoints for each subinterval
- Return type:
- get_form()[source]¶
Get the function mapping a subinterval index and a solution dictionary to a dictionary containing parts of the PDE weak form corresponding to each solution component.
Signature for the function to be returned:
` :arg index: the subinterval index :type index: :class:`int` :arg solutions: map from fields to tuples of current and previous solutions :type solutions: :class:`dict` with :class:`str` keys and :class:`tuple` values :return: map from fields to the corresponding forms :rtype: :class:`dict` with :class:`str` keys and :class:`ufl.form.Form` values `
- Returns:
the function for obtaining the form
- Return type:
see docstring above
- get_function_spaces(mesh)[source]¶
Construct the function spaces corresponding to each field, for a given mesh.
- Parameters:
mesh (
firedrake.mesh.MeshGeometry
) – the mesh to base the function spaces on- Returns:
a dictionary whose keys are field names and whose values are the corresponding function spaces
- Return type:
dict
withstr
keys andfiredrake.functionspaceimpl.FunctionSpace
values
- get_initial_condition()[source]¶
Get the initial conditions applied on the first mesh in the sequence.
- Returns:
the dictionary, whose keys are field names and whose values are the corresponding initial conditions applied
- Return type:
dict
withstr
keys andfiredrake.function.Function
values
- get_solver()[source]¶
Get the function mapping a subinterval index and an initial condition dictionary to a dictionary of solutions for the corresponding solver step.
Signature for the function to be returned: ``` :arg index: the subinterval index :type index:
int
:arg ic: map from fields to the corresponding initial condition components :type ic:dict
withstr
keys andfiredrake.function.Function
values- Returns:
map from fields to the corresponding solutions
- Return type:
dict
withstr
keys andfiredrake.function.Function
values
- Returns:
the function for obtaining the solver
- Return type:
see docstring above
- property initial_condition¶
Get the initial conditions associated with the first subinterval.
- Returns:
a dictionary whose keys are field names and whose values are the corresponding initial conditions applied on the first subinterval
- Return type:
AttrDict
withstr
keys andfiredrake.function.Function
values
- plot(fig=None, axes=None, **kwargs)[source]¶
Plot the meshes comprising a 2D
MeshSeq
.- Parameters:
fig (
matplotlib.figure.Figure
) – matplotlib figure to useaxes (
matplotlib.axes._axes.Axes
) – matplotlib axes to use
- Returns:
matplotlib figure and axes for the plots
- Rtype1:
matplotlib.figure.Figure
- Rtype2:
matplotlib.axes._axes.Axes
All keyword arguments are passed to
firedrake.pyplot.triplot()
.
- set_meshes(meshes)[source]¶
Set all meshes in the sequence and deduce various properties.
- Parameters:
meshes (
list
offiredrake.MeshGeometry
s orfiredrake.MeshGeometry
) – list of meshes to use in the sequence, or a single mesh to use for all subintervals
- property solutions¶
- Returns:
the solution data object
- Return type:
FunctionData
- solve_forward(solver_kwargs=None)[source]¶
Solve a forward problem on a sequence of subintervals.
A dictionary of solution fields is computed - see
ForwardSolutionData
for more details.
- property solver¶
See
get_solver()
.
goalie.metric module¶
Driver functions for metric-based mesh adaptation.
- enforce_variable_constraints(metrics, h_min=1e-30, h_max=1e+30, a_max=100000.0, boundary_tag=None)[source]¶
Post-process a list of metrics to enforce minimum and maximum element sizes, as well as maximum anisotropy.
- Parameters:
metrics (
list
ofRiemannianMetric
s) – the metricsh_min (
firedrake.function.Function
,float
, orint
) – minimum tolerated element sizeh_max (
firedrake.function.Function
,float
, orint
) – maximum tolerated element sizea_max (
firedrake.function.Function
,float
, orint
) – maximum tolerated element anisotropyboundary_tag (
str
orint
) – optional tag to enforce sizes on.
- ramp_complexity(base, target, iteration, num_iterations=3)[source]¶
Ramp up the target complexity over the first few iterations.
- space_time_normalise(metrics, time_partition, metric_parameters, global_factor=None, boundary=False, restrict_sizes=True, restrict_anisotropy=True)[source]¶
Apply \(L^p\) normalisation in both space and time.
Based on Equation (1) in [].
- Parameters:
metrics (
list
ofRiemannianMetric
s) – the metrics associated with each subintervaltime_partition (
TimePartition
) – temporal discretisation for the problem at handmetric_parameters (
list
ofdict
s or a singledict
to use for all subintervals) – dictionary containing the target space-time metric complexity under dm_plex_metric_target_complexity and the normalisation order under dm_plex_metric_p, or a list thereofglobal_factor (
float
) – pre-computed global normalisation factorboundary (
bool
) – ifTrue
, the normalisation to be performed over the boundaryrestrict_sizes (
bool
) – ifTrue
, minimum and maximum metric magnitudes are enforcedrestrict_anisotropy (
bool
) – ifTrue
, maximum anisotropy is enforced
- Returns:
the space-time normalised metrics
- Return type:
list
ofRiemannianMetric
s
goalie.options module¶
- class AdaptParameters(parameters=None)[source]¶
Bases:
AttrDict
A class for holding parameters associated with adaptive mesh fixed point iteration loops.
- class GoalOrientedAdaptParameters(parameters=None)[source]¶
Bases:
AdaptParameters
A class for holding parameters associated with goal-oriented adaptive mesh fixed point iteration loops.
goalie.plot module¶
Driver functions for plotting solution data.
- plot_indicator_snapshots(indicators, time_partition, field, **kwargs)[source]¶
Plot a sequence of snapshots associated with
indicators
and a givenTimePartition
Any keyword arguments are passed to
firedrake.plot.tricontourf()
.- Parameters:
indicators – list of list of indicators, indexed by mesh sequence index, then timestep
time_partition – the
TimePartition
object used to solve the problem
- plot_snapshots(solutions, time_partition, field, label, **kwargs)[source]¶
Plot a sequence of snapshots associated with
solutions.field.label
and a givenTimePartition
.Any keyword arguments are passed to
firedrake.plot.tricontourf()
.- Parameters:
solutions –
AttrDict
of solutions computed by solving a forward or adjoint problemtime_partition – the
TimePartition
object used to solve the problemfield – solution field of choice
label – choose from
'forward'
,'forward_old'
'adjoint'
and'adjoint_next'
goalie.point_seq module¶
- class PointSeq(time_partition, **kwargs)[source]¶
Bases:
MeshSeq
A simplified subset of
MeshSeq
for ODE problems.In this version, a single mesh comprised of a single vertex is shared across all subintervals.
- Parameters:
time_partition – the
TimePartition
which partitions the temporal domainget_function_spaces – a function, whose only argument is a
MeshSeq
, which constructs prognosticfiredrake.functionspaceimpl.FunctionSpace
s for each subintervalget_initial_condition – a function, whose only argument is a
MeshSeq
, which specifies initial conditions on the first meshget_form – a function, whose only argument is a
MeshSeq
, which returns a function that generates the ODE weak formget_solver – a function, whose only argument is a
MeshSeq
, which returns a function that integrates initial data over a subintervalget_bcs – a function, whose only argument is a
MeshSeq
, which returns a function that determines any Dirichlet boundary conditions
goalie.time_partition module¶
Partitioning for the temporal domain.
- class TimeInstant(field_names, **kwargs)[source]¶
Bases:
TimeInterval
A
TimePartition
for steady-state problems.Under the hood this means dividing \([0,1)\) into a single timestep.
- Parameters:
end_time (
float
orint
) – end time of the interval of interestnum_subintervals (
int
) – number of subintervals in the partitiontimesteps (
list
offloat
s orfloat
) – a list timesteps to be used on each subinterval, or a single timestep to use for all subintervalsfield_names (
list
ofstr
s orstr
) – the list of field names to considernum_timesteps_per_export (
list
of :class`int`s orint
) – a list of numbers of timesteps per export for each subinterval, or a single number to use for all subintervalsstart_time (
float
orint
) – start time of the interval of interestsubinterals – sequence of subintervals (which need not be of uniform length), or
None
to use uniform subintervals (the default)field_types (
list
ofstr
s orstr
) – a list of strings indicating whether each field is ‘unsteady’ or ‘steady’, i.e., does the corresponding equation involve time derivatives or not?
- class TimeInterval(*args, **kwargs)[source]¶
Bases:
TimePartition
A trivial
TimePartition
with a single subinterval.- Parameters:
end_time (
float
orint
) – end time of the interval of interestnum_subintervals (
int
) – number of subintervals in the partitiontimesteps (
list
offloat
s orfloat
) – a list timesteps to be used on each subinterval, or a single timestep to use for all subintervalsfield_names (
list
ofstr
s orstr
) – the list of field names to considernum_timesteps_per_export (
list
of :class`int`s orint
) – a list of numbers of timesteps per export for each subinterval, or a single number to use for all subintervalsstart_time (
float
orint
) – start time of the interval of interestsubinterals – sequence of subintervals (which need not be of uniform length), or
None
to use uniform subintervals (the default)field_types (
list
ofstr
s orstr
) – a list of strings indicating whether each field is ‘unsteady’ or ‘steady’, i.e., does the corresponding equation involve time derivatives or not?
- class TimePartition(end_time, num_subintervals, timesteps, field_names, num_timesteps_per_export=1, start_time=0.0, subintervals=None, field_types=None)[source]¶
Bases:
object
A partition of the time interval of interest into subintervals.
The subintervals are assumed to be uniform in length. However, different timestep values may be used on each subinterval.
- Parameters:
end_time (
float
orint
) – end time of the interval of interestnum_subintervals (
int
) – number of subintervals in the partitiontimesteps (
list
offloat
s orfloat
) – a list timesteps to be used on each subinterval, or a single timestep to use for all subintervalsfield_names (
list
ofstr
s orstr
) – the list of field names to considernum_timesteps_per_export (
list
of :class`int`s orint
) – a list of numbers of timesteps per export for each subinterval, or a single number to use for all subintervalsstart_time (
float
orint
) – start time of the interval of interestsubinterals – sequence of subintervals (which need not be of uniform length), or
None
to use uniform subintervals (the default)field_types (
list
ofstr
s orstr
) – a list of strings indicating whether each field is ‘unsteady’ or ‘steady’, i.e., does the corresponding equation involve time derivatives or not?
goalie.utility module¶
Utility functions and classes for mesh adaptation.
- class AttrDict(*args, **kwargs)[source]¶
Bases:
dict
Dictionary that provides both
self[key]
andself.key
access to members.Disclaimer: Copied from stackoverflow.
- create_directory(path, comm=<mpi4py.MPI.Intracomm object>)[source]¶
Create a directory on disk.
Disclaimer: Code copied from Thetis.
- effectivity_index(error_indicator, Je)[source]¶
Compute the overestimation factor of some error estimator for the QoI error.
Note that this is only typically used for simple steady-state problems with analytical solutions.
- Parameters:
error_indicator (
firedrake.function.Function
) – a \(\mathbb P0\) error indicator which localises contributions to an error estimator to individual elementsJe (
float
) – the error in the quantity of interest
- Returns:
the effectivity index
- Return type: