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));
86: PetscCall(PetscViewerDestroy(&viewer));
87: if (monitorsetup) PetscCall((*monitorsetup)(ts, vf));
88: PetscCall(TSMonitorSet(ts, (PetscErrorCode (*)(TS, PetscInt, PetscReal, Vec, void *))monitor, vf, (PetscCtxDestroyFn *)PetscViewerAndFormatDestroy));
89: }
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: /*@C
94: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
95: timestep to display the iteration's progress.
97: Logically Collective
99: Input Parameters:
100: + ts - the `TS` context obtained from `TSCreate()`
101: . monitor - monitoring routine
102: . mctx - [optional] user-defined context for private data for the monitor routine (use `NULL` if no context is desired)
103: - mdestroy - [optional] routine that frees monitor context (may be `NULL`), see `PetscCtxDestroyFn` for the calling sequence
105: Calling sequence of `monitor`:
106: + ts - the `TS` context
107: . 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)
108: . time - current time
109: . u - current iterate
110: - ctx - [optional] monitoring context
112: Level: intermediate
114: Note:
115: This routine adds an additional monitor to the list of monitors that already has been loaded.
117: Fortran Notes:
118: Only a single monitor function can be set for each `TS` object
120: .seealso: [](ch_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
121: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
122: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `PetscCtxDestroyFn`
123: @*/
124: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS ts, PetscInt steps, PetscReal time, Vec u, void *ctx), void *mctx, PetscCtxDestroyFn *mdestroy)
125: {
126: PetscInt i;
127: PetscBool identical;
129: PetscFunctionBegin;
131: for (i = 0; i < ts->numbermonitors; i++) {
132: PetscCall(PetscMonitorCompare((PetscErrorCode (*)(void))monitor, mctx, mdestroy, (PetscErrorCode (*)(void))ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical));
133: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
134: }
135: PetscCheck(ts->numbermonitors < MAXTSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set");
136: ts->monitor[ts->numbermonitors] = monitor;
137: ts->monitordestroy[ts->numbermonitors] = mdestroy;
138: ts->monitorcontext[ts->numbermonitors++] = mctx;
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: /*@C
143: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
145: Logically Collective
147: Input Parameter:
148: . ts - the `TS` context obtained from `TSCreate()`
150: Level: intermediate
152: Note:
153: There is no way to remove a single, specific monitor.
155: .seealso: [](ch_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
156: @*/
157: PetscErrorCode TSMonitorCancel(TS ts)
158: {
159: PetscInt i;
161: PetscFunctionBegin;
163: for (i = 0; i < ts->numbermonitors; i++) {
164: if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i]));
165: }
166: ts->numbermonitors = 0;
167: PetscFunctionReturn(PETSC_SUCCESS);
168: }
170: /*@C
171: TSMonitorDefault - The default monitor, prints the timestep and time for each step
173: Input Parameters:
174: + ts - the `TS` context
175: . 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)
176: . ptime - current time
177: . v - current iterate
178: - vf - the viewer and format
180: Options Database Key:
181: . -ts_monitor - monitors the time integration
183: Level: intermediate
185: Notes:
186: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
187: to be used during the `TS` integration.
189: .seealso: [](ch_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorWallClockTime()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
190: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
191: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`
192: @*/
193: PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
194: {
195: PetscViewer viewer = vf->viewer;
196: PetscBool isascii, ibinary;
198: PetscFunctionBegin;
200: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
201: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary));
202: PetscCall(PetscViewerPushFormat(viewer, vf->format));
203: if (isascii) {
204: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
205: if (step == -1) { /* this indicates it is an interpolated solution */
206: PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps));
207: } else {
208: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n"));
209: }
210: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
211: } else if (ibinary) {
212: PetscMPIInt rank;
213: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
214: if (rank == 0) {
215: PetscBool skipHeader;
216: PetscInt classid = REAL_FILE_CLASSID;
218: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
219: if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
220: PetscCall(PetscRealView(1, &ptime, viewer));
221: } else {
222: PetscCall(PetscRealView(0, &ptime, viewer));
223: }
224: }
225: PetscCall(PetscViewerPopFormat(viewer));
226: PetscFunctionReturn(PETSC_SUCCESS);
227: }
229: typedef struct {
230: PetscLogDouble time_start;
231: PetscLogDouble time_last;
232: PetscInt snes_its;
233: PetscInt ksp_its;
234: } *TSMonitorWallClockTimeContext;
236: /*@C
237: TSMonitorWallClockTimeSetUp - Setup routine passed to `TSMonitorSetFromOptions()` when using `-ts_monitor_wall_clock_time`
239: Input Parameters:
240: + ts - the `TS` context
241: - vf - the viewer and format
243: Level: intermediate
245: Note:
246: This is not called directly by users, rather one calls `TSMonitorSetFromOptions()`, with `TSMonitorWallClockTime()` and this function as arguments, to cause the monitor
247: to be used during the `TS` integration.
249: .seealso: [](ch_ts), `TSMonitorSet()`
250: @*/
251: PetscErrorCode TSMonitorWallClockTimeSetUp(TS ts, PetscViewerAndFormat *vf)
252: {
253: TSMonitorWallClockTimeContext speed;
255: PetscFunctionBegin;
256: PetscCall(PetscNew(&speed));
257: speed->time_start = PETSC_DECIDE;
258: vf->data_destroy = PetscCtxDestroyDefault;
259: vf->data = speed;
260: PetscFunctionReturn(PETSC_SUCCESS);
261: }
263: /*@C
264: TSMonitorWallClockTime - Monitor wall-clock time, KSP iterations, and SNES iterations per step.
266: Input Parameters:
267: + ts - the `TS` context
268: . 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)
269: . ptime - current time
270: . v - current solution
271: - vf - the viewer and format
273: Options Database Key:
274: . -ts_monitor_wall_clock_time - Monitor wall-clock time, KSP iterations, and SNES iterations per step.
276: Level: intermediate
278: Note:
279: This is not called directly by users, rather one calls `TSMonitorSetFromOptions()`, with this function and `TSMonitorWallClockTimeSetUp()` as arguments, to cause the monitor
280: to be used during the `TS` integration.
282: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
283: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
284: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
285: @*/
286: PetscErrorCode TSMonitorWallClockTime(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
287: {
288: PetscViewer viewer = vf->viewer;
289: TSMonitorWallClockTimeContext speed = (TSMonitorWallClockTimeContext)vf->data;
290: PetscBool isascii;
291: PetscLogDouble now;
292: PetscInt snes_its, ksp_its;
294: PetscFunctionBegin;
296: PetscCall(PetscTime(&now));
297: if (speed->time_start == PETSC_DECIDE) {
298: speed->time_start = now;
299: speed->time_last = now;
300: }
301: PetscCall(TSGetSNESIterations(ts, &snes_its));
302: PetscCall(TSGetKSPIterations(ts, &ksp_its));
303: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
304: PetscCall(PetscViewerPushFormat(viewer, vf->format));
305: if (isascii) {
306: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
307: 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,
308: now - speed->time_start, snes_its - speed->snes_its, ksp_its - speed->ksp_its));
309: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
310: }
311: PetscCall(PetscViewerPopFormat(viewer));
312: speed->time_last = now;
313: speed->snes_its = snes_its;
314: speed->ksp_its = ksp_its;
315: PetscFunctionReturn(PETSC_SUCCESS);
316: }
318: /*@C
319: TSMonitorExtreme - Prints the extreme values of the solution at each timestep
321: Input Parameters:
322: + ts - the `TS` context
323: . 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)
324: . ptime - current time
325: . v - current iterate
326: - vf - the viewer and format
328: Level: intermediate
330: Note:
331: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
332: to be used during the `TS` integration.
334: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`
335: @*/
336: PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
337: {
338: PetscViewer viewer = vf->viewer;
339: PetscBool isascii;
340: PetscReal max, min;
342: PetscFunctionBegin;
344: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
345: PetscCall(PetscViewerPushFormat(viewer, vf->format));
346: if (isascii) {
347: PetscCall(VecMax(v, NULL, &max));
348: PetscCall(VecMin(v, NULL, &min));
349: PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
350: 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));
351: PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
352: }
353: PetscCall(PetscViewerPopFormat(viewer));
354: PetscFunctionReturn(PETSC_SUCCESS);
355: }
357: /*@C
358: TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with
359: `TS` to monitor the solution process graphically in various ways
361: Collective
363: Input Parameters:
364: + comm - the MPI communicator to use
365: . host - the X display to open, or `NULL` for the local machine
366: . label - the title to put in the title bar
367: . x - the x screen coordinates of the upper left coordinate of the window
368: . y - the y screen coordinates of the upper left coordinate of the window
369: . m - the screen width in pixels
370: . n - the screen height in pixels
371: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
373: Output Parameter:
374: . ctx - the context
376: Options Database Keys:
377: + -ts_monitor_lg_timestep - automatically sets line graph monitor
378: . -ts_monitor_lg_timestep_log - automatically sets line graph monitor
379: . -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`)
380: . -ts_monitor_lg_error - monitor the error
381: . -ts_monitor_lg_ksp_iterations - monitor the number of `KSP` iterations needed for each timestep
382: . -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep
383: - -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true
385: Level: intermediate
387: Notes:
388: Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed.
390: One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()`
392: Many of the functions that control the monitoring have two forms\: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a `TS` object as the
393: 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
394: as the first argument.
396: One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()`
398: .seealso: [](ch_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
399: `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
400: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
401: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
402: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
403: @*/
404: PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
405: {
406: PetscDraw draw;
408: PetscFunctionBegin;
409: PetscCall(PetscNew(ctx));
410: PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
411: PetscCall(PetscDrawSetFromOptions(draw));
412: PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg));
413: PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg));
414: PetscCall(PetscDrawDestroy(&draw));
415: (*ctx)->howoften = howoften;
416: PetscFunctionReturn(PETSC_SUCCESS);
417: }
419: /*@C
420: TSMonitorLGTimeStep - Monitors a `TS` by printing the time-steps
422: Collective
424: Input Parameters:
425: + ts - the time integrator
426: . step - the current time step
427: . ptime - the current time
428: . v - the current state
429: - monctx - the monitor context obtained with `TSMonitorLGCtxCreate()`
431: Level: advanced
433: Note:
434: This is not called directly by users, rather one calls `TSMonitorSet()` along the `ctx` created by `TSMonitorLGCtxCreate()`
435: and `TSMonitorLGCtxDestroy()`
437: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGCtxDestroy()`
438: @*/
439: PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx)
440: {
441: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
442: PetscReal x = ptime, y;
444: PetscFunctionBegin;
445: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */
446: if (!step) {
447: PetscDrawAxis axis;
448: const char *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
449: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
450: PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel));
451: PetscCall(PetscDrawLGReset(ctx->lg));
452: }
453: PetscCall(TSGetTimeStep(ts, &y));
454: if (ctx->semilogy) y = PetscLog10Real(y);
455: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
456: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
457: PetscCall(PetscDrawLGDraw(ctx->lg));
458: PetscCall(PetscDrawLGSave(ctx->lg));
459: }
460: PetscFunctionReturn(PETSC_SUCCESS);
461: }
463: /*@C
464: TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`.
466: Collective
468: Input Parameter:
469: . ctx - the monitor context
471: Level: intermediate
473: Note:
474: Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`
476: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
477: @*/
478: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
479: {
480: PetscFunctionBegin;
481: if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)((void **)&(*ctx)->transformctx));
482: PetscCall(PetscDrawLGDestroy(&(*ctx)->lg));
483: PetscCall(PetscStrArrayDestroy(&(*ctx)->names));
484: PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames));
485: PetscCall(PetscFree((*ctx)->displayvariables));
486: PetscCall(PetscFree((*ctx)->displayvalues));
487: PetscCall(PetscFree(*ctx));
488: PetscFunctionReturn(PETSC_SUCCESS);
489: }
491: /* Creates a TSMonitorSPCtx for use with DMSwarm particle visualizations */
492: 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)
493: {
494: PetscDraw draw;
496: PetscFunctionBegin;
497: PetscCall(PetscNew(ctx));
498: PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
499: PetscCall(PetscDrawSetFromOptions(draw));
500: PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp));
501: PetscCall(PetscDrawDestroy(&draw));
502: (*ctx)->howoften = howoften;
503: (*ctx)->retain = retain;
504: (*ctx)->phase = phase;
505: (*ctx)->multispecies = multispecies;
506: PetscFunctionReturn(PETSC_SUCCESS);
507: }
509: /* Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate */
510: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
511: {
512: PetscFunctionBegin;
513: PetscCall(PetscDrawSPDestroy(&(*ctx)->sp));
514: PetscCall(PetscFree(*ctx));
515: PetscFunctionReturn(PETSC_SUCCESS);
516: }
518: /* Creates a TSMonitorHGCtx for use with DMSwarm particle visualizations */
519: 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)
520: {
521: PetscDraw draw;
522: int Nsi, Nbi;
524: PetscFunctionBegin;
525: PetscCall(PetscMPIIntCast(Ns, &Nsi));
526: PetscCall(PetscMPIIntCast(Nb, &Nbi));
527: PetscCall(PetscNew(ctx));
528: PetscCall(PetscMalloc1(Ns, &(*ctx)->hg));
529: for (int s = 0; s < Nsi; ++s) {
530: PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw));
531: PetscCall(PetscDrawSetFromOptions(draw));
532: PetscCall(PetscDrawHGCreate(draw, Nbi, &(*ctx)->hg[s]));
533: PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE));
534: PetscCall(PetscDrawDestroy(&draw));
535: }
536: (*ctx)->howoften = howoften;
537: (*ctx)->Ns = Ns;
538: (*ctx)->velocity = velocity;
539: PetscFunctionReturn(PETSC_SUCCESS);
540: }
542: /* Destroys a TSMonitorHGCtx that was created with TSMonitorHGCtxCreate */
543: PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx)
544: {
545: PetscInt s;
547: PetscFunctionBegin;
548: for (s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s]));
549: PetscCall(PetscFree((*ctx)->hg));
550: PetscCall(PetscFree(*ctx));
551: PetscFunctionReturn(PETSC_SUCCESS);
552: }
554: /*@C
555: TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling
556: `VecView()` for the solution at each timestep
558: Collective
560: Input Parameters:
561: + ts - the `TS` context
562: . step - current time-step
563: . ptime - current time
564: . u - the solution at the current time
565: - dummy - either a viewer or `NULL`
567: Options Database Keys:
568: + -ts_monitor_draw_solution - draw the solution at each time-step
569: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
571: Level: intermediate
573: Notes:
574: The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial`
575: will look bad
577: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with
578: `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration.
580: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
581: @*/
582: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
583: {
584: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
585: PetscDraw draw;
587: PetscFunctionBegin;
588: if (!step && ictx->showinitial) {
589: if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
590: PetscCall(VecCopy(u, ictx->initialsolution));
591: }
592: if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
594: if (ictx->showinitial) {
595: PetscReal pause;
596: PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
597: PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
598: PetscCall(VecView(ictx->initialsolution, ictx->viewer));
599: PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
600: PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
601: }
602: PetscCall(VecView(u, ictx->viewer));
603: if (ictx->showtimestepandtime) {
604: PetscReal xl, yl, xr, yr, h;
605: char time[32];
607: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
608: PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
609: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
610: h = yl + .95 * (yr - yl);
611: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
612: PetscCall(PetscDrawFlush(draw));
613: }
615: if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
616: PetscFunctionReturn(PETSC_SUCCESS);
617: }
619: /*@C
620: TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram
622: Collective
624: Input Parameters:
625: + ts - the `TS` context
626: . step - current time-step
627: . ptime - current time
628: . u - the solution at the current time
629: - dummy - either a viewer or `NULL`
631: Level: intermediate
633: Notes:
634: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
635: to be used during the `TS` integration.
637: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
638: @*/
639: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
640: {
641: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
642: PetscDraw draw;
643: PetscDrawAxis axis;
644: PetscInt n;
645: PetscMPIInt size;
646: PetscReal U0, U1, xl, yl, xr, yr, h;
647: char time[32];
648: const PetscScalar *U;
650: PetscFunctionBegin;
651: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
652: PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
653: PetscCall(VecGetSize(u, &n));
654: PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");
656: PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
657: PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
658: PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
659: if (!step) {
660: PetscCall(PetscDrawClear(draw));
661: PetscCall(PetscDrawAxisDraw(axis));
662: }
664: PetscCall(VecGetArrayRead(u, &U));
665: U0 = PetscRealPart(U[0]);
666: U1 = PetscRealPart(U[1]);
667: PetscCall(VecRestoreArrayRead(u, &U));
668: if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);
670: PetscDrawCollectiveBegin(draw);
671: PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
672: if (ictx->showtimestepandtime) {
673: PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
674: PetscCall(PetscSNPrintf(time, 32, "Timestep %" PetscInt_FMT " Time %g", step, (double)ptime));
675: h = yl + .95 * (yr - yl);
676: PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
677: }
678: PetscDrawCollectiveEnd(draw);
679: PetscCall(PetscDrawFlush(draw));
680: PetscCall(PetscDrawPause(draw));
681: PetscCall(PetscDrawSave(draw));
682: PetscFunctionReturn(PETSC_SUCCESS);
683: }
685: /*@C
686: TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()`
688: Collective
690: Input Parameter:
691: . ictx - the monitor context
693: Level: intermediate
695: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
696: @*/
697: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
698: {
699: PetscFunctionBegin;
700: PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
701: PetscCall(VecDestroy(&(*ictx)->initialsolution));
702: PetscCall(PetscFree(*ictx));
703: PetscFunctionReturn(PETSC_SUCCESS);
704: }
706: /*@C
707: TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`
709: Collective
711: Input Parameters:
712: + comm - the MPI communicator to use
713: . host - the X display to open, or `NULL` for the local machine
714: . label - the title to put in the title bar
715: . x - the x screen coordinates of the upper left coordinate of the window
716: . y - the y screen coordinates of the upper left coordinate of the window
717: . m - the screen width in pixels
718: . n - the screen height in pixels
719: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
721: Output Parameter:
722: . ctx - the monitor context
724: Options Database Keys:
725: + -ts_monitor_draw_solution - draw the solution at each time-step
726: - -ts_monitor_draw_solution_initial - show initial solution as well as current solution
728: Level: intermediate
730: Note:
731: The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`.
733: .seealso: [](ch_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
734: @*/
735: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
736: {
737: PetscFunctionBegin;
738: PetscCall(PetscNew(ctx));
739: PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
740: PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));
742: (*ctx)->howoften = howoften;
743: (*ctx)->showinitial = PETSC_FALSE;
744: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));
746: (*ctx)->showtimestepandtime = PETSC_FALSE;
747: PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
748: PetscFunctionReturn(PETSC_SUCCESS);
749: }
751: /*@C
752: TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling
753: `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep
755: Collective
757: Input Parameters:
758: + ts - the `TS` context
759: . step - current time-step
760: . ptime - current time
761: . u - solution at current time
762: - dummy - either a viewer or `NULL`
764: Options Database Key:
765: . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
767: Level: intermediate
769: Note:
770: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
771: to be used during the `TS` integration.
773: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
774: @*/
775: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
776: {
777: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
778: PetscViewer viewer = ctx->viewer;
779: Vec work;
781: PetscFunctionBegin;
782: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
783: PetscCall(VecDuplicate(u, &work));
784: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
785: PetscCall(VecView(work, viewer));
786: PetscCall(VecDestroy(&work));
787: PetscFunctionReturn(PETSC_SUCCESS);
788: }
790: /*@C
791: TSMonitorDrawError - Monitors progress of the `TS` solvers by calling
792: `VecView()` for the error at each timestep
794: Collective
796: Input Parameters:
797: + ts - the `TS` context
798: . step - current time-step
799: . ptime - current time
800: . u - solution at current time
801: - dummy - either a viewer or `NULL`
803: Options Database Key:
804: . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`
806: Level: intermediate
808: Notes:
809: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
810: to be used during the `TS` integration.
812: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
813: @*/
814: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
815: {
816: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
817: PetscViewer viewer = ctx->viewer;
818: Vec work;
820: PetscFunctionBegin;
821: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
822: PetscCall(VecDuplicate(u, &work));
823: PetscCall(TSComputeSolutionFunction(ts, ptime, work));
824: PetscCall(VecAXPY(work, -1.0, u));
825: PetscCall(VecView(work, viewer));
826: PetscCall(VecDestroy(&work));
827: PetscFunctionReturn(PETSC_SUCCESS);
828: }
830: /*@C
831: TSMonitorSolutionSetup - Setups the context for `TSMonitorSolution()`
833: Collective
835: Input Parameters:
836: + ts - the `TS` context
837: - vf - viewer and its format
839: Level: intermediate
841: .seealso: [](ch_ts), `TS`, `TSMonitorSolution()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSetFromOptions()`
842: @*/
843: PetscErrorCode TSMonitorSolutionSetup(TS ts, PetscViewerAndFormat *vf)
844: {
845: TSMonitorSolutionCtx ctx;
847: PetscFunctionBegin;
848: PetscCall(PetscNew(&ctx));
849: PetscCall(PetscOptionsGetBool(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_solution_skip_initial", &ctx->skip_initial, NULL));
850: vf->data = ctx;
851: vf->data_destroy = PetscCtxDestroyDefault;
852: PetscFunctionReturn(PETSC_SUCCESS);
853: }
855: /*@C
856: 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
858: Collective
860: Input Parameters:
861: + ts - the `TS` context
862: . step - current time-step
863: . ptime - current time
864: . u - current state
865: - vf - viewer and its format
867: Level: intermediate
869: Notes:
870: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
871: to be used during the `TS` integration.
873: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorSolutionSetup()`,
874: @*/
875: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
876: {
877: TSMonitorSolutionCtx ctx = (TSMonitorSolutionCtx)vf->data;
879: PetscFunctionBegin;
880: if (ctx->skip_initial && step == ts->start_step) PetscFunctionReturn(PETSC_SUCCESS);
881: if ((vf->view_interval > 0 && !(step % vf->view_interval)) || (vf->view_interval && ts->reason)) {
882: PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
883: PetscCall(VecView(u, vf->viewer));
884: PetscCall(PetscViewerPopFormat(vf->viewer));
885: }
886: PetscFunctionReturn(PETSC_SUCCESS);
887: }
889: /*@C
890: TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at selected timesteps.
892: Collective
894: Input Parameters:
895: + ts - the `TS` context
896: . step - current time-step
897: . ptime - current time
898: . u - current state
899: - ctx - monitor context obtained with `TSMonitorSolutionVTKCtxCreate()`
901: Level: developer
903: Notes:
904: 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.
905: These are named according to the file name template.
907: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
908: to be used during the `TS` integration.
910: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
911: @*/
912: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, TSMonitorVTKCtx ctx)
913: {
914: char filename[PETSC_MAX_PATH_LEN];
915: PetscViewer viewer;
917: PetscFunctionBegin;
918: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
919: if (((ctx->interval > 0) && (!(step % ctx->interval))) || (ctx->interval && ts->reason)) {
920: PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)ctx->filenametemplate, step));
921: PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
922: PetscCall(VecView(u, viewer));
923: PetscCall(PetscViewerDestroy(&viewer));
924: }
925: PetscFunctionReturn(PETSC_SUCCESS);
926: }
928: /*@C
929: TSMonitorSolutionVTKDestroy - Destroy the monitor context created with `TSMonitorSolutionVTKCtxCreate()`
931: Not Collective
933: Input Parameter:
934: . ctx - the monitor context
936: Level: developer
938: Note:
939: This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
941: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
942: @*/
943: PetscErrorCode TSMonitorSolutionVTKDestroy(TSMonitorVTKCtx *ctx)
944: {
945: PetscFunctionBegin;
946: PetscCall(PetscFree((*ctx)->filenametemplate));
947: PetscCall(PetscFree(*ctx));
948: PetscFunctionReturn(PETSC_SUCCESS);
949: }
951: /*@C
952: TSMonitorSolutionVTKCtxCreate - Create the monitor context to be used in `TSMonitorSolutionVTK()`
954: Not collective
956: Input Parameter:
957: . filenametemplate - the template file name, e.g. foo-%03d.vts
959: Output Parameter:
960: . ctx - the monitor context
962: Level: developer
964: Note:
965: This function is normally used inside `TSSetFromOptions()` to pass the context created to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.
967: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`, `TSMonitorSolutionVTKDestroy()`
968: @*/
969: PetscErrorCode TSMonitorSolutionVTKCtxCreate(const char *filenametemplate, TSMonitorVTKCtx *ctx)
970: {
971: const char *ptr = NULL, *ptr2 = NULL;
972: TSMonitorVTKCtx ictx;
974: PetscFunctionBegin;
975: PetscAssertPointer(filenametemplate, 1);
976: PetscAssertPointer(ctx, 2);
977: /* Do some cursory validation of the input. */
978: PetscCall(PetscStrstr(filenametemplate, "%", (char **)&ptr));
979: PetscCheck(ptr, PETSC_COMM_SELF, PETSC_ERR_USER, "-ts_monitor_solution_vtk requires a file template, e.g. filename-%%03" PetscInt_FMT ".vts");
980: for (ptr++; ptr && *ptr; ptr++) {
981: PetscCall(PetscStrchr("DdiouxX", *ptr, (char **)&ptr2));
982: 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");
983: if (ptr2) break;
984: }
985: PetscCall(PetscNew(&ictx));
986: PetscCall(PetscStrallocpy(filenametemplate, &ictx->filenametemplate));
987: ictx->interval = 1;
989: *ctx = ictx;
990: PetscFunctionReturn(PETSC_SUCCESS);
991: }
993: /*@C
994: TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector
995: in a time based line graph
997: Collective
999: Input Parameters:
1000: + ts - the `TS` context
1001: . step - current time-step
1002: . ptime - current time
1003: . u - current solution
1004: - dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`
1006: Options Database Key:
1007: . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables
1009: Level: intermediate
1011: Notes:
1012: Each process in a parallel run displays its component solutions in a separate window
1014: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1015: to be used during the `TS` integration.
1017: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
1018: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
1019: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
1020: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
1021: @*/
1022: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1023: {
1024: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dctx;
1025: const PetscScalar *yy;
1026: Vec v;
1028: PetscFunctionBegin;
1029: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1030: if (!step) {
1031: PetscDrawAxis axis;
1032: PetscInt dim;
1033: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1034: PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
1035: if (!ctx->names) {
1036: PetscBool flg;
1037: /* user provides names of variables to plot but no names has been set so assume names are integer values */
1038: PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
1039: if (flg) {
1040: PetscInt i, n;
1041: char **names;
1042: PetscCall(VecGetSize(u, &n));
1043: PetscCall(PetscMalloc1(n + 1, &names));
1044: for (i = 0; i < n; i++) {
1045: PetscCall(PetscMalloc1(5, &names[i]));
1046: PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
1047: }
1048: names[n] = NULL;
1049: ctx->names = names;
1050: }
1051: }
1052: if (ctx->names && !ctx->displaynames) {
1053: char **displaynames;
1054: PetscBool flg;
1055: PetscCall(VecGetLocalSize(u, &dim));
1056: PetscCall(PetscCalloc1(dim + 1, &displaynames));
1057: PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
1058: if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
1059: PetscCall(PetscStrArrayDestroy(&displaynames));
1060: }
1061: if (ctx->displaynames) {
1062: PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
1063: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
1064: } else if (ctx->names) {
1065: PetscCall(VecGetLocalSize(u, &dim));
1066: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1067: PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
1068: } else {
1069: PetscCall(VecGetLocalSize(u, &dim));
1070: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1071: }
1072: PetscCall(PetscDrawLGReset(ctx->lg));
1073: }
1075: if (!ctx->transform) v = u;
1076: else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
1077: PetscCall(VecGetArrayRead(v, &yy));
1078: if (ctx->displaynames) {
1079: PetscInt i;
1080: for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
1081: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
1082: } else {
1083: #if defined(PETSC_USE_COMPLEX)
1084: PetscInt i, n;
1085: PetscReal *yreal;
1086: PetscCall(VecGetLocalSize(v, &n));
1087: PetscCall(PetscMalloc1(n, &yreal));
1088: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1089: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1090: PetscCall(PetscFree(yreal));
1091: #else
1092: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1093: #endif
1094: }
1095: PetscCall(VecRestoreArrayRead(v, &yy));
1096: if (ctx->transform) PetscCall(VecDestroy(&v));
1098: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1099: PetscCall(PetscDrawLGDraw(ctx->lg));
1100: PetscCall(PetscDrawLGSave(ctx->lg));
1101: }
1102: PetscFunctionReturn(PETSC_SUCCESS);
1103: }
1105: /*@C
1106: TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
1108: Collective
1110: Input Parameters:
1111: + ts - the `TS` context
1112: - names - the names of the components, final string must be `NULL`
1114: Level: intermediate
1116: Notes:
1117: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1119: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
1120: @*/
1121: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
1122: {
1123: PetscInt i;
1125: PetscFunctionBegin;
1126: for (i = 0; i < ts->numbermonitors; i++) {
1127: if (ts->monitor[i] == TSMonitorLGSolution) {
1128: PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
1129: break;
1130: }
1131: }
1132: PetscFunctionReturn(PETSC_SUCCESS);
1133: }
1135: /*@C
1136: TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
1138: Collective
1140: Input Parameters:
1141: + ctx - the `TS` context
1142: - names - the names of the components, final string must be `NULL`
1144: Level: intermediate
1146: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
1147: @*/
1148: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
1149: {
1150: PetscFunctionBegin;
1151: PetscCall(PetscStrArrayDestroy(&ctx->names));
1152: PetscCall(PetscStrArrayallocpy(names, &ctx->names));
1153: PetscFunctionReturn(PETSC_SUCCESS);
1154: }
1156: /*@C
1157: TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot
1159: Collective
1161: Input Parameter:
1162: . ts - the `TS` context
1164: Output Parameter:
1165: . names - the names of the components, final string must be `NULL`
1167: Level: intermediate
1169: Note:
1170: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1172: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1173: @*/
1174: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
1175: {
1176: PetscInt i;
1178: PetscFunctionBegin;
1179: *names = NULL;
1180: for (i = 0; i < ts->numbermonitors; i++) {
1181: if (ts->monitor[i] == TSMonitorLGSolution) {
1182: TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
1183: *names = (const char *const *)ctx->names;
1184: break;
1185: }
1186: }
1187: PetscFunctionReturn(PETSC_SUCCESS);
1188: }
1190: /*@C
1191: TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor
1193: Collective
1195: Input Parameters:
1196: + ctx - the `TSMonitorLG` context
1197: - displaynames - the names of the components, final string must be `NULL`
1199: Level: intermediate
1201: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1202: @*/
1203: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
1204: {
1205: PetscInt j = 0, k;
1207: PetscFunctionBegin;
1208: if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
1209: PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
1210: PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
1211: while (displaynames[j]) j++;
1212: ctx->ndisplayvariables = j;
1213: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
1214: PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
1215: j = 0;
1216: while (displaynames[j]) {
1217: k = 0;
1218: while (ctx->names[k]) {
1219: PetscBool flg;
1220: PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
1221: if (flg) {
1222: ctx->displayvariables[j] = k;
1223: break;
1224: }
1225: k++;
1226: }
1227: j++;
1228: }
1229: PetscFunctionReturn(PETSC_SUCCESS);
1230: }
1232: /*@C
1233: TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor
1235: Collective
1237: Input Parameters:
1238: + ts - the `TS` context
1239: - displaynames - the names of the components, final string must be `NULL`
1241: Level: intermediate
1243: Note:
1244: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1246: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1247: @*/
1248: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
1249: {
1250: PetscInt i;
1252: PetscFunctionBegin;
1253: for (i = 0; i < ts->numbermonitors; i++) {
1254: if (ts->monitor[i] == TSMonitorLGSolution) {
1255: PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
1256: break;
1257: }
1258: }
1259: PetscFunctionReturn(PETSC_SUCCESS);
1260: }
1262: /*@C
1263: TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed
1265: Collective
1267: Input Parameters:
1268: + ts - the `TS` context
1269: . transform - the transform function
1270: . destroy - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1271: - tctx - optional context used by transform function
1273: Level: intermediate
1275: Note:
1276: If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored
1278: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()`, `PetscCtxDestroyFn`
1279: @*/
1280: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscCtxDestroyFn *destroy, void *tctx)
1281: {
1282: PetscInt i;
1284: PetscFunctionBegin;
1285: for (i = 0; i < ts->numbermonitors; i++) {
1286: if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
1287: }
1288: PetscFunctionReturn(PETSC_SUCCESS);
1289: }
1291: /*@C
1292: TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed
1294: Collective
1296: Input Parameters:
1297: + tctx - the `TS` context
1298: . transform - the transform function
1299: . destroy - function to destroy the optional context, see `PetscCtxDestroyFn` for its calling sequence
1300: - ctx - optional context used by transform function
1302: Level: intermediate
1304: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`, `PetscCtxDestroyFn`
1305: @*/
1306: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscCtxDestroyFn *destroy, void *tctx)
1307: {
1308: PetscFunctionBegin;
1309: ctx->transform = transform;
1310: ctx->transformdestroy = destroy;
1311: ctx->transformctx = tctx;
1312: PetscFunctionReturn(PETSC_SUCCESS);
1313: }
1315: /*@C
1316: TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error
1317: in a time based line graph
1319: Collective
1321: Input Parameters:
1322: + ts - the `TS` context
1323: . step - current time-step
1324: . ptime - current time
1325: . u - current solution
1326: - dummy - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`
1328: Options Database Key:
1329: . -ts_monitor_lg_error - create a graphical monitor of error history
1331: Level: intermediate
1333: Notes:
1334: Each process in a parallel run displays its component errors in a separate window
1336: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1338: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1339: to be used during the TS integration.
1341: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1342: @*/
1343: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
1344: {
1345: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dummy;
1346: const PetscScalar *yy;
1347: Vec y;
1349: PetscFunctionBegin;
1350: if (!step) {
1351: PetscDrawAxis axis;
1352: PetscInt dim;
1353: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1354: PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
1355: PetscCall(VecGetLocalSize(u, &dim));
1356: PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1357: PetscCall(PetscDrawLGReset(ctx->lg));
1358: }
1359: PetscCall(VecDuplicate(u, &y));
1360: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1361: PetscCall(VecAXPY(y, -1.0, u));
1362: PetscCall(VecGetArrayRead(y, &yy));
1363: #if defined(PETSC_USE_COMPLEX)
1364: {
1365: PetscReal *yreal;
1366: PetscInt i, n;
1367: PetscCall(VecGetLocalSize(y, &n));
1368: PetscCall(PetscMalloc1(n, &yreal));
1369: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1370: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1371: PetscCall(PetscFree(yreal));
1372: }
1373: #else
1374: PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1375: #endif
1376: PetscCall(VecRestoreArrayRead(y, &yy));
1377: PetscCall(VecDestroy(&y));
1378: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1379: PetscCall(PetscDrawLGDraw(ctx->lg));
1380: PetscCall(PetscDrawLGSave(ctx->lg));
1381: }
1382: PetscFunctionReturn(PETSC_SUCCESS);
1383: }
1385: /*@C
1386: TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot
1388: Input Parameters:
1389: + ts - the `TS` context
1390: . step - current time-step
1391: . ptime - current time
1392: . u - current solution
1393: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()`
1395: Options Database Keys:
1396: + -ts_monitor_sp_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1397: . -ts_monitor_sp_swarm_retain <n> - Retain n old points so we can see the history, or -1 for all points
1398: . -ts_monitor_sp_swarm_multi_species <bool> - Color each species differently
1399: - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space
1401: Level: intermediate
1403: Notes:
1404: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1405: to be used during the `TS` integration.
1407: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
1408: @*/
1409: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1410: {
1411: TSMonitorSPCtx ctx = (TSMonitorSPCtx)dctx;
1412: PetscDraw draw;
1413: DM dm, cdm;
1414: const PetscScalar *yy;
1415: PetscInt Np, p, dim = 2, *species;
1416: PetscReal species_color;
1418: PetscFunctionBegin;
1419: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1420: PetscCall(TSGetDM(ts, &dm));
1421: if (!step) {
1422: PetscDrawAxis axis;
1423: PetscReal dmboxlower[2], dmboxupper[2];
1425: PetscCall(TSGetDM(ts, &dm));
1426: PetscCall(DMGetDimension(dm, &dim));
1427: PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
1428: PetscCall(DMSwarmGetCellDM(dm, &cdm));
1429: PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
1430: PetscCall(VecGetLocalSize(u, &Np));
1431: Np /= dim * 2;
1432: PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
1433: if (ctx->phase) {
1434: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
1435: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10));
1436: } else {
1437: PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
1438: PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
1439: }
1440: PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
1441: PetscCall(PetscDrawSPReset(ctx->sp));
1442: }
1443: if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species));
1444: PetscCall(VecGetLocalSize(u, &Np));
1445: Np /= dim * 2;
1446: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1447: PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
1448: if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
1449: PetscCall(PetscDrawFlush(draw));
1450: PetscCall(PetscDrawSPReset(ctx->sp));
1451: PetscCall(VecGetArrayRead(u, &yy));
1452: for (p = 0; p < Np; ++p) {
1453: PetscReal x, y;
1455: if (ctx->phase) {
1456: x = PetscRealPart(yy[p * dim * 2]);
1457: y = PetscRealPart(yy[p * dim * 2 + dim]);
1458: } else {
1459: x = PetscRealPart(yy[p * dim * 2]);
1460: y = PetscRealPart(yy[p * dim * 2 + 1]);
1461: }
1462: if (ctx->multispecies) {
1463: species_color = species[p] + 2;
1464: PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color));
1465: } else {
1466: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1467: }
1468: PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1469: }
1470: PetscCall(VecRestoreArrayRead(u, &yy));
1471: PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
1472: PetscCall(PetscDrawSPSave(ctx->sp));
1473: if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species));
1474: }
1475: PetscFunctionReturn(PETSC_SUCCESS);
1476: }
1478: /*@C
1479: TSMonitorHGSwarmSolution - Graphically displays histograms of `DMSWARM` particles
1481: Input Parameters:
1482: + ts - the `TS` context
1483: . step - current time-step
1484: . ptime - current time
1485: . u - current solution
1486: - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorHGCtxCreate()`
1488: Options Database Keys:
1489: + -ts_monitor_hg_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1490: . -ts_monitor_hg_swarm_species <num> - Number of species to histogram
1491: . -ts_monitor_hg_swarm_bins <num> - Number of histogram bins
1492: - -ts_monitor_hg_swarm_velocity <bool> - Plot in velocity space, as opposed to coordinate space
1494: Level: intermediate
1496: Note:
1497: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1498: to be used during the `TS` integration.
1500: .seealso: `TSMonitorSet()`
1501: @*/
1502: PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1503: {
1504: TSMonitorHGCtx ctx = (TSMonitorHGCtx)dctx;
1505: PetscDraw draw;
1506: DM sw;
1507: const PetscScalar *yy;
1508: PetscInt *species;
1509: PetscInt dim, d = 0, Np, p, Ns, s;
1511: PetscFunctionBegin;
1512: if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1513: PetscCall(TSGetDM(ts, &sw));
1514: PetscCall(DMGetDimension(sw, &dim));
1515: PetscCall(DMSwarmGetNumSpecies(sw, &Ns));
1516: Ns = PetscMin(Ns, ctx->Ns);
1517: PetscCall(VecGetLocalSize(u, &Np));
1518: Np /= dim * 2;
1519: if (!step) {
1520: PetscDrawAxis axis;
1521: char title[PETSC_MAX_PATH_LEN];
1523: for (s = 0; s < Ns; ++s) {
1524: PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis));
1525: PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s));
1526: if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N"));
1527: else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N"));
1528: }
1529: }
1530: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1531: PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species));
1532: for (s = 0; s < Ns; ++s) {
1533: PetscCall(PetscDrawHGReset(ctx->hg[s]));
1534: PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw));
1535: PetscCall(PetscDrawClear(draw));
1536: PetscCall(PetscDrawFlush(draw));
1537: }
1538: PetscCall(VecGetArrayRead(u, &yy));
1539: for (p = 0; p < Np; ++p) {
1540: const PetscInt s = species[p] < Ns ? species[p] : 0;
1541: PetscReal v;
1543: if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]);
1544: else v = PetscRealPart(yy[p * dim * 2 + d]);
1545: PetscCall(PetscDrawHGAddValue(ctx->hg[s], v));
1546: }
1547: PetscCall(VecRestoreArrayRead(u, &yy));
1548: for (s = 0; s < Ns; ++s) {
1549: PetscCall(PetscDrawHGDraw(ctx->hg[s]));
1550: PetscCall(PetscDrawHGSave(ctx->hg[s]));
1551: }
1552: PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species));
1553: }
1554: PetscFunctionReturn(PETSC_SUCCESS);
1555: }
1557: /*@C
1558: TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep
1560: Collective
1562: Input Parameters:
1563: + ts - the `TS` context
1564: . step - current time-step
1565: . ptime - current time
1566: . u - current solution
1567: - vf - unused context
1569: Options Database Key:
1570: . -ts_monitor_error - create a graphical monitor of error history
1572: Level: intermediate
1574: Notes:
1575: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1576: to be used during the `TS` integration.
1578: The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.
1580: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1581: @*/
1582: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1583: {
1584: DM dm;
1585: PetscDS ds = NULL;
1586: PetscInt Nf = -1, f;
1587: PetscBool flg;
1589: PetscFunctionBegin;
1590: PetscCall(TSGetDM(ts, &dm));
1591: if (dm) PetscCall(DMGetDS(dm, &ds));
1592: if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
1593: if (Nf <= 0) {
1594: Vec y;
1595: PetscReal nrm;
1597: PetscCall(VecDuplicate(u, &y));
1598: PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1599: PetscCall(VecAXPY(y, -1.0, u));
1600: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
1601: if (flg) {
1602: PetscCall(VecNorm(y, NORM_2, &nrm));
1603: PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
1604: }
1605: PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
1606: if (flg) PetscCall(VecView(y, vf->viewer));
1607: PetscCall(VecDestroy(&y));
1608: } else {
1609: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
1610: void **ctxs;
1611: Vec v;
1612: PetscReal ferrors[1];
1614: PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
1615: for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
1616: PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
1617: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04" PetscInt_FMT " time = %-8.4g \t L_2 Error: [", step, (double)ptime));
1618: for (f = 0; f < Nf; ++f) {
1619: if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
1620: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
1621: }
1622: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));
1624: PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));
1626: PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
1627: if (flg) {
1628: PetscCall(DMGetGlobalVector(dm, &v));
1629: PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
1630: PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
1631: PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
1632: PetscCall(DMRestoreGlobalVector(dm, &v));
1633: }
1634: PetscCall(PetscFree2(exactFuncs, ctxs));
1635: }
1636: PetscFunctionReturn(PETSC_SUCCESS);
1637: }
1639: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1640: {
1641: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1642: PetscReal x = ptime, y;
1643: PetscInt its;
1645: PetscFunctionBegin;
1646: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1647: if (!n) {
1648: PetscDrawAxis axis;
1649: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1650: PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
1651: PetscCall(PetscDrawLGReset(ctx->lg));
1652: ctx->snes_its = 0;
1653: }
1654: PetscCall(TSGetSNESIterations(ts, &its));
1655: y = its - ctx->snes_its;
1656: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1657: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1658: PetscCall(PetscDrawLGDraw(ctx->lg));
1659: PetscCall(PetscDrawLGSave(ctx->lg));
1660: }
1661: ctx->snes_its = its;
1662: PetscFunctionReturn(PETSC_SUCCESS);
1663: }
1665: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1666: {
1667: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1668: PetscReal x = ptime, y;
1669: PetscInt its;
1671: PetscFunctionBegin;
1672: if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1673: if (!n) {
1674: PetscDrawAxis axis;
1675: PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1676: PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
1677: PetscCall(PetscDrawLGReset(ctx->lg));
1678: ctx->ksp_its = 0;
1679: }
1680: PetscCall(TSGetKSPIterations(ts, &its));
1681: y = its - ctx->ksp_its;
1682: PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1683: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1684: PetscCall(PetscDrawLGDraw(ctx->lg));
1685: PetscCall(PetscDrawLGSave(ctx->lg));
1686: }
1687: ctx->ksp_its = its;
1688: PetscFunctionReturn(PETSC_SUCCESS);
1689: }
1691: /*@C
1692: TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()`
1694: Collective
1696: Input Parameter:
1697: . ts - the `TS` solver object
1699: Output Parameter:
1700: . ctx - the context
1702: Level: intermediate
1704: .seealso: [](ch_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1705: @*/
1706: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1707: {
1708: PetscFunctionBegin;
1709: PetscCall(PetscNew(ctx));
1710: PetscFunctionReturn(PETSC_SUCCESS);
1711: }
1713: /*@C
1714: TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution
1716: Collective
1718: Input Parameters:
1719: + ts - the `TS` context
1720: . step - current time-step
1721: . ptime - current time
1722: . u - current solution
1723: - dctx - the envelope context
1725: Options Database Key:
1726: . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time
1728: Level: intermediate
1730: Notes:
1731: After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope
1733: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1734: to be used during the `TS` integration.
1736: .seealso: [](ch_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1737: @*/
1738: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1739: {
1740: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;
1742: PetscFunctionBegin;
1743: if (!ctx->max) {
1744: PetscCall(VecDuplicate(u, &ctx->max));
1745: PetscCall(VecDuplicate(u, &ctx->min));
1746: PetscCall(VecCopy(u, ctx->max));
1747: PetscCall(VecCopy(u, ctx->min));
1748: } else {
1749: PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
1750: PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
1751: }
1752: PetscFunctionReturn(PETSC_SUCCESS);
1753: }
1755: /*@C
1756: TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution
1758: Collective
1760: Input Parameter:
1761: . ts - the `TS` context
1763: Output Parameters:
1764: + max - the maximum values
1765: - min - the minimum values
1767: Level: intermediate
1769: Notes:
1770: If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored
1772: .seealso: [](ch_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1773: @*/
1774: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1775: {
1776: PetscInt i;
1778: PetscFunctionBegin;
1779: if (max) *max = NULL;
1780: if (min) *min = NULL;
1781: for (i = 0; i < ts->numbermonitors; i++) {
1782: if (ts->monitor[i] == TSMonitorEnvelope) {
1783: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1784: if (max) *max = ctx->max;
1785: if (min) *min = ctx->min;
1786: break;
1787: }
1788: }
1789: PetscFunctionReturn(PETSC_SUCCESS);
1790: }
1792: /*@C
1793: TSMonitorEnvelopeCtxDestroy - Destroys a context that was created with `TSMonitorEnvelopeCtxCreate()`.
1795: Collective
1797: Input Parameter:
1798: . ctx - the monitor context
1800: Level: intermediate
1802: .seealso: [](ch_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1803: @*/
1804: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1805: {
1806: PetscFunctionBegin;
1807: PetscCall(VecDestroy(&(*ctx)->min));
1808: PetscCall(VecDestroy(&(*ctx)->max));
1809: PetscCall(PetscFree(*ctx));
1810: PetscFunctionReturn(PETSC_SUCCESS);
1811: }
1813: /*@C
1814: TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS`
1816: Not Collective
1818: Input Parameters:
1819: + ts - the `TS` context
1820: . step - current timestep
1821: . t - current time
1822: . U - current solution
1823: - vf - not used
1825: Options Database Key:
1826: + -ts_dmswarm_monitor_moments - Monitor moments of particle distribution
1827: - -ts_dmswarm_monitor_moments_interval - Interval of timesteps between monitor outputs
1829: Level: intermediate
1831: Notes:
1832: This requires a `DMSWARM` be attached to the `TS`.
1834: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1835: to be used during the TS integration.
1837: .seealso: [](ch_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1838: @*/
1839: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1840: {
1841: DM sw;
1842: const PetscScalar *u;
1843: PetscReal m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1844: PetscInt dim, d, Np, p;
1845: MPI_Comm comm;
1847: PetscFunctionBeginUser;
1848: (void)t;
1849: (void)vf;
1850: PetscCall(TSGetDM(ts, &sw));
1851: if (!sw || step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS);
1852: PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
1853: PetscCall(DMGetDimension(sw, &dim));
1854: PetscCall(VecGetLocalSize(U, &Np));
1855: Np /= dim;
1856: PetscCall(VecGetArrayRead(U, &u));
1857: for (p = 0; p < Np; ++p) {
1858: for (d = 0; d < dim; ++d) {
1859: totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1860: totMom[d] += PetscRealPart(u[p * dim + d]);
1861: }
1862: }
1863: PetscCall(VecRestoreArrayRead(U, &u));
1864: for (d = 0; d < dim; ++d) totMom[d] *= m;
1865: totE *= 0.5 * m;
1866: PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
1867: for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, " Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
1868: PetscCall(PetscPrintf(comm, "\n"));
1869: PetscFunctionReturn(PETSC_SUCCESS);
1870: }