1 #include <petsc/private/dmimpl.h> 2 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 3 4 typedef struct { 5 PetscErrorCode (*boundarylocal)(DM, PetscReal, Vec, Vec, void *); 6 PetscErrorCode (*ifunctionlocal)(DM, PetscReal, Vec, Vec, Vec, void *); 7 PetscErrorCode (*ijacobianlocal)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *); 8 PetscErrorCode (*rhsfunctionlocal)(DM, PetscReal, Vec, Vec, void *); 9 void *boundarylocalctx; 10 void *ifunctionlocalctx; 11 void *ijacobianlocalctx; 12 void *rhsfunctionlocalctx; 13 Vec lumpedmassinv; 14 Mat mass; 15 KSP kspmass; 16 } DMTS_Local; 17 18 static PetscErrorCode DMTSDestroy_DMLocal(DMTS tdm) 19 { 20 PetscFunctionBegin; 21 PetscCall(PetscFree(tdm->data)); 22 PetscFunctionReturn(PETSC_SUCCESS); 23 } 24 25 static PetscErrorCode DMTSDuplicate_DMLocal(DMTS oldtdm, DMTS tdm) 26 { 27 PetscFunctionBegin; 28 PetscCall(PetscNew((DMTS_Local **)&tdm->data)); 29 if (oldtdm->data) PetscCall(PetscMemcpy(tdm->data, oldtdm->data, sizeof(DMTS_Local))); 30 PetscFunctionReturn(PETSC_SUCCESS); 31 } 32 33 static PetscErrorCode DMLocalTSGetContext(DM dm, DMTS tdm, DMTS_Local **dmlocalts) 34 { 35 PetscFunctionBegin; 36 *dmlocalts = NULL; 37 if (!tdm->data) { 38 PetscCall(PetscNew((DMTS_Local **)&tdm->data)); 39 40 tdm->ops->destroy = DMTSDestroy_DMLocal; 41 tdm->ops->duplicate = DMTSDuplicate_DMLocal; 42 } 43 *dmlocalts = (DMTS_Local *)tdm->data; 44 PetscFunctionReturn(PETSC_SUCCESS); 45 } 46 47 static PetscErrorCode TSComputeIFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, Vec F, void *ctx) 48 { 49 DM dm; 50 Vec locX, locX_t, locF; 51 DMTS_Local *dmlocalts = (DMTS_Local *)ctx; 52 53 PetscFunctionBegin; 54 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 55 PetscValidHeaderSpecific(X, VEC_CLASSID, 3); 56 PetscValidHeaderSpecific(X_t, VEC_CLASSID, 4); 57 PetscValidHeaderSpecific(F, VEC_CLASSID, 5); 58 PetscCall(TSGetDM(ts, &dm)); 59 PetscCall(DMGetLocalVector(dm, &locX)); 60 PetscCall(DMGetLocalVector(dm, &locX_t)); 61 PetscCall(DMGetLocalVector(dm, &locF)); 62 PetscCall(VecZeroEntries(locX)); 63 PetscCall(VecZeroEntries(locX_t)); 64 if (dmlocalts->boundarylocal) PetscCall((*dmlocalts->boundarylocal)(dm, time, locX, locX_t, dmlocalts->boundarylocalctx)); 65 PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX)); 66 PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX)); 67 PetscCall(DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t)); 68 PetscCall(DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t)); 69 PetscCall(VecZeroEntries(locF)); 70 CHKMEMQ; 71 PetscCall((*dmlocalts->ifunctionlocal)(dm, time, locX, locX_t, locF, dmlocalts->ifunctionlocalctx)); 72 CHKMEMQ; 73 PetscCall(VecZeroEntries(F)); 74 PetscCall(DMLocalToGlobalBegin(dm, locF, ADD_VALUES, F)); 75 PetscCall(DMLocalToGlobalEnd(dm, locF, ADD_VALUES, F)); 76 PetscCall(DMRestoreLocalVector(dm, &locX)); 77 PetscCall(DMRestoreLocalVector(dm, &locX_t)); 78 PetscCall(DMRestoreLocalVector(dm, &locF)); 79 PetscFunctionReturn(PETSC_SUCCESS); 80 } 81 82 static PetscErrorCode TSComputeRHSFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec F, void *ctx) 83 { 84 DM dm; 85 Vec locX, locF; 86 DMTS_Local *dmlocalts = (DMTS_Local *)ctx; 87 88 PetscFunctionBegin; 89 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 90 PetscValidHeaderSpecific(X, VEC_CLASSID, 3); 91 PetscValidHeaderSpecific(F, VEC_CLASSID, 4); 92 PetscCall(TSGetDM(ts, &dm)); 93 PetscCall(DMGetLocalVector(dm, &locX)); 94 PetscCall(DMGetLocalVector(dm, &locF)); 95 PetscCall(VecZeroEntries(locX)); 96 if (dmlocalts->boundarylocal) PetscCall((*dmlocalts->boundarylocal)(dm, time, locX, NULL, dmlocalts->boundarylocalctx)); 97 PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX)); 98 PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX)); 99 PetscCall(VecZeroEntries(locF)); 100 CHKMEMQ; 101 PetscCall((*dmlocalts->rhsfunctionlocal)(dm, time, locX, locF, dmlocalts->rhsfunctionlocalctx)); 102 CHKMEMQ; 103 PetscCall(VecZeroEntries(F)); 104 PetscCall(DMLocalToGlobalBegin(dm, locF, ADD_VALUES, F)); 105 PetscCall(DMLocalToGlobalEnd(dm, locF, ADD_VALUES, F)); 106 if (dmlocalts->lumpedmassinv) { 107 PetscCall(VecPointwiseMult(F, dmlocalts->lumpedmassinv, F)); 108 } else if (dmlocalts->kspmass) { 109 Vec tmp; 110 111 PetscCall(DMGetGlobalVector(dm, &tmp)); 112 PetscCall(KSPSolve(dmlocalts->kspmass, F, tmp)); 113 PetscCall(VecCopy(tmp, F)); 114 PetscCall(DMRestoreGlobalVector(dm, &tmp)); 115 } 116 PetscCall(DMRestoreLocalVector(dm, &locX)); 117 PetscCall(DMRestoreLocalVector(dm, &locF)); 118 PetscFunctionReturn(PETSC_SUCCESS); 119 } 120 121 static PetscErrorCode TSComputeIJacobian_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, PetscReal a, Mat A, Mat B, void *ctx) 122 { 123 DM dm; 124 Vec locX, locX_t; 125 DMTS_Local *dmlocalts = (DMTS_Local *)ctx; 126 127 PetscFunctionBegin; 128 PetscCall(TSGetDM(ts, &dm)); 129 if (dmlocalts->ijacobianlocal) { 130 PetscCall(DMGetLocalVector(dm, &locX)); 131 PetscCall(DMGetLocalVector(dm, &locX_t)); 132 PetscCall(VecZeroEntries(locX)); 133 PetscCall(VecZeroEntries(locX_t)); 134 if (dmlocalts->boundarylocal) PetscCall((*dmlocalts->boundarylocal)(dm, time, locX, locX_t, dmlocalts->boundarylocalctx)); 135 PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX)); 136 PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX)); 137 PetscCall(DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t)); 138 PetscCall(DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t)); 139 CHKMEMQ; 140 PetscCall((*dmlocalts->ijacobianlocal)(dm, time, locX, locX_t, a, A, B, dmlocalts->ijacobianlocalctx)); 141 CHKMEMQ; 142 PetscCall(DMRestoreLocalVector(dm, &locX)); 143 PetscCall(DMRestoreLocalVector(dm, &locX_t)); 144 } else { 145 MatFDColoring fdcoloring; 146 PetscCall(PetscObjectQuery((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject *)&fdcoloring)); 147 if (!fdcoloring) { 148 ISColoring coloring; 149 150 PetscCall(DMCreateColoring(dm, dm->coloringtype, &coloring)); 151 PetscCall(MatFDColoringCreate(B, coloring, &fdcoloring)); 152 PetscCall(ISColoringDestroy(&coloring)); 153 switch (dm->coloringtype) { 154 case IS_COLORING_GLOBAL: 155 PetscCall(MatFDColoringSetFunction(fdcoloring, (PetscErrorCode(*)(void))TSComputeIFunction_DMLocal, dmlocalts)); 156 break; 157 default: 158 SETERRQ(PetscObjectComm((PetscObject)ts), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]); 159 } 160 PetscCall(PetscObjectSetOptionsPrefix((PetscObject)fdcoloring, ((PetscObject)dm)->prefix)); 161 PetscCall(MatFDColoringSetFromOptions(fdcoloring)); 162 PetscCall(MatFDColoringSetUp(B, coloring, fdcoloring)); 163 PetscCall(PetscObjectCompose((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject)fdcoloring)); 164 PetscCall(PetscObjectDereference((PetscObject)fdcoloring)); 165 166 /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call 167 * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the 168 * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually 169 * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be 170 * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed. 171 */ 172 PetscCall(PetscObjectDereference((PetscObject)dm)); 173 } 174 PetscCall(MatFDColoringApply(B, fdcoloring, X, ts)); 175 } 176 /* This will be redundant if the user called both, but it's too common to forget. */ 177 if (A != B) { 178 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 179 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 180 } 181 PetscFunctionReturn(PETSC_SUCCESS); 182 } 183 184 /*@C 185 DMTSSetBoundaryLocal - set the function for essential boundary data for a local implicit function evaluation. 186 187 Logically Collective 188 189 Input Parameters: 190 + dm - `DM` to associate callback with 191 . func - local function evaluation 192 - ctx - context for function evaluation 193 194 Level: intermediate 195 196 Notes: 197 `func` should set the essential boundary data for the local portion of the solution, as 198 well its time derivative (if it is not `NULL`). 199 200 Vectors are initialized to zero before this function, so it is only needed for non 201 homogeneous data. 202 203 This function is somewhat optional: boundary data could potentially be inserted by a function 204 passed to `DMTSSetIFunctionLocal()`. The use case for this function is for discretizations 205 with constraints (see `DMGetDefaultConstraints()`): this function inserts boundary values 206 before constraint interpolation. 207 208 .seealso: [](ch_ts), `DM`, `TS`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 209 @*/ 210 PetscErrorCode DMTSSetBoundaryLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, void *), void *ctx) 211 { 212 DMTS tdm; 213 DMTS_Local *dmlocalts; 214 215 PetscFunctionBegin; 216 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 217 PetscCall(DMGetDMTSWrite(dm, &tdm)); 218 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 219 220 dmlocalts->boundarylocal = func; 221 dmlocalts->boundarylocalctx = ctx; 222 223 PetscFunctionReturn(PETSC_SUCCESS); 224 } 225 226 /*@C 227 DMTSGetIFunctionLocal - get the local implicit function evaluation function. This function is called with local vector 228 containing the local vector information PLUS ghost point information. It should compute a result for all local 229 elements and `DM` will automatically accumulate the overlapping values. 230 231 Logically Collective 232 233 Input Parameter: 234 . dm - `DM` to associate callback with 235 236 Output Parameters: 237 + func - local function evaluation 238 - ctx - context for function evaluation 239 240 Level: beginner 241 242 .seealso: [](ch_ts), `DM`, `DMTSSetIFunctionLocal()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 243 @*/ 244 PetscErrorCode DMTSGetIFunctionLocal(DM dm, PetscErrorCode (**func)(DM, PetscReal, Vec, Vec, Vec, void *), void **ctx) 245 { 246 DMTS tdm; 247 DMTS_Local *dmlocalts; 248 249 PetscFunctionBegin; 250 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 251 PetscCall(DMGetDMTS(dm, &tdm)); 252 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 253 if (func) { 254 PetscAssertPointer(func, 2); 255 *func = dmlocalts->ifunctionlocal; 256 } 257 if (ctx) { 258 PetscAssertPointer(ctx, 3); 259 *ctx = dmlocalts->ifunctionlocalctx; 260 } 261 PetscFunctionReturn(PETSC_SUCCESS); 262 } 263 264 /*@C 265 DMTSSetIFunctionLocal - set a local implicit function evaluation function. This function is called with local vector 266 containing the local vector information PLUS ghost point information. It should compute a result for all local 267 elements and `DM` will automatically accumulate the overlapping values. 268 269 Logically Collective 270 271 Input Parameters: 272 + dm - `DM` to associate callback with 273 . func - local function evaluation 274 - ctx - context for function evaluation 275 276 Level: beginner 277 278 .seealso: [](ch_ts), `DM`, `DMTSGetIFunctionLocal()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 279 @*/ 280 PetscErrorCode DMTSSetIFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, Vec, void *), void *ctx) 281 { 282 DMTS tdm; 283 DMTS_Local *dmlocalts; 284 285 PetscFunctionBegin; 286 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 287 PetscCall(DMGetDMTSWrite(dm, &tdm)); 288 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 289 290 dmlocalts->ifunctionlocal = func; 291 dmlocalts->ifunctionlocalctx = ctx; 292 293 PetscCall(DMTSSetIFunction(dm, TSComputeIFunction_DMLocal, dmlocalts)); 294 if (!tdm->ops->ijacobian) { /* Call us for the Jacobian too, can be overridden by the user. */ 295 PetscCall(DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts)); 296 } 297 PetscFunctionReturn(PETSC_SUCCESS); 298 } 299 300 /*@C 301 DMTSGetIJacobianLocal - get a local Jacobian evaluation function 302 303 Logically Collective 304 305 Input Parameter: 306 . dm - `DM` to associate callback with 307 308 Output Parameters: 309 + func - local Jacobian evaluation 310 - ctx - optional context for local Jacobian evaluation 311 312 Level: beginner 313 314 .seealso: [](ch_ts), `DM`, `DMTSSetIJacobianLocal()`, `DMTSSetIFunctionLocal()`, `DMTSSetIJacobian()`, `DMTSSetIFunction()` 315 @*/ 316 PetscErrorCode DMTSGetIJacobianLocal(DM dm, PetscErrorCode (**func)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void **ctx) 317 { 318 DMTS tdm; 319 DMTS_Local *dmlocalts; 320 321 PetscFunctionBegin; 322 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 323 PetscCall(DMGetDMTS(dm, &tdm)); 324 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 325 if (func) { 326 PetscAssertPointer(func, 2); 327 *func = dmlocalts->ijacobianlocal; 328 } 329 if (ctx) { 330 PetscAssertPointer(ctx, 3); 331 *ctx = dmlocalts->ijacobianlocalctx; 332 } 333 PetscFunctionReturn(PETSC_SUCCESS); 334 } 335 336 /*@C 337 DMTSSetIJacobianLocal - set a local Jacobian evaluation function 338 339 Logically Collective 340 341 Input Parameters: 342 + dm - `DM` to associate callback with 343 . func - local Jacobian evaluation 344 - ctx - optional context for local Jacobian evaluation 345 346 Level: beginner 347 348 .seealso: [](ch_ts), `DM`, `DMTSGetIJacobianLocal()`, `DMTSSetIFunctionLocal()`, `DMTSSetIJacobian()`, `DMTSSetIFunction()` 349 @*/ 350 PetscErrorCode DMTSSetIJacobianLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void *ctx) 351 { 352 DMTS tdm; 353 DMTS_Local *dmlocalts; 354 355 PetscFunctionBegin; 356 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 357 PetscCall(DMGetDMTSWrite(dm, &tdm)); 358 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 359 360 dmlocalts->ijacobianlocal = func; 361 dmlocalts->ijacobianlocalctx = ctx; 362 363 PetscCall(DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts)); 364 PetscFunctionReturn(PETSC_SUCCESS); 365 } 366 367 /*@C 368 DMTSGetRHSFunctionLocal - get a local rhs function evaluation function. This function is called with local vector 369 containing the local vector information PLUS ghost point information. It should compute a result for all local 370 elements and `DM` will automatically accumulate the overlapping values. 371 372 Logically Collective 373 374 Input Parameter: 375 . dm - `DM` to associate callback with 376 377 Output Parameters: 378 + func - local function evaluation 379 - ctx - context for function evaluation 380 381 Level: beginner 382 383 .seealso: [](ch_ts), `DM`, `DMTSSetRHSFunctionLocal()`, `DMTSSetRHSFunction()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 384 @*/ 385 PetscErrorCode DMTSGetRHSFunctionLocal(DM dm, PetscErrorCode (**func)(DM, PetscReal, Vec, Vec, void *), void **ctx) 386 { 387 DMTS tdm; 388 DMTS_Local *dmlocalts; 389 390 PetscFunctionBegin; 391 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 392 PetscCall(DMGetDMTS(dm, &tdm)); 393 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 394 if (func) { 395 PetscAssertPointer(func, 2); 396 *func = dmlocalts->rhsfunctionlocal; 397 } 398 if (ctx) { 399 PetscAssertPointer(ctx, 3); 400 *ctx = dmlocalts->rhsfunctionlocalctx; 401 } 402 PetscFunctionReturn(PETSC_SUCCESS); 403 } 404 405 /*@C 406 DMTSSetRHSFunctionLocal - set a local rhs function evaluation function. This function is called with local vector 407 containing the local vector information PLUS ghost point information. It should compute a result for all local 408 elements and `DM` will automatically accumulate the overlapping values. 409 410 Logically Collective 411 412 Input Parameters: 413 + dm - `DM` to associate callback with 414 . func - local function evaluation 415 - ctx - context for function evaluation 416 417 Level: beginner 418 419 .seealso: [](ch_ts), `DM`, `DMTSGetRHSFunctionLocal()`, `DMTSSetRHSFunction()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 420 @*/ 421 PetscErrorCode DMTSSetRHSFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, void *), void *ctx) 422 { 423 DMTS tdm; 424 DMTS_Local *dmlocalts; 425 426 PetscFunctionBegin; 427 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 428 PetscCall(DMGetDMTSWrite(dm, &tdm)); 429 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 430 431 dmlocalts->rhsfunctionlocal = func; 432 dmlocalts->rhsfunctionlocalctx = ctx; 433 434 PetscCall(DMTSSetRHSFunction(dm, TSComputeRHSFunction_DMLocal, dmlocalts)); 435 PetscFunctionReturn(PETSC_SUCCESS); 436 } 437 438 /*@C 439 DMTSCreateRHSMassMatrix - This creates the mass matrix associated with the given `DM`, and a solver to invert it, and stores them in the `DM` context. 440 441 Collective 442 443 Input Parameter: 444 . dm - `DM` providing the mass matrix 445 446 Level: developer 447 448 Note: 449 The idea here is that an explicit system can be given a mass matrix, based on the `DM`, which is inverted on the RHS at each step. 450 451 .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrixLumped()`, `DMTSDestroyRHSMassMatrix()`, `DMCreateMassMatrix()`, `DMTS` 452 @*/ 453 PetscErrorCode DMTSCreateRHSMassMatrix(DM dm) 454 { 455 DMTS tdm; 456 DMTS_Local *dmlocalts; 457 const char *prefix; 458 459 PetscFunctionBegin; 460 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 461 PetscCall(DMGetDMTSWrite(dm, &tdm)); 462 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 463 PetscCall(DMCreateMassMatrix(dm, dm, &dmlocalts->mass)); 464 PetscCall(KSPCreate(PetscObjectComm((PetscObject)dm), &dmlocalts->kspmass)); 465 PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix)); 466 PetscCall(KSPSetOptionsPrefix(dmlocalts->kspmass, prefix)); 467 PetscCall(KSPAppendOptionsPrefix(dmlocalts->kspmass, "mass_")); 468 PetscCall(KSPSetFromOptions(dmlocalts->kspmass)); 469 PetscCall(KSPSetOperators(dmlocalts->kspmass, dmlocalts->mass, dmlocalts->mass)); 470 PetscFunctionReturn(PETSC_SUCCESS); 471 } 472 473 /*@C 474 DMTSCreateRHSMassMatrixLumped - This creates the lumped mass matrix associated with the given `DM`, and a solver to invert it, and stores them in the `DM` context. 475 476 Collective 477 478 Input Parameter: 479 . dm - `DM` providing the mass matrix 480 481 Level: developer 482 483 Note: 484 The idea here is that an explicit system can be given a mass matrix, based on the `DM`, which is inverted on the RHS at each step. 485 Since the matrix is lumped, inversion is trivial. 486 487 .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrix()`, `DMTSDestroyRHSMassMatrix()`, `DMCreateMassMatrix()`, `DMTS` 488 @*/ 489 PetscErrorCode DMTSCreateRHSMassMatrixLumped(DM dm) 490 { 491 DMTS tdm; 492 DMTS_Local *dmlocalts; 493 494 PetscFunctionBegin; 495 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 496 PetscCall(DMGetDMTSWrite(dm, &tdm)); 497 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 498 PetscCall(DMCreateMassMatrixLumped(dm, &dmlocalts->lumpedmassinv)); 499 PetscCall(VecReciprocal(dmlocalts->lumpedmassinv)); 500 PetscCall(VecViewFromOptions(dmlocalts->lumpedmassinv, NULL, "-lumped_mass_inv_view")); 501 PetscFunctionReturn(PETSC_SUCCESS); 502 } 503 504 /*@C 505 DMTSDestroyRHSMassMatrix - Destroys the mass matrix and solver stored in the `DM` context, if they exist. 506 507 Logically Collective 508 509 Input Parameter: 510 . dm - `DM` providing the mass matrix 511 512 Level: developer 513 514 .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrixLumped()`, `DMCreateMassMatrix()`, `DMTS` 515 @*/ 516 PetscErrorCode DMTSDestroyRHSMassMatrix(DM dm) 517 { 518 DMTS tdm; 519 DMTS_Local *dmlocalts; 520 521 PetscFunctionBegin; 522 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 523 PetscCall(DMGetDMTSWrite(dm, &tdm)); 524 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 525 PetscCall(VecDestroy(&dmlocalts->lumpedmassinv)); 526 PetscCall(MatDestroy(&dmlocalts->mass)); 527 PetscCall(KSPDestroy(&dmlocalts->kspmass)); 528 PetscFunctionReturn(PETSC_SUCCESS); 529 } 530