xref: /petsc/src/ts/event/tsevent.c (revision 4fc747eaadbeca11629f314a99edccbc2ed7b3d3)
1 #include <petsc/private/tsimpl.h> /*I  "petscts.h" I*/
2 
3 #undef __FUNCT__
4 #define __FUNCT__ "TSEventInitialize"
5 /*
6   TSEventInitialize - Initializes TSEvent for TSSolve
7 */
8 PetscErrorCode TSEventInitialize(TSEvent event,TS ts,PetscReal t,Vec U)
9 {
10   PetscErrorCode ierr;
11 
12   PetscFunctionBegin;
13   if (!event) PetscFunctionReturn(0);
14   PetscValidPointer(event,1);
15   PetscValidHeaderSpecific(ts,TS_CLASSID,2);
16   PetscValidHeaderSpecific(U,VEC_CLASSID,4);
17   event->ptime_prev = t;
18   event->iterctr = 0;
19   ierr = (*event->eventhandler)(ts,t,U,event->fvalue_prev,event->ctx);CHKERRQ(ierr);
20   /* Initialize the event recorder */
21   event->recorder.ctr = 0;
22   PetscFunctionReturn(0);
23 }
24 
25 #undef __FUNCT__
26 #define __FUNCT__ "TSEventDestroy"
27 PetscErrorCode TSEventDestroy(TSEvent *event)
28 {
29   PetscErrorCode ierr;
30   PetscInt       i;
31 
32   PetscFunctionBegin;
33   PetscValidPointer(event,1);
34   if (!*event) PetscFunctionReturn(0);
35   ierr = PetscFree((*event)->fvalue);CHKERRQ(ierr);
36   ierr = PetscFree((*event)->fvalue_prev);CHKERRQ(ierr);
37   ierr = PetscFree((*event)->fvalue_right);CHKERRQ(ierr);
38   ierr = PetscFree((*event)->zerocrossing);CHKERRQ(ierr);
39   ierr = PetscFree((*event)->side);CHKERRQ(ierr);
40   ierr = PetscFree((*event)->direction);CHKERRQ(ierr);
41   ierr = PetscFree((*event)->terminate);CHKERRQ(ierr);
42   ierr = PetscFree((*event)->events_zero);CHKERRQ(ierr);
43   ierr = PetscFree((*event)->vtol);CHKERRQ(ierr);
44 
45   for (i=0; i < (*event)->recsize; i++) {
46     ierr = PetscFree((*event)->recorder.eventidx[i]);CHKERRQ(ierr);
47   }
48   ierr = PetscFree((*event)->recorder.eventidx);CHKERRQ(ierr);
49   ierr = PetscFree((*event)->recorder.nevents);CHKERRQ(ierr);
50   ierr = PetscFree((*event)->recorder.stepnum);CHKERRQ(ierr);
51   ierr = PetscFree((*event)->recorder.time);CHKERRQ(ierr);
52 
53   ierr = PetscViewerDestroy(&(*event)->monitor);CHKERRQ(ierr);
54   ierr = PetscFree(*event);CHKERRQ(ierr);
55   PetscFunctionReturn(0);
56 }
57 
58 #undef __FUNCT__
59 #define __FUNCT__ "TSSetEventTolerances"
60 /*@
61    TSSetEventTolerances - Set tolerances for event zero crossings when using event handler
62 
63    Logically Collective
64 
65    Input Arguments:
66 +  ts - time integration context
67 .  tol - scalar tolerance, PETSC_DECIDE to leave current value
68 -  vtol - array of tolerances or NULL, used in preference to tol if present
69 
70    Options Database Keys:
71 .  -ts_event_tol <tol> tolerance for event zero crossing
72 
73    Notes:
74    Must call TSSetEventHandler() before setting the tolerances.
75 
76    The size of vtol is equal to the number of events.
77 
78    Level: beginner
79 
80 .seealso: TS, TSEvent, TSSetEventHandler()
81 @*/
82 PetscErrorCode TSSetEventTolerances(TS ts,PetscReal tol,PetscReal vtol[])
83 {
84   TSEvent        event;
85   PetscInt       i;
86 
87   PetscFunctionBegin;
88   PetscValidHeaderSpecific(ts,TS_CLASSID,1);
89   if (vtol) PetscValidRealPointer(vtol,3);
90   if (!ts->event) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_USER,"Must set the events first by calling TSSetEventHandler()");
91 
92   event = ts->event;
93   if (vtol) {
94     for (i=0; i < event->nevents; i++) event->vtol[i] = vtol[i];
95   } else {
96     if (tol != PETSC_DECIDE || tol != PETSC_DEFAULT) {
97       for (i=0; i < event->nevents; i++) event->vtol[i] = tol;
98     }
99   }
100   PetscFunctionReturn(0);
101 }
102 
103 #undef __FUNCT__
104 #define __FUNCT__ "TSSetEventHandler"
105 /*@C
106    TSSetEventHandler - Sets a monitoring function used for detecting events
107 
108    Logically Collective on TS
109 
110    Input Parameters:
111 +  ts - the TS context obtained from TSCreate()
112 .  nevents - number of local events
113 .  direction - direction of zero crossing to be detected. -1 => Zero crossing in negative direction,
114                +1 => Zero crossing in positive direction, 0 => both ways (one for each event)
115 .  terminate - flag to indicate whether time stepping should be terminated after
116                event is detected (one for each event)
117 .  eventhandler - event monitoring routine
118 .  postevent - [optional] post-event function
119 -  ctx       - [optional] user-defined context for private data for the
120                event monitor and post event routine (use NULL if no
121                context is desired)
122 
123    Calling sequence of eventhandler:
124    PetscErrorCode EventHandler(TS ts,PetscReal t,Vec U,PetscScalar fvalue[],void* ctx)
125 
126    Input Parameters:
127 +  ts  - the TS context
128 .  t   - current time
129 .  U   - current iterate
130 -  ctx - [optional] context passed with eventhandler
131 
132    Output parameters:
133 .  fvalue    - function value of events at time t
134 
135    Calling sequence of postevent:
136    PetscErrorCode PostEvent(TS ts,PetscInt nevents_zero,PetscInt events_zero[],PetscReal t,Vec U,PetscBool forwardsolve,void* ctx)
137 
138    Input Parameters:
139 +  ts - the TS context
140 .  nevents_zero - number of local events whose event function is zero
141 .  events_zero  - indices of local events which have reached zero
142 .  t            - current time
143 .  U            - current solution
144 .  forwardsolve - Flag to indicate whether TS is doing a forward solve (1) or adjoint solve (0)
145 -  ctx          - the context passed with eventhandler
146 
147    Level: intermediate
148 
149 .keywords: TS, event, set
150 
151 .seealso: TSCreate(), TSSetTimeStep(), TSSetConvergedReason()
152 @*/
153 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)
154 {
155   PetscErrorCode ierr;
156   TSEvent        event;
157   PetscInt       i;
158   PetscBool      flg;
159 #if defined PETSC_USE_REAL_SINGLE
160   PetscReal      tol=1e-4;
161 #else
162   PetscReal      tol=1e-6;
163 #endif
164 
165   PetscFunctionBegin;
166   PetscValidHeaderSpecific(ts,TS_CLASSID,1);
167   PetscValidIntPointer(direction,2);
168   PetscValidIntPointer(terminate,3);
169 
170   ierr = PetscNewLog(ts,&event);CHKERRQ(ierr);
171   ierr = PetscMalloc1(nevents,&event->fvalue);CHKERRQ(ierr);
172   ierr = PetscMalloc1(nevents,&event->fvalue_prev);CHKERRQ(ierr);
173   ierr = PetscMalloc1(nevents,&event->fvalue_right);CHKERRQ(ierr);
174   ierr = PetscMalloc1(nevents,&event->zerocrossing);CHKERRQ(ierr);
175   ierr = PetscMalloc1(nevents,&event->side);CHKERRQ(ierr);
176   ierr = PetscMalloc1(nevents,&event->direction);CHKERRQ(ierr);
177   ierr = PetscMalloc1(nevents,&event->terminate);CHKERRQ(ierr);
178   ierr = PetscMalloc1(nevents,&event->vtol);CHKERRQ(ierr);
179   for (i=0; i < nevents; i++) {
180     event->direction[i] = direction[i];
181     event->terminate[i] = terminate[i];
182     event->zerocrossing[i] = PETSC_FALSE;
183     event->side[i] = 0;
184   }
185   ierr = PetscMalloc1(nevents,&event->events_zero);CHKERRQ(ierr);
186   event->nevents = nevents;
187   event->eventhandler = eventhandler;
188   event->postevent = postevent;
189   event->ctx = ctx;
190 
191   event->recsize = 8;  /* Initial size of the recorder */
192   ierr = PetscOptionsBegin(((PetscObject)ts)->comm,"","TS Event options","");CHKERRQ(ierr);
193   {
194     ierr = PetscOptionsReal("-ts_event_tol","Scalar event tolerance for zero crossing check","TSSetEventTolerances",tol,&tol,NULL);CHKERRQ(ierr);
195     ierr = PetscOptionsName("-ts_event_monitor","Print choices made by event handler","",&flg);CHKERRQ(ierr);
196     ierr = PetscOptionsInt("-ts_event_recorder_initial_size","Initial size of event recorder","",event->recsize,&event->recsize,NULL);CHKERRQ(ierr);
197   }
198   PetscOptionsEnd();
199 
200   ierr = PetscMalloc1(event->recsize,&event->recorder.time);CHKERRQ(ierr);
201   ierr = PetscMalloc1(event->recsize,&event->recorder.stepnum);CHKERRQ(ierr);
202   ierr = PetscMalloc1(event->recsize,&event->recorder.nevents);CHKERRQ(ierr);
203   ierr = PetscMalloc1(event->recsize,&event->recorder.eventidx);CHKERRQ(ierr);
204   for (i=0; i < event->recsize; i++) {
205     ierr = PetscMalloc1(event->nevents,&event->recorder.eventidx[i]);CHKERRQ(ierr);
206   }
207 
208   for (i=0; i < event->nevents; i++) event->vtol[i] = tol;
209   if (flg) {ierr = PetscViewerASCIIOpen(PETSC_COMM_SELF,"stdout",&event->monitor);CHKERRQ(ierr);}
210 
211   ierr = TSEventDestroy(&ts->event);CHKERRQ(ierr);
212   ts->event = event;
213   PetscFunctionReturn(0);
214 }
215 
216 #undef __FUNCT__
217 #define __FUNCT__ "TSEventRecorderResize"
218 /*
219   TSEventRecorderResize - Resizes (2X) the event recorder arrays whenever the recording limit (event->recsize)
220                           is reached.
221 */
222 static PetscErrorCode TSEventRecorderResize(TSEvent event)
223 {
224   PetscErrorCode ierr;
225   PetscReal      *time;
226   PetscInt       *stepnum;
227   PetscInt       *nevents;
228   PetscInt       **eventidx;
229   PetscInt       i,fact=2;
230 
231   PetscFunctionBegin;
232 
233   /* Create large arrays */
234   ierr = PetscMalloc1(fact*event->recsize,&time);CHKERRQ(ierr);
235   ierr = PetscMalloc1(fact*event->recsize,&stepnum);CHKERRQ(ierr);
236   ierr = PetscMalloc1(fact*event->recsize,&nevents);CHKERRQ(ierr);
237   ierr = PetscMalloc1(fact*event->recsize,&eventidx);CHKERRQ(ierr);
238   for (i=0; i < fact*event->recsize; i++) {
239     ierr = PetscMalloc1(event->nevents,&eventidx[i]);CHKERRQ(ierr);
240   }
241 
242   /* Copy over data */
243   ierr = PetscMemcpy(time,event->recorder.time,event->recsize*sizeof(PetscReal));CHKERRQ(ierr);
244   ierr = PetscMemcpy(stepnum,event->recorder.stepnum,event->recsize*sizeof(PetscInt));CHKERRQ(ierr);
245   ierr = PetscMemcpy(nevents,event->recorder.nevents,event->recsize*sizeof(PetscInt));CHKERRQ(ierr);
246   for (i=0; i < event->recsize; i++) {
247     ierr = PetscMemcpy(eventidx[i],event->recorder.eventidx[i],event->recorder.nevents[i]*sizeof(PetscInt));CHKERRQ(ierr);
248   }
249 
250   /* Destroy old arrays */
251   for (i=0; i < event->recsize; i++) {
252     ierr = PetscFree(event->recorder.eventidx[i]);CHKERRQ(ierr);
253   }
254   ierr = PetscFree(event->recorder.eventidx);CHKERRQ(ierr);
255   ierr = PetscFree(event->recorder.nevents);CHKERRQ(ierr);
256   ierr = PetscFree(event->recorder.stepnum);CHKERRQ(ierr);
257   ierr = PetscFree(event->recorder.time);CHKERRQ(ierr);
258 
259   /* Set pointers */
260   event->recorder.time = time;
261   event->recorder.stepnum = stepnum;
262   event->recorder.nevents = nevents;
263   event->recorder.eventidx = eventidx;
264 
265   /* Double size */
266   event->recsize *= fact;
267 
268   PetscFunctionReturn(0);
269 }
270 
271 /*
272    Helper rutine to handle user postenvents and recording
273 */
274 #undef __FUNCT__
275 #define __FUNCT__ "TSPostEvent"
276 static PetscErrorCode TSPostEvent(TS ts,PetscReal t,Vec U)
277 {
278   PetscErrorCode ierr;
279   TSEvent        event = ts->event;
280   PetscBool      terminate = PETSC_FALSE;
281   PetscBool      restart = PETSC_FALSE;
282   PetscInt       i,ctr,stepnum;
283   PetscBool      ts_terminate,ts_restart;
284   PetscBool      forwardsolve = PETSC_TRUE; /* Flag indicating that TS is doing a forward solve */
285 
286   PetscFunctionBegin;
287   if (event->postevent) {
288     PetscObjectState state_prev,state_post;
289     ierr = PetscObjectStateGet((PetscObject)U,&state_prev);CHKERRQ(ierr);
290     ierr = (*event->postevent)(ts,event->nevents_zero,event->events_zero,t,U,forwardsolve,event->ctx);CHKERRQ(ierr);
291     ierr = PetscObjectStateGet((PetscObject)U,&state_post);CHKERRQ(ierr);
292     if (state_prev != state_post) restart = PETSC_TRUE;
293   }
294 
295   /* Handle termination events and step restart */
296   for (i=0; i<event->nevents_zero; i++) if (event->terminate[event->events_zero[i]]) terminate = PETSC_TRUE;
297   ierr = MPIU_Allreduce(&terminate,&ts_terminate,1,MPIU_BOOL,MPI_LOR,((PetscObject)ts)->comm);CHKERRQ(ierr);
298   if (ts_terminate) {ierr = TSSetConvergedReason(ts,TS_CONVERGED_EVENT);CHKERRQ(ierr);}
299   event->status = ts_terminate ? TSEVENT_NONE : TSEVENT_RESET_NEXTSTEP;
300   ierr = MPIU_Allreduce(&restart,&ts_restart,1,MPIU_BOOL,MPI_LOR,((PetscObject)ts)->comm);CHKERRQ(ierr);
301   if (ts_restart) ts->steprestart = PETSC_TRUE;
302 
303   event->ptime_prev = t;
304   /* Reset event residual functions as states might get changed by the postevent callback */
305   if (event->postevent) {ierr = (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);CHKERRQ(ierr);}
306   /* Cache current event residual functions */
307   for (i=0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i];
308 
309   /* Record the event in the event recorder */
310   ierr = TSGetTimeStepNumber(ts,&stepnum);CHKERRQ(ierr);
311   ctr = event->recorder.ctr;
312   if (ctr == event->recsize) {
313     ierr = TSEventRecorderResize(event);CHKERRQ(ierr);
314   }
315   event->recorder.time[ctr] = t;
316   event->recorder.stepnum[ctr] = stepnum;
317   event->recorder.nevents[ctr] = event->nevents_zero;
318   for (i=0; i<event->nevents_zero; i++) event->recorder.eventidx[ctr][i] = event->events_zero[i];
319   event->recorder.ctr++;
320   PetscFunctionReturn(0);
321 }
322 
323 /* Uses Anderson-Bjorck variant of regula falsi method */
324 PETSC_STATIC_INLINE PetscReal TSEventComputeStepSize(PetscReal tleft,PetscReal t,PetscReal tright,PetscScalar fleft,PetscScalar f,PetscScalar fright,PetscInt side,PetscReal dt)
325 {
326   PetscReal new_dt, scal = 1.0;
327   if (PetscRealPart(fleft)*PetscRealPart(f) < 0) {
328     if (side == 1) {
329       scal = (PetscRealPart(fright) - PetscRealPart(f))/PetscRealPart(fright);
330       if (scal < PETSC_SMALL) scal = 0.5;
331     }
332     new_dt = (scal*PetscRealPart(fleft)*t - PetscRealPart(f)*tleft)/(scal*PetscRealPart(fleft) - PetscRealPart(f)) - tleft;
333   } else {
334     if (side == -1) {
335       scal = (PetscRealPart(fleft) - PetscRealPart(f))/PetscRealPart(fleft);
336       if (scal < PETSC_SMALL) scal = 0.5;
337     }
338     new_dt = (PetscRealPart(f)*tright - scal*PetscRealPart(fright)*t)/(PetscRealPart(f) - scal*PetscRealPart(fright)) - t;
339   }
340   return PetscMin(dt,new_dt);
341 }
342 
343 
344 #undef __FUNCT__
345 #define __FUNCT__ "TSEventHandler"
346 PetscErrorCode TSEventHandler(TS ts)
347 {
348   PetscErrorCode ierr;
349   TSEvent        event;
350   PetscReal      t;
351   Vec            U;
352   PetscInt       i;
353   PetscReal      dt,dt_min;
354   PetscInt       rollback=0,in[2],out[2];
355   PetscInt       fvalue_sign,fvalueprev_sign;
356 
357   PetscFunctionBegin;
358   PetscValidHeaderSpecific(ts,TS_CLASSID,1);
359   if (!ts->event) PetscFunctionReturn(0);
360   event = ts->event;
361 
362   ierr = TSGetTime(ts,&t);CHKERRQ(ierr);
363   ierr = TSGetTimeStep(ts,&dt);CHKERRQ(ierr);
364   ierr = TSGetSolution(ts,&U);CHKERRQ(ierr);
365 
366   if (event->status == TSEVENT_NONE) {
367     if (ts->steps == 1) /* After first successful step */
368       event->timestep_orig = ts->ptime - ts->ptime_prev;
369     event->timestep_prev = dt;
370   }
371 
372   if (event->status == TSEVENT_RESET_NEXTSTEP) {
373     /* Restore time step */
374     dt = PetscMin(event->timestep_orig,event->timestep_prev);
375     ierr = TSSetTimeStep(ts,dt);CHKERRQ(ierr);
376     event->status = TSEVENT_NONE;
377   }
378 
379   if (event->status == TSEVENT_NONE) {
380     event->ptime_end = t;
381   }
382 
383   ierr = (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);CHKERRQ(ierr);
384 
385   for (i=0; i < event->nevents; i++) {
386     if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) {
387       event->status = TSEVENT_ZERO;
388       event->fvalue_right[i] = event->fvalue[i];
389       continue;
390     }
391     fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i]));
392     fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i]));
393     if (fvalueprev_sign != 0 && (fvalue_sign != fvalueprev_sign) && (PetscAbsScalar(event->fvalue_prev[i]) > event->vtol[i])) {
394       switch (event->direction[i]) {
395       case -1:
396         if (fvalue_sign < 0) {
397           rollback = 1;
398 
399           /* Compute new time step */
400           dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],dt);
401 
402           if (event->monitor) {
403             ierr = PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Event %D interval detected [%g - %g]\n",event->iterctr,i,(double)event->ptime_prev,(double)t);CHKERRQ(ierr);
404           }
405           event->fvalue_right[i] = event->fvalue[i];
406           event->side[i] = 1;
407 
408           if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
409           event->status = TSEVENT_LOCATED_INTERVAL;
410         }
411         break;
412       case 1:
413         if (fvalue_sign > 0) {
414           rollback = 1;
415 
416           /* Compute new time step */
417           dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],dt);
418 
419           if (event->monitor) {
420             ierr = PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Event %D interval detected [%g - %g]\n",event->iterctr,i,(double)event->ptime_prev,(double)t);CHKERRQ(ierr);
421           }
422           event->fvalue_right[i] = event->fvalue[i];
423           event->side[i] = 1;
424 
425           if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
426           event->status = TSEVENT_LOCATED_INTERVAL;
427         }
428         break;
429       case 0:
430         rollback = 1;
431 
432         /* Compute new time step */
433         dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],dt);
434 
435         if (event->monitor) {
436           ierr = PetscViewerASCIIPrintf(event->monitor,"TSEvent: iter %D - Event %D interval detected [%g - %g]\n",event->iterctr,i,(double)event->ptime_prev,(double)t);CHKERRQ(ierr);
437         }
438         event->fvalue_right[i] = event->fvalue[i];
439         event->side[i] = 1;
440 
441         if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE;
442         event->status = TSEVENT_LOCATED_INTERVAL;
443 
444         break;
445       }
446     }
447   }
448 
449   in[0] = event->status; in[1] = rollback;
450   ierr = MPIU_Allreduce(in,out,2,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)ts));CHKERRQ(ierr);
451   event->status = (TSEventStatus)out[0]; rollback = out[1];
452   if (rollback) event->status = TSEVENT_LOCATED_INTERVAL;
453 
454   event->nevents_zero = 0;
455   if (event->status == TSEVENT_ZERO) {
456     for (i=0; i < event->nevents; i++) {
457       if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) {
458         event->events_zero[event->nevents_zero++] = i;
459         if (event->monitor) {
460           ierr = PetscViewerASCIIPrintf(event->monitor,"TSEvent: Event %D zero crossing at time %g located in %D iterations\n",i,(double)t,event->iterctr);CHKERRQ(ierr);
461         }
462         event->zerocrossing[i] = PETSC_FALSE;
463       }
464       event->side[i] = 0;
465     }
466     ierr = TSPostEvent(ts,t,U);CHKERRQ(ierr);
467 
468     dt = event->ptime_end - t;
469     if (PetscAbsReal(dt) < PETSC_SMALL) dt += PetscMin(event->timestep_orig,event->timestep_prev); /* XXX Should be done better */
470     ierr = TSSetTimeStep(ts,dt);CHKERRQ(ierr);
471     event->iterctr = 0;
472     PetscFunctionReturn(0);
473   }
474 
475   if (event->status == TSEVENT_LOCATED_INTERVAL) {
476     ierr = TSRollBack(ts);CHKERRQ(ierr);
477     ierr = TSSetConvergedReason(ts,TS_CONVERGED_ITERATING);CHKERRQ(ierr);
478     event->status = TSEVENT_PROCESSING;
479     event->ptime_right = t;
480   } else {
481     for (i=0; i < event->nevents; i++) {
482       if (event->status == TSEVENT_PROCESSING && event->zerocrossing[i]) {
483         /* Compute new time step */
484         dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],dt);
485         event->side[i] = -1;
486       }
487       event->fvalue_prev[i] = event->fvalue[i];
488     }
489     if (event->monitor && event->status == TSEVENT_PROCESSING) {
490       ierr = 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);CHKERRQ(ierr);
491     }
492     event->ptime_prev = t;
493   }
494 
495   if (event->status == TSEVENT_PROCESSING) event->iterctr++;
496 
497   ierr = MPIU_Allreduce(&dt,&dt_min,1,MPIU_REAL,MPIU_MIN,PetscObjectComm((PetscObject)ts));CHKERRQ(ierr);
498   ierr = TSSetTimeStep(ts,dt_min);CHKERRQ(ierr);
499   PetscFunctionReturn(0);
500 }
501 
502 #undef __FUNCT__
503 #define __FUNCT__ "TSAdjointEventHandler"
504 PetscErrorCode TSAdjointEventHandler(TS ts)
505 {
506   PetscErrorCode ierr;
507   TSEvent        event;
508   PetscReal      t;
509   Vec            U;
510   PetscInt       ctr;
511   PetscBool      forwardsolve=PETSC_FALSE; /* Flag indicating that TS is doing an adjoint solve */
512 
513   PetscFunctionBegin;
514   PetscValidHeaderSpecific(ts,TS_CLASSID,1);
515   if (!ts->event) PetscFunctionReturn(0);
516   event = ts->event;
517 
518   ierr = TSGetTime(ts,&t);CHKERRQ(ierr);
519   ierr = TSGetSolution(ts,&U);CHKERRQ(ierr);
520 
521   ctr = event->recorder.ctr-1;
522   if (ctr >= 0 && PetscAbsReal(t - event->recorder.time[ctr]) < PETSC_SMALL) {
523     /* Call the user postevent function */
524     if (event->postevent) {
525       ierr = (*event->postevent)(ts,event->recorder.nevents[ctr],event->recorder.eventidx[ctr],t,U,forwardsolve,event->ctx);CHKERRQ(ierr);
526       event->recorder.ctr--;
527     }
528   }
529 
530   PetscBarrier((PetscObject)ts);
531   PetscFunctionReturn(0);
532 }
533