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 on TS
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: Notes:
19: TSMonitor() is typically used automatically within the time stepping implementations.
20: Users would almost never call this routine directly.
22: A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions
24: Level: developer
26: @*/
27: PetscErrorCode TSMonitor(TS ts, PetscInt step, PetscReal ptime, Vec u)
28: {
29: DM dm;
30: PetscInt i, n = ts->numbermonitors;
35: TSGetDM(ts, &dm);
36: DMSetOutputSequenceNumber(dm, step, ptime);
38: VecLockReadPush(u);
39: for (i = 0; i < n; i++) (*ts->monitor[i])(ts, step, ptime, u, ts->monitorcontext[i]);
40: VecLockReadPop(u);
41: return 0;
42: }
44: /*@C
45: TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
47: Collective on TS
49: Input Parameters:
50: + ts - TS object you wish to monitor
51: . name - the monitor type one is seeking
52: . help - message indicating what monitoring is done
53: . manual - manual page for the monitor
54: . monitor - the monitor function
55: - 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
57: Level: developer
59: .seealso: `PetscOptionsGetViewer()`, `PetscOptionsGetReal()`, `PetscOptionsHasName()`, `PetscOptionsGetString()`,
60: `PetscOptionsGetIntArray()`, `PetscOptionsGetRealArray()`, `PetscOptionsBool()`
61: `PetscOptionsInt()`, `PetscOptionsString()`, `PetscOptionsReal()`, `PetscOptionsBool()`,
62: `PetscOptionsName()`, `PetscOptionsBegin()`, `PetscOptionsEnd()`, `PetscOptionsHeadBegin()`,
63: `PetscOptionsStringArray()`, `PetscOptionsRealArray()`, `PetscOptionsScalar()`,
64: `PetscOptionsBoolGroupBegin()`, `PetscOptionsBoolGroup()`, `PetscOptionsBoolGroupEnd()`,
65: `PetscOptionsFList()`, `PetscOptionsEList()`
66: @*/
67: PetscErrorCode TSMonitorSetFromOptions(TS ts, const char name[], const char help[], const char manual[], PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, PetscViewerAndFormat *), PetscErrorCode (*monitorsetup)(TS, PetscViewerAndFormat *))
68: {
69: PetscViewer viewer;
70: PetscViewerFormat format;
71: PetscBool flg;
73: PetscOptionsGetViewer(PetscObjectComm((PetscObject)ts), ((PetscObject)ts)->options, ((PetscObject)ts)->prefix, name, &viewer, &format, &flg);
74: if (flg) {
75: PetscViewerAndFormat *vf;
76: PetscViewerAndFormatCreate(viewer, format, &vf);
77: PetscObjectDereference((PetscObject)viewer);
78: if (monitorsetup) (*monitorsetup)(ts, vf);
79: TSMonitorSet(ts, (PetscErrorCode(*)(TS, PetscInt, PetscReal, Vec, void *))monitor, vf, (PetscErrorCode(*)(void **))PetscViewerAndFormatDestroy);
80: }
81: return 0;
82: }
84: /*@C
85: TSMonitorSet - Sets an ADDITIONAL function that is to be used at every
86: timestep to display the iteration's progress.
88: Logically Collective on TS
90: Input Parameters:
91: + ts - the TS context obtained from TSCreate()
92: . monitor - monitoring routine
93: . mctx - [optional] user-defined context for private data for the
94: monitor routine (use NULL if no context is desired)
95: - monitordestroy - [optional] routine that frees monitor context
96: (may be NULL)
98: Calling sequence of monitor:
99: $ PetscErrorCode monitor(TS ts,PetscInt steps,PetscReal time,Vec u,void *mctx)
101: + ts - the TS context
102: . 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)
103: . time - current time
104: . u - current iterate
105: - mctx - [optional] monitoring context
107: Notes:
108: This routine adds an additional monitor to the list of monitors that
109: already has been loaded.
111: Fortran Notes:
112: Only a single monitor function can be set for each TS object
114: Level: intermediate
116: .seealso: `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
117: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
118: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
119: @*/
120: PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, void *), void *mctx, PetscErrorCode (*mdestroy)(void **))
121: {
122: PetscInt i;
123: PetscBool identical;
126: for (i = 0; i < ts->numbermonitors; i++) {
127: PetscMonitorCompare((PetscErrorCode(*)(void))monitor, mctx, mdestroy, (PetscErrorCode(*)(void))ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical);
128: if (identical) return 0;
129: }
131: ts->monitor[ts->numbermonitors] = monitor;
132: ts->monitordestroy[ts->numbermonitors] = mdestroy;
133: ts->monitorcontext[ts->numbermonitors++] = (void *)mctx;
134: return 0;
135: }
137: /*@C
138: TSMonitorCancel - Clears all the monitors that have been set on a time-step object.
140: Logically Collective on TS
142: Input Parameters:
143: . ts - the TS context obtained from TSCreate()
145: Notes:
146: There is no way to remove a single, specific monitor.
148: Level: intermediate
150: .seealso: `TSMonitorDefault()`, `TSMonitorSet()`
151: @*/
152: PetscErrorCode TSMonitorCancel(TS ts)
153: {
154: PetscInt i;
157: for (i = 0; i < ts->numbermonitors; i++) {
158: if (ts->monitordestroy[i]) (*ts->monitordestroy[i])(&ts->monitorcontext[i]);
159: }
160: ts->numbermonitors = 0;
161: return 0;
162: }
164: /*@C
165: TSMonitorDefault - The Default monitor, prints the timestep and time for each step
167: Options Database:
168: . -ts_monitor - monitors the time integration
170: Notes:
171: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
172: to be used during the TS integration.
174: Level: intermediate
176: .seealso: `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`,
177: `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`,
178: `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()`
179: @*/
180: PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
181: {
182: PetscViewer viewer = vf->viewer;
183: PetscBool iascii, ibinary;
186: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
187: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary);
188: PetscViewerPushFormat(viewer, vf->format);
189: if (iascii) {
190: PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel);
191: if (step == -1) { /* this indicates it is an interpolated solution */
192: PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps);
193: } else {
194: PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n");
195: }
196: PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel);
197: } else if (ibinary) {
198: PetscMPIInt rank;
199: MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank);
200: if (rank == 0) {
201: PetscBool skipHeader;
202: PetscInt classid = REAL_FILE_CLASSID;
204: PetscViewerBinaryGetSkipHeader(viewer, &skipHeader);
205: if (!skipHeader) PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT);
206: PetscRealView(1, &ptime, viewer);
207: } else {
208: PetscRealView(0, &ptime, viewer);
209: }
210: }
211: PetscViewerPopFormat(viewer);
212: return 0;
213: }
215: /*@C
216: TSMonitorExtreme - Prints the extreme values of the solution at each timestep
218: Notes:
219: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
220: to be used during the TS integration.
222: Level: intermediate
224: .seealso: `TSMonitorSet()`
225: @*/
226: PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf)
227: {
228: PetscViewer viewer = vf->viewer;
229: PetscBool iascii;
230: PetscReal max, min;
233: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
234: PetscViewerPushFormat(viewer, vf->format);
235: if (iascii) {
236: VecMax(v, NULL, &max);
237: VecMin(v, NULL, &min);
238: PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel);
239: 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);
240: PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel);
241: }
242: PetscViewerPopFormat(viewer);
243: return 0;
244: }
246: /*@C
247: TSMonitorLGCtxCreate - Creates a TSMonitorLGCtx context for use with
248: TS to monitor the solution process graphically in various ways
250: Collective on TS
252: Input Parameters:
253: + host - the X display to open, or null for the local machine
254: . label - the title to put in the title bar
255: . x, y - the screen coordinates of the upper left coordinate of the window
256: . m, n - the screen width and height in pixels
257: - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time
259: Output Parameter:
260: . ctx - the context
262: Options Database Key:
263: + -ts_monitor_lg_timestep - automatically sets line graph monitor
264: + -ts_monitor_lg_timestep_log - automatically sets line graph monitor
265: . -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling TSMonitorLGSetDisplayVariables() or TSMonitorLGCtxSetDisplayVariables())
266: . -ts_monitor_lg_error - monitor the error
267: . -ts_monitor_lg_ksp_iterations - monitor the number of KSP iterations needed for each timestep
268: . -ts_monitor_lg_snes_iterations - monitor the number of SNES iterations needed for each timestep
269: - -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true
271: Notes:
272: Use TSMonitorLGCtxDestroy() to destroy.
274: One can provide a function that transforms the solution before plotting it with TSMonitorLGCtxSetTransform() or TSMonitorLGSetTransform()
276: Many of the functions that control the monitoring have two forms: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a TS object as the
277: 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
278: as the first argument.
280: One can control the names displayed for each solution or error variable with TSMonitorLGCtxSetVariableNames() or TSMonitorLGSetVariableNames()
282: Level: intermediate
284: .seealso: `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`,
285: `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
286: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
287: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
288: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
290: @*/
291: PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx)
292: {
293: PetscDraw draw;
295: PetscNew(ctx);
296: PetscDrawCreate(comm, host, label, x, y, m, n, &draw);
297: PetscDrawSetFromOptions(draw);
298: PetscDrawLGCreate(draw, 1, &(*ctx)->lg);
299: PetscDrawLGSetFromOptions((*ctx)->lg);
300: PetscDrawDestroy(&draw);
301: (*ctx)->howoften = howoften;
302: return 0;
303: }
305: PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx)
306: {
307: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
308: PetscReal x = ptime, y;
310: if (step < 0) return 0; /* -1 indicates an interpolated solution */
311: if (!step) {
312: PetscDrawAxis axis;
313: const char *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step";
314: PetscDrawLGGetAxis(ctx->lg, &axis);
315: PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel);
316: PetscDrawLGReset(ctx->lg);
317: }
318: TSGetTimeStep(ts, &y);
319: if (ctx->semilogy) y = PetscLog10Real(y);
320: PetscDrawLGAddPoint(ctx->lg, &x, &y);
321: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
322: PetscDrawLGDraw(ctx->lg);
323: PetscDrawLGSave(ctx->lg);
324: }
325: return 0;
326: }
328: /*@C
329: TSMonitorLGCtxDestroy - Destroys a line graph context that was created
330: with TSMonitorLGCtxCreate().
332: Collective on TSMonitorLGCtx
334: Input Parameter:
335: . ctx - the monitor context
337: Level: intermediate
339: .seealso: `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep();`
340: @*/
341: PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx)
342: {
343: if ((*ctx)->transformdestroy) ((*ctx)->transformdestroy)((*ctx)->transformctx);
344: PetscDrawLGDestroy(&(*ctx)->lg);
345: PetscStrArrayDestroy(&(*ctx)->names);
346: PetscStrArrayDestroy(&(*ctx)->displaynames);
347: PetscFree((*ctx)->displayvariables);
348: PetscFree((*ctx)->displayvalues);
349: PetscFree(*ctx);
350: return 0;
351: }
353: /* Creates a TS Monitor SPCtx for use with DMSwarm particle visualizations */
354: 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, TSMonitorSPCtx *ctx)
355: {
356: PetscDraw draw;
358: PetscNew(ctx);
359: PetscDrawCreate(comm, host, label, x, y, m, n, &draw);
360: PetscDrawSetFromOptions(draw);
361: PetscDrawSPCreate(draw, 1, &(*ctx)->sp);
362: PetscDrawDestroy(&draw);
363: (*ctx)->howoften = howoften;
364: (*ctx)->retain = retain;
365: (*ctx)->phase = phase;
366: return 0;
367: }
369: /*
370: Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate
371: */
372: PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx)
373: {
375: PetscDrawSPDestroy(&(*ctx)->sp);
376: PetscFree(*ctx);
378: return 0;
379: }
381: /*@C
382: TSMonitorDrawSolution - Monitors progress of the TS solvers by calling
383: VecView() for the solution at each timestep
385: Collective on TS
387: Input Parameters:
388: + ts - the TS context
389: . step - current time-step
390: . ptime - current time
391: - dummy - either a viewer or NULL
393: Options Database:
394: . -ts_monitor_draw_solution_initial - show initial solution as well as current solution
396: Notes:
397: The initial solution and current solution are not displayed with a common axis scaling so generally the option -ts_monitor_draw_solution_initial
398: will look bad
400: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
401: to be used during the TS integration.
403: Level: intermediate
405: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
406: @*/
407: PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
408: {
409: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
410: PetscDraw draw;
412: if (!step && ictx->showinitial) {
413: if (!ictx->initialsolution) VecDuplicate(u, &ictx->initialsolution);
414: VecCopy(u, ictx->initialsolution);
415: }
416: if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) return 0;
418: if (ictx->showinitial) {
419: PetscReal pause;
420: PetscViewerDrawGetPause(ictx->viewer, &pause);
421: PetscViewerDrawSetPause(ictx->viewer, 0.0);
422: VecView(ictx->initialsolution, ictx->viewer);
423: PetscViewerDrawSetPause(ictx->viewer, pause);
424: PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE);
425: }
426: VecView(u, ictx->viewer);
427: if (ictx->showtimestepandtime) {
428: PetscReal xl, yl, xr, yr, h;
429: char time[32];
431: PetscViewerDrawGetDraw(ictx->viewer, 0, &draw);
432: PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime);
433: PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr);
434: h = yl + .95 * (yr - yl);
435: PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time);
436: PetscDrawFlush(draw);
437: }
439: if (ictx->showinitial) PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE);
440: return 0;
441: }
443: /*@C
444: TSMonitorDrawSolutionPhase - Monitors progress of the TS solvers by plotting the solution as a phase diagram
446: Collective on TS
448: Input Parameters:
449: + ts - the TS context
450: . step - current time-step
451: . ptime - current time
452: - dummy - either a viewer or NULL
454: Notes:
455: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
456: to be used during the TS integration.
458: Level: intermediate
460: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
461: @*/
462: PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
463: {
464: TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy;
465: PetscDraw draw;
466: PetscDrawAxis axis;
467: PetscInt n;
468: PetscMPIInt size;
469: PetscReal U0, U1, xl, yl, xr, yr, h;
470: char time[32];
471: const PetscScalar *U;
473: MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size);
475: VecGetSize(u, &n);
478: PetscViewerDrawGetDraw(ictx->viewer, 0, &draw);
479: PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis);
480: PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr);
481: if (!step) {
482: PetscDrawClear(draw);
483: PetscDrawAxisDraw(axis);
484: }
486: VecGetArrayRead(u, &U);
487: U0 = PetscRealPart(U[0]);
488: U1 = PetscRealPart(U[1]);
489: VecRestoreArrayRead(u, &U);
490: if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) return 0;
492: PetscDrawCollectiveBegin(draw);
493: PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK);
494: if (ictx->showtimestepandtime) {
495: PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr);
496: PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime);
497: h = yl + .95 * (yr - yl);
498: PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time);
499: }
500: PetscDrawCollectiveEnd(draw);
501: PetscDrawFlush(draw);
502: PetscDrawPause(draw);
503: PetscDrawSave(draw);
504: return 0;
505: }
507: /*@C
508: TSMonitorDrawCtxDestroy - Destroys the monitor context for TSMonitorDrawSolution()
510: Collective on TS
512: Input Parameters:
513: . ctx - the monitor context
515: Level: intermediate
517: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`
518: @*/
519: PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx)
520: {
521: PetscViewerDestroy(&(*ictx)->viewer);
522: VecDestroy(&(*ictx)->initialsolution);
523: PetscFree(*ictx);
524: return 0;
525: }
527: /*@C
528: TSMonitorDrawCtxCreate - Creates the monitor context for TSMonitorDrawCtx
530: Collective on TS
532: Input Parameter:
533: . ts - time-step context
535: Output Patameter:
536: . ctx - the monitor context
538: Options Database:
539: . -ts_monitor_draw_solution_initial - show initial solution as well as current solution
541: Level: intermediate
543: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx()`
544: @*/
545: PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx)
546: {
547: PetscNew(ctx);
548: PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer);
549: PetscViewerSetFromOptions((*ctx)->viewer);
551: (*ctx)->howoften = howoften;
552: (*ctx)->showinitial = PETSC_FALSE;
553: PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL);
555: (*ctx)->showtimestepandtime = PETSC_FALSE;
556: PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL);
557: return 0;
558: }
560: /*@C
561: TSMonitorDrawSolutionFunction - Monitors progress of the TS solvers by calling
562: VecView() for the solution provided by TSSetSolutionFunction() at each timestep
564: Collective on TS
566: Input Parameters:
567: + ts - the TS context
568: . step - current time-step
569: . ptime - current time
570: - dummy - either a viewer or NULL
572: Options Database:
573: . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided TSSetSolutionFunction()
575: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
576: to be used during the TS integration.
578: Level: intermediate
580: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
581: @*/
582: PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
583: {
584: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
585: PetscViewer viewer = ctx->viewer;
586: Vec work;
588: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) return 0;
589: VecDuplicate(u, &work);
590: TSComputeSolutionFunction(ts, ptime, work);
591: VecView(work, viewer);
592: VecDestroy(&work);
593: return 0;
594: }
596: /*@C
597: TSMonitorDrawError - Monitors progress of the TS solvers by calling
598: VecView() for the error at each timestep
600: Collective on TS
602: Input Parameters:
603: + ts - the TS context
604: . step - current time-step
605: . ptime - current time
606: - dummy - either a viewer or NULL
608: Options Database:
609: . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided TSSetSolutionFunction()
611: Notes:
612: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
613: to be used during the TS integration.
615: Level: intermediate
617: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
618: @*/
619: PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
620: {
621: TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy;
622: PetscViewer viewer = ctx->viewer;
623: Vec work;
625: if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) return 0;
626: VecDuplicate(u, &work);
627: TSComputeSolutionFunction(ts, ptime, work);
628: VecAXPY(work, -1.0, u);
629: VecView(work, viewer);
630: VecDestroy(&work);
631: return 0;
632: }
634: /*@C
635: 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
637: Collective on TS
639: Input Parameters:
640: + ts - the TS context
641: . step - current time-step
642: . ptime - current time
643: . u - current state
644: - vf - viewer and its format
646: Notes:
647: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
648: to be used during the TS integration.
650: Level: intermediate
652: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
653: @*/
654: PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
655: {
656: PetscViewerPushFormat(vf->viewer, vf->format);
657: VecView(u, vf->viewer);
658: PetscViewerPopFormat(vf->viewer);
659: return 0;
660: }
662: /*@C
663: TSMonitorSolutionVTK - Monitors progress of the TS solvers by VecView() for the solution at each timestep.
665: Collective on TS
667: Input Parameters:
668: + ts - the TS context
669: . step - current time-step
670: . ptime - current time
671: . u - current state
672: - filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")
674: Level: intermediate
676: Notes:
677: 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.
678: These are named according to the file name template.
680: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
681: to be used during the TS integration.
683: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`
684: @*/
685: PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, void *filenametemplate)
686: {
687: char filename[PETSC_MAX_PATH_LEN];
688: PetscViewer viewer;
690: if (step < 0) return 0; /* -1 indicates interpolated solution */
691: PetscSNPrintf(filename, sizeof(filename), (const char *)filenametemplate, step);
692: PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer);
693: VecView(u, viewer);
694: PetscViewerDestroy(&viewer);
695: return 0;
696: }
698: /*@C
699: TSMonitorSolutionVTKDestroy - Destroy context for monitoring
701: Collective on TS
703: Input Parameters:
704: . filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ")
706: Level: intermediate
708: Note:
709: This function is normally passed to TSMonitorSet() along with TSMonitorSolutionVTK().
711: .seealso: `TSMonitorSet()`, `TSMonitorSolutionVTK()`
712: @*/
713: PetscErrorCode TSMonitorSolutionVTKDestroy(void *filenametemplate)
714: {
715: PetscFree(*(char **)filenametemplate);
716: return 0;
717: }
719: /*@C
720: TSMonitorLGSolution - Monitors progress of the TS solvers by plotting each component of the solution vector
721: in a time based line graph
723: Collective on TS
725: Input Parameters:
726: + ts - the TS context
727: . step - current time-step
728: . ptime - current time
729: . u - current solution
730: - dctx - the TSMonitorLGCtx object that contains all the options for the monitoring, this is created with TSMonitorLGCtxCreate()
732: Options Database:
733: . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables
735: Level: intermediate
737: Notes:
738: Each process in a parallel run displays its component solutions in a separate window
740: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
741: to be used during the TS integration.
743: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`,
744: `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`,
745: `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`,
746: `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()`
747: @*/
748: PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
749: {
750: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dctx;
751: const PetscScalar *yy;
752: Vec v;
754: if (step < 0) return 0; /* -1 indicates interpolated solution */
755: if (!step) {
756: PetscDrawAxis axis;
757: PetscInt dim;
758: PetscDrawLGGetAxis(ctx->lg, &axis);
759: PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution");
760: if (!ctx->names) {
761: PetscBool flg;
762: /* user provides names of variables to plot but no names has been set so assume names are integer values */
763: PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg);
764: if (flg) {
765: PetscInt i, n;
766: char **names;
767: VecGetSize(u, &n);
768: PetscMalloc1(n + 1, &names);
769: for (i = 0; i < n; i++) {
770: PetscMalloc1(5, &names[i]);
771: PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i);
772: }
773: names[n] = NULL;
774: ctx->names = names;
775: }
776: }
777: if (ctx->names && !ctx->displaynames) {
778: char **displaynames;
779: PetscBool flg;
780: VecGetLocalSize(u, &dim);
781: PetscCalloc1(dim + 1, &displaynames);
782: PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg);
783: if (flg) TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames);
784: PetscStrArrayDestroy(&displaynames);
785: }
786: if (ctx->displaynames) {
787: PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables);
788: PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames);
789: } else if (ctx->names) {
790: VecGetLocalSize(u, &dim);
791: PetscDrawLGSetDimension(ctx->lg, dim);
792: PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names);
793: } else {
794: VecGetLocalSize(u, &dim);
795: PetscDrawLGSetDimension(ctx->lg, dim);
796: }
797: PetscDrawLGReset(ctx->lg);
798: }
800: if (!ctx->transform) v = u;
801: else (*ctx->transform)(ctx->transformctx, u, &v);
802: VecGetArrayRead(v, &yy);
803: if (ctx->displaynames) {
804: PetscInt i;
805: for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]);
806: PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues);
807: } else {
808: #if defined(PETSC_USE_COMPLEX)
809: PetscInt i, n;
810: PetscReal *yreal;
811: VecGetLocalSize(v, &n);
812: PetscMalloc1(n, &yreal);
813: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
814: PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal);
815: PetscFree(yreal);
816: #else
817: PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy);
818: #endif
819: }
820: VecRestoreArrayRead(v, &yy);
821: if (ctx->transform) VecDestroy(&v);
823: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
824: PetscDrawLGDraw(ctx->lg);
825: PetscDrawLGSave(ctx->lg);
826: }
827: return 0;
828: }
830: /*@C
831: TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
833: Collective on TS
835: Input Parameters:
836: + ts - the TS context
837: - names - the names of the components, final string must be NULL
839: Level: intermediate
841: Notes:
842: If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored
844: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()`
845: @*/
846: PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names)
847: {
848: PetscInt i;
850: for (i = 0; i < ts->numbermonitors; i++) {
851: if (ts->monitor[i] == TSMonitorLGSolution) {
852: TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names);
853: break;
854: }
855: }
856: return 0;
857: }
859: /*@C
860: TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot
862: Collective on TS
864: Input Parameters:
865: + ts - the TS context
866: - names - the names of the components, final string must be NULL
868: Level: intermediate
870: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()`
871: @*/
872: PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names)
873: {
874: PetscStrArrayDestroy(&ctx->names);
875: PetscStrArrayallocpy(names, &ctx->names);
876: return 0;
877: }
879: /*@C
880: TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot
882: Collective on TS
884: Input Parameter:
885: . ts - the TS context
887: Output Parameter:
888: . names - the names of the components, final string must be NULL
890: Level: intermediate
892: Notes:
893: If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored
895: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
896: @*/
897: PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names)
898: {
899: PetscInt i;
901: *names = NULL;
902: for (i = 0; i < ts->numbermonitors; i++) {
903: if (ts->monitor[i] == TSMonitorLGSolution) {
904: TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i];
905: *names = (const char *const *)ctx->names;
906: break;
907: }
908: }
909: return 0;
910: }
912: /*@C
913: TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor
915: Collective on TS
917: Input Parameters:
918: + ctx - the TSMonitorLG context
919: - displaynames - the names of the components, final string must be NULL
921: Level: intermediate
923: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
924: @*/
925: PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames)
926: {
927: PetscInt j = 0, k;
929: if (!ctx->names) return 0;
930: PetscStrArrayDestroy(&ctx->displaynames);
931: PetscStrArrayallocpy(displaynames, &ctx->displaynames);
932: while (displaynames[j]) j++;
933: ctx->ndisplayvariables = j;
934: PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables);
935: PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues);
936: j = 0;
937: while (displaynames[j]) {
938: k = 0;
939: while (ctx->names[k]) {
940: PetscBool flg;
941: PetscStrcmp(displaynames[j], ctx->names[k], &flg);
942: if (flg) {
943: ctx->displayvariables[j] = k;
944: break;
945: }
946: k++;
947: }
948: j++;
949: }
950: return 0;
951: }
953: /*@C
954: TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor
956: Collective on TS
958: Input Parameters:
959: + ts - the TS context
960: - displaynames - the names of the components, final string must be NULL
962: Notes:
963: If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored
965: Level: intermediate
967: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`
968: @*/
969: PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames)
970: {
971: PetscInt i;
973: for (i = 0; i < ts->numbermonitors; i++) {
974: if (ts->monitor[i] == TSMonitorLGSolution) {
975: TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames);
976: break;
977: }
978: }
979: return 0;
980: }
982: /*@C
983: TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed
985: Collective on TS
987: Input Parameters:
988: + ts - the TS context
989: . transform - the transform function
990: . destroy - function to destroy the optional context
991: - ctx - optional context used by transform function
993: Notes:
994: If the TS object does not have a TSMonitorLGCtx associated with it then this function is ignored
996: Level: intermediate
998: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()`
999: @*/
1000: PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1001: {
1002: PetscInt i;
1004: for (i = 0; i < ts->numbermonitors; i++) {
1005: if (ts->monitor[i] == TSMonitorLGSolution) TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx);
1006: }
1007: return 0;
1008: }
1010: /*@C
1011: TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed
1013: Collective on TSLGCtx
1015: Input Parameters:
1016: + ts - the TS context
1017: . transform - the transform function
1018: . destroy - function to destroy the optional context
1019: - ctx - optional context used by transform function
1021: Level: intermediate
1023: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()`
1024: @*/
1025: PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx)
1026: {
1027: ctx->transform = transform;
1028: ctx->transformdestroy = destroy;
1029: ctx->transformctx = tctx;
1030: return 0;
1031: }
1033: /*@C
1034: TSMonitorLGError - Monitors progress of the TS solvers by plotting each component of the error
1035: in a time based line graph
1037: Collective on TS
1039: Input Parameters:
1040: + ts - the TS context
1041: . step - current time-step
1042: . ptime - current time
1043: . u - current solution
1044: - dctx - TSMonitorLGCtx object created with TSMonitorLGCtxCreate()
1046: Options Database Keys:
1047: . -ts_monitor_lg_error - create a graphical monitor of error history
1049: Level: intermediate
1051: Notes:
1052: Each process in a parallel run displays its component errors in a separate window
1054: The user must provide the solution using TSSetSolutionFunction() to use this monitor.
1056: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1057: to be used during the TS integration.
1059: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1060: @*/
1061: PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy)
1062: {
1063: TSMonitorLGCtx ctx = (TSMonitorLGCtx)dummy;
1064: const PetscScalar *yy;
1065: Vec y;
1067: if (!step) {
1068: PetscDrawAxis axis;
1069: PetscInt dim;
1070: PetscDrawLGGetAxis(ctx->lg, &axis);
1071: PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error");
1072: VecGetLocalSize(u, &dim);
1073: PetscDrawLGSetDimension(ctx->lg, dim);
1074: PetscDrawLGReset(ctx->lg);
1075: }
1076: VecDuplicate(u, &y);
1077: TSComputeSolutionFunction(ts, ptime, y);
1078: VecAXPY(y, -1.0, u);
1079: VecGetArrayRead(y, &yy);
1080: #if defined(PETSC_USE_COMPLEX)
1081: {
1082: PetscReal *yreal;
1083: PetscInt i, n;
1084: VecGetLocalSize(y, &n);
1085: PetscMalloc1(n, &yreal);
1086: for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]);
1087: PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal);
1088: PetscFree(yreal);
1089: }
1090: #else
1091: PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy);
1092: #endif
1093: VecRestoreArrayRead(y, &yy);
1094: VecDestroy(&y);
1095: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1096: PetscDrawLGDraw(ctx->lg);
1097: PetscDrawLGSave(ctx->lg);
1098: }
1099: return 0;
1100: }
1102: /*@C
1103: TSMonitorSPSwarmSolution - Graphically displays phase plots of DMSwarm particles on a scatter plot
1105: Input Parameters:
1106: + ts - the TS context
1107: . step - current time-step
1108: . ptime - current time
1109: . u - current solution
1110: - dctx - the TSMonitorSPCtx object that contains all the options for the monitoring, this is created with TSMonitorSPCtxCreate()
1112: Options Database:
1113: + -ts_monitor_sp_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution
1114: . -ts_monitor_sp_swarm_retain <n> - Retain n old points so we can see the history, or -1 for all points
1115: - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space
1117: Level: intermediate
1119: Notes:
1120: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1121: to be used during the TS integration.
1123: .seealso: `TSMonitoSet()`
1124: @*/
1125: PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1126: {
1127: TSMonitorSPCtx ctx = (TSMonitorSPCtx)dctx;
1128: PetscDraw draw;
1129: DM dm, cdm;
1130: const PetscScalar *yy;
1131: PetscInt Np, p, dim = 2;
1133: if (step < 0) return 0; /* -1 indicates interpolated solution */
1134: if (!step) {
1135: PetscDrawAxis axis;
1136: PetscReal dmboxlower[2], dmboxupper[2];
1138: TSGetDM(ts, &dm);
1139: DMGetDimension(dm, &dim);
1141: DMSwarmGetCellDM(dm, &cdm);
1142: DMGetBoundingBox(cdm, dmboxlower, dmboxupper);
1143: VecGetLocalSize(u, &Np);
1144: Np /= dim * 2;
1145: PetscDrawSPGetAxis(ctx->sp, &axis);
1146: if (ctx->phase) {
1147: PetscDrawAxisSetLabels(axis, "Particles", "X", "V");
1148: PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -5, 5);
1149: } else {
1150: PetscDrawAxisSetLabels(axis, "Particles", "X", "Y");
1151: PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1]);
1152: }
1153: PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE);
1154: PetscDrawSPReset(ctx->sp);
1155: }
1156: VecGetLocalSize(u, &Np);
1157: Np /= dim * 2;
1158: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) {
1159: PetscDrawSPGetDraw(ctx->sp, &draw);
1160: if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscDrawClear(draw);
1161: PetscDrawFlush(draw);
1162: PetscDrawSPReset(ctx->sp);
1163: VecGetArrayRead(u, &yy);
1164: for (p = 0; p < Np; ++p) {
1165: PetscReal x, y;
1167: if (ctx->phase) {
1168: x = PetscRealPart(yy[p * dim * 2]);
1169: y = PetscRealPart(yy[p * dim * 2 + dim]);
1170: } else {
1171: x = PetscRealPart(yy[p * dim * 2]);
1172: y = PetscRealPart(yy[p * dim * 2 + 1]);
1173: }
1174: PetscDrawSPAddPoint(ctx->sp, &x, &y);
1175: }
1176: VecRestoreArrayRead(u, &yy);
1177: PetscDrawSPDraw(ctx->sp, PETSC_FALSE);
1178: PetscDrawSPSave(ctx->sp);
1179: }
1180: return 0;
1181: }
1183: /*@C
1184: TSMonitorError - Monitors progress of the TS solvers by printing the 2 norm of the error at each timestep
1186: Collective on TS
1188: Input Parameters:
1189: + ts - the TS context
1190: . step - current time-step
1191: . ptime - current time
1192: . u - current solution
1193: - dctx - unused context
1195: Level: intermediate
1197: Notes:
1198: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1199: to be used during the TS integration.
1201: The user must provide the solution using TSSetSolutionFunction() to use this monitor.
1203: Options Database Keys:
1204: . -ts_monitor_error - create a graphical monitor of error history
1206: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()`
1207: @*/
1208: PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf)
1209: {
1210: DM dm;
1211: PetscDS ds = NULL;
1212: PetscInt Nf = -1, f;
1213: PetscBool flg;
1215: TSGetDM(ts, &dm);
1216: if (dm) DMGetDS(dm, &ds);
1217: if (ds) PetscDSGetNumFields(ds, &Nf);
1218: if (Nf <= 0) {
1219: Vec y;
1220: PetscReal nrm;
1222: VecDuplicate(u, &y);
1223: TSComputeSolutionFunction(ts, ptime, y);
1224: VecAXPY(y, -1.0, u);
1225: PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg);
1226: if (flg) {
1227: VecNorm(y, NORM_2, &nrm);
1228: PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm);
1229: }
1230: PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg);
1231: if (flg) VecView(y, vf->viewer);
1232: VecDestroy(&y);
1233: } else {
1234: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
1235: void **ctxs;
1236: Vec v;
1237: PetscReal ferrors[1];
1239: PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs);
1240: for (f = 0; f < Nf; ++f) PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f]);
1241: DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors);
1242: PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04d time = %-8.4g \t L_2 Error: [", (int)step, (double)ptime);
1243: for (f = 0; f < Nf; ++f) {
1244: if (f > 0) PetscPrintf(PETSC_COMM_WORLD, ", ");
1245: PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f]);
1246: }
1247: PetscPrintf(PETSC_COMM_WORLD, "]\n");
1249: VecViewFromOptions(u, NULL, "-sol_vec_view");
1251: PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg);
1252: if (flg) {
1253: DMGetGlobalVector(dm, &v);
1254: DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v);
1255: PetscObjectSetName((PetscObject)v, "Exact Solution");
1256: VecViewFromOptions(v, NULL, "-exact_vec_view");
1257: DMRestoreGlobalVector(dm, &v);
1258: }
1259: PetscFree2(exactFuncs, ctxs);
1260: }
1261: return 0;
1262: }
1264: PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1265: {
1266: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1267: PetscReal x = ptime, y;
1268: PetscInt its;
1270: if (n < 0) return 0; /* -1 indicates interpolated solution */
1271: if (!n) {
1272: PetscDrawAxis axis;
1273: PetscDrawLGGetAxis(ctx->lg, &axis);
1274: PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations");
1275: PetscDrawLGReset(ctx->lg);
1276: ctx->snes_its = 0;
1277: }
1278: TSGetSNESIterations(ts, &its);
1279: y = its - ctx->snes_its;
1280: PetscDrawLGAddPoint(ctx->lg, &x, &y);
1281: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1282: PetscDrawLGDraw(ctx->lg);
1283: PetscDrawLGSave(ctx->lg);
1284: }
1285: ctx->snes_its = its;
1286: return 0;
1287: }
1289: PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx)
1290: {
1291: TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx;
1292: PetscReal x = ptime, y;
1293: PetscInt its;
1295: if (n < 0) return 0; /* -1 indicates interpolated solution */
1296: if (!n) {
1297: PetscDrawAxis axis;
1298: PetscDrawLGGetAxis(ctx->lg, &axis);
1299: PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations");
1300: PetscDrawLGReset(ctx->lg);
1301: ctx->ksp_its = 0;
1302: }
1303: TSGetKSPIterations(ts, &its);
1304: y = its - ctx->ksp_its;
1305: PetscDrawLGAddPoint(ctx->lg, &x, &y);
1306: if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) {
1307: PetscDrawLGDraw(ctx->lg);
1308: PetscDrawLGSave(ctx->lg);
1309: }
1310: ctx->ksp_its = its;
1311: return 0;
1312: }
1314: /*@C
1315: TSMonitorEnvelopeCtxCreate - Creates a context for use with TSMonitorEnvelope()
1317: Collective on TS
1319: Input Parameters:
1320: . ts - the ODE solver object
1322: Output Parameter:
1323: . ctx - the context
1325: Level: intermediate
1327: .seealso: `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`
1329: @*/
1330: PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx)
1331: {
1332: PetscNew(ctx);
1333: return 0;
1334: }
1336: /*@C
1337: TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution
1339: Collective on TS
1341: Input Parameters:
1342: + ts - the TS context
1343: . step - current time-step
1344: . ptime - current time
1345: . u - current solution
1346: - dctx - the envelope context
1348: Options Database:
1349: . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time
1351: Level: intermediate
1353: Notes:
1354: After a solve you can use TSMonitorEnvelopeGetBounds() to access the envelope
1356: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1357: to be used during the TS integration.
1359: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()`
1360: @*/
1361: PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx)
1362: {
1363: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx;
1365: if (!ctx->max) {
1366: VecDuplicate(u, &ctx->max);
1367: VecDuplicate(u, &ctx->min);
1368: VecCopy(u, ctx->max);
1369: VecCopy(u, ctx->min);
1370: } else {
1371: VecPointwiseMax(ctx->max, u, ctx->max);
1372: VecPointwiseMin(ctx->min, u, ctx->min);
1373: }
1374: return 0;
1375: }
1377: /*@C
1378: TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution
1380: Collective on TS
1382: Input Parameter:
1383: . ts - the TS context
1385: Output Parameters:
1386: + max - the maximum values
1387: - min - the minimum values
1389: Notes:
1390: If the TS does not have a TSMonitorEnvelopeCtx associated with it then this function is ignored
1392: Level: intermediate
1394: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`
1395: @*/
1396: PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min)
1397: {
1398: PetscInt i;
1400: if (max) *max = NULL;
1401: if (min) *min = NULL;
1402: for (i = 0; i < ts->numbermonitors; i++) {
1403: if (ts->monitor[i] == TSMonitorEnvelope) {
1404: TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i];
1405: if (max) *max = ctx->max;
1406: if (min) *min = ctx->min;
1407: break;
1408: }
1409: }
1410: return 0;
1411: }
1413: /*@C
1414: TSMonitorEnvelopeCtxDestroy - Destroys a context that was created with TSMonitorEnvelopeCtxCreate().
1416: Collective on TSMonitorEnvelopeCtx
1418: Input Parameter:
1419: . ctx - the monitor context
1421: Level: intermediate
1423: .seealso: `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()`
1424: @*/
1425: PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx)
1426: {
1427: VecDestroy(&(*ctx)->min);
1428: VecDestroy(&(*ctx)->max);
1429: PetscFree(*ctx);
1430: return 0;
1431: }
1433: /*@C
1434: TSDMSwarmMonitorMoments - Monitors the first three moments of a DMSarm being evolved by the TS
1436: Not collective
1438: Input Parameters:
1439: + ts - the TS context
1440: . step - current timestep
1441: . t - current time
1442: . u - current solution
1443: - ctx - not used
1445: Options Database:
1446: . -ts_dmswarm_monitor_moments - Monitor moments of particle distribution
1448: Level: intermediate
1450: Notes:
1451: This requires a DMSwarm be attached to the TS.
1453: This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor
1454: to be used during the TS integration.
1456: .seealso: `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM`
1457: @*/
1458: PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf)
1459: {
1460: DM sw;
1461: const PetscScalar *u;
1462: PetscReal m = 1.0, totE = 0., totMom[3] = {0., 0., 0.};
1463: PetscInt dim, d, Np, p;
1464: MPI_Comm comm;
1467: TSGetDM(ts, &sw);
1468: if (!sw || step % ts->monitorFrequency != 0) return 0;
1469: PetscObjectGetComm((PetscObject)ts, &comm);
1470: DMGetDimension(sw, &dim);
1471: VecGetLocalSize(U, &Np);
1472: Np /= dim;
1473: VecGetArrayRead(U, &u);
1474: for (p = 0; p < Np; ++p) {
1475: for (d = 0; d < dim; ++d) {
1476: totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]);
1477: totMom[d] += PetscRealPart(u[p * dim + d]);
1478: }
1479: }
1480: VecRestoreArrayRead(U, &u);
1481: for (d = 0; d < dim; ++d) totMom[d] *= m;
1482: totE *= 0.5 * m;
1483: PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE);
1484: for (d = 0; d < dim; ++d) PetscPrintf(comm, " Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d]);
1485: PetscPrintf(comm, "\n");
1486: return 0;
1487: }