TSEventHandler#
the main function to perform a single iteration of event detection.
Synopsis#
#include "petscts.h"
PetscErrorCode TSEventHandler(TS ts)
Developer notes#
A) The ‘event->iterctr > 0’ is used as an indicator that Anderson-Bjorck refinement has started. B) If event->iterctr == 0, then justrefined_AB[i] is always false. C) The right-end quantities: ptime_right, fvalue_right[i] and fsign_right[i] are only guaranteed to be valid for event->iterctr > 0. D) If event->iterctr > 0, then event->processing is PETSC_TRUE; the opposite may not hold. E) When event->processing == PETSC_TRUE and event->iterctr == 0, the event handler iterations are complete, but the event handler continues managing the 1st and 2nd post-event steps. In this case the 1st post-event step proposed by the event handler is not checked by TSAdapt, and is always accepted (beware!). However, if the 2nd post-event step is not managed by the event handler (e.g. 1st = numerical, 2nd = PETSC_DECIDE), condition “E” does not hold, and TSAdapt may reject/adjust the 1st post-event step. F) event->side[i] may take values: 0 <=> point t is a zero-crossing for indicator function i (via vtol/dt_min criterion); -1/+1 <=> detected a bracket to the left/right of t for indicator function i; +2 <=> no brackets/zero-crossings. G) The signs event->fsign[i] (with values 0/-1/+1) are calculated for each new point. Zero sign is set if the function value is smaller than the tolerance. Besides, zero sign is enforced after marking a zero-crossing due to small bracket size criterion.
The intervals with the indicator function sign change (i.e. containing the potential zero-crossings) are called ‘brackets’. To find a zero-crossing, the algorithm first locates a bracket, and then sequentially subdivides it, generating a sequence of brackets whose length tends to zero. The bracket subdivision involves the (modified) Anderson-Bjorck method.
Apart from the comments scattered throughout the code to clarify different lines and blocks,
a few tricky aspects of the algorithm (and the underlying reasoning) are discussed in detail below#
=Sign tracking= When a zero-crossing is found, the sign variable (event->fsign[i]) is set to zero for the current point t. This happens both for zero-crossings triggered via the vtol criterion, and those triggered via the dt_min criterion. After the event, as the TS steps forward, the current sign values are handed over to event->fsign_prev[i]. The recalculation of signs is avoided if possible: e.g. if a ‘vtol’ criterion resulted in a zero-crossing at point t, but the subsequent call to postevent() handler decreased ‘vtol’, making the indicator function no longer “close to zero” at point t, the fsign[i] will still consistently keep the zero value. This allows avoiding the erroneous duplication
of events#
E.g. consider a bracket [t0, t2], where f0 < 0, f2 > 0, which resulted in a zero-crossing t1 with f1 < 0, abs(f1) < vtol. Suppose the postevent() handler changes vtol to vtol*, such that abs(f1) > vtol*. The TS makes a step t1 -> t3, where again f1 < 0, f3 > 0, and the event handler will find a new event near t1, which is actually a duplication of the original event at t1. The duplications are avoided by NOT counting the sign progressions 0 -> +1, or 0 -> -1 as brackets. Tracking (instead of recalculating) the sign values makes this procedure work more consistently.
The sign values are however recalculated if the postevent() callback has changed the current solution vector U (such a change resets everything). The sign value is also set to zero if the dt_min criterion has triggered the event. This allows the algorithm to work consistently, irrespective of the type of criterion involved (vtol/dt_min).
=Event from min bracket= When the event handler ends up with a bracket [t0, t1] with size <= dt_min, a zero crossing is reported at t1, and never at t0. If such a bracket is discovered when TS is staying at t0, one more step forward (to t1) is necessary to mark the found event. This is the situation of revisiting t1, which is described below (see Revisiting).
Why t0 is not reported as event location? Suppose it is, and let f0 < 0, f1 > 0. Also suppose that the postevent() handler has slightly changed the solution U, so the sign at t0 is recalculated: it equals -1. As the TS steps further: t0 -> t2, with sign0 == -1, and sign2 == +1, the event handler will locate the bracket [t0, t2], eventually resolving a new event near t1, i.e. finding a duplicate event. This situation is avoided by reporting the event at t1 in the first place.
=Revisiting=
When handling the situation with small bracket size, the TS solver may happen to visit the same point twice,
but with different results.
E.g. originally it discovered a bracket with sign change [t0, t10], and started resolving the zero-crossing, visiting the points t1,…,t9 : t0 < t1 < … < t9 < t10. Suppose that at t9 the algorithm discovers that [t9, t10] is a bracket with the sign change it was looking for, and that |t10 - t9| is too small. So point t10 should be revisited and marked as the zero crossing (by the minimum bracket size criterion). On re-visiting t10, via the refined sequence of steps t0,…,t10, the TS solver may arrive at a solution U* different from the solution U it found at t10 originally. Hence, the indicator functions at t10 may become different, and the condition of the sign change, which existed originally, may disappear, breaking the logic of the algorithm.
To handle such (-=unlikely=-, but possible) situations, two strategies can be considered#
[not used here] Allow the brackets with sign change to disappear during iterations. The algorithm should be able to cleanly exit the iteration and leave all the objects/variables/caches involved in a valid state.
[ADOPTED HERE!] On revisiting t10, the event handler reuses the indicator functions previously calculated for the original solution U. This U may be less precise than U*, but this trick does not allow the algorithm logic to break down. HOWEVER, the original U is not stored anywhere, it is essentially lost since the TS performed the rollback from it. On revisiting t10, the updated solution U* will inevitably be found and used everywhere EXCEPT the current indicator functions calculation, e.g. U* will be used in the postevent() handler call. Since t10 is the event location, the appropriate indicator-function-signs will be enforced to be 0 (regardless if the solution was U or U*). If the solution is then changed by the postevent(), the indicator-function-signs will be recalculated.
Whether the algorithm is revisiting a point in the current TSEventHandler() call is flagged by ‘event->revisit_right’.
See Also#
Handling of discontinuities, TS, TSEvent, TSSetEventHandler()
Level#
developer
Location#
Index of all TS routines
Table of Contents for all manual pages
Index of all manual pages