xref: /petsc/src/ts/interface/tsmon.c (revision 195e9b0206ae3287e21fe0eadbd7bac05829bf59)
1 #include <petsc/private/tsimpl.h> /*I "petscts.h"  I*/
2 #include <petscdm.h>
3 #include <petscds.h>
4 #include <petscdmswarm.h>
5 #include <petscdraw.h>
6 
7 /*@C
8    TSMonitor - Runs all user-provided monitor routines set using `TSMonitorSet()`
9 
10    Collective
11 
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
17 
18    Level: developer
19 
20    Notes:
21    `TSMonitor()` is typically used automatically within the time stepping implementations.
22    Users would almost never call this routine directly.
23 
24    A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions
25 
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;
32 
33   PetscFunctionBegin;
34   PetscValidHeaderSpecific(ts, TS_CLASSID, 1);
35   PetscValidHeaderSpecific(u, VEC_CLASSID, 4);
36 
37   PetscCall(TSGetDM(ts, &dm));
38   PetscCall(DMSetOutputSequenceNumber(dm, step, ptime));
39 
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 }
45 
46 /*@C
47    TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
48 
49    Collective
50 
51    Input Parameters:
52 +  ts - `TS` object you wish to monitor
53 .  name - the monitor type one is seeking
54 .  help - message indicating what monitoring is done
55 .  manual - manual page for the monitor
56 .  monitor - the monitor function
57 -  monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the `TS` or `PetscViewer` objects
58 
59    Level: developer
60 
61 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `PetscOptionsGetViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
62           `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
63           `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`,
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;
74 
75   PetscFunctionBegin;
76   PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)ts), ((PetscObject)ts)->options, ((PetscObject)ts)->prefix, name, &viewer, &format, &flg));
77   if (flg) {
78     PetscViewerAndFormat *vf;
79     char                  interval_key[1024];
80     PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
81     PetscCall(PetscSNPrintf(interval_key, sizeof interval_key, "%s_interval", name));
82     PetscCall(PetscOptionsGetInt(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, interval_key, &vf->view_interval, NULL));
83     PetscCall(PetscObjectDereference((PetscObject)viewer));
84     if (monitorsetup) PetscCall((*monitorsetup)(ts, vf));
85     PetscCall(TSMonitorSet(ts, (PetscErrorCode(*)(TS, PetscInt, PetscReal, Vec, void *))monitor, vf, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy));
86   }
87   PetscFunctionReturn(PETSC_SUCCESS);
88 }
89 
90 /*@C
91    TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
92    timestep to display the iteration's  progress.
93 
94    Logically Collective
95 
96    Input Parameters:
97 +  ts - the `TS` context obtained from `TSCreate()`
98 .  monitor - monitoring routine
99 .  mctx - [optional] user-defined context for private data for the monitor routine (use `NULL` if no context is desired)
100 -  monitordestroy - [optional] routine that frees monitor context (may be `NULL`)
101 
102    Calling sequence of monitor:
103 $    PetscErrorCode monitor(TS ts,PetscInt steps,PetscReal time,Vec u,void *mctx)
104 
105 +    ts - the TS context
106 .    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)
107 .    time - current time
108 .    u - current iterate
109 -    mctx - [optional] monitoring context
110 
111    Level: intermediate
112 
113    Note:
114    This routine adds an additional monitor to the list of monitors that already has been loaded.
115 
116    Fortran Note:
117    Only a single monitor function can be set for each `TS` object
118 
119 .seealso: [](chapter_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
120           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
121           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
122 @*/
123 PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, void *), void *mctx, PetscErrorCode (*mdestroy)(void **))
124 {
125   PetscInt  i;
126   PetscBool identical;
127 
128   PetscFunctionBegin;
129   PetscValidHeaderSpecific(ts, TS_CLASSID, 1);
130   for (i = 0; i < ts->numbermonitors; i++) {
131     PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))monitor, mctx, mdestroy, (PetscErrorCode(*)(void))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++] = (void *)mctx;
138   PetscFunctionReturn(PETSC_SUCCESS);
139 }
140 
141 /*@C
142    TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
143 
144    Logically Collective
145 
146    Input Parameters:
147 .  ts - the `TS` context obtained from `TSCreate()`
148 
149    Level: intermediate
150 
151    Note:
152    There is no way to remove a single, specific monitor.
153 
154 .seealso: [](chapter_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
155 @*/
156 PetscErrorCode TSMonitorCancel(TS ts)
157 {
158   PetscInt i;
159 
160   PetscFunctionBegin;
161   PetscValidHeaderSpecific(ts, TS_CLASSID, 1);
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 }
168 
169 /*@C
170    TSMonitorDefault - The default monitor, prints the timestep and time for each step
171 
172    Options Database Key:
173 .  -ts_monitor - monitors the time integration
174 
175    Level: intermediate
176 
177    Notes:
178    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
179    to be used during the `TS` integration.
180 
181 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
182           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
183           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
184 @*/
185 PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
186 {
187   PetscViewer viewer = vf->viewer;
188   PetscBool   iascii, ibinary;
189 
190   PetscFunctionBegin;
191   PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 5);
192   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
193   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary));
194   PetscCall(PetscViewerPushFormat(viewer, vf->format));
195   if (iascii) {
196     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
197     if (step == -1) { /* this indicates it is an interpolated solution */
198       PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps));
199     } else {
200       PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n"));
201     }
202     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
203   } else if (ibinary) {
204     PetscMPIInt rank;
205     PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
206     if (rank == 0) {
207       PetscBool skipHeader;
208       PetscInt  classid = REAL_FILE_CLASSID;
209 
210       PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
211       if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
212       PetscCall(PetscRealView(1, &ptime, viewer));
213     } else {
214       PetscCall(PetscRealView(0, &ptime, viewer));
215     }
216   }
217   PetscCall(PetscViewerPopFormat(viewer));
218   PetscFunctionReturn(PETSC_SUCCESS);
219 }
220 
221 /*@C
222    TSMonitorExtreme - Prints the extreme values of the solution at each timestep
223 
224    Level: intermediate
225 
226    Note:
227    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
228    to be used during the `TS` integration.
229 
230 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`
231 @*/
232 PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
233 {
234   PetscViewer viewer = vf->viewer;
235   PetscBool   iascii;
236   PetscReal   max, min;
237 
238   PetscFunctionBegin;
239   PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 5);
240   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
241   PetscCall(PetscViewerPushFormat(viewer, vf->format));
242   if (iascii) {
243     PetscCall(VecMax(v, NULL, &max));
244     PetscCall(VecMin(v, NULL, &min));
245     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
246     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));
247     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
248   }
249   PetscCall(PetscViewerPopFormat(viewer));
250   PetscFunctionReturn(PETSC_SUCCESS);
251 }
252 
253 /*@C
254    TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with
255    `TS` to monitor the solution process graphically in various ways
256 
257    Collective
258 
259    Input Parameters:
260 +  host - the X display to open, or `NULL` for the local machine
261 .  label - the title to put in the title bar
262 .  x, y - the screen coordinates of the upper left coordinate of the window
263 .  m, n - the screen width and height in pixels
264 -  howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
265 
266    Output Parameter:
267 .  ctx - the context
268 
269    Options Database Keys:
270 +  -ts_monitor_lg_timestep - automatically sets line graph monitor
271 +  -ts_monitor_lg_timestep_log - automatically sets line graph monitor
272 .  -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`)
273 .  -ts_monitor_lg_error -  monitor the error
274 .  -ts_monitor_lg_ksp_iterations - monitor the number of `KSP` iterations needed for each timestep
275 .  -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep
276 -  -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true
277 
278    Level: intermediate
279 
280    Notes:
281    Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed.
282 
283    One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()`
284 
285    Many of the functions that control the monitoring have two forms: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a `TS` object as the
286    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
287    as the first argument.
288 
289    One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()`
290 
291 .seealso: [](chapter_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
292           `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
293           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
294           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
295           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
296 @*/
297 PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
298 {
299   PetscDraw draw;
300 
301   PetscFunctionBegin;
302   PetscCall(PetscNew(ctx));
303   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
304   PetscCall(PetscDrawSetFromOptions(draw));
305   PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg));
306   PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg));
307   PetscCall(PetscDrawDestroy(&draw));
308   (*ctx)->howoften = howoften;
309   PetscFunctionReturn(PETSC_SUCCESS);
310 }
311 
312 PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx)
313 {
314   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
315   PetscReal      x   = ptime, y;
316 
317   PetscFunctionBegin;
318   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */
319   if (!step) {
320     PetscDrawAxis axis;
321     const char   *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
322     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
323     PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel));
324     PetscCall(PetscDrawLGReset(ctx->lg));
325   }
326   PetscCall(TSGetTimeStep(ts, &y));
327   if (ctx->semilogy) y = PetscLog10Real(y);
328   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
329   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
330     PetscCall(PetscDrawLGDraw(ctx->lg));
331     PetscCall(PetscDrawLGSave(ctx->lg));
332   }
333   PetscFunctionReturn(PETSC_SUCCESS);
334 }
335 
336 /*@C
337    TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`.
338 
339    Collective
340 
341    Input Parameter:
342 .  ctx - the monitor context
343 
344    Level: intermediate
345 
346    Note:
347    Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`
348 
349 .seealso: [](chapter_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep();`
350 @*/
351 PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
352 {
353   PetscFunctionBegin;
354   if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)((*ctx)->transformctx));
355   PetscCall(PetscDrawLGDestroy(&(*ctx)->lg));
356   PetscCall(PetscStrArrayDestroy(&(*ctx)->names));
357   PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames));
358   PetscCall(PetscFree((*ctx)->displayvariables));
359   PetscCall(PetscFree((*ctx)->displayvalues));
360   PetscCall(PetscFree(*ctx));
361   PetscFunctionReturn(PETSC_SUCCESS);
362 }
363 
364 /* Creates a TS Monitor SPCtx for use with DMSwarm particle visualizations */
365 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, TSMonitorSPCtx *ctx)
366 {
367   PetscDraw draw;
368 
369   PetscFunctionBegin;
370   PetscCall(PetscNew(ctx));
371   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
372   PetscCall(PetscDrawSetFromOptions(draw));
373   PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp));
374   PetscCall(PetscDrawDestroy(&draw));
375   (*ctx)->howoften = howoften;
376   (*ctx)->retain   = retain;
377   (*ctx)->phase    = phase;
378   PetscFunctionReturn(PETSC_SUCCESS);
379 }
380 
381 /*
382   Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate
383 */
384 PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
385 {
386   PetscFunctionBegin;
387 
388   PetscCall(PetscDrawSPDestroy(&(*ctx)->sp));
389   PetscCall(PetscFree(*ctx));
390 
391   PetscFunctionReturn(PETSC_SUCCESS);
392 }
393 
394 /*@C
395    TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling
396    `VecView()` for the solution at each timestep
397 
398    Collective
399 
400    Input Parameters:
401 +  ts - the `TS` context
402 .  step - current time-step
403 .  ptime - current time
404 -  dummy - either a viewer or `NULL`
405 
406    Options Database Keys:
407 +   -ts_monitor_draw_solution - draw the solution at each time-step
408 -   -ts_monitor_draw_solution_initial - show initial solution as well as current solution
409 
410    Level: intermediate
411 
412    Notes:
413    The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial`
414    will look bad
415 
416    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with
417    `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration.
418 
419 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
420 @*/
421 PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
422 {
423   TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
424   PetscDraw        draw;
425 
426   PetscFunctionBegin;
427   if (!step && ictx->showinitial) {
428     if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
429     PetscCall(VecCopy(u, ictx->initialsolution));
430   }
431   if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
432 
433   if (ictx->showinitial) {
434     PetscReal pause;
435     PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
436     PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
437     PetscCall(VecView(ictx->initialsolution, ictx->viewer));
438     PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
439     PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
440   }
441   PetscCall(VecView(u, ictx->viewer));
442   if (ictx->showtimestepandtime) {
443     PetscReal xl, yl, xr, yr, h;
444     char      time[32];
445 
446     PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
447     PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime));
448     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
449     h = yl + .95 * (yr - yl);
450     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
451     PetscCall(PetscDrawFlush(draw));
452   }
453 
454   if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
455   PetscFunctionReturn(PETSC_SUCCESS);
456 }
457 
458 /*@C
459    TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram
460 
461    Collective
462 
463    Input Parameters:
464 +  ts - the `TS` context
465 .  step - current time-step
466 .  ptime - current time
467 -  dummy - either a viewer or `NULL`
468 
469    Level: intermediate
470 
471    Notes:
472    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
473    to be used during the `TS` integration.
474 
475 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
476 @*/
477 PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
478 {
479   TSMonitorDrawCtx   ictx = (TSMonitorDrawCtx)dummy;
480   PetscDraw          draw;
481   PetscDrawAxis      axis;
482   PetscInt           n;
483   PetscMPIInt        size;
484   PetscReal          U0, U1, xl, yl, xr, yr, h;
485   char               time[32];
486   const PetscScalar *U;
487 
488   PetscFunctionBegin;
489   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
490   PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
491   PetscCall(VecGetSize(u, &n));
492   PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");
493 
494   PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
495   PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
496   PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
497   if (!step) {
498     PetscCall(PetscDrawClear(draw));
499     PetscCall(PetscDrawAxisDraw(axis));
500   }
501 
502   PetscCall(VecGetArrayRead(u, &U));
503   U0 = PetscRealPart(U[0]);
504   U1 = PetscRealPart(U[1]);
505   PetscCall(VecRestoreArrayRead(u, &U));
506   if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);
507 
508   PetscDrawCollectiveBegin(draw);
509   PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
510   if (ictx->showtimestepandtime) {
511     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
512     PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime));
513     h = yl + .95 * (yr - yl);
514     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
515   }
516   PetscDrawCollectiveEnd(draw);
517   PetscCall(PetscDrawFlush(draw));
518   PetscCall(PetscDrawPause(draw));
519   PetscCall(PetscDrawSave(draw));
520   PetscFunctionReturn(PETSC_SUCCESS);
521 }
522 
523 /*@C
524    TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()`
525 
526    Collective
527 
528    Input Parameters:
529 .    ctx - the monitor context
530 
531    Level: intermediate
532 
533 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
534 @*/
535 PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
536 {
537   PetscFunctionBegin;
538   PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
539   PetscCall(VecDestroy(&(*ictx)->initialsolution));
540   PetscCall(PetscFree(*ictx));
541   PetscFunctionReturn(PETSC_SUCCESS);
542 }
543 
544 /*@C
545    TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`
546 
547    Collective
548 
549    Input Parameter:
550 .    ts - time-step context
551 
552    Output Parameter:
553 .    ctx - the monitor context
554 
555    Options Database Keys:
556 +   -ts_monitor_draw_solution - draw the solution at each time-step
557 -   -ts_monitor_draw_solution_initial - show initial solution as well as current solution
558 
559    Level: intermediate
560 
561    Note:
562    The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`.
563 
564 .seealso: [](chapter_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
565 @*/
566 PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
567 {
568   PetscFunctionBegin;
569   PetscCall(PetscNew(ctx));
570   PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
571   PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));
572 
573   (*ctx)->howoften    = howoften;
574   (*ctx)->showinitial = PETSC_FALSE;
575   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));
576 
577   (*ctx)->showtimestepandtime = PETSC_FALSE;
578   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
579   PetscFunctionReturn(PETSC_SUCCESS);
580 }
581 
582 /*@C
583    TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling
584    `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep
585 
586    Collective
587 
588    Input Parameters:
589 +  ts - the `TS` context
590 .  step - current time-step
591 .  ptime - current time
592 -  dummy - either a viewer or `NULL`
593 
594    Options Database Key:
595 .  -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
596 
597    Level: intermediate
598 
599    Note:
600    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
601    to be used during the `TS` integration.
602 
603 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
604 @*/
605 PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
606 {
607   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
608   PetscViewer      viewer = ctx->viewer;
609   Vec              work;
610 
611   PetscFunctionBegin;
612   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
613   PetscCall(VecDuplicate(u, &work));
614   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
615   PetscCall(VecView(work, viewer));
616   PetscCall(VecDestroy(&work));
617   PetscFunctionReturn(PETSC_SUCCESS);
618 }
619 
620 /*@C
621    TSMonitorDrawError - Monitors progress of the `TS` solvers by calling
622    `VecView()` for the error at each timestep
623 
624    Collective
625 
626    Input Parameters:
627 +  ts - the `TS` context
628 .  step - current time-step
629 .  ptime - current time
630 -  dummy - either a viewer or `NULL`
631 
632    Options Database Key:
633 .  -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
634 
635    Level: intermediate
636 
637    Notes:
638    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
639    to be used during the `TS` integration.
640 
641 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
642 @*/
643 PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
644 {
645   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
646   PetscViewer      viewer = ctx->viewer;
647   Vec              work;
648 
649   PetscFunctionBegin;
650   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
651   PetscCall(VecDuplicate(u, &work));
652   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
653   PetscCall(VecAXPY(work, -1.0, u));
654   PetscCall(VecView(work, viewer));
655   PetscCall(VecDestroy(&work));
656   PetscFunctionReturn(PETSC_SUCCESS);
657 }
658 
659 /*@C
660    TSMonitorSolution - Monitors progress of the `TS` solvers by `VecView()` for the solution at each timestep. Normally the viewer is a binary file or a `PetscDraw` object
661 
662    Collective
663 
664    Input Parameters:
665 +  ts - the `TS` context
666 .  step - current time-step
667 .  ptime - current time
668 .  u - current state
669 -  vf - viewer and its format
670 
671    Level: intermediate
672 
673    Notes:
674    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
675    to be used during the `TS` integration.
676 
677 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
678 @*/
679 PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
680 {
681   PetscFunctionBegin;
682   if (vf->view_interval > 0 && !ts->reason && step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS);
683   PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
684   PetscCall(VecView(u, vf->viewer));
685   PetscCall(PetscViewerPopFormat(vf->viewer));
686   PetscFunctionReturn(PETSC_SUCCESS);
687 }
688 
689 /*@C
690    TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at each timestep.
691 
692    Collective
693 
694    Input Parameters:
695 +  ts - the `TS` context
696 .  step - current time-step
697 .  ptime - current time
698 .  u - current state
699 -  filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")
700 
701    Level: intermediate
702 
703    Notes:
704    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.
705    These are named according to the file name template.
706 
707    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
708    to be used during the `TS` integration.
709 
710 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
711 @*/
712 PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, void *filenametemplate)
713 {
714   char        filename[PETSC_MAX_PATH_LEN];
715   PetscViewer viewer;
716 
717   PetscFunctionBegin;
718   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
719   PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)filenametemplate, step));
720   PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
721   PetscCall(VecView(u, viewer));
722   PetscCall(PetscViewerDestroy(&viewer));
723   PetscFunctionReturn(PETSC_SUCCESS);
724 }
725 
726 /*@C
727    TSMonitorSolutionVTKDestroy - Destroy filename template string created for use with `TSMonitorSolutionVTK()`
728 
729    Not Collective
730 
731    Input Parameters:
732 .  filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")
733 
734    Level: intermediate
735 
736    Note:
737    This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
738 
739 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
740 @*/
741 PetscErrorCode TSMonitorSolutionVTKDestroy(void *filenametemplate)
742 {
743   PetscFunctionBegin;
744   PetscCall(PetscFree(*(char **)filenametemplate));
745   PetscFunctionReturn(PETSC_SUCCESS);
746 }
747 
748 /*@C
749    TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector
750        in a time based line graph
751 
752    Collective
753 
754    Input Parameters:
755 +  ts - the `TS` context
756 .  step - current time-step
757 .  ptime - current time
758 .  u - current solution
759 -  dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`
760 
761    Options Database Key:
762 .   -ts_monitor_lg_solution_variables - enable monitor of lg solution variables
763 
764    Level: intermediate
765 
766    Notes:
767    Each process in a parallel run displays its component solutions in a separate window
768 
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.
771 
772 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
773           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
774           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
775           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
776 @*/
777 PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
778 {
779   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)dctx;
780   const PetscScalar *yy;
781   Vec                v;
782 
783   PetscFunctionBegin;
784   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
785   if (!step) {
786     PetscDrawAxis axis;
787     PetscInt      dim;
788     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
789     PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
790     if (!ctx->names) {
791       PetscBool flg;
792       /* user provides names of variables to plot but no names has been set so assume names are integer values */
793       PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
794       if (flg) {
795         PetscInt i, n;
796         char   **names;
797         PetscCall(VecGetSize(u, &n));
798         PetscCall(PetscMalloc1(n + 1, &names));
799         for (i = 0; i < n; i++) {
800           PetscCall(PetscMalloc1(5, &names[i]));
801           PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
802         }
803         names[n]   = NULL;
804         ctx->names = names;
805       }
806     }
807     if (ctx->names && !ctx->displaynames) {
808       char    **displaynames;
809       PetscBool flg;
810       PetscCall(VecGetLocalSize(u, &dim));
811       PetscCall(PetscCalloc1(dim + 1, &displaynames));
812       PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
813       if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
814       PetscCall(PetscStrArrayDestroy(&displaynames));
815     }
816     if (ctx->displaynames) {
817       PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
818       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
819     } else if (ctx->names) {
820       PetscCall(VecGetLocalSize(u, &dim));
821       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
822       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
823     } else {
824       PetscCall(VecGetLocalSize(u, &dim));
825       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
826     }
827     PetscCall(PetscDrawLGReset(ctx->lg));
828   }
829 
830   if (!ctx->transform) v = u;
831   else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
832   PetscCall(VecGetArrayRead(v, &yy));
833   if (ctx->displaynames) {
834     PetscInt i;
835     for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
836     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
837   } else {
838 #if defined(PETSC_USE_COMPLEX)
839     PetscInt   i, n;
840     PetscReal *yreal;
841     PetscCall(VecGetLocalSize(v, &n));
842     PetscCall(PetscMalloc1(n, &yreal));
843     for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
844     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
845     PetscCall(PetscFree(yreal));
846 #else
847     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
848 #endif
849   }
850   PetscCall(VecRestoreArrayRead(v, &yy));
851   if (ctx->transform) PetscCall(VecDestroy(&v));
852 
853   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
854     PetscCall(PetscDrawLGDraw(ctx->lg));
855     PetscCall(PetscDrawLGSave(ctx->lg));
856   }
857   PetscFunctionReturn(PETSC_SUCCESS);
858 }
859 
860 /*@C
861    TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
862 
863    Collective
864 
865    Input Parameters:
866 +  ts - the `TS` context
867 -  names - the names of the components, final string must be `NULL`
868 
869    Level: intermediate
870 
871    Notes:
872     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
873 
874 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
875 @*/
876 PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
877 {
878   PetscInt i;
879 
880   PetscFunctionBegin;
881   for (i = 0; i < ts->numbermonitors; i++) {
882     if (ts->monitor[i] == TSMonitorLGSolution) {
883       PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
884       break;
885     }
886   }
887   PetscFunctionReturn(PETSC_SUCCESS);
888 }
889 
890 /*@C
891    TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
892 
893    Collective
894 
895    Input Parameters:
896 +  ts - the `TS` context
897 -  names - the names of the components, final string must be `NULL`
898 
899    Level: intermediate
900 
901 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
902 @*/
903 PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
904 {
905   PetscFunctionBegin;
906   PetscCall(PetscStrArrayDestroy(&ctx->names));
907   PetscCall(PetscStrArrayallocpy(names, &ctx->names));
908   PetscFunctionReturn(PETSC_SUCCESS);
909 }
910 
911 /*@C
912    TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot
913 
914    Collective
915 
916    Input Parameter:
917 .  ts - the `TS` context
918 
919    Output Parameter:
920 .  names - the names of the components, final string must be `NULL`
921 
922    Level: intermediate
923 
924    Note:
925     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
926 
927 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
928 @*/
929 PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
930 {
931   PetscInt i;
932 
933   PetscFunctionBegin;
934   *names = NULL;
935   for (i = 0; i < ts->numbermonitors; i++) {
936     if (ts->monitor[i] == TSMonitorLGSolution) {
937       TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
938       *names             = (const char *const *)ctx->names;
939       break;
940     }
941   }
942   PetscFunctionReturn(PETSC_SUCCESS);
943 }
944 
945 /*@C
946    TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor
947 
948    Collective
949 
950    Input Parameters:
951 +  ctx - the `TSMonitorLG` context
952 -  displaynames - the names of the components, final string must be `NULL`
953 
954    Level: intermediate
955 
956 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
957 @*/
958 PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
959 {
960   PetscInt j = 0, k;
961 
962   PetscFunctionBegin;
963   if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
964   PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
965   PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
966   while (displaynames[j]) j++;
967   ctx->ndisplayvariables = j;
968   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
969   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
970   j = 0;
971   while (displaynames[j]) {
972     k = 0;
973     while (ctx->names[k]) {
974       PetscBool flg;
975       PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
976       if (flg) {
977         ctx->displayvariables[j] = k;
978         break;
979       }
980       k++;
981     }
982     j++;
983   }
984   PetscFunctionReturn(PETSC_SUCCESS);
985 }
986 
987 /*@C
988    TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor
989 
990    Collective
991 
992    Input Parameters:
993 +  ts - the `TS` context
994 -  displaynames - the names of the components, final string must be `NULL`
995 
996    Level: intermediate
997 
998    Note:
999     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1000 
1001 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1002 @*/
1003 PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
1004 {
1005   PetscInt i;
1006 
1007   PetscFunctionBegin;
1008   for (i = 0; i < ts->numbermonitors; i++) {
1009     if (ts->monitor[i] == TSMonitorLGSolution) {
1010       PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
1011       break;
1012     }
1013   }
1014   PetscFunctionReturn(PETSC_SUCCESS);
1015 }
1016 
1017 /*@C
1018    TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed
1019 
1020    Collective
1021 
1022    Input Parameters:
1023 +  ts - the `TS` context
1024 .  transform - the transform function
1025 .  destroy - function to destroy the optional context
1026 -  ctx - optional context used by transform function
1027 
1028    Level: intermediate
1029 
1030    Note:
1031     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1032 
1033 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()`
1034 @*/
1035 PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1036 {
1037   PetscInt i;
1038 
1039   PetscFunctionBegin;
1040   for (i = 0; i < ts->numbermonitors; i++) {
1041     if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
1042   }
1043   PetscFunctionReturn(PETSC_SUCCESS);
1044 }
1045 
1046 /*@C
1047    TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed
1048 
1049    Collective
1050 
1051    Input Parameters:
1052 +  ts - the `TS` context
1053 .  transform - the transform function
1054 .  destroy - function to destroy the optional context
1055 -  ctx - optional context used by transform function
1056 
1057    Level: intermediate
1058 
1059 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`
1060 @*/
1061 PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1062 {
1063   PetscFunctionBegin;
1064   ctx->transform        = transform;
1065   ctx->transformdestroy = destroy;
1066   ctx->transformctx     = tctx;
1067   PetscFunctionReturn(PETSC_SUCCESS);
1068 }
1069 
1070 /*@C
1071    TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error
1072        in a time based line graph
1073 
1074    Collective
1075 
1076    Input Parameters:
1077 +  ts - the `TS` context
1078 .  step - current time-step
1079 .  ptime - current time
1080 .  u - current solution
1081 -  dctx - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`
1082 
1083    Options Database Key:
1084 .  -ts_monitor_lg_error - create a graphical monitor of error history
1085 
1086    Level: intermediate
1087 
1088    Notes:
1089     Each process in a parallel run displays its component errors in a separate window
1090 
1091    The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1092 
1093    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1094    to be used during the TS integration.
1095 
1096 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1097 @*/
1098 PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
1099 {
1100   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)dummy;
1101   const PetscScalar *yy;
1102   Vec                y;
1103 
1104   PetscFunctionBegin;
1105   if (!step) {
1106     PetscDrawAxis axis;
1107     PetscInt      dim;
1108     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1109     PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
1110     PetscCall(VecGetLocalSize(u, &dim));
1111     PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1112     PetscCall(PetscDrawLGReset(ctx->lg));
1113   }
1114   PetscCall(VecDuplicate(u, &y));
1115   PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1116   PetscCall(VecAXPY(y, -1.0, u));
1117   PetscCall(VecGetArrayRead(y, &yy));
1118 #if defined(PETSC_USE_COMPLEX)
1119   {
1120     PetscReal *yreal;
1121     PetscInt   i, n;
1122     PetscCall(VecGetLocalSize(y, &n));
1123     PetscCall(PetscMalloc1(n, &yreal));
1124     for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1125     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1126     PetscCall(PetscFree(yreal));
1127   }
1128 #else
1129   PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1130 #endif
1131   PetscCall(VecRestoreArrayRead(y, &yy));
1132   PetscCall(VecDestroy(&y));
1133   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1134     PetscCall(PetscDrawLGDraw(ctx->lg));
1135     PetscCall(PetscDrawLGSave(ctx->lg));
1136   }
1137   PetscFunctionReturn(PETSC_SUCCESS);
1138 }
1139 
1140 /*@C
1141    TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot
1142 
1143    Input Parameters:
1144 +  ts - the `TS` context
1145 .  step - current time-step
1146 .  ptime - current time
1147 .  u - current solution
1148 -  dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()`
1149 
1150    Options Database Keys:
1151 + -ts_monitor_sp_swarm <n>          - Monitor the solution every n steps, or -1 for plotting only the final solution
1152 . -ts_monitor_sp_swarm_retain <n>   - Retain n old points so we can see the history, or -1 for all points
1153 - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space
1154 
1155    Level: intermediate
1156 
1157    Notes:
1158    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1159    to be used during the `TS` integration.
1160 
1161 .seealso: [](chapter_ts), `TS`, `TSMonitoSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
1162 @*/
1163 PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1164 {
1165   TSMonitorSPCtx     ctx = (TSMonitorSPCtx)dctx;
1166   PetscDraw          draw;
1167   DM                 dm, cdm;
1168   const PetscScalar *yy;
1169   PetscInt           Np, p, dim = 2;
1170 
1171   PetscFunctionBegin;
1172   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1173   if (!step) {
1174     PetscDrawAxis axis;
1175     PetscReal     dmboxlower[2], dmboxupper[2];
1176 
1177     PetscCall(TSGetDM(ts, &dm));
1178     PetscCall(DMGetDimension(dm, &dim));
1179     PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
1180     PetscCall(DMSwarmGetCellDM(dm, &cdm));
1181     PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
1182     PetscCall(VecGetLocalSize(u, &Np));
1183     Np /= dim * 2;
1184     PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
1185     if (ctx->phase) {
1186       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
1187       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -5, 5));
1188     } else {
1189       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
1190       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
1191     }
1192     PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
1193     PetscCall(PetscDrawSPReset(ctx->sp));
1194   }
1195   PetscCall(VecGetLocalSize(u, &Np));
1196   Np /= dim * 2;
1197   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1198     PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
1199     if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
1200     PetscCall(PetscDrawFlush(draw));
1201     PetscCall(PetscDrawSPReset(ctx->sp));
1202     PetscCall(VecGetArrayRead(u, &yy));
1203     for (p = 0; p < Np; ++p) {
1204       PetscReal x, y;
1205 
1206       if (ctx->phase) {
1207         x = PetscRealPart(yy[p * dim * 2]);
1208         y = PetscRealPart(yy[p * dim * 2 + dim]);
1209       } else {
1210         x = PetscRealPart(yy[p * dim * 2]);
1211         y = PetscRealPart(yy[p * dim * 2 + 1]);
1212       }
1213       PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1214     }
1215     PetscCall(VecRestoreArrayRead(u, &yy));
1216     PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
1217     PetscCall(PetscDrawSPSave(ctx->sp));
1218   }
1219   PetscFunctionReturn(PETSC_SUCCESS);
1220 }
1221 
1222 /*@C
1223    TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep
1224 
1225    Collective
1226 
1227    Input Parameters:
1228 +  ts - the `TS` context
1229 .  step - current time-step
1230 .  ptime - current time
1231 .  u - current solution
1232 -  dctx - unused context
1233 
1234    Options Database Key:
1235 .  -ts_monitor_error - create a graphical monitor of error history
1236 
1237    Level: intermediate
1238 
1239    Notes:
1240    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1241    to be used during the `TS` integration.
1242 
1243    The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1244 
1245 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1246 @*/
1247 PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1248 {
1249   DM        dm;
1250   PetscDS   ds = NULL;
1251   PetscInt  Nf = -1, f;
1252   PetscBool flg;
1253 
1254   PetscFunctionBegin;
1255   PetscCall(TSGetDM(ts, &dm));
1256   if (dm) PetscCall(DMGetDS(dm, &ds));
1257   if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
1258   if (Nf <= 0) {
1259     Vec       y;
1260     PetscReal nrm;
1261 
1262     PetscCall(VecDuplicate(u, &y));
1263     PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1264     PetscCall(VecAXPY(y, -1.0, u));
1265     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
1266     if (flg) {
1267       PetscCall(VecNorm(y, NORM_2, &nrm));
1268       PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
1269     }
1270     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
1271     if (flg) PetscCall(VecView(y, vf->viewer));
1272     PetscCall(VecDestroy(&y));
1273   } else {
1274     PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
1275     void    **ctxs;
1276     Vec       v;
1277     PetscReal ferrors[1];
1278 
1279     PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
1280     for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
1281     PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
1282     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04d time = %-8.4g \t L_2 Error: [", (int)step, (double)ptime));
1283     for (f = 0; f < Nf; ++f) {
1284       if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
1285       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
1286     }
1287     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));
1288 
1289     PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));
1290 
1291     PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
1292     if (flg) {
1293       PetscCall(DMGetGlobalVector(dm, &v));
1294       PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
1295       PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
1296       PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
1297       PetscCall(DMRestoreGlobalVector(dm, &v));
1298     }
1299     PetscCall(PetscFree2(exactFuncs, ctxs));
1300   }
1301   PetscFunctionReturn(PETSC_SUCCESS);
1302 }
1303 
1304 PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1305 {
1306   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1307   PetscReal      x   = ptime, y;
1308   PetscInt       its;
1309 
1310   PetscFunctionBegin;
1311   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1312   if (!n) {
1313     PetscDrawAxis axis;
1314     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1315     PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
1316     PetscCall(PetscDrawLGReset(ctx->lg));
1317     ctx->snes_its = 0;
1318   }
1319   PetscCall(TSGetSNESIterations(ts, &its));
1320   y = its - ctx->snes_its;
1321   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1322   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1323     PetscCall(PetscDrawLGDraw(ctx->lg));
1324     PetscCall(PetscDrawLGSave(ctx->lg));
1325   }
1326   ctx->snes_its = its;
1327   PetscFunctionReturn(PETSC_SUCCESS);
1328 }
1329 
1330 PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1331 {
1332   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1333   PetscReal      x   = ptime, y;
1334   PetscInt       its;
1335 
1336   PetscFunctionBegin;
1337   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1338   if (!n) {
1339     PetscDrawAxis axis;
1340     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1341     PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
1342     PetscCall(PetscDrawLGReset(ctx->lg));
1343     ctx->ksp_its = 0;
1344   }
1345   PetscCall(TSGetKSPIterations(ts, &its));
1346   y = its - ctx->ksp_its;
1347   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1348   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1349     PetscCall(PetscDrawLGDraw(ctx->lg));
1350     PetscCall(PetscDrawLGSave(ctx->lg));
1351   }
1352   ctx->ksp_its = its;
1353   PetscFunctionReturn(PETSC_SUCCESS);
1354 }
1355 
1356 /*@C
1357    TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()`
1358 
1359    Collective
1360 
1361    Input Parameters:
1362 .  ts  - the `TS` solver object
1363 
1364    Output Parameter:
1365 .  ctx - the context
1366 
1367    Level: intermediate
1368 
1369 .seealso: [](chapter_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1370 @*/
1371 PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1372 {
1373   PetscFunctionBegin;
1374   PetscCall(PetscNew(ctx));
1375   PetscFunctionReturn(PETSC_SUCCESS);
1376 }
1377 
1378 /*@C
1379    TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution
1380 
1381    Collective
1382 
1383    Input Parameters:
1384 +  ts - the `TS` context
1385 .  step - current time-step
1386 .  ptime - current time
1387 .  u  - current solution
1388 -  dctx - the envelope context
1389 
1390    Options Database Key:
1391 .  -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time
1392 
1393    Level: intermediate
1394 
1395    Notes:
1396    After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope
1397 
1398    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1399    to be used during the `TS` integration.
1400 
1401 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1402 @*/
1403 PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1404 {
1405   TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;
1406 
1407   PetscFunctionBegin;
1408   if (!ctx->max) {
1409     PetscCall(VecDuplicate(u, &ctx->max));
1410     PetscCall(VecDuplicate(u, &ctx->min));
1411     PetscCall(VecCopy(u, ctx->max));
1412     PetscCall(VecCopy(u, ctx->min));
1413   } else {
1414     PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
1415     PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
1416   }
1417   PetscFunctionReturn(PETSC_SUCCESS);
1418 }
1419 
1420 /*@C
1421    TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution
1422 
1423    Collective
1424 
1425    Input Parameter:
1426 .  ts - the `TS` context
1427 
1428    Output Parameters:
1429 +  max - the maximum values
1430 -  min - the minimum values
1431 
1432    Level: intermediate
1433 
1434    Notes:
1435     If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored
1436 
1437 .seealso: [](chapter_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1438 @*/
1439 PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1440 {
1441   PetscInt i;
1442 
1443   PetscFunctionBegin;
1444   if (max) *max = NULL;
1445   if (min) *min = NULL;
1446   for (i = 0; i < ts->numbermonitors; i++) {
1447     if (ts->monitor[i] == TSMonitorEnvelope) {
1448       TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1449       if (max) *max = ctx->max;
1450       if (min) *min = ctx->min;
1451       break;
1452     }
1453   }
1454   PetscFunctionReturn(PETSC_SUCCESS);
1455 }
1456 
1457 /*@C
1458    TSMonitorEnvelopeCtxDestroy - Destroys a context that was created  with `TSMonitorEnvelopeCtxCreate()`.
1459 
1460    Collective
1461 
1462    Input Parameter:
1463 .  ctx - the monitor context
1464 
1465    Level: intermediate
1466 
1467 .seealso: [](chapter_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1468 @*/
1469 PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1470 {
1471   PetscFunctionBegin;
1472   PetscCall(VecDestroy(&(*ctx)->min));
1473   PetscCall(VecDestroy(&(*ctx)->max));
1474   PetscCall(PetscFree(*ctx));
1475   PetscFunctionReturn(PETSC_SUCCESS);
1476 }
1477 
1478 /*@C
1479   TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS`
1480 
1481   Not collective
1482 
1483   Input Parameters:
1484 + ts   - the `TS` context
1485 . step - current timestep
1486 . t    - current time
1487 . u    - current solution
1488 - ctx  - not used
1489 
1490   Options Database Key:
1491 . -ts_dmswarm_monitor_moments - Monitor moments of particle distribution
1492 
1493   Level: intermediate
1494 
1495   Notes:
1496   This requires a `DMSWARM` be attached to the `TS`.
1497 
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.
1500 
1501 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1502 @*/
1503 PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1504 {
1505   DM                 sw;
1506   const PetscScalar *u;
1507   PetscReal          m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1508   PetscInt           dim, d, Np, p;
1509   MPI_Comm           comm;
1510 
1511   PetscFunctionBeginUser;
1512   PetscCall(TSGetDM(ts, &sw));
1513   if (!sw || step % ts->monitorFrequency != 0) PetscFunctionReturn(PETSC_SUCCESS);
1514   PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
1515   PetscCall(DMGetDimension(sw, &dim));
1516   PetscCall(VecGetLocalSize(U, &Np));
1517   Np /= dim;
1518   PetscCall(VecGetArrayRead(U, &u));
1519   for (p = 0; p < Np; ++p) {
1520     for (d = 0; d < dim; ++d) {
1521       totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1522       totMom[d] += PetscRealPart(u[p * dim + d]);
1523     }
1524   }
1525   PetscCall(VecRestoreArrayRead(U, &u));
1526   for (d = 0; d < dim; ++d) totMom[d] *= m;
1527   totE *= 0.5 * m;
1528   PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
1529   for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, "    Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
1530   PetscCall(PetscPrintf(comm, "\n"));
1531   PetscFunctionReturn(PETSC_SUCCESS);
1532 }
1533