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