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 PetscErrorCode ierr; 249 250 PetscFunctionBegin; 251 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 252 ierr = DMGetDMTS(dm, &tdm); 253 CHKERRQ(ierr); 254 ierr = DMLocalTSGetContext(dm, tdm, &dmlocalts); 255 CHKERRQ(ierr); 256 if (func) { 257 PetscAssertPointer(func, 2); 258 *func = dmlocalts->ifunctionlocal; 259 } 260 if (ctx) { 261 PetscAssertPointer(ctx, 3); 262 *ctx = dmlocalts->ifunctionlocalctx; 263 } 264 PetscFunctionReturn(PETSC_SUCCESS); 265 } 266 267 /*@C 268 DMTSSetIFunctionLocal - set a local implicit function evaluation function. This function is called with local vector 269 containing the local vector information PLUS ghost point information. It should compute a result for all local 270 elements and `DM` will automatically accumulate the overlapping values. 271 272 Logically Collective 273 274 Input Parameters: 275 + dm - `DM` to associate callback with 276 . func - local function evaluation 277 - ctx - context for function evaluation 278 279 Level: beginner 280 281 .seealso: [](ch_ts), `DM`, `DMTSGetIFunctionLocal()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 282 @*/ 283 PetscErrorCode DMTSSetIFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, Vec, void *), void *ctx) 284 { 285 DMTS tdm; 286 DMTS_Local *dmlocalts; 287 288 PetscFunctionBegin; 289 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 290 PetscCall(DMGetDMTSWrite(dm, &tdm)); 291 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 292 293 dmlocalts->ifunctionlocal = func; 294 dmlocalts->ifunctionlocalctx = ctx; 295 296 PetscCall(DMTSSetIFunction(dm, TSComputeIFunction_DMLocal, dmlocalts)); 297 if (!tdm->ops->ijacobian) { /* Call us for the Jacobian too, can be overridden by the user. */ 298 PetscCall(DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts)); 299 } 300 PetscFunctionReturn(PETSC_SUCCESS); 301 } 302 303 /*@C 304 DMTSGetIJacobianLocal - get a local Jacobian evaluation function 305 306 Logically Collective 307 308 Input Parameter: 309 . dm - `DM` to associate callback with 310 311 Output Parameters: 312 + func - local Jacobian evaluation 313 - ctx - optional context for local Jacobian evaluation 314 315 Level: beginner 316 317 .seealso: [](ch_ts), `DM`, `DMTSSetIJacobianLocal()`, `DMTSSetIFunctionLocal()`, `DMTSSetIJacobian()`, `DMTSSetIFunction()` 318 @*/ 319 PetscErrorCode DMTSGetIJacobianLocal(DM dm, PetscErrorCode (**func)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void **ctx) 320 { 321 DMTS tdm; 322 DMTS_Local *dmlocalts; 323 PetscErrorCode ierr; 324 325 PetscFunctionBegin; 326 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 327 ierr = DMGetDMTS(dm, &tdm); 328 CHKERRQ(ierr); 329 ierr = DMLocalTSGetContext(dm, tdm, &dmlocalts); 330 CHKERRQ(ierr); 331 if (func) { 332 PetscAssertPointer(func, 2); 333 *func = dmlocalts->ijacobianlocal; 334 } 335 if (ctx) { 336 PetscAssertPointer(ctx, 3); 337 *ctx = dmlocalts->ijacobianlocalctx; 338 } 339 PetscFunctionReturn(PETSC_SUCCESS); 340 } 341 342 /*@C 343 DMTSSetIJacobianLocal - set a local Jacobian evaluation function 344 345 Logically Collective 346 347 Input Parameters: 348 + dm - `DM` to associate callback with 349 . func - local Jacobian evaluation 350 - ctx - optional context for local Jacobian evaluation 351 352 Level: beginner 353 354 .seealso: [](ch_ts), `DM`, `DMTSGetIJacobianLocal()`, `DMTSSetIFunctionLocal()`, `DMTSSetIJacobian()`, `DMTSSetIFunction()` 355 @*/ 356 PetscErrorCode DMTSSetIJacobianLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void *ctx) 357 { 358 DMTS tdm; 359 DMTS_Local *dmlocalts; 360 361 PetscFunctionBegin; 362 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 363 PetscCall(DMGetDMTSWrite(dm, &tdm)); 364 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 365 366 dmlocalts->ijacobianlocal = func; 367 dmlocalts->ijacobianlocalctx = ctx; 368 369 PetscCall(DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts)); 370 PetscFunctionReturn(PETSC_SUCCESS); 371 } 372 373 /*@C 374 DMTSGetRHSFunctionLocal - get a local rhs function evaluation function. This function is called with local vector 375 containing the local vector information PLUS ghost point information. It should compute a result for all local 376 elements and `DM` will automatically accumulate the overlapping values. 377 378 Logically Collective 379 380 Input Parameter: 381 . dm - `DM` to associate callback with 382 383 Output Parameters: 384 + func - local function evaluation 385 - ctx - context for function evaluation 386 387 Level: beginner 388 389 .seealso: [](ch_ts), `DM`, `DMTSSetRHSFunctionLocal()`, `DMTSSetRHSFunction()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 390 @*/ 391 PetscErrorCode DMTSGetRHSFunctionLocal(DM dm, PetscErrorCode (**func)(DM, PetscReal, Vec, Vec, void *), void **ctx) 392 { 393 DMTS tdm; 394 DMTS_Local *dmlocalts; 395 PetscErrorCode ierr; 396 397 PetscFunctionBegin; 398 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 399 ierr = DMGetDMTS(dm, &tdm); 400 CHKERRQ(ierr); 401 ierr = DMLocalTSGetContext(dm, tdm, &dmlocalts); 402 CHKERRQ(ierr); 403 if (func) { 404 PetscAssertPointer(func, 2); 405 *func = dmlocalts->rhsfunctionlocal; 406 } 407 if (ctx) { 408 PetscAssertPointer(ctx, 3); 409 *ctx = dmlocalts->rhsfunctionlocalctx; 410 } 411 PetscFunctionReturn(PETSC_SUCCESS); 412 } 413 414 /*@C 415 DMTSSetRHSFunctionLocal - set a local rhs function evaluation function. This function is called with local vector 416 containing the local vector information PLUS ghost point information. It should compute a result for all local 417 elements and `DM` will automatically accumulate the overlapping values. 418 419 Logically Collective 420 421 Input Parameters: 422 + dm - `DM` to associate callback with 423 . func - local function evaluation 424 - ctx - context for function evaluation 425 426 Level: beginner 427 428 .seealso: [](ch_ts), `DM`, `DMTSGetRHSFunctionLocal()`, `DMTSSetRHSFunction()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()` 429 @*/ 430 PetscErrorCode DMTSSetRHSFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, void *), void *ctx) 431 { 432 DMTS tdm; 433 DMTS_Local *dmlocalts; 434 435 PetscFunctionBegin; 436 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 437 PetscCall(DMGetDMTSWrite(dm, &tdm)); 438 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 439 440 dmlocalts->rhsfunctionlocal = func; 441 dmlocalts->rhsfunctionlocalctx = ctx; 442 443 PetscCall(DMTSSetRHSFunction(dm, TSComputeRHSFunction_DMLocal, dmlocalts)); 444 PetscFunctionReturn(PETSC_SUCCESS); 445 } 446 447 /*@C 448 DMTSCreateRHSMassMatrix - This creates the mass matrix associated with the given `DM`, and a solver to invert it, and stores them in the `DM` context. 449 450 Collective 451 452 Input Parameter: 453 . dm - `DM` providing the mass matrix 454 455 Level: developer 456 457 Note: 458 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. 459 460 .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrixLumped()`, `DMTSDestroyRHSMassMatrix()`, `DMCreateMassMatrix()`, `DMTS` 461 @*/ 462 PetscErrorCode DMTSCreateRHSMassMatrix(DM dm) 463 { 464 DMTS tdm; 465 DMTS_Local *dmlocalts; 466 const char *prefix; 467 468 PetscFunctionBegin; 469 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 470 PetscCall(DMGetDMTSWrite(dm, &tdm)); 471 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 472 PetscCall(DMCreateMassMatrix(dm, dm, &dmlocalts->mass)); 473 PetscCall(KSPCreate(PetscObjectComm((PetscObject)dm), &dmlocalts->kspmass)); 474 PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix)); 475 PetscCall(KSPSetOptionsPrefix(dmlocalts->kspmass, prefix)); 476 PetscCall(KSPAppendOptionsPrefix(dmlocalts->kspmass, "mass_")); 477 PetscCall(KSPSetFromOptions(dmlocalts->kspmass)); 478 PetscCall(KSPSetOperators(dmlocalts->kspmass, dmlocalts->mass, dmlocalts->mass)); 479 PetscFunctionReturn(PETSC_SUCCESS); 480 } 481 482 /*@C 483 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. 484 485 Collective 486 487 Input Parameter: 488 . dm - `DM` providing the mass matrix 489 490 Level: developer 491 492 Note: 493 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. 494 Since the matrix is lumped, inversion is trivial. 495 496 .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrix()`, `DMTSDestroyRHSMassMatrix()`, `DMCreateMassMatrix()`, `DMTS` 497 @*/ 498 PetscErrorCode DMTSCreateRHSMassMatrixLumped(DM dm) 499 { 500 DMTS tdm; 501 DMTS_Local *dmlocalts; 502 503 PetscFunctionBegin; 504 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 505 PetscCall(DMGetDMTSWrite(dm, &tdm)); 506 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 507 PetscCall(DMCreateMassMatrixLumped(dm, &dmlocalts->lumpedmassinv)); 508 PetscCall(VecReciprocal(dmlocalts->lumpedmassinv)); 509 PetscCall(VecViewFromOptions(dmlocalts->lumpedmassinv, NULL, "-lumped_mass_inv_view")); 510 PetscFunctionReturn(PETSC_SUCCESS); 511 } 512 513 /*@C 514 DMTSDestroyRHSMassMatrix - Destroys the mass matrix and solver stored in the `DM` context, if they exist. 515 516 Logically Collective 517 518 Input Parameter: 519 . dm - `DM` providing the mass matrix 520 521 Level: developer 522 523 .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrixLumped()`, `DMCreateMassMatrix()`, `DMTS` 524 @*/ 525 PetscErrorCode DMTSDestroyRHSMassMatrix(DM dm) 526 { 527 DMTS tdm; 528 DMTS_Local *dmlocalts; 529 530 PetscFunctionBegin; 531 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 532 PetscCall(DMGetDMTSWrite(dm, &tdm)); 533 PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts)); 534 PetscCall(VecDestroy(&dmlocalts->lumpedmassinv)); 535 PetscCall(MatDestroy(&dmlocalts->mass)); 536 PetscCall(KSPDestroy(&dmlocalts->kspmass)); 537 PetscFunctionReturn(PETSC_SUCCESS); 538 } 539