1 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 2 #include <petscdm.h> 3 #include <petscds.h> 4 #include <petscdmswarm.h> 5 #include <petscdraw.h> 6 7 /*@C 8 TSMonitor - Runs all user-provided monitor routines set using `TSMonitorSet()` 9 10 Collective 11 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 17 18 Level: developer 19 20 Notes: 21 `TSMonitor()` is typically used automatically within the time stepping implementations. 22 Users would almost never call this routine directly. 23 24 A step of -1 indicates that the monitor is being called on a solution obtained by interpolating from computed solutions 25 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; 32 33 PetscFunctionBegin; 34 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 35 PetscValidHeaderSpecific(u, VEC_CLASSID, 4); 36 37 PetscCall(TSGetDM(ts, &dm)); 38 PetscCall(DMSetOutputSequenceNumber(dm, step, ptime)); 39 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 } 45 46 /*@C 47 TSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user 48 49 Collective 50 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 58 59 Level: developer 60 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; 74 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 } 89 90 /*@C 91 TSMonitorSet - Sets an ADDITIONAL function that is to be used at every 92 timestep to display the iteration's progress. 93 94 Logically Collective 95 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`) 101 102 Calling sequence of monitor: 103 $ PetscErrorCode monitor(TS ts,PetscInt steps,PetscReal time,Vec u,void *mctx) 104 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 110 111 Level: intermediate 112 113 Note: 114 This routine adds an additional monitor to the list of monitors that already has been loaded. 115 116 Fortran Note: 117 Only a single monitor function can be set for each `TS` object 118 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; 127 128 PetscFunctionBegin; 129 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 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 } 140 141 /*@C 142 TSMonitorCancel - Clears all the monitors that have been set on a time-step object. 143 144 Logically Collective 145 146 Input Parameters: 147 . ts - the `TS` context obtained from `TSCreate()` 148 149 Level: intermediate 150 151 Note: 152 There is no way to remove a single, specific monitor. 153 154 .seealso: [](chapter_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()` 155 @*/ 156 PetscErrorCode TSMonitorCancel(TS ts) 157 { 158 PetscInt i; 159 160 PetscFunctionBegin; 161 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 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 } 168 169 /*@C 170 TSMonitorDefault - The default monitor, prints the timestep and time for each step 171 172 Options Database Key: 173 . -ts_monitor - monitors the time integration 174 175 Level: intermediate 176 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. 180 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; 189 190 PetscFunctionBegin; 191 PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 5); 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; 209 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 } 220 221 /*@C 222 TSMonitorExtreme - Prints the extreme values of the solution at each timestep 223 224 Level: intermediate 225 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. 229 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; 237 238 PetscFunctionBegin; 239 PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 5); 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 } 252 253 /*@C 254 TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with 255 `TS` to monitor the solution process graphically in various ways 256 257 Collective 258 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 265 266 Output Parameter: 267 . ctx - the context 268 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 277 278 Level: intermediate 279 280 Notes: 281 Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed. 282 283 One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()` 284 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. 288 289 One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()` 290 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; 300 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 } 311 312 PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx) 313 { 314 TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx; 315 PetscReal x = ptime, y; 316 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 } 335 336 /*@C 337 TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`. 338 339 Collective 340 341 Input Parameter: 342 . ctx - the monitor context 343 344 Level: intermediate 345 346 Note: 347 Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()` 348 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 } 363 364 /* Creates a TS Monitor SPCtx 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, TSMonitorSPCtx *ctx) 366 { 367 PetscDraw draw; 368 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 PetscFunctionReturn(PETSC_SUCCESS); 379 } 380 381 /* 382 Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate 383 */ 384 PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx) 385 { 386 PetscFunctionBegin; 387 388 PetscCall(PetscDrawSPDestroy(&(*ctx)->sp)); 389 PetscCall(PetscFree(*ctx)); 390 391 PetscFunctionReturn(PETSC_SUCCESS); 392 } 393 394 /*@C 395 TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling 396 `VecView()` for the solution at each timestep 397 398 Collective 399 400 Input Parameters: 401 + ts - the `TS` context 402 . step - current time-step 403 . ptime - current time 404 - dummy - either a viewer or `NULL` 405 406 Options Database Keys: 407 + -ts_monitor_draw_solution - draw the solution at each time-step 408 - -ts_monitor_draw_solution_initial - show initial solution as well as current solution 409 410 Level: intermediate 411 412 Notes: 413 The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial` 414 will look bad 415 416 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with 417 `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration. 418 419 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()` 420 @*/ 421 PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 422 { 423 TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy; 424 PetscDraw draw; 425 426 PetscFunctionBegin; 427 if (!step && ictx->showinitial) { 428 if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution)); 429 PetscCall(VecCopy(u, ictx->initialsolution)); 430 } 431 if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS); 432 433 if (ictx->showinitial) { 434 PetscReal pause; 435 PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause)); 436 PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0)); 437 PetscCall(VecView(ictx->initialsolution, ictx->viewer)); 438 PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause)); 439 PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE)); 440 } 441 PetscCall(VecView(u, ictx->viewer)); 442 if (ictx->showtimestepandtime) { 443 PetscReal xl, yl, xr, yr, h; 444 char time[32]; 445 446 PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw)); 447 PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime)); 448 PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr)); 449 h = yl + .95 * (yr - yl); 450 PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time)); 451 PetscCall(PetscDrawFlush(draw)); 452 } 453 454 if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE)); 455 PetscFunctionReturn(PETSC_SUCCESS); 456 } 457 458 /*@C 459 TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram 460 461 Collective 462 463 Input Parameters: 464 + ts - the `TS` context 465 . step - current time-step 466 . ptime - current time 467 - dummy - either a viewer or `NULL` 468 469 Level: intermediate 470 471 Notes: 472 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 473 to be used during the `TS` integration. 474 475 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()` 476 @*/ 477 PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 478 { 479 TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy; 480 PetscDraw draw; 481 PetscDrawAxis axis; 482 PetscInt n; 483 PetscMPIInt size; 484 PetscReal U0, U1, xl, yl, xr, yr, h; 485 char time[32]; 486 const PetscScalar *U; 487 488 PetscFunctionBegin; 489 PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size)); 490 PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs"); 491 PetscCall(VecGetSize(u, &n)); 492 PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns"); 493 494 PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw)); 495 PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis)); 496 PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr)); 497 if (!step) { 498 PetscCall(PetscDrawClear(draw)); 499 PetscCall(PetscDrawAxisDraw(axis)); 500 } 501 502 PetscCall(VecGetArrayRead(u, &U)); 503 U0 = PetscRealPart(U[0]); 504 U1 = PetscRealPart(U[1]); 505 PetscCall(VecRestoreArrayRead(u, &U)); 506 if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS); 507 508 PetscDrawCollectiveBegin(draw); 509 PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK)); 510 if (ictx->showtimestepandtime) { 511 PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr)); 512 PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime)); 513 h = yl + .95 * (yr - yl); 514 PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time)); 515 } 516 PetscDrawCollectiveEnd(draw); 517 PetscCall(PetscDrawFlush(draw)); 518 PetscCall(PetscDrawPause(draw)); 519 PetscCall(PetscDrawSave(draw)); 520 PetscFunctionReturn(PETSC_SUCCESS); 521 } 522 523 /*@C 524 TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()` 525 526 Collective 527 528 Input Parameters: 529 . ctx - the monitor context 530 531 Level: intermediate 532 533 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx` 534 @*/ 535 PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx) 536 { 537 PetscFunctionBegin; 538 PetscCall(PetscViewerDestroy(&(*ictx)->viewer)); 539 PetscCall(VecDestroy(&(*ictx)->initialsolution)); 540 PetscCall(PetscFree(*ictx)); 541 PetscFunctionReturn(PETSC_SUCCESS); 542 } 543 544 /*@C 545 TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx` 546 547 Collective 548 549 Input Parameter: 550 . ts - time-step context 551 552 Output Parameter: 553 . ctx - the monitor context 554 555 Options Database Keys: 556 + -ts_monitor_draw_solution - draw the solution at each time-step 557 - -ts_monitor_draw_solution_initial - show initial solution as well as current solution 558 559 Level: intermediate 560 561 Note: 562 The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`. 563 564 .seealso: [](chapter_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()` 565 @*/ 566 PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx) 567 { 568 PetscFunctionBegin; 569 PetscCall(PetscNew(ctx)); 570 PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer)); 571 PetscCall(PetscViewerSetFromOptions((*ctx)->viewer)); 572 573 (*ctx)->howoften = howoften; 574 (*ctx)->showinitial = PETSC_FALSE; 575 PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL)); 576 577 (*ctx)->showtimestepandtime = PETSC_FALSE; 578 PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL)); 579 PetscFunctionReturn(PETSC_SUCCESS); 580 } 581 582 /*@C 583 TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling 584 `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep 585 586 Collective 587 588 Input Parameters: 589 + ts - the `TS` context 590 . step - current time-step 591 . ptime - current time 592 - dummy - either a viewer or `NULL` 593 594 Options Database Key: 595 . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()` 596 597 Level: intermediate 598 599 Note: 600 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 601 to be used during the `TS` integration. 602 603 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 604 @*/ 605 PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 606 { 607 TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy; 608 PetscViewer viewer = ctx->viewer; 609 Vec work; 610 611 PetscFunctionBegin; 612 if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS); 613 PetscCall(VecDuplicate(u, &work)); 614 PetscCall(TSComputeSolutionFunction(ts, ptime, work)); 615 PetscCall(VecView(work, viewer)); 616 PetscCall(VecDestroy(&work)); 617 PetscFunctionReturn(PETSC_SUCCESS); 618 } 619 620 /*@C 621 TSMonitorDrawError - Monitors progress of the `TS` solvers by calling 622 `VecView()` for the error at each timestep 623 624 Collective 625 626 Input Parameters: 627 + ts - the `TS` context 628 . step - current time-step 629 . ptime - current time 630 - dummy - either a viewer or `NULL` 631 632 Options Database Key: 633 . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()` 634 635 Level: intermediate 636 637 Notes: 638 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 639 to be used during the `TS` integration. 640 641 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 642 @*/ 643 PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 644 { 645 TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy; 646 PetscViewer viewer = ctx->viewer; 647 Vec work; 648 649 PetscFunctionBegin; 650 if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS); 651 PetscCall(VecDuplicate(u, &work)); 652 PetscCall(TSComputeSolutionFunction(ts, ptime, work)); 653 PetscCall(VecAXPY(work, -1.0, u)); 654 PetscCall(VecView(work, viewer)); 655 PetscCall(VecDestroy(&work)); 656 PetscFunctionReturn(PETSC_SUCCESS); 657 } 658 659 /*@C 660 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 661 662 Collective 663 664 Input Parameters: 665 + ts - the `TS` context 666 . step - current time-step 667 . ptime - current time 668 . u - current state 669 - vf - viewer and its format 670 671 Level: intermediate 672 673 Notes: 674 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 675 to be used during the `TS` integration. 676 677 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()` 678 @*/ 679 PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf) 680 { 681 PetscFunctionBegin; 682 if (vf->view_interval > 0 && !ts->reason && step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS); 683 PetscCall(PetscViewerPushFormat(vf->viewer, vf->format)); 684 PetscCall(VecView(u, vf->viewer)); 685 PetscCall(PetscViewerPopFormat(vf->viewer)); 686 PetscFunctionReturn(PETSC_SUCCESS); 687 } 688 689 /*@C 690 TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at each timestep. 691 692 Collective 693 694 Input Parameters: 695 + ts - the `TS` context 696 . step - current time-step 697 . ptime - current time 698 . u - current state 699 - filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ") 700 701 Level: intermediate 702 703 Notes: 704 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. 705 These are named according to the file name template. 706 707 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 708 to be used during the `TS` integration. 709 710 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()` 711 @*/ 712 PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, void *filenametemplate) 713 { 714 char filename[PETSC_MAX_PATH_LEN]; 715 PetscViewer viewer; 716 717 PetscFunctionBegin; 718 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 719 PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)filenametemplate, step)); 720 PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer)); 721 PetscCall(VecView(u, viewer)); 722 PetscCall(PetscViewerDestroy(&viewer)); 723 PetscFunctionReturn(PETSC_SUCCESS); 724 } 725 726 /*@C 727 TSMonitorSolutionVTKDestroy - Destroy filename template string created for use with `TSMonitorSolutionVTK()` 728 729 Not Collective 730 731 Input Parameters: 732 . filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ") 733 734 Level: intermediate 735 736 Note: 737 This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`. 738 739 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()` 740 @*/ 741 PetscErrorCode TSMonitorSolutionVTKDestroy(void *filenametemplate) 742 { 743 PetscFunctionBegin; 744 PetscCall(PetscFree(*(char **)filenametemplate)); 745 PetscFunctionReturn(PETSC_SUCCESS); 746 } 747 748 /*@C 749 TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector 750 in a time based line graph 751 752 Collective 753 754 Input Parameters: 755 + ts - the `TS` context 756 . step - current time-step 757 . ptime - current time 758 . u - current solution 759 - dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()` 760 761 Options Database Key: 762 . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables 763 764 Level: intermediate 765 766 Notes: 767 Each process in a parallel run displays its component solutions in a separate window 768 769 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 770 to be used during the `TS` integration. 771 772 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`, 773 `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`, 774 `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`, 775 `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()` 776 @*/ 777 PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx) 778 { 779 TSMonitorLGCtx ctx = (TSMonitorLGCtx)dctx; 780 const PetscScalar *yy; 781 Vec v; 782 783 PetscFunctionBegin; 784 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 785 if (!step) { 786 PetscDrawAxis axis; 787 PetscInt dim; 788 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 789 PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution")); 790 if (!ctx->names) { 791 PetscBool flg; 792 /* user provides names of variables to plot but no names has been set so assume names are integer values */ 793 PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg)); 794 if (flg) { 795 PetscInt i, n; 796 char **names; 797 PetscCall(VecGetSize(u, &n)); 798 PetscCall(PetscMalloc1(n + 1, &names)); 799 for (i = 0; i < n; i++) { 800 PetscCall(PetscMalloc1(5, &names[i])); 801 PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i)); 802 } 803 names[n] = NULL; 804 ctx->names = names; 805 } 806 } 807 if (ctx->names && !ctx->displaynames) { 808 char **displaynames; 809 PetscBool flg; 810 PetscCall(VecGetLocalSize(u, &dim)); 811 PetscCall(PetscCalloc1(dim + 1, &displaynames)); 812 PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg)); 813 if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames)); 814 PetscCall(PetscStrArrayDestroy(&displaynames)); 815 } 816 if (ctx->displaynames) { 817 PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables)); 818 PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames)); 819 } else if (ctx->names) { 820 PetscCall(VecGetLocalSize(u, &dim)); 821 PetscCall(PetscDrawLGSetDimension(ctx->lg, dim)); 822 PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names)); 823 } else { 824 PetscCall(VecGetLocalSize(u, &dim)); 825 PetscCall(PetscDrawLGSetDimension(ctx->lg, dim)); 826 } 827 PetscCall(PetscDrawLGReset(ctx->lg)); 828 } 829 830 if (!ctx->transform) v = u; 831 else PetscCall((*ctx->transform)(ctx->transformctx, u, &v)); 832 PetscCall(VecGetArrayRead(v, &yy)); 833 if (ctx->displaynames) { 834 PetscInt i; 835 for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]); 836 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues)); 837 } else { 838 #if defined(PETSC_USE_COMPLEX) 839 PetscInt i, n; 840 PetscReal *yreal; 841 PetscCall(VecGetLocalSize(v, &n)); 842 PetscCall(PetscMalloc1(n, &yreal)); 843 for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]); 844 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal)); 845 PetscCall(PetscFree(yreal)); 846 #else 847 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy)); 848 #endif 849 } 850 PetscCall(VecRestoreArrayRead(v, &yy)); 851 if (ctx->transform) PetscCall(VecDestroy(&v)); 852 853 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 854 PetscCall(PetscDrawLGDraw(ctx->lg)); 855 PetscCall(PetscDrawLGSave(ctx->lg)); 856 } 857 PetscFunctionReturn(PETSC_SUCCESS); 858 } 859 860 /*@C 861 TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot 862 863 Collective 864 865 Input Parameters: 866 + ts - the `TS` context 867 - names - the names of the components, final string must be `NULL` 868 869 Level: intermediate 870 871 Notes: 872 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 873 874 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()` 875 @*/ 876 PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names) 877 { 878 PetscInt i; 879 880 PetscFunctionBegin; 881 for (i = 0; i < ts->numbermonitors; i++) { 882 if (ts->monitor[i] == TSMonitorLGSolution) { 883 PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names)); 884 break; 885 } 886 } 887 PetscFunctionReturn(PETSC_SUCCESS); 888 } 889 890 /*@C 891 TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot 892 893 Collective 894 895 Input Parameters: 896 + ts - the `TS` context 897 - names - the names of the components, final string must be `NULL` 898 899 Level: intermediate 900 901 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()` 902 @*/ 903 PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names) 904 { 905 PetscFunctionBegin; 906 PetscCall(PetscStrArrayDestroy(&ctx->names)); 907 PetscCall(PetscStrArrayallocpy(names, &ctx->names)); 908 PetscFunctionReturn(PETSC_SUCCESS); 909 } 910 911 /*@C 912 TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot 913 914 Collective 915 916 Input Parameter: 917 . ts - the `TS` context 918 919 Output Parameter: 920 . names - the names of the components, final string must be `NULL` 921 922 Level: intermediate 923 924 Note: 925 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 926 927 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()` 928 @*/ 929 PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names) 930 { 931 PetscInt i; 932 933 PetscFunctionBegin; 934 *names = NULL; 935 for (i = 0; i < ts->numbermonitors; i++) { 936 if (ts->monitor[i] == TSMonitorLGSolution) { 937 TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i]; 938 *names = (const char *const *)ctx->names; 939 break; 940 } 941 } 942 PetscFunctionReturn(PETSC_SUCCESS); 943 } 944 945 /*@C 946 TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor 947 948 Collective 949 950 Input Parameters: 951 + ctx - the `TSMonitorLG` context 952 - displaynames - the names of the components, final string must be `NULL` 953 954 Level: intermediate 955 956 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()` 957 @*/ 958 PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames) 959 { 960 PetscInt j = 0, k; 961 962 PetscFunctionBegin; 963 if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS); 964 PetscCall(PetscStrArrayDestroy(&ctx->displaynames)); 965 PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames)); 966 while (displaynames[j]) j++; 967 ctx->ndisplayvariables = j; 968 PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables)); 969 PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues)); 970 j = 0; 971 while (displaynames[j]) { 972 k = 0; 973 while (ctx->names[k]) { 974 PetscBool flg; 975 PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg)); 976 if (flg) { 977 ctx->displayvariables[j] = k; 978 break; 979 } 980 k++; 981 } 982 j++; 983 } 984 PetscFunctionReturn(PETSC_SUCCESS); 985 } 986 987 /*@C 988 TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor 989 990 Collective 991 992 Input Parameters: 993 + ts - the `TS` context 994 - displaynames - the names of the components, final string must be `NULL` 995 996 Level: intermediate 997 998 Note: 999 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 1000 1001 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()` 1002 @*/ 1003 PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames) 1004 { 1005 PetscInt i; 1006 1007 PetscFunctionBegin; 1008 for (i = 0; i < ts->numbermonitors; i++) { 1009 if (ts->monitor[i] == TSMonitorLGSolution) { 1010 PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames)); 1011 break; 1012 } 1013 } 1014 PetscFunctionReturn(PETSC_SUCCESS); 1015 } 1016 1017 /*@C 1018 TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed 1019 1020 Collective 1021 1022 Input Parameters: 1023 + ts - the `TS` context 1024 . transform - the transform function 1025 . destroy - function to destroy the optional context 1026 - ctx - optional context used by transform function 1027 1028 Level: intermediate 1029 1030 Note: 1031 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 1032 1033 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()` 1034 @*/ 1035 PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx) 1036 { 1037 PetscInt i; 1038 1039 PetscFunctionBegin; 1040 for (i = 0; i < ts->numbermonitors; i++) { 1041 if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx)); 1042 } 1043 PetscFunctionReturn(PETSC_SUCCESS); 1044 } 1045 1046 /*@C 1047 TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed 1048 1049 Collective 1050 1051 Input Parameters: 1052 + ts - the `TS` context 1053 . transform - the transform function 1054 . destroy - function to destroy the optional context 1055 - ctx - optional context used by transform function 1056 1057 Level: intermediate 1058 1059 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()` 1060 @*/ 1061 PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx) 1062 { 1063 PetscFunctionBegin; 1064 ctx->transform = transform; 1065 ctx->transformdestroy = destroy; 1066 ctx->transformctx = tctx; 1067 PetscFunctionReturn(PETSC_SUCCESS); 1068 } 1069 1070 /*@C 1071 TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error 1072 in a time based line graph 1073 1074 Collective 1075 1076 Input Parameters: 1077 + ts - the `TS` context 1078 . step - current time-step 1079 . ptime - current time 1080 . u - current solution 1081 - dctx - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()` 1082 1083 Options Database Key: 1084 . -ts_monitor_lg_error - create a graphical monitor of error history 1085 1086 Level: intermediate 1087 1088 Notes: 1089 Each process in a parallel run displays its component errors in a separate window 1090 1091 The user must provide the solution using `TSSetSolutionFunction()` to use this monitor. 1092 1093 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1094 to be used during the TS integration. 1095 1096 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 1097 @*/ 1098 PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 1099 { 1100 TSMonitorLGCtx ctx = (TSMonitorLGCtx)dummy; 1101 const PetscScalar *yy; 1102 Vec y; 1103 1104 PetscFunctionBegin; 1105 if (!step) { 1106 PetscDrawAxis axis; 1107 PetscInt dim; 1108 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 1109 PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error")); 1110 PetscCall(VecGetLocalSize(u, &dim)); 1111 PetscCall(PetscDrawLGSetDimension(ctx->lg, dim)); 1112 PetscCall(PetscDrawLGReset(ctx->lg)); 1113 } 1114 PetscCall(VecDuplicate(u, &y)); 1115 PetscCall(TSComputeSolutionFunction(ts, ptime, y)); 1116 PetscCall(VecAXPY(y, -1.0, u)); 1117 PetscCall(VecGetArrayRead(y, &yy)); 1118 #if defined(PETSC_USE_COMPLEX) 1119 { 1120 PetscReal *yreal; 1121 PetscInt i, n; 1122 PetscCall(VecGetLocalSize(y, &n)); 1123 PetscCall(PetscMalloc1(n, &yreal)); 1124 for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]); 1125 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal)); 1126 PetscCall(PetscFree(yreal)); 1127 } 1128 #else 1129 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy)); 1130 #endif 1131 PetscCall(VecRestoreArrayRead(y, &yy)); 1132 PetscCall(VecDestroy(&y)); 1133 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 1134 PetscCall(PetscDrawLGDraw(ctx->lg)); 1135 PetscCall(PetscDrawLGSave(ctx->lg)); 1136 } 1137 PetscFunctionReturn(PETSC_SUCCESS); 1138 } 1139 1140 /*@C 1141 TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot 1142 1143 Input Parameters: 1144 + ts - the `TS` context 1145 . step - current time-step 1146 . ptime - current time 1147 . u - current solution 1148 - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()` 1149 1150 Options Database Keys: 1151 + -ts_monitor_sp_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution 1152 . -ts_monitor_sp_swarm_retain <n> - Retain n old points so we can see the history, or -1 for all points 1153 - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space 1154 1155 Level: intermediate 1156 1157 Notes: 1158 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1159 to be used during the `TS` integration. 1160 1161 .seealso: [](chapter_ts), `TS`, `TSMonitoSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()` 1162 @*/ 1163 PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx) 1164 { 1165 TSMonitorSPCtx ctx = (TSMonitorSPCtx)dctx; 1166 PetscDraw draw; 1167 DM dm, cdm; 1168 const PetscScalar *yy; 1169 PetscInt Np, p, dim = 2; 1170 1171 PetscFunctionBegin; 1172 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 1173 if (!step) { 1174 PetscDrawAxis axis; 1175 PetscReal dmboxlower[2], dmboxupper[2]; 1176 1177 PetscCall(TSGetDM(ts, &dm)); 1178 PetscCall(DMGetDimension(dm, &dim)); 1179 PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields"); 1180 PetscCall(DMSwarmGetCellDM(dm, &cdm)); 1181 PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper)); 1182 PetscCall(VecGetLocalSize(u, &Np)); 1183 Np /= dim * 2; 1184 PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis)); 1185 if (ctx->phase) { 1186 PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V")); 1187 PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -5, 5)); 1188 } else { 1189 PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y")); 1190 PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1])); 1191 } 1192 PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE)); 1193 PetscCall(PetscDrawSPReset(ctx->sp)); 1194 } 1195 PetscCall(VecGetLocalSize(u, &Np)); 1196 Np /= dim * 2; 1197 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 1198 PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw)); 1199 if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw)); 1200 PetscCall(PetscDrawFlush(draw)); 1201 PetscCall(PetscDrawSPReset(ctx->sp)); 1202 PetscCall(VecGetArrayRead(u, &yy)); 1203 for (p = 0; p < Np; ++p) { 1204 PetscReal x, y; 1205 1206 if (ctx->phase) { 1207 x = PetscRealPart(yy[p * dim * 2]); 1208 y = PetscRealPart(yy[p * dim * 2 + dim]); 1209 } else { 1210 x = PetscRealPart(yy[p * dim * 2]); 1211 y = PetscRealPart(yy[p * dim * 2 + 1]); 1212 } 1213 PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y)); 1214 } 1215 PetscCall(VecRestoreArrayRead(u, &yy)); 1216 PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE)); 1217 PetscCall(PetscDrawSPSave(ctx->sp)); 1218 } 1219 PetscFunctionReturn(PETSC_SUCCESS); 1220 } 1221 1222 /*@C 1223 TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep 1224 1225 Collective 1226 1227 Input Parameters: 1228 + ts - the `TS` context 1229 . step - current time-step 1230 . ptime - current time 1231 . u - current solution 1232 - dctx - unused context 1233 1234 Options Database Key: 1235 . -ts_monitor_error - create a graphical monitor of error history 1236 1237 Level: intermediate 1238 1239 Notes: 1240 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1241 to be used during the `TS` integration. 1242 1243 The user must provide the solution using `TSSetSolutionFunction()` to use this monitor. 1244 1245 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 1246 @*/ 1247 PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf) 1248 { 1249 DM dm; 1250 PetscDS ds = NULL; 1251 PetscInt Nf = -1, f; 1252 PetscBool flg; 1253 1254 PetscFunctionBegin; 1255 PetscCall(TSGetDM(ts, &dm)); 1256 if (dm) PetscCall(DMGetDS(dm, &ds)); 1257 if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf)); 1258 if (Nf <= 0) { 1259 Vec y; 1260 PetscReal nrm; 1261 1262 PetscCall(VecDuplicate(u, &y)); 1263 PetscCall(TSComputeSolutionFunction(ts, ptime, y)); 1264 PetscCall(VecAXPY(y, -1.0, u)); 1265 PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg)); 1266 if (flg) { 1267 PetscCall(VecNorm(y, NORM_2, &nrm)); 1268 PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm)); 1269 } 1270 PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg)); 1271 if (flg) PetscCall(VecView(y, vf->viewer)); 1272 PetscCall(VecDestroy(&y)); 1273 } else { 1274 PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); 1275 void **ctxs; 1276 Vec v; 1277 PetscReal ferrors[1]; 1278 1279 PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs)); 1280 for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f])); 1281 PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors)); 1282 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04d time = %-8.4g \t L_2 Error: [", (int)step, (double)ptime)); 1283 for (f = 0; f < Nf; ++f) { 1284 if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", ")); 1285 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f])); 1286 } 1287 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n")); 1288 1289 PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view")); 1290 1291 PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg)); 1292 if (flg) { 1293 PetscCall(DMGetGlobalVector(dm, &v)); 1294 PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v)); 1295 PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution")); 1296 PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view")); 1297 PetscCall(DMRestoreGlobalVector(dm, &v)); 1298 } 1299 PetscCall(PetscFree2(exactFuncs, ctxs)); 1300 } 1301 PetscFunctionReturn(PETSC_SUCCESS); 1302 } 1303 1304 PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx) 1305 { 1306 TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx; 1307 PetscReal x = ptime, y; 1308 PetscInt its; 1309 1310 PetscFunctionBegin; 1311 if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 1312 if (!n) { 1313 PetscDrawAxis axis; 1314 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 1315 PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations")); 1316 PetscCall(PetscDrawLGReset(ctx->lg)); 1317 ctx->snes_its = 0; 1318 } 1319 PetscCall(TSGetSNESIterations(ts, &its)); 1320 y = its - ctx->snes_its; 1321 PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y)); 1322 if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) { 1323 PetscCall(PetscDrawLGDraw(ctx->lg)); 1324 PetscCall(PetscDrawLGSave(ctx->lg)); 1325 } 1326 ctx->snes_its = its; 1327 PetscFunctionReturn(PETSC_SUCCESS); 1328 } 1329 1330 PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx) 1331 { 1332 TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx; 1333 PetscReal x = ptime, y; 1334 PetscInt its; 1335 1336 PetscFunctionBegin; 1337 if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 1338 if (!n) { 1339 PetscDrawAxis axis; 1340 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 1341 PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations")); 1342 PetscCall(PetscDrawLGReset(ctx->lg)); 1343 ctx->ksp_its = 0; 1344 } 1345 PetscCall(TSGetKSPIterations(ts, &its)); 1346 y = its - ctx->ksp_its; 1347 PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y)); 1348 if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) { 1349 PetscCall(PetscDrawLGDraw(ctx->lg)); 1350 PetscCall(PetscDrawLGSave(ctx->lg)); 1351 } 1352 ctx->ksp_its = its; 1353 PetscFunctionReturn(PETSC_SUCCESS); 1354 } 1355 1356 /*@C 1357 TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()` 1358 1359 Collective 1360 1361 Input Parameters: 1362 . ts - the `TS` solver object 1363 1364 Output Parameter: 1365 . ctx - the context 1366 1367 Level: intermediate 1368 1369 .seealso: [](chapter_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()` 1370 @*/ 1371 PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx) 1372 { 1373 PetscFunctionBegin; 1374 PetscCall(PetscNew(ctx)); 1375 PetscFunctionReturn(PETSC_SUCCESS); 1376 } 1377 1378 /*@C 1379 TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution 1380 1381 Collective 1382 1383 Input Parameters: 1384 + ts - the `TS` context 1385 . step - current time-step 1386 . ptime - current time 1387 . u - current solution 1388 - dctx - the envelope context 1389 1390 Options Database Key: 1391 . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time 1392 1393 Level: intermediate 1394 1395 Notes: 1396 After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope 1397 1398 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1399 to be used during the `TS` integration. 1400 1401 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()` 1402 @*/ 1403 PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx) 1404 { 1405 TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx; 1406 1407 PetscFunctionBegin; 1408 if (!ctx->max) { 1409 PetscCall(VecDuplicate(u, &ctx->max)); 1410 PetscCall(VecDuplicate(u, &ctx->min)); 1411 PetscCall(VecCopy(u, ctx->max)); 1412 PetscCall(VecCopy(u, ctx->min)); 1413 } else { 1414 PetscCall(VecPointwiseMax(ctx->max, u, ctx->max)); 1415 PetscCall(VecPointwiseMin(ctx->min, u, ctx->min)); 1416 } 1417 PetscFunctionReturn(PETSC_SUCCESS); 1418 } 1419 1420 /*@C 1421 TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution 1422 1423 Collective 1424 1425 Input Parameter: 1426 . ts - the `TS` context 1427 1428 Output Parameters: 1429 + max - the maximum values 1430 - min - the minimum values 1431 1432 Level: intermediate 1433 1434 Notes: 1435 If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored 1436 1437 .seealso: [](chapter_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()` 1438 @*/ 1439 PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min) 1440 { 1441 PetscInt i; 1442 1443 PetscFunctionBegin; 1444 if (max) *max = NULL; 1445 if (min) *min = NULL; 1446 for (i = 0; i < ts->numbermonitors; i++) { 1447 if (ts->monitor[i] == TSMonitorEnvelope) { 1448 TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i]; 1449 if (max) *max = ctx->max; 1450 if (min) *min = ctx->min; 1451 break; 1452 } 1453 } 1454 PetscFunctionReturn(PETSC_SUCCESS); 1455 } 1456 1457 /*@C 1458 TSMonitorEnvelopeCtxDestroy - Destroys a context that was created with `TSMonitorEnvelopeCtxCreate()`. 1459 1460 Collective 1461 1462 Input Parameter: 1463 . ctx - the monitor context 1464 1465 Level: intermediate 1466 1467 .seealso: [](chapter_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()` 1468 @*/ 1469 PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx) 1470 { 1471 PetscFunctionBegin; 1472 PetscCall(VecDestroy(&(*ctx)->min)); 1473 PetscCall(VecDestroy(&(*ctx)->max)); 1474 PetscCall(PetscFree(*ctx)); 1475 PetscFunctionReturn(PETSC_SUCCESS); 1476 } 1477 1478 /*@C 1479 TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS` 1480 1481 Not collective 1482 1483 Input Parameters: 1484 + ts - the `TS` context 1485 . step - current timestep 1486 . t - current time 1487 . u - current solution 1488 - ctx - not used 1489 1490 Options Database Key: 1491 . -ts_dmswarm_monitor_moments - Monitor moments of particle distribution 1492 1493 Level: intermediate 1494 1495 Notes: 1496 This requires a `DMSWARM` be attached to the `TS`. 1497 1498 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1499 to be used during the TS integration. 1500 1501 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM` 1502 @*/ 1503 PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf) 1504 { 1505 DM sw; 1506 const PetscScalar *u; 1507 PetscReal m = 1.0, totE = 0., totMom[3] = {0., 0., 0.}; 1508 PetscInt dim, d, Np, p; 1509 MPI_Comm comm; 1510 1511 PetscFunctionBeginUser; 1512 PetscCall(TSGetDM(ts, &sw)); 1513 if (!sw || step % ts->monitorFrequency != 0) PetscFunctionReturn(PETSC_SUCCESS); 1514 PetscCall(PetscObjectGetComm((PetscObject)ts, &comm)); 1515 PetscCall(DMGetDimension(sw, &dim)); 1516 PetscCall(VecGetLocalSize(U, &Np)); 1517 Np /= dim; 1518 PetscCall(VecGetArrayRead(U, &u)); 1519 for (p = 0; p < Np; ++p) { 1520 for (d = 0; d < dim; ++d) { 1521 totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]); 1522 totMom[d] += PetscRealPart(u[p * dim + d]); 1523 } 1524 } 1525 PetscCall(VecRestoreArrayRead(U, &u)); 1526 for (d = 0; d < dim; ++d) totMom[d] *= m; 1527 totE *= 0.5 * m; 1528 PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE)); 1529 for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, " Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d])); 1530 PetscCall(PetscPrintf(comm, "\n")); 1531 PetscFunctionReturn(PETSC_SUCCESS); 1532 } 1533