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:   Calling sequence of `monitor`:
 60: + ts   - the `TS` to monitor
 61: . step - the current time-step
 62: . time - the current time
 63: . u    - the current solution
 64: - vf   - the `PetscViewer` and format to monitor with

 66:   Calling sequence of `monitorsetup`:
 67: + ts - the `TS` to monitor
 68: - vf - the `PetscViewer` and format to monitor with

 70:   Level: developer

 72: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `PetscOptionsCreateViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
 73:           `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`,
 74:           `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`,
 75:           `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
 76:           `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
 77:           `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
 78:           `PetscOptionsFList()`, `PetscOptionsEList()`
 79: @*/
 80: PetscErrorCode TSMonitorSetFromOptions(TS ts, const char name[], const char help[], const char manual[], PetscErrorCode (*monitor)(TS ts, PetscInt step, PetscReal time, Vec u, PetscViewerAndFormat *vf), PetscErrorCode (*monitorsetup)(TS ts, PetscViewerAndFormat *vf))
 81: {
 82:   PetscViewer       viewer;
 83:   PetscViewerFormat format;
 84:   PetscBool         flg;

 86:   PetscFunctionBegin;
 87:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)ts), ((PetscObject)ts)->options, ((PetscObject)ts)->prefix, name, &viewer, &format, &flg));
 88:   if (flg) {
 89:     PetscViewerAndFormat *vf;
 90:     char                  interval_key[1024];

 92:     PetscCall(PetscSNPrintf(interval_key, sizeof interval_key, "%s_interval", name));
 93:     PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
 94:     vf->view_interval = 1;
 95:     PetscCall(PetscOptionsGetInt(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, interval_key, &vf->view_interval, NULL));

 97:     PetscCall(PetscViewerDestroy(&viewer));
 98:     if (monitorsetup) PetscCall((*monitorsetup)(ts, vf));
 99:     PetscCall(TSMonitorSet(ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, PetscCtx))monitor, vf, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy));
100:   }
101:   PetscFunctionReturn(PETSC_SUCCESS);
102: }

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

108:   Logically Collective

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

116:   Calling sequence of `monitor`:
117: + ts    - the `TS` context
118: . 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)
119: . time  - current time
120: . u     - current iterate
121: - ctx   - [optional] monitoring context

123:   Level: intermediate

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

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

131: .seealso: [](ch_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
132:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
133:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `PetscCtxDestroyFn`
134: @*/
135: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS ts, PetscInt steps, PetscReal time, Vec u, PetscCtx ctx), PetscCtx mctx, PetscCtxDestroyFn *mdestroy)
136: {
137:   PetscFunctionBegin;
139:   for (PetscInt i = 0; i < ts->numbermonitors; i++) {
140:     PetscBool identical;

142:     PetscCall(PetscMonitorCompare((PetscErrorCode (*)(void))(PetscVoidFn *)monitor, mctx, mdestroy, (PetscErrorCode (*)(void))(PetscVoidFn *)ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical));
143:     if (identical) PetscFunctionReturn(PETSC_SUCCESS);
144:   }
145:   PetscCheck(ts->numbermonitors < MAXTSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set");
146:   ts->monitor[ts->numbermonitors]          = monitor;
147:   ts->monitordestroy[ts->numbermonitors]   = mdestroy;
148:   ts->monitorcontext[ts->numbermonitors++] = mctx;
149:   PetscFunctionReturn(PETSC_SUCCESS);
150: }

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

155:   Logically Collective

157:   Input Parameter:
158: . ts - the `TS` context obtained from `TSCreate()`

160:   Level: intermediate

162:   Note:
163:   There is no way to remove a single, specific monitor.

165: .seealso: [](ch_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
166: @*/
167: PetscErrorCode TSMonitorCancel(TS ts)
168: {
169:   PetscFunctionBegin;
171:   for (PetscInt i = 0; i < ts->numbermonitors; i++) {
172:     if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i]));
173:   }
174:   ts->numbermonitors = 0;
175:   PetscFunctionReturn(PETSC_SUCCESS);
176: }

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

181:   Input Parameters:
182: + ts    - the `TS` context
183: . 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)
184: . ptime - current time
185: . v     - current iterate
186: - vf    - the viewer and format

188:   Options Database Key:
189: . -ts_monitor - monitors the time integration

191:   Level: intermediate

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

197: .seealso: [](ch_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorWallClockTime()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
198:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
199:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`
200: @*/
201: PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
202: {
203:   PetscViewer viewer = vf->viewer;
204:   PetscBool   isascii, ibinary;

206:   PetscFunctionBegin;
208:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
209:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary));
210:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
211:   if (isascii) {
212:     const char *prefix;

214:     PetscCall(PetscObjectGetOptionsPrefix((PetscObject)ts, &prefix));
215:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
216:     if (step == -1) { /* this indicates it is an interpolated solution */
217:       PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps));
218:     } else {
219:       PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS%s%s%s dt %g time %g%s", step, prefix ? " (" : "", prefix ? prefix : "", prefix ? ")" : "", (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n"));
220:     }
221:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
222:   } else if (ibinary) {
223:     PetscMPIInt rank;
224:     PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
225:     if (rank == 0) {
226:       PetscBool skipHeader;
227:       PetscInt  classid = REAL_FILE_CLASSID;

229:       PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
230:       if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
231:       PetscCall(PetscRealView(1, &ptime, viewer));
232:     } else {
233:       PetscCall(PetscRealView(0, &ptime, viewer));
234:     }
235:   }
236:   PetscCall(PetscViewerPopFormat(viewer));
237:   PetscFunctionReturn(PETSC_SUCCESS);
238: }

240: typedef struct {
241:   PetscLogDouble time_start;
242:   PetscLogDouble time_last;
243:   PetscInt       snes_its;
244:   PetscInt       ksp_its;
245: } *TSMonitorWallClockTimeContext;

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

250:   Input Parameters:
251: + ts - the `TS` context
252: - vf - the viewer and format

254:   Level: intermediate

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

260: .seealso: [](ch_ts), `TSMonitorSet()`
261: @*/
262: PetscErrorCode TSMonitorWallClockTimeSetUp(TS ts, PetscViewerAndFormat *vf)
263: {
264:   TSMonitorWallClockTimeContext speed;

266:   PetscFunctionBegin;
267:   PetscCall(PetscNew(&speed));
268:   speed->time_start = PETSC_DECIDE;
269:   vf->data_destroy  = PetscCtxDestroyDefault;
270:   vf->data          = speed;
271:   PetscFunctionReturn(PETSC_SUCCESS);
272: }

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

277:   Input Parameters:
278: + ts    - the `TS` context
279: . 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)
280: . ptime - current time
281: . v     - current solution
282: - vf    - the viewer and format

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

287:   Level: intermediate

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

293: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
294:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
295:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
296: @*/
297: PetscErrorCode TSMonitorWallClockTime(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
298: {
299:   PetscViewer                   viewer = vf->viewer;
300:   TSMonitorWallClockTimeContext speed  = (TSMonitorWallClockTimeContext)vf->data;
301:   PetscBool                     isascii;
302:   PetscLogDouble                now;
303:   PetscInt                      snes_its, ksp_its;

305:   PetscFunctionBegin;
307:   PetscCall(PetscTime(&now));
308:   if (speed->time_start == PETSC_DECIDE) {
309:     speed->time_start = now;
310:     speed->time_last  = now;
311:   }
312:   PetscCall(TSGetSNESIterations(ts, &snes_its));
313:   PetscCall(TSGetKSPIterations(ts, &ksp_its));
314:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
315:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
316:   if (isascii) {
317:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
318:     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,
319:                                      now - speed->time_start, snes_its - speed->snes_its, ksp_its - speed->ksp_its));
320:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
321:   }
322:   PetscCall(PetscViewerPopFormat(viewer));
323:   speed->time_last = now;
324:   speed->snes_its  = snes_its;
325:   speed->ksp_its   = ksp_its;
326:   PetscFunctionReturn(PETSC_SUCCESS);
327: }

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

332:   Input Parameters:
333: + ts    - the `TS` context
334: . 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)
335: . ptime - current time
336: . v     - current iterate
337: - vf    - the viewer and format

339:   Level: intermediate

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

345: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`
346: @*/
347: PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
348: {
349:   PetscViewer viewer = vf->viewer;
350:   PetscBool   isascii;
351:   PetscReal   max, min;

353:   PetscFunctionBegin;
355:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
356:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
357:   if (isascii) {
358:     PetscCall(VecMax(v, NULL, &max));
359:     PetscCall(VecMin(v, NULL, &min));
360:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
361:     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));
362:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
363:   }
364:   PetscCall(PetscViewerPopFormat(viewer));
365:   PetscFunctionReturn(PETSC_SUCCESS);
366: }

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

372:   Collective

374:   Input Parameters:
375: + comm     - the MPI communicator to use
376: . host     - the X display to open, or `NULL` for the local machine
377: . label    - the title to put in the title bar
378: . x        - the x screen coordinates of the upper left coordinate of the window
379: . y        - the y screen coordinates of the upper left coordinate of the window
380: . m        - the screen width in pixels
381: . n        - the screen height in pixels
382: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time

384:   Output Parameter:
385: . ctx - the context

387:   Options Database Keys:
388: + -ts_monitor_lg_timestep        - automatically sets line graph monitor
389: . -ts_monitor_lg_timestep_log    - automatically sets line graph monitor
390: . -ts_monitor_lg_solution        - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`)
391: . -ts_monitor_lg_error           - monitor the error
392: . -ts_monitor_lg_ksp_iterations  - monitor the number of `KSP` iterations needed for each timestep
393: . -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep
394: - -lg_use_markers (true|false)   - mark the data points (at each time step) on the plot; default is true

396:   Level: intermediate

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

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

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

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

409: .seealso: [](ch_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
410:           `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
411:           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
412:           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
413:           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
414: @*/
415: PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
416: {
417:   PetscDraw draw;

419:   PetscFunctionBegin;
420:   PetscCall(PetscNew(ctx));
421:   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
422:   PetscCall(PetscDrawSetFromOptions(draw));
423:   PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg));
424:   PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg));
425:   PetscCall(PetscDrawDestroy(&draw));
426:   (*ctx)->howoften = howoften;
427:   PetscFunctionReturn(PETSC_SUCCESS);
428: }

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

433:   Collective

435:   Input Parameters:
436: + ts     - the time integrator
437: . step   - the current time step
438: . ptime  - the current time
439: . v      - the current state
440: - monctx - the monitor context obtained with `TSMonitorLGCtxCreate()`

442:   Level: advanced

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

448: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGCtxDestroy()`
449: @*/
450: PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscCtx monctx)
451: {
452:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
453:   PetscReal      x   = ptime, y;

455:   PetscFunctionBegin;
456:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */
457:   if (!step) {
458:     PetscDrawAxis axis;
459:     const char   *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
460:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
461:     PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel));
462:     PetscCall(PetscDrawLGReset(ctx->lg));
463:   }
464:   PetscCall(TSGetTimeStep(ts, &y));
465:   if (ctx->semilogy) y = PetscLog10Real(y);
466:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
467:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
468:     PetscCall(PetscDrawLGDraw(ctx->lg));
469:     PetscCall(PetscDrawLGSave(ctx->lg));
470:   }
471:   PetscFunctionReturn(PETSC_SUCCESS);
472: }

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

477:   Collective

479:   Input Parameter:
480: . ctx - the monitor context

482:   Level: intermediate

484:   Note:
485:   Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`

487: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
488: @*/
489: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
490: {
491:   PetscFunctionBegin;
492:   if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)(&(*ctx)->transformctx));
493:   PetscCall(PetscDrawLGDestroy(&(*ctx)->lg));
494:   PetscCall(PetscStrArrayDestroy(&(*ctx)->names));
495:   PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames));
496:   PetscCall(PetscFree((*ctx)->displayvariables));
497:   PetscCall(PetscFree((*ctx)->displayvalues));
498:   PetscCall(PetscFree(*ctx));
499:   PetscFunctionReturn(PETSC_SUCCESS);
500: }

502: /*@C
503:   TSMonitorSPCtxCreate - Creates a `TSMonitorSPCtx` scatter-plot monitor context for use with `DMSWARM` particle visualizations

505:   Collective

507:   Input Parameters:
508: + comm         - the MPI communicator to use
509: . host         - the X display to open, or `NULL` for the local machine
510: . label        - the title to put in the title bar
511: . x            - the x screen coordinates of the upper left coordinate of the window
512: . y            - the y screen coordinates of the upper left coordinate of the window
513: . m            - the screen width in pixels
514: . n            - the screen height in pixels
515: . howoften     - if positive then determines the frequency of the plotting, if -1 then only at the final time
516: . retain       - the number of old points to retain in the plot, or 0 to clear, or -1 to retain all
517: . phase        - `PETSC_TRUE` to plot in phase space rather than coordinate space
518: - multispecies - `PETSC_TRUE` to color particles by species

520:   Output Parameter:
521: . ctx - the newly created scatter plot monitor context

523:   Level: intermediate

525:   Note:
526:   Pass this context and `TSMonitorSPCtxDestroy()` to `TSMonitorSet()` with `TSMonitorSPSwarmSolution()` to display particles during the integration.

528: .seealso: [](ch_ts), `TS`, `DMSWARM`, `TSMonitorSet()`, `TSMonitorSPSwarmSolution()`, `TSMonitorSPCtxDestroy()`
529: @*/
530: 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)
531: {
532:   PetscDraw draw;

534:   PetscFunctionBegin;
535:   PetscCall(PetscNew(ctx));
536:   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
537:   PetscCall(PetscDrawSetFromOptions(draw));
538:   PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp));
539:   PetscCall(PetscDrawDestroy(&draw));
540:   (*ctx)->howoften     = howoften;
541:   (*ctx)->retain       = retain;
542:   (*ctx)->phase        = phase;
543:   (*ctx)->multispecies = multispecies;
544:   PetscFunctionReturn(PETSC_SUCCESS);
545: }

547: /*@C
548:   TSMonitorSPCtxDestroy - Destroys a `TSMonitorSPCtx` that was created with `TSMonitorSPCtxCreate()`

550:   Not Collective

552:   Input Parameter:
553: . ctx - the scatter plot monitor context

555:   Level: intermediate

557: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorSPCtxCreate()`, `TSMonitorSPSwarmSolution()`
558: @*/
559: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
560: {
561:   PetscFunctionBegin;
562:   PetscCall(PetscDrawSPDestroy(&(*ctx)->sp));
563:   PetscCall(PetscFree(*ctx));
564:   PetscFunctionReturn(PETSC_SUCCESS);
565: }

567: /*@C
568:   TSMonitorHGCtxCreate - Creates a `TSMonitorHGCtx` histogram monitor context for use with `DMSWARM` particle visualizations

570:   Collective

572:   Input Parameters:
573: + comm     - the MPI communicator to use
574: . host     - the X display to open, or `NULL` for the local machine
575: . label    - the title to put in the title bar
576: . x        - the x screen coordinates of the upper left coordinate of the window
577: . y        - the y screen coordinates of the upper left coordinate of the window
578: . m        - the screen width in pixels
579: . n        - the screen height in pixels
580: . howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
581: . Ns       - the number of species to histogram
582: . Nb       - the number of histogram bins
583: - velocity - `PETSC_TRUE` to plot histograms in velocity space, `PETSC_FALSE` for coordinate space

585:   Output Parameter:
586: . ctx - the newly created histogram monitor context

588:   Level: intermediate

590:   Note:
591:   Pass this context and `TSMonitorHGCtxDestroy()` to `TSMonitorSet()` with `TSMonitorHGSwarmSolution()` to display particle histograms during integration.

593: .seealso: [](ch_ts), `TS`, `DMSWARM`, `TSMonitorSet()`, `TSMonitorHGSwarmSolution()`, `TSMonitorHGCtxDestroy()`
594: @*/
595: 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)
596: {
597:   PetscDraw draw;
598:   int       Nsi, Nbi;

600:   PetscFunctionBegin;
601:   PetscCall(PetscMPIIntCast(Ns, &Nsi));
602:   PetscCall(PetscMPIIntCast(Nb, &Nbi));
603:   PetscCall(PetscNew(ctx));
604:   PetscCall(PetscMalloc1(Ns, &(*ctx)->hg));
605:   for (int s = 0; s < Nsi; ++s) {
606:     PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw));
607:     PetscCall(PetscDrawSetFromOptions(draw));
608:     PetscCall(PetscDrawHGCreate(draw, Nbi, &(*ctx)->hg[s]));
609:     PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE));
610:     PetscCall(PetscDrawDestroy(&draw));
611:   }
612:   (*ctx)->howoften = howoften;
613:   (*ctx)->Ns       = Ns;
614:   (*ctx)->velocity = velocity;
615:   PetscFunctionReturn(PETSC_SUCCESS);
616: }

618: /*@C
619:   TSMonitorHGCtxDestroy - Destroys a `TSMonitorHGCtx` that was created with `TSMonitorHGCtxCreate()`

621:   Not Collective

623:   Input Parameter:
624: . ctx - the histogram monitor context

626:   Level: intermediate

628: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorHGCtxCreate()`, `TSMonitorHGSwarmSolution()`
629: @*/
630: PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx)
631: {
632:   PetscFunctionBegin;
633:   for (PetscInt s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s]));
634:   PetscCall(PetscFree((*ctx)->hg));
635:   PetscCall(PetscFree(*ctx));
636:   PetscFunctionReturn(PETSC_SUCCESS);
637: }

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

643:   Collective

645:   Input Parameters:
646: + ts    - the `TS` context
647: . step  - current time-step
648: . ptime - current time
649: . u     - the solution at the current time
650: - ctx   - either a viewer or `NULL`

652:   Options Database Keys:
653: + -ts_monitor_draw_solution         - draw the solution at each time-step
654: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution

656:   Level: intermediate

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

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

665: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
666: @*/
667: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx ctx)
668: {
669:   TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)ctx;
670:   PetscDraw        draw;

672:   PetscFunctionBegin;
673:   if (!step && ictx->showinitial) {
674:     if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
675:     PetscCall(VecCopy(u, ictx->initialsolution));
676:   }
677:   if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);

679:   if (ictx->showinitial) {
680:     PetscReal pause;
681:     PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
682:     PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
683:     PetscCall(VecView(ictx->initialsolution, ictx->viewer));
684:     PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
685:     PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
686:   }
687:   PetscCall(VecView(u, ictx->viewer));
688:   if (ictx->showtimestepandtime) {
689:     PetscReal xl, yl, xr, yr, h;
690:     char      time[32];

692:     PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
693:     PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
694:     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
695:     h = yl + .95 * (yr - yl);
696:     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
697:     PetscCall(PetscDrawFlush(draw));
698:   }

700:   if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
701:   PetscFunctionReturn(PETSC_SUCCESS);
702: }

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

707:   Collective

709:   Input Parameters:
710: + ts    - the `TS` context
711: . step  - current time-step
712: . ptime - current time
713: . u     - the solution at the current time
714: - ctx   - either a viewer or `NULL`

716:   Level: intermediate

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

722: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
723: @*/
724: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx ctx)
725: {
726:   TSMonitorDrawCtx   ictx = (TSMonitorDrawCtx)ctx;
727:   PetscDraw          draw;
728:   PetscDrawAxis      axis;
729:   PetscInt           n;
730:   PetscMPIInt        size;
731:   PetscReal          U0, U1, xl, yl, xr, yr, h;
732:   char               time[32];
733:   const PetscScalar *U;

735:   PetscFunctionBegin;
736:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
737:   PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
738:   PetscCall(VecGetSize(u, &n));
739:   PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");

741:   PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
742:   PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
743:   PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
744:   if (!step) {
745:     PetscCall(PetscDrawClear(draw));
746:     PetscCall(PetscDrawAxisDraw(axis));
747:   }

749:   PetscCall(VecGetArrayRead(u, &U));
750:   U0 = PetscRealPart(U[0]);
751:   U1 = PetscRealPart(U[1]);
752:   PetscCall(VecRestoreArrayRead(u, &U));
753:   if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);

755:   PetscDrawCollectiveBegin(draw);
756:   PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
757:   if (ictx->showtimestepandtime) {
758:     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
759:     PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
760:     h = yl + .95 * (yr - yl);
761:     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
762:   }
763:   PetscDrawCollectiveEnd(draw);
764:   PetscCall(PetscDrawFlush(draw));
765:   PetscCall(PetscDrawPause(draw));
766:   PetscCall(PetscDrawSave(draw));
767:   PetscFunctionReturn(PETSC_SUCCESS);
768: }

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

773:   Collective

775:   Input Parameter:
776: . ictx - the monitor context

778:   Level: intermediate

780: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
781: @*/
782: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
783: {
784:   PetscFunctionBegin;
785:   PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
786:   PetscCall(VecDestroy(&(*ictx)->initialsolution));
787:   PetscCall(PetscFree(*ictx));
788:   PetscFunctionReturn(PETSC_SUCCESS);
789: }

791: /*@C
792:   TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`

794:   Collective

796:   Input Parameters:
797: + comm     - the MPI communicator to use
798: . host     - the X display to open, or `NULL` for the local machine
799: . label    - the title to put in the title bar
800: . x        - the x screen coordinates of the upper left coordinate of the window
801: . y        - the y screen coordinates of the upper left coordinate of the window
802: . m        - the screen width in pixels
803: . n        - the screen height in pixels
804: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time

806:   Output Parameter:
807: . ctx - the monitor context

809:   Options Database Keys:
810: + -ts_monitor_draw_solution         - draw the solution at each time-step
811: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution

813:   Level: intermediate

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

818: .seealso: [](ch_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
819: @*/
820: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
821: {
822:   PetscFunctionBegin;
823:   PetscCall(PetscNew(ctx));
824:   PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
825:   PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));

827:   (*ctx)->howoften    = howoften;
828:   (*ctx)->showinitial = PETSC_FALSE;
829:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));

831:   (*ctx)->showtimestepandtime = PETSC_FALSE;
832:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
833:   PetscFunctionReturn(PETSC_SUCCESS);
834: }

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

840:   Collective

842:   Input Parameters:
843: + ts    - the `TS` context
844: . step  - current time-step
845: . ptime - current time
846: . u     - solution at current time
847: - Ctx   - either a viewer or `NULL`

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

852:   Level: intermediate

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

858: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
859: @*/
860: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx Ctx)
861: {
862:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)Ctx;
863:   PetscViewer      viewer = ctx->viewer;
864:   Vec              work;

866:   PetscFunctionBegin;
867:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
868:   PetscCall(VecDuplicate(u, &work));
869:   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
870:   PetscCall(VecView(work, viewer));
871:   PetscCall(VecDestroy(&work));
872:   PetscFunctionReturn(PETSC_SUCCESS);
873: }

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

879:   Collective

881:   Input Parameters:
882: + ts    - the `TS` context
883: . step  - current time-step
884: . ptime - current time
885: . u     - solution at current time
886: - Ctx   - either a viewer or `NULL`

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

891:   Level: intermediate

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

897: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
898: @*/
899: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx Ctx)
900: {
901:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)Ctx;
902:   PetscViewer      viewer = ctx->viewer;
903:   Vec              work;

905:   PetscFunctionBegin;
906:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
907:   PetscCall(VecDuplicate(u, &work));
908:   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
909:   PetscCall(VecAXPY(work, -1.0, u));
910:   PetscCall(VecView(work, viewer));
911:   PetscCall(VecDestroy(&work));
912:   PetscFunctionReturn(PETSC_SUCCESS);
913: }

915: /*@C
916:   TSMonitorSolutionSetup - Setups the context for `TSMonitorSolution()`

918:   Collective

920:   Input Parameters:
921: + ts - the `TS` context
922: - vf - viewer and its format

924:   Level: intermediate

926: .seealso: [](ch_ts), `TS`, `TSMonitorSolution()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSetFromOptions()`
927: @*/
928: PetscErrorCode TSMonitorSolutionSetup(TS ts, PetscViewerAndFormat *vf)
929: {
930:   TSMonitorSolutionCtx ctx;

932:   PetscFunctionBegin;
933:   PetscCall(PetscNew(&ctx));
934:   PetscCall(PetscOptionsGetBool(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_solution_skip_initial", &ctx->skip_initial, NULL));
935:   vf->data         = ctx;
936:   vf->data_destroy = PetscCtxDestroyDefault;
937:   PetscFunctionReturn(PETSC_SUCCESS);
938: }

940: /*@C
941:   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

943:   Collective

945:   Input Parameters:
946: + ts    - the `TS` context
947: . step  - current time-step
948: . ptime - current time
949: . u     - current state
950: - vf    - viewer and its format

952:   Level: intermediate

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

958: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSolutionSetup()`
959: @*/
960: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
961: {
962:   TSMonitorSolutionCtx ctx = (TSMonitorSolutionCtx)vf->data;

964:   PetscFunctionBegin;
965:   if (ctx->skip_initial && step == ts->start_step) PetscFunctionReturn(PETSC_SUCCESS);
966:   if ((vf->view_interval > 0 && !(step % vf->view_interval)) || (vf->view_interval && ts->reason)) {
967:     PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
968:     PetscCall(VecView(u, vf->viewer));
969:     PetscCall(PetscViewerPopFormat(vf->viewer));
970:   }
971:   PetscFunctionReturn(PETSC_SUCCESS);
972: }

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

977:   Collective

979:   Input Parameters:
980: + ts    - the `TS` context
981: . step  - current time-step
982: . ptime - current time
983: . u     - current state
984: - ctx   - monitor context obtained with `TSMonitorSolutionVTKCtxCreate()`

986:   Level: developer

988:   Notes:
989:   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.
990:   These are named according to the file name template.

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

995: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
996: @*/
997: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, TSMonitorVTKCtx ctx)
998: {
999:   char        filename[PETSC_MAX_PATH_LEN];
1000:   PetscViewer viewer;

1002:   PetscFunctionBegin;
1003:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1004:   if (((ctx->interval > 0) && (!(step % ctx->interval))) || (ctx->interval && ts->reason)) {
1005:     PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)ctx->filenametemplate, step));
1006:     PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
1007:     PetscCall(VecView(u, viewer));
1008:     PetscCall(PetscViewerDestroy(&viewer));
1009:   }
1010:   PetscFunctionReturn(PETSC_SUCCESS);
1011: }

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

1016:   Not Collective

1018:   Input Parameter:
1019: . ctx - the monitor context

1021:   Level: developer

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

1026: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
1027: @*/
1028: PetscErrorCode TSMonitorSolutionVTKDestroy(TSMonitorVTKCtx *ctx)
1029: {
1030:   PetscFunctionBegin;
1031:   if (!*ctx) PetscFunctionReturn(PETSC_SUCCESS);
1032:   PetscCall(PetscFree((*ctx)->filenametemplate));
1033:   PetscCall(PetscFree(*ctx));
1034:   PetscFunctionReturn(PETSC_SUCCESS);
1035: }

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

1040:   Not collective

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

1045:   Output Parameter:
1046: . ctx - the monitor context

1048:   Level: developer

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

1053: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`, `TSMonitorSolutionVTKDestroy()`
1054: @*/
1055: PetscErrorCode TSMonitorSolutionVTKCtxCreate(const char *filenametemplate, TSMonitorVTKCtx *ctx)
1056: {
1057:   const char     *ptr = NULL, *ptr2 = NULL;
1058:   TSMonitorVTKCtx ictx;

1060:   PetscFunctionBegin;
1061:   PetscAssertPointer(filenametemplate, 1);
1062:   PetscAssertPointer(ctx, 2);
1063:   /* Do some cursory validation of the input. */
1064:   PetscCall(PetscStrstr(filenametemplate, "%", (char **)&ptr));
1065:   PetscCheck(ptr, PETSC_COMM_SELF, PETSC_ERR_USER, "-ts_monitor_solution_vtk requires a file template, e.g. filename-%%03" PetscInt_FMT ".vts");
1066:   for (ptr++; ptr && *ptr; ptr++) {
1067:     PetscCall(PetscStrchr("DdiouxX", *ptr, (char **)&ptr2));
1068:     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");
1069:     if (ptr2) break;
1070:   }
1071:   PetscCall(PetscNew(&ictx));
1072:   PetscCall(PetscStrallocpy(filenametemplate, &ictx->filenametemplate));
1073:   ictx->interval = 1;

1075:   *ctx = ictx;
1076:   PetscFunctionReturn(PETSC_SUCCESS);
1077: }

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

1083:   Collective

1085:   Input Parameters:
1086: + ts    - the `TS` context
1087: . step  - current time-step
1088: . ptime - current time
1089: . u     - current solution
1090: - dctx  - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`

1092:   Options Database Key:
1093: . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables

1095:   Level: intermediate

1097:   Notes:
1098:   Each process in a parallel run displays its component solutions in a separate window

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

1103: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
1104:           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
1105:           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
1106:           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
1107: @*/
1108: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1109: {
1110:   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)dctx;
1111:   const PetscScalar *yy;
1112:   Vec                v;

1114:   PetscFunctionBegin;
1115:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1116:   if (!step) {
1117:     PetscDrawAxis axis;
1118:     PetscInt      dim;
1119:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1120:     PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
1121:     if (!ctx->names) {
1122:       PetscBool flg;
1123:       /* user provides names of variables to plot but no names has been set so assume names are integer values */
1124:       PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
1125:       if (flg) {
1126:         PetscInt i, n;
1127:         char   **names;
1128:         PetscCall(VecGetSize(u, &n));
1129:         PetscCall(PetscMalloc1(n + 1, &names));
1130:         for (i = 0; i < n; i++) {
1131:           PetscCall(PetscMalloc1(5, &names[i]));
1132:           PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
1133:         }
1134:         names[n]   = NULL;
1135:         ctx->names = names;
1136:       }
1137:     }
1138:     if (ctx->names && !ctx->displaynames) {
1139:       char    **displaynames;
1140:       PetscBool flg;
1141:       PetscCall(VecGetLocalSize(u, &dim));
1142:       PetscCall(PetscCalloc1(dim + 1, &displaynames));
1143:       PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
1144:       if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
1145:       PetscCall(PetscStrArrayDestroy(&displaynames));
1146:     }
1147:     if (ctx->displaynames) {
1148:       PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
1149:       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
1150:     } else if (ctx->names) {
1151:       PetscCall(VecGetLocalSize(u, &dim));
1152:       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1153:       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
1154:     } else {
1155:       PetscCall(VecGetLocalSize(u, &dim));
1156:       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1157:     }
1158:     PetscCall(PetscDrawLGReset(ctx->lg));
1159:   }

1161:   if (!ctx->transform) v = u;
1162:   else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
1163:   PetscCall(VecGetArrayRead(v, &yy));
1164:   if (ctx->displaynames) {
1165:     PetscInt i;
1166:     for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
1167:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
1168:   } else {
1169: #if PetscDefined(USE_COMPLEX)
1170:     PetscInt   n;
1171:     PetscReal *yreal;
1172:     PetscCall(VecGetLocalSize(v, &n));
1173:     PetscCall(PetscMalloc1(n, &yreal));
1174:     for (PetscInt i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1175:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1176:     PetscCall(PetscFree(yreal));
1177: #else
1178:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1179: #endif
1180:   }
1181:   PetscCall(VecRestoreArrayRead(v, &yy));
1182:   if (ctx->transform) PetscCall(VecDestroy(&v));

1184:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1185:     PetscCall(PetscDrawLGDraw(ctx->lg));
1186:     PetscCall(PetscDrawLGSave(ctx->lg));
1187:   }
1188:   PetscFunctionReturn(PETSC_SUCCESS);
1189: }

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

1194:   Collective

1196:   Input Parameters:
1197: + ts    - the `TS` context
1198: - names - the names of the components, final string must be `NULL`

1200:   Level: intermediate

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

1205: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
1206: @*/
1207: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
1208: {
1209:   PetscInt i;

1211:   PetscFunctionBegin;
1212:   for (i = 0; i < ts->numbermonitors; i++) {
1213:     if (ts->monitor[i] == TSMonitorLGSolution) {
1214:       PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
1215:       break;
1216:     }
1217:   }
1218:   PetscFunctionReturn(PETSC_SUCCESS);
1219: }

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

1224:   Collective

1226:   Input Parameters:
1227: + ctx   - the `TS` context
1228: - names - the names of the components, final string must be `NULL`

1230:   Level: intermediate

1232: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
1233: @*/
1234: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
1235: {
1236:   PetscFunctionBegin;
1237:   PetscCall(PetscStrArrayDestroy(&ctx->names));
1238:   PetscCall(PetscStrArrayallocpy(names, &ctx->names));
1239:   PetscFunctionReturn(PETSC_SUCCESS);
1240: }

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

1245:   Collective

1247:   Input Parameter:
1248: . ts - the `TS` context

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

1253:   Level: intermediate

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

1258: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1259: @*/
1260: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
1261: {
1262:   PetscInt i;

1264:   PetscFunctionBegin;
1265:   *names = NULL;
1266:   for (i = 0; i < ts->numbermonitors; i++) {
1267:     if (ts->monitor[i] == TSMonitorLGSolution) {
1268:       TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
1269:       *names             = (const char *const *)ctx->names;
1270:       break;
1271:     }
1272:   }
1273:   PetscFunctionReturn(PETSC_SUCCESS);
1274: }

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

1279:   Collective

1281:   Input Parameters:
1282: + ctx          - the `TSMonitorLG` context
1283: - displaynames - the names of the components, final string must be `NULL`

1285:   Level: intermediate

1287: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1288: @*/
1289: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
1290: {
1291:   PetscInt j = 0, k;

1293:   PetscFunctionBegin;
1294:   if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
1295:   PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
1296:   PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
1297:   while (displaynames[j]) j++;
1298:   ctx->ndisplayvariables = j;
1299:   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
1300:   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
1301:   j = 0;
1302:   while (displaynames[j]) {
1303:     k = 0;
1304:     while (ctx->names[k]) {
1305:       PetscBool flg;
1306:       PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
1307:       if (flg) {
1308:         ctx->displayvariables[j] = k;
1309:         break;
1310:       }
1311:       k++;
1312:     }
1313:     j++;
1314:   }
1315:   PetscFunctionReturn(PETSC_SUCCESS);
1316: }

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

1321:   Collective

1323:   Input Parameters:
1324: + ts           - the `TS` context
1325: - displaynames - the names of the components, final string must be `NULL`

1327:   Level: intermediate

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

1332: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1333: @*/
1334: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
1335: {
1336:   PetscInt i;

1338:   PetscFunctionBegin;
1339:   for (i = 0; i < ts->numbermonitors; i++) {
1340:     if (ts->monitor[i] == TSMonitorLGSolution) {
1341:       PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
1342:       break;
1343:     }
1344:   }
1345:   PetscFunctionReturn(PETSC_SUCCESS);
1346: }

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

1351:   Collective

1353:   Input Parameters:
1354: + ts        - the `TS` context
1355: . transform - the transform function
1356: . destroy   - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1357: - tctx      - optional context used by transform function

1359:   Calling sequence of `transform`:
1360: + tctx - context used by the transform function
1361: . u    - the input solution vector
1362: - w    - the output transformed vector

1364:   Level: intermediate

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

1369: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorLGCtxSetTransform()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `PetscCtxDestroyFn`
1370: @*/
1371: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(PetscCtx tctx, Vec u, Vec *w), PetscCtxDestroyFn *destroy, PetscCtx tctx)
1372: {
1373:   PetscFunctionBegin;
1374:   for (PetscInt i = 0; i < ts->numbermonitors; i++) {
1375:     if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
1376:   }
1377:   PetscFunctionReturn(PETSC_SUCCESS);
1378: }

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

1383:   Collective

1385:   Input Parameters:
1386: + tctx      - the `TS` context
1387: . transform - the transform function
1388: . destroy   - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1389: - ctx       - optional context used by transform function

1391:   Calling sequence of `transform`:
1392: + tctx - context used by the transform function
1393: . u    - the input solution vector
1394: - w    - the output transformed vector

1396:   Level: intermediate

1398: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`, `PetscCtxDestroyFn`
1399: @*/
1400: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(PetscCtx tctx, Vec u, Vec *w), PetscCtxDestroyFn *destroy, PetscCtx tctx)
1401: {
1402:   PetscFunctionBegin;
1403:   ctx->transform        = transform;
1404:   ctx->transformdestroy = destroy;
1405:   ctx->transformctx     = tctx;
1406:   PetscFunctionReturn(PETSC_SUCCESS);
1407: }

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

1413:   Collective

1415:   Input Parameters:
1416: + ts    - the `TS` context
1417: . step  - current time-step
1418: . ptime - current time
1419: . u     - current solution
1420: - Ctx   - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`

1422:   Options Database Key:
1423: . -ts_monitor_lg_error - create a graphical monitor of error history

1425:   Level: intermediate

1427:   Notes:
1428:   Each process in a parallel run displays its component errors in a separate window

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

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

1435: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1436: @*/
1437: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx Ctx)
1438: {
1439:   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)Ctx;
1440:   const PetscScalar *yy;
1441:   Vec                y;

1443:   PetscFunctionBegin;
1444:   if (!step) {
1445:     PetscDrawAxis axis;
1446:     PetscInt      dim;
1447:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1448:     PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
1449:     PetscCall(VecGetLocalSize(u, &dim));
1450:     PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1451:     PetscCall(PetscDrawLGReset(ctx->lg));
1452:   }
1453:   PetscCall(VecDuplicate(u, &y));
1454:   PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1455:   PetscCall(VecAXPY(y, -1.0, u));
1456:   PetscCall(VecGetArrayRead(y, &yy));
1457: #if PetscDefined(USE_COMPLEX)
1458:   {
1459:     PetscReal *yreal;
1460:     PetscInt   n;
1461:     PetscCall(VecGetLocalSize(y, &n));
1462:     PetscCall(PetscMalloc1(n, &yreal));
1463:     for (PetscInt i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1464:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1465:     PetscCall(PetscFree(yreal));
1466:   }
1467: #else
1468:   PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1469: #endif
1470:   PetscCall(VecRestoreArrayRead(y, &yy));
1471:   PetscCall(VecDestroy(&y));
1472:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1473:     PetscCall(PetscDrawLGDraw(ctx->lg));
1474:     PetscCall(PetscDrawLGSave(ctx->lg));
1475:   }
1476:   PetscFunctionReturn(PETSC_SUCCESS);
1477: }

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

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

1489:   Options Database Keys:
1490: + -ts_monitor_sp_swarm n                          - Monitor the solution every n steps, or -1 for plotting only the final solution
1491: . -ts_monitor_sp_swarm_retain n                   - Retain n old points so we can see the history, or -1 for all points
1492: . -ts_monitor_sp_swarm_multi_species (true|false) - Color each species differently
1493: - -ts_monitor_sp_swarm_phase (true|false)         - Plot in phase space, as opposed to coordinate space

1495:   Level: intermediate

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

1501: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
1502: @*/
1503: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx dctx)
1504: {
1505:   TSMonitorSPCtx     ctx = (TSMonitorSPCtx)dctx;
1506:   PetscDraw          draw;
1507:   DM                 dm, cdm;
1508:   const PetscScalar *yy;
1509:   PetscInt           Np, p, dim = 2, *species;
1510:   PetscReal          species_color;

1512:   PetscFunctionBegin;
1513:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1514:   PetscCall(TSGetDM(ts, &dm));
1515:   if (!step) {
1516:     PetscDrawAxis axis;
1517:     PetscReal     dmboxlower[2], dmboxupper[2];

1519:     PetscCall(TSGetDM(ts, &dm));
1520:     PetscCall(DMGetDimension(dm, &dim));
1521:     PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
1522:     PetscCall(DMSwarmGetCellDM(dm, &cdm));
1523:     PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
1524:     PetscCall(VecGetLocalSize(u, &Np));
1525:     Np /= dim * 2;
1526:     PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
1527:     if (ctx->phase) {
1528:       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
1529:       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10));
1530:     } else {
1531:       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
1532:       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
1533:     }
1534:     PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
1535:     PetscCall(PetscDrawSPReset(ctx->sp));
1536:   }
1537:   if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species));
1538:   PetscCall(VecGetLocalSize(u, &Np));
1539:   Np /= dim * 2;
1540:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1541:     PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
1542:     if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
1543:     PetscCall(PetscDrawFlush(draw));
1544:     PetscCall(PetscDrawSPReset(ctx->sp));
1545:     PetscCall(VecGetArrayRead(u, &yy));
1546:     for (p = 0; p < Np; ++p) {
1547:       PetscReal x, y;

1549:       if (ctx->phase) {
1550:         x = PetscRealPart(yy[p * dim * 2]);
1551:         y = PetscRealPart(yy[p * dim * 2 + dim]);
1552:       } else {
1553:         x = PetscRealPart(yy[p * dim * 2]);
1554:         y = PetscRealPart(yy[p * dim * 2 + 1]);
1555:       }
1556:       if (ctx->multispecies) {
1557:         species_color = species[p] + 2;
1558:         PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color));
1559:       } else {
1560:         PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1561:       }
1562:       PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1563:     }
1564:     PetscCall(VecRestoreArrayRead(u, &yy));
1565:     PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
1566:     PetscCall(PetscDrawSPSave(ctx->sp));
1567:     if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species));
1568:   }
1569:   PetscFunctionReturn(PETSC_SUCCESS);
1570: }

1572: /*@C
1573:   TSMonitorHGSwarmSolution - Graphically displays histograms of `DMSWARM` particles

1575:   Input Parameters:
1576: + ts    - the `TS` context
1577: . step  - current time-step
1578: . ptime - current time
1579: . u     - current solution
1580: - dctx  - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorHGCtxCreate()`

1582:   Options Database Keys:
1583: + -ts_monitor_hg_swarm n                     - Monitor the solution every n steps, or -1 for plotting only the final solution
1584: . -ts_monitor_hg_swarm_species num           - Number of species to histogram
1585: . -ts_monitor_hg_swarm_bins num              - Number of histogram bins
1586: - -ts_monitor_hg_swarm_velocity (true|false) - Plot in velocity space, as opposed to coordinate space

1588:   Level: intermediate

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

1594: .seealso: `TSMonitorSet()`
1595: @*/
1596: PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx dctx)
1597: {
1598:   TSMonitorHGCtx     ctx = (TSMonitorHGCtx)dctx;
1599:   PetscDraw          draw;
1600:   DM                 sw;
1601:   const PetscScalar *yy;
1602:   PetscInt          *species;
1603:   PetscInt           dim, d = 0, Np, p, Ns, s;

1605:   PetscFunctionBegin;
1606:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1607:   PetscCall(TSGetDM(ts, &sw));
1608:   PetscCall(DMGetDimension(sw, &dim));
1609:   PetscCall(DMSwarmGetNumSpecies(sw, &Ns));
1610:   Ns = PetscMin(Ns, ctx->Ns);
1611:   PetscCall(VecGetLocalSize(u, &Np));
1612:   Np /= dim * 2;
1613:   if (!step) {
1614:     PetscDrawAxis axis;
1615:     char          title[PETSC_MAX_PATH_LEN];

1617:     for (s = 0; s < Ns; ++s) {
1618:       PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis));
1619:       PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s));
1620:       if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N"));
1621:       else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N"));
1622:     }
1623:   }
1624:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1625:     PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species));
1626:     for (s = 0; s < Ns; ++s) {
1627:       PetscCall(PetscDrawHGReset(ctx->hg[s]));
1628:       PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw));
1629:       PetscCall(PetscDrawClear(draw));
1630:       PetscCall(PetscDrawFlush(draw));
1631:     }
1632:     PetscCall(VecGetArrayRead(u, &yy));
1633:     for (p = 0; p < Np; ++p) {
1634:       const PetscInt s = species[p] < Ns ? species[p] : 0;
1635:       PetscReal      v;

1637:       if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]);
1638:       else v = PetscRealPart(yy[p * dim * 2 + d]);
1639:       PetscCall(PetscDrawHGAddValue(ctx->hg[s], v));
1640:     }
1641:     PetscCall(VecRestoreArrayRead(u, &yy));
1642:     for (s = 0; s < Ns; ++s) {
1643:       PetscCall(PetscDrawHGDraw(ctx->hg[s]));
1644:       PetscCall(PetscDrawHGSave(ctx->hg[s]));
1645:     }
1646:     PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species));
1647:   }
1648:   PetscFunctionReturn(PETSC_SUCCESS);
1649: }

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

1654:   Collective

1656:   Input Parameters:
1657: + ts    - the `TS` context
1658: . step  - current time-step
1659: . ptime - current time
1660: . u     - current solution
1661: - vf    - unused context

1663:   Options Database Key:
1664: . -ts_monitor_error - create a graphical monitor of error history

1666:   Level: intermediate

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

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

1674: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1675: @*/
1676: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1677: {
1678:   DM        dm;
1679:   PetscDS   ds = NULL;
1680:   PetscInt  Nf = -1, f;
1681:   PetscBool flg;

1683:   PetscFunctionBegin;
1684:   PetscCall(TSGetDM(ts, &dm));
1685:   if (dm) PetscCall(DMGetDS(dm, &ds));
1686:   if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
1687:   if (Nf <= 0) {
1688:     Vec       y;
1689:     PetscReal nrm;

1691:     PetscCall(VecDuplicate(u, &y));
1692:     PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1693:     PetscCall(VecAXPY(y, -1.0, u));
1694:     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
1695:     if (flg) {
1696:       PetscCall(VecNorm(y, NORM_2, &nrm));
1697:       PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
1698:     }
1699:     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
1700:     if (flg) PetscCall(VecView(y, vf->viewer));
1701:     PetscCall(VecDestroy(&y));
1702:   } else {
1703:     PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, PetscCtx ctx);
1704:     void    **ctxs;
1705:     Vec       v;
1706:     PetscReal ferrors[1];

1708:     PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
1709:     for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
1710:     PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
1711:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04" PetscInt_FMT " time = %-8.4g \t L_2 Error: [", step, (double)ptime));
1712:     for (f = 0; f < Nf; ++f) {
1713:       if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
1714:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
1715:     }
1716:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));

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

1720:     PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
1721:     if (flg) {
1722:       PetscCall(DMGetGlobalVector(dm, &v));
1723:       PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
1724:       PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
1725:       PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
1726:       PetscCall(DMRestoreGlobalVector(dm, &v));
1727:     }
1728:     PetscCall(PetscFree2(exactFuncs, ctxs));
1729:   }
1730:   PetscFunctionReturn(PETSC_SUCCESS);
1731: }

1733: /*@C
1734:   TSMonitorLGSNESIterations - Monitors the number of nonlinear (`SNES`) iterations used per time step in a line-graph plot

1736:   Collective

1738:   Input Parameters:
1739: + ts     - the `TS` context
1740: . n      - iteration number (a negative value indicates an interpolated solution and is ignored)
1741: . ptime  - current time
1742: . v      - current solution
1743: - monctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, created with `TSMonitorLGCtxCreate()`

1745:   Level: intermediate

1747:   Note:
1748:   This is not called directly by users; pass this function to `TSMonitorSet()` along with the context created by `TSMonitorLGCtxCreate()` and `TSMonitorLGCtxDestroy()`.

1750: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGKSPIterations()`
1751: @*/
1752: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, PetscCtx monctx)
1753: {
1754:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1755:   PetscReal      x   = ptime, y;
1756:   PetscInt       its;

1758:   PetscFunctionBegin;
1759:   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1760:   if (!n) {
1761:     PetscDrawAxis axis;
1762:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1763:     PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
1764:     PetscCall(PetscDrawLGReset(ctx->lg));
1765:     ctx->snes_its = 0;
1766:   }
1767:   PetscCall(TSGetSNESIterations(ts, &its));
1768:   y = its - ctx->snes_its;
1769:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1770:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1771:     PetscCall(PetscDrawLGDraw(ctx->lg));
1772:     PetscCall(PetscDrawLGSave(ctx->lg));
1773:   }
1774:   ctx->snes_its = its;
1775:   PetscFunctionReturn(PETSC_SUCCESS);
1776: }

1778: /*@C
1779:   TSMonitorLGKSPIterations - Monitors the number of linear (`KSP`) iterations used per time step in a line-graph plot

1781:   Collective

1783:   Input Parameters:
1784: + ts     - the `TS` context
1785: . n      - iteration number (a negative value indicates an interpolated solution and is ignored)
1786: . ptime  - current time
1787: . v      - current solution
1788: - monctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, created with `TSMonitorLGCtxCreate()`

1790:   Level: intermediate

1792:   Note:
1793:   This is not called directly by users; pass this function to `TSMonitorSet()` along with the context created by `TSMonitorLGCtxCreate()` and `TSMonitorLGCtxDestroy()`.

1795: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGSNESIterations()`
1796: @*/
1797: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, PetscCtx monctx)
1798: {
1799:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1800:   PetscReal      x   = ptime, y;
1801:   PetscInt       its;

1803:   PetscFunctionBegin;
1804:   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1805:   if (!n) {
1806:     PetscDrawAxis axis;
1807:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1808:     PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
1809:     PetscCall(PetscDrawLGReset(ctx->lg));
1810:     ctx->ksp_its = 0;
1811:   }
1812:   PetscCall(TSGetKSPIterations(ts, &its));
1813:   y = its - ctx->ksp_its;
1814:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1815:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1816:     PetscCall(PetscDrawLGDraw(ctx->lg));
1817:     PetscCall(PetscDrawLGSave(ctx->lg));
1818:   }
1819:   ctx->ksp_its = its;
1820:   PetscFunctionReturn(PETSC_SUCCESS);
1821: }

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

1826:   Collective

1828:   Input Parameter:
1829: . ts - the `TS` solver object

1831:   Output Parameter:
1832: . ctx - the context

1834:   Level: intermediate

1836: .seealso: [](ch_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1837: @*/
1838: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1839: {
1840:   PetscFunctionBegin;
1841:   PetscCall(PetscNew(ctx));
1842:   PetscFunctionReturn(PETSC_SUCCESS);
1843: }

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

1848:   Collective

1850:   Input Parameters:
1851: + ts    - the `TS` context
1852: . step  - current time-step
1853: . ptime - current time
1854: . u     - current solution
1855: - dctx  - the envelope context

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

1860:   Level: intermediate

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

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

1868: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1869: @*/
1870: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx dctx)
1871: {
1872:   TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;

1874:   PetscFunctionBegin;
1875:   if (!ctx->max) {
1876:     PetscCall(VecDuplicate(u, &ctx->max));
1877:     PetscCall(VecDuplicate(u, &ctx->min));
1878:     PetscCall(VecCopy(u, ctx->max));
1879:     PetscCall(VecCopy(u, ctx->min));
1880:   } else {
1881:     PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
1882:     PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
1883:   }
1884:   PetscFunctionReturn(PETSC_SUCCESS);
1885: }

1887: /*@C
1888:   TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution

1890:   Collective

1892:   Input Parameter:
1893: . ts - the `TS` context

1895:   Output Parameters:
1896: + max - the maximum values
1897: - min - the minimum values

1899:   Level: intermediate

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

1904: .seealso: [](ch_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1905: @*/
1906: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1907: {
1908:   PetscFunctionBegin;
1909:   if (max) *max = NULL;
1910:   if (min) *min = NULL;
1911:   for (PetscInt i = 0; i < ts->numbermonitors; i++) {
1912:     if (ts->monitor[i] == TSMonitorEnvelope) {
1913:       TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1914:       if (max) *max = ctx->max;
1915:       if (min) *min = ctx->min;
1916:       break;
1917:     }
1918:   }
1919:   PetscFunctionReturn(PETSC_SUCCESS);
1920: }

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

1925:   Collective

1927:   Input Parameter:
1928: . ctx - the monitor context

1930:   Level: intermediate

1932: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1933: @*/
1934: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1935: {
1936:   PetscFunctionBegin;
1937:   PetscCall(VecDestroy(&(*ctx)->min));
1938:   PetscCall(VecDestroy(&(*ctx)->max));
1939:   PetscCall(PetscFree(*ctx));
1940:   PetscFunctionReturn(PETSC_SUCCESS);
1941: }

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

1946:   Not Collective

1948:   Input Parameters:
1949: + ts   - the `TS` context
1950: . step - current timestep
1951: . t    - current time
1952: . U    - current solution
1953: - vf   - not used

1955:   Options Database Key:
1956: + -ts_dmswarm_monitor_moments          - Monitor moments of particle distribution
1957: - -ts_dmswarm_monitor_moments_interval - Interval of timesteps between monitor outputs

1959:   Level: intermediate

1961:   Notes:
1962:   This requires a `DMSWARM` be attached to the `TS`.

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

1967: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1968: @*/
1969: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1970: {
1971:   DM                 sw;
1972:   const PetscScalar *u;
1973:   PetscReal          m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1974:   PetscInt           dim, d, Np, p;
1975:   MPI_Comm           comm;

1977:   PetscFunctionBeginUser;
1978:   (void)t;
1979:   (void)vf;
1980:   PetscCall(TSGetDM(ts, &sw));
1981:   if (!sw || step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS);
1982:   PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
1983:   PetscCall(DMGetDimension(sw, &dim));
1984:   PetscCall(VecGetLocalSize(U, &Np));
1985:   Np /= dim;
1986:   PetscCall(VecGetArrayRead(U, &u));
1987:   for (p = 0; p < Np; ++p) {
1988:     for (d = 0; d < dim; ++d) {
1989:       totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1990:       totMom[d] += PetscRealPart(u[p * dim + d]);
1991:     }
1992:   }
1993:   PetscCall(VecRestoreArrayRead(U, &u));
1994:   for (d = 0; d < dim; ++d) totMom[d] *= m;
1995:   totE *= 0.5 * m;
1996:   PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
1997:   for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, "    Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
1998:   PetscCall(PetscPrintf(comm, "\n"));
1999:   PetscFunctionReturn(PETSC_SUCCESS);
2000: }