Actual source code: tsevent.c
1: #include <petsc/private/tsimpl.h>
3: /*
4: TSEventInitialize - Initializes TSEvent for TSSolve
5: */
6: PetscErrorCode TSEventInitialize(TSEvent event, TS ts, PetscReal t, Vec U)
7: {
8: if (!event) return 0;
12: event->ptime_prev = t;
13: event->iterctr = 0;
14: (*event->eventhandler)(ts, t, U, event->fvalue_prev, event->ctx);
15: return 0;
16: }
18: PetscErrorCode TSEventDestroy(TSEvent *event)
19: {
20: PetscInt i;
23: if (!*event) return 0;
24: if (--(*event)->refct > 0) {
25: *event = NULL;
26: return 0;
27: }
29: PetscFree((*event)->fvalue);
30: PetscFree((*event)->fvalue_prev);
31: PetscFree((*event)->fvalue_right);
32: PetscFree((*event)->zerocrossing);
33: PetscFree((*event)->side);
34: PetscFree((*event)->direction);
35: PetscFree((*event)->terminate);
36: PetscFree((*event)->events_zero);
37: PetscFree((*event)->vtol);
39: for (i = 0; i < (*event)->recsize; i++) PetscFree((*event)->recorder.eventidx[i]);
40: PetscFree((*event)->recorder.eventidx);
41: PetscFree((*event)->recorder.nevents);
42: PetscFree((*event)->recorder.stepnum);
43: PetscFree((*event)->recorder.time);
45: PetscViewerDestroy(&(*event)->monitor);
46: PetscFree(*event);
47: return 0;
48: }
50: /*@
51: TSSetPostEventIntervalStep - Set the time-step used immediately following the event interval
53: Logically Collective
55: Input Parameters:
56: + ts - time integration context
57: - dt - post event interval step
59: Options Database Key:
60: . -ts_event_post_eventinterval_step <dt> time-step after event interval
62: Level: advanced
64: Notes:
65: `TSSetPostEventIntervalStep()` allows one to set a time-step that is used immediately following an event interval.
67: This function should be called from the postevent function set with `TSSetEventHandler()`.
69: The post event interval time-step should be selected based on the dynamics following the event.
70: If the dynamics are stiff, a conservative (small) step should be used.
71: If not, then a larger time-step can be used.
73: .seealso: [](chapter_ts), `TS`, `TSEvent`, `TSSetEventHandler()`
74: @*/
75: PetscErrorCode TSSetPostEventIntervalStep(TS ts, PetscReal dt)
76: {
77: ts->event->timestep_posteventinterval = dt;
78: return 0;
79: }
81: /*@
82: TSSetEventTolerances - Set tolerances for event zero crossings when using event handler
84: Logically Collective
86: Input Parameters:
87: + ts - time integration context
88: . tol - scalar tolerance, `PETSC_DECIDE` to leave current value
89: - vtol - array of tolerances or NULL, used in preference to tol if present
91: Options Database Key:
92: . -ts_event_tol <tol> - tolerance for event zero crossing
94: Level: beginner
96: Notes:
97: Must call `TSSetEventHandler(`) before setting the tolerances.
99: The size of vtol is equal to the number of events.
101: .seealso: [](chapter_ts), `TS`, `TSEvent`, `TSSetEventHandler()`
102: @*/
103: PetscErrorCode TSSetEventTolerances(TS ts, PetscReal tol, PetscReal vtol[])
104: {
105: TSEvent event;
106: PetscInt i;
112: event = ts->event;
113: if (vtol) {
114: for (i = 0; i < event->nevents; i++) event->vtol[i] = vtol[i];
115: } else {
116: if (tol != PETSC_DECIDE || tol != PETSC_DEFAULT) {
117: for (i = 0; i < event->nevents; i++) event->vtol[i] = tol;
118: }
119: }
120: return 0;
121: }
123: /*@C
124: TSSetEventHandler - Sets a function used for detecting events
126: Logically Collective
128: Input Parameters:
129: + ts - the `TS` context obtained from `TSCreate()`
130: . nevents - number of local events
131: . direction - direction of zero crossing to be detected. -1 => Zero crossing in negative direction,
132: +1 => Zero crossing in positive direction, 0 => both ways (one for each event)
133: . terminate - flag to indicate whether time stepping should be terminated after
134: event is detected (one for each event)
135: . eventhandler - event monitoring routine
136: . postevent - [optional] post-event function
137: - ctx - [optional] user-defined context for private data for the
138: event monitor and post event routine (use NULL if no
139: context is desired)
141: Calling sequence of eventhandler:
142: PetscErrorCode PetscEventHandler(TS ts,PetscReal t,Vec U,PetscScalar fvalue[],void* ctx)
144: Input Parameters:
145: + ts - the TS context
146: . t - current time
147: . U - current iterate
148: - ctx - [optional] context passed with eventhandler
150: Output parameters:
151: . fvalue - function value of events at time t
153: Calling sequence of postevent:
154: PetscErrorCode PostEvent(TS ts,PetscInt nevents_zero,PetscInt events_zero[],PetscReal t,Vec U,PetscBool forwardsolve,void* ctx)
156: Input Parameters:
157: + ts - the TS context
158: . nevents_zero - number of local events whose event function is zero
159: . events_zero - indices of local events which have reached zero
160: . t - current time
161: . U - current solution
162: . forwardsolve - Flag to indicate whether TS is doing a forward solve (1) or adjoint solve (0)
163: - ctx - the context passed with eventhandler
165: Level: intermediate
167: .seealso: [](chapter_ts), `TSEvent`, `TSCreate()`, `TSSetTimeStep()`, `TSSetConvergedReason()`
168: @*/
169: PetscErrorCode TSSetEventHandler(TS ts, PetscInt nevents, PetscInt direction[], PetscBool terminate[], PetscErrorCode (*eventhandler)(TS, PetscReal, Vec, PetscScalar[], void *), PetscErrorCode (*postevent)(TS, PetscInt, PetscInt[], PetscReal, Vec, PetscBool, void *), void *ctx)
170: {
171: TSAdapt adapt;
172: PetscReal hmin;
173: TSEvent event;
174: PetscInt i;
175: PetscBool flg;
176: #if defined PETSC_USE_REAL_SINGLE
177: PetscReal tol = 1e-4;
178: #else
179: PetscReal tol = 1e-6;
180: #endif
183: if (nevents) {
186: }
188: PetscNew(&event);
189: PetscMalloc1(nevents, &event->fvalue);
190: PetscMalloc1(nevents, &event->fvalue_prev);
191: PetscMalloc1(nevents, &event->fvalue_right);
192: PetscMalloc1(nevents, &event->zerocrossing);
193: PetscMalloc1(nevents, &event->side);
194: PetscMalloc1(nevents, &event->direction);
195: PetscMalloc1(nevents, &event->terminate);
196: PetscMalloc1(nevents, &event->vtol);
197: for (i = 0; i < nevents; i++) {
198: event->direction[i] = direction[i];
199: event->terminate[i] = terminate[i];
200: event->zerocrossing[i] = PETSC_FALSE;
201: event->side[i] = 0;
202: }
203: PetscMalloc1(nevents, &event->events_zero);
204: event->nevents = nevents;
205: event->eventhandler = eventhandler;
206: event->postevent = postevent;
207: event->ctx = ctx;
208: event->timestep_posteventinterval = ts->time_step;
209: TSGetAdapt(ts, &adapt);
210: TSAdaptGetStepLimits(adapt, &hmin, NULL);
211: event->timestep_min = hmin;
213: event->recsize = 8; /* Initial size of the recorder */
214: PetscOptionsBegin(((PetscObject)ts)->comm, ((PetscObject)ts)->prefix, "TS Event options", "TS");
215: {
216: PetscOptionsReal("-ts_event_tol", "Scalar event tolerance for zero crossing check", "TSSetEventTolerances", tol, &tol, NULL);
217: PetscOptionsName("-ts_event_monitor", "Print choices made by event handler", "", &flg);
218: PetscOptionsInt("-ts_event_recorder_initial_size", "Initial size of event recorder", "", event->recsize, &event->recsize, NULL);
219: PetscOptionsReal("-ts_event_post_eventinterval_step", "Time step after event interval", "", event->timestep_posteventinterval, &event->timestep_posteventinterval, NULL);
220: PetscOptionsReal("-ts_event_post_event_step", "Time step after event", "", event->timestep_postevent, &event->timestep_postevent, NULL);
221: PetscOptionsReal("-ts_event_dt_min", "Minimum time step considered for TSEvent", "", event->timestep_min, &event->timestep_min, NULL);
222: }
223: PetscOptionsEnd();
225: PetscMalloc1(event->recsize, &event->recorder.time);
226: PetscMalloc1(event->recsize, &event->recorder.stepnum);
227: PetscMalloc1(event->recsize, &event->recorder.nevents);
228: PetscMalloc1(event->recsize, &event->recorder.eventidx);
229: for (i = 0; i < event->recsize; i++) PetscMalloc1(event->nevents, &event->recorder.eventidx[i]);
230: /* Initialize the event recorder */
231: event->recorder.ctr = 0;
233: for (i = 0; i < event->nevents; i++) event->vtol[i] = tol;
234: if (flg) PetscViewerASCIIOpen(PETSC_COMM_SELF, "stdout", &event->monitor);
236: TSEventDestroy(&ts->event);
237: ts->event = event;
238: ts->event->refct = 1;
239: return 0;
240: }
242: /*
243: TSEventRecorderResize - Resizes (2X) the event recorder arrays whenever the recording limit (event->recsize)
244: is reached.
245: */
246: static PetscErrorCode TSEventRecorderResize(TSEvent event)
247: {
248: PetscReal *time;
249: PetscInt *stepnum;
250: PetscInt *nevents;
251: PetscInt **eventidx;
252: PetscInt i, fact = 2;
255: /* Create large arrays */
256: PetscMalloc1(fact * event->recsize, &time);
257: PetscMalloc1(fact * event->recsize, &stepnum);
258: PetscMalloc1(fact * event->recsize, &nevents);
259: PetscMalloc1(fact * event->recsize, &eventidx);
260: for (i = 0; i < fact * event->recsize; i++) PetscMalloc1(event->nevents, &eventidx[i]);
262: /* Copy over data */
263: PetscArraycpy(time, event->recorder.time, event->recsize);
264: PetscArraycpy(stepnum, event->recorder.stepnum, event->recsize);
265: PetscArraycpy(nevents, event->recorder.nevents, event->recsize);
266: for (i = 0; i < event->recsize; i++) PetscArraycpy(eventidx[i], event->recorder.eventidx[i], event->recorder.nevents[i]);
268: /* Destroy old arrays */
269: for (i = 0; i < event->recsize; i++) PetscFree(event->recorder.eventidx[i]);
270: PetscFree(event->recorder.eventidx);
271: PetscFree(event->recorder.nevents);
272: PetscFree(event->recorder.stepnum);
273: PetscFree(event->recorder.time);
275: /* Set pointers */
276: event->recorder.time = time;
277: event->recorder.stepnum = stepnum;
278: event->recorder.nevents = nevents;
279: event->recorder.eventidx = eventidx;
281: /* Double size */
282: event->recsize *= fact;
284: return 0;
285: }
287: /*
288: Helper routine to handle user postevents and recording
289: */
290: static PetscErrorCode TSPostEvent(TS ts, PetscReal t, Vec U)
291: {
292: TSEvent event = ts->event;
293: PetscBool terminate = PETSC_FALSE;
294: PetscBool restart = PETSC_FALSE;
295: PetscInt i, ctr, stepnum;
296: PetscBool inflag[2], outflag[2];
297: PetscBool forwardsolve = PETSC_TRUE; /* Flag indicating that TS is doing a forward solve */
299: if (event->postevent) {
300: PetscObjectState state_prev, state_post;
301: PetscObjectStateGet((PetscObject)U, &state_prev);
302: (*event->postevent)(ts, event->nevents_zero, event->events_zero, t, U, forwardsolve, event->ctx);
303: PetscObjectStateGet((PetscObject)U, &state_post);
304: if (state_prev != state_post) restart = PETSC_TRUE;
305: }
307: /* Handle termination events and step restart */
308: for (i = 0; i < event->nevents_zero; i++)
309: if (event->terminate[event->events_zero[i]]) terminate = PETSC_TRUE;
310: inflag[0] = restart;
311: inflag[1] = terminate;
312: MPIU_Allreduce(inflag, outflag, 2, MPIU_BOOL, MPI_LOR, ((PetscObject)ts)->comm);
313: restart = outflag[0];
314: terminate = outflag[1];
315: if (restart) TSRestartStep(ts);
316: if (terminate) TSSetConvergedReason(ts, TS_CONVERGED_EVENT);
317: event->status = terminate ? TSEVENT_NONE : TSEVENT_RESET_NEXTSTEP;
319: /* Reset event residual functions as states might get changed by the postevent callback */
320: if (event->postevent) {
321: VecLockReadPush(U);
322: (*event->eventhandler)(ts, t, U, event->fvalue, event->ctx);
323: VecLockReadPop(U);
324: }
326: /* Cache current time and event residual functions */
327: event->ptime_prev = t;
328: for (i = 0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i];
330: /* Record the event in the event recorder */
331: TSGetStepNumber(ts, &stepnum);
332: ctr = event->recorder.ctr;
333: if (ctr == event->recsize) TSEventRecorderResize(event);
334: event->recorder.time[ctr] = t;
335: event->recorder.stepnum[ctr] = stepnum;
336: event->recorder.nevents[ctr] = event->nevents_zero;
337: for (i = 0; i < event->nevents_zero; i++) event->recorder.eventidx[ctr][i] = event->events_zero[i];
338: event->recorder.ctr++;
339: return 0;
340: }
342: /* Uses Anderson-Bjorck variant of regula falsi method */
343: static inline PetscReal TSEventComputeStepSize(PetscReal tleft, PetscReal t, PetscReal tright, PetscScalar fleft, PetscScalar f, PetscScalar fright, PetscInt side, PetscReal dt)
344: {
345: PetscReal new_dt, scal = 1.0;
346: if (PetscRealPart(fleft) * PetscRealPart(f) < 0) {
347: if (side == 1) {
348: scal = (PetscRealPart(fright) - PetscRealPart(f)) / PetscRealPart(fright);
349: if (scal < PETSC_SMALL) scal = 0.5;
350: }
351: new_dt = (scal * PetscRealPart(fleft) * t - PetscRealPart(f) * tleft) / (scal * PetscRealPart(fleft) - PetscRealPart(f)) - tleft;
352: } else {
353: if (side == -1) {
354: scal = (PetscRealPart(fleft) - PetscRealPart(f)) / PetscRealPart(fleft);
355: if (scal < PETSC_SMALL) scal = 0.5;
356: }
357: new_dt = (PetscRealPart(f) * tright - scal * PetscRealPart(fright) * t) / (PetscRealPart(f) - scal * PetscRealPart(fright)) - t;
358: }
359: return PetscMin(dt, new_dt);
360: }
362: static PetscErrorCode TSEventDetection(TS ts)
363: {
364: TSEvent event = ts->event;
365: PetscReal t;
366: PetscInt i;
367: PetscInt fvalue_sign, fvalueprev_sign;
368: PetscInt in, out;
370: TSGetTime(ts, &t);
371: for (i = 0; i < event->nevents; i++) {
372: if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) {
373: if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
374: event->status = TSEVENT_LOCATED_INTERVAL;
375: if (event->monitor) {
376: PetscViewerASCIIPrintf(event->monitor, "TSEvent: iter %" PetscInt_FMT " - Event %" PetscInt_FMT " interval detected due to zero value (tol=%g) [%g - %g]\n", event->iterctr, i, (double)event->vtol[i], (double)event->ptime_prev, (double)t);
377: }
378: continue;
379: }
380: if (PetscAbsScalar(event->fvalue_prev[i]) < event->vtol[i]) continue; /* avoid duplicative detection if the previous endpoint is an event location */
381: fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
382: fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
383: if (fvalueprev_sign != 0 && (fvalue_sign != fvalueprev_sign)) {
384: if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
385: event->status = TSEVENT_LOCATED_INTERVAL;
386: if (event->monitor) PetscViewerASCIIPrintf(event->monitor, "TSEvent: iter %" PetscInt_FMT " - Event %" PetscInt_FMT " interval detected due to sign change [%g - %g]\n", event->iterctr, i, (double)event->ptime_prev, (double)t);
387: }
388: }
389: in = (PetscInt)event->status;
390: MPIU_Allreduce(&in, &out, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)ts));
391: event->status = (TSEventStatus)out;
392: return 0;
393: }
395: static PetscErrorCode TSEventLocation(TS ts, PetscReal *dt)
396: {
397: TSEvent event = ts->event;
398: PetscInt i;
399: PetscReal t;
400: PetscInt fvalue_sign, fvalueprev_sign;
401: PetscInt rollback = 0, in[2], out[2];
403: TSGetTime(ts, &t);
404: event->nevents_zero = 0;
405: for (i = 0; i < event->nevents; i++) {
406: if (event->zerocrossing[i]) {
407: if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i] || *dt < event->timestep_min || PetscAbsReal((*dt) / ((event->ptime_right - event->ptime_prev) / 2)) < event->vtol[i]) { /* stopping criteria */
408: event->status = TSEVENT_ZERO;
409: event->fvalue_right[i] = event->fvalue[i];
410: continue;
411: }
412: /* Compute new time step */
413: *dt = TSEventComputeStepSize(event->ptime_prev, t, event->ptime_right, event->fvalue_prev[i], event->fvalue[i], event->fvalue_right[i], event->side[i], *dt);
414: fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
415: fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
416: switch (event->direction[i]) {
417: case -1:
418: if (fvalue_sign < 0) {
419: rollback = 1;
420: event->fvalue_right[i] = event->fvalue[i];
421: event->side[i] = 1;
422: }
423: break;
424: case 1:
425: if (fvalue_sign > 0) {
426: rollback = 1;
427: event->fvalue_right[i] = event->fvalue[i];
428: event->side[i] = 1;
429: }
430: break;
431: case 0:
432: if (fvalue_sign != fvalueprev_sign) { /* trigger rollback only when there is a sign change */
433: rollback = 1;
434: event->fvalue_right[i] = event->fvalue[i];
435: event->side[i] = 1;
436: }
437: break;
438: }
439: if (event->status == TSEVENT_PROCESSING) event->side[i] = -1;
440: }
441: }
442: in[0] = (PetscInt)event->status;
443: in[1] = rollback;
444: MPIU_Allreduce(in, out, 2, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)ts));
445: event->status = (TSEventStatus)out[0];
446: rollback = out[1];
447: /* If rollback is true, the status will be overwritten so that an event at the endtime of current time step will be postponed to guarantee correct order */
448: if (rollback) event->status = TSEVENT_LOCATED_INTERVAL;
449: if (event->status == TSEVENT_ZERO) {
450: for (i = 0; i < event->nevents; i++) {
451: if (event->zerocrossing[i]) {
452: if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i] || *dt < event->timestep_min || PetscAbsReal((*dt) / ((event->ptime_right - event->ptime_prev) / 2)) < event->vtol[i]) { /* stopping criteria */
453: event->events_zero[event->nevents_zero++] = i;
454: if (event->monitor) PetscViewerASCIIPrintf(event->monitor, "TSEvent: iter %" PetscInt_FMT " - Event %" PetscInt_FMT " zero crossing located at time %g\n", event->iterctr, i, (double)t);
455: event->zerocrossing[i] = PETSC_FALSE;
456: }
457: }
458: event->side[i] = 0;
459: }
460: }
461: return 0;
462: }
464: PetscErrorCode TSEventHandler(TS ts)
465: {
466: TSEvent event;
467: PetscReal t;
468: Vec U;
469: PetscInt i;
470: PetscReal dt, dt_min, dt_reset = 0.0;
473: if (!ts->event) return 0;
474: event = ts->event;
476: TSGetTime(ts, &t);
477: TSGetTimeStep(ts, &dt);
478: TSGetSolution(ts, &U);
480: if (event->status == TSEVENT_NONE) {
481: event->timestep_prev = dt;
482: event->ptime_end = t;
483: }
484: if (event->status == TSEVENT_RESET_NEXTSTEP) {
485: /* user has specified a PostEventInterval dt */
486: dt = event->timestep_posteventinterval;
487: if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
488: PetscReal maxdt = ts->max_time - t;
489: dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt, maxdt, 10 * PETSC_MACHINE_EPSILON, 0) ? maxdt : dt);
490: }
491: TSSetTimeStep(ts, dt);
492: event->status = TSEVENT_NONE;
493: }
495: VecLockReadPush(U);
496: (*event->eventhandler)(ts, t, U, event->fvalue, event->ctx);
497: VecLockReadPop(U);
499: /* Detect the events */
500: TSEventDetection(ts);
502: /* Locate the events */
503: if (event->status == TSEVENT_LOCATED_INTERVAL || event->status == TSEVENT_PROCESSING) {
504: /* Approach the zero crosing by setting a new step size */
505: TSEventLocation(ts, &dt);
506: /* Roll back when new events are detected */
507: if (event->status == TSEVENT_LOCATED_INTERVAL) {
508: TSRollBack(ts);
509: TSSetConvergedReason(ts, TS_CONVERGED_ITERATING);
510: event->iterctr++;
511: }
512: MPIU_Allreduce(&dt, &dt_min, 1, MPIU_REAL, MPIU_MIN, PetscObjectComm((PetscObject)ts));
513: if (dt_reset > 0.0 && dt_reset < dt_min) dt_min = dt_reset;
514: TSSetTimeStep(ts, dt_min);
515: /* Found the zero crossing */
516: if (event->status == TSEVENT_ZERO) {
517: TSPostEvent(ts, t, U);
519: dt = event->ptime_end - t;
520: if (PetscAbsReal(dt) < PETSC_SMALL) { /* we hit the event, continue with the candidate time step */
521: dt = event->timestep_prev;
522: event->status = TSEVENT_NONE;
523: }
524: if (event->timestep_postevent) { /* user has specified a PostEvent dt*/
525: dt = event->timestep_postevent;
526: }
527: if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
528: PetscReal maxdt = ts->max_time - t;
529: dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt, maxdt, 10 * PETSC_MACHINE_EPSILON, 0) ? maxdt : dt);
530: }
531: TSSetTimeStep(ts, dt);
532: event->iterctr = 0;
533: }
534: /* Have not found the zero crosing yet */
535: if (event->status == TSEVENT_PROCESSING) {
536: if (event->monitor) PetscViewerASCIIPrintf(event->monitor, "TSEvent: iter %" PetscInt_FMT " - Stepping forward as no event detected in interval [%g - %g]\n", event->iterctr, (double)event->ptime_prev, (double)t);
537: event->iterctr++;
538: }
539: }
540: if (event->status == TSEVENT_LOCATED_INTERVAL) { /* The step has been rolled back */
541: event->status = TSEVENT_PROCESSING;
542: event->ptime_right = t;
543: } else {
544: for (i = 0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i];
545: event->ptime_prev = t;
546: }
547: return 0;
548: }
550: PetscErrorCode TSAdjointEventHandler(TS ts)
551: {
552: TSEvent event;
553: PetscReal t;
554: Vec U;
555: PetscInt ctr;
556: PetscBool forwardsolve = PETSC_FALSE; /* Flag indicating that TS is doing an adjoint solve */
559: if (!ts->event) return 0;
560: event = ts->event;
562: TSGetTime(ts, &t);
563: TSGetSolution(ts, &U);
565: ctr = event->recorder.ctr - 1;
566: if (ctr >= 0 && PetscAbsReal(t - event->recorder.time[ctr]) < PETSC_SMALL) {
567: /* Call the user postevent function */
568: if (event->postevent) {
569: (*event->postevent)(ts, event->recorder.nevents[ctr], event->recorder.eventidx[ctr], t, U, forwardsolve, event->ctx);
570: event->recorder.ctr--;
571: }
572: }
574: PetscBarrier((PetscObject)ts);
575: return 0;
576: }
578: /*@
579: TSGetNumEvents - Get the numbers of events set
581: Logically Collective
583: Input Parameter:
584: . ts - the `TS` context
586: Output Parameter:
587: . nevents - number of events
589: Level: intermediate
591: .seealso: [](chapter_ts), `TSEvent`, `TSSetEventHandler()`
592: @*/
593: PetscErrorCode TSGetNumEvents(TS ts, PetscInt *nevents)
594: {
595: *nevents = ts->event->nevents;
596: return 0;
597: }