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 + ts - the `TS` context 105 . 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) 106 . time - current time 107 . u - current iterate 108 - mctx - [optional] monitoring context 109 110 Level: intermediate 111 112 Note: 113 This routine adds an additional monitor to the list of monitors that already has been loaded. 114 115 Fortran Note: 116 Only a single monitor function can be set for each `TS` object 117 118 .seealso: [](chapter_ts), `TSMonitorDefault()`, `TSMonitorCancel()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`, 119 `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`, 120 `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()` 121 @*/ 122 PetscErrorCode TSMonitorSet(TS ts, PetscErrorCode (*monitor)(TS, PetscInt, PetscReal, Vec, void *), void *mctx, PetscErrorCode (*mdestroy)(void **)) 123 { 124 PetscInt i; 125 PetscBool identical; 126 127 PetscFunctionBegin; 128 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 129 for (i = 0; i < ts->numbermonitors; i++) { 130 PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))monitor, mctx, mdestroy, (PetscErrorCode(*)(void))ts->monitor[i], ts->monitorcontext[i], ts->monitordestroy[i], &identical)); 131 if (identical) PetscFunctionReturn(PETSC_SUCCESS); 132 } 133 PetscCheck(ts->numbermonitors < MAXTSMONITORS, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Too many monitors set"); 134 ts->monitor[ts->numbermonitors] = monitor; 135 ts->monitordestroy[ts->numbermonitors] = mdestroy; 136 ts->monitorcontext[ts->numbermonitors++] = (void *)mctx; 137 PetscFunctionReturn(PETSC_SUCCESS); 138 } 139 140 /*@C 141 TSMonitorCancel - Clears all the monitors that have been set on a time-step object. 142 143 Logically Collective 144 145 Input Parameters: 146 . ts - the `TS` context obtained from `TSCreate()` 147 148 Level: intermediate 149 150 Note: 151 There is no way to remove a single, specific monitor. 152 153 .seealso: [](chapter_ts), `TS`, `TSMonitorDefault()`, `TSMonitorSet()` 154 @*/ 155 PetscErrorCode TSMonitorCancel(TS ts) 156 { 157 PetscInt i; 158 159 PetscFunctionBegin; 160 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 161 for (i = 0; i < ts->numbermonitors; i++) { 162 if (ts->monitordestroy[i]) PetscCall((*ts->monitordestroy[i])(&ts->monitorcontext[i])); 163 } 164 ts->numbermonitors = 0; 165 PetscFunctionReturn(PETSC_SUCCESS); 166 } 167 168 /*@C 169 TSMonitorDefault - The default monitor, prints the timestep and time for each step 170 171 Options Database Key: 172 . -ts_monitor - monitors the time integration 173 174 Level: intermediate 175 176 Notes: 177 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 178 to be used during the `TS` integration. 179 180 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSDMSwarmMonitorMoments()`, `TSMonitorExtreme()`, `TSMonitorDrawSolution()`, 181 `TSMonitorDrawSolutionPhase()`, `TSMonitorDrawSolutionFunction()`, `TSMonitorDrawError()`, `TSMonitorSolution()`, `TSMonitorSolutionVTK()`, 182 `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorSPSwarmSolution()`, `TSMonitorError()`, `TSMonitorEnvelope()`, `TSDMSwarmMonitorMoments()` 183 @*/ 184 PetscErrorCode TSMonitorDefault(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf) 185 { 186 PetscViewer viewer = vf->viewer; 187 PetscBool iascii, ibinary; 188 189 PetscFunctionBegin; 190 PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 5); 191 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 192 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &ibinary)); 193 PetscCall(PetscViewerPushFormat(viewer, vf->format)); 194 if (iascii) { 195 PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel)); 196 if (step == -1) { /* this indicates it is an interpolated solution */ 197 PetscCall(PetscViewerASCIIPrintf(viewer, "Interpolated solution at time %g between steps %" PetscInt_FMT " and %" PetscInt_FMT "\n", (double)ptime, ts->steps - 1, ts->steps)); 198 } else { 199 PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " TS dt %g time %g%s", step, (double)ts->time_step, (double)ptime, ts->steprollback ? " (r)\n" : "\n")); 200 } 201 PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel)); 202 } else if (ibinary) { 203 PetscMPIInt rank; 204 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)viewer), &rank)); 205 if (rank == 0) { 206 PetscBool skipHeader; 207 PetscInt classid = REAL_FILE_CLASSID; 208 209 PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &skipHeader)); 210 if (!skipHeader) PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT)); 211 PetscCall(PetscRealView(1, &ptime, viewer)); 212 } else { 213 PetscCall(PetscRealView(0, &ptime, viewer)); 214 } 215 } 216 PetscCall(PetscViewerPopFormat(viewer)); 217 PetscFunctionReturn(PETSC_SUCCESS); 218 } 219 220 /*@C 221 TSMonitorExtreme - Prints the extreme values of the solution at each timestep 222 223 Level: intermediate 224 225 Note: 226 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 227 to be used during the `TS` integration. 228 229 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()` 230 @*/ 231 PetscErrorCode TSMonitorExtreme(TS ts, PetscInt step, PetscReal ptime, Vec v, PetscViewerAndFormat *vf) 232 { 233 PetscViewer viewer = vf->viewer; 234 PetscBool iascii; 235 PetscReal max, min; 236 237 PetscFunctionBegin; 238 PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 5); 239 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 240 PetscCall(PetscViewerPushFormat(viewer, vf->format)); 241 if (iascii) { 242 PetscCall(VecMax(v, NULL, &max)); 243 PetscCall(VecMin(v, NULL, &min)); 244 PetscCall(PetscViewerASCIIAddTab(viewer, ((PetscObject)ts)->tablevel)); 245 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)); 246 PetscCall(PetscViewerASCIISubtractTab(viewer, ((PetscObject)ts)->tablevel)); 247 } 248 PetscCall(PetscViewerPopFormat(viewer)); 249 PetscFunctionReturn(PETSC_SUCCESS); 250 } 251 252 /*@C 253 TSMonitorLGCtxCreate - Creates a `TSMonitorLGCtx` context for use with 254 `TS` to monitor the solution process graphically in various ways 255 256 Collective 257 258 Input Parameters: 259 + host - the X display to open, or `NULL` for the local machine 260 . label - the title to put in the title bar 261 . x, y - the screen coordinates of the upper left coordinate of the window 262 . m, n - the screen width and height in pixels 263 - howoften - if positive then determines the frequency of the plotting, if -1 then only at the final time 264 265 Output Parameter: 266 . ctx - the context 267 268 Options Database Keys: 269 + -ts_monitor_lg_timestep - automatically sets line graph monitor 270 + -ts_monitor_lg_timestep_log - automatically sets line graph monitor 271 . -ts_monitor_lg_solution - monitor the solution (or certain values of the solution by calling `TSMonitorLGSetDisplayVariables()` or `TSMonitorLGCtxSetDisplayVariables()`) 272 . -ts_monitor_lg_error - monitor the error 273 . -ts_monitor_lg_ksp_iterations - monitor the number of `KSP` iterations needed for each timestep 274 . -ts_monitor_lg_snes_iterations - monitor the number of `SNES` iterations needed for each timestep 275 - -lg_use_markers <true,false> - mark the data points (at each time step) on the plot; default is true 276 277 Level: intermediate 278 279 Notes: 280 Pass the context and `TSMonitorLGCtxDestroy()` to `TSMonitorSet()` to have the context destroyed when no longer needed. 281 282 One can provide a function that transforms the solution before plotting it with `TSMonitorLGCtxSetTransform()` or `TSMonitorLGSetTransform()` 283 284 Many of the functions that control the monitoring have two forms: TSMonitorLGSet/GetXXXX() and TSMonitorLGCtxSet/GetXXXX() the first take a `TS` object as the 285 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 286 as the first argument. 287 288 One can control the names displayed for each solution or error variable with `TSMonitorLGCtxSetVariableNames()` or `TSMonitorLGSetVariableNames()` 289 290 .seealso: [](chapter_ts), `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()`, `TSMonitorDefault()`, `VecView()`, 291 `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`, 292 `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`, 293 `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`, 294 `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()` 295 @*/ 296 PetscErrorCode TSMonitorLGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtx *ctx) 297 { 298 PetscDraw draw; 299 300 PetscFunctionBegin; 301 PetscCall(PetscNew(ctx)); 302 PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw)); 303 PetscCall(PetscDrawSetFromOptions(draw)); 304 PetscCall(PetscDrawLGCreate(draw, 1, &(*ctx)->lg)); 305 PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg)); 306 PetscCall(PetscDrawDestroy(&draw)); 307 (*ctx)->howoften = howoften; 308 PetscFunctionReturn(PETSC_SUCCESS); 309 } 310 311 PetscErrorCode TSMonitorLGTimeStep(TS ts, PetscInt step, PetscReal ptime, Vec v, void *monctx) 312 { 313 TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx; 314 PetscReal x = ptime, y; 315 316 PetscFunctionBegin; 317 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates an interpolated solution */ 318 if (!step) { 319 PetscDrawAxis axis; 320 const char *ylabel = ctx->semilogy ? "Log Time Step" : "Time Step"; 321 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 322 PetscCall(PetscDrawAxisSetLabels(axis, "Timestep as function of time", "Time", ylabel)); 323 PetscCall(PetscDrawLGReset(ctx->lg)); 324 } 325 PetscCall(TSGetTimeStep(ts, &y)); 326 if (ctx->semilogy) y = PetscLog10Real(y); 327 PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y)); 328 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 329 PetscCall(PetscDrawLGDraw(ctx->lg)); 330 PetscCall(PetscDrawLGSave(ctx->lg)); 331 } 332 PetscFunctionReturn(PETSC_SUCCESS); 333 } 334 335 /*@C 336 TSMonitorLGCtxDestroy - Destroys a line graph context that was created with `TSMonitorLGCtxCreate()`. 337 338 Collective 339 340 Input Parameter: 341 . ctx - the monitor context 342 343 Level: intermediate 344 345 Note: 346 Pass to `TSMonitorSet()` along with the context and `TSMonitorLGTimeStep()` 347 348 .seealso: [](chapter_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep();` 349 @*/ 350 PetscErrorCode TSMonitorLGCtxDestroy(TSMonitorLGCtx *ctx) 351 { 352 PetscFunctionBegin; 353 if ((*ctx)->transformdestroy) PetscCall(((*ctx)->transformdestroy)((*ctx)->transformctx)); 354 PetscCall(PetscDrawLGDestroy(&(*ctx)->lg)); 355 PetscCall(PetscStrArrayDestroy(&(*ctx)->names)); 356 PetscCall(PetscStrArrayDestroy(&(*ctx)->displaynames)); 357 PetscCall(PetscFree((*ctx)->displayvariables)); 358 PetscCall(PetscFree((*ctx)->displayvalues)); 359 PetscCall(PetscFree(*ctx)); 360 PetscFunctionReturn(PETSC_SUCCESS); 361 } 362 363 /* Creates a TSMonitorSPCtx for use with DMSwarm particle visualizations */ 364 PetscErrorCode TSMonitorSPCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, PetscInt retain, PetscBool phase, PetscBool multispecies, TSMonitorSPCtx *ctx) 365 { 366 PetscDraw draw; 367 368 PetscFunctionBegin; 369 PetscCall(PetscNew(ctx)); 370 PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw)); 371 PetscCall(PetscDrawSetFromOptions(draw)); 372 PetscCall(PetscDrawSPCreate(draw, 1, &(*ctx)->sp)); 373 PetscCall(PetscDrawDestroy(&draw)); 374 (*ctx)->howoften = howoften; 375 (*ctx)->retain = retain; 376 (*ctx)->phase = phase; 377 (*ctx)->multispecies = multispecies; 378 PetscFunctionReturn(PETSC_SUCCESS); 379 } 380 381 /* Destroys a TSMonitorSPCtx that was created with TSMonitorSPCtxCreate */ 382 PetscErrorCode TSMonitorSPCtxDestroy(TSMonitorSPCtx *ctx) 383 { 384 PetscFunctionBegin; 385 386 PetscCall(PetscDrawSPDestroy(&(*ctx)->sp)); 387 PetscCall(PetscFree(*ctx)); 388 PetscFunctionReturn(PETSC_SUCCESS); 389 } 390 391 /* Creates a TSMonitorHGCtx for use with DMSwarm particle visualizations */ 392 PetscErrorCode TSMonitorHGCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, PetscInt Ns, PetscInt Nb, PetscBool velocity, TSMonitorHGCtx *ctx) 393 { 394 PetscDraw draw; 395 PetscInt s; 396 397 PetscFunctionBegin; 398 PetscCall(PetscNew(ctx)); 399 PetscCall(PetscMalloc1(Ns, &(*ctx)->hg)); 400 for (s = 0; s < Ns; ++s) { 401 PetscCall(PetscDrawCreate(comm, host, label, x + s * m, y, m, n, &draw)); 402 PetscCall(PetscDrawSetFromOptions(draw)); 403 PetscCall(PetscDrawHGCreate(draw, Nb, &(*ctx)->hg[s])); 404 PetscCall(PetscDrawHGCalcStats((*ctx)->hg[s], PETSC_TRUE)); 405 PetscCall(PetscDrawDestroy(&draw)); 406 } 407 (*ctx)->howoften = howoften; 408 (*ctx)->Ns = Ns; 409 (*ctx)->velocity = velocity; 410 PetscFunctionReturn(PETSC_SUCCESS); 411 } 412 413 /* Destroys a TSMonitorHGCtx that was created with TSMonitorHGCtxCreate */ 414 PetscErrorCode TSMonitorHGCtxDestroy(TSMonitorHGCtx *ctx) 415 { 416 PetscInt s; 417 418 PetscFunctionBegin; 419 for (s = 0; s < (*ctx)->Ns; ++s) PetscCall(PetscDrawHGDestroy(&(*ctx)->hg[s])); 420 PetscCall(PetscFree((*ctx)->hg)); 421 PetscCall(PetscFree(*ctx)); 422 PetscFunctionReturn(PETSC_SUCCESS); 423 } 424 425 /*@C 426 TSMonitorDrawSolution - Monitors progress of the `TS` solvers by calling 427 `VecView()` for the solution at each timestep 428 429 Collective 430 431 Input Parameters: 432 + ts - the `TS` context 433 . step - current time-step 434 . ptime - current time 435 - dummy - either a viewer or `NULL` 436 437 Options Database Keys: 438 + -ts_monitor_draw_solution - draw the solution at each time-step 439 - -ts_monitor_draw_solution_initial - show initial solution as well as current solution 440 441 Level: intermediate 442 443 Notes: 444 The initial solution and current solution are not displayed with a common axis scaling so generally the option `-ts_monitor_draw_solution_initial` 445 will look bad 446 447 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, as well as the context created with 448 `TSMonitorDrawCtxCreate()` and the function `TSMonitorDrawCtxDestroy()` to cause the monitor to be used during the `TS` integration. 449 450 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtxCreate()`, `TSMonitorDrawCtxDestroy()` 451 @*/ 452 PetscErrorCode TSMonitorDrawSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 453 { 454 TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy; 455 PetscDraw draw; 456 457 PetscFunctionBegin; 458 if (!step && ictx->showinitial) { 459 if (!ictx->initialsolution) PetscCall(VecDuplicate(u, &ictx->initialsolution)); 460 PetscCall(VecCopy(u, ictx->initialsolution)); 461 } 462 if (!(((ictx->howoften > 0) && (!(step % ictx->howoften))) || ((ictx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS); 463 464 if (ictx->showinitial) { 465 PetscReal pause; 466 PetscCall(PetscViewerDrawGetPause(ictx->viewer, &pause)); 467 PetscCall(PetscViewerDrawSetPause(ictx->viewer, 0.0)); 468 PetscCall(VecView(ictx->initialsolution, ictx->viewer)); 469 PetscCall(PetscViewerDrawSetPause(ictx->viewer, pause)); 470 PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_TRUE)); 471 } 472 PetscCall(VecView(u, ictx->viewer)); 473 if (ictx->showtimestepandtime) { 474 PetscReal xl, yl, xr, yr, h; 475 char time[32]; 476 477 PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw)); 478 PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime)); 479 PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr)); 480 h = yl + .95 * (yr - yl); 481 PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time)); 482 PetscCall(PetscDrawFlush(draw)); 483 } 484 485 if (ictx->showinitial) PetscCall(PetscViewerDrawSetHold(ictx->viewer, PETSC_FALSE)); 486 PetscFunctionReturn(PETSC_SUCCESS); 487 } 488 489 /*@C 490 TSMonitorDrawSolutionPhase - Monitors progress of the `TS` solvers by plotting the solution as a phase diagram 491 492 Collective 493 494 Input Parameters: 495 + ts - the `TS` context 496 . step - current time-step 497 . ptime - current time 498 - dummy - either a viewer or `NULL` 499 500 Level: intermediate 501 502 Notes: 503 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 504 to be used during the `TS` integration. 505 506 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()` 507 @*/ 508 PetscErrorCode TSMonitorDrawSolutionPhase(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 509 { 510 TSMonitorDrawCtx ictx = (TSMonitorDrawCtx)dummy; 511 PetscDraw draw; 512 PetscDrawAxis axis; 513 PetscInt n; 514 PetscMPIInt size; 515 PetscReal U0, U1, xl, yl, xr, yr, h; 516 char time[32]; 517 const PetscScalar *U; 518 519 PetscFunctionBegin; 520 PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ts), &size)); 521 PetscCheck(size == 1, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only allowed for sequential runs"); 522 PetscCall(VecGetSize(u, &n)); 523 PetscCheck(n == 2, PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "Only for ODEs with two unknowns"); 524 525 PetscCall(PetscViewerDrawGetDraw(ictx->viewer, 0, &draw)); 526 PetscCall(PetscViewerDrawGetDrawAxis(ictx->viewer, 0, &axis)); 527 PetscCall(PetscDrawAxisGetLimits(axis, &xl, &xr, &yl, &yr)); 528 if (!step) { 529 PetscCall(PetscDrawClear(draw)); 530 PetscCall(PetscDrawAxisDraw(axis)); 531 } 532 533 PetscCall(VecGetArrayRead(u, &U)); 534 U0 = PetscRealPart(U[0]); 535 U1 = PetscRealPart(U[1]); 536 PetscCall(VecRestoreArrayRead(u, &U)); 537 if ((U0 < xl) || (U1 < yl) || (U0 > xr) || (U1 > yr)) PetscFunctionReturn(PETSC_SUCCESS); 538 539 PetscDrawCollectiveBegin(draw); 540 PetscCall(PetscDrawPoint(draw, U0, U1, PETSC_DRAW_BLACK)); 541 if (ictx->showtimestepandtime) { 542 PetscCall(PetscDrawGetCoordinates(draw, &xl, &yl, &xr, &yr)); 543 PetscCall(PetscSNPrintf(time, 32, "Timestep %d Time %g", (int)step, (double)ptime)); 544 h = yl + .95 * (yr - yl); 545 PetscCall(PetscDrawStringCentered(draw, .5 * (xl + xr), h, PETSC_DRAW_BLACK, time)); 546 } 547 PetscDrawCollectiveEnd(draw); 548 PetscCall(PetscDrawFlush(draw)); 549 PetscCall(PetscDrawPause(draw)); 550 PetscCall(PetscDrawSave(draw)); 551 PetscFunctionReturn(PETSC_SUCCESS); 552 } 553 554 /*@C 555 TSMonitorDrawCtxDestroy - Destroys the monitor context for `TSMonitorDrawSolution()` 556 557 Collective 558 559 Input Parameters: 560 . ctx - the monitor context 561 562 Level: intermediate 563 564 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawSolution()`, `TSMonitorDrawError()`, `TSMonitorDrawCtx` 565 @*/ 566 PetscErrorCode TSMonitorDrawCtxDestroy(TSMonitorDrawCtx *ictx) 567 { 568 PetscFunctionBegin; 569 PetscCall(PetscViewerDestroy(&(*ictx)->viewer)); 570 PetscCall(VecDestroy(&(*ictx)->initialsolution)); 571 PetscCall(PetscFree(*ictx)); 572 PetscFunctionReturn(PETSC_SUCCESS); 573 } 574 575 /*@C 576 TSMonitorDrawCtxCreate - Creates the monitor context for `TSMonitorDrawCtx` 577 578 Collective 579 580 Input Parameter: 581 . ts - time-step context 582 583 Output Parameter: 584 . ctx - the monitor context 585 586 Options Database Keys: 587 + -ts_monitor_draw_solution - draw the solution at each time-step 588 - -ts_monitor_draw_solution_initial - show initial solution as well as current solution 589 590 Level: intermediate 591 592 Note: 593 The context created by this function, `PetscMonitorDrawSolution()`, and `TSMonitorDrawCtxDestroy()` should be passed together to `TSMonitorSet()`. 594 595 .seealso: [](chapter_ts), `TS`, `TSMonitorDrawCtxDestroy()`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorDrawCtx`, `PetscMonitorDrawSolution()` 596 @*/ 597 PetscErrorCode TSMonitorDrawCtxCreate(MPI_Comm comm, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorDrawCtx *ctx) 598 { 599 PetscFunctionBegin; 600 PetscCall(PetscNew(ctx)); 601 PetscCall(PetscViewerDrawOpen(comm, host, label, x, y, m, n, &(*ctx)->viewer)); 602 PetscCall(PetscViewerSetFromOptions((*ctx)->viewer)); 603 604 (*ctx)->howoften = howoften; 605 (*ctx)->showinitial = PETSC_FALSE; 606 PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_initial", &(*ctx)->showinitial, NULL)); 607 608 (*ctx)->showtimestepandtime = PETSC_FALSE; 609 PetscCall(PetscOptionsGetBool(NULL, NULL, "-ts_monitor_draw_solution_show_time", &(*ctx)->showtimestepandtime, NULL)); 610 PetscFunctionReturn(PETSC_SUCCESS); 611 } 612 613 /*@C 614 TSMonitorDrawSolutionFunction - Monitors progress of the `TS` solvers by calling 615 `VecView()` for the solution provided by `TSSetSolutionFunction()` at each timestep 616 617 Collective 618 619 Input Parameters: 620 + ts - the `TS` context 621 . step - current time-step 622 . ptime - current time 623 - dummy - either a viewer or `NULL` 624 625 Options Database Key: 626 . -ts_monitor_draw_solution_function - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()` 627 628 Level: intermediate 629 630 Note: 631 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 632 to be used during the `TS` integration. 633 634 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 635 @*/ 636 PetscErrorCode TSMonitorDrawSolutionFunction(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 637 { 638 TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy; 639 PetscViewer viewer = ctx->viewer; 640 Vec work; 641 642 PetscFunctionBegin; 643 if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS); 644 PetscCall(VecDuplicate(u, &work)); 645 PetscCall(TSComputeSolutionFunction(ts, ptime, work)); 646 PetscCall(VecView(work, viewer)); 647 PetscCall(VecDestroy(&work)); 648 PetscFunctionReturn(PETSC_SUCCESS); 649 } 650 651 /*@C 652 TSMonitorDrawError - Monitors progress of the `TS` solvers by calling 653 `VecView()` for the error at each timestep 654 655 Collective 656 657 Input Parameters: 658 + ts - the `TS` context 659 . step - current time-step 660 . ptime - current time 661 - dummy - either a viewer or `NULL` 662 663 Options Database Key: 664 . -ts_monitor_draw_error - Monitor error graphically, requires user to have provided `TSSetSolutionFunction()` 665 666 Level: intermediate 667 668 Notes: 669 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 670 to be used during the `TS` integration. 671 672 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 673 @*/ 674 PetscErrorCode TSMonitorDrawError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 675 { 676 TSMonitorDrawCtx ctx = (TSMonitorDrawCtx)dummy; 677 PetscViewer viewer = ctx->viewer; 678 Vec work; 679 680 PetscFunctionBegin; 681 if (!(((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason))) PetscFunctionReturn(PETSC_SUCCESS); 682 PetscCall(VecDuplicate(u, &work)); 683 PetscCall(TSComputeSolutionFunction(ts, ptime, work)); 684 PetscCall(VecAXPY(work, -1.0, u)); 685 PetscCall(VecView(work, viewer)); 686 PetscCall(VecDestroy(&work)); 687 PetscFunctionReturn(PETSC_SUCCESS); 688 } 689 690 /*@C 691 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 692 693 Collective 694 695 Input Parameters: 696 + ts - the `TS` context 697 . step - current time-step 698 . ptime - current time 699 . u - current state 700 - vf - viewer and its format 701 702 Level: intermediate 703 704 Notes: 705 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 706 to be used during the `TS` integration. 707 708 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()` 709 @*/ 710 PetscErrorCode TSMonitorSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf) 711 { 712 PetscFunctionBegin; 713 if (vf->view_interval > 0 && !ts->reason && step % vf->view_interval != 0) PetscFunctionReturn(PETSC_SUCCESS); 714 PetscCall(PetscViewerPushFormat(vf->viewer, vf->format)); 715 PetscCall(VecView(u, vf->viewer)); 716 PetscCall(PetscViewerPopFormat(vf->viewer)); 717 PetscFunctionReturn(PETSC_SUCCESS); 718 } 719 720 /*@C 721 TSMonitorSolutionVTK - Monitors progress of the `TS` solvers by `VecView()` for the solution at each timestep. 722 723 Collective 724 725 Input Parameters: 726 + ts - the `TS` context 727 . step - current time-step 728 . ptime - current time 729 . u - current state 730 - filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ") 731 732 Level: intermediate 733 734 Notes: 735 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. 736 These are named according to the file name template. 737 738 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 739 to be used during the `TS` integration. 740 741 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()` 742 @*/ 743 PetscErrorCode TSMonitorSolutionVTK(TS ts, PetscInt step, PetscReal ptime, Vec u, void *filenametemplate) 744 { 745 char filename[PETSC_MAX_PATH_LEN]; 746 PetscViewer viewer; 747 748 PetscFunctionBegin; 749 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 750 PetscCall(PetscSNPrintf(filename, sizeof(filename), (const char *)filenametemplate, step)); 751 PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)ts), filename, FILE_MODE_WRITE, &viewer)); 752 PetscCall(VecView(u, viewer)); 753 PetscCall(PetscViewerDestroy(&viewer)); 754 PetscFunctionReturn(PETSC_SUCCESS); 755 } 756 757 /*@C 758 TSMonitorSolutionVTKDestroy - Destroy filename template string created for use with `TSMonitorSolutionVTK()` 759 760 Not Collective 761 762 Input Parameters: 763 . filenametemplate - string containing a format specifier for the integer time step (e.g. %03" PetscInt_FMT ") 764 765 Level: intermediate 766 767 Note: 768 This function is normally passed to `TSMonitorSet()` along with `TSMonitorSolutionVTK()`. 769 770 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorSolutionVTK()` 771 @*/ 772 PetscErrorCode TSMonitorSolutionVTKDestroy(void *filenametemplate) 773 { 774 PetscFunctionBegin; 775 PetscCall(PetscFree(*(char **)filenametemplate)); 776 PetscFunctionReturn(PETSC_SUCCESS); 777 } 778 779 /*@C 780 TSMonitorLGSolution - Monitors progress of the `TS` solvers by plotting each component of the solution vector 781 in a time based line graph 782 783 Collective 784 785 Input Parameters: 786 + ts - the `TS` context 787 . step - current time-step 788 . ptime - current time 789 . u - current solution 790 - dctx - the `TSMonitorLGCtx` object that contains all the options for the monitoring, this is created with `TSMonitorLGCtxCreate()` 791 792 Options Database Key: 793 . -ts_monitor_lg_solution_variables - enable monitor of lg solution variables 794 795 Level: intermediate 796 797 Notes: 798 Each process in a parallel run displays its component solutions in a separate window 799 800 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 801 to be used during the `TS` integration. 802 803 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGCtxCreate()`, `TSMonitorLGCtxSetVariableNames()`, `TSMonitorLGCtxGetVariableNames()`, 804 `TSMonitorLGSetVariableNames()`, `TSMonitorLGGetVariableNames()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetDisplayVariables()`, 805 `TSMonitorLGCtxSetTransform()`, `TSMonitorLGSetTransform()`, `TSMonitorLGError()`, `TSMonitorLGSNESIterations()`, `TSMonitorLGKSPIterations()`, 806 `TSMonitorEnvelopeCtxCreate()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxDestroy()`, `TSMonitorEnvelop()` 807 @*/ 808 PetscErrorCode TSMonitorLGSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx) 809 { 810 TSMonitorLGCtx ctx = (TSMonitorLGCtx)dctx; 811 const PetscScalar *yy; 812 Vec v; 813 814 PetscFunctionBegin; 815 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 816 if (!step) { 817 PetscDrawAxis axis; 818 PetscInt dim; 819 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 820 PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution")); 821 if (!ctx->names) { 822 PetscBool flg; 823 /* user provides names of variables to plot but no names has been set so assume names are integer values */ 824 PetscCall(PetscOptionsHasName(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", &flg)); 825 if (flg) { 826 PetscInt i, n; 827 char **names; 828 PetscCall(VecGetSize(u, &n)); 829 PetscCall(PetscMalloc1(n + 1, &names)); 830 for (i = 0; i < n; i++) { 831 PetscCall(PetscMalloc1(5, &names[i])); 832 PetscCall(PetscSNPrintf(names[i], 5, "%" PetscInt_FMT, i)); 833 } 834 names[n] = NULL; 835 ctx->names = names; 836 } 837 } 838 if (ctx->names && !ctx->displaynames) { 839 char **displaynames; 840 PetscBool flg; 841 PetscCall(VecGetLocalSize(u, &dim)); 842 PetscCall(PetscCalloc1(dim + 1, &displaynames)); 843 PetscCall(PetscOptionsGetStringArray(((PetscObject)ts)->options, ((PetscObject)ts)->prefix, "-ts_monitor_lg_solution_variables", displaynames, &dim, &flg)); 844 if (flg) PetscCall(TSMonitorLGCtxSetDisplayVariables(ctx, (const char *const *)displaynames)); 845 PetscCall(PetscStrArrayDestroy(&displaynames)); 846 } 847 if (ctx->displaynames) { 848 PetscCall(PetscDrawLGSetDimension(ctx->lg, ctx->ndisplayvariables)); 849 PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->displaynames)); 850 } else if (ctx->names) { 851 PetscCall(VecGetLocalSize(u, &dim)); 852 PetscCall(PetscDrawLGSetDimension(ctx->lg, dim)); 853 PetscCall(PetscDrawLGSetLegend(ctx->lg, (const char *const *)ctx->names)); 854 } else { 855 PetscCall(VecGetLocalSize(u, &dim)); 856 PetscCall(PetscDrawLGSetDimension(ctx->lg, dim)); 857 } 858 PetscCall(PetscDrawLGReset(ctx->lg)); 859 } 860 861 if (!ctx->transform) v = u; 862 else PetscCall((*ctx->transform)(ctx->transformctx, u, &v)); 863 PetscCall(VecGetArrayRead(v, &yy)); 864 if (ctx->displaynames) { 865 PetscInt i; 866 for (i = 0; i < ctx->ndisplayvariables; i++) ctx->displayvalues[i] = PetscRealPart(yy[ctx->displayvariables[i]]); 867 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, ctx->displayvalues)); 868 } else { 869 #if defined(PETSC_USE_COMPLEX) 870 PetscInt i, n; 871 PetscReal *yreal; 872 PetscCall(VecGetLocalSize(v, &n)); 873 PetscCall(PetscMalloc1(n, &yreal)); 874 for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]); 875 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal)); 876 PetscCall(PetscFree(yreal)); 877 #else 878 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy)); 879 #endif 880 } 881 PetscCall(VecRestoreArrayRead(v, &yy)); 882 if (ctx->transform) PetscCall(VecDestroy(&v)); 883 884 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 885 PetscCall(PetscDrawLGDraw(ctx->lg)); 886 PetscCall(PetscDrawLGSave(ctx->lg)); 887 } 888 PetscFunctionReturn(PETSC_SUCCESS); 889 } 890 891 /*@C 892 TSMonitorLGSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot 893 894 Collective 895 896 Input Parameters: 897 + ts - the `TS` context 898 - names - the names of the components, final string must be `NULL` 899 900 Level: intermediate 901 902 Notes: 903 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 904 905 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGCtxSetVariableNames()` 906 @*/ 907 PetscErrorCode TSMonitorLGSetVariableNames(TS ts, const char *const *names) 908 { 909 PetscInt i; 910 911 PetscFunctionBegin; 912 for (i = 0; i < ts->numbermonitors; i++) { 913 if (ts->monitor[i] == TSMonitorLGSolution) { 914 PetscCall(TSMonitorLGCtxSetVariableNames((TSMonitorLGCtx)ts->monitorcontext[i], names)); 915 break; 916 } 917 } 918 PetscFunctionReturn(PETSC_SUCCESS); 919 } 920 921 /*@C 922 TSMonitorLGCtxSetVariableNames - Sets the name of each component in the solution vector so that it may be displayed in the plot 923 924 Collective 925 926 Input Parameters: 927 + ts - the `TS` context 928 - names - the names of the components, final string must be `NULL` 929 930 Level: intermediate 931 932 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()`, `TSMonitorLGSetVariableNames()` 933 @*/ 934 PetscErrorCode TSMonitorLGCtxSetVariableNames(TSMonitorLGCtx ctx, const char *const *names) 935 { 936 PetscFunctionBegin; 937 PetscCall(PetscStrArrayDestroy(&ctx->names)); 938 PetscCall(PetscStrArrayallocpy(names, &ctx->names)); 939 PetscFunctionReturn(PETSC_SUCCESS); 940 } 941 942 /*@C 943 TSMonitorLGGetVariableNames - Gets the name of each component in the solution vector so that it may be displayed in the plot 944 945 Collective 946 947 Input Parameter: 948 . ts - the `TS` context 949 950 Output Parameter: 951 . names - the names of the components, final string must be `NULL` 952 953 Level: intermediate 954 955 Note: 956 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 957 958 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()` 959 @*/ 960 PetscErrorCode TSMonitorLGGetVariableNames(TS ts, const char *const **names) 961 { 962 PetscInt i; 963 964 PetscFunctionBegin; 965 *names = NULL; 966 for (i = 0; i < ts->numbermonitors; i++) { 967 if (ts->monitor[i] == TSMonitorLGSolution) { 968 TSMonitorLGCtx ctx = (TSMonitorLGCtx)ts->monitorcontext[i]; 969 *names = (const char *const *)ctx->names; 970 break; 971 } 972 } 973 PetscFunctionReturn(PETSC_SUCCESS); 974 } 975 976 /*@C 977 TSMonitorLGCtxSetDisplayVariables - Sets the variables that are to be display in the monitor 978 979 Collective 980 981 Input Parameters: 982 + ctx - the `TSMonitorLG` context 983 - displaynames - the names of the components, final string must be `NULL` 984 985 Level: intermediate 986 987 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()` 988 @*/ 989 PetscErrorCode TSMonitorLGCtxSetDisplayVariables(TSMonitorLGCtx ctx, const char *const *displaynames) 990 { 991 PetscInt j = 0, k; 992 993 PetscFunctionBegin; 994 if (!ctx->names) PetscFunctionReturn(PETSC_SUCCESS); 995 PetscCall(PetscStrArrayDestroy(&ctx->displaynames)); 996 PetscCall(PetscStrArrayallocpy(displaynames, &ctx->displaynames)); 997 while (displaynames[j]) j++; 998 ctx->ndisplayvariables = j; 999 PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvariables)); 1000 PetscCall(PetscMalloc1(ctx->ndisplayvariables, &ctx->displayvalues)); 1001 j = 0; 1002 while (displaynames[j]) { 1003 k = 0; 1004 while (ctx->names[k]) { 1005 PetscBool flg; 1006 PetscCall(PetscStrcmp(displaynames[j], ctx->names[k], &flg)); 1007 if (flg) { 1008 ctx->displayvariables[j] = k; 1009 break; 1010 } 1011 k++; 1012 } 1013 j++; 1014 } 1015 PetscFunctionReturn(PETSC_SUCCESS); 1016 } 1017 1018 /*@C 1019 TSMonitorLGSetDisplayVariables - Sets the variables that are to be display in the monitor 1020 1021 Collective 1022 1023 Input Parameters: 1024 + ts - the `TS` context 1025 - displaynames - the names of the components, final string must be `NULL` 1026 1027 Level: intermediate 1028 1029 Note: 1030 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 1031 1032 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()` 1033 @*/ 1034 PetscErrorCode TSMonitorLGSetDisplayVariables(TS ts, const char *const *displaynames) 1035 { 1036 PetscInt i; 1037 1038 PetscFunctionBegin; 1039 for (i = 0; i < ts->numbermonitors; i++) { 1040 if (ts->monitor[i] == TSMonitorLGSolution) { 1041 PetscCall(TSMonitorLGCtxSetDisplayVariables((TSMonitorLGCtx)ts->monitorcontext[i], displaynames)); 1042 break; 1043 } 1044 } 1045 PetscFunctionReturn(PETSC_SUCCESS); 1046 } 1047 1048 /*@C 1049 TSMonitorLGSetTransform - Solution vector will be transformed by provided function before being displayed 1050 1051 Collective 1052 1053 Input Parameters: 1054 + ts - the `TS` context 1055 . transform - the transform function 1056 . destroy - function to destroy the optional context 1057 - ctx - optional context used by transform function 1058 1059 Level: intermediate 1060 1061 Note: 1062 If the `TS` object does not have a `TSMonitorLGCtx` associated with it then this function is ignored 1063 1064 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGCtxSetTransform()` 1065 @*/ 1066 PetscErrorCode TSMonitorLGSetTransform(TS ts, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx) 1067 { 1068 PetscInt i; 1069 1070 PetscFunctionBegin; 1071 for (i = 0; i < ts->numbermonitors; i++) { 1072 if (ts->monitor[i] == TSMonitorLGSolution) PetscCall(TSMonitorLGCtxSetTransform((TSMonitorLGCtx)ts->monitorcontext[i], transform, destroy, tctx)); 1073 } 1074 PetscFunctionReturn(PETSC_SUCCESS); 1075 } 1076 1077 /*@C 1078 TSMonitorLGCtxSetTransform - Solution vector will be transformed by provided function before being displayed 1079 1080 Collective 1081 1082 Input Parameters: 1083 + ts - the `TS` context 1084 . transform - the transform function 1085 . destroy - function to destroy the optional context 1086 - ctx - optional context used by transform function 1087 1088 Level: intermediate 1089 1090 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetVariableNames()`, `TSMonitorLGSetTransform()` 1091 @*/ 1092 PetscErrorCode TSMonitorLGCtxSetTransform(TSMonitorLGCtx ctx, PetscErrorCode (*transform)(void *, Vec, Vec *), PetscErrorCode (*destroy)(void *), void *tctx) 1093 { 1094 PetscFunctionBegin; 1095 ctx->transform = transform; 1096 ctx->transformdestroy = destroy; 1097 ctx->transformctx = tctx; 1098 PetscFunctionReturn(PETSC_SUCCESS); 1099 } 1100 1101 /*@C 1102 TSMonitorLGError - Monitors progress of the `TS` solvers by plotting each component of the error 1103 in a time based line graph 1104 1105 Collective 1106 1107 Input Parameters: 1108 + ts - the `TS` context 1109 . step - current time-step 1110 . ptime - current time 1111 . u - current solution 1112 - dctx - `TSMonitorLGCtx` object created with `TSMonitorLGCtxCreate()` 1113 1114 Options Database Key: 1115 . -ts_monitor_lg_error - create a graphical monitor of error history 1116 1117 Level: intermediate 1118 1119 Notes: 1120 Each process in a parallel run displays its component errors in a separate window 1121 1122 The user must provide the solution using `TSSetSolutionFunction()` to use this monitor. 1123 1124 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1125 to be used during the TS integration. 1126 1127 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 1128 @*/ 1129 PetscErrorCode TSMonitorLGError(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dummy) 1130 { 1131 TSMonitorLGCtx ctx = (TSMonitorLGCtx)dummy; 1132 const PetscScalar *yy; 1133 Vec y; 1134 1135 PetscFunctionBegin; 1136 if (!step) { 1137 PetscDrawAxis axis; 1138 PetscInt dim; 1139 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 1140 PetscCall(PetscDrawAxisSetLabels(axis, "Error in solution as function of time", "Time", "Error")); 1141 PetscCall(VecGetLocalSize(u, &dim)); 1142 PetscCall(PetscDrawLGSetDimension(ctx->lg, dim)); 1143 PetscCall(PetscDrawLGReset(ctx->lg)); 1144 } 1145 PetscCall(VecDuplicate(u, &y)); 1146 PetscCall(TSComputeSolutionFunction(ts, ptime, y)); 1147 PetscCall(VecAXPY(y, -1.0, u)); 1148 PetscCall(VecGetArrayRead(y, &yy)); 1149 #if defined(PETSC_USE_COMPLEX) 1150 { 1151 PetscReal *yreal; 1152 PetscInt i, n; 1153 PetscCall(VecGetLocalSize(y, &n)); 1154 PetscCall(PetscMalloc1(n, &yreal)); 1155 for (i = 0; i < n; i++) yreal[i] = PetscRealPart(yy[i]); 1156 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yreal)); 1157 PetscCall(PetscFree(yreal)); 1158 } 1159 #else 1160 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg, ptime, yy)); 1161 #endif 1162 PetscCall(VecRestoreArrayRead(y, &yy)); 1163 PetscCall(VecDestroy(&y)); 1164 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 1165 PetscCall(PetscDrawLGDraw(ctx->lg)); 1166 PetscCall(PetscDrawLGSave(ctx->lg)); 1167 } 1168 PetscFunctionReturn(PETSC_SUCCESS); 1169 } 1170 1171 /*@C 1172 TSMonitorSPSwarmSolution - Graphically displays phase plots of `DMSWARM` particles on a scatter plot 1173 1174 Input Parameters: 1175 + ts - the `TS` context 1176 . step - current time-step 1177 . ptime - current time 1178 . u - current solution 1179 - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorSPCtxCreate()` 1180 1181 Options Database Keys: 1182 + -ts_monitor_sp_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution 1183 . -ts_monitor_sp_swarm_retain <n> - Retain n old points so we can see the history, or -1 for all points 1184 . -ts_monitor_sp_swarm_multi_species <bool> - Color each species differently 1185 - -ts_monitor_sp_swarm_phase <bool> - Plot in phase space, as opposed to coordinate space 1186 1187 Level: intermediate 1188 1189 Notes: 1190 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1191 to be used during the `TS` integration. 1192 1193 .seealso: [](chapter_ts), `TS`, `TSMonitoSet()`, `DMSWARM`, `TSMonitorSPCtxCreate()` 1194 @*/ 1195 PetscErrorCode TSMonitorSPSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx) 1196 { 1197 TSMonitorSPCtx ctx = (TSMonitorSPCtx)dctx; 1198 PetscDraw draw; 1199 DM dm, cdm; 1200 const PetscScalar *yy; 1201 PetscInt Np, p, dim = 2, *species; 1202 PetscReal species_color; 1203 1204 PetscFunctionBegin; 1205 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 1206 PetscCall(TSGetDM(ts, &dm)); 1207 if (!step) { 1208 PetscDrawAxis axis; 1209 PetscReal dmboxlower[2], dmboxupper[2]; 1210 1211 PetscCall(TSGetDM(ts, &dm)); 1212 PetscCall(DMGetDimension(dm, &dim)); 1213 PetscCheck(dim == 2, PETSC_COMM_SELF, PETSC_ERR_SUP, "Monitor only supports two dimensional fields"); 1214 PetscCall(DMSwarmGetCellDM(dm, &cdm)); 1215 PetscCall(DMGetBoundingBox(cdm, dmboxlower, dmboxupper)); 1216 PetscCall(VecGetLocalSize(u, &Np)); 1217 Np /= dim * 2; 1218 PetscCall(PetscDrawSPGetAxis(ctx->sp, &axis)); 1219 if (ctx->phase) { 1220 PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "V")); 1221 PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], -10, 10)); 1222 } else { 1223 PetscCall(PetscDrawAxisSetLabels(axis, "Particles", "X", "Y")); 1224 PetscCall(PetscDrawAxisSetLimits(axis, dmboxlower[0], dmboxupper[0], dmboxlower[1], dmboxupper[1])); 1225 } 1226 PetscCall(PetscDrawAxisSetHoldLimits(axis, PETSC_TRUE)); 1227 PetscCall(PetscDrawSPReset(ctx->sp)); 1228 } 1229 if (ctx->multispecies) PetscCall(DMSwarmGetField(dm, "species", NULL, NULL, (void **)&species)); 1230 PetscCall(VecGetLocalSize(u, &Np)); 1231 Np /= dim * 2; 1232 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 1233 PetscCall(PetscDrawSPGetDraw(ctx->sp, &draw)); 1234 if ((ctx->retain == 0) || (ctx->retain > 0 && !(step % ctx->retain))) PetscCall(PetscDrawClear(draw)); 1235 PetscCall(PetscDrawFlush(draw)); 1236 PetscCall(PetscDrawSPReset(ctx->sp)); 1237 PetscCall(VecGetArrayRead(u, &yy)); 1238 for (p = 0; p < Np; ++p) { 1239 PetscReal x, y; 1240 1241 if (ctx->phase) { 1242 x = PetscRealPart(yy[p * dim * 2]); 1243 y = PetscRealPart(yy[p * dim * 2 + dim]); 1244 } else { 1245 x = PetscRealPart(yy[p * dim * 2]); 1246 y = PetscRealPart(yy[p * dim * 2 + 1]); 1247 } 1248 if (ctx->multispecies) { 1249 species_color = species[p] + 2; 1250 PetscCall(PetscDrawSPAddPointColorized(ctx->sp, &x, &y, &species_color)); 1251 } else { 1252 PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y)); 1253 } 1254 PetscCall(PetscDrawSPAddPoint(ctx->sp, &x, &y)); 1255 } 1256 PetscCall(VecRestoreArrayRead(u, &yy)); 1257 PetscCall(PetscDrawSPDraw(ctx->sp, PETSC_FALSE)); 1258 PetscCall(PetscDrawSPSave(ctx->sp)); 1259 if (ctx->multispecies) PetscCall(DMSwarmRestoreField(dm, "species", NULL, NULL, (void **)&species)); 1260 } 1261 PetscFunctionReturn(PETSC_SUCCESS); 1262 } 1263 1264 /*@C 1265 TSMonitorHGSwarmSolution - Graphically displays histograms of `DMSWARM` particles 1266 1267 Input Parameters: 1268 + ts - the `TS` context 1269 . step - current time-step 1270 . ptime - current time 1271 . u - current solution 1272 - dctx - the `TSMonitorSPCtx` object that contains all the options for the monitoring, this is created with `TSMonitorHGCtxCreate()` 1273 1274 Options Database Keys: 1275 + -ts_monitor_hg_swarm <n> - Monitor the solution every n steps, or -1 for plotting only the final solution 1276 . -ts_monitor_hg_swarm_species <num> - Number of species to histogram 1277 . -ts_monitor_hg_swarm_bins <num> - Number of histogram bins 1278 - -ts_monitor_hg_swarm_velocity <bool> - Plot in velocity space, as opposed to coordinate space 1279 1280 Level: intermediate 1281 1282 Note: 1283 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1284 to be used during the `TS` integration. 1285 1286 .seealso: `TSMonitoSet()` 1287 @*/ 1288 PetscErrorCode TSMonitorHGSwarmSolution(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx) 1289 { 1290 TSMonitorHGCtx ctx = (TSMonitorHGCtx)dctx; 1291 PetscDraw draw; 1292 DM sw; 1293 const PetscScalar *yy; 1294 PetscInt *species; 1295 PetscInt dim, d = 0, Np, p, Ns, s; 1296 1297 PetscFunctionBegin; 1298 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 1299 PetscCall(TSGetDM(ts, &sw)); 1300 PetscCall(DMGetDimension(sw, &dim)); 1301 PetscCall(DMSwarmGetNumSpecies(sw, &Ns)); 1302 Ns = PetscMin(Ns, ctx->Ns); 1303 PetscCall(VecGetLocalSize(u, &Np)); 1304 Np /= dim * 2; 1305 if (!step) { 1306 PetscDrawAxis axis; 1307 char title[PETSC_MAX_PATH_LEN]; 1308 1309 for (s = 0; s < Ns; ++s) { 1310 PetscCall(PetscDrawHGGetAxis(ctx->hg[s], &axis)); 1311 PetscCall(PetscSNPrintf(title, PETSC_MAX_PATH_LEN, "Species %" PetscInt_FMT, s)); 1312 if (ctx->velocity) PetscCall(PetscDrawAxisSetLabels(axis, title, "V", "N")); 1313 else PetscCall(PetscDrawAxisSetLabels(axis, title, "X", "N")); 1314 } 1315 } 1316 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && ts->reason)) { 1317 PetscCall(DMSwarmGetField(sw, "species", NULL, NULL, (void **)&species)); 1318 for (s = 0; s < Ns; ++s) { 1319 PetscCall(PetscDrawHGReset(ctx->hg[s])); 1320 PetscCall(PetscDrawHGGetDraw(ctx->hg[s], &draw)); 1321 PetscCall(PetscDrawClear(draw)); 1322 PetscCall(PetscDrawFlush(draw)); 1323 } 1324 PetscCall(VecGetArrayRead(u, &yy)); 1325 for (p = 0; p < Np; ++p) { 1326 const PetscInt s = species[p] < Ns ? species[p] : 0; 1327 PetscReal v; 1328 1329 if (ctx->velocity) v = PetscRealPart(yy[p * dim * 2 + d + dim]); 1330 else v = PetscRealPart(yy[p * dim * 2 + d]); 1331 PetscCall(PetscDrawHGAddValue(ctx->hg[s], v)); 1332 } 1333 PetscCall(VecRestoreArrayRead(u, &yy)); 1334 for (s = 0; s < Ns; ++s) { 1335 PetscCall(PetscDrawHGDraw(ctx->hg[s])); 1336 PetscCall(PetscDrawHGSave(ctx->hg[s])); 1337 } 1338 PetscCall(DMSwarmRestoreField(sw, "species", NULL, NULL, (void **)&species)); 1339 } 1340 PetscFunctionReturn(PETSC_SUCCESS); 1341 } 1342 1343 /*@C 1344 TSMonitorError - Monitors progress of the `TS` solvers by printing the 2 norm of the error at each timestep 1345 1346 Collective 1347 1348 Input Parameters: 1349 + ts - the `TS` context 1350 . step - current time-step 1351 . ptime - current time 1352 . u - current solution 1353 - dctx - unused context 1354 1355 Options Database Key: 1356 . -ts_monitor_error - create a graphical monitor of error history 1357 1358 Level: intermediate 1359 1360 Notes: 1361 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1362 to be used during the `TS` integration. 1363 1364 The user must provide the solution using `TSSetSolutionFunction()` to use this monitor. 1365 1366 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSSetSolutionFunction()` 1367 @*/ 1368 PetscErrorCode TSMonitorError(TS ts, PetscInt step, PetscReal ptime, Vec u, PetscViewerAndFormat *vf) 1369 { 1370 DM dm; 1371 PetscDS ds = NULL; 1372 PetscInt Nf = -1, f; 1373 PetscBool flg; 1374 1375 PetscFunctionBegin; 1376 PetscCall(TSGetDM(ts, &dm)); 1377 if (dm) PetscCall(DMGetDS(dm, &ds)); 1378 if (ds) PetscCall(PetscDSGetNumFields(ds, &Nf)); 1379 if (Nf <= 0) { 1380 Vec y; 1381 PetscReal nrm; 1382 1383 PetscCall(VecDuplicate(u, &y)); 1384 PetscCall(TSComputeSolutionFunction(ts, ptime, y)); 1385 PetscCall(VecAXPY(y, -1.0, u)); 1386 PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERASCII, &flg)); 1387 if (flg) { 1388 PetscCall(VecNorm(y, NORM_2, &nrm)); 1389 PetscCall(PetscViewerASCIIPrintf(vf->viewer, "2-norm of error %g\n", (double)nrm)); 1390 } 1391 PetscCall(PetscObjectTypeCompare((PetscObject)vf->viewer, PETSCVIEWERDRAW, &flg)); 1392 if (flg) PetscCall(VecView(y, vf->viewer)); 1393 PetscCall(VecDestroy(&y)); 1394 } else { 1395 PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx); 1396 void **ctxs; 1397 Vec v; 1398 PetscReal ferrors[1]; 1399 1400 PetscCall(PetscMalloc2(Nf, &exactFuncs, Nf, &ctxs)); 1401 for (f = 0; f < Nf; ++f) PetscCall(PetscDSGetExactSolution(ds, f, &exactFuncs[f], &ctxs[f])); 1402 PetscCall(DMComputeL2FieldDiff(dm, ptime, exactFuncs, ctxs, u, ferrors)); 1403 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Timestep: %04d time = %-8.4g \t L_2 Error: [", (int)step, (double)ptime)); 1404 for (f = 0; f < Nf; ++f) { 1405 if (f > 0) PetscCall(PetscPrintf(PETSC_COMM_WORLD, ", ")); 1406 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%2.3g", (double)ferrors[f])); 1407 } 1408 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "]\n")); 1409 1410 PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view")); 1411 1412 PetscCall(PetscOptionsHasName(NULL, NULL, "-exact_vec_view", &flg)); 1413 if (flg) { 1414 PetscCall(DMGetGlobalVector(dm, &v)); 1415 PetscCall(DMProjectFunction(dm, ptime, exactFuncs, ctxs, INSERT_ALL_VALUES, v)); 1416 PetscCall(PetscObjectSetName((PetscObject)v, "Exact Solution")); 1417 PetscCall(VecViewFromOptions(v, NULL, "-exact_vec_view")); 1418 PetscCall(DMRestoreGlobalVector(dm, &v)); 1419 } 1420 PetscCall(PetscFree2(exactFuncs, ctxs)); 1421 } 1422 PetscFunctionReturn(PETSC_SUCCESS); 1423 } 1424 1425 PetscErrorCode TSMonitorLGSNESIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx) 1426 { 1427 TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx; 1428 PetscReal x = ptime, y; 1429 PetscInt its; 1430 1431 PetscFunctionBegin; 1432 if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 1433 if (!n) { 1434 PetscDrawAxis axis; 1435 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 1436 PetscCall(PetscDrawAxisSetLabels(axis, "Nonlinear iterations as function of time", "Time", "SNES Iterations")); 1437 PetscCall(PetscDrawLGReset(ctx->lg)); 1438 ctx->snes_its = 0; 1439 } 1440 PetscCall(TSGetSNESIterations(ts, &its)); 1441 y = its - ctx->snes_its; 1442 PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y)); 1443 if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) { 1444 PetscCall(PetscDrawLGDraw(ctx->lg)); 1445 PetscCall(PetscDrawLGSave(ctx->lg)); 1446 } 1447 ctx->snes_its = its; 1448 PetscFunctionReturn(PETSC_SUCCESS); 1449 } 1450 1451 PetscErrorCode TSMonitorLGKSPIterations(TS ts, PetscInt n, PetscReal ptime, Vec v, void *monctx) 1452 { 1453 TSMonitorLGCtx ctx = (TSMonitorLGCtx)monctx; 1454 PetscReal x = ptime, y; 1455 PetscInt its; 1456 1457 PetscFunctionBegin; 1458 if (n < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 1459 if (!n) { 1460 PetscDrawAxis axis; 1461 PetscCall(PetscDrawLGGetAxis(ctx->lg, &axis)); 1462 PetscCall(PetscDrawAxisSetLabels(axis, "Linear iterations as function of time", "Time", "KSP Iterations")); 1463 PetscCall(PetscDrawLGReset(ctx->lg)); 1464 ctx->ksp_its = 0; 1465 } 1466 PetscCall(TSGetKSPIterations(ts, &its)); 1467 y = its - ctx->ksp_its; 1468 PetscCall(PetscDrawLGAddPoint(ctx->lg, &x, &y)); 1469 if (((ctx->howoften > 0) && (!(n % ctx->howoften)) && (n > -1)) || ((ctx->howoften == -1) && (n == -1))) { 1470 PetscCall(PetscDrawLGDraw(ctx->lg)); 1471 PetscCall(PetscDrawLGSave(ctx->lg)); 1472 } 1473 ctx->ksp_its = its; 1474 PetscFunctionReturn(PETSC_SUCCESS); 1475 } 1476 1477 /*@C 1478 TSMonitorEnvelopeCtxCreate - Creates a context for use with `TSMonitorEnvelope()` 1479 1480 Collective 1481 1482 Input Parameters: 1483 . ts - the `TS` solver object 1484 1485 Output Parameter: 1486 . ctx - the context 1487 1488 Level: intermediate 1489 1490 .seealso: [](chapter_ts), `TS`, `TSMonitorLGTimeStep()`, `TSMonitorSet()`, `TSMonitorLGSolution()`, `TSMonitorLGError()` 1491 @*/ 1492 PetscErrorCode TSMonitorEnvelopeCtxCreate(TS ts, TSMonitorEnvelopeCtx *ctx) 1493 { 1494 PetscFunctionBegin; 1495 PetscCall(PetscNew(ctx)); 1496 PetscFunctionReturn(PETSC_SUCCESS); 1497 } 1498 1499 /*@C 1500 TSMonitorEnvelope - Monitors the maximum and minimum value of each component of the solution 1501 1502 Collective 1503 1504 Input Parameters: 1505 + ts - the `TS` context 1506 . step - current time-step 1507 . ptime - current time 1508 . u - current solution 1509 - dctx - the envelope context 1510 1511 Options Database Key: 1512 . -ts_monitor_envelope - determine maximum and minimum value of each component of the solution over the solution time 1513 1514 Level: intermediate 1515 1516 Notes: 1517 After a solve you can use `TSMonitorEnvelopeGetBounds()` to access the envelope 1518 1519 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1520 to be used during the `TS` integration. 1521 1522 .seealso: [](chapter_ts), `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorEnvelopeGetBounds()`, `TSMonitorEnvelopeCtxCreate()` 1523 @*/ 1524 PetscErrorCode TSMonitorEnvelope(TS ts, PetscInt step, PetscReal ptime, Vec u, void *dctx) 1525 { 1526 TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)dctx; 1527 1528 PetscFunctionBegin; 1529 if (!ctx->max) { 1530 PetscCall(VecDuplicate(u, &ctx->max)); 1531 PetscCall(VecDuplicate(u, &ctx->min)); 1532 PetscCall(VecCopy(u, ctx->max)); 1533 PetscCall(VecCopy(u, ctx->min)); 1534 } else { 1535 PetscCall(VecPointwiseMax(ctx->max, u, ctx->max)); 1536 PetscCall(VecPointwiseMin(ctx->min, u, ctx->min)); 1537 } 1538 PetscFunctionReturn(PETSC_SUCCESS); 1539 } 1540 1541 /*@C 1542 TSMonitorEnvelopeGetBounds - Gets the bounds for the components of the solution 1543 1544 Collective 1545 1546 Input Parameter: 1547 . ts - the `TS` context 1548 1549 Output Parameters: 1550 + max - the maximum values 1551 - min - the minimum values 1552 1553 Level: intermediate 1554 1555 Notes: 1556 If the `TS` does not have a `TSMonitorEnvelopeCtx` associated with it then this function is ignored 1557 1558 .seealso: [](chapter_ts), `TSMonitorEnvelopeCtx`, `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `VecView()`, `TSMonitorLGSetDisplayVariables()` 1559 @*/ 1560 PetscErrorCode TSMonitorEnvelopeGetBounds(TS ts, Vec *max, Vec *min) 1561 { 1562 PetscInt i; 1563 1564 PetscFunctionBegin; 1565 if (max) *max = NULL; 1566 if (min) *min = NULL; 1567 for (i = 0; i < ts->numbermonitors; i++) { 1568 if (ts->monitor[i] == TSMonitorEnvelope) { 1569 TSMonitorEnvelopeCtx ctx = (TSMonitorEnvelopeCtx)ts->monitorcontext[i]; 1570 if (max) *max = ctx->max; 1571 if (min) *min = ctx->min; 1572 break; 1573 } 1574 } 1575 PetscFunctionReturn(PETSC_SUCCESS); 1576 } 1577 1578 /*@C 1579 TSMonitorEnvelopeCtxDestroy - Destroys a context that was created with `TSMonitorEnvelopeCtxCreate()`. 1580 1581 Collective 1582 1583 Input Parameter: 1584 . ctx - the monitor context 1585 1586 Level: intermediate 1587 1588 .seealso: [](chapter_ts), `TS`, `TSMonitorLGCtxCreate()`, `TSMonitorSet()`, `TSMonitorLGTimeStep()` 1589 @*/ 1590 PetscErrorCode TSMonitorEnvelopeCtxDestroy(TSMonitorEnvelopeCtx *ctx) 1591 { 1592 PetscFunctionBegin; 1593 PetscCall(VecDestroy(&(*ctx)->min)); 1594 PetscCall(VecDestroy(&(*ctx)->max)); 1595 PetscCall(PetscFree(*ctx)); 1596 PetscFunctionReturn(PETSC_SUCCESS); 1597 } 1598 1599 /*@C 1600 TSDMSwarmMonitorMoments - Monitors the first three moments of a `DMSWARM` being evolved by the `TS` 1601 1602 Not Collective 1603 1604 Input Parameters: 1605 + ts - the `TS` context 1606 . step - current timestep 1607 . t - current time 1608 . u - current solution 1609 - ctx - not used 1610 1611 Options Database Key: 1612 . -ts_dmswarm_monitor_moments - Monitor moments of particle distribution 1613 1614 Level: intermediate 1615 1616 Notes: 1617 This requires a `DMSWARM` be attached to the `TS`. 1618 1619 This is not called directly by users, rather one calls `TSMonitorSet()`, with this function as an argument, to cause the monitor 1620 to be used during the TS integration. 1621 1622 .seealso: [](chapter_ts), `TS`, `TSMonitorSet()`, `TSMonitorDefault()`, `DMSWARM` 1623 @*/ 1624 PetscErrorCode TSDMSwarmMonitorMoments(TS ts, PetscInt step, PetscReal t, Vec U, PetscViewerAndFormat *vf) 1625 { 1626 DM sw; 1627 const PetscScalar *u; 1628 PetscReal m = 1.0, totE = 0., totMom[3] = {0., 0., 0.}; 1629 PetscInt dim, d, Np, p; 1630 MPI_Comm comm; 1631 1632 PetscFunctionBeginUser; 1633 PetscCall(TSGetDM(ts, &sw)); 1634 if (!sw || step % ts->monitorFrequency != 0) PetscFunctionReturn(PETSC_SUCCESS); 1635 PetscCall(PetscObjectGetComm((PetscObject)ts, &comm)); 1636 PetscCall(DMGetDimension(sw, &dim)); 1637 PetscCall(VecGetLocalSize(U, &Np)); 1638 Np /= dim; 1639 PetscCall(VecGetArrayRead(U, &u)); 1640 for (p = 0; p < Np; ++p) { 1641 for (d = 0; d < dim; ++d) { 1642 totE += PetscRealPart(u[p * dim + d] * u[p * dim + d]); 1643 totMom[d] += PetscRealPart(u[p * dim + d]); 1644 } 1645 } 1646 PetscCall(VecRestoreArrayRead(U, &u)); 1647 for (d = 0; d < dim; ++d) totMom[d] *= m; 1648 totE *= 0.5 * m; 1649 PetscCall(PetscPrintf(comm, "Step %4" PetscInt_FMT " Total Energy: %10.8lf", step, (double)totE)); 1650 for (d = 0; d < dim; ++d) PetscCall(PetscPrintf(comm, " Total Momentum %c: %10.8lf", (char)('x' + d), (double)totMom[d])); 1651 PetscCall(PetscPrintf(comm, "\n")); 1652 PetscFunctionReturn(PETSC_SUCCESS); 1653 } 1654