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