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 Keys:
60: . -ts_event_post_eventinterval_step <dt> time-step after event interval
62: Notes:
63: TSSetPostEventIntervalStep allows one to set a time-step that is used immediately following an event interval.
65: This function should be called from the postevent function set with TSSetEventHandler().
67: The post event interval time-step should be selected based on the dynamics following the event.
68: If the dynamics are stiff, a conservative (small) step should be used.
69: If not, then a larger time-step can be used.
71: Level: Advanced
72: .seealso: `TS`, `TSEvent`, `TSSetEventHandler()`
73: @*/
74: PetscErrorCode TSSetPostEventIntervalStep(TS ts, PetscReal dt)
75: {
76: ts->event->timestep_posteventinterval = dt;
77: return 0;
78: }
80: /*@
81: TSSetEventTolerances - Set tolerances for event zero crossings when using event handler
83: Logically Collective
85: Input Parameters:
86: + ts - time integration context
87: . tol - scalar tolerance, PETSC_DECIDE to leave current value
88: - vtol - array of tolerances or NULL, used in preference to tol if present
90: Options Database Keys:
91: . -ts_event_tol <tol> - tolerance for event zero crossing
93: Notes:
94: Must call TSSetEventHandler() before setting the tolerances.
96: The size of vtol is equal to the number of events.
98: Level: beginner
100: .seealso: `TS`, `TSEvent`, `TSSetEventHandler()`
101: @*/
102: PetscErrorCode TSSetEventTolerances(TS ts, PetscReal tol, PetscReal vtol[])
103: {
104: TSEvent event;
105: PetscInt i;
111: event = ts->event;
112: if (vtol) {
113: for (i = 0; i < event->nevents; i++) event->vtol[i] = vtol[i];
114: } else {
115: if (tol != PETSC_DECIDE || tol != PETSC_DEFAULT) {
116: for (i = 0; i < event->nevents; i++) event->vtol[i] = tol;
117: }
118: }
119: return 0;
120: }
122: /*@C
123: TSSetEventHandler - Sets a function used for detecting events
125: Logically Collective on TS
127: Input Parameters:
128: + ts - the TS context obtained from TSCreate()
129: . nevents - number of local events
130: . direction - direction of zero crossing to be detected. -1 => Zero crossing in negative direction,
131: +1 => Zero crossing in positive direction, 0 => both ways (one for each event)
132: . terminate - flag to indicate whether time stepping should be terminated after
133: event is detected (one for each event)
134: . eventhandler - event monitoring routine
135: . postevent - [optional] post-event function
136: - ctx - [optional] user-defined context for private data for the
137: event monitor and post event routine (use NULL if no
138: context is desired)
140: Calling sequence of eventhandler:
141: PetscErrorCode PetscEventHandler(TS ts,PetscReal t,Vec U,PetscScalar fvalue[],void* ctx)
143: Input Parameters:
144: + ts - the TS context
145: . t - current time
146: . U - current iterate
147: - ctx - [optional] context passed with eventhandler
149: Output parameters:
150: . fvalue - function value of events at time t
152: Calling sequence of postevent:
153: PetscErrorCode PostEvent(TS ts,PetscInt nevents_zero,PetscInt events_zero[],PetscReal t,Vec U,PetscBool forwardsolve,void* ctx)
155: Input Parameters:
156: + ts - the TS context
157: . nevents_zero - number of local events whose event function is zero
158: . events_zero - indices of local events which have reached zero
159: . t - current time
160: . U - current solution
161: . forwardsolve - Flag to indicate whether TS is doing a forward solve (1) or adjoint solve (0)
162: - ctx - the context passed with eventhandler
164: Level: intermediate
166: .seealso: `TSCreate()`, `TSSetTimeStep()`, `TSSetConvergedReason()`
167: @*/
168: 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)
169: {
170: TSAdapt adapt;
171: PetscReal hmin;
172: TSEvent event;
173: PetscInt i;
174: PetscBool flg;
175: #if defined PETSC_USE_REAL_SINGLE
176: PetscReal tol = 1e-4;
177: #else
178: PetscReal tol = 1e-6;
179: #endif
182: if (nevents) {
185: }
187: PetscNew(&event);
188: PetscMalloc1(nevents, &event->fvalue);
189: PetscMalloc1(nevents, &event->fvalue_prev);
190: PetscMalloc1(nevents, &event->fvalue_right);
191: PetscMalloc1(nevents, &event->zerocrossing);
192: PetscMalloc1(nevents, &event->side);
193: PetscMalloc1(nevents, &event->direction);
194: PetscMalloc1(nevents, &event->terminate);
195: PetscMalloc1(nevents, &event->vtol);
196: for (i = 0; i < nevents; i++) {
197: event->direction[i] = direction[i];
198: event->terminate[i] = terminate[i];
199: event->zerocrossing[i] = PETSC_FALSE;
200: event->side[i] = 0;
201: }
202: PetscMalloc1(nevents, &event->events_zero);
203: event->nevents = nevents;
204: event->eventhandler = eventhandler;
205: event->postevent = postevent;
206: event->ctx = ctx;
207: event->timestep_posteventinterval = ts->time_step;
208: TSGetAdapt(ts, &adapt);
209: TSAdaptGetStepLimits(adapt, &hmin, NULL);
210: event->timestep_min = hmin;
212: event->recsize = 8; /* Initial size of the recorder */
213: PetscOptionsBegin(((PetscObject)ts)->comm, ((PetscObject)ts)->prefix, "TS Event options", "TS");
214: {
215: PetscOptionsReal("-ts_event_tol", "Scalar event tolerance for zero crossing check", "TSSetEventTolerances", tol, &tol, NULL);
216: PetscOptionsName("-ts_event_monitor", "Print choices made by event handler", "", &flg);
217: PetscOptionsInt("-ts_event_recorder_initial_size", "Initial size of event recorder", "", event->recsize, &event->recsize, NULL);
218: PetscOptionsReal("-ts_event_post_eventinterval_step", "Time step after event interval", "", event->timestep_posteventinterval, &event->timestep_posteventinterval, NULL);
219: PetscOptionsReal("-ts_event_post_event_step", "Time step after event", "", event->timestep_postevent, &event->timestep_postevent, NULL);
220: PetscOptionsReal("-ts_event_dt_min", "Minimum time step considered for TSEvent", "", event->timestep_min, &event->timestep_min, NULL);
221: }
222: PetscOptionsEnd();
224: PetscMalloc1(event->recsize, &event->recorder.time);
225: PetscMalloc1(event->recsize, &event->recorder.stepnum);
226: PetscMalloc1(event->recsize, &event->recorder.nevents);
227: PetscMalloc1(event->recsize, &event->recorder.eventidx);
228: for (i = 0; i < event->recsize; i++) PetscMalloc1(event->nevents, &event->recorder.eventidx[i]);
229: /* Initialize the event recorder */
230: event->recorder.ctr = 0;
232: for (i = 0; i < event->nevents; i++) event->vtol[i] = tol;
233: if (flg) PetscViewerASCIIOpen(PETSC_COMM_SELF, "stdout", &event->monitor);
235: TSEventDestroy(&ts->event);
236: ts->event = event;
237: ts->event->refct = 1;
238: return 0;
239: }
241: /*
242: TSEventRecorderResize - Resizes (2X) the event recorder arrays whenever the recording limit (event->recsize)
243: is reached.
244: */
245: static PetscErrorCode TSEventRecorderResize(TSEvent event)
246: {
247: PetscReal *time;
248: PetscInt *stepnum;
249: PetscInt *nevents;
250: PetscInt **eventidx;
251: PetscInt i, fact = 2;
254: /* Create large arrays */
255: PetscMalloc1(fact * event->recsize, &time);
256: PetscMalloc1(fact * event->recsize, &stepnum);
257: PetscMalloc1(fact * event->recsize, &nevents);
258: PetscMalloc1(fact * event->recsize, &eventidx);
259: for (i = 0; i < fact * event->recsize; i++) PetscMalloc1(event->nevents, &eventidx[i]);
261: /* Copy over data */
262: PetscArraycpy(time, event->recorder.time, event->recsize);
263: PetscArraycpy(stepnum, event->recorder.stepnum, event->recsize);
264: PetscArraycpy(nevents, event->recorder.nevents, event->recsize);
265: for (i = 0; i < event->recsize; i++) PetscArraycpy(eventidx[i], event->recorder.eventidx[i], event->recorder.nevents[i]);
267: /* Destroy old arrays */
268: for (i = 0; i < event->recsize; i++) PetscFree(event->recorder.eventidx[i]);
269: PetscFree(event->recorder.eventidx);
270: PetscFree(event->recorder.nevents);
271: PetscFree(event->recorder.stepnum);
272: PetscFree(event->recorder.time);
274: /* Set pointers */
275: event->recorder.time = time;
276: event->recorder.stepnum = stepnum;
277: event->recorder.nevents = nevents;
278: event->recorder.eventidx = eventidx;
280: /* Double size */
281: event->recsize *= fact;
283: return 0;
284: }
286: /*
287: Helper routine to handle user postevents and recording
288: */
289: static PetscErrorCode TSPostEvent(TS ts, PetscReal t, Vec U)
290: {
291: TSEvent event = ts->event;
292: PetscBool terminate = PETSC_FALSE;
293: PetscBool restart = PETSC_FALSE;
294: PetscInt i, ctr, stepnum;
295: PetscBool inflag[2], outflag[2];
296: PetscBool forwardsolve = PETSC_TRUE; /* Flag indicating that TS is doing a forward solve */
298: if (event->postevent) {
299: PetscObjectState state_prev, state_post;
300: PetscObjectStateGet((PetscObject)U, &state_prev);
301: (*event->postevent)(ts, event->nevents_zero, event->events_zero, t, U, forwardsolve, event->ctx);
302: PetscObjectStateGet((PetscObject)U, &state_post);
303: if (state_prev != state_post) restart = PETSC_TRUE;
304: }
306: /* Handle termination events and step restart */
307: for (i = 0; i < event->nevents_zero; i++)
308: if (event->terminate[event->events_zero[i]]) terminate = PETSC_TRUE;
309: inflag[0] = restart;
310: inflag[1] = terminate;
311: MPIU_Allreduce(inflag, outflag, 2, MPIU_BOOL, MPI_LOR, ((PetscObject)ts)->comm);
312: restart = outflag[0];
313: terminate = outflag[1];
314: if (restart) TSRestartStep(ts);
315: if (terminate) TSSetConvergedReason(ts, TS_CONVERGED_EVENT);
316: event->status = terminate ? TSEVENT_NONE : TSEVENT_RESET_NEXTSTEP;
318: /* Reset event residual functions as states might get changed by the postevent callback */
319: if (event->postevent) {
320: VecLockReadPush(U);
321: (*event->eventhandler)(ts, t, U, event->fvalue, event->ctx);
322: VecLockReadPop(U);
323: }
325: /* Cache current time and event residual functions */
326: event->ptime_prev = t;
327: for (i = 0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i];
329: /* Record the event in the event recorder */
330: TSGetStepNumber(ts, &stepnum);
331: ctr = event->recorder.ctr;
332: if (ctr == event->recsize) TSEventRecorderResize(event);
333: event->recorder.time[ctr] = t;
334: event->recorder.stepnum[ctr] = stepnum;
335: event->recorder.nevents[ctr] = event->nevents_zero;
336: for (i = 0; i < event->nevents_zero; i++) event->recorder.eventidx[ctr][i] = event->events_zero[i];
337: event->recorder.ctr++;
338: return 0;
339: }
341: /* Uses Anderson-Bjorck variant of regula falsi method */
342: static inline PetscReal TSEventComputeStepSize(PetscReal tleft, PetscReal t, PetscReal tright, PetscScalar fleft, PetscScalar f, PetscScalar fright, PetscInt side, PetscReal dt)
343: {
344: PetscReal new_dt, scal = 1.0;
345: if (PetscRealPart(fleft) * PetscRealPart(f) < 0) {
346: if (side == 1) {
347: scal = (PetscRealPart(fright) - PetscRealPart(f)) / PetscRealPart(fright);
348: if (scal < PETSC_SMALL) scal = 0.5;
349: }
350: new_dt = (scal * PetscRealPart(fleft) * t - PetscRealPart(f) * tleft) / (scal * PetscRealPart(fleft) - PetscRealPart(f)) - tleft;
351: } else {
352: if (side == -1) {
353: scal = (PetscRealPart(fleft) - PetscRealPart(f)) / PetscRealPart(fleft);
354: if (scal < PETSC_SMALL) scal = 0.5;
355: }
356: new_dt = (PetscRealPart(f) * tright - scal * PetscRealPart(fright) * t) / (PetscRealPart(f) - scal * PetscRealPart(fright)) - t;
357: }
358: return PetscMin(dt, new_dt);
359: }
361: static PetscErrorCode TSEventDetection(TS ts)
362: {
363: TSEvent event = ts->event;
364: PetscReal t;
365: PetscInt i;
366: PetscInt fvalue_sign, fvalueprev_sign;
367: PetscInt in, out;
369: TSGetTime(ts, &t);
370: for (i = 0; i < event->nevents; i++) {
371: if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) {
372: if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
373: event->status = TSEVENT_LOCATED_INTERVAL;
374: if (event->monitor) {
375: 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);
376: }
377: continue;
378: }
379: if (PetscAbsScalar(event->fvalue_prev[i]) < event->vtol[i]) continue; /* avoid duplicative detection if the previous endpoint is an event location */
380: fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
381: fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
382: if (fvalueprev_sign != 0 && (fvalue_sign != fvalueprev_sign)) {
383: if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
384: event->status = TSEVENT_LOCATED_INTERVAL;
385: 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);
386: }
387: }
388: in = (PetscInt)event->status;
389: MPIU_Allreduce(&in, &out, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)ts));
390: event->status = (TSEventStatus)out;
391: return 0;
392: }
394: static PetscErrorCode TSEventLocation(TS ts, PetscReal *dt)
395: {
396: TSEvent event = ts->event;
397: PetscInt i;
398: PetscReal t;
399: PetscInt fvalue_sign, fvalueprev_sign;
400: PetscInt rollback = 0, in[2], out[2];
402: TSGetTime(ts, &t);
403: event->nevents_zero = 0;
404: for (i = 0; i < event->nevents; i++) {
405: if (event->zerocrossing[i]) {
406: 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 */
407: event->status = TSEVENT_ZERO;
408: event->fvalue_right[i] = event->fvalue[i];
409: continue;
410: }
411: /* Compute new time step */
412: *dt = TSEventComputeStepSize(event->ptime_prev, t, event->ptime_right, event->fvalue_prev[i], event->fvalue[i], event->fvalue_right[i], event->side[i], *dt);
413: fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
414: fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
415: switch (event->direction[i]) {
416: case -1:
417: if (fvalue_sign < 0) {
418: rollback = 1;
419: event->fvalue_right[i] = event->fvalue[i];
420: event->side[i] = 1;
421: }
422: break;
423: case 1:
424: if (fvalue_sign > 0) {
425: rollback = 1;
426: event->fvalue_right[i] = event->fvalue[i];
427: event->side[i] = 1;
428: }
429: break;
430: case 0:
431: if (fvalue_sign != fvalueprev_sign) { /* trigger rollback only when there is a sign change */
432: rollback = 1;
433: event->fvalue_right[i] = event->fvalue[i];
434: event->side[i] = 1;
435: }
436: break;
437: }
438: if (event->status == TSEVENT_PROCESSING) event->side[i] = -1;
439: }
440: }
441: in[0] = (PetscInt)event->status;
442: in[1] = rollback;
443: MPIU_Allreduce(in, out, 2, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)ts));
444: event->status = (TSEventStatus)out[0];
445: rollback = out[1];
446: /* 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 corret order */
447: if (rollback) event->status = TSEVENT_LOCATED_INTERVAL;
448: if (event->status == TSEVENT_ZERO) {
449: for (i = 0; i < event->nevents; i++) {
450: if (event->zerocrossing[i]) {
451: 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 */
452: event->events_zero[event->nevents_zero++] = i;
453: 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);
454: event->zerocrossing[i] = PETSC_FALSE;
455: }
456: }
457: event->side[i] = 0;
458: }
459: }
460: return 0;
461: }
463: PetscErrorCode TSEventHandler(TS ts)
464: {
465: TSEvent event;
466: PetscReal t;
467: Vec U;
468: PetscInt i;
469: PetscReal dt, dt_min, dt_reset = 0.0;
472: if (!ts->event) return 0;
473: event = ts->event;
475: TSGetTime(ts, &t);
476: TSGetTimeStep(ts, &dt);
477: TSGetSolution(ts, &U);
479: if (event->status == TSEVENT_NONE) {
480: event->timestep_prev = dt;
481: event->ptime_end = t;
482: }
483: if (event->status == TSEVENT_RESET_NEXTSTEP) {
484: /* user has specified a PostEventInterval dt */
485: dt = event->timestep_posteventinterval;
486: if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
487: PetscReal maxdt = ts->max_time - t;
488: dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt, maxdt, 10 * PETSC_MACHINE_EPSILON, 0) ? maxdt : dt);
489: }
490: TSSetTimeStep(ts, dt);
491: event->status = TSEVENT_NONE;
492: }
494: VecLockReadPush(U);
495: (*event->eventhandler)(ts, t, U, event->fvalue, event->ctx);
496: VecLockReadPop(U);
498: /* Detect the events */
499: TSEventDetection(ts);
501: /* Locate the events */
502: if (event->status == TSEVENT_LOCATED_INTERVAL || event->status == TSEVENT_PROCESSING) {
503: /* Approach the zero crosing by setting a new step size */
504: TSEventLocation(ts, &dt);
505: /* Roll back when new events are detected */
506: if (event->status == TSEVENT_LOCATED_INTERVAL) {
507: TSRollBack(ts);
508: TSSetConvergedReason(ts, TS_CONVERGED_ITERATING);
509: event->iterctr++;
510: }
511: MPIU_Allreduce(&dt, &dt_min, 1, MPIU_REAL, MPIU_MIN, PetscObjectComm((PetscObject)ts));
512: if (dt_reset > 0.0 && dt_reset < dt_min) dt_min = dt_reset;
513: TSSetTimeStep(ts, dt_min);
514: /* Found the zero crossing */
515: if (event->status == TSEVENT_ZERO) {
516: TSPostEvent(ts, t, U);
518: dt = event->ptime_end - t;
519: if (PetscAbsReal(dt) < PETSC_SMALL) { /* we hit the event, continue with the candidate time step */
520: dt = event->timestep_prev;
521: event->status = TSEVENT_NONE;
522: }
523: if (event->timestep_postevent) { /* user has specified a PostEvent dt*/
524: dt = event->timestep_postevent;
525: }
526: if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
527: PetscReal maxdt = ts->max_time - t;
528: dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt, maxdt, 10 * PETSC_MACHINE_EPSILON, 0) ? maxdt : dt);
529: }
530: TSSetTimeStep(ts, dt);
531: event->iterctr = 0;
532: }
533: /* Have not found the zero crosing yet */
534: if (event->status == TSEVENT_PROCESSING) {
535: 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);
536: event->iterctr++;
537: }
538: }
539: if (event->status == TSEVENT_LOCATED_INTERVAL) { /* The step has been rolled back */
540: event->status = TSEVENT_PROCESSING;
541: event->ptime_right = t;
542: } else {
543: for (i = 0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i];
544: event->ptime_prev = t;
545: }
546: return 0;
547: }
549: PetscErrorCode TSAdjointEventHandler(TS ts)
550: {
551: TSEvent event;
552: PetscReal t;
553: Vec U;
554: PetscInt ctr;
555: PetscBool forwardsolve = PETSC_FALSE; /* Flag indicating that TS is doing an adjoint solve */
558: if (!ts->event) return 0;
559: event = ts->event;
561: TSGetTime(ts, &t);
562: TSGetSolution(ts, &U);
564: ctr = event->recorder.ctr - 1;
565: if (ctr >= 0 && PetscAbsReal(t - event->recorder.time[ctr]) < PETSC_SMALL) {
566: /* Call the user postevent function */
567: if (event->postevent) {
568: (*event->postevent)(ts, event->recorder.nevents[ctr], event->recorder.eventidx[ctr], t, U, forwardsolve, event->ctx);
569: event->recorder.ctr--;
570: }
571: }
573: PetscBarrier((PetscObject)ts);
574: return 0;
575: }
577: /*@
578: TSGetNumEvents - Get the numbers of events set
580: Logically Collective
582: Input Parameter:
583: . ts - the TS context
585: Output Parameter:
586: . nevents - number of events
588: Level: intermediate
590: .seealso: `TSSetEventHandler()`
592: @*/
593: PetscErrorCode TSGetNumEvents(TS ts, PetscInt *nevents)
594: {
595: *nevents = ts->event->nevents;
596: return 0;
597: }