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 PetscFunctionReturn(0); 21 } 22 23 #undef __FUNCT__ 24 #define __FUNCT__ "TSEventDestroy" 25 PetscErrorCode TSEventDestroy(TSEvent *event) 26 { 27 PetscErrorCode ierr; 28 PetscInt i; 29 30 PetscFunctionBegin; 31 PetscValidPointer(event,1); 32 if (!*event) PetscFunctionReturn(0); 33 if (--(*event)->refct > 0) {*event = 0; PetscFunctionReturn(0);} 34 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 PetscEventHandler(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 if(nevents) { 168 PetscValidIntPointer(direction,2); 169 PetscValidIntPointer(terminate,3); 170 } 171 172 ierr = PetscNewLog(ts,&event);CHKERRQ(ierr); 173 ierr = PetscMalloc1(nevents,&event->fvalue);CHKERRQ(ierr); 174 ierr = PetscMalloc1(nevents,&event->fvalue_prev);CHKERRQ(ierr); 175 ierr = PetscMalloc1(nevents,&event->fvalue_right);CHKERRQ(ierr); 176 ierr = PetscMalloc1(nevents,&event->zerocrossing);CHKERRQ(ierr); 177 ierr = PetscMalloc1(nevents,&event->side);CHKERRQ(ierr); 178 ierr = PetscMalloc1(nevents,&event->direction);CHKERRQ(ierr); 179 ierr = PetscMalloc1(nevents,&event->terminate);CHKERRQ(ierr); 180 ierr = PetscMalloc1(nevents,&event->vtol);CHKERRQ(ierr); 181 for (i=0; i < nevents; i++) { 182 event->direction[i] = direction[i]; 183 event->terminate[i] = terminate[i]; 184 event->zerocrossing[i] = PETSC_FALSE; 185 event->side[i] = 0; 186 } 187 ierr = PetscMalloc1(nevents,&event->events_zero);CHKERRQ(ierr); 188 event->nevents = nevents; 189 event->eventhandler = eventhandler; 190 event->postevent = postevent; 191 event->ctx = ctx; 192 193 event->recsize = 8; /* Initial size of the recorder */ 194 ierr = PetscOptionsBegin(((PetscObject)ts)->comm,"","TS Event options","");CHKERRQ(ierr); 195 { 196 ierr = PetscOptionsReal("-ts_event_tol","Scalar event tolerance for zero crossing check","TSSetEventTolerances",tol,&tol,NULL);CHKERRQ(ierr); 197 ierr = PetscOptionsName("-ts_event_monitor","Print choices made by event handler","",&flg);CHKERRQ(ierr); 198 ierr = PetscOptionsInt("-ts_event_recorder_initial_size","Initial size of event recorder","",event->recsize,&event->recsize,NULL);CHKERRQ(ierr); 199 } 200 PetscOptionsEnd(); 201 202 ierr = PetscMalloc1(event->recsize,&event->recorder.time);CHKERRQ(ierr); 203 ierr = PetscMalloc1(event->recsize,&event->recorder.stepnum);CHKERRQ(ierr); 204 ierr = PetscMalloc1(event->recsize,&event->recorder.nevents);CHKERRQ(ierr); 205 ierr = PetscMalloc1(event->recsize,&event->recorder.eventidx);CHKERRQ(ierr); 206 for (i=0; i < event->recsize; i++) { 207 ierr = PetscMalloc1(event->nevents,&event->recorder.eventidx[i]);CHKERRQ(ierr); 208 } 209 /* Initialize the event recorder */ 210 event->recorder.ctr = 0; 211 212 for (i=0; i < event->nevents; i++) event->vtol[i] = tol; 213 if (flg) {ierr = PetscViewerASCIIOpen(PETSC_COMM_SELF,"stdout",&event->monitor);CHKERRQ(ierr);} 214 215 ierr = TSEventDestroy(&ts->event);CHKERRQ(ierr); 216 ts->event = event; 217 ts->event->refct = 1; 218 PetscFunctionReturn(0); 219 } 220 221 #undef __FUNCT__ 222 #define __FUNCT__ "TSEventRecorderResize" 223 /* 224 TSEventRecorderResize - Resizes (2X) the event recorder arrays whenever the recording limit (event->recsize) 225 is reached. 226 */ 227 static PetscErrorCode TSEventRecorderResize(TSEvent event) 228 { 229 PetscErrorCode ierr; 230 PetscReal *time; 231 PetscInt *stepnum; 232 PetscInt *nevents; 233 PetscInt **eventidx; 234 PetscInt i,fact=2; 235 236 PetscFunctionBegin; 237 238 /* Create large arrays */ 239 ierr = PetscMalloc1(fact*event->recsize,&time);CHKERRQ(ierr); 240 ierr = PetscMalloc1(fact*event->recsize,&stepnum);CHKERRQ(ierr); 241 ierr = PetscMalloc1(fact*event->recsize,&nevents);CHKERRQ(ierr); 242 ierr = PetscMalloc1(fact*event->recsize,&eventidx);CHKERRQ(ierr); 243 for (i=0; i < fact*event->recsize; i++) { 244 ierr = PetscMalloc1(event->nevents,&eventidx[i]);CHKERRQ(ierr); 245 } 246 247 /* Copy over data */ 248 ierr = PetscMemcpy(time,event->recorder.time,event->recsize*sizeof(PetscReal));CHKERRQ(ierr); 249 ierr = PetscMemcpy(stepnum,event->recorder.stepnum,event->recsize*sizeof(PetscInt));CHKERRQ(ierr); 250 ierr = PetscMemcpy(nevents,event->recorder.nevents,event->recsize*sizeof(PetscInt));CHKERRQ(ierr); 251 for (i=0; i < event->recsize; i++) { 252 ierr = PetscMemcpy(eventidx[i],event->recorder.eventidx[i],event->recorder.nevents[i]*sizeof(PetscInt));CHKERRQ(ierr); 253 } 254 255 /* Destroy old arrays */ 256 for (i=0; i < event->recsize; i++) { 257 ierr = PetscFree(event->recorder.eventidx[i]);CHKERRQ(ierr); 258 } 259 ierr = PetscFree(event->recorder.eventidx);CHKERRQ(ierr); 260 ierr = PetscFree(event->recorder.nevents);CHKERRQ(ierr); 261 ierr = PetscFree(event->recorder.stepnum);CHKERRQ(ierr); 262 ierr = PetscFree(event->recorder.time);CHKERRQ(ierr); 263 264 /* Set pointers */ 265 event->recorder.time = time; 266 event->recorder.stepnum = stepnum; 267 event->recorder.nevents = nevents; 268 event->recorder.eventidx = eventidx; 269 270 /* Double size */ 271 event->recsize *= fact; 272 273 PetscFunctionReturn(0); 274 } 275 276 /* 277 Helper rutine to handle user postenvents and recording 278 */ 279 #undef __FUNCT__ 280 #define __FUNCT__ "TSPostEvent" 281 static PetscErrorCode TSPostEvent(TS ts,PetscReal t,Vec U) 282 { 283 PetscErrorCode ierr; 284 TSEvent event = ts->event; 285 PetscBool terminate = PETSC_FALSE; 286 PetscBool restart = PETSC_FALSE; 287 PetscInt i,ctr,stepnum; 288 PetscBool ts_terminate,ts_restart; 289 PetscBool forwardsolve = PETSC_TRUE; /* Flag indicating that TS is doing a forward solve */ 290 291 PetscFunctionBegin; 292 if (event->postevent) { 293 PetscObjectState state_prev,state_post; 294 ierr = PetscObjectStateGet((PetscObject)U,&state_prev);CHKERRQ(ierr); 295 ierr = (*event->postevent)(ts,event->nevents_zero,event->events_zero,t,U,forwardsolve,event->ctx);CHKERRQ(ierr); 296 ierr = PetscObjectStateGet((PetscObject)U,&state_post);CHKERRQ(ierr); 297 if (state_prev != state_post) restart = PETSC_TRUE; 298 } 299 300 /* Handle termination events and step restart */ 301 for (i=0; i<event->nevents_zero; i++) if (event->terminate[event->events_zero[i]]) terminate = PETSC_TRUE; 302 ierr = MPIU_Allreduce(&terminate,&ts_terminate,1,MPIU_BOOL,MPI_LOR,((PetscObject)ts)->comm);CHKERRQ(ierr); 303 if (ts_terminate) {ierr = TSSetConvergedReason(ts,TS_CONVERGED_EVENT);CHKERRQ(ierr);} 304 event->status = ts_terminate ? TSEVENT_NONE : TSEVENT_RESET_NEXTSTEP; 305 ierr = MPIU_Allreduce(&restart,&ts_restart,1,MPIU_BOOL,MPI_LOR,((PetscObject)ts)->comm);CHKERRQ(ierr); 306 if (ts_restart) ts->steprestart = PETSC_TRUE; 307 308 event->ptime_prev = t; 309 /* Reset event residual functions as states might get changed by the postevent callback */ 310 if (event->postevent) {ierr = (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);CHKERRQ(ierr);} 311 /* Cache current event residual functions */ 312 for (i=0; i < event->nevents; i++) event->fvalue_prev[i] = event->fvalue[i]; 313 314 /* Record the event in the event recorder */ 315 ierr = TSGetTimeStepNumber(ts,&stepnum);CHKERRQ(ierr); 316 ctr = event->recorder.ctr; 317 if (ctr == event->recsize) { 318 ierr = TSEventRecorderResize(event);CHKERRQ(ierr); 319 } 320 event->recorder.time[ctr] = t; 321 event->recorder.stepnum[ctr] = stepnum; 322 event->recorder.nevents[ctr] = event->nevents_zero; 323 for (i=0; i<event->nevents_zero; i++) event->recorder.eventidx[ctr][i] = event->events_zero[i]; 324 event->recorder.ctr++; 325 PetscFunctionReturn(0); 326 } 327 328 /* Uses Anderson-Bjorck variant of regula falsi method */ 329 PETSC_STATIC_INLINE PetscReal TSEventComputeStepSize(PetscReal tleft,PetscReal t,PetscReal tright,PetscScalar fleft,PetscScalar f,PetscScalar fright,PetscInt side,PetscReal dt) 330 { 331 PetscReal new_dt, scal = 1.0; 332 if (PetscRealPart(fleft)*PetscRealPart(f) < 0) { 333 if (side == 1) { 334 scal = (PetscRealPart(fright) - PetscRealPart(f))/PetscRealPart(fright); 335 if (scal < PETSC_SMALL) scal = 0.5; 336 } 337 new_dt = (scal*PetscRealPart(fleft)*t - PetscRealPart(f)*tleft)/(scal*PetscRealPart(fleft) - PetscRealPart(f)) - tleft; 338 } else { 339 if (side == -1) { 340 scal = (PetscRealPart(fleft) - PetscRealPart(f))/PetscRealPart(fleft); 341 if (scal < PETSC_SMALL) scal = 0.5; 342 } 343 new_dt = (PetscRealPart(f)*tright - scal*PetscRealPart(fright)*t)/(PetscRealPart(f) - scal*PetscRealPart(fright)) - t; 344 } 345 return PetscMin(dt,new_dt); 346 } 347 348 349 #undef __FUNCT__ 350 #define __FUNCT__ "TSEventHandler" 351 PetscErrorCode TSEventHandler(TS ts) 352 { 353 PetscErrorCode ierr; 354 TSEvent event; 355 PetscReal t; 356 Vec U; 357 PetscInt i; 358 PetscReal dt,dt_min; 359 PetscInt rollback=0,in[2],out[2]; 360 PetscInt fvalue_sign,fvalueprev_sign; 361 362 PetscFunctionBegin; 363 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 364 if (!ts->event) PetscFunctionReturn(0); 365 event = ts->event; 366 367 ierr = TSGetTime(ts,&t);CHKERRQ(ierr); 368 ierr = TSGetTimeStep(ts,&dt);CHKERRQ(ierr); 369 ierr = TSGetSolution(ts,&U);CHKERRQ(ierr); 370 371 if (event->status == TSEVENT_NONE) { 372 if (ts->steps == 1) /* After first successful step */ 373 event->timestep_orig = ts->ptime - ts->ptime_prev; 374 event->timestep_prev = dt; 375 } 376 377 if (event->status == TSEVENT_RESET_NEXTSTEP) { 378 /* Restore time step */ 379 dt = PetscMin(event->timestep_orig,event->timestep_prev); 380 ierr = TSSetTimeStep(ts,dt);CHKERRQ(ierr); 381 event->status = TSEVENT_NONE; 382 } 383 384 if (event->status == TSEVENT_NONE) { 385 event->ptime_end = t; 386 } 387 388 ierr = (*event->eventhandler)(ts,t,U,event->fvalue,event->ctx);CHKERRQ(ierr); 389 390 for (i=0; i < event->nevents; i++) { 391 if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) { 392 event->status = TSEVENT_ZERO; 393 event->fvalue_right[i] = event->fvalue[i]; 394 continue; 395 } 396 fvalue_sign = PetscSign(PetscRealPart(event->fvalue[i])); 397 fvalueprev_sign = PetscSign(PetscRealPart(event->fvalue_prev[i])); 398 if (fvalueprev_sign != 0 && (fvalue_sign != fvalueprev_sign) && (PetscAbsScalar(event->fvalue_prev[i]) > event->vtol[i])) { 399 switch (event->direction[i]) { 400 case -1: 401 if (fvalue_sign < 0) { 402 rollback = 1; 403 404 /* Compute new time step */ 405 dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],dt); 406 407 if (event->monitor) { 408 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); 409 } 410 event->fvalue_right[i] = event->fvalue[i]; 411 event->side[i] = 1; 412 413 if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE; 414 event->status = TSEVENT_LOCATED_INTERVAL; 415 } 416 break; 417 case 1: 418 if (fvalue_sign > 0) { 419 rollback = 1; 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 424 if (event->monitor) { 425 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); 426 } 427 event->fvalue_right[i] = event->fvalue[i]; 428 event->side[i] = 1; 429 430 if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE; 431 event->status = TSEVENT_LOCATED_INTERVAL; 432 } 433 break; 434 case 0: 435 rollback = 1; 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 440 if (event->monitor) { 441 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); 442 } 443 event->fvalue_right[i] = event->fvalue[i]; 444 event->side[i] = 1; 445 446 if (!event->iterctr) event->zerocrossing[i] = PETSC_TRUE; 447 event->status = TSEVENT_LOCATED_INTERVAL; 448 449 break; 450 } 451 } 452 } 453 454 in[0] = event->status; in[1] = rollback; 455 ierr = MPIU_Allreduce(in,out,2,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)ts));CHKERRQ(ierr); 456 event->status = (TSEventStatus)out[0]; rollback = out[1]; 457 if (rollback) event->status = TSEVENT_LOCATED_INTERVAL; 458 459 event->nevents_zero = 0; 460 if (event->status == TSEVENT_ZERO) { 461 for (i=0; i < event->nevents; i++) { 462 if (PetscAbsScalar(event->fvalue[i]) < event->vtol[i]) { 463 event->events_zero[event->nevents_zero++] = i; 464 if (event->monitor) { 465 ierr = PetscViewerASCIIPrintf(event->monitor,"TSEvent: Event %D zero crossing at time %g located in %D iterations\n",i,(double)t,event->iterctr);CHKERRQ(ierr); 466 } 467 event->zerocrossing[i] = PETSC_FALSE; 468 } 469 event->side[i] = 0; 470 } 471 ierr = TSPostEvent(ts,t,U);CHKERRQ(ierr); 472 473 dt = event->ptime_end - t; 474 if (PetscAbsReal(dt) < PETSC_SMALL) dt += PetscMin(event->timestep_orig,event->timestep_prev); /* XXX Should be done better */ 475 ierr = TSSetTimeStep(ts,dt);CHKERRQ(ierr); 476 event->iterctr = 0; 477 PetscFunctionReturn(0); 478 } 479 480 if (event->status == TSEVENT_LOCATED_INTERVAL) { 481 ierr = TSRollBack(ts);CHKERRQ(ierr); 482 ierr = TSSetConvergedReason(ts,TS_CONVERGED_ITERATING);CHKERRQ(ierr); 483 event->status = TSEVENT_PROCESSING; 484 event->ptime_right = t; 485 } else { 486 for (i=0; i < event->nevents; i++) { 487 if (event->status == TSEVENT_PROCESSING && event->zerocrossing[i]) { 488 /* Compute new time step */ 489 dt = TSEventComputeStepSize(event->ptime_prev,t,event->ptime_right,event->fvalue_prev[i],event->fvalue[i],event->fvalue_right[i],event->side[i],dt); 490 event->side[i] = -1; 491 } 492 event->fvalue_prev[i] = event->fvalue[i]; 493 } 494 if (event->monitor && event->status == TSEVENT_PROCESSING) { 495 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); 496 } 497 event->ptime_prev = t; 498 } 499 500 if (event->status == TSEVENT_PROCESSING) event->iterctr++; 501 502 ierr = MPIU_Allreduce(&dt,&dt_min,1,MPIU_REAL,MPIU_MIN,PetscObjectComm((PetscObject)ts));CHKERRQ(ierr); 503 ierr = TSSetTimeStep(ts,dt_min);CHKERRQ(ierr); 504 PetscFunctionReturn(0); 505 } 506 507 #undef __FUNCT__ 508 #define __FUNCT__ "TSAdjointEventHandler" 509 PetscErrorCode TSAdjointEventHandler(TS ts) 510 { 511 PetscErrorCode ierr; 512 TSEvent event; 513 PetscReal t; 514 Vec U; 515 PetscInt ctr; 516 PetscBool forwardsolve=PETSC_FALSE; /* Flag indicating that TS is doing an adjoint solve */ 517 518 PetscFunctionBegin; 519 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 520 if (!ts->event) PetscFunctionReturn(0); 521 event = ts->event; 522 523 ierr = TSGetTime(ts,&t);CHKERRQ(ierr); 524 ierr = TSGetSolution(ts,&U);CHKERRQ(ierr); 525 526 ctr = event->recorder.ctr-1; 527 if (ctr >= 0 && PetscAbsReal(t - event->recorder.time[ctr]) < PETSC_SMALL) { 528 /* Call the user postevent function */ 529 if (event->postevent) { 530 ierr = (*event->postevent)(ts,event->recorder.nevents[ctr],event->recorder.eventidx[ctr],t,U,forwardsolve,event->ctx);CHKERRQ(ierr); 531 event->recorder.ctr--; 532 } 533 } 534 535 PetscBarrier((PetscObject)ts); 536 PetscFunctionReturn(0); 537 } 538