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