Actual source code: tsmon.c
1: #include <petsc/private/tsimpl.h>
2: #include <petscdm.h>
3: #include <petscds.h>
4: #include <petscdmswarm.h>
5: #include <petscdraw.h>
7: /*@C
8: TSMonitor - Runs all user-provided monitor routines set using `TSMonitorSet()`
10: Collective
12: Input Parameters:
13: + ts - time stepping context obtained from `TSCreate()`
14: . step - step number that has just completed
15: . ptime - model time of the state
16: - u - state at the current model time
18: Level: developer
20: Notes:
21: `TSMonitor()` is typically used automatically within the time stepping implementations.
22: Users would almost never call this routine directly.
24: A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions
26: .seealso: `TS`, `TSMonitorSet()`, `TSMonitorSetFromOptions()`
27: @*/
28: PetscErrorCode TSMonitor(TS ts, PetscInt step, PetscReal ptime, Vec u)
29: {
30: DM dm;
31: PetscInt i, n = ts->numbermonitors;
33: PetscFunctionBegin;
37: PetscCall(TSGetDM(ts, &dm));
38: PetscCall(DMSetOutputSequenceNumber(dm, step, ptime));
40: PetscCall(VecLockReadPush(u));
41: for (i = 0; i < n; i++) PetscCall((*ts->monitor[i])(ts, step, ptime, u, ts->monitorcontext[i]));
42: PetscCall(VecLockReadPop(u));
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: /*@C
47: TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
49: Collective
51: Input Parameters:
52: + ts - `TS` object you wish to monitor
53: . name - the monitor type one is seeking
54: . help - message indicating what monitoring is done
55: . manual - manual page for the monitor
56: . monitor - the monitor function
57: - monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the `TS` or `PetscViewer` objects
59: Level: developer
61: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `PetscOptionsGetViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
62: `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
63: `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`,
64: `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
65: `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
66: `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
67: `PetscOptionsFList()`, `PetscOptionsEList()`
68: @*/
69: PetscErrorCode TSMonitorSetFromOptions(TS ts, const char name[], const char help[], const char manual[], PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *), PetscErrorCode (*monitorsetup)(TS, PetscViewerAndFormat *))
70: {
71: PetscViewer viewer;
72: PetscViewerFormat format;
73: PetscBool flg;
75: PetscFunctionBegin;
76: PetscCall(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: }
90: /*@C
91: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
92: timestep to display the iteration's progress.
94: Logically Collective
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: - mdestroy - [optional] routine that frees monitor context (may be `NULL`)
102: Calling sequence of `monitor`:
103: + ts - the `TS` context
104: . 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)
105: . time - current time
106: . u - current iterate
107: - ctx - [optional] monitoring context
109: Level: intermediate
111: Note:
112: This routine adds an additional monitor to the list of monitors that already has been loaded.
114: Fortran Notes:
115: Only a single monitor function can be set for each `TS` object
117: .seealso: [](ch_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
118: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
119: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`
120: @*/
121: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS ts, PetscInt steps, PetscReal time, Vec u, void *ctx), void *mctx, PetscErrorCode (*mdestroy)(void **))
122: {
123: PetscInt i;
124: PetscBool identical;
126: PetscFunctionBegin;
128: for (i = 0; i < ts->numbermonitors; i++) {
129: PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))monitor, mctx, mdestroy, (PetscErrorCode(*)(void))ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical));
130: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
131: }
132: PetscCheck(ts->numbermonitors < MAXTSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set");
133: ts->monitor[ts->numbermonitors] = monitor;
134: ts->monitordestroy[ts->numbermonitors] = mdestroy;
135: ts->monitorcontext[ts->numbermonitors++] = (void *)mctx;
136: PetscFunctionReturn(PETSC_SUCCESS);
137: }
139: /*@C
140: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
142: Logically Collective
144: Input Parameter:
145: . ts - the `TS` context obtained from `TSCreate()`
147: Level: intermediate
149: Note:
150: There is no way to remove a single, specific monitor.
152: .seealso: [](ch_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
153: @*/
154: PetscErrorCode TSMonitorCancel(TS ts)
155: {
156: PetscInt i;
158: PetscFunctionBegin;
160: for (i = 0; i < ts->numbermonitors; i++) {
161: if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i]));
162: }
163: ts->numbermonitors = 0;
164: PetscFunctionReturn(PETSC_SUCCESS);
165: }
167: /*@C
168: TSMonitorDefault - The default monitor, prints the timestep and time for each step
170: Input Parameters:
171: + ts - the `TS` context
172: . 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)
173: . ptime - current time
174: . v - current iterate
175: - vf - the viewer and format
177: Options Database Key:
178: . -ts_monitor - monitors the time integration
180: Level: intermediate
182: Notes:
183: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
184: to be used during the `TS` integration.
186: .seealso: [](ch_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
187: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
188: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`
189: @*/
190: PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
191: {
192: PetscViewer viewer = vf->viewer;
193: PetscBool iascii, ibinary;
195: PetscFunctionBegin;
197: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
198: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary));
199: PetscCall(PetscViewerPushFormat(viewer, vf->format));
200: if (iascii) {
201: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
202: if (step == -1) { /* this indicates it is an interpolated solution */
203: PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps));
204: } else {
205: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n"));
206: }
207: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
208: } else if (ibinary) {
209: PetscMPIInt rank;
210: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
211: if (rank == 0) {
212: PetscBool skipHeader;
213: PetscInt classid = REAL_FILE_CLASSID;
215: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
216: if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
217: PetscCall(PetscRealView(1, &ptime, viewer));
218: } else {
219: PetscCall(PetscRealView(0, &ptime, viewer));
220: }
221: }
222: PetscCall(PetscViewerPopFormat(viewer));
223: PetscFunctionReturn(PETSC_SUCCESS);
224: }
226: /*@C
227: TSMonitorExtreme - Prints the extreme values of the solution at each timestep
229: Input Parameters:
230: + ts - the `TS` context
231: . 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)
232: . ptime - current time
233: . v - current iterate
234: - vf - the viewer and format
236: Level: intermediate
238: Note:
239: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
240: to be used during the `TS` integration.
242: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`
243: @*/
244: PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
245: {
246: PetscViewer viewer = vf->viewer;
247: PetscBool iascii;
248: PetscReal max, min;
250: PetscFunctionBegin;
252: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
253: PetscCall(PetscViewerPushFormat(viewer, vf->format));
254: if (iascii) {
255: PetscCall(VecMax(v, NULL, &max));
256: PetscCall(VecMin(v, NULL, &min));
257: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
258: 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));
259: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
260: }
261: PetscCall(PetscViewerPopFormat(viewer));
262: PetscFunctionReturn(PETSC_SUCCESS);
263: }
265: /*@C
266: TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with
267: `TS` to monitor the solution process graphically in various ways
269: Collective
271: Input Parameters:
272: + comm - the MPI communicator to use
273: . host - the X display to open, or `NULL` for the local machine
274: . label - the title to put in the title bar
275: . x - the x screen coordinates of the upper left coordinate of the window
276: . y - the y screen coordinates of the upper left coordinate of the window
277: . m - the screen width in pixels
278: . n - the screen height in pixels
279: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
281: Output Parameter:
282: . ctx - the context
284: Options Database Keys:
285: + -ts_monitor_lg_timestep - automatically sets line graph monitor
286: . -ts_monitor_lg_timestep_log - automatically sets line graph monitor
287: . -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`)
288: . -ts_monitor_lg_error - monitor the error
289: . -ts_monitor_lg_ksp_iterations - monitor the number of `KSP` iterations needed for each timestep
290: . -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep
291: - -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true
293: Level: intermediate
295: Notes:
296: Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed.
298: One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()`
300: Many of the functions that control the monitoring have two forms\: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a `TS` object as the
301: 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
302: as the first argument.
304: One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()`
306: .seealso: [](ch_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
307: `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
308: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
309: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
310: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
311: @*/
312: PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
313: {
314: PetscDraw draw;
316: PetscFunctionBegin;
317: PetscCall(PetscNew(ctx));
318: PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
319: PetscCall(PetscDrawSetFromOptions(draw));
320: PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg));
321: PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg));
322: PetscCall(PetscDrawDestroy(&draw));
323: (*ctx)->howoften = howoften;
324: PetscFunctionReturn(PETSC_SUCCESS);
325: }
327: PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx)
328: {
329: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
330: PetscReal x = ptime, y;
332: PetscFunctionBegin;
333: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */
334: if (!step) {
335: PetscDrawAxis axis;
336: const char *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
337: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
338: PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel));
339: PetscCall(PetscDrawLGReset(ctx->lg));
340: }
341: PetscCall(TSGetTimeStep(ts, &y));
342: if (ctx->semilogy) y = PetscLog10Real(y);
343: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
344: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
345: PetscCall(PetscDrawLGDraw(ctx->lg));
346: PetscCall(PetscDrawLGSave(ctx->lg));
347: }
348: PetscFunctionReturn(PETSC_SUCCESS);
349: }
351: /*@C
352: TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`.
354: Collective
356: Input Parameter:
357: . ctx - the monitor context
359: Level: intermediate
361: Note:
362: Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`
364: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep();`
365: @*/
366: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
367: {
368: PetscFunctionBegin;
369: if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)((*ctx)->transformctx));
370: PetscCall(PetscDrawLGDestroy(&(*ctx)->lg));
371: PetscCall(PetscStrArrayDestroy(&(*ctx)->names));
372: PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames));
373: PetscCall(PetscFree((*ctx)->displayvariables));
374: PetscCall(PetscFree((*ctx)->displayvalues));
375: PetscCall(PetscFree(*ctx));
376: PetscFunctionReturn(PETSC_SUCCESS);
377: }
379: /* Creates a TSMonitorSPCtx for use with DMSwarm particle visualizations */
380: 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)
381: {
382: PetscDraw draw;
384: PetscFunctionBegin;
385: PetscCall(PetscNew(ctx));
386: PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
387: PetscCall(PetscDrawSetFromOptions(draw));
388: PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp));
389: PetscCall(PetscDrawDestroy(&draw));
390: (*ctx)->howoften = howoften;
391: (*ctx)->retain = retain;
392: (*ctx)->phase = phase;
393: (*ctx)->multispecies = multispecies;
394: PetscFunctionReturn(PETSC_SUCCESS);
395: }
397: /* Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate */
398: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
399: {
400: PetscFunctionBegin;
402: PetscCall(PetscDrawSPDestroy(&(*ctx)->sp));
403: PetscCall(PetscFree(*ctx));
404: PetscFunctionReturn(PETSC_SUCCESS);
405: }
407: /* Creates a TSMonitorHGCtx for use with DMSwarm particle visualizations */
408: 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)
409: {
410: PetscDraw draw;
411: PetscInt s;
413: PetscFunctionBegin;
414: PetscCall(PetscNew(ctx));
415: PetscCall(PetscMalloc1(Ns, &(*ctx)->hg));
416: for (s = 0; s < Ns; ++s) {
417: PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw));
418: PetscCall(PetscDrawSetFromOptions(draw));
419: PetscCall(PetscDrawHGCreate(draw, Nb, &(*ctx)->hg[s]));
420: PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE));
421: PetscCall(PetscDrawDestroy(&draw));
422: }
423: (*ctx)->howoften = howoften;
424: (*ctx)->Ns = Ns;
425: (*ctx)->velocity = velocity;
426: PetscFunctionReturn(PETSC_SUCCESS);
427: }
429: /* Destroys a TSMonitorHGCtx that was created with TSMonitorHGCtxCreate */
430: PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx)
431: {
432: PetscInt s;
434: PetscFunctionBegin;
435: for (s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s]));
436: PetscCall(PetscFree((*ctx)->hg));
437: PetscCall(PetscFree(*ctx));
438: PetscFunctionReturn(PETSC_SUCCESS);
439: }
441: /*@C
442: TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling
443: `VecView()` for the solution at each timestep
445: Collective
447: Input Parameters:
448: + ts - the `TS` context
449: . step - current time-step
450: . ptime - current time
451: . u - the solution at the current time
452: - dummy - either a viewer or `NULL`
454: Options Database Keys:
455: + -ts_monitor_draw_solution - draw the solution at each time-step
456: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
458: Level: intermediate
460: Notes:
461: The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial`
462: will look bad
464: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with
465: `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration.
467: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
468: @*/
469: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
470: {
471: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
472: PetscDraw draw;
474: PetscFunctionBegin;
475: if (!step && ictx->showinitial) {
476: if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
477: PetscCall(VecCopy(u, ictx->initialsolution));
478: }
479: if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
481: if (ictx->showinitial) {
482: PetscReal pause;
483: PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
484: PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
485: PetscCall(VecView(ictx->initialsolution, ictx->viewer));
486: PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
487: PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
488: }
489: PetscCall(VecView(u, ictx->viewer));
490: if (ictx->showtimestepandtime) {
491: PetscReal xl, yl, xr, yr, h;
492: char time[32];
494: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
495: PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime));
496: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
497: h = yl + .95 * (yr - yl);
498: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
499: PetscCall(PetscDrawFlush(draw));
500: }
502: if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
503: PetscFunctionReturn(PETSC_SUCCESS);
504: }
506: /*@C
507: TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram
509: Collective
511: Input Parameters:
512: + ts - the `TS` context
513: . step - current time-step
514: . ptime - current time
515: . u - the solution at the current time
516: - dummy - either a viewer or `NULL`
518: Level: intermediate
520: Notes:
521: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
522: to be used during the `TS` integration.
524: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
525: @*/
526: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
527: {
528: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
529: PetscDraw draw;
530: PetscDrawAxis axis;
531: PetscInt n;
532: PetscMPIInt size;
533: PetscReal U0, U1, xl, yl, xr, yr, h;
534: char time[32];
535: const PetscScalar *U;
537: PetscFunctionBegin;
538: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
539: PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
540: PetscCall(VecGetSize(u, &n));
541: PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");
543: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
544: PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
545: PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
546: if (!step) {
547: PetscCall(PetscDrawClear(draw));
548: PetscCall(PetscDrawAxisDraw(axis));
549: }
551: PetscCall(VecGetArrayRead(u, &U));
552: U0 = PetscRealPart(U[0]);
553: U1 = PetscRealPart(U[1]);
554: PetscCall(VecRestoreArrayRead(u, &U));
555: if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);
557: PetscDrawCollectiveBegin(draw);
558: PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
559: if (ictx->showtimestepandtime) {
560: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
561: PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime));
562: h = yl + .95 * (yr - yl);
563: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
564: }
565: PetscDrawCollectiveEnd(draw);
566: PetscCall(PetscDrawFlush(draw));
567: PetscCall(PetscDrawPause(draw));
568: PetscCall(PetscDrawSave(draw));
569: PetscFunctionReturn(PETSC_SUCCESS);
570: }
572: /*@C
573: TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()`
575: Collective
577: Input Parameter:
578: . ictx - the monitor context
580: Level: intermediate
582: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
583: @*/
584: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
585: {
586: PetscFunctionBegin;
587: PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
588: PetscCall(VecDestroy(&(*ictx)->initialsolution));
589: PetscCall(PetscFree(*ictx));
590: PetscFunctionReturn(PETSC_SUCCESS);
591: }
593: /*@C
594: TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`
596: Collective
598: Input Parameters:
599: + comm - the MPI communicator to use
600: . host - the X display to open, or `NULL` for the local machine
601: . label - the title to put in the title bar
602: . x - the x screen coordinates of the upper left coordinate of the window
603: . y - the y screen coordinates of the upper left coordinate of the window
604: . m - the screen width in pixels
605: . n - the screen height in pixels
606: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
608: Output Parameter:
609: . ctx - the monitor context
611: Options Database Keys:
612: + -ts_monitor_draw_solution - draw the solution at each time-step
613: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
615: Level: intermediate
617: Note:
618: The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`.
620: .seealso: [](ch_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
621: @*/
622: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
623: {
624: PetscFunctionBegin;
625: PetscCall(PetscNew(ctx));
626: PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
627: PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));
629: (*ctx)->howoften = howoften;
630: (*ctx)->showinitial = PETSC_FALSE;
631: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));
633: (*ctx)->showtimestepandtime = PETSC_FALSE;
634: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
635: PetscFunctionReturn(PETSC_SUCCESS);
636: }
638: /*@C
639: TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling
640: `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep
642: Collective
644: Input Parameters:
645: + ts - the `TS` context
646: . step - current time-step
647: . ptime - current time
648: . u - solution at current time
649: - dummy - either a viewer or `NULL`
651: Options Database Key:
652: . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
654: Level: intermediate
656: Note:
657: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
658: to be used during the `TS` integration.
660: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
661: @*/
662: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
663: {
664: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
665: PetscViewer viewer = ctx->viewer;
666: Vec work;
668: PetscFunctionBegin;
669: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
670: PetscCall(VecDuplicate(u, &work));
671: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
672: PetscCall(VecView(work, viewer));
673: PetscCall(VecDestroy(&work));
674: PetscFunctionReturn(PETSC_SUCCESS);
675: }
677: /*@C
678: TSMonitorDrawError - Monitors progress of the `TS` solvers by calling
679: `VecView()` for the error at each timestep
681: Collective
683: Input Parameters:
684: + ts - the `TS` context
685: . step - current time-step
686: . ptime - current time
687: . u - solution at current time
688: - dummy - either a viewer or `NULL`
690: Options Database Key:
691: . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
693: Level: intermediate
695: Notes:
696: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
697: to be used during the `TS` integration.
699: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
700: @*/
701: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
702: {
703: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
704: PetscViewer viewer = ctx->viewer;
705: Vec work;
707: PetscFunctionBegin;
708: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
709: PetscCall(VecDuplicate(u, &work));
710: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
711: PetscCall(VecAXPY(work, -1.0, u));
712: PetscCall(VecView(work, viewer));
713: PetscCall(VecDestroy(&work));
714: PetscFunctionReturn(PETSC_SUCCESS);
715: }
717: /*@C
718: 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
720: Collective
722: Input Parameters:
723: + ts - the `TS` context
724: . step - current time-step
725: . ptime - current time
726: . u - current state
727: - vf - viewer and its format
729: Level: intermediate
731: Notes:
732: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
733: to be used during the `TS` integration.
735: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
736: @*/
737: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
738: {
739: PetscFunctionBegin;
740: if (vf->view_interval > 0 && !ts->reason && step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS);
741: PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
742: PetscCall(VecView(u, vf->viewer));
743: PetscCall(PetscViewerPopFormat(vf->viewer));
744: PetscFunctionReturn(PETSC_SUCCESS);
745: }
747: /*@C
748: TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at each timestep.
750: Collective
752: Input Parameters:
753: + ts - the `TS` context
754: . step - current time-step
755: . ptime - current time
756: . u - current state
757: - filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")
759: Level: intermediate
761: Notes:
762: 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.
763: These are named according to the file name template.
765: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
766: to be used during the `TS` integration.
768: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
769: @*/
770: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, void *filenametemplate)
771: {
772: char filename[PETSC_MAX_PATH_LEN];
773: PetscViewer viewer;
775: PetscFunctionBegin;
776: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
777: PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)filenametemplate, step));
778: PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
779: PetscCall(VecView(u, viewer));
780: PetscCall(PetscViewerDestroy(&viewer));
781: PetscFunctionReturn(PETSC_SUCCESS);
782: }
784: /*@C
785: TSMonitorSolutionVTKDestroy - Destroy filename template string created for use with `TSMonitorSolutionVTK()`
787: Not Collective
789: Input Parameter:
790: . filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")
792: Level: intermediate
794: Note:
795: This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
797: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
798: @*/
799: PetscErrorCode TSMonitorSolutionVTKDestroy(void *filenametemplate)
800: {
801: PetscFunctionBegin;
802: PetscCall(PetscFree(*(char **)filenametemplate));
803: PetscFunctionReturn(PETSC_SUCCESS);
804: }
806: /*@C
807: TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector
808: in a time based line graph
810: Collective
812: Input Parameters:
813: + ts - the `TS` context
814: . step - current time-step
815: . ptime - current time
816: . u - current solution
817: - dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`
819: Options Database Key:
820: . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables
822: Level: intermediate
824: Notes:
825: Each process in a parallel run displays its component solutions in a separate window
827: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
828: to be used during the `TS` integration.
830: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
831: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
832: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
833: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
834: @*/
835: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
836: {
837: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dctx;
838: const PetscScalar *yy;
839: Vec v;
841: PetscFunctionBegin;
842: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
843: if (!step) {
844: PetscDrawAxis axis;
845: PetscInt dim;
846: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
847: PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
848: if (!ctx->names) {
849: PetscBool flg;
850: /* user provides names of variables to plot but no names has been set so assume names are integer values */
851: PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
852: if (flg) {
853: PetscInt i, n;
854: char **names;
855: PetscCall(VecGetSize(u, &n));
856: PetscCall(PetscMalloc1(n + 1, &names));
857: for (i = 0; i < n; i++) {
858: PetscCall(PetscMalloc1(5, &names[i]));
859: PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
860: }
861: names[n] = NULL;
862: ctx->names = names;
863: }
864: }
865: if (ctx->names && !ctx->displaynames) {
866: char **displaynames;
867: PetscBool flg;
868: PetscCall(VecGetLocalSize(u, &dim));
869: PetscCall(PetscCalloc1(dim + 1, &displaynames));
870: PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
871: if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
872: PetscCall(PetscStrArrayDestroy(&displaynames));
873: }
874: if (ctx->displaynames) {
875: PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
876: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
877: } else if (ctx->names) {
878: PetscCall(VecGetLocalSize(u, &dim));
879: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
880: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
881: } else {
882: PetscCall(VecGetLocalSize(u, &dim));
883: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
884: }
885: PetscCall(PetscDrawLGReset(ctx->lg));
886: }
888: if (!ctx->transform) v = u;
889: else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
890: PetscCall(VecGetArrayRead(v, &yy));
891: if (ctx->displaynames) {
892: PetscInt i;
893: for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
894: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
895: } else {
896: #if defined(PETSC_USE_COMPLEX)
897: PetscInt i, n;
898: PetscReal *yreal;
899: PetscCall(VecGetLocalSize(v, &n));
900: PetscCall(PetscMalloc1(n, &yreal));
901: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
902: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
903: PetscCall(PetscFree(yreal));
904: #else
905: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
906: #endif
907: }
908: PetscCall(VecRestoreArrayRead(v, &yy));
909: if (ctx->transform) PetscCall(VecDestroy(&v));
911: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
912: PetscCall(PetscDrawLGDraw(ctx->lg));
913: PetscCall(PetscDrawLGSave(ctx->lg));
914: }
915: PetscFunctionReturn(PETSC_SUCCESS);
916: }
918: /*@C
919: TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
921: Collective
923: Input Parameters:
924: + ts - the `TS` context
925: - names - the names of the components, final string must be `NULL`
927: Level: intermediate
929: Notes:
930: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
932: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
933: @*/
934: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
935: {
936: PetscInt i;
938: PetscFunctionBegin;
939: for (i = 0; i < ts->numbermonitors; i++) {
940: if (ts->monitor[i] == TSMonitorLGSolution) {
941: PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
942: break;
943: }
944: }
945: PetscFunctionReturn(PETSC_SUCCESS);
946: }
948: /*@C
949: TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
951: Collective
953: Input Parameters:
954: + ctx - the `TS` context
955: - names - the names of the components, final string must be `NULL`
957: Level: intermediate
959: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
960: @*/
961: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
962: {
963: PetscFunctionBegin;
964: PetscCall(PetscStrArrayDestroy(&ctx->names));
965: PetscCall(PetscStrArrayallocpy(names, &ctx->names));
966: PetscFunctionReturn(PETSC_SUCCESS);
967: }
969: /*@C
970: TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot
972: Collective
974: Input Parameter:
975: . ts - the `TS` context
977: Output Parameter:
978: . names - the names of the components, final string must be `NULL`
980: Level: intermediate
982: Note:
983: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
985: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
986: @*/
987: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
988: {
989: PetscInt i;
991: PetscFunctionBegin;
992: *names = NULL;
993: for (i = 0; i < ts->numbermonitors; i++) {
994: if (ts->monitor[i] == TSMonitorLGSolution) {
995: TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
996: *names = (const char *const *)ctx->names;
997: break;
998: }
999: }
1000: PetscFunctionReturn(PETSC_SUCCESS);
1001: }
1003: /*@C
1004: TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor
1006: Collective
1008: Input Parameters:
1009: + ctx - the `TSMonitorLG` context
1010: - displaynames - the names of the components, final string must be `NULL`
1012: Level: intermediate
1014: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1015: @*/
1016: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
1017: {
1018: PetscInt j = 0, k;
1020: PetscFunctionBegin;
1021: if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
1022: PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
1023: PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
1024: while (displaynames[j]) j++;
1025: ctx->ndisplayvariables = j;
1026: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
1027: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
1028: j = 0;
1029: while (displaynames[j]) {
1030: k = 0;
1031: while (ctx->names[k]) {
1032: PetscBool flg;
1033: PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
1034: if (flg) {
1035: ctx->displayvariables[j] = k;
1036: break;
1037: }
1038: k++;
1039: }
1040: j++;
1041: }
1042: PetscFunctionReturn(PETSC_SUCCESS);
1043: }
1045: /*@C
1046: TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor
1048: Collective
1050: Input Parameters:
1051: + ts - the `TS` context
1052: - displaynames - the names of the components, final string must be `NULL`
1054: Level: intermediate
1056: Note:
1057: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1059: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1060: @*/
1061: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
1062: {
1063: PetscInt i;
1065: PetscFunctionBegin;
1066: for (i = 0; i < ts->numbermonitors; i++) {
1067: if (ts->monitor[i] == TSMonitorLGSolution) {
1068: PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
1069: break;
1070: }
1071: }
1072: PetscFunctionReturn(PETSC_SUCCESS);
1073: }
1075: /*@C
1076: TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed
1078: Collective
1080: Input Parameters:
1081: + ts - the `TS` context
1082: . transform - the transform function
1083: . destroy - function to destroy the optional context
1084: - tctx - optional context used by transform function
1086: Level: intermediate
1088: Note:
1089: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1091: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()`
1092: @*/
1093: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1094: {
1095: PetscInt i;
1097: PetscFunctionBegin;
1098: for (i = 0; i < ts->numbermonitors; i++) {
1099: if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
1100: }
1101: PetscFunctionReturn(PETSC_SUCCESS);
1102: }
1104: /*@C
1105: TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed
1107: Collective
1109: Input Parameters:
1110: + tctx - the `TS` context
1111: . transform - the transform function
1112: . destroy - function to destroy the optional context
1113: - ctx - optional context used by transform function
1115: Level: intermediate
1117: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`
1118: @*/
1119: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1120: {
1121: PetscFunctionBegin;
1122: ctx->transform = transform;
1123: ctx->transformdestroy = destroy;
1124: ctx->transformctx = tctx;
1125: PetscFunctionReturn(PETSC_SUCCESS);
1126: }
1128: /*@C
1129: TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error
1130: in a time based line graph
1132: Collective
1134: Input Parameters:
1135: + ts - the `TS` context
1136: . step - current time-step
1137: . ptime - current time
1138: . u - current solution
1139: - dummy - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`
1141: Options Database Key:
1142: . -ts_monitor_lg_error - create a graphical monitor of error history
1144: Level: intermediate
1146: Notes:
1147: Each process in a parallel run displays its component errors in a separate window
1149: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1151: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1152: to be used during the TS integration.
1154: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1155: @*/
1156: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
1157: {
1158: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dummy;
1159: const PetscScalar *yy;
1160: Vec y;
1162: PetscFunctionBegin;
1163: if (!step) {
1164: PetscDrawAxis axis;
1165: PetscInt dim;
1166: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1167: PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
1168: PetscCall(VecGetLocalSize(u, &dim));
1169: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1170: PetscCall(PetscDrawLGReset(ctx->lg));
1171: }
1172: PetscCall(VecDuplicate(u, &y));
1173: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1174: PetscCall(VecAXPY(y, -1.0, u));
1175: PetscCall(VecGetArrayRead(y, &yy));
1176: #if defined(PETSC_USE_COMPLEX)
1177: {
1178: PetscReal *yreal;
1179: PetscInt i, n;
1180: PetscCall(VecGetLocalSize(y, &n));
1181: PetscCall(PetscMalloc1(n, &yreal));
1182: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1183: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1184: PetscCall(PetscFree(yreal));
1185: }
1186: #else
1187: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1188: #endif
1189: PetscCall(VecRestoreArrayRead(y, &yy));
1190: PetscCall(VecDestroy(&y));
1191: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1192: PetscCall(PetscDrawLGDraw(ctx->lg));
1193: PetscCall(PetscDrawLGSave(ctx->lg));
1194: }
1195: PetscFunctionReturn(PETSC_SUCCESS);
1196: }
1198: /*@C
1199: TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot
1201: Input Parameters:
1202: + ts - the `TS` context
1203: . step - current time-step
1204: . ptime - current time
1205: . u - current solution
1206: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()`
1208: Options Database Keys:
1209: + -ts_monitor_sp_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1210: . -ts_monitor_sp_swarm_retain <n> - Retain n old points so we can see the history, or -1 for all points
1211: . -ts_monitor_sp_swarm_multi_species <bool> - Color each species differently
1212: - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space
1214: Level: intermediate
1216: Notes:
1217: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1218: to be used during the `TS` integration.
1220: .seealso: [](ch_ts), `TS`, `TSMonitoSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
1221: @*/
1222: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1223: {
1224: TSMonitorSPCtx ctx = (TSMonitorSPCtx)dctx;
1225: PetscDraw draw;
1226: DM dm, cdm;
1227: const PetscScalar *yy;
1228: PetscInt Np, p, dim = 2, *species;
1229: PetscReal species_color;
1231: PetscFunctionBegin;
1232: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1233: PetscCall(TSGetDM(ts, &dm));
1234: if (!step) {
1235: PetscDrawAxis axis;
1236: PetscReal dmboxlower[2], dmboxupper[2];
1238: PetscCall(TSGetDM(ts, &dm));
1239: PetscCall(DMGetDimension(dm, &dim));
1240: PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
1241: PetscCall(DMSwarmGetCellDM(dm, &cdm));
1242: PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
1243: PetscCall(VecGetLocalSize(u, &Np));
1244: Np /= dim * 2;
1245: PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
1246: if (ctx->phase) {
1247: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
1248: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10));
1249: } else {
1250: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
1251: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
1252: }
1253: PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
1254: PetscCall(PetscDrawSPReset(ctx->sp));
1255: }
1256: if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species));
1257: PetscCall(VecGetLocalSize(u, &Np));
1258: Np /= dim * 2;
1259: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1260: PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
1261: if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
1262: PetscCall(PetscDrawFlush(draw));
1263: PetscCall(PetscDrawSPReset(ctx->sp));
1264: PetscCall(VecGetArrayRead(u, &yy));
1265: for (p = 0; p < Np; ++p) {
1266: PetscReal x, y;
1268: if (ctx->phase) {
1269: x = PetscRealPart(yy[p * dim * 2]);
1270: y = PetscRealPart(yy[p * dim * 2 + dim]);
1271: } else {
1272: x = PetscRealPart(yy[p * dim * 2]);
1273: y = PetscRealPart(yy[p * dim * 2 + 1]);
1274: }
1275: if (ctx->multispecies) {
1276: species_color = species[p] + 2;
1277: PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color));
1278: } else {
1279: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1280: }
1281: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1282: }
1283: PetscCall(VecRestoreArrayRead(u, &yy));
1284: PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
1285: PetscCall(PetscDrawSPSave(ctx->sp));
1286: if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species));
1287: }
1288: PetscFunctionReturn(PETSC_SUCCESS);
1289: }
1291: /*@C
1292: TSMonitorHGSwarmSolution - Graphically displays histograms of `DMSWARM` particles
1294: Input Parameters:
1295: + ts - the `TS` context
1296: . step - current time-step
1297: . ptime - current time
1298: . u - current solution
1299: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorHGCtxCreate()`
1301: Options Database Keys:
1302: + -ts_monitor_hg_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1303: . -ts_monitor_hg_swarm_species <num> - Number of species to histogram
1304: . -ts_monitor_hg_swarm_bins <num> - Number of histogram bins
1305: - -ts_monitor_hg_swarm_velocity <bool> - Plot in velocity space, as opposed to coordinate space
1307: Level: intermediate
1309: Note:
1310: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1311: to be used during the `TS` integration.
1313: .seealso: `TSMonitoSet()`
1314: @*/
1315: PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1316: {
1317: TSMonitorHGCtx ctx = (TSMonitorHGCtx)dctx;
1318: PetscDraw draw;
1319: DM sw;
1320: const PetscScalar *yy;
1321: PetscInt *species;
1322: PetscInt dim, d = 0, Np, p, Ns, s;
1324: PetscFunctionBegin;
1325: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1326: PetscCall(TSGetDM(ts, &sw));
1327: PetscCall(DMGetDimension(sw, &dim));
1328: PetscCall(DMSwarmGetNumSpecies(sw, &Ns));
1329: Ns = PetscMin(Ns, ctx->Ns);
1330: PetscCall(VecGetLocalSize(u, &Np));
1331: Np /= dim * 2;
1332: if (!step) {
1333: PetscDrawAxis axis;
1334: char title[PETSC_MAX_PATH_LEN];
1336: for (s = 0; s < Ns; ++s) {
1337: PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis));
1338: PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s));
1339: if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N"));
1340: else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N"));
1341: }
1342: }
1343: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1344: PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species));
1345: for (s = 0; s < Ns; ++s) {
1346: PetscCall(PetscDrawHGReset(ctx->hg[s]));
1347: PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw));
1348: PetscCall(PetscDrawClear(draw));
1349: PetscCall(PetscDrawFlush(draw));
1350: }
1351: PetscCall(VecGetArrayRead(u, &yy));
1352: for (p = 0; p < Np; ++p) {
1353: const PetscInt s = species[p] < Ns ? species[p] : 0;
1354: PetscReal v;
1356: if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]);
1357: else v = PetscRealPart(yy[p * dim * 2 + d]);
1358: PetscCall(PetscDrawHGAddValue(ctx->hg[s], v));
1359: }
1360: PetscCall(VecRestoreArrayRead(u, &yy));
1361: for (s = 0; s < Ns; ++s) {
1362: PetscCall(PetscDrawHGDraw(ctx->hg[s]));
1363: PetscCall(PetscDrawHGSave(ctx->hg[s]));
1364: }
1365: PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species));
1366: }
1367: PetscFunctionReturn(PETSC_SUCCESS);
1368: }
1370: /*@C
1371: TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep
1373: Collective
1375: Input Parameters:
1376: + ts - the `TS` context
1377: . step - current time-step
1378: . ptime - current time
1379: . u - current solution
1380: - vf - unused context
1382: Options Database Key:
1383: . -ts_monitor_error - create a graphical monitor of error history
1385: Level: intermediate
1387: Notes:
1388: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1389: to be used during the `TS` integration.
1391: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1393: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1394: @*/
1395: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1396: {
1397: DM dm;
1398: PetscDS ds = NULL;
1399: PetscInt Nf = -1, f;
1400: PetscBool flg;
1402: PetscFunctionBegin;
1403: PetscCall(TSGetDM(ts, &dm));
1404: if (dm) PetscCall(DMGetDS(dm, &ds));
1405: if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
1406: if (Nf <= 0) {
1407: Vec y;
1408: PetscReal nrm;
1410: PetscCall(VecDuplicate(u, &y));
1411: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1412: PetscCall(VecAXPY(y, -1.0, u));
1413: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
1414: if (flg) {
1415: PetscCall(VecNorm(y, NORM_2, &nrm));
1416: PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
1417: }
1418: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
1419: if (flg) PetscCall(VecView(y, vf->viewer));
1420: PetscCall(VecDestroy(&y));
1421: } else {
1422: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
1423: void **ctxs;
1424: Vec v;
1425: PetscReal ferrors[1];
1427: PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
1428: for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
1429: PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
1430: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04d time = %-8.4g \t L_2 Error: [", (int)step, (double)ptime));
1431: for (f = 0; f < Nf; ++f) {
1432: if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
1433: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
1434: }
1435: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));
1437: PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));
1439: PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
1440: if (flg) {
1441: PetscCall(DMGetGlobalVector(dm, &v));
1442: PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
1443: PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
1444: PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
1445: PetscCall(DMRestoreGlobalVector(dm, &v));
1446: }
1447: PetscCall(PetscFree2(exactFuncs, ctxs));
1448: }
1449: PetscFunctionReturn(PETSC_SUCCESS);
1450: }
1452: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1453: {
1454: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1455: PetscReal x = ptime, y;
1456: PetscInt its;
1458: PetscFunctionBegin;
1459: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1460: if (!n) {
1461: PetscDrawAxis axis;
1462: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1463: PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
1464: PetscCall(PetscDrawLGReset(ctx->lg));
1465: ctx->snes_its = 0;
1466: }
1467: PetscCall(TSGetSNESIterations(ts, &its));
1468: y = its - ctx->snes_its;
1469: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1470: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1471: PetscCall(PetscDrawLGDraw(ctx->lg));
1472: PetscCall(PetscDrawLGSave(ctx->lg));
1473: }
1474: ctx->snes_its = its;
1475: PetscFunctionReturn(PETSC_SUCCESS);
1476: }
1478: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1479: {
1480: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1481: PetscReal x = ptime, y;
1482: PetscInt its;
1484: PetscFunctionBegin;
1485: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1486: if (!n) {
1487: PetscDrawAxis axis;
1488: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1489: PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
1490: PetscCall(PetscDrawLGReset(ctx->lg));
1491: ctx->ksp_its = 0;
1492: }
1493: PetscCall(TSGetKSPIterations(ts, &its));
1494: y = its - ctx->ksp_its;
1495: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1496: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1497: PetscCall(PetscDrawLGDraw(ctx->lg));
1498: PetscCall(PetscDrawLGSave(ctx->lg));
1499: }
1500: ctx->ksp_its = its;
1501: PetscFunctionReturn(PETSC_SUCCESS);
1502: }
1504: /*@C
1505: TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()`
1507: Collective
1509: Input Parameter:
1510: . ts - the `TS` solver object
1512: Output Parameter:
1513: . ctx - the context
1515: Level: intermediate
1517: .seealso: [](ch_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1518: @*/
1519: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1520: {
1521: PetscFunctionBegin;
1522: PetscCall(PetscNew(ctx));
1523: PetscFunctionReturn(PETSC_SUCCESS);
1524: }
1526: /*@C
1527: TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution
1529: Collective
1531: Input Parameters:
1532: + ts - the `TS` context
1533: . step - current time-step
1534: . ptime - current time
1535: . u - current solution
1536: - dctx - the envelope context
1538: Options Database Key:
1539: . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time
1541: Level: intermediate
1543: Notes:
1544: After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope
1546: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1547: to be used during the `TS` integration.
1549: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1550: @*/
1551: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1552: {
1553: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;
1555: PetscFunctionBegin;
1556: if (!ctx->max) {
1557: PetscCall(VecDuplicate(u, &ctx->max));
1558: PetscCall(VecDuplicate(u, &ctx->min));
1559: PetscCall(VecCopy(u, ctx->max));
1560: PetscCall(VecCopy(u, ctx->min));
1561: } else {
1562: PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
1563: PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
1564: }
1565: PetscFunctionReturn(PETSC_SUCCESS);
1566: }
1568: /*@C
1569: TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution
1571: Collective
1573: Input Parameter:
1574: . ts - the `TS` context
1576: Output Parameters:
1577: + max - the maximum values
1578: - min - the minimum values
1580: Level: intermediate
1582: Notes:
1583: If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored
1585: .seealso: [](ch_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1586: @*/
1587: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1588: {
1589: PetscInt i;
1591: PetscFunctionBegin;
1592: if (max) *max = NULL;
1593: if (min) *min = NULL;
1594: for (i = 0; i < ts->numbermonitors; i++) {
1595: if (ts->monitor[i] == TSMonitorEnvelope) {
1596: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1597: if (max) *max = ctx->max;
1598: if (min) *min = ctx->min;
1599: break;
1600: }
1601: }
1602: PetscFunctionReturn(PETSC_SUCCESS);
1603: }
1605: /*@C
1606: TSMonitorEnvelopeCtxDestroy - Destroys a context that was created with `TSMonitorEnvelopeCtxCreate()`.
1608: Collective
1610: Input Parameter:
1611: . ctx - the monitor context
1613: Level: intermediate
1615: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1616: @*/
1617: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1618: {
1619: PetscFunctionBegin;
1620: PetscCall(VecDestroy(&(*ctx)->min));
1621: PetscCall(VecDestroy(&(*ctx)->max));
1622: PetscCall(PetscFree(*ctx));
1623: PetscFunctionReturn(PETSC_SUCCESS);
1624: }
1626: /*@C
1627: TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS`
1629: Not Collective
1631: Input Parameters:
1632: + ts - the `TS` context
1633: . step - current timestep
1634: . t - current time
1635: . U - current solution
1636: - vf - not used
1638: Options Database Key:
1639: . -ts_dmswarm_monitor_moments - Monitor moments of particle distribution
1641: Level: intermediate
1643: Notes:
1644: This requires a `DMSWARM` be attached to the `TS`.
1646: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1647: to be used during the TS integration.
1649: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1650: @*/
1651: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1652: {
1653: DM sw;
1654: const PetscScalar *u;
1655: PetscReal m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1656: PetscInt dim, d, Np, p;
1657: MPI_Comm comm;
1659: PetscFunctionBeginUser;
1660: (void)t;
1661: (void)vf;
1662: PetscCall(TSGetDM(ts, &sw));
1663: if (!sw || step % ts->monitorFrequency != 0) PetscFunctionReturn(PETSC_SUCCESS);
1664: PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
1665: PetscCall(DMGetDimension(sw, &dim));
1666: PetscCall(VecGetLocalSize(U, &Np));
1667: Np /= dim;
1668: PetscCall(VecGetArrayRead(U, &u));
1669: for (p = 0; p < Np; ++p) {
1670: for (d = 0; d < dim; ++d) {
1671: totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1672: totMom[d] += PetscRealPart(u[p * dim + d]);
1673: }
1674: }
1675: PetscCall(VecRestoreArrayRead(U, &u));
1676: for (d = 0; d < dim; ++d) totMom[d] *= m;
1677: totE *= 0.5 * m;
1678: PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
1679: for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, " Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
1680: PetscCall(PetscPrintf(comm, "\n"));
1681: PetscFunctionReturn(PETSC_SUCCESS);
1682: }