Better syntax and functionality for functionals

Registered by David Ham

Libadjoint functionals should be sums of terms, each of which is either a form integrated over a time interval, or a form evaluated at a particular time interval.

Blueprint information

Status:
Complete
Approver:
None
Priority:
High
Drafter:
None
Direction:
Needs approval
Assignee:
None
Definition:
Approved
Series goal:
None
Implementation:
Implemented
Milestone target:
None
Started by
Patrick Farrell
Completed by
Patrick Farrell

Sprints

Whiteboard

There are two distinct cases:

1. How to express time in current DOLFIN models with no time abstraction.

2. Whatever happens in the future with a DOLFIN time abstraction.

1. is what we should do now. We should handle when 2. someone does
the DOLFIN time abstraction.

TimeForms
=========

We proceeded from the ansatz that functionals will be forms integrated
over time, forms evaluated at a particular time (which might be the
start or end), and sums of the preceding two cases.

We can follow UFL in handling sums via a container object containing a
list of summed terms. For instance, we can define a TimeForm object
containing a list of (form, time) tuples. This easily allows addition
to be defined and closed (ie the sum of two TimeForms is a TimeForm).
The Functional object would take a TimeForm and would loop over its
list as required to evaluate all the terms.

In the case of a form evaluated at one time, a possible constructor might be:
at_time(form, time). Time would ,in this case, be a float or a named
constant for start or end.

In the case of a form integrated over a time interval, a dt time
measure object would be a suitable constructor. So, for example,
x*y*dx*dt would indicate an integral over all time, while
x*y*dx*dt[a:b] would indicate an integral from time a to time b. Time,
in this case, would be a slice object.

Some other time adjustments
=======================

In order to allow for Functional evaluation on the first execution of
a model (this is needed for optimisation), it should be possible to
register a Functional before we start. This requires a little more
timestep-related annotation than we have at the moment
(unfortunately). It would also be a really good idea not to assume
that time progresses in constant timesteps from zero. I think this
means we need the following:

adj_time_start(time) : start set the start time.
adj_inc_timestep(time): adj_inc_timestep now records the time at which
the timestep ends.
adj_time_finish(): indicate that time has ended and any final
functionals can be computed.

(incidentally, using prefixes is unpythonic. Should we define an
adjoint or a libadjoint namespace and have adjoint.time_start, for
example? Actually one way we could do it is to recognise the adjointer
as the object embodying the recorded simulation and therefore have
adjointer.time_start() or even adjointer.time.start() opinions?)

Evaluating functionals automatically during forward execution
================================================

A registered functional should be automatically evaluated during forward evaluation using the trapezoidal rule for integrals, and a pointwise interpolation for at_time functionals.

The functional evaluation should occur inside the adj_inc_timestep, before any values are deleted from libadjoint's in-core memory by the checkpointer. A final evaluation will need to occur inside adj_time_finish.

Evaluating time integrals
---------------------------------

Assume we have a TimeForm of the form F(u)*dt[a:b]

The support of u(T^n) in a time integral is bounded by the time interval (T^(n-1), T^(n+1)). At time T(n+1) we should therefore calculate the contribution w*(F(u(T^n)). It will be possible to calculate w by intersecting (a,b) with (T^(n-1), T^(n+1)) and libadjoint cannot possibly have forgotten u(T^n) yet since we have not called the checkpointer . Since we always calculate the previous integral contribution, we need to calculate the final contribution for the interval (T^(-2),T^(-1)) at time n.

Evaluating at_time terms
---------------------------------

For at_time terms, the situation depends on whether we prefer:

F(u(T^(n+theta))) ~= F( theta(u(T^n)) + (1-theta)(u(T^(n+1))) )

or

F(u(T^(n+theta))) ~= theta(F(u(T^n))) + (1-theta)(F(u(T^(n+1)))

The former is more accurate, the latter requires fewer functional evaluations. For sufficiently regular u and F, both have the same order error term.

(?)

Work Items

Work items:
Implement TimeForms : INPROGRESS
Improve time annotations : TODO
Re-implement Functional : TODO

This blueprint contains Public information 
Everyone can see this information.