Actual source code: tsmon.c
1: #include <petsc/private/tsimpl.h>
2: #include <petscdm.h>
3: #include <petscds.h>
4: #include <petscdmswarm.h>
5: #include <petscdraw.h>
7: /*@C
8: TSMonitor - Runs all user-provided monitor routines set using `TSMonitorSet()`
10: Collective
12: Input Parameters:
13: + ts - time stepping context obtained from `TSCreate()`
14: . step - step number that has just completed
15: . ptime - model time of the state
16: - u - state at the current model time
18: Level: developer
20: Notes:
21: `TSMonitor()` is typically used automatically within the time stepping implementations.
22: Users would almost never call this routine directly.
24: A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions
26: .seealso: `TS`, `TSMonitorSet()`, `TSMonitorSetFromOptions()`
27: @*/
28: PetscErrorCode TSMonitor(TS ts, PetscInt step, PetscReal ptime, Vec u)
29: {
30: DM dm;
31: PetscInt i, n = ts->numbermonitors;
33: PetscFunctionBegin;
37: PetscCall(TSGetDM(ts, &dm));
38: PetscCall(DMSetOutputSequenceNumber(dm, step, ptime));
40: PetscCall(VecLockReadPush(u));
41: for (i = 0; i < n; i++) PetscCall((*ts->monitor[i])(ts, step, ptime, u, ts->monitorcontext[i]));
42: PetscCall(VecLockReadPop(u));
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: /*@C
47: TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
49: Collective
51: Input Parameters:
52: + ts - `TS` object you wish to monitor
53: . name - the monitor type one is seeking
54: . help - message indicating what monitoring is done
55: . manual - manual page for the monitor
56: . monitor - the monitor function, this must use a `PetscViewerFormat` as its context
57: - monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the `TS` or `PetscViewer` objects
59: Level: developer
61: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `PetscOptionsCreateViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
62: `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
63: `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`,
64: `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
65: `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
66: `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
67: `PetscOptionsFList()`, `PetscOptionsEList()`
68: @*/
69: PetscErrorCode TSMonitorSetFromOptions(TS ts, const char name[], const char help[], const char manual[], PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *), PetscErrorCode (*monitorsetup)(TS, PetscViewerAndFormat *))
70: {
71: PetscViewer viewer;
72: PetscViewerFormat format;
73: PetscBool flg;
75: PetscFunctionBegin;
76: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)ts), ((PetscObject)ts)->options, ((PetscObject)ts)->prefix, name, &viewer, &format, &flg));
77: if (flg) {
78: PetscViewerAndFormat *vf;
79: char interval_key[1024];
81: PetscCall(PetscSNPrintf(interval_key, sizeof interval_key, "%s_interval", name));
82: PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
83: vf->view_interval = 1;
84: PetscCall(PetscOptionsGetInt(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, interval_key, &vf->view_interval, NULL));
85: PetscCall(PetscViewerDestroy(&viewer));
86: if (monitorsetup) PetscCall((*monitorsetup)(ts, vf));
87: PetscCall(TSMonitorSet(ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, void *))monitor, vf, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy));
88: }
89: PetscFunctionReturn(PETSC_SUCCESS);
90: }
92: /*@C
93: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
94: timestep to display the iteration's progress.
96: Logically Collective
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`), see `PetscCtxDestroyFn` for the calling sequence
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
111: Level: intermediate
113: Note:
114: This routine adds an additional monitor to the list of monitors that already has been loaded.
116: Fortran Notes:
117: Only a single monitor function can be set for each `TS` object
119: .seealso: [](ch_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
120: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
121: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `PetscCtxDestroyFn`
122: @*/
123: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS ts, PetscInt steps, PetscReal time, Vec u, void *ctx), void *mctx, PetscCtxDestroyFn *mdestroy)
124: {
125: PetscInt i;
126: PetscBool identical;
128: PetscFunctionBegin;
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++] = mctx;
138: PetscFunctionReturn(PETSC_SUCCESS);
139: }
141: /*@C
142: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
144: Logically Collective
146: Input Parameter:
147: . ts - the `TS` context obtained from `TSCreate()`
149: Level: intermediate
151: Note:
152: There is no way to remove a single, specific monitor.
154: .seealso: [](ch_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
155: @*/
156: PetscErrorCode TSMonitorCancel(TS ts)
157: {
158: PetscInt i;
160: PetscFunctionBegin;
162: for (i = 0; i < ts->numbermonitors; i++) {
163: if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i]));
164: }
165: ts->numbermonitors = 0;
166: PetscFunctionReturn(PETSC_SUCCESS);
167: }
169: /*@C
170: TSMonitorDefault - The default monitor, prints the timestep and time for each step
172: Input Parameters:
173: + ts - the `TS` context
174: . step - iteration number (after the final time step the monitor routine may be called with a step of -1, this indicates the solution has been interpolated to this time)
175: . ptime - current time
176: . v - current iterate
177: - vf - the viewer and format
179: Options Database Key:
180: . -ts_monitor - monitors the time integration
182: Level: intermediate
184: Notes:
185: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
186: to be used during the `TS` integration.
188: .seealso: [](ch_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `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;
197: PetscFunctionBegin;
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;
217: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
218: if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
219: PetscCall(PetscRealView(1, &ptime, viewer));
220: } else {
221: PetscCall(PetscRealView(0, &ptime, viewer));
222: }
223: }
224: PetscCall(PetscViewerPopFormat(viewer));
225: PetscFunctionReturn(PETSC_SUCCESS);
226: }
228: /*@C
229: TSMonitorExtreme - Prints the extreme values of the solution at each timestep
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
238: Level: intermediate
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.
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;
252: PetscFunctionBegin;
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: }
267: /*@C
268: TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with
269: `TS` to monitor the solution process graphically in various ways
271: Collective
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
283: Output Parameter:
284: . ctx - the context
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
295: Level: intermediate
297: Notes:
298: Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed.
300: One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()`
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.
306: One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()`
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;
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: }
329: /*@C
330: TSMonitorLGTimeStep - Monitors a `TS` by printing the time-steps
332: Collective
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()`
341: Level: advanced
343: Note:
344: This is not called directly by users, rather one calls `TSMonitorSet()` along the `ctx` created by `TSMonitorLGCtxCreate()`
345: and `TSMonitorLGCtxDestroy()`
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;
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: }
373: /*@C
374: TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`.
376: Collective
378: Input Parameter:
379: . ctx - the monitor context
381: Level: intermediate
383: Note:
384: Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`
386: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
387: @*/
388: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
389: {
390: PetscFunctionBegin;
391: if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)((void **)&(*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: }
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;
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: }
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: }
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: int Nsi, Nbi;
434: PetscFunctionBegin;
435: PetscCall(PetscMPIIntCast(Ns, &Nsi));
436: PetscCall(PetscMPIIntCast(Nb, &Nbi));
437: PetscCall(PetscNew(ctx));
438: PetscCall(PetscMalloc1(Ns, &(*ctx)->hg));
439: for (int s = 0; s < Nsi; ++s) {
440: PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw));
441: PetscCall(PetscDrawSetFromOptions(draw));
442: PetscCall(PetscDrawHGCreate(draw, Nbi, &(*ctx)->hg[s]));
443: PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE));
444: PetscCall(PetscDrawDestroy(&draw));
445: }
446: (*ctx)->howoften = howoften;
447: (*ctx)->Ns = Ns;
448: (*ctx)->velocity = velocity;
449: PetscFunctionReturn(PETSC_SUCCESS);
450: }
452: /* Destroys a TSMonitorHGCtx that was created with TSMonitorHGCtxCreate */
453: PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx)
454: {
455: PetscInt s;
457: PetscFunctionBegin;
458: for (s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s]));
459: PetscCall(PetscFree((*ctx)->hg));
460: PetscCall(PetscFree(*ctx));
461: PetscFunctionReturn(PETSC_SUCCESS);
462: }
464: /*@C
465: TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling
466: `VecView()` for the solution at each timestep
468: Collective
470: Input Parameters:
471: + ts - the `TS` context
472: . step - current time-step
473: . ptime - current time
474: . u - the solution at the current time
475: - dummy - either a viewer or `NULL`
477: Options Database Keys:
478: + -ts_monitor_draw_solution - draw the solution at each time-step
479: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
481: Level: intermediate
483: Notes:
484: The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial`
485: will look bad
487: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with
488: `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration.
490: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
491: @*/
492: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
493: {
494: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
495: PetscDraw draw;
497: PetscFunctionBegin;
498: if (!step && ictx->showinitial) {
499: if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
500: PetscCall(VecCopy(u, ictx->initialsolution));
501: }
502: if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
504: if (ictx->showinitial) {
505: PetscReal pause;
506: PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
507: PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
508: PetscCall(VecView(ictx->initialsolution, ictx->viewer));
509: PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
510: PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
511: }
512: PetscCall(VecView(u, ictx->viewer));
513: if (ictx->showtimestepandtime) {
514: PetscReal xl, yl, xr, yr, h;
515: char time[32];
517: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
518: PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
519: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
520: h = yl + .95 * (yr - yl);
521: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
522: PetscCall(PetscDrawFlush(draw));
523: }
525: if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
526: PetscFunctionReturn(PETSC_SUCCESS);
527: }
529: /*@C
530: TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram
532: Collective
534: Input Parameters:
535: + ts - the `TS` context
536: . step - current time-step
537: . ptime - current time
538: . u - the solution at the current time
539: - dummy - either a viewer or `NULL`
541: Level: intermediate
543: Notes:
544: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
545: to be used during the `TS` integration.
547: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
548: @*/
549: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
550: {
551: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
552: PetscDraw draw;
553: PetscDrawAxis axis;
554: PetscInt n;
555: PetscMPIInt size;
556: PetscReal U0, U1, xl, yl, xr, yr, h;
557: char time[32];
558: const PetscScalar *U;
560: PetscFunctionBegin;
561: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
562: PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
563: PetscCall(VecGetSize(u, &n));
564: PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");
566: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
567: PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
568: PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
569: if (!step) {
570: PetscCall(PetscDrawClear(draw));
571: PetscCall(PetscDrawAxisDraw(axis));
572: }
574: PetscCall(VecGetArrayRead(u, &U));
575: U0 = PetscRealPart(U[0]);
576: U1 = PetscRealPart(U[1]);
577: PetscCall(VecRestoreArrayRead(u, &U));
578: if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);
580: PetscDrawCollectiveBegin(draw);
581: PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
582: if (ictx->showtimestepandtime) {
583: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
584: PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
585: h = yl + .95 * (yr - yl);
586: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
587: }
588: PetscDrawCollectiveEnd(draw);
589: PetscCall(PetscDrawFlush(draw));
590: PetscCall(PetscDrawPause(draw));
591: PetscCall(PetscDrawSave(draw));
592: PetscFunctionReturn(PETSC_SUCCESS);
593: }
595: /*@C
596: TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()`
598: Collective
600: Input Parameter:
601: . ictx - the monitor context
603: Level: intermediate
605: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
606: @*/
607: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
608: {
609: PetscFunctionBegin;
610: PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
611: PetscCall(VecDestroy(&(*ictx)->initialsolution));
612: PetscCall(PetscFree(*ictx));
613: PetscFunctionReturn(PETSC_SUCCESS);
614: }
616: /*@C
617: TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`
619: Collective
621: Input Parameters:
622: + comm - the MPI communicator to use
623: . host - the X display to open, or `NULL` for the local machine
624: . label - the title to put in the title bar
625: . x - the x screen coordinates of the upper left coordinate of the window
626: . y - the y screen coordinates of the upper left coordinate of the window
627: . m - the screen width in pixels
628: . n - the screen height in pixels
629: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
631: Output Parameter:
632: . ctx - the monitor context
634: Options Database Keys:
635: + -ts_monitor_draw_solution - draw the solution at each time-step
636: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
638: Level: intermediate
640: Note:
641: The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`.
643: .seealso: [](ch_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
644: @*/
645: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
646: {
647: PetscFunctionBegin;
648: PetscCall(PetscNew(ctx));
649: PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
650: PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));
652: (*ctx)->howoften = howoften;
653: (*ctx)->showinitial = PETSC_FALSE;
654: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));
656: (*ctx)->showtimestepandtime = PETSC_FALSE;
657: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
658: PetscFunctionReturn(PETSC_SUCCESS);
659: }
661: /*@C
662: TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling
663: `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep
665: Collective
667: Input Parameters:
668: + ts - the `TS` context
669: . step - current time-step
670: . ptime - current time
671: . u - solution at current time
672: - dummy - either a viewer or `NULL`
674: Options Database Key:
675: . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
677: Level: intermediate
679: Note:
680: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
681: to be used during the `TS` integration.
683: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
684: @*/
685: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
686: {
687: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
688: PetscViewer viewer = ctx->viewer;
689: Vec work;
691: PetscFunctionBegin;
692: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
693: PetscCall(VecDuplicate(u, &work));
694: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
695: PetscCall(VecView(work, viewer));
696: PetscCall(VecDestroy(&work));
697: PetscFunctionReturn(PETSC_SUCCESS);
698: }
700: /*@C
701: TSMonitorDrawError - Monitors progress of the `TS` solvers by calling
702: `VecView()` for the error at each timestep
704: Collective
706: Input Parameters:
707: + ts - the `TS` context
708: . step - current time-step
709: . ptime - current time
710: . u - solution at current time
711: - dummy - either a viewer or `NULL`
713: Options Database Key:
714: . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
716: Level: intermediate
718: Notes:
719: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
720: to be used during the `TS` integration.
722: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
723: @*/
724: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
725: {
726: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
727: PetscViewer viewer = ctx->viewer;
728: Vec work;
730: PetscFunctionBegin;
731: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
732: PetscCall(VecDuplicate(u, &work));
733: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
734: PetscCall(VecAXPY(work, -1.0, u));
735: PetscCall(VecView(work, viewer));
736: PetscCall(VecDestroy(&work));
737: PetscFunctionReturn(PETSC_SUCCESS);
738: }
740: /*@C
741: 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
743: Collective
745: Input Parameters:
746: + ts - the `TS` context
747: . step - current time-step
748: . ptime - current time
749: . u - current state
750: - vf - viewer and its format
752: Level: intermediate
754: Notes:
755: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
756: to be used during the `TS` integration.
758: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
759: @*/
760: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
761: {
762: PetscFunctionBegin;
763: if ((vf->view_interval > 0 && !(step % vf->view_interval)) || (vf->view_interval && ts->reason)) {
764: PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
765: PetscCall(VecView(u, vf->viewer));
766: PetscCall(PetscViewerPopFormat(vf->viewer));
767: }
768: PetscFunctionReturn(PETSC_SUCCESS);
769: }
771: /*@C
772: TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at selected timesteps.
774: Collective
776: Input Parameters:
777: + ts - the `TS` context
778: . step - current time-step
779: . ptime - current time
780: . u - current state
781: - ctx - monitor context obtained with `TSMonitorSolutionVTKCtxCreate()`
783: Level: developer
785: Notes:
786: 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.
787: These are named according to the file name template.
789: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
790: to be used during the `TS` integration.
792: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
793: @*/
794: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, TSMonitorVTKCtx ctx)
795: {
796: char filename[PETSC_MAX_PATH_LEN];
797: PetscViewer viewer;
799: PetscFunctionBegin;
800: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
801: if (((ctx->interval > 0) && (!(step % ctx->interval))) || (ctx->interval && ts->reason)) {
802: PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)ctx->filenametemplate, step));
803: PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
804: PetscCall(VecView(u, viewer));
805: PetscCall(PetscViewerDestroy(&viewer));
806: }
807: PetscFunctionReturn(PETSC_SUCCESS);
808: }
810: /*@C
811: TSMonitorSolutionVTKDestroy - Destroy the monitor context created with `TSMonitorSolutionVTKCtxCreate()`
813: Not Collective
815: Input Parameter:
816: . ctx - the monitor context
818: Level: developer
820: Note:
821: This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
823: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
824: @*/
825: PetscErrorCode TSMonitorSolutionVTKDestroy(TSMonitorVTKCtx *ctx)
826: {
827: PetscFunctionBegin;
828: PetscCall(PetscFree((*ctx)->filenametemplate));
829: PetscCall(PetscFree(*ctx));
830: PetscFunctionReturn(PETSC_SUCCESS);
831: }
833: /*@C
834: TSMonitorSolutionVTKCtxCreate - Create the monitor context to be used in `TSMonitorSolutionVTK()`
836: Not collective
838: Input Parameter:
839: . filenametemplate - the template file name, e.g. foo-%03d.vts
841: Output Parameter:
842: . ctx - the monitor context
844: Level: developer
846: Note:
847: This function is normally used inside `TSSetFromOptions()` to pass the context created to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
849: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`, `TSMonitorSolutionVTKDestroy()`
850: @*/
851: PetscErrorCode TSMonitorSolutionVTKCtxCreate(const char *filenametemplate, TSMonitorVTKCtx *ctx)
852: {
853: const char *ptr = NULL, *ptr2 = NULL;
854: TSMonitorVTKCtx ictx;
856: PetscFunctionBegin;
857: PetscAssertPointer(filenametemplate, 1);
858: PetscAssertPointer(ctx, 2);
859: /* Do some cursory validation of the input. */
860: PetscCall(PetscStrstr(filenametemplate, "%", (char **)&ptr));
861: PetscCheck(ptr, PETSC_COMM_SELF, PETSC_ERR_USER, "-ts_monitor_solution_vtk requires a file template, e.g. filename-%%03" PetscInt_FMT ".vts");
862: for (ptr++; ptr && *ptr; ptr++) {
863: PetscCall(PetscStrchr("DdiouxX", *ptr, (char **)&ptr2));
864: 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");
865: if (ptr2) break;
866: }
867: PetscCall(PetscNew(&ictx));
868: PetscCall(PetscStrallocpy(filenametemplate, &ictx->filenametemplate));
869: ictx->interval = 1;
871: *ctx = ictx;
872: PetscFunctionReturn(PETSC_SUCCESS);
873: }
875: /*@C
876: TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector
877: in a time based line graph
879: Collective
881: Input Parameters:
882: + ts - the `TS` context
883: . step - current time-step
884: . ptime - current time
885: . u - current solution
886: - dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`
888: Options Database Key:
889: . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables
891: Level: intermediate
893: Notes:
894: Each process in a parallel run displays its component solutions in a separate window
896: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
897: to be used during the `TS` integration.
899: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
900: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
901: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
902: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
903: @*/
904: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
905: {
906: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dctx;
907: const PetscScalar *yy;
908: Vec v;
910: PetscFunctionBegin;
911: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
912: if (!step) {
913: PetscDrawAxis axis;
914: PetscInt dim;
915: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
916: PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
917: if (!ctx->names) {
918: PetscBool flg;
919: /* user provides names of variables to plot but no names has been set so assume names are integer values */
920: PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
921: if (flg) {
922: PetscInt i, n;
923: char **names;
924: PetscCall(VecGetSize(u, &n));
925: PetscCall(PetscMalloc1(n + 1, &names));
926: for (i = 0; i < n; i++) {
927: PetscCall(PetscMalloc1(5, &names[i]));
928: PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
929: }
930: names[n] = NULL;
931: ctx->names = names;
932: }
933: }
934: if (ctx->names && !ctx->displaynames) {
935: char **displaynames;
936: PetscBool flg;
937: PetscCall(VecGetLocalSize(u, &dim));
938: PetscCall(PetscCalloc1(dim + 1, &displaynames));
939: PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
940: if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
941: PetscCall(PetscStrArrayDestroy(&displaynames));
942: }
943: if (ctx->displaynames) {
944: PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
945: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
946: } else if (ctx->names) {
947: PetscCall(VecGetLocalSize(u, &dim));
948: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
949: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
950: } else {
951: PetscCall(VecGetLocalSize(u, &dim));
952: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
953: }
954: PetscCall(PetscDrawLGReset(ctx->lg));
955: }
957: if (!ctx->transform) v = u;
958: else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
959: PetscCall(VecGetArrayRead(v, &yy));
960: if (ctx->displaynames) {
961: PetscInt i;
962: for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
963: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
964: } else {
965: #if defined(PETSC_USE_COMPLEX)
966: PetscInt i, n;
967: PetscReal *yreal;
968: PetscCall(VecGetLocalSize(v, &n));
969: PetscCall(PetscMalloc1(n, &yreal));
970: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
971: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
972: PetscCall(PetscFree(yreal));
973: #else
974: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
975: #endif
976: }
977: PetscCall(VecRestoreArrayRead(v, &yy));
978: if (ctx->transform) PetscCall(VecDestroy(&v));
980: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
981: PetscCall(PetscDrawLGDraw(ctx->lg));
982: PetscCall(PetscDrawLGSave(ctx->lg));
983: }
984: PetscFunctionReturn(PETSC_SUCCESS);
985: }
987: /*@C
988: TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
990: Collective
992: Input Parameters:
993: + ts - the `TS` context
994: - names - the names of the components, final string must be `NULL`
996: Level: intermediate
998: Notes:
999: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1001: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
1002: @*/
1003: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
1004: {
1005: PetscInt i;
1007: PetscFunctionBegin;
1008: for (i = 0; i < ts->numbermonitors; i++) {
1009: if (ts->monitor[i] == TSMonitorLGSolution) {
1010: PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
1011: break;
1012: }
1013: }
1014: PetscFunctionReturn(PETSC_SUCCESS);
1015: }
1017: /*@C
1018: TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
1020: Collective
1022: Input Parameters:
1023: + ctx - the `TS` context
1024: - names - the names of the components, final string must be `NULL`
1026: Level: intermediate
1028: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
1029: @*/
1030: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
1031: {
1032: PetscFunctionBegin;
1033: PetscCall(PetscStrArrayDestroy(&ctx->names));
1034: PetscCall(PetscStrArrayallocpy(names, &ctx->names));
1035: PetscFunctionReturn(PETSC_SUCCESS);
1036: }
1038: /*@C
1039: TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot
1041: Collective
1043: Input Parameter:
1044: . ts - the `TS` context
1046: Output Parameter:
1047: . names - the names of the components, final string must be `NULL`
1049: Level: intermediate
1051: Note:
1052: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1054: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1055: @*/
1056: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
1057: {
1058: PetscInt i;
1060: PetscFunctionBegin;
1061: *names = NULL;
1062: for (i = 0; i < ts->numbermonitors; i++) {
1063: if (ts->monitor[i] == TSMonitorLGSolution) {
1064: TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
1065: *names = (const char *const *)ctx->names;
1066: break;
1067: }
1068: }
1069: PetscFunctionReturn(PETSC_SUCCESS);
1070: }
1072: /*@C
1073: TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor
1075: Collective
1077: Input Parameters:
1078: + ctx - the `TSMonitorLG` context
1079: - displaynames - the names of the components, final string must be `NULL`
1081: Level: intermediate
1083: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1084: @*/
1085: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
1086: {
1087: PetscInt j = 0, k;
1089: PetscFunctionBegin;
1090: if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
1091: PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
1092: PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
1093: while (displaynames[j]) j++;
1094: ctx->ndisplayvariables = j;
1095: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
1096: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
1097: j = 0;
1098: while (displaynames[j]) {
1099: k = 0;
1100: while (ctx->names[k]) {
1101: PetscBool flg;
1102: PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
1103: if (flg) {
1104: ctx->displayvariables[j] = k;
1105: break;
1106: }
1107: k++;
1108: }
1109: j++;
1110: }
1111: PetscFunctionReturn(PETSC_SUCCESS);
1112: }
1114: /*@C
1115: TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor
1117: Collective
1119: Input Parameters:
1120: + ts - the `TS` context
1121: - displaynames - the names of the components, final string must be `NULL`
1123: Level: intermediate
1125: Note:
1126: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1128: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1129: @*/
1130: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
1131: {
1132: PetscInt i;
1134: PetscFunctionBegin;
1135: for (i = 0; i < ts->numbermonitors; i++) {
1136: if (ts->monitor[i] == TSMonitorLGSolution) {
1137: PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
1138: break;
1139: }
1140: }
1141: PetscFunctionReturn(PETSC_SUCCESS);
1142: }
1144: /*@C
1145: TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed
1147: Collective
1149: Input Parameters:
1150: + ts - the `TS` context
1151: . transform - the transform function
1152: . destroy - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1153: - tctx - optional context used by transform function
1155: Level: intermediate
1157: Note:
1158: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1160: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()`, `PetscCtxDestroyFn`
1161: @*/
1162: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscCtxDestroyFn *destroy, void *tctx)
1163: {
1164: PetscInt i;
1166: PetscFunctionBegin;
1167: for (i = 0; i < ts->numbermonitors; i++) {
1168: if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
1169: }
1170: PetscFunctionReturn(PETSC_SUCCESS);
1171: }
1173: /*@C
1174: TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed
1176: Collective
1178: Input Parameters:
1179: + tctx - the `TS` context
1180: . transform - the transform function
1181: . destroy - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1182: - ctx - optional context used by transform function
1184: Level: intermediate
1186: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`, `PetscCtxDestroyFn`
1187: @*/
1188: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscCtxDestroyFn *destroy, void *tctx)
1189: {
1190: PetscFunctionBegin;
1191: ctx->transform = transform;
1192: ctx->transformdestroy = destroy;
1193: ctx->transformctx = tctx;
1194: PetscFunctionReturn(PETSC_SUCCESS);
1195: }
1197: /*@C
1198: TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error
1199: in a time based line graph
1201: Collective
1203: Input Parameters:
1204: + ts - the `TS` context
1205: . step - current time-step
1206: . ptime - current time
1207: . u - current solution
1208: - dummy - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`
1210: Options Database Key:
1211: . -ts_monitor_lg_error - create a graphical monitor of error history
1213: Level: intermediate
1215: Notes:
1216: Each process in a parallel run displays its component errors in a separate window
1218: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1220: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1221: to be used during the TS integration.
1223: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1224: @*/
1225: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
1226: {
1227: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dummy;
1228: const PetscScalar *yy;
1229: Vec y;
1231: PetscFunctionBegin;
1232: if (!step) {
1233: PetscDrawAxis axis;
1234: PetscInt dim;
1235: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1236: PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
1237: PetscCall(VecGetLocalSize(u, &dim));
1238: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1239: PetscCall(PetscDrawLGReset(ctx->lg));
1240: }
1241: PetscCall(VecDuplicate(u, &y));
1242: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1243: PetscCall(VecAXPY(y, -1.0, u));
1244: PetscCall(VecGetArrayRead(y, &yy));
1245: #if defined(PETSC_USE_COMPLEX)
1246: {
1247: PetscReal *yreal;
1248: PetscInt i, n;
1249: PetscCall(VecGetLocalSize(y, &n));
1250: PetscCall(PetscMalloc1(n, &yreal));
1251: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1252: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1253: PetscCall(PetscFree(yreal));
1254: }
1255: #else
1256: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1257: #endif
1258: PetscCall(VecRestoreArrayRead(y, &yy));
1259: PetscCall(VecDestroy(&y));
1260: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1261: PetscCall(PetscDrawLGDraw(ctx->lg));
1262: PetscCall(PetscDrawLGSave(ctx->lg));
1263: }
1264: PetscFunctionReturn(PETSC_SUCCESS);
1265: }
1267: /*@C
1268: TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot
1270: Input Parameters:
1271: + ts - the `TS` context
1272: . step - current time-step
1273: . ptime - current time
1274: . u - current solution
1275: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()`
1277: Options Database Keys:
1278: + -ts_monitor_sp_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1279: . -ts_monitor_sp_swarm_retain <n> - Retain n old points so we can see the history, or -1 for all points
1280: . -ts_monitor_sp_swarm_multi_species <bool> - Color each species differently
1281: - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space
1283: Level: intermediate
1285: Notes:
1286: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1287: to be used during the `TS` integration.
1289: .seealso: [](ch_ts), `TS`, `TSMonitoSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
1290: @*/
1291: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1292: {
1293: TSMonitorSPCtx ctx = (TSMonitorSPCtx)dctx;
1294: PetscDraw draw;
1295: DM dm, cdm;
1296: const PetscScalar *yy;
1297: PetscInt Np, p, dim = 2, *species;
1298: PetscReal species_color;
1300: PetscFunctionBegin;
1301: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1302: PetscCall(TSGetDM(ts, &dm));
1303: if (!step) {
1304: PetscDrawAxis axis;
1305: PetscReal dmboxlower[2], dmboxupper[2];
1307: PetscCall(TSGetDM(ts, &dm));
1308: PetscCall(DMGetDimension(dm, &dim));
1309: PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
1310: PetscCall(DMSwarmGetCellDM(dm, &cdm));
1311: PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
1312: PetscCall(VecGetLocalSize(u, &Np));
1313: Np /= dim * 2;
1314: PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
1315: if (ctx->phase) {
1316: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
1317: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10));
1318: } else {
1319: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
1320: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
1321: }
1322: PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
1323: PetscCall(PetscDrawSPReset(ctx->sp));
1324: }
1325: if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species));
1326: PetscCall(VecGetLocalSize(u, &Np));
1327: Np /= dim * 2;
1328: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1329: PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
1330: if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
1331: PetscCall(PetscDrawFlush(draw));
1332: PetscCall(PetscDrawSPReset(ctx->sp));
1333: PetscCall(VecGetArrayRead(u, &yy));
1334: for (p = 0; p < Np; ++p) {
1335: PetscReal x, y;
1337: if (ctx->phase) {
1338: x = PetscRealPart(yy[p * dim * 2]);
1339: y = PetscRealPart(yy[p * dim * 2 + dim]);
1340: } else {
1341: x = PetscRealPart(yy[p * dim * 2]);
1342: y = PetscRealPart(yy[p * dim * 2 + 1]);
1343: }
1344: if (ctx->multispecies) {
1345: species_color = species[p] + 2;
1346: PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color));
1347: } else {
1348: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1349: }
1350: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1351: }
1352: PetscCall(VecRestoreArrayRead(u, &yy));
1353: PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
1354: PetscCall(PetscDrawSPSave(ctx->sp));
1355: if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species));
1356: }
1357: PetscFunctionReturn(PETSC_SUCCESS);
1358: }
1360: /*@C
1361: TSMonitorHGSwarmSolution - Graphically displays histograms of `DMSWARM` particles
1363: Input Parameters:
1364: + ts - the `TS` context
1365: . step - current time-step
1366: . ptime - current time
1367: . u - current solution
1368: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorHGCtxCreate()`
1370: Options Database Keys:
1371: + -ts_monitor_hg_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1372: . -ts_monitor_hg_swarm_species <num> - Number of species to histogram
1373: . -ts_monitor_hg_swarm_bins <num> - Number of histogram bins
1374: - -ts_monitor_hg_swarm_velocity <bool> - Plot in velocity space, as opposed to coordinate space
1376: Level: intermediate
1378: Note:
1379: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1380: to be used during the `TS` integration.
1382: .seealso: `TSMonitoSet()`
1383: @*/
1384: PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1385: {
1386: TSMonitorHGCtx ctx = (TSMonitorHGCtx)dctx;
1387: PetscDraw draw;
1388: DM sw;
1389: const PetscScalar *yy;
1390: PetscInt *species;
1391: PetscInt dim, d = 0, Np, p, Ns, s;
1393: PetscFunctionBegin;
1394: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1395: PetscCall(TSGetDM(ts, &sw));
1396: PetscCall(DMGetDimension(sw, &dim));
1397: PetscCall(DMSwarmGetNumSpecies(sw, &Ns));
1398: Ns = PetscMin(Ns, ctx->Ns);
1399: PetscCall(VecGetLocalSize(u, &Np));
1400: Np /= dim * 2;
1401: if (!step) {
1402: PetscDrawAxis axis;
1403: char title[PETSC_MAX_PATH_LEN];
1405: for (s = 0; s < Ns; ++s) {
1406: PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis));
1407: PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s));
1408: if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N"));
1409: else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N"));
1410: }
1411: }
1412: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1413: PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species));
1414: for (s = 0; s < Ns; ++s) {
1415: PetscCall(PetscDrawHGReset(ctx->hg[s]));
1416: PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw));
1417: PetscCall(PetscDrawClear(draw));
1418: PetscCall(PetscDrawFlush(draw));
1419: }
1420: PetscCall(VecGetArrayRead(u, &yy));
1421: for (p = 0; p < Np; ++p) {
1422: const PetscInt s = species[p] < Ns ? species[p] : 0;
1423: PetscReal v;
1425: if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]);
1426: else v = PetscRealPart(yy[p * dim * 2 + d]);
1427: PetscCall(PetscDrawHGAddValue(ctx->hg[s], v));
1428: }
1429: PetscCall(VecRestoreArrayRead(u, &yy));
1430: for (s = 0; s < Ns; ++s) {
1431: PetscCall(PetscDrawHGDraw(ctx->hg[s]));
1432: PetscCall(PetscDrawHGSave(ctx->hg[s]));
1433: }
1434: PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species));
1435: }
1436: PetscFunctionReturn(PETSC_SUCCESS);
1437: }
1439: /*@C
1440: TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep
1442: Collective
1444: Input Parameters:
1445: + ts - the `TS` context
1446: . step - current time-step
1447: . ptime - current time
1448: . u - current solution
1449: - vf - unused context
1451: Options Database Key:
1452: . -ts_monitor_error - create a graphical monitor of error history
1454: Level: intermediate
1456: Notes:
1457: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1458: to be used during the `TS` integration.
1460: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1462: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1463: @*/
1464: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1465: {
1466: DM dm;
1467: PetscDS ds = NULL;
1468: PetscInt Nf = -1, f;
1469: PetscBool flg;
1471: PetscFunctionBegin;
1472: PetscCall(TSGetDM(ts, &dm));
1473: if (dm) PetscCall(DMGetDS(dm, &ds));
1474: if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
1475: if (Nf <= 0) {
1476: Vec y;
1477: PetscReal nrm;
1479: PetscCall(VecDuplicate(u, &y));
1480: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1481: PetscCall(VecAXPY(y, -1.0, u));
1482: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
1483: if (flg) {
1484: PetscCall(VecNorm(y, NORM_2, &nrm));
1485: PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
1486: }
1487: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
1488: if (flg) PetscCall(VecView(y, vf->viewer));
1489: PetscCall(VecDestroy(&y));
1490: } else {
1491: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
1492: void **ctxs;
1493: Vec v;
1494: PetscReal ferrors[1];
1496: PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
1497: for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
1498: PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
1499: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04" PetscInt_FMT " time = %-8.4g \t L_2 Error: [", step, (double)ptime));
1500: for (f = 0; f < Nf; ++f) {
1501: if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
1502: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
1503: }
1504: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));
1506: PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));
1508: PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
1509: if (flg) {
1510: PetscCall(DMGetGlobalVector(dm, &v));
1511: PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
1512: PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
1513: PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
1514: PetscCall(DMRestoreGlobalVector(dm, &v));
1515: }
1516: PetscCall(PetscFree2(exactFuncs, ctxs));
1517: }
1518: PetscFunctionReturn(PETSC_SUCCESS);
1519: }
1521: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1522: {
1523: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1524: PetscReal x = ptime, y;
1525: PetscInt its;
1527: PetscFunctionBegin;
1528: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1529: if (!n) {
1530: PetscDrawAxis axis;
1531: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1532: PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
1533: PetscCall(PetscDrawLGReset(ctx->lg));
1534: ctx->snes_its = 0;
1535: }
1536: PetscCall(TSGetSNESIterations(ts, &its));
1537: y = its - ctx->snes_its;
1538: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1539: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1540: PetscCall(PetscDrawLGDraw(ctx->lg));
1541: PetscCall(PetscDrawLGSave(ctx->lg));
1542: }
1543: ctx->snes_its = its;
1544: PetscFunctionReturn(PETSC_SUCCESS);
1545: }
1547: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1548: {
1549: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1550: PetscReal x = ptime, y;
1551: PetscInt its;
1553: PetscFunctionBegin;
1554: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1555: if (!n) {
1556: PetscDrawAxis axis;
1557: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1558: PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
1559: PetscCall(PetscDrawLGReset(ctx->lg));
1560: ctx->ksp_its = 0;
1561: }
1562: PetscCall(TSGetKSPIterations(ts, &its));
1563: y = its - ctx->ksp_its;
1564: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1565: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1566: PetscCall(PetscDrawLGDraw(ctx->lg));
1567: PetscCall(PetscDrawLGSave(ctx->lg));
1568: }
1569: ctx->ksp_its = its;
1570: PetscFunctionReturn(PETSC_SUCCESS);
1571: }
1573: /*@C
1574: TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()`
1576: Collective
1578: Input Parameter:
1579: . ts - the `TS` solver object
1581: Output Parameter:
1582: . ctx - the context
1584: Level: intermediate
1586: .seealso: [](ch_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1587: @*/
1588: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1589: {
1590: PetscFunctionBegin;
1591: PetscCall(PetscNew(ctx));
1592: PetscFunctionReturn(PETSC_SUCCESS);
1593: }
1595: /*@C
1596: TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution
1598: Collective
1600: Input Parameters:
1601: + ts - the `TS` context
1602: . step - current time-step
1603: . ptime - current time
1604: . u - current solution
1605: - dctx - the envelope context
1607: Options Database Key:
1608: . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time
1610: Level: intermediate
1612: Notes:
1613: After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope
1615: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1616: to be used during the `TS` integration.
1618: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1619: @*/
1620: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1621: {
1622: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;
1624: PetscFunctionBegin;
1625: if (!ctx->max) {
1626: PetscCall(VecDuplicate(u, &ctx->max));
1627: PetscCall(VecDuplicate(u, &ctx->min));
1628: PetscCall(VecCopy(u, ctx->max));
1629: PetscCall(VecCopy(u, ctx->min));
1630: } else {
1631: PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
1632: PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
1633: }
1634: PetscFunctionReturn(PETSC_SUCCESS);
1635: }
1637: /*@C
1638: TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution
1640: Collective
1642: Input Parameter:
1643: . ts - the `TS` context
1645: Output Parameters:
1646: + max - the maximum values
1647: - min - the minimum values
1649: Level: intermediate
1651: Notes:
1652: If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored
1654: .seealso: [](ch_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1655: @*/
1656: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1657: {
1658: PetscInt i;
1660: PetscFunctionBegin;
1661: if (max) *max = NULL;
1662: if (min) *min = NULL;
1663: for (i = 0; i < ts->numbermonitors; i++) {
1664: if (ts->monitor[i] == TSMonitorEnvelope) {
1665: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1666: if (max) *max = ctx->max;
1667: if (min) *min = ctx->min;
1668: break;
1669: }
1670: }
1671: PetscFunctionReturn(PETSC_SUCCESS);
1672: }
1674: /*@C
1675: TSMonitorEnvelopeCtxDestroy - Destroys a context that was created with `TSMonitorEnvelopeCtxCreate()`.
1677: Collective
1679: Input Parameter:
1680: . ctx - the monitor context
1682: Level: intermediate
1684: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1685: @*/
1686: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1687: {
1688: PetscFunctionBegin;
1689: PetscCall(VecDestroy(&(*ctx)->min));
1690: PetscCall(VecDestroy(&(*ctx)->max));
1691: PetscCall(PetscFree(*ctx));
1692: PetscFunctionReturn(PETSC_SUCCESS);
1693: }
1695: /*@C
1696: TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS`
1698: Not Collective
1700: Input Parameters:
1701: + ts - the `TS` context
1702: . step - current timestep
1703: . t - current time
1704: . U - current solution
1705: - vf - not used
1707: Options Database Key:
1708: . -ts_dmswarm_monitor_moments - Monitor moments of particle distribution
1710: Level: intermediate
1712: Notes:
1713: This requires a `DMSWARM` be attached to the `TS`.
1715: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1716: to be used during the TS integration.
1718: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1719: @*/
1720: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1721: {
1722: DM sw;
1723: const PetscScalar *u;
1724: PetscReal m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1725: PetscInt dim, d, Np, p;
1726: MPI_Comm comm;
1728: PetscFunctionBeginUser;
1729: (void)t;
1730: (void)vf;
1731: PetscCall(TSGetDM(ts, &sw));
1732: if (!sw || step % ts->monitorFrequency != 0) PetscFunctionReturn(PETSC_SUCCESS);
1733: PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
1734: PetscCall(DMGetDimension(sw, &dim));
1735: PetscCall(VecGetLocalSize(U, &Np));
1736: Np /= dim;
1737: PetscCall(VecGetArrayRead(U, &u));
1738: for (p = 0; p < Np; ++p) {
1739: for (d = 0; d < dim; ++d) {
1740: totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1741: totMom[d] += PetscRealPart(u[p * dim + d]);
1742: }
1743: }
1744: PetscCall(VecRestoreArrayRead(U, &u));
1745: for (d = 0; d < dim; ++d) totMom[d] *= m;
1746: totE *= 0.5 * m;
1747: PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
1748: for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, " Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
1749: PetscCall(PetscPrintf(comm, "\n"));
1750: PetscFunctionReturn(PETSC_SUCCESS);
1751: }