Actual source code: tsmon.c

  1: #include <petsc/private/tsimpl.h>
  2: #include <petscdm.h>
  3: #include <petscds.h>
  4: #include <petscdmswarm.h>
  5: #include <petscdraw.h>

  7: /*@C
  8:    TSMonitor - Runs all user-provided monitor routines set using `TSMonitorSet()`

 10:    Collective

 12:    Input Parameters:
 13: +  ts - time stepping context obtained from `TSCreate()`
 14: .  step - step number that has just completed
 15: .  ptime - model time of the state
 16: -  u - state at the current model time

 18:    Level: developer

 20:    Notes:
 21:    `TSMonitor()` is typically used automatically within the time stepping implementations.
 22:    Users would almost never call this routine directly.

 24:    A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions

 26: .seealso: `TS`, `TSMonitorSet()`, `TSMonitorSetFromOptions()`
 27: @*/
 28: PetscErrorCode TSMonitor(TS ts, PetscInt step, PetscReal ptime, Vec u)
 29: {
 30:   DM       dm;
 31:   PetscInt i, n = ts->numbermonitors;

 33:   PetscFunctionBegin;

 37:   PetscCall(TSGetDM(ts, &dm));
 38:   PetscCall(DMSetOutputSequenceNumber(dm, step, ptime));

 40:   PetscCall(VecLockReadPush(u));
 41:   for (i = 0; i < n; i++) PetscCall((*ts->monitor[i])(ts, step, ptime, u, ts->monitorcontext[i]));
 42:   PetscCall(VecLockReadPop(u));
 43:   PetscFunctionReturn(PETSC_SUCCESS);
 44: }

 46: /*@C
 47:    TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user

 49:    Collective

 51:    Input Parameters:
 52: +  ts - `TS` object you wish to monitor
 53: .  name - the monitor type one is seeking
 54: .  help - message indicating what monitoring is done
 55: .  manual - manual page for the monitor
 56: .  monitor - the monitor function
 57: -  monitorsetup - a function that is called once ONLY if the user selected this monitor that may set additional features of the `TS` or `PetscViewer` objects

 59:    Level: developer

 61: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `PetscOptionsGetViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
 62:           `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
 63:           `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`,
 64:           `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
 65:           `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
 66:           `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
 67:           `PetscOptionsFList()`, `PetscOptionsEList()`
 68: @*/
 69: PetscErrorCode TSMonitorSetFromOptions(TS ts, const char name[], const char help[], const char manual[], PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *), PetscErrorCode (*monitorsetup)(TS, PetscViewerAndFormat *))
 70: {
 71:   PetscViewer       viewer;
 72:   PetscViewerFormat format;
 73:   PetscBool         flg;

 75:   PetscFunctionBegin;
 76:   PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)ts), ((PetscObject)ts)->options, ((PetscObject)ts)->prefix, name, &viewer, &format, &flg));
 77:   if (flg) {
 78:     PetscViewerAndFormat *vf;
 79:     char                  interval_key[1024];
 80:     PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
 81:     PetscCall(PetscSNPrintf(interval_key, sizeof interval_key, "%s_interval", name));
 82:     PetscCall(PetscOptionsGetInt(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, interval_key, &vf->view_interval, NULL));
 83:     PetscCall(PetscObjectDereference((PetscObject)viewer));
 84:     if (monitorsetup) PetscCall((*monitorsetup)(ts, vf));
 85:     PetscCall(TSMonitorSet(ts, (PetscErrorCode(*)(TS, PetscInt, PetscReal, Vec, void *))monitor, vf, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy));
 86:   }
 87:   PetscFunctionReturn(PETSC_SUCCESS);
 88: }

 90: /*@C
 91:    TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
 92:    timestep to display the iteration's  progress.

 94:    Logically Collective

 96:    Input Parameters:
 97: +  ts - the `TS` context obtained from `TSCreate()`
 98: .  monitor - monitoring routine
 99: .  mctx - [optional] user-defined context for private data for the monitor routine (use `NULL` if no context is desired)
100: -  monitordestroy - [optional] routine that frees monitor context (may be `NULL`)

102:    Calling sequence of monitor:
103: $    PetscErrorCode monitor(TS ts,PetscInt steps,PetscReal time,Vec u,void *mctx)

105: +    ts - the TS context
106: .    steps - iteration number (after the final time step the monitor routine may be called with a step of -1, this indicates the solution has been interpolated to this time)
107: .    time - current time
108: .    u - current iterate
109: -    mctx - [optional] monitoring context

111:    Level: intermediate

113:    Note:
114:    This routine adds an additional monitor to the list of monitors that already has been loaded.

116:    Fortran Note:
117:    Only a single monitor function can be set for each `TS` object

119: .seealso: [](chapter_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
120:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
121:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
122: @*/
123: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, void *), void *mctx, PetscErrorCode (*mdestroy)(void **))
124: {
125:   PetscInt  i;
126:   PetscBool identical;

128:   PetscFunctionBegin;
130:   for (i = 0; i < ts->numbermonitors; i++) {
131:     PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))monitor, mctx, mdestroy, (PetscErrorCode(*)(void))ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical));
132:     if (identical) PetscFunctionReturn(PETSC_SUCCESS);
133:   }
134:   PetscCheck(ts->numbermonitors < MAXTSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set");
135:   ts->monitor[ts->numbermonitors]          = monitor;
136:   ts->monitordestroy[ts->numbermonitors]   = mdestroy;
137:   ts->monitorcontext[ts->numbermonitors++] = (void *)mctx;
138:   PetscFunctionReturn(PETSC_SUCCESS);
139: }

141: /*@C
142:    TSMonitorCancel - Clears all the monitors that have been set on a time-step object.

144:    Logically Collective

146:    Input Parameters:
147: .  ts - the `TS` context obtained from `TSCreate()`

149:    Level: intermediate

151:    Note:
152:    There is no way to remove a single, specific monitor.

154: .seealso: [](chapter_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()`
155: @*/
156: PetscErrorCode TSMonitorCancel(TS ts)
157: {
158:   PetscInt i;

160:   PetscFunctionBegin;
162:   for (i = 0; i < ts->numbermonitors; i++) {
163:     if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i]));
164:   }
165:   ts->numbermonitors = 0;
166:   PetscFunctionReturn(PETSC_SUCCESS);
167: }

169: /*@C
170:    TSMonitorDefault - The default monitor, prints the timestep and time for each step

172:    Options Database Key:
173: .  -ts_monitor - monitors the time integration

175:    Level: intermediate

177:    Notes:
178:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
179:    to be used during the `TS` integration.

181: .seealso: [](chapter_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
182:           `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
183:           `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
184: @*/
185: PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
186: {
187:   PetscViewer viewer = vf->viewer;
188:   PetscBool   iascii, ibinary;

190:   PetscFunctionBegin;
192:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
193:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary));
194:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
195:   if (iascii) {
196:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
197:     if (step == -1) { /* this indicates it is an interpolated solution */
198:       PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps));
199:     } else {
200:       PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n"));
201:     }
202:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
203:   } else if (ibinary) {
204:     PetscMPIInt rank;
205:     PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank));
206:     if (rank == 0) {
207:       PetscBool skipHeader;
208:       PetscInt  classid = REAL_FILE_CLASSID;

210:       PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader));
211:       if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
212:       PetscCall(PetscRealView(1, &ptime, viewer));
213:     } else {
214:       PetscCall(PetscRealView(0, &ptime, viewer));
215:     }
216:   }
217:   PetscCall(PetscViewerPopFormat(viewer));
218:   PetscFunctionReturn(PETSC_SUCCESS);
219: }

221: /*@C
222:    TSMonitorExtreme - Prints the extreme values of the solution at each timestep

224:    Level: intermediate

226:    Note:
227:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
228:    to be used during the `TS` integration.

230: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`
231: @*/
232: PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
233: {
234:   PetscViewer viewer = vf->viewer;
235:   PetscBool   iascii;
236:   PetscReal   max, min;

238:   PetscFunctionBegin;
240:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
241:   PetscCall(PetscViewerPushFormat(viewer, vf->format));
242:   if (iascii) {
243:     PetscCall(VecMax(v, NULL, &max));
244:     PetscCall(VecMin(v, NULL, &min));
245:     PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel));
246:     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));
247:     PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel));
248:   }
249:   PetscCall(PetscViewerPopFormat(viewer));
250:   PetscFunctionReturn(PETSC_SUCCESS);
251: }

253: /*@C
254:    TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with
255:    `TS` to monitor the solution process graphically in various ways

257:    Collective

259:    Input Parameters:
260: +  host - the X display to open, or `NULL` for the local machine
261: .  label - the title to put in the title bar
262: .  x, y - the screen coordinates of the upper left coordinate of the window
263: .  m, n - the screen width and height in pixels
264: -  howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time

266:    Output Parameter:
267: .  ctx - the context

269:    Options Database Keys:
270: +  -ts_monitor_lg_timestep - automatically sets line graph monitor
271: +  -ts_monitor_lg_timestep_log - automatically sets line graph monitor
272: .  -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`)
273: .  -ts_monitor_lg_error -  monitor the error
274: .  -ts_monitor_lg_ksp_iterations - monitor the number of `KSP` iterations needed for each timestep
275: .  -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep
276: -  -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true

278:    Level: intermediate

280:    Notes:
281:    Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed.

283:    One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()`

285:    Many of the functions that control the monitoring have two forms: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a `TS` object as the
286:    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
287:    as the first argument.

289:    One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()`

291: .seealso: [](chapter_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
292:           `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
293:           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
294:           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
295:           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
296: @*/
297: PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
298: {
299:   PetscDraw draw;

301:   PetscFunctionBegin;
302:   PetscCall(PetscNew(ctx));
303:   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
304:   PetscCall(PetscDrawSetFromOptions(draw));
305:   PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg));
306:   PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg));
307:   PetscCall(PetscDrawDestroy(&draw));
308:   (*ctx)->howoften = howoften;
309:   PetscFunctionReturn(PETSC_SUCCESS);
310: }

312: PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx)
313: {
314:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
315:   PetscReal      x   = ptime, y;

317:   PetscFunctionBegin;
318:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */
319:   if (!step) {
320:     PetscDrawAxis axis;
321:     const char   *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
322:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
323:     PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel));
324:     PetscCall(PetscDrawLGReset(ctx->lg));
325:   }
326:   PetscCall(TSGetTimeStep(ts, &y));
327:   if (ctx->semilogy) y = PetscLog10Real(y);
328:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
329:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
330:     PetscCall(PetscDrawLGDraw(ctx->lg));
331:     PetscCall(PetscDrawLGSave(ctx->lg));
332:   }
333:   PetscFunctionReturn(PETSC_SUCCESS);
334: }

336: /*@C
337:    TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`.

339:    Collective

341:    Input Parameter:
342: .  ctx - the monitor context

344:    Level: intermediate

346:    Note:
347:    Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()`

349: .seealso: [](chapter_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep();`
350: @*/
351: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
352: {
353:   PetscFunctionBegin;
354:   if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)((*ctx)->transformctx));
355:   PetscCall(PetscDrawLGDestroy(&(*ctx)->lg));
356:   PetscCall(PetscStrArrayDestroy(&(*ctx)->names));
357:   PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames));
358:   PetscCall(PetscFree((*ctx)->displayvariables));
359:   PetscCall(PetscFree((*ctx)->displayvalues));
360:   PetscCall(PetscFree(*ctx));
361:   PetscFunctionReturn(PETSC_SUCCESS);
362: }

364: /* Creates a TSMonitorSPCtx for use with DMSwarm particle visualizations */
365: 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)
366: {
367:   PetscDraw draw;

369:   PetscFunctionBegin;
370:   PetscCall(PetscNew(ctx));
371:   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
372:   PetscCall(PetscDrawSetFromOptions(draw));
373:   PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp));
374:   PetscCall(PetscDrawDestroy(&draw));
375:   (*ctx)->howoften     = howoften;
376:   (*ctx)->retain       = retain;
377:   (*ctx)->phase        = phase;
378:   (*ctx)->multispecies = multispecies;
379:   PetscFunctionReturn(PETSC_SUCCESS);
380: }

382: /* Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate */
383: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
384: {
385:   PetscFunctionBegin;

387:   PetscCall(PetscDrawSPDestroy(&(*ctx)->sp));
388:   PetscCall(PetscFree(*ctx));
389:   PetscFunctionReturn(PETSC_SUCCESS);
390: }

392: /* Creates a TSMonitorHGCtx for use with DMSwarm particle visualizations */
393: 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)
394: {
395:   PetscDraw draw;
396:   PetscInt  s;

398:   PetscFunctionBegin;
399:   PetscCall(PetscNew(ctx));
400:   PetscCall(PetscMalloc1(Ns, &(*ctx)->hg));
401:   for (s = 0; s < Ns; ++s) {
402:     PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw));
403:     PetscCall(PetscDrawSetFromOptions(draw));
404:     PetscCall(PetscDrawHGCreate(draw, Nb, &(*ctx)->hg[s]));
405:     PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE));
406:     PetscCall(PetscDrawDestroy(&draw));
407:   }
408:   (*ctx)->howoften = howoften;
409:   (*ctx)->Ns       = Ns;
410:   (*ctx)->velocity = velocity;
411:   PetscFunctionReturn(PETSC_SUCCESS);
412: }

414: /* Destroys a TSMonitorHGCtx that was created with TSMonitorHGCtxCreate */
415: PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx)
416: {
417:   PetscInt s;

419:   PetscFunctionBegin;
420:   for (s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s]));
421:   PetscCall(PetscFree((*ctx)->hg));
422:   PetscCall(PetscFree(*ctx));
423:   PetscFunctionReturn(PETSC_SUCCESS);
424: }

426: /*@C
427:    TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling
428:    `VecView()` for the solution at each timestep

430:    Collective

432:    Input Parameters:
433: +  ts - the `TS` context
434: .  step - current time-step
435: .  ptime - current time
436: -  dummy - either a viewer or `NULL`

438:    Options Database Keys:
439: +   -ts_monitor_draw_solution - draw the solution at each time-step
440: -   -ts_monitor_draw_solution_initial - show initial solution as well as current solution

442:    Level: intermediate

444:    Notes:
445:    The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial`
446:    will look bad

448:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with
449:    `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration.

451: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()`
452: @*/
453: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
454: {
455:   TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
456:   PetscDraw        draw;

458:   PetscFunctionBegin;
459:   if (!step && ictx->showinitial) {
460:     if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution));
461:     PetscCall(VecCopy(u, ictx->initialsolution));
462:   }
463:   if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);

465:   if (ictx->showinitial) {
466:     PetscReal pause;
467:     PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause));
468:     PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0));
469:     PetscCall(VecView(ictx->initialsolution, ictx->viewer));
470:     PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause));
471:     PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE));
472:   }
473:   PetscCall(VecView(u, ictx->viewer));
474:   if (ictx->showtimestepandtime) {
475:     PetscReal xl, yl, xr, yr, h;
476:     char      time[32];

478:     PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
479:     PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime));
480:     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
481:     h = yl + .95 * (yr - yl);
482:     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
483:     PetscCall(PetscDrawFlush(draw));
484:   }

486:   if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE));
487:   PetscFunctionReturn(PETSC_SUCCESS);
488: }

490: /*@C
491:    TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram

493:    Collective

495:    Input Parameters:
496: +  ts - the `TS` context
497: .  step - current time-step
498: .  ptime - current time
499: -  dummy - either a viewer or `NULL`

501:    Level: intermediate

503:    Notes:
504:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
505:    to be used during the `TS` integration.

507: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
508: @*/
509: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
510: {
511:   TSMonitorDrawCtx   ictx = (TSMonitorDrawCtx)dummy;
512:   PetscDraw          draw;
513:   PetscDrawAxis      axis;
514:   PetscInt           n;
515:   PetscMPIInt        size;
516:   PetscReal          U0, U1, xl, yl, xr, yr, h;
517:   char               time[32];
518:   const PetscScalar *U;

520:   PetscFunctionBegin;
521:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size));
522:   PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs");
523:   PetscCall(VecGetSize(u, &n));
524:   PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns");

526:   PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw));
527:   PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis));
528:   PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr));
529:   if (!step) {
530:     PetscCall(PetscDrawClear(draw));
531:     PetscCall(PetscDrawAxisDraw(axis));
532:   }

534:   PetscCall(VecGetArrayRead(u, &U));
535:   U0 = PetscRealPart(U[0]);
536:   U1 = PetscRealPart(U[1]);
537:   PetscCall(VecRestoreArrayRead(u, &U));
538:   if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS);

540:   PetscDrawCollectiveBegin(draw);
541:   PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK));
542:   if (ictx->showtimestepandtime) {
543:     PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr));
544:     PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime));
545:     h = yl + .95 * (yr - yl);
546:     PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time));
547:   }
548:   PetscDrawCollectiveEnd(draw);
549:   PetscCall(PetscDrawFlush(draw));
550:   PetscCall(PetscDrawPause(draw));
551:   PetscCall(PetscDrawSave(draw));
552:   PetscFunctionReturn(PETSC_SUCCESS);
553: }

555: /*@C
556:    TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()`

558:    Collective

560:    Input Parameters:
561: .    ctx - the monitor context

563:    Level: intermediate

565: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx`
566: @*/
567: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
568: {
569:   PetscFunctionBegin;
570:   PetscCall(PetscViewerDestroy(&(*ictx)->viewer));
571:   PetscCall(VecDestroy(&(*ictx)->initialsolution));
572:   PetscCall(PetscFree(*ictx));
573:   PetscFunctionReturn(PETSC_SUCCESS);
574: }

576: /*@C
577:    TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx`

579:    Collective

581:    Input Parameter:
582: .    ts - time-step context

584:    Output Parameter:
585: .    ctx - the monitor context

587:    Options Database Keys:
588: +   -ts_monitor_draw_solution - draw the solution at each time-step
589: -   -ts_monitor_draw_solution_initial - show initial solution as well as current solution

591:    Level: intermediate

593:    Note:
594:    The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`.

596: .seealso: [](chapter_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()`
597: @*/
598: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
599: {
600:   PetscFunctionBegin;
601:   PetscCall(PetscNew(ctx));
602:   PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer));
603:   PetscCall(PetscViewerSetFromOptions((*ctx)->viewer));

605:   (*ctx)->howoften    = howoften;
606:   (*ctx)->showinitial = PETSC_FALSE;
607:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL));

609:   (*ctx)->showtimestepandtime = PETSC_FALSE;
610:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL));
611:   PetscFunctionReturn(PETSC_SUCCESS);
612: }

614: /*@C
615:    TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling
616:    `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep

618:    Collective

620:    Input Parameters:
621: +  ts - the `TS` context
622: .  step - current time-step
623: .  ptime - current time
624: -  dummy - either a viewer or `NULL`

626:    Options Database Key:
627: .  -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`

629:    Level: intermediate

631:    Note:
632:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
633:    to be used during the `TS` integration.

635: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
636: @*/
637: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
638: {
639:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
640:   PetscViewer      viewer = ctx->viewer;
641:   Vec              work;

643:   PetscFunctionBegin;
644:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
645:   PetscCall(VecDuplicate(u, &work));
646:   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
647:   PetscCall(VecView(work, viewer));
648:   PetscCall(VecDestroy(&work));
649:   PetscFunctionReturn(PETSC_SUCCESS);
650: }

652: /*@C
653:    TSMonitorDrawError - Monitors progress of the `TS` solvers by calling
654:    `VecView()` for the error at each timestep

656:    Collective

658:    Input Parameters:
659: +  ts - the `TS` context
660: .  step - current time-step
661: .  ptime - current time
662: -  dummy - either a viewer or `NULL`

664:    Options Database Key:
665: .  -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()`

667:    Level: intermediate

669:    Notes:
670:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
671:    to be used during the `TS` integration.

673: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
674: @*/
675: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
676: {
677:   TSMonitorDrawCtx ctx    = (TSMonitorDrawCtx)dummy;
678:   PetscViewer      viewer = ctx->viewer;
679:   Vec              work;

681:   PetscFunctionBegin;
682:   if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS);
683:   PetscCall(VecDuplicate(u, &work));
684:   PetscCall(TSComputeSolutionFunction(ts, ptime, work));
685:   PetscCall(VecAXPY(work, -1.0, u));
686:   PetscCall(VecView(work, viewer));
687:   PetscCall(VecDestroy(&work));
688:   PetscFunctionReturn(PETSC_SUCCESS);
689: }

691: /*@C
692:    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

694:    Collective

696:    Input Parameters:
697: +  ts - the `TS` context
698: .  step - current time-step
699: .  ptime - current time
700: .  u - current state
701: -  vf - viewer and its format

703:    Level: intermediate

705:    Notes:
706:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
707:    to be used during the `TS` integration.

709: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
710: @*/
711: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
712: {
713:   PetscFunctionBegin;
714:   if (vf->view_interval > 0 && !ts->reason && step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS);
715:   PetscCall(PetscViewerPushFormat(vf->viewer, vf->format));
716:   PetscCall(VecView(u, vf->viewer));
717:   PetscCall(PetscViewerPopFormat(vf->viewer));
718:   PetscFunctionReturn(PETSC_SUCCESS);
719: }

721: /*@C
722:    TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at each timestep.

724:    Collective

726:    Input Parameters:
727: +  ts - the `TS` context
728: .  step - current time-step
729: .  ptime - current time
730: .  u - current state
731: -  filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")

733:    Level: intermediate

735:    Notes:
736:    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.
737:    These are named according to the file name template.

739:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
740:    to be used during the `TS` integration.

742: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
743: @*/
744: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, void *filenametemplate)
745: {
746:   char        filename[PETSC_MAX_PATH_LEN];
747:   PetscViewer viewer;

749:   PetscFunctionBegin;
750:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
751:   PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)filenametemplate, step));
752:   PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer));
753:   PetscCall(VecView(u, viewer));
754:   PetscCall(PetscViewerDestroy(&viewer));
755:   PetscFunctionReturn(PETSC_SUCCESS);
756: }

758: /*@C
759:    TSMonitorSolutionVTKDestroy - Destroy filename template string created for use with `TSMonitorSolutionVTK()`

761:    Not Collective

763:    Input Parameters:
764: .  filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")

766:    Level: intermediate

768:    Note:
769:    This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`.

771: .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()`
772: @*/
773: PetscErrorCode TSMonitorSolutionVTKDestroy(void *filenametemplate)
774: {
775:   PetscFunctionBegin;
776:   PetscCall(PetscFree(*(char **)filenametemplate));
777:   PetscFunctionReturn(PETSC_SUCCESS);
778: }

780: /*@C
781:    TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector
782:        in a time based line graph

784:    Collective

786:    Input Parameters:
787: +  ts - the `TS` context
788: .  step - current time-step
789: .  ptime - current time
790: .  u - current solution
791: -  dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()`

793:    Options Database Key:
794: .   -ts_monitor_lg_solution_variables - enable monitor of lg solution variables

796:    Level: intermediate

798:    Notes:
799:    Each process in a parallel run displays its component solutions in a separate window

801:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
802:    to be used during the `TS` integration.

804: .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
805:           `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
806:           `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
807:           `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
808: @*/
809: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
810: {
811:   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)dctx;
812:   const PetscScalar *yy;
813:   Vec                v;

815:   PetscFunctionBegin;
816:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
817:   if (!step) {
818:     PetscDrawAxis axis;
819:     PetscInt      dim;
820:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
821:     PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution"));
822:     if (!ctx->names) {
823:       PetscBool flg;
824:       /* user provides names of variables to plot but no names has been set so assume names are integer values */
825:       PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg));
826:       if (flg) {
827:         PetscInt i, n;
828:         char   **names;
829:         PetscCall(VecGetSize(u, &n));
830:         PetscCall(PetscMalloc1(n + 1, &names));
831:         for (i = 0; i < n; i++) {
832:           PetscCall(PetscMalloc1(5, &names[i]));
833:           PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i));
834:         }
835:         names[n]   = NULL;
836:         ctx->names = names;
837:       }
838:     }
839:     if (ctx->names && !ctx->displaynames) {
840:       char    **displaynames;
841:       PetscBool flg;
842:       PetscCall(VecGetLocalSize(u, &dim));
843:       PetscCall(PetscCalloc1(dim + 1, &displaynames));
844:       PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg));
845:       if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames));
846:       PetscCall(PetscStrArrayDestroy(&displaynames));
847:     }
848:     if (ctx->displaynames) {
849:       PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables));
850:       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames));
851:     } else if (ctx->names) {
852:       PetscCall(VecGetLocalSize(u, &dim));
853:       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
854:       PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names));
855:     } else {
856:       PetscCall(VecGetLocalSize(u, &dim));
857:       PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
858:     }
859:     PetscCall(PetscDrawLGReset(ctx->lg));
860:   }

862:   if (!ctx->transform) v = u;
863:   else PetscCall((*ctx->transform)(ctx->transformctx, u, &v));
864:   PetscCall(VecGetArrayRead(v, &yy));
865:   if (ctx->displaynames) {
866:     PetscInt i;
867:     for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
868:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues));
869:   } else {
870: #if defined(PETSC_USE_COMPLEX)
871:     PetscInt   i, n;
872:     PetscReal *yreal;
873:     PetscCall(VecGetLocalSize(v, &n));
874:     PetscCall(PetscMalloc1(n, &yreal));
875:     for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
876:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
877:     PetscCall(PetscFree(yreal));
878: #else
879:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
880: #endif
881:   }
882:   PetscCall(VecRestoreArrayRead(v, &yy));
883:   if (ctx->transform) PetscCall(VecDestroy(&v));

885:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
886:     PetscCall(PetscDrawLGDraw(ctx->lg));
887:     PetscCall(PetscDrawLGSave(ctx->lg));
888:   }
889:   PetscFunctionReturn(PETSC_SUCCESS);
890: }

892: /*@C
893:    TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot

895:    Collective

897:    Input Parameters:
898: +  ts - the `TS` context
899: -  names - the names of the components, final string must be `NULL`

901:    Level: intermediate

903:    Notes:
904:     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

906: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
907: @*/
908: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
909: {
910:   PetscInt i;

912:   PetscFunctionBegin;
913:   for (i = 0; i < ts->numbermonitors; i++) {
914:     if (ts->monitor[i] == TSMonitorLGSolution) {
915:       PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names));
916:       break;
917:     }
918:   }
919:   PetscFunctionReturn(PETSC_SUCCESS);
920: }

922: /*@C
923:    TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot

925:    Collective

927:    Input Parameters:
928: +  ts - the `TS` context
929: -  names - the names of the components, final string must be `NULL`

931:    Level: intermediate

933: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
934: @*/
935: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
936: {
937:   PetscFunctionBegin;
938:   PetscCall(PetscStrArrayDestroy(&ctx->names));
939:   PetscCall(PetscStrArrayallocpy(names, &ctx->names));
940:   PetscFunctionReturn(PETSC_SUCCESS);
941: }

943: /*@C
944:    TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot

946:    Collective

948:    Input Parameter:
949: .  ts - the `TS` context

951:    Output Parameter:
952: .  names - the names of the components, final string must be `NULL`

954:    Level: intermediate

956:    Note:
957:     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

959: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
960: @*/
961: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
962: {
963:   PetscInt i;

965:   PetscFunctionBegin;
966:   *names = NULL;
967:   for (i = 0; i < ts->numbermonitors; i++) {
968:     if (ts->monitor[i] == TSMonitorLGSolution) {
969:       TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
970:       *names             = (const char *const *)ctx->names;
971:       break;
972:     }
973:   }
974:   PetscFunctionReturn(PETSC_SUCCESS);
975: }

977: /*@C
978:    TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor

980:    Collective

982:    Input Parameters:
983: +  ctx - the `TSMonitorLG` context
984: -  displaynames - the names of the components, final string must be `NULL`

986:    Level: intermediate

988: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
989: @*/
990: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
991: {
992:   PetscInt j = 0, k;

994:   PetscFunctionBegin;
995:   if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS);
996:   PetscCall(PetscStrArrayDestroy(&ctx->displaynames));
997:   PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames));
998:   while (displaynames[j]) j++;
999:   ctx->ndisplayvariables = j;
1000:   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables));
1001:   PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues));
1002:   j = 0;
1003:   while (displaynames[j]) {
1004:     k = 0;
1005:     while (ctx->names[k]) {
1006:       PetscBool flg;
1007:       PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg));
1008:       if (flg) {
1009:         ctx->displayvariables[j] = k;
1010:         break;
1011:       }
1012:       k++;
1013:     }
1014:     j++;
1015:   }
1016:   PetscFunctionReturn(PETSC_SUCCESS);
1017: }

1019: /*@C
1020:    TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor

1022:    Collective

1024:    Input Parameters:
1025: +  ts - the `TS` context
1026: -  displaynames - the names of the components, final string must be `NULL`

1028:    Level: intermediate

1030:    Note:
1031:     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

1033: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
1034: @*/
1035: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
1036: {
1037:   PetscInt i;

1039:   PetscFunctionBegin;
1040:   for (i = 0; i < ts->numbermonitors; i++) {
1041:     if (ts->monitor[i] == TSMonitorLGSolution) {
1042:       PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames));
1043:       break;
1044:     }
1045:   }
1046:   PetscFunctionReturn(PETSC_SUCCESS);
1047: }

1049: /*@C
1050:    TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed

1052:    Collective

1054:    Input Parameters:
1055: +  ts - the `TS` context
1056: .  transform - the transform function
1057: .  destroy - function to destroy the optional context
1058: -  ctx - optional context used by transform function

1060:    Level: intermediate

1062:    Note:
1063:     If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored

1065: .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()`
1066: @*/
1067: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1068: {
1069:   PetscInt i;

1071:   PetscFunctionBegin;
1072:   for (i = 0; i < ts->numbermonitors; i++) {
1073:     if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx));
1074:   }
1075:   PetscFunctionReturn(PETSC_SUCCESS);
1076: }

1078: /*@C
1079:    TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed

1081:    Collective

1083:    Input Parameters:
1084: +  ts - the `TS` context
1085: .  transform - the transform function
1086: .  destroy - function to destroy the optional context
1087: -  ctx - optional context used by transform function

1089:    Level: intermediate

1091: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`
1092: @*/
1093: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1094: {
1095:   PetscFunctionBegin;
1096:   ctx->transform        = transform;
1097:   ctx->transformdestroy = destroy;
1098:   ctx->transformctx     = tctx;
1099:   PetscFunctionReturn(PETSC_SUCCESS);
1100: }

1102: /*@C
1103:    TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error
1104:        in a time based line graph

1106:    Collective

1108:    Input Parameters:
1109: +  ts - the `TS` context
1110: .  step - current time-step
1111: .  ptime - current time
1112: .  u - current solution
1113: -  dctx - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()`

1115:    Options Database Key:
1116: .  -ts_monitor_lg_error - create a graphical monitor of error history

1118:    Level: intermediate

1120:    Notes:
1121:     Each process in a parallel run displays its component errors in a separate window

1123:    The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.

1125:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1126:    to be used during the TS integration.

1128: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1129: @*/
1130: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
1131: {
1132:   TSMonitorLGCtx     ctx = (TSMonitorLGCtx)dummy;
1133:   const PetscScalar *yy;
1134:   Vec                y;

1136:   PetscFunctionBegin;
1137:   if (!step) {
1138:     PetscDrawAxis axis;
1139:     PetscInt      dim;
1140:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1141:     PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error"));
1142:     PetscCall(VecGetLocalSize(u, &dim));
1143:     PetscCall(PetscDrawLGSetDimension(ctx->lg, dim));
1144:     PetscCall(PetscDrawLGReset(ctx->lg));
1145:   }
1146:   PetscCall(VecDuplicate(u, &y));
1147:   PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1148:   PetscCall(VecAXPY(y, -1.0, u));
1149:   PetscCall(VecGetArrayRead(y, &yy));
1150: #if defined(PETSC_USE_COMPLEX)
1151:   {
1152:     PetscReal *yreal;
1153:     PetscInt   i, n;
1154:     PetscCall(VecGetLocalSize(y, &n));
1155:     PetscCall(PetscMalloc1(n, &yreal));
1156:     for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1157:     PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal));
1158:     PetscCall(PetscFree(yreal));
1159:   }
1160: #else
1161:   PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy));
1162: #endif
1163:   PetscCall(VecRestoreArrayRead(y, &yy));
1164:   PetscCall(VecDestroy(&y));
1165:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1166:     PetscCall(PetscDrawLGDraw(ctx->lg));
1167:     PetscCall(PetscDrawLGSave(ctx->lg));
1168:   }
1169:   PetscFunctionReturn(PETSC_SUCCESS);
1170: }

1172: /*@C
1173:    TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot

1175:    Input Parameters:
1176: +  ts - the `TS` context
1177: .  step - current time-step
1178: .  ptime - current time
1179: .  u - current solution
1180: -  dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()`

1182:    Options Database:
1183: + -ts_monitor_sp_swarm <n>                  - Monitor the solution every n steps, or -1 for plotting only the final solution
1184: . -ts_monitor_sp_swarm_retain <n>           - Retain n old points so we can see the history, or -1 for all points
1185: . -ts_monitor_sp_swarm_multi_species <bool> - Color each species differently
1186: - -ts_monitor_sp_swarm_phase <bool>         - Plot in phase space, as opposed to coordinate space

1188:    Options Database Keys:
1189: + -ts_monitor_sp_swarm <n>          - Monitor the solution every n steps, or -1 for plotting only the final solution
1190: . -ts_monitor_sp_swarm_retain <n>   - Retain n old points so we can see the history, or -1 for all points
1191: - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space

1193:    Level: intermediate

1195:    Notes:
1196:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1197:    to be used during the `TS` integration.

1199: .seealso: [](chapter_ts), `TS`, `TSMonitoSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()`
1200: @*/
1201: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1202: {
1203:   TSMonitorSPCtx     ctx = (TSMonitorSPCtx)dctx;
1204:   PetscDraw          draw;
1205:   DM                 dm, cdm;
1206:   const PetscScalar *yy;
1207:   PetscInt           Np, p, dim = 2, *species;
1208:   PetscReal          species_color;

1210:   PetscFunctionBegin;
1211:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1212:   PetscCall(TSGetDM(ts, &dm));
1213:   if (!step) {
1214:     PetscDrawAxis axis;
1215:     PetscReal     dmboxlower[2], dmboxupper[2];

1217:     PetscCall(TSGetDM(ts, &dm));
1218:     PetscCall(DMGetDimension(dm, &dim));
1219:     PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields");
1220:     PetscCall(DMSwarmGetCellDM(dm, &cdm));
1221:     PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper));
1222:     PetscCall(VecGetLocalSize(u, &Np));
1223:     Np /= dim * 2;
1224:     PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis));
1225:     if (ctx->phase) {
1226:       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V"));
1227:       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10));
1228:     } else {
1229:       PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y"));
1230:       PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]));
1231:     }
1232:     PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE));
1233:     PetscCall(PetscDrawSPReset(ctx->sp));
1234:   }
1235:   if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species));
1236:   PetscCall(VecGetLocalSize(u, &Np));
1237:   Np /= dim * 2;
1238:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1239:     PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw));
1240:     if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw));
1241:     PetscCall(PetscDrawFlush(draw));
1242:     PetscCall(PetscDrawSPReset(ctx->sp));
1243:     PetscCall(VecGetArrayRead(u, &yy));
1244:     for (p = 0; p < Np; ++p) {
1245:       PetscReal x, y;

1247:       if (ctx->phase) {
1248:         x = PetscRealPart(yy[p * dim * 2]);
1249:         y = PetscRealPart(yy[p * dim * 2 + dim]);
1250:       } else {
1251:         x = PetscRealPart(yy[p * dim * 2]);
1252:         y = PetscRealPart(yy[p * dim * 2 + 1]);
1253:       }
1254:       if (ctx->multispecies) {
1255:         species_color = species[p] + 2;
1256:         PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color));
1257:       } else {
1258:         PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1259:       }
1260:       PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y));
1261:     }
1262:     PetscCall(VecRestoreArrayRead(u, &yy));
1263:     PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE));
1264:     PetscCall(PetscDrawSPSave(ctx->sp));
1265:     if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species));
1266:   }
1267:   PetscFunctionReturn(PETSC_SUCCESS);
1268: }

1270: /*@C
1271:    TSMonitorHGSwarmSolution - Graphically displays histograms of DMSwarm particles

1273:    Input Parameters:
1274: +  ts - the TS context
1275: .  step - current time-step
1276: .  ptime - current time
1277: .  u - current solution
1278: -  dctx - the TSMonitorSPCtx object that contains all the options for the monitoring, this is created with TSMonitorHGCtxCreate()

1280:    Options Database:
1281: + -ts_monitor_hg_swarm <n>             - Monitor the solution every n steps, or -1 for plotting only the final solution
1282: . -ts_monitor_hg_swarm_species <num>   - Number of species to histogram
1283: . -ts_monitor_hg_swarm_bins <num>      - Number of histogram bins
1284: - -ts_monitor_hg_swarm_velocity <bool> - Plot in velocity space, as opposed to coordinate space

1286:    Level: intermediate

1288:    Notes:
1289:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1290:    to be used during the `TS` integration.

1292: .seealso: `TSMonitoSet()`
1293: @*/
1294: PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1295: {
1296:   TSMonitorHGCtx     ctx = (TSMonitorHGCtx)dctx;
1297:   PetscDraw          draw;
1298:   DM                 sw;
1299:   const PetscScalar *yy;
1300:   PetscInt          *species;
1301:   PetscInt           dim, d = 0, Np, p, Ns, s;

1303:   PetscFunctionBegin;
1304:   if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1305:   PetscCall(TSGetDM(ts, &sw));
1306:   PetscCall(DMGetDimension(sw, &dim));
1307:   PetscCall(DMSwarmGetNumSpecies(sw, &Ns));
1308:   Ns = PetscMin(Ns, ctx->Ns);
1309:   PetscCall(VecGetLocalSize(u, &Np));
1310:   Np /= dim * 2;
1311:   if (!step) {
1312:     PetscDrawAxis axis;
1313:     char          title[PETSC_MAX_PATH_LEN];

1315:     for (s = 0; s < Ns; ++s) {
1316:       PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis));
1317:       PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s));
1318:       if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N"));
1319:       else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N"));
1320:     }
1321:   }
1322:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1323:     PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species));
1324:     for (s = 0; s < Ns; ++s) {
1325:       PetscCall(PetscDrawHGReset(ctx->hg[s]));
1326:       PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw));
1327:       PetscCall(PetscDrawClear(draw));
1328:       PetscCall(PetscDrawFlush(draw));
1329:     }
1330:     PetscCall(VecGetArrayRead(u, &yy));
1331:     for (p = 0; p < Np; ++p) {
1332:       const PetscInt s = species[p] < Ns ? species[p] : 0;
1333:       PetscReal      v;

1335:       if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]);
1336:       else v = PetscRealPart(yy[p * dim * 2 + d]);
1337:       PetscCall(PetscDrawHGAddValue(ctx->hg[s], v));
1338:     }
1339:     PetscCall(VecRestoreArrayRead(u, &yy));
1340:     for (s = 0; s < Ns; ++s) {
1341:       PetscCall(PetscDrawHGDraw(ctx->hg[s]));
1342:       PetscCall(PetscDrawHGSave(ctx->hg[s]));
1343:     }
1344:     PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species));
1345:   }
1346:   PetscFunctionReturn(PETSC_SUCCESS);
1347: }

1349: /*@C
1350:    TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep

1352:    Collective

1354:    Input Parameters:
1355: +  ts - the `TS` context
1356: .  step - current time-step
1357: .  ptime - current time
1358: .  u - current solution
1359: -  dctx - unused context

1361:    Options Database Key:
1362: .  -ts_monitor_error - create a graphical monitor of error history

1364:    Level: intermediate

1366:    Notes:
1367:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1368:    to be used during the `TS` integration.

1370:    The user must provide the solution using `TSSetSolutionFunction()` to use this monitor.

1372: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1373: @*/
1374: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1375: {
1376:   DM        dm;
1377:   PetscDS   ds = NULL;
1378:   PetscInt  Nf = -1, f;
1379:   PetscBool flg;

1381:   PetscFunctionBegin;
1382:   PetscCall(TSGetDM(ts, &dm));
1383:   if (dm) PetscCall(DMGetDS(dm, &ds));
1384:   if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf));
1385:   if (Nf <= 0) {
1386:     Vec       y;
1387:     PetscReal nrm;

1389:     PetscCall(VecDuplicate(u, &y));
1390:     PetscCall(TSComputeSolutionFunction(ts, ptime, y));
1391:     PetscCall(VecAXPY(y, -1.0, u));
1392:     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg));
1393:     if (flg) {
1394:       PetscCall(VecNorm(y, NORM_2, &nrm));
1395:       PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm));
1396:     }
1397:     PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg));
1398:     if (flg) PetscCall(VecView(y, vf->viewer));
1399:     PetscCall(VecDestroy(&y));
1400:   } else {
1401:     PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
1402:     void    **ctxs;
1403:     Vec       v;
1404:     PetscReal ferrors[1];

1406:     PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs));
1407:     for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]));
1408:     PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors));
1409:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04d time = %-8.4g \t L_2 Error: [", (int)step, (double)ptime));
1410:     for (f = 0; f < Nf; ++f) {
1411:       if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", "));
1412:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]));
1413:     }
1414:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n"));

1416:     PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));

1418:     PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg));
1419:     if (flg) {
1420:       PetscCall(DMGetGlobalVector(dm, &v));
1421:       PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v));
1422:       PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution"));
1423:       PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view"));
1424:       PetscCall(DMRestoreGlobalVector(dm, &v));
1425:     }
1426:     PetscCall(PetscFree2(exactFuncs, ctxs));
1427:   }
1428:   PetscFunctionReturn(PETSC_SUCCESS);
1429: }

1431: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1432: {
1433:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1434:   PetscReal      x   = ptime, y;
1435:   PetscInt       its;

1437:   PetscFunctionBegin;
1438:   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1439:   if (!n) {
1440:     PetscDrawAxis axis;
1441:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1442:     PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations"));
1443:     PetscCall(PetscDrawLGReset(ctx->lg));
1444:     ctx->snes_its = 0;
1445:   }
1446:   PetscCall(TSGetSNESIterations(ts, &its));
1447:   y = its - ctx->snes_its;
1448:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1449:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1450:     PetscCall(PetscDrawLGDraw(ctx->lg));
1451:     PetscCall(PetscDrawLGSave(ctx->lg));
1452:   }
1453:   ctx->snes_its = its;
1454:   PetscFunctionReturn(PETSC_SUCCESS);
1455: }

1457: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1458: {
1459:   TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1460:   PetscReal      x   = ptime, y;
1461:   PetscInt       its;

1463:   PetscFunctionBegin;
1464:   if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */
1465:   if (!n) {
1466:     PetscDrawAxis axis;
1467:     PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis));
1468:     PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations"));
1469:     PetscCall(PetscDrawLGReset(ctx->lg));
1470:     ctx->ksp_its = 0;
1471:   }
1472:   PetscCall(TSGetKSPIterations(ts, &its));
1473:   y = its - ctx->ksp_its;
1474:   PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y));
1475:   if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1476:     PetscCall(PetscDrawLGDraw(ctx->lg));
1477:     PetscCall(PetscDrawLGSave(ctx->lg));
1478:   }
1479:   ctx->ksp_its = its;
1480:   PetscFunctionReturn(PETSC_SUCCESS);
1481: }

1483: /*@C
1484:    TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()`

1486:    Collective

1488:    Input Parameters:
1489: .  ts  - the `TS` solver object

1491:    Output Parameter:
1492: .  ctx - the context

1494:    Level: intermediate

1496: .seealso: [](chapter_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1497: @*/
1498: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1499: {
1500:   PetscFunctionBegin;
1501:   PetscCall(PetscNew(ctx));
1502:   PetscFunctionReturn(PETSC_SUCCESS);
1503: }

1505: /*@C
1506:    TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution

1508:    Collective

1510:    Input Parameters:
1511: +  ts - the `TS` context
1512: .  step - current time-step
1513: .  ptime - current time
1514: .  u  - current solution
1515: -  dctx - the envelope context

1517:    Options Database Key:
1518: .  -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time

1520:    Level: intermediate

1522:    Notes:
1523:    After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope

1525:    This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1526:    to be used during the `TS` integration.

1528: .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1529: @*/
1530: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1531: {
1532:   TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;

1534:   PetscFunctionBegin;
1535:   if (!ctx->max) {
1536:     PetscCall(VecDuplicate(u, &ctx->max));
1537:     PetscCall(VecDuplicate(u, &ctx->min));
1538:     PetscCall(VecCopy(u, ctx->max));
1539:     PetscCall(VecCopy(u, ctx->min));
1540:   } else {
1541:     PetscCall(VecPointwiseMax(ctx->max, u, ctx->max));
1542:     PetscCall(VecPointwiseMin(ctx->min, u, ctx->min));
1543:   }
1544:   PetscFunctionReturn(PETSC_SUCCESS);
1545: }

1547: /*@C
1548:    TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution

1550:    Collective

1552:    Input Parameter:
1553: .  ts - the `TS` context

1555:    Output Parameters:
1556: +  max - the maximum values
1557: -  min - the minimum values

1559:    Level: intermediate

1561:    Notes:
1562:     If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored

1564: .seealso: [](chapter_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1565: @*/
1566: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1567: {
1568:   PetscInt i;

1570:   PetscFunctionBegin;
1571:   if (max) *max = NULL;
1572:   if (min) *min = NULL;
1573:   for (i = 0; i < ts->numbermonitors; i++) {
1574:     if (ts->monitor[i] == TSMonitorEnvelope) {
1575:       TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1576:       if (max) *max = ctx->max;
1577:       if (min) *min = ctx->min;
1578:       break;
1579:     }
1580:   }
1581:   PetscFunctionReturn(PETSC_SUCCESS);
1582: }

1584: /*@C
1585:    TSMonitorEnvelopeCtxDestroy - Destroys a context that was created  with `TSMonitorEnvelopeCtxCreate()`.

1587:    Collective

1589:    Input Parameter:
1590: .  ctx - the monitor context

1592:    Level: intermediate

1594: .seealso: [](chapter_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1595: @*/
1596: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1597: {
1598:   PetscFunctionBegin;
1599:   PetscCall(VecDestroy(&(*ctx)->min));
1600:   PetscCall(VecDestroy(&(*ctx)->max));
1601:   PetscCall(PetscFree(*ctx));
1602:   PetscFunctionReturn(PETSC_SUCCESS);
1603: }

1605: /*@C
1606:   TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS`

1608:   Not collective

1610:   Input Parameters:
1611: + ts   - the `TS` context
1612: . step - current timestep
1613: . t    - current time
1614: . u    - current solution
1615: - ctx  - not used

1617:   Options Database Key:
1618: . -ts_dmswarm_monitor_moments - Monitor moments of particle distribution

1620:   Level: intermediate

1622:   Notes:
1623:   This requires a `DMSWARM` be attached to the `TS`.

1625:   This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1626:   to be used during the TS integration.

1628: .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1629: @*/
1630: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1631: {
1632:   DM                 sw;
1633:   const PetscScalar *u;
1634:   PetscReal          m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1635:   PetscInt           dim, d, Np, p;
1636:   MPI_Comm           comm;

1638:   PetscFunctionBeginUser;
1639:   PetscCall(TSGetDM(ts, &sw));
1640:   if (!sw || step % ts->monitorFrequency != 0) PetscFunctionReturn(PETSC_SUCCESS);
1641:   PetscCall(PetscObjectGetComm((PetscObject)ts, &comm));
1642:   PetscCall(DMGetDimension(sw, &dim));
1643:   PetscCall(VecGetLocalSize(U, &Np));
1644:   Np /= dim;
1645:   PetscCall(VecGetArrayRead(U, &u));
1646:   for (p = 0; p < Np; ++p) {
1647:     for (d = 0; d < dim; ++d) {
1648:       totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1649:       totMom[d] += PetscRealPart(u[p * dim + d]);
1650:     }
1651:   }
1652:   PetscCall(VecRestoreArrayRead(U, &u));
1653:   for (d = 0; d < dim; ++d) totMom[d] *= m;
1654:   totE *= 0.5 * m;
1655:   PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE));
1656:   for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, "    Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]));
1657:   PetscCall(PetscPrintf(comm, "\n"));
1658:   PetscFunctionReturn(PETSC_SUCCESS);
1659: }