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