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) {*event = NULL; return 0;}

 26:   PetscFree((*event)->fvalue);
 27:   PetscFree((*event)->fvalue_prev);
 28:   PetscFree((*event)->fvalue_right);
 29:   PetscFree((*event)->zerocrossing);
 30:   PetscFree((*event)->side);
 31:   PetscFree((*event)->direction);
 32:   PetscFree((*event)->terminate);
 33:   PetscFree((*event)->events_zero);
 34:   PetscFree((*event)->vtol);

 36:   for (i=0; i < (*event)->recsize; i++) {
 37:     PetscFree((*event)->recorder.eventidx[i]);
 38:   }
 39:   PetscFree((*event)->recorder.eventidx);
 40:   PetscFree((*event)->recorder.nevents);
 41:   PetscFree((*event)->recorder.stepnum);
 42:   PetscFree((*event)->recorder.time);

 44:   PetscViewerDestroy(&(*event)->monitor);
 45:   PetscFree(*event);
 46:   return 0;
 47: }

 49: /*@
 50:   TSSetPostEventIntervalStep - Set the time-step used immediately following the event interval

 52:   Logically Collective

 54:   Input Parameters:
 55: + ts - time integration context
 56: - dt - post event interval step

 58:   Options Database Keys:
 59: . -ts_event_post_eventinterval_step <dt> time-step after event interval

 61:   Notes:
 62:   TSSetPostEventIntervalStep allows one to set a time-step that is used immediately following an event interval.

 64:   This function should be called from the postevent function set with TSSetEventHandler().

 66:   The post event interval time-step should be selected based on the dynamics following the event.
 67:   If the dynamics are stiff, a conservative (small) step should be used.
 68:   If not, then a larger time-step can be used.

 70:   Level: Advanced
 71:   .seealso: TS, TSEvent, TSSetEventHandler()
 72: @*/
 73: PetscErrorCode TSSetPostEventIntervalStep(TS ts,PetscReal dt)
 74: {
 75:   ts->event->timestep_posteventinterval = dt;
 76:   return 0;
 77: }

 79: /*@
 80:    TSSetEventTolerances - Set tolerances for event zero crossings when using event handler

 82:    Logically Collective

 84:    Input Parameters:
 85: +  ts - time integration context
 86: .  tol - scalar tolerance, PETSC_DECIDE to leave current value
 87: -  vtol - array of tolerances or NULL, used in preference to tol if present

 89:    Options Database Keys:
 90: .  -ts_event_tol <tol> - tolerance for event zero crossing

 92:    Notes:
 93:    Must call TSSetEventHandler() before setting the tolerances.

 95:    The size of vtol is equal to the number of events.

 97:    Level: beginner

 99: .seealso: TS, TSEvent, TSSetEventHandler()
100: @*/
101: PetscErrorCode TSSetEventTolerances(TS ts,PetscReal tol,PetscReal vtol[])
102: {
103:   TSEvent        event;
104:   PetscInt       i;


110:   event = ts->event;
111:   if (vtol) {
112:     for (i=0; i < event->nevents; i++) event->vtol[i] = vtol[i];
113:   } else {
114:     if (tol != PETSC_DECIDE || tol != PETSC_DEFAULT) {
115:       for (i=0; i < event->nevents; i++) event->vtol[i] = tol;
116:     }
117:   }
118:   return 0;
119: }

121: /*@C
122:    TSSetEventHandler - Sets a function used for detecting events

124:    Logically Collective on TS

126:    Input Parameters:
127: +  ts - the TS context obtained from TSCreate()
128: .  nevents - number of local events
129: .  direction - direction of zero crossing to be detected. -1 => Zero crossing in negative direction,
130:                +1 => Zero crossing in positive direction, 0 => both ways (one for each event)
131: .  terminate - flag to indicate whether time stepping should be terminated after
132:                event is detected (one for each event)
133: .  eventhandler - event monitoring routine
134: .  postevent - [optional] post-event function
135: -  ctx       - [optional] user-defined context for private data for the
136:                event monitor and post event routine (use NULL if no
137:                context is desired)

139:    Calling sequence of eventhandler:
140:    PetscErrorCode PetscEventHandler(TS ts,PetscReal t,Vec U,PetscScalar fvalue[],void* ctx)

142:    Input Parameters:
143: +  ts  - the TS context
144: .  t   - current time
145: .  U   - current iterate
146: -  ctx - [optional] context passed with eventhandler

148:    Output parameters:
149: .  fvalue    - function value of events at time t

151:    Calling sequence of postevent:
152:    PetscErrorCode PostEvent(TS ts,PetscInt nevents_zero,PetscInt events_zero[],PetscReal t,Vec U,PetscBool forwardsolve,void* ctx)

154:    Input Parameters:
155: +  ts - the TS context
156: .  nevents_zero - number of local events whose event function is zero
157: .  events_zero  - indices of local events which have reached zero
158: .  t            - current time
159: .  U            - current solution
160: .  forwardsolve - Flag to indicate whether TS is doing a forward solve (1) or adjoint solve (0)
161: -  ctx          - the context passed with eventhandler

163:    Level: intermediate

165: .seealso: TSCreate(), TSSetTimeStep(), TSSetConvergedReason()
166: @*/
167: 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)
168: {
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:   PetscNewLog(ts,&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++) {
229:     PetscMalloc1(event->nevents,&event->recorder.eventidx[i]);
230:   }
231:   /* Initialize the event recorder */
232:   event->recorder.ctr = 0;

234:   for (i=0; i < event->nevents; i++) event->vtol[i] = tol;
235:   if (flg) PetscViewerASCIIOpen(PETSC_COMM_SELF,"stdout",&event->monitor);

237:   TSEventDestroy(&ts->event);
238:   ts->event = event;
239:   ts->event->refct = 1;
240:   return 0;
241: }

243: /*
244:   TSEventRecorderResize - Resizes (2X) the event recorder arrays whenever the recording limit (event->recsize)
245:                           is reached.
246: */
247: static PetscErrorCode TSEventRecorderResize(TSEvent event)
248: {
249:   PetscReal      *time;
250:   PetscInt       *stepnum;
251:   PetscInt       *nevents;
252:   PetscInt       **eventidx;
253:   PetscInt       i,fact=2;


256:   /* Create large arrays */
257:   PetscMalloc1(fact*event->recsize,&time);
258:   PetscMalloc1(fact*event->recsize,&stepnum);
259:   PetscMalloc1(fact*event->recsize,&nevents);
260:   PetscMalloc1(fact*event->recsize,&eventidx);
261:   for (i=0; i < fact*event->recsize; i++) {
262:     PetscMalloc1(event->nevents,&eventidx[i]);
263:   }

265:   /* Copy over data */
266:   PetscArraycpy(time,event->recorder.time,event->recsize);
267:   PetscArraycpy(stepnum,event->recorder.stepnum,event->recsize);
268:   PetscArraycpy(nevents,event->recorder.nevents,event->recsize);
269:   for (i=0; i < event->recsize; i++) {
270:     PetscArraycpy(eventidx[i],event->recorder.eventidx[i],event->recorder.nevents[i]);
271:   }

273:   /* Destroy old arrays */
274:   for (i=0; i < event->recsize; i++) {
275:     PetscFree(event->recorder.eventidx[i]);
276:   }
277:   PetscFree(event->recorder.eventidx);
278:   PetscFree(event->recorder.nevents);
279:   PetscFree(event->recorder.stepnum);
280:   PetscFree(event->recorder.time);

282:   /* Set pointers */
283:   event->recorder.time = time;
284:   event->recorder.stepnum = stepnum;
285:   event->recorder.nevents = nevents;
286:   event->recorder.eventidx = eventidx;

288:   /* Double size */
289:   event->recsize *= fact;

291:   return 0;
292: }

294: /*
295:    Helper routine to handle user postevents and recording
296: */
297: static PetscErrorCode TSPostEvent(TS ts,PetscReal t,Vec U)
298: {
299:   TSEvent        event = ts->event;
300:   PetscBool      terminate = PETSC_FALSE;
301:   PetscBool      restart = PETSC_FALSE;
302:   PetscInt       i,ctr,stepnum;
303:   PetscBool      inflag[2],outflag[2];
304:   PetscBool      forwardsolve = PETSC_TRUE; /* Flag indicating that TS is doing a forward solve */

306:   if (event->postevent) {
307:     PetscObjectState state_prev,state_post;
308:     PetscObjectStateGet((PetscObject)U,&state_prev);
309:     (*event->postevent)(ts,event->nevents_zero,event->events_zero,t,U,forwardsolve,event->ctx);
310:     PetscObjectStateGet((PetscObject)U,&state_post);
311:     if (state_prev != state_post) restart = PETSC_TRUE;
312:   }

314:   /* Handle termination events and step restart */
315:   for (i=0; i<event->nevents_zero; i++) if (event->terminate[event->events_zero[i]]) terminate = PETSC_TRUE;
316:   inflag[0] = restart; inflag[1] = terminate;
317:   MPIU_Allreduce(inflag,outflag,2,MPIU_BOOL,MPI_LOR,((PetscObject)ts)->comm);
318:   restart = outflag[0]; terminate = outflag[1];
319:   if (restart) TSRestartStep(ts);
320:   if (terminate) TSSetConvergedReason(ts,TS_CONVERGED_EVENT);
321:   event->status = terminate ? TSEVENT_NONE : TSEVENT_RESET_NEXTSTEP;

323:   /* Reset event residual functions as states might get changed by the postevent callback */
324:   if (event->postevent) {
325:     VecLockReadPush(U);
326:     (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);
327:     VecLockReadPop(U);
328:   }

330:   /* Cache current time and event residual functions */
331:   event->ptime_prev = t;
332:   for (i=0; i<event->nevents; i++)
333:     event->fvalue_prev[i] = event->fvalue[i];

335:   /* Record the event in the event recorder */
336:   TSGetStepNumber(ts,&stepnum);
337:   ctr = event->recorder.ctr;
338:   if (ctr == event->recsize) {
339:     TSEventRecorderResize(event);
340:   }
341:   event->recorder.time[ctr] = t;
342:   event->recorder.stepnum[ctr] = stepnum;
343:   event->recorder.nevents[ctr] = event->nevents_zero;
344:   for (i=0; i<event->nevents_zero; i++) event->recorder.eventidx[ctr][i] = event->events_zero[i];
345:   event->recorder.ctr++;
346:   return 0;
347: }

349: /* Uses Anderson-Bjorck variant of regula falsi method */
350: static inline PetscReal TSEventComputeStepSize(PetscReal tleft,PetscReal t,PetscReal tright,PetscScalar fleft,PetscScalar f,PetscScalar fright,PetscInt side,PetscReal dt)
351: {
352:   PetscReal new_dt, scal = 1.0;
353:   if (PetscRealPart(fleft)*PetscRealPart(f) < 0) {
354:     if (side == 1) {
355:       scal = (PetscRealPart(fright) - PetscRealPart(f))/PetscRealPart(fright);
356:       if (scal < PETSC_SMALL) scal = 0.5;
357:     }
358:     new_dt = (scal*PetscRealPart(fleft)*t - PetscRealPart(f)*tleft)/(scal*PetscRealPart(fleft) - PetscRealPart(f)) - tleft;
359:   } else {
360:     if (side == -1) {
361:       scal = (PetscRealPart(fleft) - PetscRealPart(f))/PetscRealPart(fleft);
362:       if (scal < PETSC_SMALL) scal = 0.5;
363:     }
364:     new_dt = (PetscRealPart(f)*tright - scal*PetscRealPart(fright)*t)/(PetscRealPart(f) - scal*PetscRealPart(fright)) - t;
365:   }
366:   return PetscMin(dt,new_dt);
367: }

369: static PetscErrorCode TSEventDetection(TS ts)
370: {
371:   TSEvent        event = ts->event;
372:   PetscReal      t;
373:   PetscInt       i;
374:   PetscInt       fvalue_sign,fvalueprev_sign;
375:   PetscInt       in,out;

377:   TSGetTime(ts,&t);
378:   for (i=0; i < event->nevents; i++) {
379:     if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) {
380:       if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
381:       event->status = TSEVENT_LOCATED_INTERVAL;
382:       if (event->monitor) {
383:         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);
384:       }
385:       continue;
386:     }
387:     if (PetscAbsScalar(event->fvalue_prev[i]) < event->vtol[i]) continue; /* avoid duplicative detection if the previous endpoint is an event location */
388:     fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
389:     fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
390:     if (fvalueprev_sign != 0 && (fvalue_sign != fvalueprev_sign)) {
391:       if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
392:       event->status = TSEVENT_LOCATED_INTERVAL;
393:       if (event->monitor) {
394:         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);
395:       }
396:     }
397:   }
398:   in = (PetscInt)event->status;
399:   MPIU_Allreduce(&in,&out,1,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)ts));
400:   event->status = (TSEventStatus)out;
401:   return 0;
402: }

404: static PetscErrorCode TSEventLocation(TS ts,PetscReal *dt)
405: {
406:   TSEvent        event = ts->event;
407:   PetscInt       i;
408:   PetscReal      t;
409:   PetscInt       fvalue_sign,fvalueprev_sign;
410:   PetscInt       rollback=0,in[2],out[2];

412:   TSGetTime(ts,&t);
413:   event->nevents_zero = 0;
414:   for (i=0; i < event->nevents; i++) {
415:     if (event->zerocrossing[i]) {
416:       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 */
417:         event->status = TSEVENT_ZERO;
418:         event->fvalue_right[i] = event->fvalue[i];
419:         continue;
420:       }
421:       /* Compute new time step */
422:       *dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],*dt);
423:       fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
424:       fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
425:       switch (event->direction[i]) {
426:       case -1:
427:         if (fvalue_sign < 0) {
428:           rollback = 1;
429:           event->fvalue_right[i] = event->fvalue[i];
430:           event->side[i] = 1;
431:         }
432:         break;
433:       case 1:
434:         if (fvalue_sign > 0) {
435:           rollback = 1;
436:           event->fvalue_right[i] = event->fvalue[i];
437:           event->side[i] = 1;
438:         }
439:         break;
440:       case 0:
441:         if (fvalue_sign != fvalueprev_sign) { /* trigger rollback only when there is a sign change */
442:           rollback = 1;
443:           event->fvalue_right[i] = event->fvalue[i];
444:           event->side[i] = 1;
445:         }
446:         break;
447:       }
448:       if (event->status == TSEVENT_PROCESSING) event->side[i] = -1;
449:     }
450:   }
451:   in[0] = (PetscInt)event->status; in[1] = rollback;
452:   MPIU_Allreduce(in,out,2,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)ts));
453:   event->status = (TSEventStatus)out[0]; rollback = out[1];
454:   /* 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 */
455:   if (rollback) event->status = TSEVENT_LOCATED_INTERVAL;
456:   if (event->status == TSEVENT_ZERO) {
457:     for (i=0; i < event->nevents; i++) {
458:       if (event->zerocrossing[i]) {
459:         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 */
460:           event->events_zero[event->nevents_zero++] = i;
461:           if (event->monitor) {
462:             PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Event %D zero crossing located at time %g\n",event->iterctr,i,(double)t);
463:           }
464:           event->zerocrossing[i] = PETSC_FALSE;
465:         }
466:       }
467:       event->side[i] = 0;
468:     }
469:   }
470:   return 0;
471: }

473: PetscErrorCode TSEventHandler(TS ts)
474: {
475:   TSEvent        event;
476:   PetscReal      t;
477:   Vec            U;
478:   PetscInt       i;
479:   PetscReal      dt,dt_min,dt_reset = 0.0;

482:   if (!ts->event) return 0;
483:   event = ts->event;

485:   TSGetTime(ts,&t);
486:   TSGetTimeStep(ts,&dt);
487:   TSGetSolution(ts,&U);

489:   if (event->status == TSEVENT_NONE) {
490:     event->timestep_prev = dt;
491:     event->ptime_end = t;
492:   }
493:   if (event->status == TSEVENT_RESET_NEXTSTEP) {
494:     /* user has specified a PostEventInterval dt */
495:     dt = event->timestep_posteventinterval;
496:     if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
497:       PetscReal maxdt = ts->max_time-t;
498:       dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt,maxdt,10*PETSC_MACHINE_EPSILON,0) ? maxdt : dt);
499:     }
500:     TSSetTimeStep(ts,dt);
501:     event->status = TSEVENT_NONE;
502:   }

504:   VecLockReadPush(U);
505:   (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);
506:   VecLockReadPop(U);

508:   /* Detect the events */
509:   TSEventDetection(ts);

511:   /* Locate the events */
512:   if (event->status == TSEVENT_LOCATED_INTERVAL || event->status == TSEVENT_PROCESSING) {
513:     /* Approach the zero crosing by setting a new step size */
514:     TSEventLocation(ts,&dt);
515:     /* Roll back when new events are detected */
516:     if (event->status == TSEVENT_LOCATED_INTERVAL) {
517:       TSRollBack(ts);
518:       TSSetConvergedReason(ts,TS_CONVERGED_ITERATING);
519:       event->iterctr++;
520:     }
521:     MPIU_Allreduce(&dt,&dt_min,1,MPIU_REAL,MPIU_MIN,PetscObjectComm((PetscObject)ts));
522:     if (dt_reset > 0.0 && dt_reset < dt_min) dt_min = dt_reset;
523:     TSSetTimeStep(ts,dt_min);
524:     /* Found the zero crossing */
525:     if (event->status == TSEVENT_ZERO) {
526:       TSPostEvent(ts,t,U);

528:       dt = event->ptime_end - t;
529:       if (PetscAbsReal(dt) < PETSC_SMALL) { /* we hit the event, continue with the candidate time step */
530:         dt = event->timestep_prev;
531:         event->status = TSEVENT_NONE;
532:       }
533:       if (event->timestep_postevent) { /* user has specified a PostEvent dt*/
534:         dt = event->timestep_postevent;
535:       }
536:       if (ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
537:         PetscReal maxdt = ts->max_time-t;
538:         dt = dt > maxdt ? maxdt : (PetscIsCloseAtTol(dt,maxdt,10*PETSC_MACHINE_EPSILON,0) ? maxdt : dt);
539:       }
540:       TSSetTimeStep(ts,dt);
541:       event->iterctr = 0;
542:     }
543:     /* Have not found the zero crosing yet */
544:     if (event->status == TSEVENT_PROCESSING) {
545:       if (event->monitor) {
546:         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);
547:       }
548:       event->iterctr++;
549:     }
550:   }
551:   if (event->status == TSEVENT_LOCATED_INTERVAL) { /* The step has been rolled back */
552:     event->status = TSEVENT_PROCESSING;
553:     event->ptime_right = t;
554:   } else {
555:     for (i=0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i];
556:     event->ptime_prev = t;
557:   }
558:   return 0;
559: }

561: PetscErrorCode TSAdjointEventHandler(TS ts)
562: {
563:   TSEvent        event;
564:   PetscReal      t;
565:   Vec            U;
566:   PetscInt       ctr;
567:   PetscBool      forwardsolve=PETSC_FALSE; /* Flag indicating that TS is doing an adjoint solve */

570:   if (!ts->event) return 0;
571:   event = ts->event;

573:   TSGetTime(ts,&t);
574:   TSGetSolution(ts,&U);

576:   ctr = event->recorder.ctr-1;
577:   if (ctr >= 0 && PetscAbsReal(t - event->recorder.time[ctr]) < PETSC_SMALL) {
578:     /* Call the user postevent function */
579:     if (event->postevent) {
580:       (*event->postevent)(ts,event->recorder.nevents[ctr],event->recorder.eventidx[ctr],t,U,forwardsolve,event->ctx);
581:       event->recorder.ctr--;
582:     }
583:   }

585:   PetscBarrier((PetscObject)ts);
586:   return 0;
587: }

589: /*@
590:   TSGetNumEvents - Get the numbers of events set

592:   Logically Collective

594:   Input Parameter:
595: . ts - the TS context

597:   Output Parameter:
598: . nevents - number of events

600:   Level: intermediate

602: .seealso: TSSetEventHandler()

604: @*/
605: PetscErrorCode TSGetNumEvents(TS ts,PetscInt * nevents)
606: {
607:   *nevents = ts->event->nevents;
608:   return 0;
609: }