Actual source code: tsmon.c
1: #include <petsc/private/tsimpl.h>
2: #include <petscdm.h>
3: #include <petscds.h>
4: #include <petscdmswarm.h>
5: #include <petscdraw.h>
7: /*@C
8: TSMonitor - Runs all user-provided monitor routines set using `TSMonitorSet()`
10: Collective
12: Input Parameters:
13: + ts - time stepping context obtained from `TSCreate()`
14: . step - step number that has just completed
15: . ptime - model time of the state
16: - u - state at the current model time
18: Level: developer
20: Notes:
21: `TSMonitor()` is typically used automatically within the time stepping implementations.
22: Users would almost never call this routine directly.
24: A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions
26: .seealso: `TS`, `TSMonitorSet()`, `TSMonitorSetFromOptions()`
27: @*/
28: PetscErrorCode TSMonitor(TS ts, PetscInt step, PetscReal ptime, Vec u)
29: {
30: DM dm;
31: PetscInt i, n = ts->numbermonitors;
33: PetscFunctionBegin;
37: PetscCall(TSGetDM(ts, &dm));
38: PetscCall(DMSetOutputSequenceNumber(dm, step, ptime));
40: PetscCall(VecLockReadPush(u));
41: for (i = 0; i < n; i++) PetscCall((*ts->monitor[i])(ts, step, ptime, u, ts->monitorcontext[i]));
42: PetscCall(VecLockReadPop(u));
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: /*@C
47: TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
49: Collective
51: Input Parameters:
52: + ts - `TS` object you wish to monitor
53: . name - the monitor type one is seeking
54: . help - message indicating what monitoring is done
55: . manual - manual page for the monitor
56: . monitor - the monitor function, this must use a `PetscViewerFormat` as its context
57: - monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the `TS` or `PetscViewer` objects
59: Calling sequence of `monitor`:
60: + ts - the `TS` to monitor
61: . step - the current time-step
62: . time - the current time
63: . u - the current solution
64: - vf - the `PetscViewer` and format to monitor with
66: Calling sequence of `monitorsetup`:
67: + ts - the `TS` to monitor
68: - vf - the `PetscViewer` and format to monitor with
70: Level: developer
72: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `PetscOptionsCreateViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
73: `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
74: `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`,
75: `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
76: `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
77: `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
78: `PetscOptionsFList()`, `PetscOptionsEList()`
79: @*/
80: PetscErrorCode TSMonitorSetFromOptions(TS ts, const char name[], const char help[], const char manual[], PetscErrorCode (*monitor)(TS ts, PetscInt step, PetscReal time, Vec u, PetscViewerAndFormat *vf), PetscErrorCode (*monitorsetup)(TS ts, PetscViewerAndFormat *vf))
81: {
82: PetscViewer viewer;
83: PetscViewerFormat format;
84: PetscBool flg;
86: PetscFunctionBegin;
87: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)ts), ((PetscObject)ts)->options, ((PetscObject)ts)->prefix, name, &viewer, &format, &flg));
88: if (flg) {
89: PetscViewerAndFormat *vf;
90: char interval_key[1024];
92: PetscCall(PetscSNPrintf(interval_key, sizeof interval_key, "%s_interval", name));
93: PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
94: vf->view_interval = 1;
95: PetscCall(PetscOptionsGetInt(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, interval_key, &vf->view_interval, NULL));
97: PetscCall(PetscViewerDestroy(&viewer));
98: if (monitorsetup) PetscCall((*monitorsetup)(ts, vf));
99: PetscCall(TSMonitorSet(ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, PetscCtx))monitor, vf, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy));
100: }
101: PetscFunctionReturn(PETSC_SUCCESS);
102: }
104: /*@C
105: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
106: timestep to display the iteration's progress.
108: Logically Collective
110: Input Parameters:
111: + ts - the `TS` context obtained from `TSCreate()`
112: . monitor - monitoring routine
113: . mctx - [optional] user-defined context for private data for the monitor routine (use `NULL` if no context is desired)
114: - mdestroy - [optional] routine that frees monitor context (may be `NULL`), see `PetscCtxDestroyFn` for the calling sequence
116: Calling sequence of `monitor`:
117: + ts - the `TS` context
118: . steps - iteration number (after the final time step the monitor routine may be called with a step of -1, this indicates the solution has been interpolated to this time)
119: . time - current time
120: . u - current iterate
121: - ctx - [optional] monitoring context
123: Level: intermediate
125: Note:
126: This routine adds an additional monitor to the list of monitors that already has been loaded.
128: Fortran Notes:
129: Only a single monitor function can be set for each `TS` object
131: .seealso: [](ch_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
132: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
133: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `PetscCtxDestroyFn`
134: @*/
135: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS ts, PetscInt steps, PetscReal time, Vec u, PetscCtx ctx), PetscCtx mctx, PetscCtxDestroyFn *mdestroy)
136: {
137: PetscFunctionBegin;
139: for (PetscInt i = 0; i < ts->numbermonitors; i++) {
140: PetscBool identical;
142: PetscCall(PetscMonitorCompare((PetscErrorCode (*)(void))(PetscVoidFn *)monitor, mctx, mdestroy, (PetscErrorCode (*)(void))(PetscVoidFn *)ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical));
143: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
144: }
145: PetscCheck(ts->numbermonitors < MAXTSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set");
146: ts->monitor[ts->numbermonitors] = monitor;
147: ts->monitordestroy[ts->numbermonitors] = mdestroy;
148: ts->monitorcontext[ts->numbermonitors++] = mctx;
149: PetscFunctionReturn(PETSC_SUCCESS);
150: }
152: /*@C
153: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
155: Logically Collective
157: Input Parameter:
158: . ts - the `TS` context obtained from `TSCreate()`
160: Level: intermediate
162: Note:
163: There is no way to remove a single, specific monitor.
165: .seealso: [](ch_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
166: @*/
167: PetscErrorCode TSMonitorCancel(TS ts)
168: {
169: PetscInt i;
171: PetscFunctionBegin;
173: for (i = 0; i < ts->numbermonitors; i++) {
174: if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i]));
175: }
176: ts->numbermonitors = 0;
177: PetscFunctionReturn(PETSC_SUCCESS);
178: }
180: /*@C
181: TSMonitorDefault - The default monitor, prints the timestep and time for each step
183: Input Parameters:
184: + ts - the `TS` context
185: . 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)
186: . ptime - current time
187: . v - current iterate
188: - vf - the viewer and format
190: Options Database Key:
191: . -ts_monitor - monitors the time integration
193: Level: intermediate
195: Notes:
196: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
197: to be used during the `TS` integration.
199: .seealso: [](ch_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorWallClockTime()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
200: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
201: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`
202: @*/
203: PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
204: {
205: PetscViewer viewer = vf->viewer;
206: PetscBool isascii, ibinary;
208: PetscFunctionBegin;
210: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
211: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary));
212: PetscCall(PetscViewerPushFormat(viewer, vf->format));
213: if (isascii) {
214: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
215: if (step == -1) { /* this indicates it is an interpolated solution */
216: PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps));
217: } else {
218: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n"));
219: }
220: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
221: } else if (ibinary) {
222: PetscMPIInt rank;
223: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
224: if (rank == 0) {
225: PetscBool skipHeader;
226: PetscInt classid = REAL_FILE_CLASSID;
228: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
229: if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
230: PetscCall(PetscRealView(1, &ptime, viewer));
231: } else {
232: PetscCall(PetscRealView(0, &ptime, viewer));
233: }
234: }
235: PetscCall(PetscViewerPopFormat(viewer));
236: PetscFunctionReturn(PETSC_SUCCESS);
237: }
239: typedef struct {
240: PetscLogDouble time_start;
241: PetscLogDouble time_last;
242: PetscInt snes_its;
243: PetscInt ksp_its;
244: } *TSMonitorWallClockTimeContext;
246: /*@C
247: TSMonitorWallClockTimeSetUp - Setup routine passed to `TSMonitorSetFromOptions()` when using `-ts_monitor_wall_clock_time`
249: Input Parameters:
250: + ts - the `TS` context
251: - vf - the viewer and format
253: Level: intermediate
255: Note:
256: This is not called directly by users, rather one calls `TSMonitorSetFromOptions()`, with `TSMonitorWallClockTime()` and this function as arguments, to cause the monitor
257: to be used during the `TS` integration.
259: .seealso: [](ch_ts), `TSMonitorSet()`
260: @*/
261: PetscErrorCode TSMonitorWallClockTimeSetUp(TS ts, PetscViewerAndFormat *vf)
262: {
263: TSMonitorWallClockTimeContext speed;
265: PetscFunctionBegin;
266: PetscCall(PetscNew(&speed));
267: speed->time_start = PETSC_DECIDE;
268: vf->data_destroy = PetscCtxDestroyDefault;
269: vf->data = speed;
270: PetscFunctionReturn(PETSC_SUCCESS);
271: }
273: /*@C
274: TSMonitorWallClockTime - Monitor wall-clock time, KSP iterations, and SNES iterations per step.
276: Input Parameters:
277: + ts - the `TS` context
278: . 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)
279: . ptime - current time
280: . v - current solution
281: - vf - the viewer and format
283: Options Database Key:
284: . -ts_monitor_wall_clock_time - Monitor wall-clock time, KSP iterations, and SNES iterations per step.
286: Level: intermediate
288: Note:
289: This is not called directly by users, rather one calls `TSMonitorSetFromOptions()`, with this function and `TSMonitorWallClockTimeSetUp()` as arguments, to cause the monitor
290: to be used during the `TS` integration.
292: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
293: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
294: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
295: @*/
296: PetscErrorCode TSMonitorWallClockTime(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
297: {
298: PetscViewer viewer = vf->viewer;
299: TSMonitorWallClockTimeContext speed = (TSMonitorWallClockTimeContext)vf->data;
300: PetscBool isascii;
301: PetscLogDouble now;
302: PetscInt snes_its, ksp_its;
304: PetscFunctionBegin;
306: PetscCall(PetscTime(&now));
307: if (speed->time_start == PETSC_DECIDE) {
308: speed->time_start = now;
309: speed->time_last = now;
310: }
311: PetscCall(TSGetSNESIterations(ts, &snes_its));
312: PetscCall(TSGetKSPIterations(ts, &ksp_its));
313: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
314: PetscCall(PetscViewerPushFormat(viewer, vf->format));
315: if (isascii) {
316: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
317: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s elapsed %.6f of %.6f snes %" PetscInt_FMT " ksp %" PetscInt_FMT "\n", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)" : "", now - speed->time_last,
318: now - speed->time_start, snes_its - speed->snes_its, ksp_its - speed->ksp_its));
319: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
320: }
321: PetscCall(PetscViewerPopFormat(viewer));
322: speed->time_last = now;
323: speed->snes_its = snes_its;
324: speed->ksp_its = ksp_its;
325: PetscFunctionReturn(PETSC_SUCCESS);
326: }
328: /*@C
329: TSMonitorExtreme - Prints the extreme values of the solution at each timestep
331: Input Parameters:
332: + ts - the `TS` context
333: . 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)
334: . ptime - current time
335: . v - current iterate
336: - vf - the viewer and format
338: Level: intermediate
340: Note:
341: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
342: to be used during the `TS` integration.
344: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`
345: @*/
346: PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
347: {
348: PetscViewer viewer = vf->viewer;
349: PetscBool isascii;
350: PetscReal max, min;
352: PetscFunctionBegin;
354: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
355: PetscCall(PetscViewerPushFormat(viewer, vf->format));
356: if (isascii) {
357: PetscCall(VecMax(v, NULL, &max));
358: PetscCall(VecMin(v, NULL, &min));
359: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
360: 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));
361: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
362: }
363: PetscCall(PetscViewerPopFormat(viewer));
364: PetscFunctionReturn(PETSC_SUCCESS);
365: }
367: /*@C
368: TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with
369: `TS` to monitor the solution process graphically in various ways
371: Collective
373: Input Parameters:
374: + comm - the MPI communicator to use
375: . host - the X display to open, or `NULL` for the local machine
376: . label - the title to put in the title bar
377: . x - the x screen coordinates of the upper left coordinate of the window
378: . y - the y screen coordinates of the upper left coordinate of the window
379: . m - the screen width in pixels
380: . n - the screen height in pixels
381: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
383: Output Parameter:
384: . ctx - the context
386: Options Database Keys:
387: + -ts_monitor_lg_timestep - automatically sets line graph monitor
388: . -ts_monitor_lg_timestep_log - automatically sets line graph monitor
389: . -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`)
390: . -ts_monitor_lg_error - monitor the error
391: . -ts_monitor_lg_ksp_iterations - monitor the number of `KSP` iterations needed for each timestep
392: . -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep
393: - -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true
395: Level: intermediate
397: Notes:
398: Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed.
400: One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()`
402: Many of the functions that control the monitoring have two forms\: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a `TS` object as the
403: 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
404: as the first argument.
406: One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()`
408: .seealso: [](ch_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
409: `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
410: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
411: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
412: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
413: @*/
414: PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
415: {
416: PetscDraw draw;
418: PetscFunctionBegin;
419: PetscCall(PetscNew(ctx));
420: PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
421: PetscCall(PetscDrawSetFromOptions(draw));
422: PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg));
423: PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg));
424: PetscCall(PetscDrawDestroy(&draw));
425: (*ctx)->howoften = howoften;
426: PetscFunctionReturn(PETSC_SUCCESS);
427: }
429: /*@C
430: TSMonitorLGTimeStep - Monitors a `TS` by printing the time-steps
432: Collective
434: Input Parameters:
435: + ts - the time integrator
436: . step - the current time step
437: . ptime - the current time
438: . v - the current state
439: - monctx - the monitor context obtained with `TSMonitorLGCtxCreate()`
441: Level: advanced
443: Note:
444: This is not called directly by users, rather one calls `TSMonitorSet()` along the `ctx` created by `TSMonitorLGCtxCreate()`
445: and `TSMonitorLGCtxDestroy()`
447: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGCtxDestroy()`
448: @*/
449: PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscCtx monctx)
450: {
451: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
452: PetscReal x = ptime, y;
454: PetscFunctionBegin;
455: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */
456: if (!step) {
457: PetscDrawAxis axis;
458: const char *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
459: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
460: PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel));
461: PetscCall(PetscDrawLGReset(ctx->lg));
462: }
463: PetscCall(TSGetTimeStep(ts, &y));
464: if (ctx->semilogy) y = PetscLog10Real(y);
465: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
466: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
467: PetscCall(PetscDrawLGDraw(ctx->lg));
468: PetscCall(PetscDrawLGSave(ctx->lg));
469: }
470: PetscFunctionReturn(PETSC_SUCCESS);
471: }
473: /*@C
474: TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`.
476: Collective
478: Input Parameter:
479: . ctx - the monitor context
481: Level: intermediate
483: Note:
484: Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`
486: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
487: @*/
488: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
489: {
490: PetscFunctionBegin;
491: if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)(&(*ctx)->transformctx));
492: PetscCall(PetscDrawLGDestroy(&(*ctx)->lg));
493: PetscCall(PetscStrArrayDestroy(&(*ctx)->names));
494: PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames));
495: PetscCall(PetscFree((*ctx)->displayvariables));
496: PetscCall(PetscFree((*ctx)->displayvalues));
497: PetscCall(PetscFree(*ctx));
498: PetscFunctionReturn(PETSC_SUCCESS);
499: }
501: /* Creates a TSMonitorSPCtx for use with DMSwarm particle visualizations */
502: 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)
503: {
504: PetscDraw draw;
506: PetscFunctionBegin;
507: PetscCall(PetscNew(ctx));
508: PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
509: PetscCall(PetscDrawSetFromOptions(draw));
510: PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp));
511: PetscCall(PetscDrawDestroy(&draw));
512: (*ctx)->howoften = howoften;
513: (*ctx)->retain = retain;
514: (*ctx)->phase = phase;
515: (*ctx)->multispecies = multispecies;
516: PetscFunctionReturn(PETSC_SUCCESS);
517: }
519: /* Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate */
520: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
521: {
522: PetscFunctionBegin;
523: PetscCall(PetscDrawSPDestroy(&(*ctx)->sp));
524: PetscCall(PetscFree(*ctx));
525: PetscFunctionReturn(PETSC_SUCCESS);
526: }
528: /* Creates a TSMonitorHGCtx for use with DMSwarm particle visualizations */
529: 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)
530: {
531: PetscDraw draw;
532: int Nsi, Nbi;
534: PetscFunctionBegin;
535: PetscCall(PetscMPIIntCast(Ns, &Nsi));
536: PetscCall(PetscMPIIntCast(Nb, &Nbi));
537: PetscCall(PetscNew(ctx));
538: PetscCall(PetscMalloc1(Ns, &(*ctx)->hg));
539: for (int s = 0; s < Nsi; ++s) {
540: PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw));
541: PetscCall(PetscDrawSetFromOptions(draw));
542: PetscCall(PetscDrawHGCreate(draw, Nbi, &(*ctx)->hg[s]));
543: PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE));
544: PetscCall(PetscDrawDestroy(&draw));
545: }
546: (*ctx)->howoften = howoften;
547: (*ctx)->Ns = Ns;
548: (*ctx)->velocity = velocity;
549: PetscFunctionReturn(PETSC_SUCCESS);
550: }
552: /* Destroys a TSMonitorHGCtx that was created with TSMonitorHGCtxCreate */
553: PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx)
554: {
555: PetscInt s;
557: PetscFunctionBegin;
558: for (s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s]));
559: PetscCall(PetscFree((*ctx)->hg));
560: PetscCall(PetscFree(*ctx));
561: PetscFunctionReturn(PETSC_SUCCESS);
562: }
564: /*@C
565: TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling
566: `VecView()` for the solution at each timestep
568: Collective
570: Input Parameters:
571: + ts - the `TS` context
572: . step - current time-step
573: . ptime - current time
574: . u - the solution at the current time
575: - ctx - either a viewer or `NULL`
577: Options Database Keys:
578: + -ts_monitor_draw_solution - draw the solution at each time-step
579: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
581: Level: intermediate
583: Notes:
584: The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial`
585: will look bad
587: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with
588: `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration.
590: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
591: @*/
592: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx ctx)
593: {
594: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)ctx;
595: PetscDraw draw;
597: PetscFunctionBegin;
598: if (!step && ictx->showinitial) {
599: if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
600: PetscCall(VecCopy(u, ictx->initialsolution));
601: }
602: if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
604: if (ictx->showinitial) {
605: PetscReal pause;
606: PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
607: PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
608: PetscCall(VecView(ictx->initialsolution, ictx->viewer));
609: PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
610: PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
611: }
612: PetscCall(VecView(u, ictx->viewer));
613: if (ictx->showtimestepandtime) {
614: PetscReal xl, yl, xr, yr, h;
615: char time[32];
617: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
618: PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
619: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
620: h = yl + .95 * (yr - yl);
621: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
622: PetscCall(PetscDrawFlush(draw));
623: }
625: if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
626: PetscFunctionReturn(PETSC_SUCCESS);
627: }
629: /*@C
630: TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram
632: Collective
634: Input Parameters:
635: + ts - the `TS` context
636: . step - current time-step
637: . ptime - current time
638: . u - the solution at the current time
639: - ctx - either a viewer or `NULL`
641: Level: intermediate
643: Notes:
644: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
645: to be used during the `TS` integration.
647: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
648: @*/
649: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx ctx)
650: {
651: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)ctx;
652: PetscDraw draw;
653: PetscDrawAxis axis;
654: PetscInt n;
655: PetscMPIInt size;
656: PetscReal U0, U1, xl, yl, xr, yr, h;
657: char time[32];
658: const PetscScalar *U;
660: PetscFunctionBegin;
661: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
662: PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
663: PetscCall(VecGetSize(u, &n));
664: PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");
666: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
667: PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
668: PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
669: if (!step) {
670: PetscCall(PetscDrawClear(draw));
671: PetscCall(PetscDrawAxisDraw(axis));
672: }
674: PetscCall(VecGetArrayRead(u, &U));
675: U0 = PetscRealPart(U[0]);
676: U1 = PetscRealPart(U[1]);
677: PetscCall(VecRestoreArrayRead(u, &U));
678: if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);
680: PetscDrawCollectiveBegin(draw);
681: PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
682: if (ictx->showtimestepandtime) {
683: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
684: PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
685: h = yl + .95 * (yr - yl);
686: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
687: }
688: PetscDrawCollectiveEnd(draw);
689: PetscCall(PetscDrawFlush(draw));
690: PetscCall(PetscDrawPause(draw));
691: PetscCall(PetscDrawSave(draw));
692: PetscFunctionReturn(PETSC_SUCCESS);
693: }
695: /*@C
696: TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()`
698: Collective
700: Input Parameter:
701: . ictx - the monitor context
703: Level: intermediate
705: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
706: @*/
707: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
708: {
709: PetscFunctionBegin;
710: PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
711: PetscCall(VecDestroy(&(*ictx)->initialsolution));
712: PetscCall(PetscFree(*ictx));
713: PetscFunctionReturn(PETSC_SUCCESS);
714: }
716: /*@C
717: TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`
719: Collective
721: Input Parameters:
722: + comm - the MPI communicator to use
723: . host - the X display to open, or `NULL` for the local machine
724: . label - the title to put in the title bar
725: . x - the x screen coordinates of the upper left coordinate of the window
726: . y - the y screen coordinates of the upper left coordinate of the window
727: . m - the screen width in pixels
728: . n - the screen height in pixels
729: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
731: Output Parameter:
732: . ctx - the monitor context
734: Options Database Keys:
735: + -ts_monitor_draw_solution - draw the solution at each time-step
736: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
738: Level: intermediate
740: Note:
741: The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`.
743: .seealso: [](ch_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
744: @*/
745: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
746: {
747: PetscFunctionBegin;
748: PetscCall(PetscNew(ctx));
749: PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
750: PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));
752: (*ctx)->howoften = howoften;
753: (*ctx)->showinitial = PETSC_FALSE;
754: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));
756: (*ctx)->showtimestepandtime = PETSC_FALSE;
757: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
758: PetscFunctionReturn(PETSC_SUCCESS);
759: }
761: /*@C
762: TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling
763: `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep
765: Collective
767: Input Parameters:
768: + ts - the `TS` context
769: . step - current time-step
770: . ptime - current time
771: . u - solution at current time
772: - Ctx - either a viewer or `NULL`
774: Options Database Key:
775: . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
777: Level: intermediate
779: Note:
780: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
781: to be used during the `TS` integration.
783: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
784: @*/
785: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx Ctx)
786: {
787: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)Ctx;
788: PetscViewer viewer = ctx->viewer;
789: Vec work;
791: PetscFunctionBegin;
792: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
793: PetscCall(VecDuplicate(u, &work));
794: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
795: PetscCall(VecView(work, viewer));
796: PetscCall(VecDestroy(&work));
797: PetscFunctionReturn(PETSC_SUCCESS);
798: }
800: /*@C
801: TSMonitorDrawError - Monitors progress of the `TS` solvers by calling
802: `VecView()` for the error at each timestep
804: Collective
806: Input Parameters:
807: + ts - the `TS` context
808: . step - current time-step
809: . ptime - current time
810: . u - solution at current time
811: - Ctx - either a viewer or `NULL`
813: Options Database Key:
814: . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
816: Level: intermediate
818: Notes:
819: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
820: to be used during the `TS` integration.
822: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
823: @*/
824: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx Ctx)
825: {
826: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)Ctx;
827: PetscViewer viewer = ctx->viewer;
828: Vec work;
830: PetscFunctionBegin;
831: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
832: PetscCall(VecDuplicate(u, &work));
833: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
834: PetscCall(VecAXPY(work, -1.0, u));
835: PetscCall(VecView(work, viewer));
836: PetscCall(VecDestroy(&work));
837: PetscFunctionReturn(PETSC_SUCCESS);
838: }
840: /*@C
841: TSMonitorSolutionSetup - Setups the context for `TSMonitorSolution()`
843: Collective
845: Input Parameters:
846: + ts - the `TS` context
847: - vf - viewer and its format
849: Level: intermediate
851: .seealso: [](ch_ts), `TS`, `TSMonitorSolution()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSetFromOptions()`
852: @*/
853: PetscErrorCode TSMonitorSolutionSetup(TS ts, PetscViewerAndFormat *vf)
854: {
855: TSMonitorSolutionCtx ctx;
857: PetscFunctionBegin;
858: PetscCall(PetscNew(&ctx));
859: PetscCall(PetscOptionsGetBool(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_solution_skip_initial", &ctx->skip_initial, NULL));
860: vf->data = ctx;
861: vf->data_destroy = PetscCtxDestroyDefault;
862: PetscFunctionReturn(PETSC_SUCCESS);
863: }
865: /*@C
866: 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
868: Collective
870: Input Parameters:
871: + ts - the `TS` context
872: . step - current time-step
873: . ptime - current time
874: . u - current state
875: - vf - viewer and its format
877: Level: intermediate
879: Notes:
880: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
881: to be used during the `TS` integration.
883: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSolutionSetup()`,
884: @*/
885: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
886: {
887: TSMonitorSolutionCtx ctx = (TSMonitorSolutionCtx)vf->data;
889: PetscFunctionBegin;
890: if (ctx->skip_initial && step == ts->start_step) PetscFunctionReturn(PETSC_SUCCESS);
891: if ((vf->view_interval > 0 && !(step % vf->view_interval)) || (vf->view_interval && ts->reason)) {
892: PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
893: PetscCall(VecView(u, vf->viewer));
894: PetscCall(PetscViewerPopFormat(vf->viewer));
895: }
896: PetscFunctionReturn(PETSC_SUCCESS);
897: }
899: /*@C
900: TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at selected timesteps.
902: Collective
904: Input Parameters:
905: + ts - the `TS` context
906: . step - current time-step
907: . ptime - current time
908: . u - current state
909: - ctx - monitor context obtained with `TSMonitorSolutionVTKCtxCreate()`
911: Level: developer
913: Notes:
914: 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.
915: These are named according to the file name template.
917: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
918: to be used during the `TS` integration.
920: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
921: @*/
922: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, TSMonitorVTKCtx ctx)
923: {
924: char filename[PETSC_MAX_PATH_LEN];
925: PetscViewer viewer;
927: PetscFunctionBegin;
928: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
929: if (((ctx->interval > 0) && (!(step % ctx->interval))) || (ctx->interval && ts->reason)) {
930: PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)ctx->filenametemplate, step));
931: PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
932: PetscCall(VecView(u, viewer));
933: PetscCall(PetscViewerDestroy(&viewer));
934: }
935: PetscFunctionReturn(PETSC_SUCCESS);
936: }
938: /*@C
939: TSMonitorSolutionVTKDestroy - Destroy the monitor context created with `TSMonitorSolutionVTKCtxCreate()`
941: Not Collective
943: Input Parameter:
944: . ctx - the monitor context
946: Level: developer
948: Note:
949: This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
951: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
952: @*/
953: PetscErrorCode TSMonitorSolutionVTKDestroy(TSMonitorVTKCtx *ctx)
954: {
955: PetscFunctionBegin;
956: PetscCall(PetscFree((*ctx)->filenametemplate));
957: PetscCall(PetscFree(*ctx));
958: PetscFunctionReturn(PETSC_SUCCESS);
959: }
961: /*@C
962: TSMonitorSolutionVTKCtxCreate - Create the monitor context to be used in `TSMonitorSolutionVTK()`
964: Not collective
966: Input Parameter:
967: . filenametemplate - the template file name, e.g. foo-%03d.vts
969: Output Parameter:
970: . ctx - the monitor context
972: Level: developer
974: Note:
975: This function is normally used inside `TSSetFromOptions()` to pass the context created to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
977: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`, `TSMonitorSolutionVTKDestroy()`
978: @*/
979: PetscErrorCode TSMonitorSolutionVTKCtxCreate(const char *filenametemplate, TSMonitorVTKCtx *ctx)
980: {
981: const char *ptr = NULL, *ptr2 = NULL;
982: TSMonitorVTKCtx ictx;
984: PetscFunctionBegin;
985: PetscAssertPointer(filenametemplate, 1);
986: PetscAssertPointer(ctx, 2);
987: /* Do some cursory validation of the input. */
988: PetscCall(PetscStrstr(filenametemplate, "%", (char **)&ptr));
989: PetscCheck(ptr, PETSC_COMM_SELF, PETSC_ERR_USER, "-ts_monitor_solution_vtk requires a file template, e.g. filename-%%03" PetscInt_FMT ".vts");
990: for (ptr++; ptr && *ptr; ptr++) {
991: PetscCall(PetscStrchr("DdiouxX", *ptr, (char **)&ptr2));
992: 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");
993: if (ptr2) break;
994: }
995: PetscCall(PetscNew(&ictx));
996: PetscCall(PetscStrallocpy(filenametemplate, &ictx->filenametemplate));
997: ictx->interval = 1;
999: *ctx = ictx;
1000: PetscFunctionReturn(PETSC_SUCCESS);
1001: }
1003: /*@C
1004: TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector
1005: in a time based line graph
1007: Collective
1009: Input Parameters:
1010: + ts - the `TS` context
1011: . step - current time-step
1012: . ptime - current time
1013: . u - current solution
1014: - dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`
1016: Options Database Key:
1017: . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables
1019: Level: intermediate
1021: Notes:
1022: Each process in a parallel run displays its component solutions in a separate window
1024: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1025: to be used during the `TS` integration.
1027: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
1028: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
1029: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
1030: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
1031: @*/
1032: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1033: {
1034: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dctx;
1035: const PetscScalar *yy;
1036: Vec v;
1038: PetscFunctionBegin;
1039: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1040: if (!step) {
1041: PetscDrawAxis axis;
1042: PetscInt dim;
1043: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1044: PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
1045: if (!ctx->names) {
1046: PetscBool flg;
1047: /* user provides names of variables to plot but no names has been set so assume names are integer values */
1048: PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
1049: if (flg) {
1050: PetscInt i, n;
1051: char **names;
1052: PetscCall(VecGetSize(u, &n));
1053: PetscCall(PetscMalloc1(n + 1, &names));
1054: for (i = 0; i < n; i++) {
1055: PetscCall(PetscMalloc1(5, &names[i]));
1056: PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
1057: }
1058: names[n] = NULL;
1059: ctx->names = names;
1060: }
1061: }
1062: if (ctx->names && !ctx->displaynames) {
1063: char **displaynames;
1064: PetscBool flg;
1065: PetscCall(VecGetLocalSize(u, &dim));
1066: PetscCall(PetscCalloc1(dim + 1, &displaynames));
1067: PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
1068: if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
1069: PetscCall(PetscStrArrayDestroy(&displaynames));
1070: }
1071: if (ctx->displaynames) {
1072: PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
1073: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
1074: } else if (ctx->names) {
1075: PetscCall(VecGetLocalSize(u, &dim));
1076: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1077: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
1078: } else {
1079: PetscCall(VecGetLocalSize(u, &dim));
1080: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1081: }
1082: PetscCall(PetscDrawLGReset(ctx->lg));
1083: }
1085: if (!ctx->transform) v = u;
1086: else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
1087: PetscCall(VecGetArrayRead(v, &yy));
1088: if (ctx->displaynames) {
1089: PetscInt i;
1090: for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
1091: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
1092: } else {
1093: #if defined(PETSC_USE_COMPLEX)
1094: PetscInt i, n;
1095: PetscReal *yreal;
1096: PetscCall(VecGetLocalSize(v, &n));
1097: PetscCall(PetscMalloc1(n, &yreal));
1098: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1099: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1100: PetscCall(PetscFree(yreal));
1101: #else
1102: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1103: #endif
1104: }
1105: PetscCall(VecRestoreArrayRead(v, &yy));
1106: if (ctx->transform) PetscCall(VecDestroy(&v));
1108: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1109: PetscCall(PetscDrawLGDraw(ctx->lg));
1110: PetscCall(PetscDrawLGSave(ctx->lg));
1111: }
1112: PetscFunctionReturn(PETSC_SUCCESS);
1113: }
1115: /*@C
1116: TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
1118: Collective
1120: Input Parameters:
1121: + ts - the `TS` context
1122: - names - the names of the components, final string must be `NULL`
1124: Level: intermediate
1126: Notes:
1127: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1129: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
1130: @*/
1131: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
1132: {
1133: PetscInt i;
1135: PetscFunctionBegin;
1136: for (i = 0; i < ts->numbermonitors; i++) {
1137: if (ts->monitor[i] == TSMonitorLGSolution) {
1138: PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
1139: break;
1140: }
1141: }
1142: PetscFunctionReturn(PETSC_SUCCESS);
1143: }
1145: /*@C
1146: TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
1148: Collective
1150: Input Parameters:
1151: + ctx - the `TS` context
1152: - names - the names of the components, final string must be `NULL`
1154: Level: intermediate
1156: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
1157: @*/
1158: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
1159: {
1160: PetscFunctionBegin;
1161: PetscCall(PetscStrArrayDestroy(&ctx->names));
1162: PetscCall(PetscStrArrayallocpy(names, &ctx->names));
1163: PetscFunctionReturn(PETSC_SUCCESS);
1164: }
1166: /*@C
1167: TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot
1169: Collective
1171: Input Parameter:
1172: . ts - the `TS` context
1174: Output Parameter:
1175: . names - the names of the components, final string must be `NULL`
1177: Level: intermediate
1179: Note:
1180: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1182: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1183: @*/
1184: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
1185: {
1186: PetscInt i;
1188: PetscFunctionBegin;
1189: *names = NULL;
1190: for (i = 0; i < ts->numbermonitors; i++) {
1191: if (ts->monitor[i] == TSMonitorLGSolution) {
1192: TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
1193: *names = (const char *const *)ctx->names;
1194: break;
1195: }
1196: }
1197: PetscFunctionReturn(PETSC_SUCCESS);
1198: }
1200: /*@C
1201: TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor
1203: Collective
1205: Input Parameters:
1206: + ctx - the `TSMonitorLG` context
1207: - displaynames - the names of the components, final string must be `NULL`
1209: Level: intermediate
1211: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1212: @*/
1213: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
1214: {
1215: PetscInt j = 0, k;
1217: PetscFunctionBegin;
1218: if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
1219: PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
1220: PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
1221: while (displaynames[j]) j++;
1222: ctx->ndisplayvariables = j;
1223: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
1224: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
1225: j = 0;
1226: while (displaynames[j]) {
1227: k = 0;
1228: while (ctx->names[k]) {
1229: PetscBool flg;
1230: PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
1231: if (flg) {
1232: ctx->displayvariables[j] = k;
1233: break;
1234: }
1235: k++;
1236: }
1237: j++;
1238: }
1239: PetscFunctionReturn(PETSC_SUCCESS);
1240: }
1242: /*@C
1243: TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor
1245: Collective
1247: Input Parameters:
1248: + ts - the `TS` context
1249: - displaynames - the names of the components, final string must be `NULL`
1251: Level: intermediate
1253: Note:
1254: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1256: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1257: @*/
1258: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
1259: {
1260: PetscInt i;
1262: PetscFunctionBegin;
1263: for (i = 0; i < ts->numbermonitors; i++) {
1264: if (ts->monitor[i] == TSMonitorLGSolution) {
1265: PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
1266: break;
1267: }
1268: }
1269: PetscFunctionReturn(PETSC_SUCCESS);
1270: }
1272: /*@C
1273: TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed
1275: Collective
1277: Input Parameters:
1278: + ts - the `TS` context
1279: . transform - the transform function
1280: . destroy - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1281: - tctx - optional context used by transform function
1283: Calling sequence of `transform`:
1284: + tctx - context used by the transform function
1285: . u - the input solution vector
1286: - w - the output transformed vector
1288: Level: intermediate
1290: Note:
1291: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1293: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorLGCtxSetTransform()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `PetscCtxDestroyFn`
1294: @*/
1295: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(PetscCtx tctx, Vec u, Vec *w), PetscCtxDestroyFn *destroy, PetscCtx tctx)
1296: {
1297: PetscInt i;
1299: PetscFunctionBegin;
1300: for (i = 0; i < ts->numbermonitors; i++) {
1301: if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
1302: }
1303: PetscFunctionReturn(PETSC_SUCCESS);
1304: }
1306: /*@C
1307: TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed
1309: Collective
1311: Input Parameters:
1312: + tctx - the `TS` context
1313: . transform - the transform function
1314: . destroy - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1315: - ctx - optional context used by transform function
1317: Calling sequence of `transform`:
1318: + tctx - context used by the transform function
1319: . u - the input solution vector
1320: - w - the output transformed vector
1322: Level: intermediate
1324: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`, `PetscCtxDestroyFn`
1325: @*/
1326: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(PetscCtx tctx, Vec u, Vec *w), PetscCtxDestroyFn *destroy, PetscCtx tctx)
1327: {
1328: PetscFunctionBegin;
1329: ctx->transform = transform;
1330: ctx->transformdestroy = destroy;
1331: ctx->transformctx = tctx;
1332: PetscFunctionReturn(PETSC_SUCCESS);
1333: }
1335: /*@C
1336: TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error
1337: in a time based line graph
1339: Collective
1341: Input Parameters:
1342: + ts - the `TS` context
1343: . step - current time-step
1344: . ptime - current time
1345: . u - current solution
1346: - Ctx - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`
1348: Options Database Key:
1349: . -ts_monitor_lg_error - create a graphical monitor of error history
1351: Level: intermediate
1353: Notes:
1354: Each process in a parallel run displays its component errors in a separate window
1356: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1358: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1359: to be used during the TS integration.
1361: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1362: @*/
1363: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx Ctx)
1364: {
1365: TSMonitorLGCtx ctx = (TSMonitorLGCtx)Ctx;
1366: const PetscScalar *yy;
1367: Vec y;
1369: PetscFunctionBegin;
1370: if (!step) {
1371: PetscDrawAxis axis;
1372: PetscInt dim;
1373: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1374: PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
1375: PetscCall(VecGetLocalSize(u, &dim));
1376: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1377: PetscCall(PetscDrawLGReset(ctx->lg));
1378: }
1379: PetscCall(VecDuplicate(u, &y));
1380: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1381: PetscCall(VecAXPY(y, -1.0, u));
1382: PetscCall(VecGetArrayRead(y, &yy));
1383: #if defined(PETSC_USE_COMPLEX)
1384: {
1385: PetscReal *yreal;
1386: PetscInt i, n;
1387: PetscCall(VecGetLocalSize(y, &n));
1388: PetscCall(PetscMalloc1(n, &yreal));
1389: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1390: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1391: PetscCall(PetscFree(yreal));
1392: }
1393: #else
1394: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1395: #endif
1396: PetscCall(VecRestoreArrayRead(y, &yy));
1397: PetscCall(VecDestroy(&y));
1398: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1399: PetscCall(PetscDrawLGDraw(ctx->lg));
1400: PetscCall(PetscDrawLGSave(ctx->lg));
1401: }
1402: PetscFunctionReturn(PETSC_SUCCESS);
1403: }
1405: /*@C
1406: TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot
1408: Input Parameters:
1409: + ts - the `TS` context
1410: . step - current time-step
1411: . ptime - current time
1412: . u - current solution
1413: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()`
1415: Options Database Keys:
1416: + -ts_monitor_sp_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1417: . -ts_monitor_sp_swarm_retain <n> - Retain n old points so we can see the history, or -1 for all points
1418: . -ts_monitor_sp_swarm_multi_species <bool> - Color each species differently
1419: - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space
1421: Level: intermediate
1423: Notes:
1424: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1425: to be used during the `TS` integration.
1427: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
1428: @*/
1429: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx dctx)
1430: {
1431: TSMonitorSPCtx ctx = (TSMonitorSPCtx)dctx;
1432: PetscDraw draw;
1433: DM dm, cdm;
1434: const PetscScalar *yy;
1435: PetscInt Np, p, dim = 2, *species;
1436: PetscReal species_color;
1438: PetscFunctionBegin;
1439: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1440: PetscCall(TSGetDM(ts, &dm));
1441: if (!step) {
1442: PetscDrawAxis axis;
1443: PetscReal dmboxlower[2], dmboxupper[2];
1445: PetscCall(TSGetDM(ts, &dm));
1446: PetscCall(DMGetDimension(dm, &dim));
1447: PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
1448: PetscCall(DMSwarmGetCellDM(dm, &cdm));
1449: PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
1450: PetscCall(VecGetLocalSize(u, &Np));
1451: Np /= dim * 2;
1452: PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
1453: if (ctx->phase) {
1454: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
1455: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10));
1456: } else {
1457: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
1458: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
1459: }
1460: PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
1461: PetscCall(PetscDrawSPReset(ctx->sp));
1462: }
1463: if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species));
1464: PetscCall(VecGetLocalSize(u, &Np));
1465: Np /= dim * 2;
1466: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1467: PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
1468: if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
1469: PetscCall(PetscDrawFlush(draw));
1470: PetscCall(PetscDrawSPReset(ctx->sp));
1471: PetscCall(VecGetArrayRead(u, &yy));
1472: for (p = 0; p < Np; ++p) {
1473: PetscReal x, y;
1475: if (ctx->phase) {
1476: x = PetscRealPart(yy[p * dim * 2]);
1477: y = PetscRealPart(yy[p * dim * 2 + dim]);
1478: } else {
1479: x = PetscRealPart(yy[p * dim * 2]);
1480: y = PetscRealPart(yy[p * dim * 2 + 1]);
1481: }
1482: if (ctx->multispecies) {
1483: species_color = species[p] + 2;
1484: PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color));
1485: } else {
1486: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1487: }
1488: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1489: }
1490: PetscCall(VecRestoreArrayRead(u, &yy));
1491: PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
1492: PetscCall(PetscDrawSPSave(ctx->sp));
1493: if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species));
1494: }
1495: PetscFunctionReturn(PETSC_SUCCESS);
1496: }
1498: /*@C
1499: TSMonitorHGSwarmSolution - Graphically displays histograms of `DMSWARM` particles
1501: Input Parameters:
1502: + ts - the `TS` context
1503: . step - current time-step
1504: . ptime - current time
1505: . u - current solution
1506: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorHGCtxCreate()`
1508: Options Database Keys:
1509: + -ts_monitor_hg_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1510: . -ts_monitor_hg_swarm_species <num> - Number of species to histogram
1511: . -ts_monitor_hg_swarm_bins <num> - Number of histogram bins
1512: - -ts_monitor_hg_swarm_velocity <bool> - Plot in velocity space, as opposed to coordinate space
1514: Level: intermediate
1516: Note:
1517: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1518: to be used during the `TS` integration.
1520: .seealso: `TSMonitorSet()`
1521: @*/
1522: PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx dctx)
1523: {
1524: TSMonitorHGCtx ctx = (TSMonitorHGCtx)dctx;
1525: PetscDraw draw;
1526: DM sw;
1527: const PetscScalar *yy;
1528: PetscInt *species;
1529: PetscInt dim, d = 0, Np, p, Ns, s;
1531: PetscFunctionBegin;
1532: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1533: PetscCall(TSGetDM(ts, &sw));
1534: PetscCall(DMGetDimension(sw, &dim));
1535: PetscCall(DMSwarmGetNumSpecies(sw, &Ns));
1536: Ns = PetscMin(Ns, ctx->Ns);
1537: PetscCall(VecGetLocalSize(u, &Np));
1538: Np /= dim * 2;
1539: if (!step) {
1540: PetscDrawAxis axis;
1541: char title[PETSC_MAX_PATH_LEN];
1543: for (s = 0; s < Ns; ++s) {
1544: PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis));
1545: PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s));
1546: if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N"));
1547: else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N"));
1548: }
1549: }
1550: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1551: PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species));
1552: for (s = 0; s < Ns; ++s) {
1553: PetscCall(PetscDrawHGReset(ctx->hg[s]));
1554: PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw));
1555: PetscCall(PetscDrawClear(draw));
1556: PetscCall(PetscDrawFlush(draw));
1557: }
1558: PetscCall(VecGetArrayRead(u, &yy));
1559: for (p = 0; p < Np; ++p) {
1560: const PetscInt s = species[p] < Ns ? species[p] : 0;
1561: PetscReal v;
1563: if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]);
1564: else v = PetscRealPart(yy[p * dim * 2 + d]);
1565: PetscCall(PetscDrawHGAddValue(ctx->hg[s], v));
1566: }
1567: PetscCall(VecRestoreArrayRead(u, &yy));
1568: for (s = 0; s < Ns; ++s) {
1569: PetscCall(PetscDrawHGDraw(ctx->hg[s]));
1570: PetscCall(PetscDrawHGSave(ctx->hg[s]));
1571: }
1572: PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species));
1573: }
1574: PetscFunctionReturn(PETSC_SUCCESS);
1575: }
1577: /*@C
1578: TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep
1580: Collective
1582: Input Parameters:
1583: + ts - the `TS` context
1584: . step - current time-step
1585: . ptime - current time
1586: . u - current solution
1587: - vf - unused context
1589: Options Database Key:
1590: . -ts_monitor_error - create a graphical monitor of error history
1592: Level: intermediate
1594: Notes:
1595: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1596: to be used during the `TS` integration.
1598: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1600: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1601: @*/
1602: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1603: {
1604: DM dm;
1605: PetscDS ds = NULL;
1606: PetscInt Nf = -1, f;
1607: PetscBool flg;
1609: PetscFunctionBegin;
1610: PetscCall(TSGetDM(ts, &dm));
1611: if (dm) PetscCall(DMGetDS(dm, &ds));
1612: if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
1613: if (Nf <= 0) {
1614: Vec y;
1615: PetscReal nrm;
1617: PetscCall(VecDuplicate(u, &y));
1618: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1619: PetscCall(VecAXPY(y, -1.0, u));
1620: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
1621: if (flg) {
1622: PetscCall(VecNorm(y, NORM_2, &nrm));
1623: PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
1624: }
1625: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
1626: if (flg) PetscCall(VecView(y, vf->viewer));
1627: PetscCall(VecDestroy(&y));
1628: } else {
1629: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, PetscCtx ctx);
1630: void **ctxs;
1631: Vec v;
1632: PetscReal ferrors[1];
1634: PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
1635: for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
1636: PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
1637: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04" PetscInt_FMT " time = %-8.4g \t L_2 Error: [", step, (double)ptime));
1638: for (f = 0; f < Nf; ++f) {
1639: if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
1640: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
1641: }
1642: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));
1644: PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));
1646: PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
1647: if (flg) {
1648: PetscCall(DMGetGlobalVector(dm, &v));
1649: PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
1650: PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
1651: PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
1652: PetscCall(DMRestoreGlobalVector(dm, &v));
1653: }
1654: PetscCall(PetscFree2(exactFuncs, ctxs));
1655: }
1656: PetscFunctionReturn(PETSC_SUCCESS);
1657: }
1659: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, PetscCtx monctx)
1660: {
1661: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1662: PetscReal x = ptime, y;
1663: PetscInt its;
1665: PetscFunctionBegin;
1666: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1667: if (!n) {
1668: PetscDrawAxis axis;
1669: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1670: PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
1671: PetscCall(PetscDrawLGReset(ctx->lg));
1672: ctx->snes_its = 0;
1673: }
1674: PetscCall(TSGetSNESIterations(ts, &its));
1675: y = its - ctx->snes_its;
1676: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1677: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1678: PetscCall(PetscDrawLGDraw(ctx->lg));
1679: PetscCall(PetscDrawLGSave(ctx->lg));
1680: }
1681: ctx->snes_its = its;
1682: PetscFunctionReturn(PETSC_SUCCESS);
1683: }
1685: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, PetscCtx monctx)
1686: {
1687: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1688: PetscReal x = ptime, y;
1689: PetscInt its;
1691: PetscFunctionBegin;
1692: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1693: if (!n) {
1694: PetscDrawAxis axis;
1695: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1696: PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
1697: PetscCall(PetscDrawLGReset(ctx->lg));
1698: ctx->ksp_its = 0;
1699: }
1700: PetscCall(TSGetKSPIterations(ts, &its));
1701: y = its - ctx->ksp_its;
1702: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1703: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1704: PetscCall(PetscDrawLGDraw(ctx->lg));
1705: PetscCall(PetscDrawLGSave(ctx->lg));
1706: }
1707: ctx->ksp_its = its;
1708: PetscFunctionReturn(PETSC_SUCCESS);
1709: }
1711: /*@C
1712: TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()`
1714: Collective
1716: Input Parameter:
1717: . ts - the `TS` solver object
1719: Output Parameter:
1720: . ctx - the context
1722: Level: intermediate
1724: .seealso: [](ch_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1725: @*/
1726: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1727: {
1728: PetscFunctionBegin;
1729: PetscCall(PetscNew(ctx));
1730: PetscFunctionReturn(PETSC_SUCCESS);
1731: }
1733: /*@C
1734: TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution
1736: Collective
1738: Input Parameters:
1739: + ts - the `TS` context
1740: . step - current time-step
1741: . ptime - current time
1742: . u - current solution
1743: - dctx - the envelope context
1745: Options Database Key:
1746: . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time
1748: Level: intermediate
1750: Notes:
1751: After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope
1753: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1754: to be used during the `TS` integration.
1756: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1757: @*/
1758: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscCtx dctx)
1759: {
1760: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;
1762: PetscFunctionBegin;
1763: if (!ctx->max) {
1764: PetscCall(VecDuplicate(u, &ctx->max));
1765: PetscCall(VecDuplicate(u, &ctx->min));
1766: PetscCall(VecCopy(u, ctx->max));
1767: PetscCall(VecCopy(u, ctx->min));
1768: } else {
1769: PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
1770: PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
1771: }
1772: PetscFunctionReturn(PETSC_SUCCESS);
1773: }
1775: /*@C
1776: TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution
1778: Collective
1780: Input Parameter:
1781: . ts - the `TS` context
1783: Output Parameters:
1784: + max - the maximum values
1785: - min - the minimum values
1787: Level: intermediate
1789: Notes:
1790: If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored
1792: .seealso: [](ch_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1793: @*/
1794: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1795: {
1796: PetscInt i;
1798: PetscFunctionBegin;
1799: if (max) *max = NULL;
1800: if (min) *min = NULL;
1801: for (i = 0; i < ts->numbermonitors; i++) {
1802: if (ts->monitor[i] == TSMonitorEnvelope) {
1803: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1804: if (max) *max = ctx->max;
1805: if (min) *min = ctx->min;
1806: break;
1807: }
1808: }
1809: PetscFunctionReturn(PETSC_SUCCESS);
1810: }
1812: /*@C
1813: TSMonitorEnvelopeCtxDestroy - Destroys a context that was created with `TSMonitorEnvelopeCtxCreate()`.
1815: Collective
1817: Input Parameter:
1818: . ctx - the monitor context
1820: Level: intermediate
1822: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1823: @*/
1824: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1825: {
1826: PetscFunctionBegin;
1827: PetscCall(VecDestroy(&(*ctx)->min));
1828: PetscCall(VecDestroy(&(*ctx)->max));
1829: PetscCall(PetscFree(*ctx));
1830: PetscFunctionReturn(PETSC_SUCCESS);
1831: }
1833: /*@C
1834: TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS`
1836: Not Collective
1838: Input Parameters:
1839: + ts - the `TS` context
1840: . step - current timestep
1841: . t - current time
1842: . U - current solution
1843: - vf - not used
1845: Options Database Key:
1846: + -ts_dmswarm_monitor_moments - Monitor moments of particle distribution
1847: - -ts_dmswarm_monitor_moments_interval - Interval of timesteps between monitor outputs
1849: Level: intermediate
1851: Notes:
1852: This requires a `DMSWARM` be attached to the `TS`.
1854: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1855: to be used during the TS integration.
1857: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1858: @*/
1859: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1860: {
1861: DM sw;
1862: const PetscScalar *u;
1863: PetscReal m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1864: PetscInt dim, d, Np, p;
1865: MPI_Comm comm;
1867: PetscFunctionBeginUser;
1868: (void)t;
1869: (void)vf;
1870: PetscCall(TSGetDM(ts, &sw));
1871: if (!sw || step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS);
1872: PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
1873: PetscCall(DMGetDimension(sw, &dim));
1874: PetscCall(VecGetLocalSize(U, &Np));
1875: Np /= dim;
1876: PetscCall(VecGetArrayRead(U, &u));
1877: for (p = 0; p < Np; ++p) {
1878: for (d = 0; d < dim; ++d) {
1879: totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1880: totMom[d] += PetscRealPart(u[p * dim + d]);
1881: }
1882: }
1883: PetscCall(VecRestoreArrayRead(U, &u));
1884: for (d = 0; d < dim; ++d) totMom[d] *= m;
1885: totE *= 0.5 * m;
1886: PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
1887: for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, " Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
1888: PetscCall(PetscPrintf(comm, "\n"));
1889: PetscFunctionReturn(PETSC_SUCCESS);
1890: }