Actual source code: tsmon.c

      1: #include <petsc/private/tsimpl.h>
      2: #include <petscdm.h>
      3: #include <petscds.h>
      4: #include <petscdmswarm.h>
      5: #include <petscdraw.h>

      7: /*@C
      8:   TSMonitor - Runs all user-provided monitor routines set using `TSMonitorSet()`

     10:   Collective

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

     18:   Level: developer

     20:   Notes:
     21:   `TSMonitor()` is typically used automatically within the time stepping implementations.
     22:   Users would almost never call this routine directly.

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

     26: .seealso: `TS`, `TSMonitorSet()`, `TSMonitorSetFromOptions()`
     27: @*/
     28: PetscErrorCode TSMonitor(TS ts, PetscInt step, PetscReal ptime, Vec u)
     29: {
     30:   DM       dm;
     31:   PetscInt i, n = ts->numbermonitors;

     33:   PetscFunctionBegin;

     37:   PetscCall(TSGetDM(ts, &dm));
     38:   PetscCall(DMSetOutputSequenceNumber(dm, step, ptime));

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

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

     49:   Collective

     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, this must use a `PetscViewerFormat` as its context
     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: [](ch_ts), `TS`, `TSMonitorSet()`, `PetscOptionsCreateViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
     62:           `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
     63:           `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`,
     64:           `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
     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:   PetscViewer       viewer;
     72:   PetscViewerFormat format;
     73:   PetscBool         flg;

     75:   PetscFunctionBegin;
     76:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)ts), ((PetscObject)ts)->options, ((PetscObject)ts)->prefix, name, &viewer, &format, &flg));
     77:   if (flg) {
     78:     PetscViewerAndFormat *vf;
     79:     char                  interval_key[1024];

     81:     PetscCall(PetscSNPrintf(interval_key, sizeof interval_key, "%s_interval", name));
     82:     PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
     83:     vf->view_interval = 1;
     84:     PetscCall(PetscOptionsGetInt(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, interval_key, &vf->view_interval, NULL));

     86:     PetscCall(PetscViewerDestroy(&viewer));
     87:     if (monitorsetup) PetscCall((*monitorsetup)(ts, vf));
     88:     PetscCall(TSMonitorSet(ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, void *))monitor, vf, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy));
     89:   }
     90:   PetscFunctionReturn(PETSC_SUCCESS);
     91: }

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

     97:   Logically Collective

     99:   Input Parameters:
    100: + ts       - the `TS` context obtained from `TSCreate()`
    101: . monitor  - monitoring routine
    102: . mctx     - [optional] user-defined context for private data for the monitor routine (use `NULL` if no context is desired)
    103: - mdestroy - [optional] routine that frees monitor context (may be `NULL`), see `PetscCtxDestroyFn` for the calling sequence

    105:   Calling sequence of `monitor`:
    106: + ts    - the `TS` context
    107: . 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)
    108: . time  - current time
    109: . u     - current iterate
    110: - ctx   - [optional] monitoring context

    112:   Level: intermediate

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

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

    120: .seealso: [](ch_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
    121:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
    122:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`,  `PetscCtxDestroyFn`
    123: @*/
    124: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS ts, PetscInt steps, PetscReal time, Vec u, void *ctx), void *mctx, PetscCtxDestroyFn *mdestroy)
    125: {
    126:   PetscFunctionBegin;
    128:   for (PetscInt i = 0; i < ts->numbermonitors; i++) {
    129:     PetscBool identical;

    131:     PetscCall(PetscMonitorCompare((PetscErrorCode (*)(void))(PetscVoidFn *)monitor, mctx, mdestroy, (PetscErrorCode (*)(void))(PetscVoidFn *)ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical));
    132:     if (identical) PetscFunctionReturn(PETSC_SUCCESS);
    133:   }
    134:   PetscCheck(ts->numbermonitors < MAXTSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set");
    135:   ts->monitor[ts->numbermonitors]          = monitor;
    136:   ts->monitordestroy[ts->numbermonitors]   = mdestroy;
    137:   ts->monitorcontext[ts->numbermonitors++] = mctx;
    138:   PetscFunctionReturn(PETSC_SUCCESS);
    139: }

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

    144:   Logically Collective

    146:   Input Parameter:
    147: . ts - the `TS` context obtained from `TSCreate()`

    149:   Level: intermediate

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

    154: .seealso: [](ch_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
    155: @*/
    156: PetscErrorCode TSMonitorCancel(TS ts)
    157: {
    158:   PetscInt i;

    160:   PetscFunctionBegin;
    162:   for (i = 0; i < ts->numbermonitors; i++) {
    163:     if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i]));
    164:   }
    165:   ts->numbermonitors = 0;
    166:   PetscFunctionReturn(PETSC_SUCCESS);
    167: }

    169: /*@C
    170:   TSMonitorDefault - The default monitor, prints the timestep and time for each step

    172:   Input Parameters:
    173: + ts    - the `TS` context
    174: . step  - 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)
    175: . ptime - current time
    176: . v     - current iterate
    177: - vf    - the viewer and format

    179:   Options Database Key:
    180: . -ts_monitor - monitors the time integration

    182:   Level: intermediate

    184:   Notes:
    185:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    186:   to be used during the `TS` integration.

    188: .seealso: [](ch_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorWallClockTime()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
    189:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
    190:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`
    191: @*/
    192: PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
    193: {
    194:   PetscViewer viewer = vf->viewer;
    195:   PetscBool   isascii, ibinary;

    197:   PetscFunctionBegin;
    199:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
    200:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary));
    201:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
    202:   if (isascii) {
    203:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
    204:     if (step == -1) { /* this indicates it is an interpolated solution */
    205:       PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps));
    206:     } else {
    207:       PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n"));
    208:     }
    209:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
    210:   } else if (ibinary) {
    211:     PetscMPIInt rank;
    212:     PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
    213:     if (rank == 0) {
    214:       PetscBool skipHeader;
    215:       PetscInt  classid = REAL_FILE_CLASSID;

    217:       PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
    218:       if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
    219:       PetscCall(PetscRealView(1, &ptime, viewer));
    220:     } else {
    221:       PetscCall(PetscRealView(0, &ptime, viewer));
    222:     }
    223:   }
    224:   PetscCall(PetscViewerPopFormat(viewer));
    225:   PetscFunctionReturn(PETSC_SUCCESS);
    226: }

    228: typedef struct {
    229:   PetscLogDouble time_start;
    230:   PetscLogDouble time_last;
    231:   PetscInt       snes_its;
    232:   PetscInt       ksp_its;
    233: } *TSMonitorWallClockTimeContext;

    235: /*@C
    236:   TSMonitorWallClockTimeSetUp - Setup routine passed to `TSMonitorSetFromOptions()` when using `-ts_monitor_wall_clock_time`

    238:   Input Parameters:
    239: + ts - the `TS` context
    240: - vf - the viewer and format

    242:   Level: intermediate

    244:   Note:
    245:   This is not called directly by users, rather one calls `TSMonitorSetFromOptions()`, with `TSMonitorWallClockTime()` and this function as arguments, to cause the monitor
    246:   to be used during the `TS` integration.

    248: .seealso: [](ch_ts), `TSMonitorSet()`
    249: @*/
    250: PetscErrorCode TSMonitorWallClockTimeSetUp(TS ts, PetscViewerAndFormat *vf)
    251: {
    252:   TSMonitorWallClockTimeContext speed;

    254:   PetscFunctionBegin;
    255:   PetscCall(PetscNew(&speed));
    256:   speed->time_start = PETSC_DECIDE;
    257:   vf->data_destroy  = PetscCtxDestroyDefault;
    258:   vf->data          = speed;
    259:   PetscFunctionReturn(PETSC_SUCCESS);
    260: }

    262: /*@C
    263:   TSMonitorWallClockTime - Monitor wall-clock time, KSP iterations, and SNES iterations per step.

    265:   Input Parameters:
    266: + ts    - the `TS` context
    267: . step  - 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)
    268: . ptime - current time
    269: . v     - current solution
    270: - vf    - the viewer and format

    272:   Options Database Key:
    273: . -ts_monitor_wall_clock_time - Monitor wall-clock time, KSP iterations, and SNES iterations per step.

    275:   Level: intermediate

    277:   Note:
    278:   This is not called directly by users, rather one calls `TSMonitorSetFromOptions()`, with this function and `TSMonitorWallClockTimeSetUp()` as arguments, to cause the monitor
    279:   to be used during the `TS` integration.

    281: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
    282:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
    283:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
    284: @*/
    285: PetscErrorCode TSMonitorWallClockTime(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
    286: {
    287:   PetscViewer                   viewer = vf->viewer;
    288:   TSMonitorWallClockTimeContext speed  = (TSMonitorWallClockTimeContext)vf->data;
    289:   PetscBool                     isascii;
    290:   PetscLogDouble                now;
    291:   PetscInt                      snes_its, ksp_its;

    293:   PetscFunctionBegin;
    295:   PetscCall(PetscTime(&now));
    296:   if (speed->time_start == PETSC_DECIDE) {
    297:     speed->time_start = now;
    298:     speed->time_last  = now;
    299:   }
    300:   PetscCall(TSGetSNESIterations(ts, &snes_its));
    301:   PetscCall(TSGetKSPIterations(ts, &ksp_its));
    302:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
    303:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
    304:   if (isascii) {
    305:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
    306:     PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s elapsed %.6f of %.6f snes %" PetscInt_FMT " ksp %" PetscInt_FMT "\n", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)" : "", now - speed->time_last,
    307:                                      now - speed->time_start, snes_its - speed->snes_its, ksp_its - speed->ksp_its));
    308:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
    309:   }
    310:   PetscCall(PetscViewerPopFormat(viewer));
    311:   speed->time_last = now;
    312:   speed->snes_its  = snes_its;
    313:   speed->ksp_its   = ksp_its;
    314:   PetscFunctionReturn(PETSC_SUCCESS);
    315: }

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

    320:   Input Parameters:
    321: + ts    - the `TS` context
    322: . step  - 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)
    323: . ptime - current time
    324: . v     - current iterate
    325: - vf    - the viewer and format

    327:   Level: intermediate

    329:   Note:
    330:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    331:   to be used during the `TS` integration.

    333: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`
    334: @*/
    335: PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
    336: {
    337:   PetscViewer viewer = vf->viewer;
    338:   PetscBool   isascii;
    339:   PetscReal   max, min;

    341:   PetscFunctionBegin;
    343:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
    344:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
    345:   if (isascii) {
    346:     PetscCall(VecMax(v, NULL, &max));
    347:     PetscCall(VecMin(v, NULL, &min));
    348:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
    349:     PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " 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));
    350:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
    351:   }
    352:   PetscCall(PetscViewerPopFormat(viewer));
    353:   PetscFunctionReturn(PETSC_SUCCESS);
    354: }

    356: /*@C
    357:   TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with
    358:   `TS` to monitor the solution process graphically in various ways

    360:   Collective

    362:   Input Parameters:
    363: + comm     - the MPI communicator to use
    364: . host     - the X display to open, or `NULL` for the local machine
    365: . label    - the title to put in the title bar
    366: . x        - the x screen coordinates of the upper left coordinate of the window
    367: . y        - the y screen coordinates of the upper left coordinate of the window
    368: . m        - the screen width in pixels
    369: . n        - the screen height in pixels
    370: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time

    372:   Output Parameter:
    373: . ctx - the context

    375:   Options Database Keys:
    376: + -ts_monitor_lg_timestep        - automatically sets line graph monitor
    377: . -ts_monitor_lg_timestep_log    - automatically sets line graph monitor
    378: . -ts_monitor_lg_solution        - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`)
    379: . -ts_monitor_lg_error           - monitor the error
    380: . -ts_monitor_lg_ksp_iterations  - monitor the number of `KSP` iterations needed for each timestep
    381: . -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep
    382: - -lg_use_markers <true,false>   - mark the data points (at each time step) on the plot; default is true

    384:   Level: intermediate

    386:   Notes:
    387:   Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed.

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

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

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

    397: .seealso: [](ch_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
    398:           `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
    399:           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
    400:           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
    401:           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
    402: @*/
    403: PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
    404: {
    405:   PetscDraw draw;

    407:   PetscFunctionBegin;
    408:   PetscCall(PetscNew(ctx));
    409:   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
    410:   PetscCall(PetscDrawSetFromOptions(draw));
    411:   PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg));
    412:   PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg));
    413:   PetscCall(PetscDrawDestroy(&draw));
    414:   (*ctx)->howoften = howoften;
    415:   PetscFunctionReturn(PETSC_SUCCESS);
    416: }

    418: /*@C
    419:   TSMonitorLGTimeStep - Monitors a `TS` by printing the time-steps

    421:   Collective

    423:   Input Parameters:
    424: + ts     - the time integrator
    425: . step   - the current time step
    426: . ptime  - the current time
    427: . v      - the current state
    428: - monctx - the monitor context obtained with `TSMonitorLGCtxCreate()`

    430:   Level: advanced

    432:   Note:
    433:   This is not called directly by users, rather one calls `TSMonitorSet()` along the `ctx` created by `TSMonitorLGCtxCreate()`
    434:   and `TSMonitorLGCtxDestroy()`

    436: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGCtxDestroy()`
    437: @*/
    438: PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx)
    439: {
    440:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
    441:   PetscReal      x   = ptime, y;

    443:   PetscFunctionBegin;
    444:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */
    445:   if (!step) {
    446:     PetscDrawAxis axis;
    447:     const char   *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
    448:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
    449:     PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel));
    450:     PetscCall(PetscDrawLGReset(ctx->lg));
    451:   }
    452:   PetscCall(TSGetTimeStep(ts, &y));
    453:   if (ctx->semilogy) y = PetscLog10Real(y);
    454:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
    455:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
    456:     PetscCall(PetscDrawLGDraw(ctx->lg));
    457:     PetscCall(PetscDrawLGSave(ctx->lg));
    458:   }
    459:   PetscFunctionReturn(PETSC_SUCCESS);
    460: }

    462: /*@C
    463:   TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`.

    465:   Collective

    467:   Input Parameter:
    468: . ctx - the monitor context

    470:   Level: intermediate

    472:   Note:
    473:   Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`

    475: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
    476: @*/
    477: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
    478: {
    479:   PetscFunctionBegin;
    480:   if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)((void **)&(*ctx)->transformctx));
    481:   PetscCall(PetscDrawLGDestroy(&(*ctx)->lg));
    482:   PetscCall(PetscStrArrayDestroy(&(*ctx)->names));
    483:   PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames));
    484:   PetscCall(PetscFree((*ctx)->displayvariables));
    485:   PetscCall(PetscFree((*ctx)->displayvalues));
    486:   PetscCall(PetscFree(*ctx));
    487:   PetscFunctionReturn(PETSC_SUCCESS);
    488: }

    490: /* Creates a TSMonitorSPCtx for use with DMSwarm particle visualizations */
    491: PetscErrorCode TSMonitorSPCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, PetscInt retain, PetscBool phase, PetscBool multispecies, TSMonitorSPCtx *ctx)
    492: {
    493:   PetscDraw draw;

    495:   PetscFunctionBegin;
    496:   PetscCall(PetscNew(ctx));
    497:   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
    498:   PetscCall(PetscDrawSetFromOptions(draw));
    499:   PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp));
    500:   PetscCall(PetscDrawDestroy(&draw));
    501:   (*ctx)->howoften     = howoften;
    502:   (*ctx)->retain       = retain;
    503:   (*ctx)->phase        = phase;
    504:   (*ctx)->multispecies = multispecies;
    505:   PetscFunctionReturn(PETSC_SUCCESS);
    506: }

    508: /* Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate */
    509: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
    510: {
    511:   PetscFunctionBegin;
    512:   PetscCall(PetscDrawSPDestroy(&(*ctx)->sp));
    513:   PetscCall(PetscFree(*ctx));
    514:   PetscFunctionReturn(PETSC_SUCCESS);
    515: }

    517: /* Creates a TSMonitorHGCtx for use with DMSwarm particle visualizations */
    518: PetscErrorCode TSMonitorHGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, PetscInt Ns, PetscInt Nb, PetscBool velocity, TSMonitorHGCtx *ctx)
    519: {
    520:   PetscDraw draw;
    521:   int       Nsi, Nbi;

    523:   PetscFunctionBegin;
    524:   PetscCall(PetscMPIIntCast(Ns, &Nsi));
    525:   PetscCall(PetscMPIIntCast(Nb, &Nbi));
    526:   PetscCall(PetscNew(ctx));
    527:   PetscCall(PetscMalloc1(Ns, &(*ctx)->hg));
    528:   for (int s = 0; s < Nsi; ++s) {
    529:     PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw));
    530:     PetscCall(PetscDrawSetFromOptions(draw));
    531:     PetscCall(PetscDrawHGCreate(draw, Nbi, &(*ctx)->hg[s]));
    532:     PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE));
    533:     PetscCall(PetscDrawDestroy(&draw));
    534:   }
    535:   (*ctx)->howoften = howoften;
    536:   (*ctx)->Ns       = Ns;
    537:   (*ctx)->velocity = velocity;
    538:   PetscFunctionReturn(PETSC_SUCCESS);
    539: }

    541: /* Destroys a TSMonitorHGCtx that was created with TSMonitorHGCtxCreate */
    542: PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx)
    543: {
    544:   PetscInt s;

    546:   PetscFunctionBegin;
    547:   for (s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s]));
    548:   PetscCall(PetscFree((*ctx)->hg));
    549:   PetscCall(PetscFree(*ctx));
    550:   PetscFunctionReturn(PETSC_SUCCESS);
    551: }

    553: /*@C
    554:   TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling
    555:   `VecView()` for the solution at each timestep

    557:   Collective

    559:   Input Parameters:
    560: + ts    - the `TS` context
    561: . step  - current time-step
    562: . ptime - current time
    563: . u     - the solution at the current time
    564: - dummy - either a viewer or `NULL`

    566:   Options Database Keys:
    567: + -ts_monitor_draw_solution         - draw the solution at each time-step
    568: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution

    570:   Level: intermediate

    572:   Notes:
    573:   The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial`
    574:   will look bad

    576:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with
    577:   `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration.

    579: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
    580: @*/
    581: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
    582: {
    583:   TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
    584:   PetscDraw        draw;

    586:   PetscFunctionBegin;
    587:   if (!step && ictx->showinitial) {
    588:     if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
    589:     PetscCall(VecCopy(u, ictx->initialsolution));
    590:   }
    591:   if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);

    593:   if (ictx->showinitial) {
    594:     PetscReal pause;
    595:     PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
    596:     PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
    597:     PetscCall(VecView(ictx->initialsolution, ictx->viewer));
    598:     PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
    599:     PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
    600:   }
    601:   PetscCall(VecView(u, ictx->viewer));
    602:   if (ictx->showtimestepandtime) {
    603:     PetscReal xl, yl, xr, yr, h;
    604:     char      time[32];

    606:     PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
    607:     PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
    608:     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
    609:     h = yl + .95 * (yr - yl);
    610:     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
    611:     PetscCall(PetscDrawFlush(draw));
    612:   }

    614:   if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
    615:   PetscFunctionReturn(PETSC_SUCCESS);
    616: }

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

    621:   Collective

    623:   Input Parameters:
    624: + ts    - the `TS` context
    625: . step  - current time-step
    626: . ptime - current time
    627: . u     - the solution at the current time
    628: - dummy - either a viewer or `NULL`

    630:   Level: intermediate

    632:   Notes:
    633:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    634:   to be used during the `TS` integration.

    636: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
    637: @*/
    638: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
    639: {
    640:   TSMonitorDrawCtx   ictx = (TSMonitorDrawCtx)dummy;
    641:   PetscDraw          draw;
    642:   PetscDrawAxis      axis;
    643:   PetscInt           n;
    644:   PetscMPIInt        size;
    645:   PetscReal          U0, U1, xl, yl, xr, yr, h;
    646:   char               time[32];
    647:   const PetscScalar *U;

    649:   PetscFunctionBegin;
    650:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
    651:   PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
    652:   PetscCall(VecGetSize(u, &n));
    653:   PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");

    655:   PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
    656:   PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
    657:   PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
    658:   if (!step) {
    659:     PetscCall(PetscDrawClear(draw));
    660:     PetscCall(PetscDrawAxisDraw(axis));
    661:   }

    663:   PetscCall(VecGetArrayRead(u, &U));
    664:   U0 = PetscRealPart(U[0]);
    665:   U1 = PetscRealPart(U[1]);
    666:   PetscCall(VecRestoreArrayRead(u, &U));
    667:   if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);

    669:   PetscDrawCollectiveBegin(draw);
    670:   PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
    671:   if (ictx->showtimestepandtime) {
    672:     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
    673:     PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
    674:     h = yl + .95 * (yr - yl);
    675:     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
    676:   }
    677:   PetscDrawCollectiveEnd(draw);
    678:   PetscCall(PetscDrawFlush(draw));
    679:   PetscCall(PetscDrawPause(draw));
    680:   PetscCall(PetscDrawSave(draw));
    681:   PetscFunctionReturn(PETSC_SUCCESS);
    682: }

    684: /*@C
    685:   TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()`

    687:   Collective

    689:   Input Parameter:
    690: . ictx - the monitor context

    692:   Level: intermediate

    694: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
    695: @*/
    696: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
    697: {
    698:   PetscFunctionBegin;
    699:   PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
    700:   PetscCall(VecDestroy(&(*ictx)->initialsolution));
    701:   PetscCall(PetscFree(*ictx));
    702:   PetscFunctionReturn(PETSC_SUCCESS);
    703: }

    705: /*@C
    706:   TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`

    708:   Collective

    710:   Input Parameters:
    711: + comm     - the MPI communicator to use
    712: . host     - the X display to open, or `NULL` for the local machine
    713: . label    - the title to put in the title bar
    714: . x        - the x screen coordinates of the upper left coordinate of the window
    715: . y        - the y screen coordinates of the upper left coordinate of the window
    716: . m        - the screen width in pixels
    717: . n        - the screen height in pixels
    718: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time

    720:   Output Parameter:
    721: . ctx - the monitor context

    723:   Options Database Keys:
    724: + -ts_monitor_draw_solution         - draw the solution at each time-step
    725: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution

    727:   Level: intermediate

    729:   Note:
    730:   The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`.

    732: .seealso: [](ch_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
    733: @*/
    734: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
    735: {
    736:   PetscFunctionBegin;
    737:   PetscCall(PetscNew(ctx));
    738:   PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
    739:   PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));

    741:   (*ctx)->howoften    = howoften;
    742:   (*ctx)->showinitial = PETSC_FALSE;
    743:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));

    745:   (*ctx)->showtimestepandtime = PETSC_FALSE;
    746:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
    747:   PetscFunctionReturn(PETSC_SUCCESS);
    748: }

    750: /*@C
    751:   TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling
    752:   `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep

    754:   Collective

    756:   Input Parameters:
    757: + ts    - the `TS` context
    758: . step  - current time-step
    759: . ptime - current time
    760: . u     - solution at current time
    761: - dummy - either a viewer or `NULL`

    763:   Options Database Key:
    764: . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`

    766:   Level: intermediate

    768:   Note:
    769:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    770:   to be used during the `TS` integration.

    772: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
    773: @*/
    774: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
    775: {
    776:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
    777:   PetscViewer      viewer = ctx->viewer;
    778:   Vec              work;

    780:   PetscFunctionBegin;
    781:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
    782:   PetscCall(VecDuplicate(u, &work));
    783:   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
    784:   PetscCall(VecView(work, viewer));
    785:   PetscCall(VecDestroy(&work));
    786:   PetscFunctionReturn(PETSC_SUCCESS);
    787: }

    789: /*@C
    790:   TSMonitorDrawError - Monitors progress of the `TS` solvers by calling
    791:   `VecView()` for the error at each timestep

    793:   Collective

    795:   Input Parameters:
    796: + ts    - the `TS` context
    797: . step  - current time-step
    798: . ptime - current time
    799: . u     - solution at current time
    800: - dummy - either a viewer or `NULL`

    802:   Options Database Key:
    803: . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`

    805:   Level: intermediate

    807:   Notes:
    808:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    809:   to be used during the `TS` integration.

    811: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
    812: @*/
    813: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
    814: {
    815:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
    816:   PetscViewer      viewer = ctx->viewer;
    817:   Vec              work;

    819:   PetscFunctionBegin;
    820:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
    821:   PetscCall(VecDuplicate(u, &work));
    822:   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
    823:   PetscCall(VecAXPY(work, -1.0, u));
    824:   PetscCall(VecView(work, viewer));
    825:   PetscCall(VecDestroy(&work));
    826:   PetscFunctionReturn(PETSC_SUCCESS);
    827: }

    829: /*@C
    830:   TSMonitorSolutionSetup - Setups the context for `TSMonitorSolution()`

    832:   Collective

    834:   Input Parameters:
    835: + ts - the `TS` context
    836: - vf - viewer and its format

    838:   Level: intermediate

    840: .seealso: [](ch_ts), `TS`, `TSMonitorSolution()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSetFromOptions()`
    841: @*/
    842: PetscErrorCode TSMonitorSolutionSetup(TS ts, PetscViewerAndFormat *vf)
    843: {
    844:   TSMonitorSolutionCtx ctx;

    846:   PetscFunctionBegin;
    847:   PetscCall(PetscNew(&ctx));
    848:   PetscCall(PetscOptionsGetBool(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_solution_skip_initial", &ctx->skip_initial, NULL));
    849:   vf->data         = ctx;
    850:   vf->data_destroy = PetscCtxDestroyDefault;
    851:   PetscFunctionReturn(PETSC_SUCCESS);
    852: }

    854: /*@C
    855:   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

    857:   Collective

    859:   Input Parameters:
    860: + ts    - the `TS` context
    861: . step  - current time-step
    862: . ptime - current time
    863: . u     - current state
    864: - vf    - viewer and its format

    866:   Level: intermediate

    868:   Notes:
    869:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    870:   to be used during the `TS` integration.

    872: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSolutionSetup()`,
    873: @*/
    874: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
    875: {
    876:   TSMonitorSolutionCtx ctx = (TSMonitorSolutionCtx)vf->data;

    878:   PetscFunctionBegin;
    879:   if (ctx->skip_initial && step == ts->start_step) PetscFunctionReturn(PETSC_SUCCESS);
    880:   if ((vf->view_interval > 0 && !(step % vf->view_interval)) || (vf->view_interval && ts->reason)) {
    881:     PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
    882:     PetscCall(VecView(u, vf->viewer));
    883:     PetscCall(PetscViewerPopFormat(vf->viewer));
    884:   }
    885:   PetscFunctionReturn(PETSC_SUCCESS);
    886: }

    888: /*@C
    889:   TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at selected timesteps.

    891:   Collective

    893:   Input Parameters:
    894: + ts    - the `TS` context
    895: . step  - current time-step
    896: . ptime - current time
    897: . u     - current state
    898: - ctx   - monitor context obtained with `TSMonitorSolutionVTKCtxCreate()`

    900:   Level: developer

    902:   Notes:
    903:   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.
    904:   These are named according to the file name template.

    906:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    907:   to be used during the `TS` integration.

    909: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
    910: @*/
    911: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, TSMonitorVTKCtx ctx)
    912: {
    913:   char        filename[PETSC_MAX_PATH_LEN];
    914:   PetscViewer viewer;

    916:   PetscFunctionBegin;
    917:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
    918:   if (((ctx->interval > 0) && (!(step % ctx->interval))) || (ctx->interval && ts->reason)) {
    919:     PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)ctx->filenametemplate, step));
    920:     PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
    921:     PetscCall(VecView(u, viewer));
    922:     PetscCall(PetscViewerDestroy(&viewer));
    923:   }
    924:   PetscFunctionReturn(PETSC_SUCCESS);
    925: }

    927: /*@C
    928:   TSMonitorSolutionVTKDestroy - Destroy the monitor context created with `TSMonitorSolutionVTKCtxCreate()`

    930:   Not Collective

    932:   Input Parameter:
    933: . ctx - the monitor context

    935:   Level: developer

    937:   Note:
    938:   This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.

    940: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
    941: @*/
    942: PetscErrorCode TSMonitorSolutionVTKDestroy(TSMonitorVTKCtx *ctx)
    943: {
    944:   PetscFunctionBegin;
    945:   PetscCall(PetscFree((*ctx)->filenametemplate));
    946:   PetscCall(PetscFree(*ctx));
    947:   PetscFunctionReturn(PETSC_SUCCESS);
    948: }

    950: /*@C
    951:   TSMonitorSolutionVTKCtxCreate - Create the monitor context to be used in `TSMonitorSolutionVTK()`

    953:   Not collective

    955:   Input Parameter:
    956: . filenametemplate - the template file name, e.g. foo-%03d.vts

    958:   Output Parameter:
    959: . ctx - the monitor context

    961:   Level: developer

    963:   Note:
    964:   This function is normally used inside `TSSetFromOptions()` to pass the context created to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.

    966: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`, `TSMonitorSolutionVTKDestroy()`
    967: @*/
    968: PetscErrorCode TSMonitorSolutionVTKCtxCreate(const char *filenametemplate, TSMonitorVTKCtx *ctx)
    969: {
    970:   const char     *ptr = NULL, *ptr2 = NULL;
    971:   TSMonitorVTKCtx ictx;

    973:   PetscFunctionBegin;
    974:   PetscAssertPointer(filenametemplate, 1);
    975:   PetscAssertPointer(ctx, 2);
    976:   /* Do some cursory validation of the input. */
    977:   PetscCall(PetscStrstr(filenametemplate, "%", (char **)&ptr));
    978:   PetscCheck(ptr, PETSC_COMM_SELF, PETSC_ERR_USER, "-ts_monitor_solution_vtk requires a file template, e.g. filename-%%03" PetscInt_FMT ".vts");
    979:   for (ptr++; ptr && *ptr; ptr++) {
    980:     PetscCall(PetscStrchr("DdiouxX", *ptr, (char **)&ptr2));
    981:     PetscCheck(ptr2 || (*ptr >= '0' && *ptr <= '9'), PETSC_COMM_SELF, PETSC_ERR_USER, "Invalid file template argument to -ts_monitor_solution_vtk, should look like filename-%%03" PetscInt_FMT ".vts");
    982:     if (ptr2) break;
    983:   }
    984:   PetscCall(PetscNew(&ictx));
    985:   PetscCall(PetscStrallocpy(filenametemplate, &ictx->filenametemplate));
    986:   ictx->interval = 1;

    988:   *ctx = ictx;
    989:   PetscFunctionReturn(PETSC_SUCCESS);
    990: }

    992: /*@C
    993:   TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector
    994:   in a time based line graph

    996:   Collective

    998:   Input Parameters:
    999: + ts    - the `TS` context
    1000: . step  - current time-step
    1001: . ptime - current time
    1002: . u     - current solution
    1003: - dctx  - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`

    1005:   Options Database Key:
    1006: . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables

    1008:   Level: intermediate

    1010:   Notes:
    1011:   Each process in a parallel run displays its component solutions in a separate window

    1013:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    1014:   to be used during the `TS` integration.

    1016: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
    1017:           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
    1018:           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
    1019:           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
    1020: @*/
    1021: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
    1022: {
    1023:   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)dctx;
    1024:   const PetscScalar *yy;
    1025:   Vec                v;

    1027:   PetscFunctionBegin;
    1028:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
    1029:   if (!step) {
    1030:     PetscDrawAxis axis;
    1031:     PetscInt      dim;
    1032:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
    1033:     PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
    1034:     if (!ctx->names) {
    1035:       PetscBool flg;
    1036:       /* user provides names of variables to plot but no names has been set so assume names are integer values */
    1037:       PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
    1038:       if (flg) {
    1039:         PetscInt i, n;
    1040:         char   **names;
    1041:         PetscCall(VecGetSize(u, &n));
    1042:         PetscCall(PetscMalloc1(n + 1, &names));
    1043:         for (i = 0; i < n; i++) {
    1044:           PetscCall(PetscMalloc1(5, &names[i]));
    1045:           PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
    1046:         }
    1047:         names[n]   = NULL;
    1048:         ctx->names = names;
    1049:       }
    1050:     }
    1051:     if (ctx->names && !ctx->displaynames) {
    1052:       char    **displaynames;
    1053:       PetscBool flg;
    1054:       PetscCall(VecGetLocalSize(u, &dim));
    1055:       PetscCall(PetscCalloc1(dim + 1, &displaynames));
    1056:       PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
    1057:       if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
    1058:       PetscCall(PetscStrArrayDestroy(&displaynames));
    1059:     }
    1060:     if (ctx->displaynames) {
    1061:       PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
    1062:       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
    1063:     } else if (ctx->names) {
    1064:       PetscCall(VecGetLocalSize(u, &dim));
    1065:       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
    1066:       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
    1067:     } else {
    1068:       PetscCall(VecGetLocalSize(u, &dim));
    1069:       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
    1070:     }
    1071:     PetscCall(PetscDrawLGReset(ctx->lg));
    1072:   }

    1074:   if (!ctx->transform) v = u;
    1075:   else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
    1076:   PetscCall(VecGetArrayRead(v, &yy));
    1077:   if (ctx->displaynames) {
    1078:     PetscInt i;
    1079:     for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
    1080:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
    1081:   } else {
    1082: #if defined(PETSC_USE_COMPLEX)
    1083:     PetscInt   i, n;
    1084:     PetscReal *yreal;
    1085:     PetscCall(VecGetLocalSize(v, &n));
    1086:     PetscCall(PetscMalloc1(n, &yreal));
    1087:     for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
    1088:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
    1089:     PetscCall(PetscFree(yreal));
    1090: #else
    1091:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
    1092: #endif
    1093:   }
    1094:   PetscCall(VecRestoreArrayRead(v, &yy));
    1095:   if (ctx->transform) PetscCall(VecDestroy(&v));

    1097:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
    1098:     PetscCall(PetscDrawLGDraw(ctx->lg));
    1099:     PetscCall(PetscDrawLGSave(ctx->lg));
    1100:   }
    1101:   PetscFunctionReturn(PETSC_SUCCESS);
    1102: }

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

    1107:   Collective

    1109:   Input Parameters:
    1110: + ts    - the `TS` context
    1111: - names - the names of the components, final string must be `NULL`

    1113:   Level: intermediate

    1115:   Notes:
    1116:   If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

    1118: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
    1119: @*/
    1120: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
    1121: {
    1122:   PetscInt i;

    1124:   PetscFunctionBegin;
    1125:   for (i = 0; i < ts->numbermonitors; i++) {
    1126:     if (ts->monitor[i] == TSMonitorLGSolution) {
    1127:       PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
    1128:       break;
    1129:     }
    1130:   }
    1131:   PetscFunctionReturn(PETSC_SUCCESS);
    1132: }

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

    1137:   Collective

    1139:   Input Parameters:
    1140: + ctx   - the `TS` context
    1141: - names - the names of the components, final string must be `NULL`

    1143:   Level: intermediate

    1145: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
    1146: @*/
    1147: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
    1148: {
    1149:   PetscFunctionBegin;
    1150:   PetscCall(PetscStrArrayDestroy(&ctx->names));
    1151:   PetscCall(PetscStrArrayallocpy(names, &ctx->names));
    1152:   PetscFunctionReturn(PETSC_SUCCESS);
    1153: }

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

    1158:   Collective

    1160:   Input Parameter:
    1161: . ts - the `TS` context

    1163:   Output Parameter:
    1164: . names - the names of the components, final string must be `NULL`

    1166:   Level: intermediate

    1168:   Note:
    1169:   If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

    1171: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
    1172: @*/
    1173: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
    1174: {
    1175:   PetscInt i;

    1177:   PetscFunctionBegin;
    1178:   *names = NULL;
    1179:   for (i = 0; i < ts->numbermonitors; i++) {
    1180:     if (ts->monitor[i] == TSMonitorLGSolution) {
    1181:       TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
    1182:       *names             = (const char *const *)ctx->names;
    1183:       break;
    1184:     }
    1185:   }
    1186:   PetscFunctionReturn(PETSC_SUCCESS);
    1187: }

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

    1192:   Collective

    1194:   Input Parameters:
    1195: + ctx          - the `TSMonitorLG` context
    1196: - displaynames - the names of the components, final string must be `NULL`

    1198:   Level: intermediate

    1200: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
    1201: @*/
    1202: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
    1203: {
    1204:   PetscInt j = 0, k;

    1206:   PetscFunctionBegin;
    1207:   if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
    1208:   PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
    1209:   PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
    1210:   while (displaynames[j]) j++;
    1211:   ctx->ndisplayvariables = j;
    1212:   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
    1213:   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
    1214:   j = 0;
    1215:   while (displaynames[j]) {
    1216:     k = 0;
    1217:     while (ctx->names[k]) {
    1218:       PetscBool flg;
    1219:       PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
    1220:       if (flg) {
    1221:         ctx->displayvariables[j] = k;
    1222:         break;
    1223:       }
    1224:       k++;
    1225:     }
    1226:     j++;
    1227:   }
    1228:   PetscFunctionReturn(PETSC_SUCCESS);
    1229: }

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

    1234:   Collective

    1236:   Input Parameters:
    1237: + ts           - the `TS` context
    1238: - displaynames - the names of the components, final string must be `NULL`

    1240:   Level: intermediate

    1242:   Note:
    1243:   If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

    1245: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
    1246: @*/
    1247: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
    1248: {
    1249:   PetscInt i;

    1251:   PetscFunctionBegin;
    1252:   for (i = 0; i < ts->numbermonitors; i++) {
    1253:     if (ts->monitor[i] == TSMonitorLGSolution) {
    1254:       PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
    1255:       break;
    1256:     }
    1257:   }
    1258:   PetscFunctionReturn(PETSC_SUCCESS);
    1259: }

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

    1264:   Collective

    1266:   Input Parameters:
    1267: + ts        - the `TS` context
    1268: . transform - the transform function
    1269: . destroy   - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
    1270: - tctx      - optional context used by transform function

    1272:   Level: intermediate

    1274:   Note:
    1275:   If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

    1277: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()`, `PetscCtxDestroyFn`
    1278: @*/
    1279: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscCtxDestroyFn *destroy, void *tctx)
    1280: {
    1281:   PetscInt i;

    1283:   PetscFunctionBegin;
    1284:   for (i = 0; i < ts->numbermonitors; i++) {
    1285:     if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
    1286:   }
    1287:   PetscFunctionReturn(PETSC_SUCCESS);
    1288: }

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

    1293:   Collective

    1295:   Input Parameters:
    1296: + tctx      - the `TS` context
    1297: . transform - the transform function
    1298: . destroy   - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
    1299: - ctx       - optional context used by transform function

    1301:   Level: intermediate

    1303: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`, `PetscCtxDestroyFn`
    1304: @*/
    1305: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscCtxDestroyFn *destroy, void *tctx)
    1306: {
    1307:   PetscFunctionBegin;
    1308:   ctx->transform        = transform;
    1309:   ctx->transformdestroy = destroy;
    1310:   ctx->transformctx     = tctx;
    1311:   PetscFunctionReturn(PETSC_SUCCESS);
    1312: }

    1314: /*@C
    1315:   TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error
    1316:   in a time based line graph

    1318:   Collective

    1320:   Input Parameters:
    1321: + ts    - the `TS` context
    1322: . step  - current time-step
    1323: . ptime - current time
    1324: . u     - current solution
    1325: - dummy - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`

    1327:   Options Database Key:
    1328: . -ts_monitor_lg_error - create a graphical monitor of error history

    1330:   Level: intermediate

    1332:   Notes:
    1333:   Each process in a parallel run displays its component errors in a separate window

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

    1337:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    1338:   to be used during the TS integration.

    1340: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
    1341: @*/
    1342: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
    1343: {
    1344:   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)dummy;
    1345:   const PetscScalar *yy;
    1346:   Vec                y;

    1348:   PetscFunctionBegin;
    1349:   if (!step) {
    1350:     PetscDrawAxis axis;
    1351:     PetscInt      dim;
    1352:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
    1353:     PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
    1354:     PetscCall(VecGetLocalSize(u, &dim));
    1355:     PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
    1356:     PetscCall(PetscDrawLGReset(ctx->lg));
    1357:   }
    1358:   PetscCall(VecDuplicate(u, &y));
    1359:   PetscCall(TSComputeSolutionFunction(ts, ptime, y));
    1360:   PetscCall(VecAXPY(y, -1.0, u));
    1361:   PetscCall(VecGetArrayRead(y, &yy));
    1362: #if defined(PETSC_USE_COMPLEX)
    1363:   {
    1364:     PetscReal *yreal;
    1365:     PetscInt   i, n;
    1366:     PetscCall(VecGetLocalSize(y, &n));
    1367:     PetscCall(PetscMalloc1(n, &yreal));
    1368:     for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
    1369:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
    1370:     PetscCall(PetscFree(yreal));
    1371:   }
    1372: #else
    1373:   PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
    1374: #endif
    1375:   PetscCall(VecRestoreArrayRead(y, &yy));
    1376:   PetscCall(VecDestroy(&y));
    1377:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
    1378:     PetscCall(PetscDrawLGDraw(ctx->lg));
    1379:     PetscCall(PetscDrawLGSave(ctx->lg));
    1380:   }
    1381:   PetscFunctionReturn(PETSC_SUCCESS);
    1382: }

    1384: /*@C
    1385:   TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot

    1387:   Input Parameters:
    1388: + ts    - the `TS` context
    1389: . step  - current time-step
    1390: . ptime - current time
    1391: . u     - current solution
    1392: - dctx  - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()`

    1394:   Options Database Keys:
    1395: + -ts_monitor_sp_swarm <n>                  - Monitor the solution every n steps, or -1 for plotting only the final solution
    1396: . -ts_monitor_sp_swarm_retain <n>           - Retain n old points so we can see the history, or -1 for all points
    1397: . -ts_monitor_sp_swarm_multi_species <bool> - Color each species differently
    1398: - -ts_monitor_sp_swarm_phase <bool>         - Plot in phase space, as opposed to coordinate space

    1400:   Level: intermediate

    1402:   Notes:
    1403:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    1404:   to be used during the `TS` integration.

    1406: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
    1407: @*/
    1408: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
    1409: {
    1410:   TSMonitorSPCtx     ctx = (TSMonitorSPCtx)dctx;
    1411:   PetscDraw          draw;
    1412:   DM                 dm, cdm;
    1413:   const PetscScalar *yy;
    1414:   PetscInt           Np, p, dim = 2, *species;
    1415:   PetscReal          species_color;

    1417:   PetscFunctionBegin;
    1418:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
    1419:   PetscCall(TSGetDM(ts, &dm));
    1420:   if (!step) {
    1421:     PetscDrawAxis axis;
    1422:     PetscReal     dmboxlower[2], dmboxupper[2];

    1424:     PetscCall(TSGetDM(ts, &dm));
    1425:     PetscCall(DMGetDimension(dm, &dim));
    1426:     PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
    1427:     PetscCall(DMSwarmGetCellDM(dm, &cdm));
    1428:     PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
    1429:     PetscCall(VecGetLocalSize(u, &Np));
    1430:     Np /= dim * 2;
    1431:     PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
    1432:     if (ctx->phase) {
    1433:       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
    1434:       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10));
    1435:     } else {
    1436:       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
    1437:       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
    1438:     }
    1439:     PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
    1440:     PetscCall(PetscDrawSPReset(ctx->sp));
    1441:   }
    1442:   if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species));
    1443:   PetscCall(VecGetLocalSize(u, &Np));
    1444:   Np /= dim * 2;
    1445:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
    1446:     PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
    1447:     if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
    1448:     PetscCall(PetscDrawFlush(draw));
    1449:     PetscCall(PetscDrawSPReset(ctx->sp));
    1450:     PetscCall(VecGetArrayRead(u, &yy));
    1451:     for (p = 0; p < Np; ++p) {
    1452:       PetscReal x, y;

    1454:       if (ctx->phase) {
    1455:         x = PetscRealPart(yy[p * dim * 2]);
    1456:         y = PetscRealPart(yy[p * dim * 2 + dim]);
    1457:       } else {
    1458:         x = PetscRealPart(yy[p * dim * 2]);
    1459:         y = PetscRealPart(yy[p * dim * 2 + 1]);
    1460:       }
    1461:       if (ctx->multispecies) {
    1462:         species_color = species[p] + 2;
    1463:         PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color));
    1464:       } else {
    1465:         PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
    1466:       }
    1467:       PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
    1468:     }
    1469:     PetscCall(VecRestoreArrayRead(u, &yy));
    1470:     PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
    1471:     PetscCall(PetscDrawSPSave(ctx->sp));
    1472:     if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species));
    1473:   }
    1474:   PetscFunctionReturn(PETSC_SUCCESS);
    1475: }

    1477: /*@C
    1478:   TSMonitorHGSwarmSolution - Graphically displays histograms of `DMSWARM` particles

    1480:   Input Parameters:
    1481: + ts    - the `TS` context
    1482: . step  - current time-step
    1483: . ptime - current time
    1484: . u     - current solution
    1485: - dctx  - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorHGCtxCreate()`

    1487:   Options Database Keys:
    1488: + -ts_monitor_hg_swarm <n>             - Monitor the solution every n steps, or -1 for plotting only the final solution
    1489: . -ts_monitor_hg_swarm_species <num>   - Number of species to histogram
    1490: . -ts_monitor_hg_swarm_bins <num>      - Number of histogram bins
    1491: - -ts_monitor_hg_swarm_velocity <bool> - Plot in velocity space, as opposed to coordinate space

    1493:   Level: intermediate

    1495:   Note:
    1496:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    1497:   to be used during the `TS` integration.

    1499: .seealso: `TSMonitorSet()`
    1500: @*/
    1501: PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
    1502: {
    1503:   TSMonitorHGCtx     ctx = (TSMonitorHGCtx)dctx;
    1504:   PetscDraw          draw;
    1505:   DM                 sw;
    1506:   const PetscScalar *yy;
    1507:   PetscInt          *species;
    1508:   PetscInt           dim, d = 0, Np, p, Ns, s;

    1510:   PetscFunctionBegin;
    1511:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
    1512:   PetscCall(TSGetDM(ts, &sw));
    1513:   PetscCall(DMGetDimension(sw, &dim));
    1514:   PetscCall(DMSwarmGetNumSpecies(sw, &Ns));
    1515:   Ns = PetscMin(Ns, ctx->Ns);
    1516:   PetscCall(VecGetLocalSize(u, &Np));
    1517:   Np /= dim * 2;
    1518:   if (!step) {
    1519:     PetscDrawAxis axis;
    1520:     char          title[PETSC_MAX_PATH_LEN];

    1522:     for (s = 0; s < Ns; ++s) {
    1523:       PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis));
    1524:       PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s));
    1525:       if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N"));
    1526:       else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N"));
    1527:     }
    1528:   }
    1529:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
    1530:     PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species));
    1531:     for (s = 0; s < Ns; ++s) {
    1532:       PetscCall(PetscDrawHGReset(ctx->hg[s]));
    1533:       PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw));
    1534:       PetscCall(PetscDrawClear(draw));
    1535:       PetscCall(PetscDrawFlush(draw));
    1536:     }
    1537:     PetscCall(VecGetArrayRead(u, &yy));
    1538:     for (p = 0; p < Np; ++p) {
    1539:       const PetscInt s = species[p] < Ns ? species[p] : 0;
    1540:       PetscReal      v;

    1542:       if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]);
    1543:       else v = PetscRealPart(yy[p * dim * 2 + d]);
    1544:       PetscCall(PetscDrawHGAddValue(ctx->hg[s], v));
    1545:     }
    1546:     PetscCall(VecRestoreArrayRead(u, &yy));
    1547:     for (s = 0; s < Ns; ++s) {
    1548:       PetscCall(PetscDrawHGDraw(ctx->hg[s]));
    1549:       PetscCall(PetscDrawHGSave(ctx->hg[s]));
    1550:     }
    1551:     PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species));
    1552:   }
    1553:   PetscFunctionReturn(PETSC_SUCCESS);
    1554: }

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

    1559:   Collective

    1561:   Input Parameters:
    1562: + ts    - the `TS` context
    1563: . step  - current time-step
    1564: . ptime - current time
    1565: . u     - current solution
    1566: - vf    - unused context

    1568:   Options Database Key:
    1569: . -ts_monitor_error - create a graphical monitor of error history

    1571:   Level: intermediate

    1573:   Notes:
    1574:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    1575:   to be used during the `TS` integration.

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

    1579: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
    1580: @*/
    1581: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
    1582: {
    1583:   DM        dm;
    1584:   PetscDS   ds = NULL;
    1585:   PetscInt  Nf = -1, f;
    1586:   PetscBool flg;

    1588:   PetscFunctionBegin;
    1589:   PetscCall(TSGetDM(ts, &dm));
    1590:   if (dm) PetscCall(DMGetDS(dm, &ds));
    1591:   if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
    1592:   if (Nf <= 0) {
    1593:     Vec       y;
    1594:     PetscReal nrm;

    1596:     PetscCall(VecDuplicate(u, &y));
    1597:     PetscCall(TSComputeSolutionFunction(ts, ptime, y));
    1598:     PetscCall(VecAXPY(y, -1.0, u));
    1599:     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
    1600:     if (flg) {
    1601:       PetscCall(VecNorm(y, NORM_2, &nrm));
    1602:       PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
    1603:     }
    1604:     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
    1605:     if (flg) PetscCall(VecView(y, vf->viewer));
    1606:     PetscCall(VecDestroy(&y));
    1607:   } else {
    1608:     PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
    1609:     void    **ctxs;
    1610:     Vec       v;
    1611:     PetscReal ferrors[1];

    1613:     PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
    1614:     for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
    1615:     PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
    1616:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04" PetscInt_FMT " time = %-8.4g \t L_2 Error: [", step, (double)ptime));
    1617:     for (f = 0; f < Nf; ++f) {
    1618:       if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
    1619:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
    1620:     }
    1621:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));

    1623:     PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));

    1625:     PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
    1626:     if (flg) {
    1627:       PetscCall(DMGetGlobalVector(dm, &v));
    1628:       PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
    1629:       PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
    1630:       PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
    1631:       PetscCall(DMRestoreGlobalVector(dm, &v));
    1632:     }
    1633:     PetscCall(PetscFree2(exactFuncs, ctxs));
    1634:   }
    1635:   PetscFunctionReturn(PETSC_SUCCESS);
    1636: }

    1638: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
    1639: {
    1640:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
    1641:   PetscReal      x   = ptime, y;
    1642:   PetscInt       its;

    1644:   PetscFunctionBegin;
    1645:   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
    1646:   if (!n) {
    1647:     PetscDrawAxis axis;
    1648:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
    1649:     PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
    1650:     PetscCall(PetscDrawLGReset(ctx->lg));
    1651:     ctx->snes_its = 0;
    1652:   }
    1653:   PetscCall(TSGetSNESIterations(ts, &its));
    1654:   y = its - ctx->snes_its;
    1655:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
    1656:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
    1657:     PetscCall(PetscDrawLGDraw(ctx->lg));
    1658:     PetscCall(PetscDrawLGSave(ctx->lg));
    1659:   }
    1660:   ctx->snes_its = its;
    1661:   PetscFunctionReturn(PETSC_SUCCESS);
    1662: }

    1664: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
    1665: {
    1666:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
    1667:   PetscReal      x   = ptime, y;
    1668:   PetscInt       its;

    1670:   PetscFunctionBegin;
    1671:   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
    1672:   if (!n) {
    1673:     PetscDrawAxis axis;
    1674:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
    1675:     PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
    1676:     PetscCall(PetscDrawLGReset(ctx->lg));
    1677:     ctx->ksp_its = 0;
    1678:   }
    1679:   PetscCall(TSGetKSPIterations(ts, &its));
    1680:   y = its - ctx->ksp_its;
    1681:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
    1682:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
    1683:     PetscCall(PetscDrawLGDraw(ctx->lg));
    1684:     PetscCall(PetscDrawLGSave(ctx->lg));
    1685:   }
    1686:   ctx->ksp_its = its;
    1687:   PetscFunctionReturn(PETSC_SUCCESS);
    1688: }

    1690: /*@C
    1691:   TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()`

    1693:   Collective

    1695:   Input Parameter:
    1696: . ts - the `TS` solver object

    1698:   Output Parameter:
    1699: . ctx - the context

    1701:   Level: intermediate

    1703: .seealso: [](ch_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
    1704: @*/
    1705: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
    1706: {
    1707:   PetscFunctionBegin;
    1708:   PetscCall(PetscNew(ctx));
    1709:   PetscFunctionReturn(PETSC_SUCCESS);
    1710: }

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

    1715:   Collective

    1717:   Input Parameters:
    1718: + ts    - the `TS` context
    1719: . step  - current time-step
    1720: . ptime - current time
    1721: . u     - current solution
    1722: - dctx  - the envelope context

    1724:   Options Database Key:
    1725: . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time

    1727:   Level: intermediate

    1729:   Notes:
    1730:   After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope

    1732:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    1733:   to be used during the `TS` integration.

    1735: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
    1736: @*/
    1737: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
    1738: {
    1739:   TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;

    1741:   PetscFunctionBegin;
    1742:   if (!ctx->max) {
    1743:     PetscCall(VecDuplicate(u, &ctx->max));
    1744:     PetscCall(VecDuplicate(u, &ctx->min));
    1745:     PetscCall(VecCopy(u, ctx->max));
    1746:     PetscCall(VecCopy(u, ctx->min));
    1747:   } else {
    1748:     PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
    1749:     PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
    1750:   }
    1751:   PetscFunctionReturn(PETSC_SUCCESS);
    1752: }

    1754: /*@C
    1755:   TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution

    1757:   Collective

    1759:   Input Parameter:
    1760: . ts - the `TS` context

    1762:   Output Parameters:
    1763: + max - the maximum values
    1764: - min - the minimum values

    1766:   Level: intermediate

    1768:   Notes:
    1769:   If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored

    1771: .seealso: [](ch_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
    1772: @*/
    1773: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
    1774: {
    1775:   PetscInt i;

    1777:   PetscFunctionBegin;
    1778:   if (max) *max = NULL;
    1779:   if (min) *min = NULL;
    1780:   for (i = 0; i < ts->numbermonitors; i++) {
    1781:     if (ts->monitor[i] == TSMonitorEnvelope) {
    1782:       TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
    1783:       if (max) *max = ctx->max;
    1784:       if (min) *min = ctx->min;
    1785:       break;
    1786:     }
    1787:   }
    1788:   PetscFunctionReturn(PETSC_SUCCESS);
    1789: }

    1791: /*@C
    1792:   TSMonitorEnvelopeCtxDestroy - Destroys a context that was created  with `TSMonitorEnvelopeCtxCreate()`.

    1794:   Collective

    1796:   Input Parameter:
    1797: . ctx - the monitor context

    1799:   Level: intermediate

    1801: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
    1802: @*/
    1803: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
    1804: {
    1805:   PetscFunctionBegin;
    1806:   PetscCall(VecDestroy(&(*ctx)->min));
    1807:   PetscCall(VecDestroy(&(*ctx)->max));
    1808:   PetscCall(PetscFree(*ctx));
    1809:   PetscFunctionReturn(PETSC_SUCCESS);
    1810: }

    1812: /*@C
    1813:   TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS`

    1815:   Not Collective

    1817:   Input Parameters:
    1818: + ts   - the `TS` context
    1819: . step - current timestep
    1820: . t    - current time
    1821: . U    - current solution
    1822: - vf   - not used

    1824:   Options Database Key:
    1825: + -ts_dmswarm_monitor_moments          - Monitor moments of particle distribution
    1826: - -ts_dmswarm_monitor_moments_interval - Interval of timesteps between monitor outputs

    1828:   Level: intermediate

    1830:   Notes:
    1831:   This requires a `DMSWARM` be attached to the `TS`.

    1833:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
    1834:   to be used during the TS integration.

    1836: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
    1837: @*/
    1838: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
    1839: {
    1840:   DM                 sw;
    1841:   const PetscScalar *u;
    1842:   PetscReal          m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
    1843:   PetscInt           dim, d, Np, p;
    1844:   MPI_Comm           comm;

    1846:   PetscFunctionBeginUser;
    1847:   (void)t;
    1848:   (void)vf;
    1849:   PetscCall(TSGetDM(ts, &sw));
    1850:   if (!sw || step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS);
    1851:   PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
    1852:   PetscCall(DMGetDimension(sw, &dim));
    1853:   PetscCall(VecGetLocalSize(U, &Np));
    1854:   Np /= dim;
    1855:   PetscCall(VecGetArrayRead(U, &u));
    1856:   for (p = 0; p < Np; ++p) {
    1857:     for (d = 0; d < dim; ++d) {
    1858:       totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
    1859:       totMom[d] += PetscRealPart(u[p * dim + d]);
    1860:     }
    1861:   }
    1862:   PetscCall(VecRestoreArrayRead(U, &u));
    1863:   for (d = 0; d < dim; ++d) totMom[d] *= m;
    1864:   totE *= 0.5 * m;
    1865:   PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
    1866:   for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, "    Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
    1867:   PetscCall(PetscPrintf(comm, "\n"));
    1868:   PetscFunctionReturn(PETSC_SUCCESS);
    1869: }