Actual source code: tsmon.c

petsc-main 2021-04-20
Report Typos and Errors
  1: #include <petsc/private/tsimpl.h>
  2: #include <petscdm.h>
  3: #include <petscdraw.h>

  5: /*@C
  6:    TSMonitor - Runs all user-provided monitor routines set using TSMonitorSet()

  8:    Collective on TS

 10:    Input Parameters:
 11: +  ts - time stepping context obtained from TSCreate()
 12: .  step - step number that has just completed
 13: .  ptime - model time of the state
 14: -  u - state at the current model time

 16:    Notes:
 17:    TSMonitor() is typically used automatically within the time stepping implementations.
 18:    Users would almost never call this routine directly.

 20:    A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions

 22:    Level: developer

 24: @*/
 25: PetscErrorCode TSMonitor(TS ts,PetscInt step,PetscReal ptime,Vec u)
 26: {
 27:   DM             dm;
 28:   PetscInt       i,n = ts->numbermonitors;


 35:   TSGetDM(ts,&dm);
 36:   DMSetOutputSequenceNumber(dm,step,ptime);

 38:   VecLockReadPush(u);
 39:   for (i=0; i<n; i++) {
 40:     (*ts->monitor[i])(ts,step,ptime,u,ts->monitorcontext[i]);
 41:   }
 42:   VecLockReadPop(u);
 43:   return(0);
 44: }

 46: /*@C
 47:    TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user

 49:    Collective on TS

 51:    Input Parameters:
 52: +  ts - TS object you wish to monitor
 53: .  name - the monitor type one is seeking
 54: .  help - message indicating what monitoring is done
 55: .  manual - manual page for the monitor
 56: .  monitor - the monitor function
 57: -  monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the TS or PetscViewer objects

 59:    Level: developer

 61: .seealso: PetscOptionsGetViewer(), PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(),
 62:           PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
 63:           PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(),
 64:           PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
 65:           PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
 66:           PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
 67:           PetscOptionsFList(), PetscOptionsEList()
 68: @*/
 69: PetscErrorCode  TSMonitorSetFromOptions(TS ts,const char name[],const char help[], const char manual[],PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,PetscViewerAndFormat*),PetscErrorCode (*monitorsetup)(TS,PetscViewerAndFormat*))
 70: {
 71:   PetscErrorCode    ierr;
 72:   PetscViewer       viewer;
 73:   PetscViewerFormat format;
 74:   PetscBool         flg;

 77:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)ts),((PetscObject) ts)->options,((PetscObject)ts)->prefix,name,&viewer,&format,&flg);
 78:   if (flg) {
 79:     PetscViewerAndFormat *vf;
 80:     PetscViewerAndFormatCreate(viewer,format,&vf);
 81:     PetscObjectDereference((PetscObject)viewer);
 82:     if (monitorsetup) {
 83:       (*monitorsetup)(ts,vf);
 84:     }
 85:     TSMonitorSet(ts,(PetscErrorCode (*)(TS,PetscInt,PetscReal,Vec,void*))monitor,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
 86:   }
 87:   return(0);
 88: }

 90: /*@C
 91:    TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
 92:    timestep to display the iteration's  progress.

 94:    Logically Collective on TS

 96:    Input Parameters:
 97: +  ts - the TS context obtained from TSCreate()
 98: .  monitor - monitoring routine
 99: .  mctx - [optional] user-defined context for private data for the
100:              monitor routine (use NULL if no context is desired)
101: -  monitordestroy - [optional] routine that frees monitor context
102:           (may be NULL)

104:    Calling sequence of monitor:
105: $    PetscErrorCode monitor(TS ts,PetscInt steps,PetscReal time,Vec u,void *mctx)

107: +    ts - the TS context
108: .    steps - iteration number (after the final time step the monitor routine may be called with a step of -1, this indicates the solution has been interpolated to this time)
109: .    time - current time
110: .    u - current iterate
111: -    mctx - [optional] monitoring context

113:    Notes:
114:    This routine adds an additional monitor to the list of monitors that
115:    already has been loaded.

117:    Fortran Notes:
118:     Only a single monitor function can be set for each TS object

120:    Level: intermediate

122: .seealso: TSMonitorDefault(), TSMonitorCancel()
123: @*/
124: PetscErrorCode  TSMonitorSet(TS ts,PetscErrorCode (*monitor)(TS,PetscInt,PetscReal,Vec,void*),void *mctx,PetscErrorCode (*mdestroy)(void**))
125: {
127:   PetscInt       i;
128:   PetscBool      identical;

132:   for (i=0; i<ts->numbermonitors;i++) {
133:     PetscMonitorCompare((PetscErrorCode (*)(void))monitor,mctx,mdestroy,(PetscErrorCode (*)(void))ts->monitor[i],ts->monitorcontext[i],ts->monitordestroy[i],&identical);
134:     if (identical) return(0);
135:   }
136:   if (ts->numbermonitors >= MAXTSMONITORS) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Too many monitors set");
137:   ts->monitor[ts->numbermonitors]          = monitor;
138:   ts->monitordestroy[ts->numbermonitors]   = mdestroy;
139:   ts->monitorcontext[ts->numbermonitors++] = (void*)mctx;
140:   return(0);
141: }

143: /*@C
144:    TSMonitorCancel - Clears all the monitors that have been set on a time-step object.

146:    Logically Collective on TS

148:    Input Parameters:
149: .  ts - the TS context obtained from TSCreate()

151:    Notes:
152:    There is no way to remove a single, specific monitor.

154:    Level: intermediate

156: .seealso: TSMonitorDefault(), TSMonitorSet()
157: @*/
158: PetscErrorCode  TSMonitorCancel(TS ts)
159: {
161:   PetscInt       i;

165:   for (i=0; i<ts->numbermonitors; i++) {
166:     if (ts->monitordestroy[i]) {
167:       (*ts->monitordestroy[i])(&ts->monitorcontext[i]);
168:     }
169:   }
170:   ts->numbermonitors = 0;
171:   return(0);
172: }

174: /*@C
175:    TSMonitorDefault - The Default monitor, prints the timestep and time for each step

177:    Level: intermediate

179: .seealso:  TSMonitorSet()
180: @*/
181: PetscErrorCode TSMonitorDefault(TS ts,PetscInt step,PetscReal ptime,Vec v,PetscViewerAndFormat *vf)
182: {
184:   PetscViewer    viewer =  vf->viewer;
185:   PetscBool      iascii,ibinary;

189:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
190:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&ibinary);
191:   PetscViewerPushFormat(viewer,vf->format);
192:   if (iascii) {
193:     PetscViewerASCIIAddTab(viewer,((PetscObject)ts)->tablevel);
194:     if (step == -1){ /* this indicates it is an interpolated solution */
195:       PetscViewerASCIIPrintf(viewer,"Interpolated solution at time %g between steps %D and %D\n",(double)ptime,ts->steps-1,ts->steps);
196:     } else {
197:       PetscViewerASCIIPrintf(viewer,"%D TS dt %g time %g%s",step,(double)ts->time_step,(double)ptime,ts->steprollback ? " (r)\n" : "\n");
198:     }
199:     PetscViewerASCIISubtractTab(viewer,((PetscObject)ts)->tablevel);
200:   } else if (ibinary) {
201:     PetscMPIInt rank;
202:     MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
203:     if (!rank) {
204:       PetscBool skipHeader;
205:       PetscInt  classid = REAL_FILE_CLASSID;

207:       PetscViewerBinaryGetSkipHeader(viewer,&skipHeader);
208:       if (!skipHeader) {
209:          PetscViewerBinaryWrite(viewer,&classid,1,PETSC_INT);
210:        }
211:       PetscRealView(1,&ptime,viewer);
212:     } else {
213:       PetscRealView(0,&ptime,viewer);
214:     }
215:   }
216:   PetscViewerPopFormat(viewer);
217:   return(0);
218: }

220: /*@C
221:    TSMonitorExtreme - Prints the extreme values of the solution at each timestep

223:    Level: intermediate

225: .seealso:  TSMonitorSet()
226: @*/
227: PetscErrorCode TSMonitorExtreme(TS ts,PetscInt step,PetscReal ptime,Vec v,PetscViewerAndFormat *vf)
228: {
230:   PetscViewer    viewer =  vf->viewer;
231:   PetscBool      iascii;
232:   PetscReal      max,min;


237:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
238:   PetscViewerPushFormat(viewer,vf->format);
239:   if (iascii) {
240:     VecMax(v,NULL,&max);
241:     VecMin(v,NULL,&min);
242:     PetscViewerASCIIAddTab(viewer,((PetscObject)ts)->tablevel);
243:     PetscViewerASCIIPrintf(viewer,"%D TS dt %g time %g%s max %g min %g\n",step,(double)ts->time_step,(double)ptime,ts->steprollback ? " (r)" : "",(double)max,(double)min);
244:     PetscViewerASCIISubtractTab(viewer,((PetscObject)ts)->tablevel);
245:   }
246:   PetscViewerPopFormat(viewer);
247:   return(0);
248: }

250: /*@C
251:    TSMonitorLGCtxCreate - Creates a TSMonitorLGCtx context for use with
252:    TS to monitor the solution process graphically in various ways

254:    Collective on TS

256:    Input Parameters:
257: +  host - the X display to open, or null for the local machine
258: .  label - the title to put in the title bar
259: .  x, y - the screen coordinates of the upper left coordinate of the window
260: .  m, n - the screen width and height in pixels
261: -  howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time

263:    Output Parameter:
264: .  ctx - the context

266:    Options Database Key:
267: +  -ts_monitor_lg_timestep - automatically sets line graph monitor
268: +  -ts_monitor_lg_timestep_log - automatically sets line graph monitor
269: .  -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling TSMonitorLGSetDisplayVariables() or TSMonitorLGCtxSetDisplayVariables())
270: .  -ts_monitor_lg_error -  monitor the error
271: .  -ts_monitor_lg_ksp_iterations - monitor the number of KSP iterations needed for each timestep
272: .  -ts_monitor_lg_snes_iterations - monitor the number of SNES iterations needed for each timestep
273: -  -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true

275:    Notes:
276:    Use TSMonitorLGCtxDestroy() to destroy.

278:    One can provide a function that transforms the solution before plotting it with TSMonitorLGCtxSetTransform() or TSMonitorLGSetTransform()

280:    Many of the functions that control the monitoring have two forms: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a TS object as the
281:    first argument (if that TS object does not have a TSMonitorLGCtx associated with it the function call is ignored) and the second takes a TSMonitorLGCtx object
282:    as the first argument.

284:    One can control the names displayed for each solution or error variable with TSMonitorLGCtxSetVariableNames() or TSMonitorLGSetVariableNames()

286:    Level: intermediate

288: .seealso: TSMonitorLGTimeStep(), TSMonitorSet(), TSMonitorLGSolution(), TSMonitorLGError(), TSMonitorDefault(), VecView(),
289:            TSMonitorLGCtxCreate(), TSMonitorLGCtxSetVariableNames(), TSMonitorLGCtxGetVariableNames(),
290:            TSMonitorLGSetVariableNames(), TSMonitorLGGetVariableNames(), TSMonitorLGSetDisplayVariables(), TSMonitorLGCtxSetDisplayVariables(),
291:            TSMonitorLGCtxSetTransform(), TSMonitorLGSetTransform(), TSMonitorLGError(), TSMonitorLGSNESIterations(), TSMonitorLGKSPIterations(),
292:            TSMonitorEnvelopeCtxCreate(), TSMonitorEnvelopeGetBounds(), TSMonitorEnvelopeCtxDestroy(), TSMonitorEnvelop()

294: @*/
295: PetscErrorCode  TSMonitorLGCtxCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscInt howoften,TSMonitorLGCtx *ctx)
296: {
297:   PetscDraw      draw;

301:   PetscNew(ctx);
302:   PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
303:   PetscDrawSetFromOptions(draw);
304:   PetscDrawLGCreate(draw,1,&(*ctx)->lg);
305:   PetscDrawLGSetFromOptions((*ctx)->lg);
306:   PetscDrawDestroy(&draw);
307:   (*ctx)->howoften = howoften;
308:   return(0);
309: }

311: PetscErrorCode TSMonitorLGTimeStep(TS ts,PetscInt step,PetscReal ptime,Vec v,void *monctx)
312: {
313:   TSMonitorLGCtx ctx = (TSMonitorLGCtx) monctx;
314:   PetscReal      x   = ptime,y;

318:   if (step < 0) return(0); /* -1 indicates an interpolated solution */
319:   if (!step) {
320:     PetscDrawAxis axis;
321:     const char *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
322:     PetscDrawLGGetAxis(ctx->lg,&axis);
323:     PetscDrawAxisSetLabels(axis,"Timestep as function of time","Time",ylabel);
324:     PetscDrawLGReset(ctx->lg);
325:   }
326:   TSGetTimeStep(ts,&y);
327:   if (ctx->semilogy) y = PetscLog10Real(y);
328:   PetscDrawLGAddPoint(ctx->lg,&x,&y);
329:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
330:     PetscDrawLGDraw(ctx->lg);
331:     PetscDrawLGSave(ctx->lg);
332:   }
333:   return(0);
334: }

336: /*@C
337:    TSMonitorLGCtxDestroy - Destroys a line graph context that was created
338:    with TSMonitorLGCtxCreate().

340:    Collective on TSMonitorLGCtx

342:    Input Parameter:
343: .  ctx - the monitor context

345:    Level: intermediate

347: .seealso: TSMonitorLGCtxCreate(),  TSMonitorSet(), TSMonitorLGTimeStep();
348: @*/
349: PetscErrorCode  TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
350: {

354:   if ((*ctx)->transformdestroy) {
355:     ((*ctx)->transformdestroy)((*ctx)->transformctx);
356:   }
357:   PetscDrawLGDestroy(&(*ctx)->lg);
358:   PetscStrArrayDestroy(&(*ctx)->names);
359:   PetscStrArrayDestroy(&(*ctx)->displaynames);
360:   PetscFree((*ctx)->displayvariables);
361:   PetscFree((*ctx)->displayvalues);
362:   PetscFree(*ctx);
363:   return(0);
364: }

366: /*

368:   Creates a TS Monitor SPCtx for use with DM Swarm particle visualizations

370: */
371: PetscErrorCode TSMonitorSPCtxCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscInt howoften,TSMonitorSPCtx *ctx)
372: {
373:   PetscDraw      draw;

377:   PetscNew(ctx);
378:   PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
379:   PetscDrawSetFromOptions(draw);
380:   PetscDrawSPCreate(draw,1,&(*ctx)->sp);
381:   PetscDrawDestroy(&draw);
382:   (*ctx)->howoften = howoften;
383:   return(0);

385: }

387: /*
388:   Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate
389: */
390: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
391: {


396:   PetscDrawSPDestroy(&(*ctx)->sp);
397:   PetscFree(*ctx);

399:   return(0);

401: }

403: /*@C
404:    TSMonitorDrawSolution - Monitors progress of the TS solvers by calling
405:    VecView() for the solution at each timestep

407:    Collective on TS

409:    Input Parameters:
410: +  ts - the TS context
411: .  step - current time-step
412: .  ptime - current time
413: -  dummy - either a viewer or NULL

415:    Options Database:
416: .   -ts_monitor_draw_solution_initial - show initial solution as well as current solution

418:    Notes:
419:     the initial solution and current solution are not display with a common axis scaling so generally the option -ts_monitor_draw_solution_initial
420:        will look bad

422:    Level: intermediate

424: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView()
425: @*/
426: PetscErrorCode  TSMonitorDrawSolution(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dummy)
427: {
428:   PetscErrorCode   ierr;
429:   TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
430:   PetscDraw        draw;

433:   if (!step && ictx->showinitial) {
434:     if (!ictx->initialsolution) {
435:       VecDuplicate(u,&ictx->initialsolution);
436:     }
437:     VecCopy(u,ictx->initialsolution);
438:   }
439:   if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) return(0);

441:   if (ictx->showinitial) {
442:     PetscReal pause;
443:     PetscViewerDrawGetPause(ictx->viewer,&pause);
444:     PetscViewerDrawSetPause(ictx->viewer,0.0);
445:     VecView(ictx->initialsolution,ictx->viewer);
446:     PetscViewerDrawSetPause(ictx->viewer,pause);
447:     PetscViewerDrawSetHold(ictx->viewer,PETSC_TRUE);
448:   }
449:   VecView(u,ictx->viewer);
450:   if (ictx->showtimestepandtime) {
451:     PetscReal xl,yl,xr,yr,h;
452:     char      time[32];

454:     PetscViewerDrawGetDraw(ictx->viewer,0,&draw);
455:     PetscSNPrintf(time,32,"Timestep %d Time %g",(int)step,(double)ptime);
456:     PetscDrawGetCoordinates(draw,&xl,&yl,&xr,&yr);
457:     h    = yl + .95*(yr - yl);
458:     PetscDrawStringCentered(draw,.5*(xl+xr),h,PETSC_DRAW_BLACK,time);
459:     PetscDrawFlush(draw);
460:   }

462:   if (ictx->showinitial) {
463:     PetscViewerDrawSetHold(ictx->viewer,PETSC_FALSE);
464:   }
465:   return(0);
466: }

468: /*@C
469:    TSMonitorDrawSolutionPhase - Monitors progress of the TS solvers by plotting the solution as a phase diagram

471:    Collective on TS

473:    Input Parameters:
474: +  ts - the TS context
475: .  step - current time-step
476: .  ptime - current time
477: -  dummy - either a viewer or NULL

479:    Level: intermediate

481: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView()
482: @*/
483: PetscErrorCode  TSMonitorDrawSolutionPhase(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dummy)
484: {
485:   PetscErrorCode    ierr;
486:   TSMonitorDrawCtx  ictx = (TSMonitorDrawCtx)dummy;
487:   PetscDraw         draw;
488:   PetscDrawAxis     axis;
489:   PetscInt          n;
490:   PetscMPIInt       size;
491:   PetscReal         U0,U1,xl,yl,xr,yr,h;
492:   char              time[32];
493:   const PetscScalar *U;

496:   MPI_Comm_size(PetscObjectComm((PetscObject)ts),&size);
497:   if (size != 1) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_SUP,"Only allowed for sequential runs");
498:   VecGetSize(u,&n);
499:   if (n != 2) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_SUP,"Only for ODEs with two unknowns");

501:   PetscViewerDrawGetDraw(ictx->viewer,0,&draw);
502:   PetscViewerDrawGetDrawAxis(ictx->viewer,0,&axis);
503:   PetscDrawAxisGetLimits(axis,&xl,&xr,&yl,&yr);
504:   if (!step) {
505:     PetscDrawClear(draw);
506:     PetscDrawAxisDraw(axis);
507:   }

509:   VecGetArrayRead(u,&U);
510:   U0 = PetscRealPart(U[0]);
511:   U1 = PetscRealPart(U[1]);
512:   VecRestoreArrayRead(u,&U);
513:   if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) return(0);

515:   PetscDrawCollectiveBegin(draw);
516:   PetscDrawPoint(draw,U0,U1,PETSC_DRAW_BLACK);
517:   if (ictx->showtimestepandtime) {
518:     PetscDrawGetCoordinates(draw,&xl,&yl,&xr,&yr);
519:     PetscSNPrintf(time,32,"Timestep %d Time %g",(int)step,(double)ptime);
520:     h    = yl + .95*(yr - yl);
521:     PetscDrawStringCentered(draw,.5*(xl+xr),h,PETSC_DRAW_BLACK,time);
522:   }
523:   PetscDrawCollectiveEnd(draw);
524:   PetscDrawFlush(draw);
525:   PetscDrawPause(draw);
526:   PetscDrawSave(draw);
527:   return(0);
528: }

530: /*@C
531:    TSMonitorDrawCtxDestroy - Destroys the monitor context for TSMonitorDrawSolution()

533:    Collective on TS

535:    Input Parameters:
536: .    ctx - the monitor context

538:    Level: intermediate

540: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorDrawSolution(), TSMonitorDrawError()
541: @*/
542: PetscErrorCode  TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
543: {

547:   PetscViewerDestroy(&(*ictx)->viewer);
548:   VecDestroy(&(*ictx)->initialsolution);
549:   PetscFree(*ictx);
550:   return(0);
551: }

553: /*@C
554:    TSMonitorDrawCtxCreate - Creates the monitor context for TSMonitorDrawCtx

556:    Collective on TS

558:    Input Parameter:
559: .    ts - time-step context

561:    Output Patameter:
562: .    ctx - the monitor context

564:    Options Database:
565: .   -ts_monitor_draw_solution_initial - show initial solution as well as current solution

567:    Level: intermediate

569: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorDrawCtx()
570: @*/
571: PetscErrorCode  TSMonitorDrawCtxCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscInt howoften,TSMonitorDrawCtx *ctx)
572: {
573:   PetscErrorCode   ierr;

576:   PetscNew(ctx);
577:   PetscViewerDrawOpen(comm,host,label,x,y,m,n,&(*ctx)->viewer);
578:   PetscViewerSetFromOptions((*ctx)->viewer);

580:   (*ctx)->howoften    = howoften;
581:   (*ctx)->showinitial = PETSC_FALSE;
582:   PetscOptionsGetBool(NULL,NULL,"-ts_monitor_draw_solution_initial",&(*ctx)->showinitial,NULL);

584:   (*ctx)->showtimestepandtime = PETSC_FALSE;
585:   PetscOptionsGetBool(NULL,NULL,"-ts_monitor_draw_solution_show_time",&(*ctx)->showtimestepandtime,NULL);
586:   return(0);
587: }

589: /*@C
590:    TSMonitorDrawSolutionFunction - Monitors progress of the TS solvers by calling
591:    VecView() for the solution provided by TSSetSolutionFunction() at each timestep

593:    Collective on TS

595:    Input Parameters:
596: +  ts - the TS context
597: .  step - current time-step
598: .  ptime - current time
599: -  dummy - either a viewer or NULL

601:    Options Database:
602: .  -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided TSSetSolutionFunction()

604:    Level: intermediate

606: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSSetSolutionFunction()
607: @*/
608: PetscErrorCode  TSMonitorDrawSolutionFunction(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dummy)
609: {
610:   PetscErrorCode   ierr;
611:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
612:   PetscViewer      viewer = ctx->viewer;
613:   Vec              work;

616:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) return(0);
617:   VecDuplicate(u,&work);
618:   TSComputeSolutionFunction(ts,ptime,work);
619:   VecView(work,viewer);
620:   VecDestroy(&work);
621:   return(0);
622: }

624: /*@C
625:    TSMonitorDrawError - Monitors progress of the TS solvers by calling
626:    VecView() for the error at each timestep

628:    Collective on TS

630:    Input Parameters:
631: +  ts - the TS context
632: .  step - current time-step
633: .  ptime - current time
634: -  dummy - either a viewer or NULL

636:    Options Database:
637: .  -ts_monitor_draw_error - Monitor error graphically, requires user to have provided TSSetSolutionFunction()

639:    Level: intermediate

641: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSSetSolutionFunction()
642: @*/
643: PetscErrorCode  TSMonitorDrawError(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dummy)
644: {
645:   PetscErrorCode   ierr;
646:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
647:   PetscViewer      viewer = ctx->viewer;
648:   Vec              work;

651:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) return(0);
652:   VecDuplicate(u,&work);
653:   TSComputeSolutionFunction(ts,ptime,work);
654:   VecAXPY(work,-1.0,u);
655:   VecView(work,viewer);
656:   VecDestroy(&work);
657:   return(0);
658: }

660: /*@C
661:    TSMonitorSolution - Monitors progress of the TS solvers by VecView() for the solution at each timestep. Normally the viewer is a binary file or a PetscDraw object

663:    Collective on TS

665:    Input Parameters:
666: +  ts - the TS context
667: .  step - current time-step
668: .  ptime - current time
669: .  u - current state
670: -  vf - viewer and its format

672:    Level: intermediate

674: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView()
675: @*/
676: PetscErrorCode  TSMonitorSolution(TS ts,PetscInt step,PetscReal ptime,Vec u,PetscViewerAndFormat *vf)
677: {

681:   PetscViewerPushFormat(vf->viewer,vf->format);
682:   VecView(u,vf->viewer);
683:   PetscViewerPopFormat(vf->viewer);
684:   return(0);
685: }

687: /*@C
688:    TSMonitorSolutionVTK - Monitors progress of the TS solvers by VecView() for the solution at each timestep.

690:    Collective on TS

692:    Input Parameters:
693: +  ts - the TS context
694: .  step - current time-step
695: .  ptime - current time
696: .  u - current state
697: -  filenametemplate - string containing a format specifier for the integer time step (e.g. %03D)

699:    Level: intermediate

701:    Notes:
702:    The VTK format does not allow writing multiple time steps in the same file, therefore a different file will be written for each time step.
703:    These are named according to the file name template.

705:    This function is normally passed as an argument to TSMonitorSet() along with TSMonitorSolutionVTKDestroy().

707: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView()
708: @*/
709: PetscErrorCode TSMonitorSolutionVTK(TS ts,PetscInt step,PetscReal ptime,Vec u,void *filenametemplate)
710: {
712:   char           filename[PETSC_MAX_PATH_LEN];
713:   PetscViewer    viewer;

716:   if (step < 0) return(0); /* -1 indicates interpolated solution */
717:   PetscSNPrintf(filename,sizeof(filename),(const char*)filenametemplate,step);
718:   PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts),filename,FILE_MODE_WRITE,&viewer);
719:   VecView(u,viewer);
720:   PetscViewerDestroy(&viewer);
721:   return(0);
722: }

724: /*@C
725:    TSMonitorSolutionVTKDestroy - Destroy context for monitoring

727:    Collective on TS

729:    Input Parameters:
730: .  filenametemplate - string containing a format specifier for the integer time step (e.g. %03D)

732:    Level: intermediate

734:    Note:
735:    This function is normally passed to TSMonitorSet() along with TSMonitorSolutionVTK().

737: .seealso: TSMonitorSet(), TSMonitorSolutionVTK()
738: @*/
739: PetscErrorCode TSMonitorSolutionVTKDestroy(void *filenametemplate)
740: {

744:   PetscFree(*(char**)filenametemplate);
745:   return(0);
746: }

748: /*@C
749:    TSMonitorLGSolution - Monitors progress of the TS solvers by plotting each component of the solution vector
750:        in a time based line graph

752:    Collective on TS

754:    Input Parameters:
755: +  ts - the TS context
756: .  step - current time-step
757: .  ptime - current time
758: .  u - current solution
759: -  dctx - the TSMonitorLGCtx object that contains all the options for the monitoring, this is created with TSMonitorLGCtxCreate()

761:    Options Database:
762: .   -ts_monitor_lg_solution_variables

764:    Level: intermediate

766:    Notes:
767:     Each process in a parallel run displays its component solutions in a separate window

769: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGCtxCreate(), TSMonitorLGCtxSetVariableNames(), TSMonitorLGCtxGetVariableNames(),
770:            TSMonitorLGSetVariableNames(), TSMonitorLGGetVariableNames(), TSMonitorLGSetDisplayVariables(), TSMonitorLGCtxSetDisplayVariables(),
771:            TSMonitorLGCtxSetTransform(), TSMonitorLGSetTransform(), TSMonitorLGError(), TSMonitorLGSNESIterations(), TSMonitorLGKSPIterations(),
772:            TSMonitorEnvelopeCtxCreate(), TSMonitorEnvelopeGetBounds(), TSMonitorEnvelopeCtxDestroy(), TSMonitorEnvelop()
773: @*/
774: PetscErrorCode  TSMonitorLGSolution(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dctx)
775: {
776:   PetscErrorCode    ierr;
777:   TSMonitorLGCtx    ctx = (TSMonitorLGCtx)dctx;
778:   const PetscScalar *yy;
779:   Vec               v;

782:   if (step < 0) return(0); /* -1 indicates interpolated solution */
783:   if (!step) {
784:     PetscDrawAxis axis;
785:     PetscInt      dim;
786:     PetscDrawLGGetAxis(ctx->lg,&axis);
787:     PetscDrawAxisSetLabels(axis,"Solution as function of time","Time","Solution");
788:     if (!ctx->names) {
789:       PetscBool flg;
790:       /* user provides names of variables to plot but no names has been set so assume names are integer values */
791:       PetscOptionsHasName(((PetscObject)ts)->options,((PetscObject)ts)->prefix,"-ts_monitor_lg_solution_variables",&flg);
792:       if (flg) {
793:         PetscInt i,n;
794:         char     **names;
795:         VecGetSize(u,&n);
796:         PetscMalloc1(n+1,&names);
797:         for (i=0; i<n; i++) {
798:           PetscMalloc1(5,&names[i]);
799:           PetscSNPrintf(names[i],5,"%D",i);
800:         }
801:         names[n] = NULL;
802:         ctx->names = names;
803:       }
804:     }
805:     if (ctx->names && !ctx->displaynames) {
806:       char      **displaynames;
807:       PetscBool flg;
808:       VecGetLocalSize(u,&dim);
809:       PetscCalloc1(dim+1,&displaynames);
810:       PetscOptionsGetStringArray(((PetscObject)ts)->options,((PetscObject)ts)->prefix,"-ts_monitor_lg_solution_variables",displaynames,&dim,&flg);
811:       if (flg) {
812:         TSMonitorLGCtxSetDisplayVariables(ctx,(const char *const *)displaynames);
813:       }
814:       PetscStrArrayDestroy(&displaynames);
815:     }
816:     if (ctx->displaynames) {
817:       PetscDrawLGSetDimension(ctx->lg,ctx->ndisplayvariables);
818:       PetscDrawLGSetLegend(ctx->lg,(const char *const *)ctx->displaynames);
819:     } else if (ctx->names) {
820:       VecGetLocalSize(u,&dim);
821:       PetscDrawLGSetDimension(ctx->lg,dim);
822:       PetscDrawLGSetLegend(ctx->lg,(const char *const *)ctx->names);
823:     } else {
824:       VecGetLocalSize(u,&dim);
825:       PetscDrawLGSetDimension(ctx->lg,dim);
826:     }
827:     PetscDrawLGReset(ctx->lg);
828:   }

830:   if (!ctx->transform) v = u;
831:   else {(*ctx->transform)(ctx->transformctx,u,&v);}
832:   VecGetArrayRead(v,&yy);
833:   if (ctx->displaynames) {
834:     PetscInt i;
835:     for (i=0; i<ctx->ndisplayvariables; i++)
836:       ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
837:     PetscDrawLGAddCommonPoint(ctx->lg,ptime,ctx->displayvalues);
838:   } else {
839: #if defined(PETSC_USE_COMPLEX)
840:     PetscInt  i,n;
841:     PetscReal *yreal;
842:     VecGetLocalSize(v,&n);
843:     PetscMalloc1(n,&yreal);
844:     for (i=0; i<n; i++) yreal[i] = PetscRealPart(yy[i]);
845:     PetscDrawLGAddCommonPoint(ctx->lg,ptime,yreal);
846:     PetscFree(yreal);
847: #else
848:     PetscDrawLGAddCommonPoint(ctx->lg,ptime,yy);
849: #endif
850:   }
851:   VecRestoreArrayRead(v,&yy);
852:   if (ctx->transform) {VecDestroy(&v);}

854:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
855:     PetscDrawLGDraw(ctx->lg);
856:     PetscDrawLGSave(ctx->lg);
857:   }
858:   return(0);
859: }

861: /*@C
862:    TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot

864:    Collective on TS

866:    Input Parameters:
867: +  ts - the TS context
868: -  names - the names of the components, final string must be NULL

870:    Level: intermediate

872:    Notes:
873:     If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored

875: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetDisplayVariables(), TSMonitorLGCtxSetVariableNames()
876: @*/
877: PetscErrorCode  TSMonitorLGSetVariableNames(TS ts,const char * const *names)
878: {
879:   PetscErrorCode    ierr;
880:   PetscInt          i;

883:   for (i=0; i<ts->numbermonitors; i++) {
884:     if (ts->monitor[i] == TSMonitorLGSolution) {
885:       TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i],names);
886:       break;
887:     }
888:   }
889:   return(0);
890: }

892: /*@C
893:    TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot

895:    Collective on TS

897:    Input Parameters:
898: +  ts - the TS context
899: -  names - the names of the components, final string must be NULL

901:    Level: intermediate

903: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetDisplayVariables(), TSMonitorLGSetVariableNames()
904: @*/
905: PetscErrorCode  TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx,const char * const *names)
906: {
907:   PetscErrorCode    ierr;

910:   PetscStrArrayDestroy(&ctx->names);
911:   PetscStrArrayallocpy(names,&ctx->names);
912:   return(0);
913: }

915: /*@C
916:    TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot

918:    Collective on TS

920:    Input Parameter:
921: .  ts - the TS context

923:    Output Parameter:
924: .  names - the names of the components, final string must be NULL

926:    Level: intermediate

928:    Notes:
929:     If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored

931: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetDisplayVariables()
932: @*/
933: PetscErrorCode  TSMonitorLGGetVariableNames(TS ts,const char *const **names)
934: {
935:   PetscInt       i;

938:   *names = NULL;
939:   for (i=0; i<ts->numbermonitors; i++) {
940:     if (ts->monitor[i] == TSMonitorLGSolution) {
941:       TSMonitorLGCtx  ctx = (TSMonitorLGCtx) ts->monitorcontext[i];
942:       *names = (const char *const *)ctx->names;
943:       break;
944:     }
945:   }
946:   return(0);
947: }

949: /*@C
950:    TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor

952:    Collective on TS

954:    Input Parameters:
955: +  ctx - the TSMonitorLG context
956: -  displaynames - the names of the components, final string must be NULL

958:    Level: intermediate

960: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetVariableNames()
961: @*/
962: PetscErrorCode  TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx,const char * const *displaynames)
963: {
964:   PetscInt          j = 0,k;
965:   PetscErrorCode    ierr;

968:   if (!ctx->names) return(0);
969:   PetscStrArrayDestroy(&ctx->displaynames);
970:   PetscStrArrayallocpy(displaynames,&ctx->displaynames);
971:   while (displaynames[j]) j++;
972:   ctx->ndisplayvariables = j;
973:   PetscMalloc1(ctx->ndisplayvariables,&ctx->displayvariables);
974:   PetscMalloc1(ctx->ndisplayvariables,&ctx->displayvalues);
975:   j = 0;
976:   while (displaynames[j]) {
977:     k = 0;
978:     while (ctx->names[k]) {
979:       PetscBool flg;
980:       PetscStrcmp(displaynames[j],ctx->names[k],&flg);
981:       if (flg) {
982:         ctx->displayvariables[j] = k;
983:         break;
984:       }
985:       k++;
986:     }
987:     j++;
988:   }
989:   return(0);
990: }

992: /*@C
993:    TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor

995:    Collective on TS

997:    Input Parameters:
998: +  ts - the TS context
999: -  displaynames - the names of the components, final string must be NULL

1001:    Notes:
1002:     If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored

1004:    Level: intermediate

1006: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetVariableNames()
1007: @*/
1008: PetscErrorCode  TSMonitorLGSetDisplayVariables(TS ts,const char * const *displaynames)
1009: {
1010:   PetscInt          i;
1011:   PetscErrorCode    ierr;

1014:   for (i=0; i<ts->numbermonitors; i++) {
1015:     if (ts->monitor[i] == TSMonitorLGSolution) {
1016:       TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i],displaynames);
1017:       break;
1018:     }
1019:   }
1020:   return(0);
1021: }

1023: /*@C
1024:    TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed

1026:    Collective on TS

1028:    Input Parameters:
1029: +  ts - the TS context
1030: .  transform - the transform function
1031: .  destroy - function to destroy the optional context
1032: -  ctx - optional context used by transform function

1034:    Notes:
1035:     If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored

1037:    Level: intermediate

1039: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetVariableNames(), TSMonitorLGCtxSetTransform()
1040: @*/
1041: PetscErrorCode  TSMonitorLGSetTransform(TS ts,PetscErrorCode (*transform)(void*,Vec,Vec*),PetscErrorCode (*destroy)(void*),void *tctx)
1042: {
1043:   PetscInt          i;
1044:   PetscErrorCode    ierr;

1047:   for (i=0; i<ts->numbermonitors; i++) {
1048:     if (ts->monitor[i] == TSMonitorLGSolution) {
1049:       TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i],transform,destroy,tctx);
1050:     }
1051:   }
1052:   return(0);
1053: }

1055: /*@C
1056:    TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed

1058:    Collective on TSLGCtx

1060:    Input Parameters:
1061: +  ts - the TS context
1062: .  transform - the transform function
1063: .  destroy - function to destroy the optional context
1064: -  ctx - optional context used by transform function

1066:    Level: intermediate

1068: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetVariableNames(), TSMonitorLGSetTransform()
1069: @*/
1070: PetscErrorCode  TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx,PetscErrorCode (*transform)(void*,Vec,Vec*),PetscErrorCode (*destroy)(void*),void *tctx)
1071: {
1073:   ctx->transform    = transform;
1074:   ctx->transformdestroy = destroy;
1075:   ctx->transformctx = tctx;
1076:   return(0);
1077: }

1079: /*@C
1080:    TSMonitorLGError - Monitors progress of the TS solvers by plotting each component of the error
1081:        in a time based line graph

1083:    Collective on TS

1085:    Input Parameters:
1086: +  ts - the TS context
1087: .  step - current time-step
1088: .  ptime - current time
1089: .  u - current solution
1090: -  dctx - TSMonitorLGCtx object created with TSMonitorLGCtxCreate()

1092:    Level: intermediate

1094:    Notes:
1095:     Each process in a parallel run displays its component errors in a separate window

1097:    The user must provide the solution using TSSetSolutionFunction() to use this monitor.

1099:    Options Database Keys:
1100: .  -ts_monitor_lg_error - create a graphical monitor of error history

1102: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSSetSolutionFunction()
1103: @*/
1104: PetscErrorCode  TSMonitorLGError(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dummy)
1105: {
1106:   PetscErrorCode    ierr;
1107:   TSMonitorLGCtx    ctx = (TSMonitorLGCtx)dummy;
1108:   const PetscScalar *yy;
1109:   Vec               y;

1112:   if (!step) {
1113:     PetscDrawAxis axis;
1114:     PetscInt      dim;
1115:     PetscDrawLGGetAxis(ctx->lg,&axis);
1116:     PetscDrawAxisSetLabels(axis,"Error in solution as function of time","Time","Error");
1117:     VecGetLocalSize(u,&dim);
1118:     PetscDrawLGSetDimension(ctx->lg,dim);
1119:     PetscDrawLGReset(ctx->lg);
1120:   }
1121:   VecDuplicate(u,&y);
1122:   TSComputeSolutionFunction(ts,ptime,y);
1123:   VecAXPY(y,-1.0,u);
1124:   VecGetArrayRead(y,&yy);
1125: #if defined(PETSC_USE_COMPLEX)
1126:   {
1127:     PetscReal *yreal;
1128:     PetscInt  i,n;
1129:     VecGetLocalSize(y,&n);
1130:     PetscMalloc1(n,&yreal);
1131:     for (i=0; i<n; i++) yreal[i] = PetscRealPart(yy[i]);
1132:     PetscDrawLGAddCommonPoint(ctx->lg,ptime,yreal);
1133:     PetscFree(yreal);
1134:   }
1135: #else
1136:   PetscDrawLGAddCommonPoint(ctx->lg,ptime,yy);
1137: #endif
1138:   VecRestoreArrayRead(y,&yy);
1139:   VecDestroy(&y);
1140:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1141:     PetscDrawLGDraw(ctx->lg);
1142:     PetscDrawLGSave(ctx->lg);
1143:   }
1144:   return(0);
1145: }

1147: /*@C
1148:    TSMonitorSPSwarmSolution - Graphically displays phase plots of DMSwarm particles on a scatter plot

1150:    Input Parameters:
1151: +  ts - the TS context
1152: .  step - current time-step
1153: .  ptime - current time
1154: .  u - current solution
1155: -  dctx - the TSMonitorSPCtx object that contains all the options for the monitoring, this is created with TSMonitorSPCtxCreate()

1157:    Options Database:
1158: .   -ts_monitor_sp_swarm

1160:    Level: intermediate

1162: @*/
1163: PetscErrorCode TSMonitorSPSwarmSolution(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dctx)
1164: {
1165:   PetscErrorCode    ierr;
1166:   TSMonitorSPCtx    ctx = (TSMonitorSPCtx)dctx;
1167:   const PetscScalar *yy;
1168:   PetscReal       *y,*x;
1169:   PetscInt          Np, p, dim=2;
1170:   DM                dm;

1173:   if (step < 0) return(0); /* -1 indicates interpolated solution */
1174:   if (!step) {
1175:     PetscDrawAxis axis;
1176:     PetscDrawSPGetAxis(ctx->sp,&axis);
1177:     PetscDrawAxisSetLabels(axis,"Particles","X","Y");
1178:     PetscDrawAxisSetLimits(axis, -5, 5, -5, 5);
1179:     PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE);
1180:     TSGetDM(ts, &dm);
1181:     DMGetDimension(dm, &dim);
1182:     if (dim!=2) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Dimensions improper for monitor arguments! Current support: two dimensions.");
1183:     VecGetLocalSize(u, &Np);
1184:     Np /= 2*dim;
1185:     PetscDrawSPSetDimension(ctx->sp, Np);
1186:     PetscDrawSPReset(ctx->sp);
1187:   }
1188:   VecGetLocalSize(u, &Np);
1189:   Np /= 2*dim;
1190:   VecGetArrayRead(u,&yy);
1191:   PetscMalloc2(Np, &x, Np, &y);
1192:   /* get points from solution vector */
1193:   for (p=0; p<Np; ++p){
1194:     x[p] = PetscRealPart(yy[2*dim*p]);
1195:     y[p] = PetscRealPart(yy[2*dim*p+1]);
1196:   }
1197:   VecRestoreArrayRead(u,&yy);
1198:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1199:     PetscDrawSPAddPoint(ctx->sp,x,y);
1200:     PetscDrawSPDraw(ctx->sp,PETSC_FALSE);
1201:     PetscDrawSPSave(ctx->sp);
1202:   }
1203:   PetscFree2(x, y);
1204:   return(0);
1205: }



1209: /*@C
1210:    TSMonitorError - Monitors progress of the TS solvers by printing the 2 norm of the error at each timestep

1212:    Collective on TS

1214:    Input Parameters:
1215: +  ts - the TS context
1216: .  step - current time-step
1217: .  ptime - current time
1218: .  u - current solution
1219: -  dctx - unused context

1221:    Level: intermediate

1223:    The user must provide the solution using TSSetSolutionFunction() to use this monitor.

1225:    Options Database Keys:
1226: .  -ts_monitor_error - create a graphical monitor of error history

1228: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSSetSolutionFunction()
1229: @*/
1230: PetscErrorCode  TSMonitorError(TS ts,PetscInt step,PetscReal ptime,Vec u,PetscViewerAndFormat *vf)
1231: {
1232:   PetscErrorCode    ierr;
1233:   Vec               y;
1234:   PetscReal         nrm;
1235:   PetscBool         flg;

1238:   VecDuplicate(u,&y);
1239:   TSComputeSolutionFunction(ts,ptime,y);
1240:   VecAXPY(y,-1.0,u);
1241:   PetscObjectTypeCompare((PetscObject)vf->viewer,PETSCVIEWERASCII,&flg);
1242:   if (flg) {
1243:     VecNorm(y,NORM_2,&nrm);
1244:     PetscViewerASCIIPrintf(vf->viewer,"2-norm of error %g\n",(double)nrm);
1245:   }
1246:   PetscObjectTypeCompare((PetscObject)vf->viewer,PETSCVIEWERDRAW,&flg);
1247:   if (flg) {
1248:     VecView(y,vf->viewer);
1249:   }
1250:   VecDestroy(&y);
1251:   return(0);
1252: }

1254: PetscErrorCode TSMonitorLGSNESIterations(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1255: {
1256:   TSMonitorLGCtx ctx = (TSMonitorLGCtx) monctx;
1257:   PetscReal      x   = ptime,y;
1259:   PetscInt       its;

1262:   if (n < 0) return(0); /* -1 indicates interpolated solution */
1263:   if (!n) {
1264:     PetscDrawAxis axis;
1265:     PetscDrawLGGetAxis(ctx->lg,&axis);
1266:     PetscDrawAxisSetLabels(axis,"Nonlinear iterations as function of time","Time","SNES Iterations");
1267:     PetscDrawLGReset(ctx->lg);
1268:     ctx->snes_its = 0;
1269:   }
1270:   TSGetSNESIterations(ts,&its);
1271:   y    = its - ctx->snes_its;
1272:   PetscDrawLGAddPoint(ctx->lg,&x,&y);
1273:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1274:     PetscDrawLGDraw(ctx->lg);
1275:     PetscDrawLGSave(ctx->lg);
1276:   }
1277:   ctx->snes_its = its;
1278:   return(0);
1279: }

1281: PetscErrorCode TSMonitorLGKSPIterations(TS ts,PetscInt n,PetscReal ptime,Vec v,void *monctx)
1282: {
1283:   TSMonitorLGCtx ctx = (TSMonitorLGCtx) monctx;
1284:   PetscReal      x   = ptime,y;
1286:   PetscInt       its;

1289:   if (n < 0) return(0); /* -1 indicates interpolated solution */
1290:   if (!n) {
1291:     PetscDrawAxis axis;
1292:     PetscDrawLGGetAxis(ctx->lg,&axis);
1293:     PetscDrawAxisSetLabels(axis,"Linear iterations as function of time","Time","KSP Iterations");
1294:     PetscDrawLGReset(ctx->lg);
1295:     ctx->ksp_its = 0;
1296:   }
1297:   TSGetKSPIterations(ts,&its);
1298:   y    = its - ctx->ksp_its;
1299:   PetscDrawLGAddPoint(ctx->lg,&x,&y);
1300:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1301:     PetscDrawLGDraw(ctx->lg);
1302:     PetscDrawLGSave(ctx->lg);
1303:   }
1304:   ctx->ksp_its = its;
1305:   return(0);
1306: }

1308: /*@C
1309:    TSMonitorEnvelopeCtxCreate - Creates a context for use with TSMonitorEnvelope()

1311:    Collective on TS

1313:    Input Parameters:
1314: .  ts  - the ODE solver object

1316:    Output Parameter:
1317: .  ctx - the context

1319:    Level: intermediate

1321: .seealso: TSMonitorLGTimeStep(), TSMonitorSet(), TSMonitorLGSolution(), TSMonitorLGError()

1323: @*/
1324: PetscErrorCode  TSMonitorEnvelopeCtxCreate(TS ts,TSMonitorEnvelopeCtx *ctx)
1325: {

1329:   PetscNew(ctx);
1330:   return(0);
1331: }

1333: /*@C
1334:    TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution

1336:    Collective on TS

1338:    Input Parameters:
1339: +  ts - the TS context
1340: .  step - current time-step
1341: .  ptime - current time
1342: .  u  - current solution
1343: -  dctx - the envelope context

1345:    Options Database:
1346: .  -ts_monitor_envelope

1348:    Level: intermediate

1350:    Notes:
1351:     after a solve you can use TSMonitorEnvelopeGetBounds() to access the envelope

1353: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorEnvelopeGetBounds(), TSMonitorEnvelopeCtxCreate()
1354: @*/
1355: PetscErrorCode  TSMonitorEnvelope(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dctx)
1356: {
1357:   PetscErrorCode       ierr;
1358:   TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;

1361:   if (!ctx->max) {
1362:     VecDuplicate(u,&ctx->max);
1363:     VecDuplicate(u,&ctx->min);
1364:     VecCopy(u,ctx->max);
1365:     VecCopy(u,ctx->min);
1366:   } else {
1367:     VecPointwiseMax(ctx->max,u,ctx->max);
1368:     VecPointwiseMin(ctx->min,u,ctx->min);
1369:   }
1370:   return(0);
1371: }

1373: /*@C
1374:    TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution

1376:    Collective on TS

1378:    Input Parameter:
1379: .  ts - the TS context

1381:    Output Parameter:
1382: +  max - the maximum values
1383: -  min - the minimum values

1385:    Notes:
1386:     If the TS does not have a TSMonitorEnvelopeCtx associated with it then this function is ignored

1388:    Level: intermediate

1390: .seealso: TSMonitorSet(), TSMonitorDefault(), VecView(), TSMonitorLGSetDisplayVariables()
1391: @*/
1392: PetscErrorCode  TSMonitorEnvelopeGetBounds(TS ts,Vec *max,Vec *min)
1393: {
1394:   PetscInt i;

1397:   if (max) *max = NULL;
1398:   if (min) *min = NULL;
1399:   for (i=0; i<ts->numbermonitors; i++) {
1400:     if (ts->monitor[i] == TSMonitorEnvelope) {
1401:       TSMonitorEnvelopeCtx  ctx = (TSMonitorEnvelopeCtx) ts->monitorcontext[i];
1402:       if (max) *max = ctx->max;
1403:       if (min) *min = ctx->min;
1404:       break;
1405:     }
1406:   }
1407:   return(0);
1408: }

1410: /*@C
1411:    TSMonitorEnvelopeCtxDestroy - Destroys a context that was created  with TSMonitorEnvelopeCtxCreate().

1413:    Collective on TSMonitorEnvelopeCtx

1415:    Input Parameter:
1416: .  ctx - the monitor context

1418:    Level: intermediate

1420: .seealso: TSMonitorLGCtxCreate(),  TSMonitorSet(), TSMonitorLGTimeStep()
1421: @*/
1422: PetscErrorCode  TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1423: {

1427:   VecDestroy(&(*ctx)->min);
1428:   VecDestroy(&(*ctx)->max);
1429:   PetscFree(*ctx);
1430:   return(0);
1431: }

1433: /*@C
1434:   TSDMSwarmMonitorMoments - Monitors the first three moments of a DMSarm being evolved by the TS

1436:   Not collective

1438:   Input Parameters:
1439: + ts   - the TS context
1440: . step - current timestep
1441: . t    - current time
1442: . u    - current solution
1443: - ctx  - not used

1445:   Options Database:
1446: . -ts_dmswarm_monitor_moments

1448:   Level: intermediate

1450:   Notes:
1451:   This requires a DMSwarm be attached to the TS.

1453: .seealso: TSMonitorSet(), TSMonitorDefault(), DMSWARM
1454: @*/
1455: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1456: {
1457:   DM                 sw;
1458:   const PetscScalar *u;
1459:   PetscReal          m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1460:   PetscInt           dim, d, Np, p;
1461:   MPI_Comm           comm;
1462:   PetscErrorCode     ierr;

1465:   TSGetDM(ts, &sw);
1466:   if (!sw || step%ts->monitorFrequency != 0) return(0);
1467:   PetscObjectGetComm((PetscObject) ts, &comm);
1468:   DMGetDimension(sw, &dim);
1469:   VecGetLocalSize(U, &Np);
1470:   Np  /= dim;
1471:   VecGetArrayRead(U, &u);
1472:   for (p = 0; p < Np; ++p) {
1473:     for (d = 0; d < dim; ++d) {
1474:       totE      += PetscRealPart(u[p*dim+d]*u[p*dim+d]);
1475:       totMom[d] += PetscRealPart(u[p*dim+d]);
1476:     }
1477:   }
1478:   VecRestoreArrayRead(U, &u);
1479:   for (d = 0; d < dim; ++d) totMom[d] *= m;
1480:   totE *= 0.5*m;
1481:   PetscPrintf(comm, "Step %4D Total Energy: %10.8lf", step, (double) totE);
1482:   for (d = 0; d < dim; ++d) {PetscPrintf(comm, "    Total Momentum %c: %10.8lf", 'x'+d, (double) totMom[d]);}
1483:   PetscPrintf(comm, "\n");
1484:   return(0);
1485: }