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: {
11: if (!event) return(0);
15: event->ptime_prev = t;
16: event->iterctr = 0;
17: (*event->eventhandler)(ts,t,U,event->fvalue_prev,event->ctx);
18: return(0);
19: }
21: PetscErrorCode TSEventDestroy(TSEvent *event)
22: {
24: PetscInt i;
28: if (!*event) return(0);
29: if (--(*event)->refct > 0) {*event = NULL; return(0);}
31: PetscFree((*event)->fvalue);
32: PetscFree((*event)->fvalue_prev);
33: PetscFree((*event)->fvalue_right);
34: PetscFree((*event)->zerocrossing);
35: PetscFree((*event)->side);
36: PetscFree((*event)->direction);
37: PetscFree((*event)->terminate);
38: PetscFree((*event)->events_zero);
39: PetscFree((*event)->vtol);
41: for (i=0; i < (*event)->recsize; i++) {
42: PetscFree((*event)->recorder.eventidx[i]);
43: }
44: PetscFree((*event)->recorder.eventidx);
45: PetscFree((*event)->recorder.nevents);
46: PetscFree((*event)->recorder.stepnum);
47: PetscFree((*event)->recorder.time);
49: PetscViewerDestroy(&(*event)->monitor);
50: PetscFree(*event);
51: return(0);
52: }
54: /*@
55: TSSetPostEventIntervalStep - Set the time-step used immediately following the event interval
57: Logically Collective
59: Input Parameters:
60: + ts - time integration context
61: - dt - post event interval step
63: Options Database Keys:
64: . -ts_event_post_eventinterval_step <dt> time-step after event interval
66: Notes:
67: TSSetPostEventIntervalStep allows one to set a time-step that is used immediately following an event interval.
69: This function should be called from the postevent function set with TSSetEventHandler().
71: The post event interval time-step should be selected based on the dynamics following the event.
72: If the dynamics are stiff, a conservative (small) step should be used.
73: If not, then a larger time-step can be used.
75: Level: Advanced
76: .seealso: TS, TSEvent, TSSetEventHandler()
77: @*/
78: PetscErrorCode TSSetPostEventIntervalStep(TS ts,PetscReal dt)
79: {
81: ts->event->timestep_posteventinterval = dt;
82: return(0);
83: }
85: /*@
86: TSSetEventTolerances - Set tolerances for event zero crossings when using event handler
88: Logically Collective
90: Input Parameters:
91: + ts - time integration context
92: . tol - scalar tolerance, PETSC_DECIDE to leave current value
93: - vtol - array of tolerances or NULL, used in preference to tol if present
95: Options Database Keys:
96: . -ts_event_tol <tol> tolerance for event zero crossing
98: Notes:
99: Must call TSSetEventHandler() before setting the tolerances.
101: The size of vtol is equal to the number of events.
103: Level: beginner
105: .seealso: TS, TSEvent, TSSetEventHandler()
106: @*/
107: PetscErrorCode TSSetEventTolerances(TS ts,PetscReal tol,PetscReal vtol[])
108: {
109: TSEvent event;
110: PetscInt i;
115: if (!ts->event) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_USER,"Must set the events first by calling TSSetEventHandler()");
117: event = ts->event;
118: if (vtol) {
119: for (i=0; i < event->nevents; i++) event->vtol[i] = vtol[i];
120: } else {
121: if (tol != PETSC_DECIDE || tol != PETSC_DEFAULT) {
122: for (i=0; i < event->nevents; i++) event->vtol[i] = tol;
123: }
124: }
125: return(0);
126: }
128: /*@C
129: TSSetEventHandler - Sets a function used for detecting events
131: Logically Collective on TS
133: Input Parameters:
134: + ts - the TS context obtained from TSCreate()
135: . nevents - number of local events
136: . direction - direction of zero crossing to be detected. -1 => Zero crossing in negative direction,
137: +1 => Zero crossing in positive direction, 0 => both ways (one for each event)
138: . terminate - flag to indicate whether time stepping should be terminated after
139: event is detected (one for each event)
140: . eventhandler - event monitoring routine
141: . postevent - [optional] post-event function
142: - ctx - [optional] user-defined context for private data for the
143: event monitor and post event routine (use NULL if no
144: context is desired)
146: Calling sequence of eventhandler:
147: PetscErrorCode PetscEventHandler(TS ts,PetscReal t,Vec U,PetscScalar fvalue[],void* ctx)
149: Input Parameters:
150: + ts - the TS context
151: . t - current time
152: . U - current iterate
153: - ctx - [optional] context passed with eventhandler
155: Output parameters:
156: . fvalue - function value of events at time t
158: Calling sequence of postevent:
159: PetscErrorCode PostEvent(TS ts,PetscInt nevents_zero,PetscInt events_zero[],PetscReal t,Vec U,PetscBool forwardsolve,void* ctx)
161: Input Parameters:
162: + ts - the TS context
163: . nevents_zero - number of local events whose event function is zero
164: . events_zero - indices of local events which have reached zero
165: . t - current time
166: . U - current solution
167: . forwardsolve - Flag to indicate whether TS is doing a forward solve (1) or adjoint solve (0)
168: - ctx - the context passed with eventhandler
170: Level: intermediate
172: .seealso: TSCreate(), TSSetTimeStep(), TSSetConvergedReason()
173: @*/
174: 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)
175: {
177: TSAdapt adapt;
178: PetscReal hmin;
179: TSEvent event;
180: PetscInt i;
181: PetscBool flg;
182: #if defined PETSC_USE_REAL_SINGLE
183: PetscReal tol=1e-4;
184: #else
185: PetscReal tol=1e-6;
186: #endif
190: if (nevents) {
193: }
195: PetscNewLog(ts,&event);
196: PetscMalloc1(nevents,&event->fvalue);
197: PetscMalloc1(nevents,&event->fvalue_prev);
198: PetscMalloc1(nevents,&event->fvalue_right);
199: PetscMalloc1(nevents,&event->zerocrossing);
200: PetscMalloc1(nevents,&event->side);
201: PetscMalloc1(nevents,&event->direction);
202: PetscMalloc1(nevents,&event->terminate);
203: PetscMalloc1(nevents,&event->vtol);
204: for (i=0; i < nevents; i++) {
205: event->direction[i] = direction[i];
206: event->terminate[i] = terminate[i];
207: event->zerocrossing[i] = PETSC_FALSE;
208: event->side[i] = 0;
209: }
210: PetscMalloc1(nevents,&event->events_zero);
211: event->nevents = nevents;
212: event->eventhandler = eventhandler;
213: event->postevent = postevent;
214: event->ctx = ctx;
215: event->timestep_posteventinterval = ts->time_step;
216: TSGetAdapt(ts,&adapt);
217: TSAdaptGetStepLimits(adapt,&hmin,NULL);
218: event->timestep_min = hmin;
220: event->recsize = 8; /* Initial size of the recorder */
221: PetscOptionsBegin(((PetscObject)ts)->comm,((PetscObject)ts)->prefix,"TS Event options","TS");
222: {
223: PetscOptionsReal("-ts_event_tol","Scalar event tolerance for zero crossing check","TSSetEventTolerances",tol,&tol,NULL);
224: PetscOptionsName("-ts_event_monitor","Print choices made by event handler","",&flg);
225: PetscOptionsInt("-ts_event_recorder_initial_size","Initial size of event recorder","",event->recsize,&event->recsize,NULL);
226: PetscOptionsReal("-ts_event_post_eventinterval_step","Time step after event interval","",event->timestep_posteventinterval,&event->timestep_posteventinterval,NULL);
227: PetscOptionsReal("-ts_event_post_event_step","Time step after event","",event->timestep_postevent,&event->timestep_postevent,NULL);
228: PetscOptionsReal("-ts_event_dt_min","Minimum time step considered for TSEvent","",event->timestep_min,&event->timestep_min,NULL);
229: }
230: PetscOptionsEnd();
232: PetscMalloc1(event->recsize,&event->recorder.time);
233: PetscMalloc1(event->recsize,&event->recorder.stepnum);
234: PetscMalloc1(event->recsize,&event->recorder.nevents);
235: PetscMalloc1(event->recsize,&event->recorder.eventidx);
236: for (i=0; i < event->recsize; i++) {
237: PetscMalloc1(event->nevents,&event->recorder.eventidx[i]);
238: }
239: /* Initialize the event recorder */
240: event->recorder.ctr = 0;
242: for (i=0; i < event->nevents; i++) event->vtol[i] = tol;
243: if (flg) {PetscViewerASCIIOpen(PETSC_COMM_SELF,"stdout",&event->monitor);}
245: TSEventDestroy(&ts->event);
246: ts->event = event;
247: ts->event->refct = 1;
248: return(0);
249: }
251: /*
252: TSEventRecorderResize - Resizes (2X) the event recorder arrays whenever the recording limit (event->recsize)
253: is reached.
254: */
255: static PetscErrorCode TSEventRecorderResize(TSEvent event)
256: {
258: PetscReal *time;
259: PetscInt *stepnum;
260: PetscInt *nevents;
261: PetscInt **eventidx;
262: PetscInt i,fact=2;
266: /* Create large arrays */
267: PetscMalloc1(fact*event->recsize,&time);
268: PetscMalloc1(fact*event->recsize,&stepnum);
269: PetscMalloc1(fact*event->recsize,&nevents);
270: PetscMalloc1(fact*event->recsize,&eventidx);
271: for (i=0; i < fact*event->recsize; i++) {
272: PetscMalloc1(event->nevents,&eventidx[i]);
273: }
275: /* Copy over data */
276: PetscArraycpy(time,event->recorder.time,event->recsize);
277: PetscArraycpy(stepnum,event->recorder.stepnum,event->recsize);
278: PetscArraycpy(nevents,event->recorder.nevents,event->recsize);
279: for (i=0; i < event->recsize; i++) {
280: PetscArraycpy(eventidx[i],event->recorder.eventidx[i],event->recorder.nevents[i]);
281: }
283: /* Destroy old arrays */
284: for (i=0; i < event->recsize; i++) {
285: PetscFree(event->recorder.eventidx[i]);
286: }
287: PetscFree(event->recorder.eventidx);
288: PetscFree(event->recorder.nevents);
289: PetscFree(event->recorder.stepnum);
290: PetscFree(event->recorder.time);
292: /* Set pointers */
293: event->recorder.time = time;
294: event->recorder.stepnum = stepnum;
295: event->recorder.nevents = nevents;
296: event->recorder.eventidx = eventidx;
298: /* Double size */
299: event->recsize *= fact;
301: return(0);
302: }
304: /*
305: Helper routine to handle user postevents and recording
306: */
307: static PetscErrorCode TSPostEvent(TS ts,PetscReal t,Vec U)
308: {
310: TSEvent event = ts->event;
311: PetscBool terminate = PETSC_FALSE;
312: PetscBool restart = PETSC_FALSE;
313: PetscInt i,ctr,stepnum;
314: PetscBool inflag[2],outflag[2];
315: PetscBool forwardsolve = PETSC_TRUE; /* Flag indicating that TS is doing a forward solve */
318: if (event->postevent) {
319: PetscObjectState state_prev,state_post;
320: PetscObjectStateGet((PetscObject)U,&state_prev);
321: (*event->postevent)(ts,event->nevents_zero,event->events_zero,t,U,forwardsolve,event->ctx);
322: PetscObjectStateGet((PetscObject)U,&state_post);
323: if (state_prev != state_post) restart = PETSC_TRUE;
324: }
326: /* Handle termination events and step restart */
327: for (i=0; i<event->nevents_zero; i++) if (event->terminate[event->events_zero[i]]) terminate = PETSC_TRUE;
328: inflag[0] = restart; inflag[1] = terminate;
329: MPIU_Allreduce(inflag,outflag,2,MPIU_BOOL,MPI_LOR,((PetscObject)ts)->comm);
330: restart = outflag[0]; terminate = outflag[1];
331: if (restart) {TSRestartStep(ts);}
332: if (terminate) {TSSetConvergedReason(ts,TS_CONVERGED_EVENT);}
333: event->status = terminate ? TSEVENT_NONE : TSEVENT_RESET_NEXTSTEP;
335: /* Reset event residual functions as states might get changed by the postevent callback */
336: if (event->postevent) {
337: VecLockReadPush(U);
338: (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);
339: VecLockReadPop(U);
340: }
342: /* Cache current time and event residual functions */
343: event->ptime_prev = t;
344: for (i=0; i<event->nevents; i++)
345: event->fvalue_prev[i] = event->fvalue[i];
347: /* Record the event in the event recorder */
348: TSGetStepNumber(ts,&stepnum);
349: ctr = event->recorder.ctr;
350: if (ctr == event->recsize) {
351: TSEventRecorderResize(event);
352: }
353: event->recorder.time[ctr] = t;
354: event->recorder.stepnum[ctr] = stepnum;
355: event->recorder.nevents[ctr] = event->nevents_zero;
356: for (i=0; i<event->nevents_zero; i++) event->recorder.eventidx[ctr][i] = event->events_zero[i];
357: event->recorder.ctr++;
358: return(0);
359: }
361: /* Uses Anderson-Bjorck variant of regula falsi method */
362: PETSC_STATIC_INLINE PetscReal TSEventComputeStepSize(PetscReal tleft,PetscReal t,PetscReal tright,PetscScalar fleft,PetscScalar f,PetscScalar fright,PetscInt side,PetscReal dt)
363: {
364: PetscReal new_dt, scal = 1.0;
365: if (PetscRealPart(fleft)*PetscRealPart(f) < 0) {
366: if (side == 1) {
367: scal = (PetscRealPart(fright) - PetscRealPart(f))/PetscRealPart(fright);
368: if (scal < PETSC_SMALL) scal = 0.5;
369: }
370: new_dt = (scal*PetscRealPart(fleft)*t - PetscRealPart(f)*tleft)/(scal*PetscRealPart(fleft) - PetscRealPart(f)) - tleft;
371: } else {
372: if (side == -1) {
373: scal = (PetscRealPart(fleft) - PetscRealPart(f))/PetscRealPart(fleft);
374: if (scal < PETSC_SMALL) scal = 0.5;
375: }
376: new_dt = (PetscRealPart(f)*tright - scal*PetscRealPart(fright)*t)/(PetscRealPart(f) - scal*PetscRealPart(fright)) - t;
377: }
378: return PetscMin(dt,new_dt);
379: }
381: static PetscErrorCode TSEventDetection(TS ts)
382: {
384: TSEvent event = ts->event;
385: PetscReal t;
386: PetscInt i;
387: PetscInt fvalue_sign,fvalueprev_sign;
388: PetscInt in,out;
391: TSGetTime(ts,&t);
392: for (i=0; i < event->nevents; i++) {
393: if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) {
394: if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
395: event->status = TSEVENT_LOCATED_INTERVAL;
396: if (event->monitor) {
397: PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Event %D interval detected due to zero value (tol=%g) [%g - %g]\n",event->iterctr,i,(double)event->vtol[i],(double)event->ptime_prev,(double)t);
398: }
399: continue;
400: }
401: if (PetscAbsScalar(event->fvalue_prev[i]) < event->vtol[i]) continue; /* avoid duplicative detection if the previous endpoint is an event location */
402: fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
403: fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
404: if (fvalueprev_sign != 0 && (fvalue_sign != fvalueprev_sign)) {
405: if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
406: event->status = TSEVENT_LOCATED_INTERVAL;
407: if (event->monitor) {
408: PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Event %D interval detected due to sign change [%g - %g]\n",event->iterctr,i,(double)event->ptime_prev,(double)t);
409: }
410: }
411: }
412: in = event->status;
413: MPIU_Allreduce(&in,&out,1,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)ts));
414: event->status = (TSEventStatus)out;
415: return(0);
416: }
418: static PetscErrorCode TSEventLocation(TS ts,PetscReal *dt)
419: {
421: TSEvent event = ts->event;
422: PetscInt i;
423: PetscReal t;
424: PetscInt fvalue_sign,fvalueprev_sign;
425: PetscInt rollback=0,in[2],out[2];
428: TSGetTime(ts,&t);
429: event->nevents_zero = 0;
430: for (i=0; i < event->nevents; i++) {
431: if (event->zerocrossing[i]) {
432: 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 */
433: event->status = TSEVENT_ZERO;
434: event->fvalue_right[i] = event->fvalue[i];
435: continue;
436: }
437: /* Compute new time step */
438: *dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],*dt);
439: fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
440: fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
441: switch (event->direction[i]) {
442: case -1:
443: if (fvalue_sign < 0) {
444: rollback = 1;
445: event->fvalue_right[i] = event->fvalue[i];
446: event->side[i] = 1;
447: }
448: break;
449: case 1:
450: if (fvalue_sign > 0) {
451: rollback = 1;
452: event->fvalue_right[i] = event->fvalue[i];
453: event->side[i] = 1;
454: }
455: break;
456: case 0:
457: if (fvalue_sign != fvalueprev_sign) { /* trigger rollback only when there is a sign change */
458: rollback = 1;
459: event->fvalue_right[i] = event->fvalue[i];
460: event->side[i] = 1;
461: }
462: break;
463: }
464: if (event->status == TSEVENT_PROCESSING) event->side[i] = -1;
465: }
466: }
467: in[0] = event->status; in[1] = rollback;
468: MPIU_Allreduce(in,out,2,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)ts));
469: event->status = (TSEventStatus)out[0]; rollback = out[1];
470: /* 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 */
471: if (rollback) event->status = TSEVENT_LOCATED_INTERVAL;
472: if (event->status == TSEVENT_ZERO) {
473: for (i=0; i < event->nevents; i++) {
474: if (event->zerocrossing[i]) {
475: 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 */
476: event->events_zero[event->nevents_zero++] = i;
477: if (event->monitor) {
478: PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Event %D zero crossing located at time %g\n",event->iterctr,i,(double)t);
479: }
480: event->zerocrossing[i] = PETSC_FALSE;
481: }
482: }
483: event->side[i] = 0;
484: }
485: }
486: return(0);
487: }
489: PetscErrorCode TSEventHandler(TS ts)
490: {
492: TSEvent event;
493: PetscReal t;
494: Vec U;
495: PetscInt i;
496: PetscReal dt,dt_min,dt_reset = 0.0;
500: if (!ts->event) return(0);
501: event = ts->event;
503: TSGetTime(ts,&t);
504: TSGetTimeStep(ts,&dt);
505: TSGetSolution(ts,&U);
507: if (event->status == TSEVENT_NONE) {
508: event->timestep_prev = dt;
509: event->ptime_end = t;
510: }
511: if (event->status == TSEVENT_RESET_NEXTSTEP) {
512: /* user has specified a PostEventInterval dt */
513: dt = event->timestep_posteventinterval;
514: if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
515: PetscReal maxdt = ts->max_time-t;
516: dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt,maxdt,10*PETSC_MACHINE_EPSILON,0) ? maxdt : dt);
517: }
518: TSSetTimeStep(ts,dt);
519: event->status = TSEVENT_NONE;
520: }
522: VecLockReadPush(U);
523: (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);
524: VecLockReadPop(U);
526: /* Detect the events */
527: TSEventDetection(ts);
529: /* Locate the events */
530: if (event->status == TSEVENT_LOCATED_INTERVAL || event->status == TSEVENT_PROCESSING) {
531: /* Approach the zero crosing by setting a new step size */
532: TSEventLocation(ts,&dt);
533: /* Roll back when new events are detected */
534: if (event->status == TSEVENT_LOCATED_INTERVAL) {
535: TSRollBack(ts);
536: TSSetConvergedReason(ts,TS_CONVERGED_ITERATING);
537: event->iterctr++;
538: }
539: MPIU_Allreduce(&dt,&dt_min,1,MPIU_REAL,MPIU_MIN,PetscObjectComm((PetscObject)ts));
540: if (dt_reset > 0.0 && dt_reset < dt_min) dt_min = dt_reset;
541: TSSetTimeStep(ts,dt_min);
542: /* Found the zero crossing */
543: if (event->status == TSEVENT_ZERO) {
544: TSPostEvent(ts,t,U);
546: dt = event->ptime_end - t;
547: if (PetscAbsReal(dt) < PETSC_SMALL) { /* we hit the event, continue with the candidate time step */
548: dt = event->timestep_prev;
549: event->status = TSEVENT_NONE;
550: }
551: if (event->timestep_postevent) { /* user has specified a PostEvent dt*/
552: dt = event->timestep_postevent;
553: }
554: if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
555: PetscReal maxdt = ts->max_time-t;
556: dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt,maxdt,10*PETSC_MACHINE_EPSILON,0) ? maxdt : dt);
557: }
558: TSSetTimeStep(ts,dt);
559: event->iterctr = 0;
560: }
561: /* Have not found the zero crosing yet */
562: if (event->status == TSEVENT_PROCESSING) {
563: if (event->monitor) {
564: PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Stepping forward as no event detected in interval [%g - %g]\n",event->iterctr,(double)event->ptime_prev,(double)t);
565: }
566: event->iterctr++;
567: }
568: }
569: if (event->status == TSEVENT_LOCATED_INTERVAL) { /* The step has been rolled back */
570: event->status = TSEVENT_PROCESSING;
571: event->ptime_right = t;
572: } else {
573: for (i=0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i];
574: event->ptime_prev = t;
575: }
576: return(0);
577: }
579: PetscErrorCode TSAdjointEventHandler(TS ts)
580: {
582: TSEvent event;
583: PetscReal t;
584: Vec U;
585: PetscInt ctr;
586: PetscBool forwardsolve=PETSC_FALSE; /* Flag indicating that TS is doing an adjoint solve */
590: if (!ts->event) return(0);
591: event = ts->event;
593: TSGetTime(ts,&t);
594: TSGetSolution(ts,&U);
596: ctr = event->recorder.ctr-1;
597: if (ctr >= 0 && PetscAbsReal(t - event->recorder.time[ctr]) < PETSC_SMALL) {
598: /* Call the user postevent function */
599: if (event->postevent) {
600: (*event->postevent)(ts,event->recorder.nevents[ctr],event->recorder.eventidx[ctr],t,U,forwardsolve,event->ctx);
601: event->recorder.ctr--;
602: }
603: }
605: PetscBarrier((PetscObject)ts);
606: return(0);
607: }
609: /*@
610: TSGetNumEvents - Get the numbers of events set
612: Logically Collective
614: Input Parameter:
615: . ts - the TS context
617: Output Parameter:
618: . nevents - number of events
620: Level: intermediate
622: .seealso: TSSetEventHandler()
624: @*/
625: PetscErrorCode TSGetNumEvents(TS ts,PetscInt * nevents)
626: {
628: *nevents = ts->event->nevents;
629: return(0);
630: }