Actual source code: tsmon.c

  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;

236:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
237:   PetscViewerPushFormat(viewer,vf->format);
238:   if (iascii) {
239:     VecMax(v,NULL,&max);
240:     VecMin(v,NULL,&min);
241:     PetscViewerASCIIAddTab(viewer,((PetscObject)ts)->tablevel);
242:     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);
243:     PetscViewerASCIISubtractTab(viewer,((PetscObject)ts)->tablevel);
244:   }
245:   PetscViewerPopFormat(viewer);
246:   return(0);
247: }

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

253:    Collective on TS

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

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

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

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

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

279:    Many of the functions that control the monitoring have two forms: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a TS object as the
280:    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
281:    as the first argument.

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

285:    Level: intermediate

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

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

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

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

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

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

339:    Collective on TSMonitorLGCtx

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

344:    Level: intermediate

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

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

365: /*

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

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

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

384: }

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


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

398:   return(0);

400: }

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

406:    Collective on TS

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

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

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

421:    Level: intermediate

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

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

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

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

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

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

470:    Collective on TS

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

478:    Level: intermediate

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

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

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

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

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

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

532:    Collective on TS

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

537:    Level: intermediate

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

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

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

555:    Collective on TS

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

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

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

566:    Level: intermediate

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

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

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

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

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

592:    Collective on TS

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

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

603:    Level: intermediate

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

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

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

627:    Collective on TS

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

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

638:    Level: intermediate

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

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

659: /*@C
660:    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

662:    Collective on TS

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

671:    Level: intermediate

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

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

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

689:    Collective on TS

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

698:    Level: intermediate

700:    Notes:
701:    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.
702:    These are named according to the file name template.

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

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

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

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

726:    Collective on TS

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

731:    Level: intermediate

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

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

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

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

751:    Collective on TS

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

760:    Options Database:
761: .   -ts_monitor_lg_solution_variables

763:    Level: intermediate

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

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

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

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

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

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

863:    Collective on TS

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

869:    Level: intermediate

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

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

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

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

894:    Collective on TS

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

900:    Level: intermediate

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

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

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

917:    Collective on TS

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

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

925:    Level: intermediate

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

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

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

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

951:    Collective on TS

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

957:    Level: intermediate

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

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

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

994:    Collective on TS

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

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

1003:    Level: intermediate

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

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

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

1025:    Collective on TS

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

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

1036:    Level: intermediate

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

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

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

1057:    Collective on TSLGCtx

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

1065:    Level: intermediate

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

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

1082:    Collective on TS

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

1091:    Level: intermediate

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

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

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

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

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

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

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

1156:    Options Database:
1157: .   -ts_monitor_sp_swarm

1159:    Level: intermediate

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

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

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

1209:    Collective on TS

1211:    Input Parameters:
1212: +  ts - the TS context
1213: .  step - current time-step
1214: .  ptime - current time
1215: .  u - current solution
1216: -  dctx - unused context

1218:    Level: intermediate

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

1222:    Options Database Keys:
1223: .  -ts_monitor_error - create a graphical monitor of error history

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

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

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

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

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

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

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

1308:    Collective on TS

1310:    Input Parameters:
1311: .  ts  - the ODE solver object

1313:    Output Parameter:
1314: .  ctx - the context

1316:    Level: intermediate

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

1320: @*/
1321: PetscErrorCode  TSMonitorEnvelopeCtxCreate(TS ts,TSMonitorEnvelopeCtx *ctx)
1322: {

1326:   PetscNew(ctx);
1327:   return(0);
1328: }

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

1333:    Collective on TS

1335:    Input Parameters:
1336: +  ts - the TS context
1337: .  step - current time-step
1338: .  ptime - current time
1339: .  u  - current solution
1340: -  dctx - the envelope context

1342:    Options Database:
1343: .  -ts_monitor_envelope

1345:    Level: intermediate

1347:    Notes:
1348:     after a solve you can use TSMonitorEnvelopeGetBounds() to access the envelope

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

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

1370: /*@C
1371:    TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution

1373:    Collective on TS

1375:    Input Parameter:
1376: .  ts - the TS context

1378:    Output Parameters:
1379: +  max - the maximum values
1380: -  min - the minimum values

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

1385:    Level: intermediate

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

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

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

1410:    Collective on TSMonitorEnvelopeCtx

1412:    Input Parameter:
1413: .  ctx - the monitor context

1415:    Level: intermediate

1417: .seealso: TSMonitorLGCtxCreate(),  TSMonitorSet(), TSMonitorLGTimeStep()
1418: @*/
1419: PetscErrorCode  TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1420: {

1424:   VecDestroy(&(*ctx)->min);
1425:   VecDestroy(&(*ctx)->max);
1426:   PetscFree(*ctx);
1427:   return(0);
1428: }

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

1433:   Not collective

1435:   Input Parameters:
1436: + ts   - the TS context
1437: . step - current timestep
1438: . t    - current time
1439: . u    - current solution
1440: - ctx  - not used

1442:   Options Database:
1443: . -ts_dmswarm_monitor_moments

1445:   Level: intermediate

1447:   Notes:
1448:   This requires a DMSwarm be attached to the TS.

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

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