Actual source code: precon.c

  1: /*
  2:     The PC (preconditioner) interface routines, callable by users.
  3: */
  4: #include <petsc/private/pcimpl.h>
  5: #include <petscdm.h>

  7: /* Logging support */
  8: PetscClassId  PC_CLASSID;
  9: PetscLogEvent PC_SetUp, PC_SetUpOnBlocks, PC_Apply, PC_MatApply, PC_ApplyCoarse, PC_ApplyMultiple, PC_ApplySymmetricLeft;
 10: PetscLogEvent PC_ApplySymmetricRight, PC_ModifySubMatrices, PC_ApplyOnBlocks, PC_ApplyTransposeOnBlocks;
 11: PetscInt      PetscMGLevelId;
 12: PetscLogStage PCMPIStage;

 14: PETSC_INTERN PetscErrorCode PCGetDefaultType_Private(PC pc, const char *type[])
 15: {
 16:   PetscMPIInt size;
 17:   PetscBool   hasopblock, hasopsolve, flg1, flg2, set, flg3, isnormal;

 19:   PetscFunctionBegin;
 20:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)pc), &size));
 21:   if (pc->pmat) {
 22:     PetscCall(MatHasOperation(pc->pmat, MATOP_GET_DIAGONAL_BLOCK, &hasopblock));
 23:     PetscCall(MatHasOperation(pc->pmat, MATOP_SOLVE, &hasopsolve));
 24:     if (size == 1) {
 25:       PetscCall(MatGetFactorAvailable(pc->pmat, "petsc", MAT_FACTOR_ICC, &flg1));
 26:       PetscCall(MatGetFactorAvailable(pc->pmat, "petsc", MAT_FACTOR_ILU, &flg2));
 27:       PetscCall(MatIsSymmetricKnown(pc->pmat, &set, &flg3));
 28:       PetscCall(PetscObjectTypeCompareAny((PetscObject)pc->pmat, &isnormal, MATNORMAL, MATNORMALHERMITIAN, NULL));
 29:       if (flg1 && (!flg2 || (set && flg3))) {
 30:         *type = PCICC;
 31:       } else if (flg2) {
 32:         *type = PCILU;
 33:       } else if (isnormal) {
 34:         *type = PCNONE;
 35:       } else if (hasopblock) { /* likely is a parallel matrix run on one processor */
 36:         *type = PCBJACOBI;
 37:       } else if (hasopsolve) {
 38:         *type = PCMAT;
 39:       } else {
 40:         *type = PCNONE;
 41:       }
 42:     } else {
 43:       if (hasopblock) {
 44:         *type = PCBJACOBI;
 45:       } else if (hasopsolve) {
 46:         *type = PCMAT;
 47:       } else {
 48:         *type = PCNONE;
 49:       }
 50:     }
 51:   } else *type = NULL;
 52:   PetscFunctionReturn(PETSC_SUCCESS);
 53: }

 55: /*@
 56:   PCReset - Resets a `PC` context to the pcsetupcalled = 0 state and removes any allocated `Vec`s and `Mat`s

 58:   Collective

 60:   Input Parameter:
 61: . pc - the preconditioner context

 63:   Level: developer

 65:   Note:
 66:   This allows a `PC` to be reused for a different sized linear system but using the same options that have been previously set in `pc`

 68: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCSetUp()`
 69: @*/
 70: PetscErrorCode PCReset(PC pc)
 71: {
 72:   PetscFunctionBegin;
 74:   PetscTryTypeMethod(pc, reset);
 75:   PetscCall(VecDestroy(&pc->diagonalscaleright));
 76:   PetscCall(VecDestroy(&pc->diagonalscaleleft));
 77:   PetscCall(MatDestroy(&pc->pmat));
 78:   PetscCall(MatDestroy(&pc->mat));

 80:   pc->setupcalled = 0;
 81:   PetscFunctionReturn(PETSC_SUCCESS);
 82: }

 84: /*@C
 85:   PCDestroy - Destroys `PC` context that was created with `PCCreate()`.

 87:   Collective

 89:   Input Parameter:
 90: . pc - the preconditioner context

 92:   Level: developer

 94: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCSetUp()`
 95: @*/
 96: PetscErrorCode PCDestroy(PC *pc)
 97: {
 98:   PetscFunctionBegin;
 99:   if (!*pc) PetscFunctionReturn(PETSC_SUCCESS);
101:   if (--((PetscObject)*pc)->refct > 0) {
102:     *pc = NULL;
103:     PetscFunctionReturn(PETSC_SUCCESS);
104:   }

106:   PetscCall(PCReset(*pc));

108:   /* if memory was published with SAWs then destroy it */
109:   PetscCall(PetscObjectSAWsViewOff((PetscObject)*pc));
110:   PetscTryTypeMethod(*pc, destroy);
111:   PetscCall(DMDestroy(&(*pc)->dm));
112:   PetscCall(PetscHeaderDestroy(pc));
113:   PetscFunctionReturn(PETSC_SUCCESS);
114: }

116: /*@
117:   PCGetDiagonalScale - Indicates if the preconditioner applies an additional left and right
118:   scaling as needed by certain time-stepping codes.

120:   Logically Collective

122:   Input Parameter:
123: . pc - the preconditioner context

125:   Output Parameter:
126: . flag - `PETSC_TRUE` if it applies the scaling

128:   Level: developer

130:   Note:
131:   If this returns `PETSC_TRUE` then the system solved via the Krylov method is, for left and right preconditioning,

133:   $$
134:   \begin{align*}
135:   D M A D^{-1} y = D M b  \\
136:   D A M D^{-1} z = D b.
137:   \end{align*}
138:   $$

140: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCSetUp()`, `PCDiagonalScaleLeft()`, `PCDiagonalScaleRight()`, `PCSetDiagonalScale()`
141: @*/
142: PetscErrorCode PCGetDiagonalScale(PC pc, PetscBool *flag)
143: {
144:   PetscFunctionBegin;
146:   PetscAssertPointer(flag, 2);
147:   *flag = pc->diagonalscale;
148:   PetscFunctionReturn(PETSC_SUCCESS);
149: }

151: /*@
152:   PCSetDiagonalScale - Indicates the left scaling to use to apply an additional left and right
153:   scaling as needed by certain time-stepping codes.

155:   Logically Collective

157:   Input Parameters:
158: + pc - the preconditioner context
159: - s  - scaling vector

161:   Level: intermediate

163:   Notes:
164:   The system solved via the Krylov method is, for left and right preconditioning,
165:   $$
166:   \begin{align*}
167:   D M A D^{-1} y = D M b \\
168:   D A M D^{-1} z = D b.
169:   \end{align*}
170:   $$

172:   `PCDiagonalScaleLeft()` scales a vector by $D$. `PCDiagonalScaleRight()` scales a vector by $D^{-1}$.

174: .seealso: [](ch_ksp), `PCCreate()`, `PCSetUp()`, `PCDiagonalScaleLeft()`, `PCDiagonalScaleRight()`, `PCGetDiagonalScale()`
175: @*/
176: PetscErrorCode PCSetDiagonalScale(PC pc, Vec s)
177: {
178:   PetscFunctionBegin;
181:   pc->diagonalscale = PETSC_TRUE;

183:   PetscCall(PetscObjectReference((PetscObject)s));
184:   PetscCall(VecDestroy(&pc->diagonalscaleleft));

186:   pc->diagonalscaleleft = s;

188:   PetscCall(VecDuplicate(s, &pc->diagonalscaleright));
189:   PetscCall(VecCopy(s, pc->diagonalscaleright));
190:   PetscCall(VecReciprocal(pc->diagonalscaleright));
191:   PetscFunctionReturn(PETSC_SUCCESS);
192: }

194: /*@
195:   PCDiagonalScaleLeft - Scales a vector by the left scaling as needed by certain time-stepping codes.

197:   Logically Collective

199:   Input Parameters:
200: + pc  - the preconditioner context
201: . in  - input vector
202: - out - scaled vector (maybe the same as in)

204:   Level: intermediate

206:   Notes:
207:   The system solved via the Krylov method is, for left and right preconditioning,

209:   $$
210:   \begin{align*}
211:   D M A D^{-1} y = D M b  \\
212:   D A M D^{-1} z = D b.
213:   \end{align*}
214:   $$

216:   `PCDiagonalScaleLeft()` scales a vector by $D$. `PCDiagonalScaleRight()` scales a vector by $D^{-1}$.

218:   If diagonal scaling is turned off and `in` is not `out` then `in` is copied to `out`

220: .seealso: [](ch_ksp), `PCCreate()`, `PCSetUp()`, `PCSetDiagonalScale()`, `PCDiagonalScaleRight()`, `PCDiagonalScale()`
221: @*/
222: PetscErrorCode PCDiagonalScaleLeft(PC pc, Vec in, Vec out)
223: {
224:   PetscFunctionBegin;
228:   if (pc->diagonalscale) {
229:     PetscCall(VecPointwiseMult(out, pc->diagonalscaleleft, in));
230:   } else if (in != out) {
231:     PetscCall(VecCopy(in, out));
232:   }
233:   PetscFunctionReturn(PETSC_SUCCESS);
234: }

236: /*@
237:   PCDiagonalScaleRight - Scales a vector by the right scaling as needed by certain time-stepping codes.

239:   Logically Collective

241:   Input Parameters:
242: + pc  - the preconditioner context
243: . in  - input vector
244: - out - scaled vector (maybe the same as in)

246:   Level: intermediate

248:   Notes:
249:   The system solved via the Krylov method is, for left and right preconditioning,

251:   $$
252:   \begin{align*}
253:   D M A D^{-1} y = D M b  \\
254:   D A M D^{-1} z = D b.
255:   \end{align*}
256:   $$

258:   `PCDiagonalScaleLeft()` scales a vector by $D$. `PCDiagonalScaleRight()` scales a vector by $D^{-1}$.

260:   If diagonal scaling is turned off and `in` is not `out` then `in` is copied to `out`

262: .seealso: [](ch_ksp), `PCCreate()`, `PCSetUp()`, `PCDiagonalScaleLeft()`, `PCSetDiagonalScale()`, `PCDiagonalScale()`
263: @*/
264: PetscErrorCode PCDiagonalScaleRight(PC pc, Vec in, Vec out)
265: {
266:   PetscFunctionBegin;
270:   if (pc->diagonalscale) {
271:     PetscCall(VecPointwiseMult(out, pc->diagonalscaleright, in));
272:   } else if (in != out) {
273:     PetscCall(VecCopy(in, out));
274:   }
275:   PetscFunctionReturn(PETSC_SUCCESS);
276: }

278: /*@
279:   PCSetUseAmat - Sets a flag to indicate that when the preconditioner needs to apply (part of) the
280:   operator during the preconditioning process it applies the Amat provided to `TSSetRHSJacobian()`,
281:   `TSSetIJacobian()`, `SNESSetJacobian()`, `KSPSetOperators()` or `PCSetOperators()` not the Pmat.

283:   Logically Collective

285:   Input Parameters:
286: + pc  - the preconditioner context
287: - flg - `PETSC_TRUE` to use the Amat, `PETSC_FALSE` to use the Pmat (default is false)

289:   Options Database Key:
290: . -pc_use_amat <true,false> - use the amat argument to `KSPSetOperators()` or `PCSetOperators()` to apply the operator

292:   Level: intermediate

294:   Note:
295:   For the common case in which the linear system matrix and the matrix used to construct the
296:   preconditioner are identical, this routine has no affect.

298: .seealso: [](ch_ksp), `PC`, `PCGetUseAmat()`, `PCBJACOBI`, `PGMG`, `PCFIELDSPLIT`, `PCCOMPOSITE`,
299:           `KSPSetOperators()`, `PCSetOperators()`
300: @*/
301: PetscErrorCode PCSetUseAmat(PC pc, PetscBool flg)
302: {
303:   PetscFunctionBegin;
305:   pc->useAmat = flg;
306:   PetscFunctionReturn(PETSC_SUCCESS);
307: }

309: /*@
310:   PCSetErrorIfFailure - Causes `PC` to generate an error if a floating point exception, for example a zero pivot, is detected.

312:   Logically Collective

314:   Input Parameters:
315: + pc  - iterative context obtained from `PCCreate()`
316: - flg - `PETSC_TRUE` indicates you want the error generated

318:   Level: advanced

320:   Notes:
321:   Normally PETSc continues if a linear solver fails due to a failed setup of a preconditioner, you can call `KSPGetConvergedReason()` after a `KSPSolve()`
322:   to determine if it has converged or failed. Or use -ksp_error_if_not_converged to cause the program to terminate as soon as lack of convergence is
323:   detected.

325:   This is propagated into `KSP`s used by this `PC`, which then propagate it into `PC`s used by those `KSP`s

327: .seealso: [](ch_ksp), `PC`, `KSPSetErrorIfNotConverged()`, `PCGetInitialGuessNonzero()`, `PCSetInitialGuessKnoll()`, `PCGetInitialGuessKnoll()`
328: @*/
329: PetscErrorCode PCSetErrorIfFailure(PC pc, PetscBool flg)
330: {
331:   PetscFunctionBegin;
334:   pc->erroriffailure = flg;
335:   PetscFunctionReturn(PETSC_SUCCESS);
336: }

338: /*@
339:   PCGetUseAmat - Gets a flag to indicate that when the preconditioner needs to apply (part of) the
340:   operator during the preconditioning process it applies the Amat provided to `TSSetRHSJacobian()`,
341:   `TSSetIJacobian()`, `SNESSetJacobian()`, `KSPSetOperators()` or `PCSetOperators()` not the Pmat.

343:   Logically Collective

345:   Input Parameter:
346: . pc - the preconditioner context

348:   Output Parameter:
349: . flg - `PETSC_TRUE` to use the Amat, `PETSC_FALSE` to use the Pmat (default is false)

351:   Level: intermediate

353:   Note:
354:   For the common case in which the linear system matrix and the matrix used to construct the
355:   preconditioner are identical, this routine is does nothing.

357: .seealso: [](ch_ksp), `PC`, `PCSetUseAmat()`, `PCBJACOBI`, `PGMG`, `PCFIELDSPLIT`, `PCCOMPOSITE`
358: @*/
359: PetscErrorCode PCGetUseAmat(PC pc, PetscBool *flg)
360: {
361:   PetscFunctionBegin;
363:   *flg = pc->useAmat;
364:   PetscFunctionReturn(PETSC_SUCCESS);
365: }

367: /*@
368:   PCSetKSPNestLevel - sets the amount of nesting the `KSP` that contains this `PC` has

370:   Collective

372:   Input Parameters:
373: + pc    - the `PC`
374: - level - the nest level

376:   Level: developer

378: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPGMRES`, `KSPType`, `KSPGetNestLevel()`, `PCGetKSPNestLevel()`, `KSPSetNestLevel()`
379: @*/
380: PetscErrorCode PCSetKSPNestLevel(PC pc, PetscInt level)
381: {
382:   PetscFunctionBegin;
385:   pc->kspnestlevel = level;
386:   PetscFunctionReturn(PETSC_SUCCESS);
387: }

389: /*@
390:   PCGetKSPNestLevel - gets the amount of nesting the `KSP` that contains this `PC` has

392:   Not Collective

394:   Input Parameter:
395: . pc - the `PC`

397:   Output Parameter:
398: . level - the nest level

400:   Level: developer

402: .seealso: [](ch_ksp), `KSPSetUp()`, `KSPSolve()`, `KSPDestroy()`, `KSP`, `KSPGMRES`, `KSPType`, `KSPSetNestLevel()`, `PCSetKSPNestLevel()`, `KSPGetNestLevel()`
403: @*/
404: PetscErrorCode PCGetKSPNestLevel(PC pc, PetscInt *level)
405: {
406:   PetscFunctionBegin;
408:   PetscAssertPointer(level, 2);
409:   *level = pc->kspnestlevel;
410:   PetscFunctionReturn(PETSC_SUCCESS);
411: }

413: /*@
414:   PCCreate - Creates a preconditioner context, `PC`

416:   Collective

418:   Input Parameter:
419: . comm - MPI communicator

421:   Output Parameter:
422: . newpc - location to put the preconditioner context

424:   Level: developer

426:   Note:
427:   The default preconditioner for sparse matrices is `PCILU` or `PCICC` with 0 fill on one process and block Jacobi (`PCBJACOBI`) with `PCILU` or `PCICC`
428:   in parallel. For dense matrices it is always `PCNONE`.

430: .seealso: [](ch_ksp), `PC`, `PCSetUp()`, `PCApply()`, `PCDestroy()`
431: @*/
432: PetscErrorCode PCCreate(MPI_Comm comm, PC *newpc)
433: {
434:   PC pc;

436:   PetscFunctionBegin;
437:   PetscAssertPointer(newpc, 2);
438:   *newpc = NULL;
439:   PetscCall(PCInitializePackage());

441:   PetscCall(PetscHeaderCreate(pc, PC_CLASSID, "PC", "Preconditioner", "PC", comm, PCDestroy, PCView));

443:   pc->mat                  = NULL;
444:   pc->pmat                 = NULL;
445:   pc->setupcalled          = 0;
446:   pc->setfromoptionscalled = 0;
447:   pc->data                 = NULL;
448:   pc->diagonalscale        = PETSC_FALSE;
449:   pc->diagonalscaleleft    = NULL;
450:   pc->diagonalscaleright   = NULL;

452:   pc->modifysubmatrices  = NULL;
453:   pc->modifysubmatricesP = NULL;

455:   *newpc = pc;
456:   PetscFunctionReturn(PETSC_SUCCESS);
457: }

459: /*@
460:   PCApply - Applies the preconditioner to a vector.

462:   Collective

464:   Input Parameters:
465: + pc - the preconditioner context
466: - x  - input vector

468:   Output Parameter:
469: . y - output vector

471:   Level: developer

473: .seealso: [](ch_ksp), `PC`, `PCApplyTranspose()`, `PCApplyBAorAB()`
474: @*/
475: PetscErrorCode PCApply(PC pc, Vec x, Vec y)
476: {
477:   PetscInt m, n, mv, nv;

479:   PetscFunctionBegin;
483:   PetscCheck(x != y, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "x and y must be different vectors");
484:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(x, 2, PETSC_TRUE));
485:   /* use pmat to check vector sizes since for KSPLSQR the pmat may be of a different size than mat */
486:   PetscCall(MatGetLocalSize(pc->pmat, &m, &n));
487:   PetscCall(VecGetLocalSize(x, &mv));
488:   PetscCall(VecGetLocalSize(y, &nv));
489:   /* check pmat * y = x is feasible */
490:   PetscCheck(mv == m, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Preconditioner number of local rows %" PetscInt_FMT " does not equal input vector size %" PetscInt_FMT, m, mv);
491:   PetscCheck(nv == n, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Preconditioner number of local columns %" PetscInt_FMT " does not equal output vector size %" PetscInt_FMT, n, nv);
492:   PetscCall(VecSetErrorIfLocked(y, 3));

494:   PetscCall(PCSetUp(pc));
495:   PetscCall(VecLockReadPush(x));
496:   PetscCall(PetscLogEventBegin(PC_Apply, pc, x, y, 0));
497:   PetscUseTypeMethod(pc, apply, x, y);
498:   PetscCall(PetscLogEventEnd(PC_Apply, pc, x, y, 0));
499:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(y, 3, PETSC_FALSE));
500:   PetscCall(VecLockReadPop(x));
501:   PetscFunctionReturn(PETSC_SUCCESS);
502: }

504: /*@
505:   PCMatApply - Applies the preconditioner to multiple vectors stored as a `MATDENSE`. Like `PCApply()`, `Y` and `X` must be different matrices.

507:   Collective

509:   Input Parameters:
510: + pc - the preconditioner context
511: - X  - block of input vectors

513:   Output Parameter:
514: . Y - block of output vectors

516:   Level: developer

518: .seealso: [](ch_ksp), `PC`, `PCApply()`, `KSPMatSolve()`
519: @*/
520: PetscErrorCode PCMatApply(PC pc, Mat X, Mat Y)
521: {
522:   Mat       A;
523:   Vec       cy, cx;
524:   PetscInt  m1, M1, m2, M2, n1, N1, n2, N2, m3, M3, n3, N3;
525:   PetscBool match;

527:   PetscFunctionBegin;
531:   PetscCheckSameComm(pc, 1, X, 2);
532:   PetscCheckSameComm(pc, 1, Y, 3);
533:   PetscCheck(Y != X, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "Y and X must be different matrices");
534:   PetscCall(PCGetOperators(pc, NULL, &A));
535:   PetscCall(MatGetLocalSize(A, &m3, &n3));
536:   PetscCall(MatGetLocalSize(X, &m2, &n2));
537:   PetscCall(MatGetLocalSize(Y, &m1, &n1));
538:   PetscCall(MatGetSize(A, &M3, &N3));
539:   PetscCall(MatGetSize(X, &M2, &N2));
540:   PetscCall(MatGetSize(Y, &M1, &N1));
541:   PetscCheck(n1 == n2 && N1 == N2, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Incompatible number of columns between block of input vectors (n,N) = (%" PetscInt_FMT ",%" PetscInt_FMT ") and block of output vectors (n,N) = (%" PetscInt_FMT ",%" PetscInt_FMT ")", n2, N2, n1, N1);
542:   PetscCheck(m2 == m3 && M2 == M3, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Incompatible layout between block of input vectors (m,M) = (%" PetscInt_FMT ",%" PetscInt_FMT ") and Pmat (m,M)x(n,N) = (%" PetscInt_FMT ",%" PetscInt_FMT ")x(%" PetscInt_FMT ",%" PetscInt_FMT ")", m2, M2, m3, M3, n3, N3);
543:   PetscCheck(m1 == n3 && M1 == N3, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Incompatible layout between block of output vectors (m,M) = (%" PetscInt_FMT ",%" PetscInt_FMT ") and Pmat (m,M)x(n,N) = (%" PetscInt_FMT ",%" PetscInt_FMT ")x(%" PetscInt_FMT ",%" PetscInt_FMT ")", m1, M1, m3, M3, n3, N3);
544:   PetscCall(PetscObjectBaseTypeCompareAny((PetscObject)Y, &match, MATSEQDENSE, MATMPIDENSE, ""));
545:   PetscCheck(match, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Provided block of output vectors not stored in a dense Mat");
546:   PetscCall(PetscObjectBaseTypeCompareAny((PetscObject)X, &match, MATSEQDENSE, MATMPIDENSE, ""));
547:   PetscCheck(match, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Provided block of input vectors not stored in a dense Mat");
548:   PetscCall(PCSetUp(pc));
549:   if (pc->ops->matapply) {
550:     PetscCall(PetscLogEventBegin(PC_MatApply, pc, X, Y, 0));
551:     PetscUseTypeMethod(pc, matapply, X, Y);
552:     PetscCall(PetscLogEventEnd(PC_MatApply, pc, X, Y, 0));
553:   } else {
554:     PetscCall(PetscInfo(pc, "PC type %s applying column by column\n", ((PetscObject)pc)->type_name));
555:     for (n1 = 0; n1 < N1; ++n1) {
556:       PetscCall(MatDenseGetColumnVecRead(X, n1, &cx));
557:       PetscCall(MatDenseGetColumnVecWrite(Y, n1, &cy));
558:       PetscCall(PCApply(pc, cx, cy));
559:       PetscCall(MatDenseRestoreColumnVecWrite(Y, n1, &cy));
560:       PetscCall(MatDenseRestoreColumnVecRead(X, n1, &cx));
561:     }
562:   }
563:   PetscFunctionReturn(PETSC_SUCCESS);
564: }

566: /*@
567:   PCApplySymmetricLeft - Applies the left part of a symmetric preconditioner to a vector.

569:   Collective

571:   Input Parameters:
572: + pc - the preconditioner context
573: - x  - input vector

575:   Output Parameter:
576: . y - output vector

578:   Level: developer

580:   Note:
581:   Currently, this routine is implemented only for `PCICC` and `PCJACOBI` preconditioners.

583: .seealso: [](ch_ksp), `PC`, `PCApply()`, `PCApplySymmetricRight()`
584: @*/
585: PetscErrorCode PCApplySymmetricLeft(PC pc, Vec x, Vec y)
586: {
587:   PetscFunctionBegin;
591:   PetscCheck(x != y, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "x and y must be different vectors");
592:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(x, 2, PETSC_TRUE));
593:   PetscCall(PCSetUp(pc));
594:   PetscCall(VecLockReadPush(x));
595:   PetscCall(PetscLogEventBegin(PC_ApplySymmetricLeft, pc, x, y, 0));
596:   PetscUseTypeMethod(pc, applysymmetricleft, x, y);
597:   PetscCall(PetscLogEventEnd(PC_ApplySymmetricLeft, pc, x, y, 0));
598:   PetscCall(VecLockReadPop(x));
599:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(y, 3, PETSC_FALSE));
600:   PetscFunctionReturn(PETSC_SUCCESS);
601: }

603: /*@
604:   PCApplySymmetricRight - Applies the right part of a symmetric preconditioner to a vector.

606:   Collective

608:   Input Parameters:
609: + pc - the preconditioner context
610: - x  - input vector

612:   Output Parameter:
613: . y - output vector

615:   Level: developer

617:   Note:
618:   Currently, this routine is implemented only for `PCICC` and `PCJACOBI` preconditioners.

620: .seealso: [](ch_ksp), `PC`, `PCApply()`, `PCApplySymmetricLeft()`
621: @*/
622: PetscErrorCode PCApplySymmetricRight(PC pc, Vec x, Vec y)
623: {
624:   PetscFunctionBegin;
628:   PetscCheck(x != y, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "x and y must be different vectors");
629:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(x, 2, PETSC_TRUE));
630:   PetscCall(PCSetUp(pc));
631:   PetscCall(VecLockReadPush(x));
632:   PetscCall(PetscLogEventBegin(PC_ApplySymmetricRight, pc, x, y, 0));
633:   PetscUseTypeMethod(pc, applysymmetricright, x, y);
634:   PetscCall(PetscLogEventEnd(PC_ApplySymmetricRight, pc, x, y, 0));
635:   PetscCall(VecLockReadPop(x));
636:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(y, 3, PETSC_FALSE));
637:   PetscFunctionReturn(PETSC_SUCCESS);
638: }

640: /*@
641:   PCApplyTranspose - Applies the transpose of preconditioner to a vector.

643:   Collective

645:   Input Parameters:
646: + pc - the preconditioner context
647: - x  - input vector

649:   Output Parameter:
650: . y - output vector

652:   Level: developer

654:   Note:
655:   For complex numbers this applies the non-Hermitian transpose.

657:   Developer Note:
658:   We need to implement a `PCApplyHermitianTranspose()`

660: .seealso: [](ch_ksp), `PC`, `PCApply()`, `PCApplyBAorAB()`, `PCApplyBAorABTranspose()`, `PCApplyTransposeExists()`
661: @*/
662: PetscErrorCode PCApplyTranspose(PC pc, Vec x, Vec y)
663: {
664:   PetscFunctionBegin;
668:   PetscCheck(x != y, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "x and y must be different vectors");
669:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(x, 2, PETSC_TRUE));
670:   PetscCall(PCSetUp(pc));
671:   PetscCall(VecLockReadPush(x));
672:   PetscCall(PetscLogEventBegin(PC_Apply, pc, x, y, 0));
673:   PetscUseTypeMethod(pc, applytranspose, x, y);
674:   PetscCall(PetscLogEventEnd(PC_Apply, pc, x, y, 0));
675:   PetscCall(VecLockReadPop(x));
676:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(y, 3, PETSC_FALSE));
677:   PetscFunctionReturn(PETSC_SUCCESS);
678: }

680: /*@
681:   PCApplyTransposeExists - Test whether the preconditioner has a transpose apply operation

683:   Collective

685:   Input Parameter:
686: . pc - the preconditioner context

688:   Output Parameter:
689: . flg - `PETSC_TRUE` if a transpose operation is defined

691:   Level: developer

693: .seealso: [](ch_ksp), `PC`, `PCApplyTranspose()`
694: @*/
695: PetscErrorCode PCApplyTransposeExists(PC pc, PetscBool *flg)
696: {
697:   PetscFunctionBegin;
699:   PetscAssertPointer(flg, 2);
700:   if (pc->ops->applytranspose) *flg = PETSC_TRUE;
701:   else *flg = PETSC_FALSE;
702:   PetscFunctionReturn(PETSC_SUCCESS);
703: }

705: /*@
706:   PCApplyBAorAB - Applies the preconditioner and operator to a vector. $y = B*A*x $ or $ y = A*B*x$.

708:   Collective

710:   Input Parameters:
711: + pc   - the preconditioner context
712: . side - indicates the preconditioner side, one of `PC_LEFT`, `PC_RIGHT`, or `PC_SYMMETRIC`
713: . x    - input vector
714: - work - work vector

716:   Output Parameter:
717: . y - output vector

719:   Level: developer

721:   Note:
722:   If the `PC` has had `PCSetDiagonalScale()` set then $ D M A D^{-1} $ for left preconditioning or $ D A M D^{-1} $ is actually applied.
723:   The specific `KSPSolve()` method must also be written to handle the post-solve "correction" for the diagonal scaling.

725: .seealso: [](ch_ksp), `PC`, `PCApply()`, `PCApplyTranspose()`, `PCApplyBAorABTranspose()`
726: @*/
727: PetscErrorCode PCApplyBAorAB(PC pc, PCSide side, Vec x, Vec y, Vec work)
728: {
729:   PetscFunctionBegin;
735:   PetscCheckSameComm(pc, 1, x, 3);
736:   PetscCheckSameComm(pc, 1, y, 4);
737:   PetscCheckSameComm(pc, 1, work, 5);
738:   PetscCheck(x != y, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "x and y must be different vectors");
739:   PetscCheck(side == PC_LEFT || side == PC_SYMMETRIC || side == PC_RIGHT, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Side must be right, left, or symmetric");
740:   PetscCheck(!pc->diagonalscale || side != PC_SYMMETRIC, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Cannot include diagonal scaling with symmetric preconditioner application");
741:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(x, 3, PETSC_TRUE));

743:   PetscCall(PCSetUp(pc));
744:   if (pc->diagonalscale) {
745:     if (pc->ops->applyBA) {
746:       Vec work2; /* this is expensive, but to fix requires a second work vector argument to PCApplyBAorAB() */
747:       PetscCall(VecDuplicate(x, &work2));
748:       PetscCall(PCDiagonalScaleRight(pc, x, work2));
749:       PetscUseTypeMethod(pc, applyBA, side, work2, y, work);
750:       PetscCall(PCDiagonalScaleLeft(pc, y, y));
751:       PetscCall(VecDestroy(&work2));
752:     } else if (side == PC_RIGHT) {
753:       PetscCall(PCDiagonalScaleRight(pc, x, y));
754:       PetscCall(PCApply(pc, y, work));
755:       PetscCall(MatMult(pc->mat, work, y));
756:       PetscCall(PCDiagonalScaleLeft(pc, y, y));
757:     } else if (side == PC_LEFT) {
758:       PetscCall(PCDiagonalScaleRight(pc, x, y));
759:       PetscCall(MatMult(pc->mat, y, work));
760:       PetscCall(PCApply(pc, work, y));
761:       PetscCall(PCDiagonalScaleLeft(pc, y, y));
762:     } else PetscCheck(side != PC_SYMMETRIC, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Cannot provide diagonal scaling with symmetric application of preconditioner");
763:   } else {
764:     if (pc->ops->applyBA) {
765:       PetscUseTypeMethod(pc, applyBA, side, x, y, work);
766:     } else if (side == PC_RIGHT) {
767:       PetscCall(PCApply(pc, x, work));
768:       PetscCall(MatMult(pc->mat, work, y));
769:     } else if (side == PC_LEFT) {
770:       PetscCall(MatMult(pc->mat, x, work));
771:       PetscCall(PCApply(pc, work, y));
772:     } else if (side == PC_SYMMETRIC) {
773:       /* There's an extra copy here; maybe should provide 2 work vectors instead? */
774:       PetscCall(PCApplySymmetricRight(pc, x, work));
775:       PetscCall(MatMult(pc->mat, work, y));
776:       PetscCall(VecCopy(y, work));
777:       PetscCall(PCApplySymmetricLeft(pc, work, y));
778:     }
779:   }
780:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(y, 4, PETSC_FALSE));
781:   PetscFunctionReturn(PETSC_SUCCESS);
782: }

784: /*@
785:   PCApplyBAorABTranspose - Applies the transpose of the preconditioner
786:   and operator to a vector. That is, applies $B^T * A^T$ with left preconditioning,
787:   NOT $(B*A)^T = A^T*B^T$.

789:   Collective

791:   Input Parameters:
792: + pc   - the preconditioner context
793: . side - indicates the preconditioner side, one of `PC_LEFT`, `PC_RIGHT`, or `PC_SYMMETRIC`
794: . x    - input vector
795: - work - work vector

797:   Output Parameter:
798: . y - output vector

800:   Level: developer

802:   Note:
803:   This routine is used internally so that the same Krylov code can be used to solve $A x = b$ and $A^T x = b$, with a preconditioner
804:   defined by $B^T$. This is why this has the funny form that it computes $B^T * A^T$

806: .seealso: [](ch_ksp), `PC`, `PCApply()`, `PCApplyTranspose()`, `PCApplyBAorAB()`
807: @*/
808: PetscErrorCode PCApplyBAorABTranspose(PC pc, PCSide side, Vec x, Vec y, Vec work)
809: {
810:   PetscFunctionBegin;
815:   PetscCheck(x != y, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "x and y must be different vectors");
816:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(x, 3, PETSC_TRUE));
817:   if (pc->ops->applyBAtranspose) {
818:     PetscUseTypeMethod(pc, applyBAtranspose, side, x, y, work);
819:     if (pc->erroriffailure) PetscCall(VecValidValues_Internal(y, 4, PETSC_FALSE));
820:     PetscFunctionReturn(PETSC_SUCCESS);
821:   }
822:   PetscCheck(side == PC_LEFT || side == PC_RIGHT, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_OUTOFRANGE, "Side must be right or left");

824:   PetscCall(PCSetUp(pc));
825:   if (side == PC_RIGHT) {
826:     PetscCall(PCApplyTranspose(pc, x, work));
827:     PetscCall(MatMultTranspose(pc->mat, work, y));
828:   } else if (side == PC_LEFT) {
829:     PetscCall(MatMultTranspose(pc->mat, x, work));
830:     PetscCall(PCApplyTranspose(pc, work, y));
831:   }
832:   /* add support for PC_SYMMETRIC */
833:   if (pc->erroriffailure) PetscCall(VecValidValues_Internal(y, 4, PETSC_FALSE));
834:   PetscFunctionReturn(PETSC_SUCCESS);
835: }

837: /*@
838:   PCApplyRichardsonExists - Determines whether a particular preconditioner has a
839:   built-in fast application of Richardson's method.

841:   Not Collective

843:   Input Parameter:
844: . pc - the preconditioner

846:   Output Parameter:
847: . exists - `PETSC_TRUE` or `PETSC_FALSE`

849:   Level: developer

851: .seealso: [](ch_ksp), `PC`, `KSPRICHARDSON`, `PCApplyRichardson()`
852: @*/
853: PetscErrorCode PCApplyRichardsonExists(PC pc, PetscBool *exists)
854: {
855:   PetscFunctionBegin;
857:   PetscAssertPointer(exists, 2);
858:   if (pc->ops->applyrichardson) *exists = PETSC_TRUE;
859:   else *exists = PETSC_FALSE;
860:   PetscFunctionReturn(PETSC_SUCCESS);
861: }

863: /*@
864:   PCApplyRichardson - Applies several steps of Richardson iteration with
865:   the particular preconditioner. This routine is usually used by the
866:   Krylov solvers and not the application code directly.

868:   Collective

870:   Input Parameters:
871: + pc        - the preconditioner context
872: . b         - the right-hand side
873: . w         - one work vector
874: . rtol      - relative decrease in residual norm convergence criteria
875: . abstol    - absolute residual norm convergence criteria
876: . dtol      - divergence residual norm increase criteria
877: . its       - the number of iterations to apply.
878: - guesszero - if the input x contains nonzero initial guess

880:   Output Parameters:
881: + outits - number of iterations actually used (for SOR this always equals its)
882: . reason - the reason the apply terminated
883: - y      - the solution (also contains initial guess if guesszero is `PETSC_FALSE`

885:   Level: developer

887:   Notes:
888:   Most preconditioners do not support this function. Use the command
889:   `PCApplyRichardsonExists()` to determine if one does.

891:   Except for the `PCMG` this routine ignores the convergence tolerances
892:   and always runs for the number of iterations

894: .seealso: [](ch_ksp), `PC`, `PCApplyRichardsonExists()`
895: @*/
896: PetscErrorCode PCApplyRichardson(PC pc, Vec b, Vec y, Vec w, PetscReal rtol, PetscReal abstol, PetscReal dtol, PetscInt its, PetscBool guesszero, PetscInt *outits, PCRichardsonConvergedReason *reason)
897: {
898:   PetscFunctionBegin;
903:   PetscCheck(b != y, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "b and y must be different vectors");
904:   PetscCall(PCSetUp(pc));
905:   PetscUseTypeMethod(pc, applyrichardson, b, y, w, rtol, abstol, dtol, its, guesszero, outits, reason);
906:   PetscFunctionReturn(PETSC_SUCCESS);
907: }

909: /*@
910:   PCSetFailedReason - Sets the reason a `PCSetUp()` failed or `PC_NOERROR` if it did not fail

912:   Logically Collective

914:   Input Parameters:
915: + pc     - the preconditioner context
916: - reason - the reason it failedx

918:   Level: advanced

920: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCApply()`, `PCDestroy()`, `PCFailedReason`
921: @*/
922: PetscErrorCode PCSetFailedReason(PC pc, PCFailedReason reason)
923: {
924:   PetscFunctionBegin;
926:   pc->failedreason = reason;
927:   PetscFunctionReturn(PETSC_SUCCESS);
928: }

930: /*@
931:   PCGetFailedReason - Gets the reason a `PCSetUp()` failed or `PC_NOERROR` if it did not fail

933:   Logically Collective

935:   Input Parameter:
936: . pc - the preconditioner context

938:   Output Parameter:
939: . reason - the reason it failed

941:   Level: advanced

943:   Note:
944:   This is the maximum over reason over all ranks in the PC communicator. It is only valid after
945:   a call `KSPCheckDot()` or  `KSPCheckNorm()` inside a `KSPSolve()` or `PCReduceFailedReason()`.
946:   It is not valid immediately after a `PCSetUp()` or `PCApply()`, then use `PCGetFailedReasonRank()`

948: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCApply()`, `PCDestroy()`, `PCGetFailedReasonRank()`, `PCSetFailedReason()`
949: @*/
950: PetscErrorCode PCGetFailedReason(PC pc, PCFailedReason *reason)
951: {
952:   PetscFunctionBegin;
954:   if (pc->setupcalled < 0) *reason = (PCFailedReason)pc->setupcalled;
955:   else *reason = pc->failedreason;
956:   PetscFunctionReturn(PETSC_SUCCESS);
957: }

959: /*@
960:   PCGetFailedReasonRank - Gets the reason a `PCSetUp()` failed or `PC_NOERROR` if it did not fail on this MPI rank

962:   Not Collective

964:   Input Parameter:
965: . pc - the preconditioner context

967:   Output Parameter:
968: . reason - the reason it failed

970:   Level: advanced

972:   Note:
973:   Different processes may have different reasons or no reason, see `PCGetFailedReason()`

975: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCApply()`, `PCDestroy()`, `PCGetFailedReason()`, `PCSetFailedReason()`, `PCReduceFailedReason()`
976: @*/
977: PetscErrorCode PCGetFailedReasonRank(PC pc, PCFailedReason *reason)
978: {
979:   PetscFunctionBegin;
981:   if (pc->setupcalled < 0) *reason = (PCFailedReason)pc->setupcalled;
982:   else *reason = pc->failedreason;
983:   PetscFunctionReturn(PETSC_SUCCESS);
984: }

986: /*@
987:   PCReduceFailedReason - Reduce the failed reason among the MPI processes that share the `PC`

989:   Collective

991:   Input Parameter:
992: . pc - the preconditioner context

994:   Level: advanced

996:   Note:
997:   Different MPI processes may have different reasons or no reason, see `PCGetFailedReason()`. This routine
998:   makes them have a common value (failure if any MPI process had a failure).

1000: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCApply()`, `PCDestroy()`, `PCGetFailedReason()`, `PCSetFailedReason()`
1001: @*/
1002: PetscErrorCode PCReduceFailedReason(PC pc)
1003: {
1004:   PetscInt buf;

1006:   PetscFunctionBegin;
1008:   buf = (PetscInt)pc->failedreason;
1009:   PetscCall(MPIU_Allreduce(MPI_IN_PLACE, &buf, 1, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject)pc)));
1010:   pc->failedreason = (PCFailedReason)buf;
1011:   PetscFunctionReturn(PETSC_SUCCESS);
1012: }

1014: /*  Next line needed to deactivate KSP_Solve logging */
1015: #include <petsc/private/kspimpl.h>

1017: /*
1018:       a setupcall of 0 indicates never setup,
1019:                      1 indicates has been previously setup
1020:                     -1 indicates a PCSetUp() was attempted and failed
1021: */
1022: /*@
1023:   PCSetUp - Prepares for the use of a preconditioner.

1025:   Collective

1027:   Input Parameter:
1028: . pc - the preconditioner context

1030:   Level: developer

1032: .seealso: [](ch_ksp), `PC`, `PCCreate()`, `PCApply()`, `PCDestroy()`
1033: @*/
1034: PetscErrorCode PCSetUp(PC pc)
1035: {
1036:   const char      *def;
1037:   PetscObjectState matstate, matnonzerostate;

1039:   PetscFunctionBegin;
1041:   PetscCheck(pc->mat, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "Matrix must be set first");

1043:   if (pc->setupcalled && pc->reusepreconditioner) {
1044:     PetscCall(PetscInfo(pc, "Leaving PC with identical preconditioner since reuse preconditioner is set\n"));
1045:     PetscFunctionReturn(PETSC_SUCCESS);
1046:   }

1048:   PetscCall(PetscObjectStateGet((PetscObject)pc->pmat, &matstate));
1049:   PetscCall(MatGetNonzeroState(pc->pmat, &matnonzerostate));
1050:   if (!pc->setupcalled) {
1051:     //PetscCall(PetscInfo(pc, "Setting up PC for first time\n"));
1052:     pc->flag = DIFFERENT_NONZERO_PATTERN;
1053:   } else if (matstate == pc->matstate) PetscFunctionReturn(PETSC_SUCCESS);
1054:   else {
1055:     if (matnonzerostate != pc->matnonzerostate) {
1056:       PetscCall(PetscInfo(pc, "Setting up PC with different nonzero pattern\n"));
1057:       pc->flag = DIFFERENT_NONZERO_PATTERN;
1058:     } else {
1059:       //PetscCall(PetscInfo(pc, "Setting up PC with same nonzero pattern\n"));
1060:       pc->flag = SAME_NONZERO_PATTERN;
1061:     }
1062:   }
1063:   pc->matstate        = matstate;
1064:   pc->matnonzerostate = matnonzerostate;

1066:   if (!((PetscObject)pc)->type_name) {
1067:     PetscCall(PCGetDefaultType_Private(pc, &def));
1068:     PetscCall(PCSetType(pc, def));
1069:   }

1071:   PetscCall(MatSetErrorIfFailure(pc->pmat, pc->erroriffailure));
1072:   PetscCall(MatSetErrorIfFailure(pc->mat, pc->erroriffailure));
1073:   PetscCall(PetscLogEventBegin(PC_SetUp, pc, 0, 0, 0));
1074:   if (pc->ops->setup) {
1075:     /* do not log solves and applications of preconditioners while constructing preconditioners; perhaps they should be logged separately from the regular solves */
1076:     PetscCall(KSPInitializePackage());
1077:     PetscCall(PetscLogEventDeactivatePush(KSP_Solve));
1078:     PetscCall(PetscLogEventDeactivatePush(PC_Apply));
1079:     PetscUseTypeMethod(pc, setup);
1080:     PetscCall(PetscLogEventDeactivatePop(KSP_Solve));
1081:     PetscCall(PetscLogEventDeactivatePop(PC_Apply));
1082:   }
1083:   PetscCall(PetscLogEventEnd(PC_SetUp, pc, 0, 0, 0));
1084:   if (!pc->setupcalled) pc->setupcalled = 1;
1085:   PetscFunctionReturn(PETSC_SUCCESS);
1086: }

1088: /*@
1089:   PCSetUpOnBlocks - Sets up the preconditioner for each block in
1090:   the block Jacobi, block Gauss-Seidel, and overlapping Schwarz
1091:   methods.

1093:   Collective

1095:   Input Parameter:
1096: . pc - the preconditioner context

1098:   Level: developer

1100:   Note:
1101:   For nested preconditioners such as `PCBJACOBI` `PCSetUp()` is not called on each sub-`KSP` when `PCSetUp()` is
1102:   called on the outer `PC`, this routine ensures it is called.

1104: .seealso: [](ch_ksp), `PC`, `PCSetUp()`, `PCCreate()`, `PCApply()`, `PCDestroy()`
1105: @*/
1106: PetscErrorCode PCSetUpOnBlocks(PC pc)
1107: {
1108:   PetscFunctionBegin;
1110:   if (!pc->ops->setuponblocks) PetscFunctionReturn(PETSC_SUCCESS);
1111:   PetscCall(PetscLogEventBegin(PC_SetUpOnBlocks, pc, 0, 0, 0));
1112:   PetscUseTypeMethod(pc, setuponblocks);
1113:   PetscCall(PetscLogEventEnd(PC_SetUpOnBlocks, pc, 0, 0, 0));
1114:   PetscFunctionReturn(PETSC_SUCCESS);
1115: }

1117: /*@C
1118:   PCSetModifySubMatrices - Sets a user-defined routine for modifying the
1119:   submatrices that arise within certain subdomain-based preconditioners such as `PCASM`

1121:   Logically Collective

1123:   Input Parameters:
1124: + pc   - the preconditioner context
1125: . func - routine for modifying the submatrices
1126: - ctx  - optional user-defined context (may be `NULL`)

1128:   Calling sequence of `func`:
1129: + pc     - the preconditioner context
1130: . nsub   - number of index sets
1131: . row    - an array of index sets that contain the global row numbers
1132:          that comprise each local submatrix
1133: . col    - an array of index sets that contain the global column numbers
1134:          that comprise each local submatrix
1135: . submat - array of local submatrices
1136: - ctx    - optional user-defined context for private data for the
1137:          user-defined func routine (may be `NULL`)

1139:   Level: advanced

1141:   Notes:
1142:   The basic submatrices are extracted from the preconditioner matrix as
1143:   usual; the user can then alter these (for example, to set different boundary
1144:   conditions for each submatrix) before they are used for the local solves.

1146:   `PCSetModifySubMatrices()` MUST be called before `KSPSetUp()` and
1147:   `KSPSolve()`.

1149:   A routine set by `PCSetModifySubMatrices()` is currently called within
1150:   the block Jacobi (`PCBJACOBI`) and additive Schwarz (`PCASM`)
1151:   preconditioners.  All other preconditioners ignore this routine.

1153: .seealso: [](ch_ksp), `PC`, `PCBJACOBI`, `PCASM`, `PCModifySubMatrices()`
1154: @*/
1155: PetscErrorCode PCSetModifySubMatrices(PC pc, PetscErrorCode (*func)(PC pc, PetscInt nsub, const IS row[], const IS col[], Mat submat[], void *ctx), void *ctx)
1156: {
1157:   PetscFunctionBegin;
1159:   pc->modifysubmatrices  = func;
1160:   pc->modifysubmatricesP = ctx;
1161:   PetscFunctionReturn(PETSC_SUCCESS);
1162: }

1164: /*@C
1165:   PCModifySubMatrices - Calls an optional user-defined routine within
1166:   certain preconditioners if one has been set with `PCSetModifySubMatrices()`.

1168:   Collective

1170:   Input Parameters:
1171: + pc     - the preconditioner context
1172: . nsub   - the number of local submatrices
1173: . row    - an array of index sets that contain the global row numbers
1174:          that comprise each local submatrix
1175: . col    - an array of index sets that contain the global column numbers
1176:          that comprise each local submatrix
1177: . submat - array of local submatrices
1178: - ctx    - optional user-defined context for private data for the
1179:          user-defined routine (may be `NULL`)

1181:   Output Parameter:
1182: . submat - array of local submatrices (the entries of which may
1183:             have been modified)

1185:   Level: developer

1187:   Note:
1188:   The user should NOT generally call this routine, as it will
1189:   automatically be called within certain preconditioners.

1191: .seealso: [](ch_ksp), `PC`, `PCSetModifySubMatrices()`
1192: @*/
1193: PetscErrorCode PCModifySubMatrices(PC pc, PetscInt nsub, const IS row[], const IS col[], Mat submat[], void *ctx)
1194: {
1195:   PetscFunctionBegin;
1197:   if (!pc->modifysubmatrices) PetscFunctionReturn(PETSC_SUCCESS);
1198:   PetscCall(PetscLogEventBegin(PC_ModifySubMatrices, pc, 0, 0, 0));
1199:   PetscCall((*pc->modifysubmatrices)(pc, nsub, row, col, submat, ctx));
1200:   PetscCall(PetscLogEventEnd(PC_ModifySubMatrices, pc, 0, 0, 0));
1201:   PetscFunctionReturn(PETSC_SUCCESS);
1202: }

1204: /*@
1205:   PCSetOperators - Sets the matrix associated with the linear system and
1206:   a (possibly) different one associated with the preconditioner.

1208:   Logically Collective

1210:   Input Parameters:
1211: + pc   - the preconditioner context
1212: . Amat - the matrix that defines the linear system
1213: - Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat.

1215:   Level: intermediate

1217:   Notes:
1218:   Passing a `NULL` for `Amat` or `Pmat` removes the matrix that is currently used.

1220:   If you wish to replace either `Amat` or `Pmat` but leave the other one untouched then
1221:   first call `KSPGetOperators()` to get the one you wish to keep, call `PetscObjectReference()`
1222:   on it and then pass it back in in your call to `KSPSetOperators()`.

1224:   More Notes about Repeated Solution of Linear Systems:
1225:   PETSc does NOT reset the matrix entries of either `Amat` or `Pmat`
1226:   to zero after a linear solve; the user is completely responsible for
1227:   matrix assembly.  See the routine `MatZeroEntries()` if desiring to
1228:   zero all elements of a matrix.

1230: .seealso: [](ch_ksp), `PC`, `PCGetOperators()`, `MatZeroEntries()`
1231:  @*/
1232: PetscErrorCode PCSetOperators(PC pc, Mat Amat, Mat Pmat)
1233: {
1234:   PetscInt m1, n1, m2, n2;

1236:   PetscFunctionBegin;
1240:   if (Amat) PetscCheckSameComm(pc, 1, Amat, 2);
1241:   if (Pmat) PetscCheckSameComm(pc, 1, Pmat, 3);
1242:   if (pc->setupcalled && pc->mat && pc->pmat && Amat && Pmat) {
1243:     PetscCall(MatGetLocalSize(Amat, &m1, &n1));
1244:     PetscCall(MatGetLocalSize(pc->mat, &m2, &n2));
1245:     PetscCheck(m1 == m2 && n1 == n2, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Cannot change local size of Amat after use old sizes %" PetscInt_FMT " %" PetscInt_FMT " new sizes %" PetscInt_FMT " %" PetscInt_FMT, m2, n2, m1, n1);
1246:     PetscCall(MatGetLocalSize(Pmat, &m1, &n1));
1247:     PetscCall(MatGetLocalSize(pc->pmat, &m2, &n2));
1248:     PetscCheck(m1 == m2 && n1 == n2, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Cannot change local size of Pmat after use old sizes %" PetscInt_FMT " %" PetscInt_FMT " new sizes %" PetscInt_FMT " %" PetscInt_FMT, m2, n2, m1, n1);
1249:   }

1251:   if (Pmat != pc->pmat) {
1252:     /* changing the operator that defines the preconditioner thus reneed to clear current states so new preconditioner is built */
1253:     pc->matnonzerostate = -1;
1254:     pc->matstate        = -1;
1255:   }

1257:   /* reference first in case the matrices are the same */
1258:   if (Amat) PetscCall(PetscObjectReference((PetscObject)Amat));
1259:   PetscCall(MatDestroy(&pc->mat));
1260:   if (Pmat) PetscCall(PetscObjectReference((PetscObject)Pmat));
1261:   PetscCall(MatDestroy(&pc->pmat));
1262:   pc->mat  = Amat;
1263:   pc->pmat = Pmat;
1264:   PetscFunctionReturn(PETSC_SUCCESS);
1265: }

1267: /*@
1268:   PCSetReusePreconditioner - reuse the current preconditioner even if the operator in the preconditioner has changed.

1270:   Logically Collective

1272:   Input Parameters:
1273: + pc   - the preconditioner context
1274: - flag - `PETSC_TRUE` do not compute a new preconditioner, `PETSC_FALSE` do compute a new preconditioner

1276:   Level: intermediate

1278:   Note:
1279:   Normally if a matrix inside a `PC` changes the `PC` automatically updates itself using information from the changed matrix. This option
1280:   prevents this.

1282: .seealso: [](ch_ksp), `PC`, `PCGetOperators()`, `MatZeroEntries()`, `PCGetReusePreconditioner()`, `KSPSetReusePreconditioner()`
1283:  @*/
1284: PetscErrorCode PCSetReusePreconditioner(PC pc, PetscBool flag)
1285: {
1286:   PetscFunctionBegin;
1289:   pc->reusepreconditioner = flag;
1290:   PetscFunctionReturn(PETSC_SUCCESS);
1291: }

1293: /*@
1294:   PCGetReusePreconditioner - Determines if the `PC` reuses the current preconditioner even if the operator in the preconditioner has changed.

1296:   Not Collective

1298:   Input Parameter:
1299: . pc - the preconditioner context

1301:   Output Parameter:
1302: . flag - `PETSC_TRUE` do not compute a new preconditioner, `PETSC_FALSE` do compute a new preconditioner

1304:   Level: intermediate

1306: .seealso: [](ch_ksp), `PC`, `PCGetOperators()`, `MatZeroEntries()`, `PCSetReusePreconditioner()`
1307:  @*/
1308: PetscErrorCode PCGetReusePreconditioner(PC pc, PetscBool *flag)
1309: {
1310:   PetscFunctionBegin;
1312:   PetscAssertPointer(flag, 2);
1313:   *flag = pc->reusepreconditioner;
1314:   PetscFunctionReturn(PETSC_SUCCESS);
1315: }

1317: /*@
1318:   PCGetOperators - Gets the matrix associated with the linear system and
1319:   possibly a different one associated with the preconditioner.

1321:   Not Collective, though parallel `Mat`s are returned if `pc` is parallel

1323:   Input Parameter:
1324: . pc - the preconditioner context

1326:   Output Parameters:
1327: + Amat - the matrix defining the linear system
1328: - Pmat - the matrix from which the preconditioner is constructed, usually the same as Amat.

1330:   Level: intermediate

1332:   Note:
1333:   Does not increase the reference count of the matrices, so you should not destroy them

1335:   Alternative usage: If the operators have NOT been set with `KSPSetOperators()`/`PCSetOperators()` then the operators
1336:   are created in `PC` and returned to the user. In this case, if both operators
1337:   mat and pmat are requested, two DIFFERENT operators will be returned. If
1338:   only one is requested both operators in the PC will be the same (i.e. as
1339:   if one had called `KSPSetOperators()`/`PCSetOperators()` with the same argument for both Mats).
1340:   The user must set the sizes of the returned matrices and their type etc just
1341:   as if the user created them with `MatCreate()`. For example,

1343: .vb
1344:          KSP/PCGetOperators(ksp/pc,&Amat,NULL); is equivalent to
1345:            set size, type, etc of Amat

1347:          MatCreate(comm,&mat);
1348:          KSP/PCSetOperators(ksp/pc,Amat,Amat);
1349:          PetscObjectDereference((PetscObject)mat);
1350:            set size, type, etc of Amat
1351: .ve

1353:   and

1355: .vb
1356:          KSP/PCGetOperators(ksp/pc,&Amat,&Pmat); is equivalent to
1357:            set size, type, etc of Amat and Pmat

1359:          MatCreate(comm,&Amat);
1360:          MatCreate(comm,&Pmat);
1361:          KSP/PCSetOperators(ksp/pc,Amat,Pmat);
1362:          PetscObjectDereference((PetscObject)Amat);
1363:          PetscObjectDereference((PetscObject)Pmat);
1364:            set size, type, etc of Amat and Pmat
1365: .ve

1367:   The rationale for this support is so that when creating a `TS`, `SNES`, or `KSP` the hierarchy
1368:   of underlying objects (i.e. `SNES`, `KSP`, `PC`, `Mat`) and their lifespans can be completely
1369:   managed by the top most level object (i.e. the `TS`, `SNES`, or `KSP`). Another way to look
1370:   at this is when you create a `SNES` you do not NEED to create a `KSP` and attach it to
1371:   the `SNES` object (the `SNES` object manages it for you). Similarly when you create a KSP
1372:   you do not need to attach a `PC` to it (the `KSP` object manages the `PC` object for you).
1373:   Thus, why should YOU have to create the `Mat` and attach it to the `SNES`/`KSP`/`PC`, when
1374:   it can be created for you?

1376: .seealso: [](ch_ksp), `PC`, `PCSetOperators()`, `KSPGetOperators()`, `KSPSetOperators()`, `PCGetOperatorsSet()`
1377: @*/
1378: PetscErrorCode PCGetOperators(PC pc, Mat *Amat, Mat *Pmat)
1379: {
1380:   PetscFunctionBegin;
1382:   if (Amat) {
1383:     if (!pc->mat) {
1384:       if (pc->pmat && !Pmat) { /* Pmat has been set, but user did not request it, so use for Amat */
1385:         pc->mat = pc->pmat;
1386:         PetscCall(PetscObjectReference((PetscObject)pc->mat));
1387:       } else { /* both Amat and Pmat are empty */
1388:         PetscCall(MatCreate(PetscObjectComm((PetscObject)pc), &pc->mat));
1389:         if (!Pmat) { /* user did NOT request Pmat, so make same as Amat */
1390:           pc->pmat = pc->mat;
1391:           PetscCall(PetscObjectReference((PetscObject)pc->pmat));
1392:         }
1393:       }
1394:     }
1395:     *Amat = pc->mat;
1396:   }
1397:   if (Pmat) {
1398:     if (!pc->pmat) {
1399:       if (pc->mat && !Amat) { /* Amat has been set but was not requested, so use for pmat */
1400:         pc->pmat = pc->mat;
1401:         PetscCall(PetscObjectReference((PetscObject)pc->pmat));
1402:       } else {
1403:         PetscCall(MatCreate(PetscObjectComm((PetscObject)pc), &pc->pmat));
1404:         if (!Amat) { /* user did NOT request Amat, so make same as Pmat */
1405:           pc->mat = pc->pmat;
1406:           PetscCall(PetscObjectReference((PetscObject)pc->mat));
1407:         }
1408:       }
1409:     }
1410:     *Pmat = pc->pmat;
1411:   }
1412:   PetscFunctionReturn(PETSC_SUCCESS);
1413: }

1415: /*@C
1416:   PCGetOperatorsSet - Determines if the matrix associated with the linear system and
1417:   possibly a different one associated with the preconditioner have been set in the `PC`.

1419:   Not Collective, though the results on all processes should be the same

1421:   Input Parameter:
1422: . pc - the preconditioner context

1424:   Output Parameters:
1425: + mat  - the matrix associated with the linear system was set
1426: - pmat - matrix associated with the preconditioner was set, usually the same

1428:   Level: intermediate

1430: .seealso: [](ch_ksp), `PC`, `PCSetOperators()`, `KSPGetOperators()`, `KSPSetOperators()`, `PCGetOperators()`
1431: @*/
1432: PetscErrorCode PCGetOperatorsSet(PC pc, PetscBool *mat, PetscBool *pmat)
1433: {
1434:   PetscFunctionBegin;
1436:   if (mat) *mat = (pc->mat) ? PETSC_TRUE : PETSC_FALSE;
1437:   if (pmat) *pmat = (pc->pmat) ? PETSC_TRUE : PETSC_FALSE;
1438:   PetscFunctionReturn(PETSC_SUCCESS);
1439: }

1441: /*@
1442:   PCFactorGetMatrix - Gets the factored matrix from the
1443:   preconditioner context.  This routine is valid only for the `PCLU`,
1444:   `PCILU`, `PCCHOLESKY`, and `PCICC` methods.

1446:   Not Collective though `mat` is parallel if `pc` is parallel

1448:   Input Parameter:
1449: . pc - the preconditioner context

1451:   Output Parameters:
1452: . mat - the factored matrix

1454:   Level: advanced

1456:   Note:
1457:   Does not increase the reference count for `mat` so DO NOT destroy it

1459: .seealso: [](ch_ksp), `PC`, `PCLU`, `PCILU`, `PCCHOLESKY`, `PCICC`
1460: @*/
1461: PetscErrorCode PCFactorGetMatrix(PC pc, Mat *mat)
1462: {
1463:   PetscFunctionBegin;
1465:   PetscAssertPointer(mat, 2);
1466:   PetscCall(PCFactorSetUpMatSolverType(pc));
1467:   PetscUseTypeMethod(pc, getfactoredmatrix, mat);
1468:   PetscFunctionReturn(PETSC_SUCCESS);
1469: }

1471: /*@
1472:   PCSetOptionsPrefix - Sets the prefix used for searching for all
1473:   `PC` options in the database.

1475:   Logically Collective

1477:   Input Parameters:
1478: + pc     - the preconditioner context
1479: - prefix - the prefix string to prepend to all `PC` option requests

1481:   Note:
1482:   A hyphen (-) must NOT be given at the beginning of the prefix name.
1483:   The first character of all runtime options is AUTOMATICALLY the
1484:   hyphen.

1486:   Level: advanced

1488: .seealso: [](ch_ksp), `PC`, `PCSetFromOptions`, `PCAppendOptionsPrefix()`, `PCGetOptionsPrefix()`
1489: @*/
1490: PetscErrorCode PCSetOptionsPrefix(PC pc, const char prefix[])
1491: {
1492:   PetscFunctionBegin;
1494:   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)pc, prefix));
1495:   PetscFunctionReturn(PETSC_SUCCESS);
1496: }

1498: /*@
1499:   PCAppendOptionsPrefix - Appends to the prefix used for searching for all
1500:   `PC` options in the database.

1502:   Logically Collective

1504:   Input Parameters:
1505: + pc     - the preconditioner context
1506: - prefix - the prefix string to prepend to all `PC` option requests

1508:   Note:
1509:   A hyphen (-) must NOT be given at the beginning of the prefix name.
1510:   The first character of all runtime options is AUTOMATICALLY the
1511:   hyphen.

1513:   Level: advanced

1515: .seealso: [](ch_ksp), `PC`, `PCSetFromOptions`, `PCSetOptionsPrefix()`, `PCGetOptionsPrefix()`
1516: @*/
1517: PetscErrorCode PCAppendOptionsPrefix(PC pc, const char prefix[])
1518: {
1519:   PetscFunctionBegin;
1521:   PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)pc, prefix));
1522:   PetscFunctionReturn(PETSC_SUCCESS);
1523: }

1525: /*@
1526:   PCGetOptionsPrefix - Gets the prefix used for searching for all
1527:   PC options in the database.

1529:   Not Collective

1531:   Input Parameter:
1532: . pc - the preconditioner context

1534:   Output Parameter:
1535: . prefix - pointer to the prefix string used, is returned

1537:   Level: advanced

1539:   Fortran Note:
1540:   The user should pass in a string `prefix` of
1541:   sufficient length to hold the prefix.

1543: .seealso: [](ch_ksp), `PC`, `PCSetFromOptions`, `PCSetOptionsPrefix()`, `PCAppendOptionsPrefix()`
1544: @*/
1545: PetscErrorCode PCGetOptionsPrefix(PC pc, const char *prefix[])
1546: {
1547:   PetscFunctionBegin;
1549:   PetscAssertPointer(prefix, 2);
1550:   PetscCall(PetscObjectGetOptionsPrefix((PetscObject)pc, prefix));
1551:   PetscFunctionReturn(PETSC_SUCCESS);
1552: }

1554: /*
1555:    Indicates the right-hand side will be changed by KSPSolve(), this occurs for a few
1556:   preconditioners including BDDC and Eisentat that transform the equations before applying
1557:   the Krylov methods
1558: */
1559: PETSC_INTERN PetscErrorCode PCPreSolveChangeRHS(PC pc, PetscBool *change)
1560: {
1561:   PetscFunctionBegin;
1563:   PetscAssertPointer(change, 2);
1564:   *change = PETSC_FALSE;
1565:   PetscTryMethod(pc, "PCPreSolveChangeRHS_C", (PC, PetscBool *), (pc, change));
1566:   PetscFunctionReturn(PETSC_SUCCESS);
1567: }

1569: /*@
1570:   PCPreSolve - Optional pre-solve phase, intended for any
1571:   preconditioner-specific actions that must be performed before
1572:   the iterative solve itself.

1574:   Collective

1576:   Input Parameters:
1577: + pc  - the preconditioner context
1578: - ksp - the Krylov subspace context

1580:   Level: developer

1582:   Example Usage:
1583: .vb
1584:     PCPreSolve(pc,ksp);
1585:     KSPSolve(ksp,b,x);
1586:     PCPostSolve(pc,ksp);
1587: .ve

1589:   Notes:
1590:   The pre-solve phase is distinct from the `PCSetUp()` phase.

1592:   `KSPSolve()` calls this directly, so is rarely called by the user.

1594: .seealso: [](ch_ksp), `PC`, `PCPostSolve()`
1595: @*/
1596: PetscErrorCode PCPreSolve(PC pc, KSP ksp)
1597: {
1598:   Vec x, rhs;

1600:   PetscFunctionBegin;
1603:   pc->presolvedone++;
1604:   PetscCheck(pc->presolvedone <= 2, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Cannot embed PCPreSolve() more than twice");
1605:   PetscCall(KSPGetSolution(ksp, &x));
1606:   PetscCall(KSPGetRhs(ksp, &rhs));

1608:   if (pc->ops->presolve) PetscUseTypeMethod(pc, presolve, ksp, rhs, x);
1609:   else if (pc->presolve) PetscCall(pc->presolve(pc, ksp));
1610:   PetscFunctionReturn(PETSC_SUCCESS);
1611: }

1613: /*@C
1614:   PCSetPreSolve - Sets function used by `PCPreSolve()` which is intended for any
1615:   preconditioner-specific actions that must be performed before
1616:   the iterative solve itself.

1618:   Logically Collective

1620:   Input Parameters:
1621: + pc       - the preconditioner object
1622: - presolve - the function to call before the solve

1624:   Calling sequence of `presolve`:
1625: + pc  - the `PC` context
1626: - ksp - the `KSP` context

1628:   Level: developer

1630: .seealso: [](ch_ksp), `PC`, `PCSetUp()`, `PCPreSolve()`
1631: @*/
1632: PetscErrorCode PCSetPreSolve(PC pc, PetscErrorCode (*presolve)(PC pc, KSP ksp))
1633: {
1634:   PetscFunctionBegin;
1636:   pc->presolve = presolve;
1637:   PetscFunctionReturn(PETSC_SUCCESS);
1638: }

1640: /*@
1641:   PCPostSolve - Optional post-solve phase, intended for any
1642:   preconditioner-specific actions that must be performed after
1643:   the iterative solve itself.

1645:   Collective

1647:   Input Parameters:
1648: + pc  - the preconditioner context
1649: - ksp - the Krylov subspace context

1651:   Example Usage:
1652: .vb
1653:     PCPreSolve(pc,ksp);
1654:     KSPSolve(ksp,b,x);
1655:     PCPostSolve(pc,ksp);
1656: .ve

1658:   Level: developer

1660:   Note:
1661:   `KSPSolve()` calls this routine directly, so it is rarely called by the user.

1663: .seealso: [](ch_ksp), `PC`, `PCSetPostSolve()`, `PCSetPresolve()`, `PCPreSolve()`, `KSPSolve()`
1664: @*/
1665: PetscErrorCode PCPostSolve(PC pc, KSP ksp)
1666: {
1667:   Vec x, rhs;

1669:   PetscFunctionBegin;
1672:   pc->presolvedone--;
1673:   PetscCall(KSPGetSolution(ksp, &x));
1674:   PetscCall(KSPGetRhs(ksp, &rhs));
1675:   PetscTryTypeMethod(pc, postsolve, ksp, rhs, x);
1676:   PetscFunctionReturn(PETSC_SUCCESS);
1677: }

1679: /*@
1680:   PCLoad - Loads a `PC` that has been stored in binary  with `PCView()`.

1682:   Collective

1684:   Input Parameters:
1685: + newdm  - the newly loaded `PC`, this needs to have been created with `PCCreate()` or
1686:            some related function before a call to `PCLoad()`.
1687: - viewer - binary file viewer, obtained from `PetscViewerBinaryOpen()`

1689:   Level: intermediate

1691:   Note:
1692:   The type is determined by the data in the file, any `PCType` set into the `PC` before this call is ignored.

1694: .seealso: [](ch_ksp), `PC`, `PetscViewerBinaryOpen()`, `PCView()`, `MatLoad()`, `VecLoad()`
1695: @*/
1696: PetscErrorCode PCLoad(PC newdm, PetscViewer viewer)
1697: {
1698:   PetscBool isbinary;
1699:   PetscInt  classid;
1700:   char      type[256];

1702:   PetscFunctionBegin;
1705:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
1706:   PetscCheck(isbinary, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid viewer; open viewer with PetscViewerBinaryOpen()");

1708:   PetscCall(PetscViewerBinaryRead(viewer, &classid, 1, NULL, PETSC_INT));
1709:   PetscCheck(classid == PC_FILE_CLASSID, PetscObjectComm((PetscObject)newdm), PETSC_ERR_ARG_WRONG, "Not PC next in file");
1710:   PetscCall(PetscViewerBinaryRead(viewer, type, 256, NULL, PETSC_CHAR));
1711:   PetscCall(PCSetType(newdm, type));
1712:   PetscTryTypeMethod(newdm, load, viewer);
1713:   PetscFunctionReturn(PETSC_SUCCESS);
1714: }

1716: #include <petscdraw.h>
1717: #if defined(PETSC_HAVE_SAWS)
1718: #include <petscviewersaws.h>
1719: #endif

1721: /*@
1722:   PCViewFromOptions - View from the `PC` based on options in the options database

1724:   Collective

1726:   Input Parameters:
1727: + A    - the `PC` context
1728: . obj  - Optional object that provides the options prefix
1729: - name - command line option

1731:   Level: intermediate

1733: .seealso: [](ch_ksp), `PC`, `PCView`, `PetscObjectViewFromOptions()`, `PCCreate()`
1734: @*/
1735: PetscErrorCode PCViewFromOptions(PC A, PetscObject obj, const char name[])
1736: {
1737:   PetscFunctionBegin;
1739:   PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
1740:   PetscFunctionReturn(PETSC_SUCCESS);
1741: }

1743: /*@
1744:   PCView - Prints information about the `PC`

1746:   Collective

1748:   Input Parameters:
1749: + pc     - the `PC` context
1750: - viewer - optional visualization context

1752:   Level: developer

1754:   Notes:
1755:   The available visualization contexts include
1756: +     `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
1757: -     `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard
1758:   output where only the first processor opens
1759:   the file. All other processors send their
1760:   data to the first processor to print.

1762:   The user can open an alternative visualization contexts with
1763:   `PetscViewerASCIIOpen()` (output to a specified file).

1765: .seealso: [](ch_ksp), `PC`, `PetscViewer`, `KSPView()`, `PetscViewerASCIIOpen()`
1766: @*/
1767: PetscErrorCode PCView(PC pc, PetscViewer viewer)
1768: {
1769:   PCType            cstr;
1770:   PetscViewerFormat format;
1771:   PetscBool         iascii, isstring, isbinary, isdraw, pop = PETSC_FALSE;
1772: #if defined(PETSC_HAVE_SAWS)
1773:   PetscBool issaws;
1774: #endif

1776:   PetscFunctionBegin;
1778:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pc), &viewer));
1780:   PetscCheckSameComm(pc, 1, viewer, 2);

1782:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
1783:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring));
1784:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
1785:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
1786: #if defined(PETSC_HAVE_SAWS)
1787:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSAWS, &issaws));
1788: #endif

1790:   if (iascii) {
1791:     PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)pc, viewer));
1792:     if (!pc->setupcalled) PetscCall(PetscViewerASCIIPrintf(viewer, "  PC has not been set up so information may be incomplete\n"));
1793:     PetscCall(PetscViewerASCIIPushTab(viewer));
1794:     PetscTryTypeMethod(pc, view, viewer);
1795:     PetscCall(PetscViewerASCIIPopTab(viewer));
1796:     if (pc->mat) {
1797:       PetscCall(PetscViewerGetFormat(viewer, &format));
1798:       if (format != PETSC_VIEWER_ASCII_INFO_DETAIL) {
1799:         PetscCall(PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_INFO));
1800:         pop = PETSC_TRUE;
1801:       }
1802:       if (pc->pmat == pc->mat) {
1803:         PetscCall(PetscViewerASCIIPrintf(viewer, "  linear system matrix = precond matrix:\n"));
1804:         PetscCall(PetscViewerASCIIPushTab(viewer));
1805:         PetscCall(MatView(pc->mat, viewer));
1806:         PetscCall(PetscViewerASCIIPopTab(viewer));
1807:       } else {
1808:         if (pc->pmat) {
1809:           PetscCall(PetscViewerASCIIPrintf(viewer, "  linear system matrix followed by preconditioner matrix:\n"));
1810:         } else {
1811:           PetscCall(PetscViewerASCIIPrintf(viewer, "  linear system matrix:\n"));
1812:         }
1813:         PetscCall(PetscViewerASCIIPushTab(viewer));
1814:         PetscCall(MatView(pc->mat, viewer));
1815:         if (pc->pmat) PetscCall(MatView(pc->pmat, viewer));
1816:         PetscCall(PetscViewerASCIIPopTab(viewer));
1817:       }
1818:       if (pop) PetscCall(PetscViewerPopFormat(viewer));
1819:     }
1820:   } else if (isstring) {
1821:     PetscCall(PCGetType(pc, &cstr));
1822:     PetscCall(PetscViewerStringSPrintf(viewer, " PCType: %-7.7s", cstr));
1823:     PetscTryTypeMethod(pc, view, viewer);
1824:     if (pc->mat) PetscCall(MatView(pc->mat, viewer));
1825:     if (pc->pmat && pc->pmat != pc->mat) PetscCall(MatView(pc->pmat, viewer));
1826:   } else if (isbinary) {
1827:     PetscInt    classid = PC_FILE_CLASSID;
1828:     MPI_Comm    comm;
1829:     PetscMPIInt rank;
1830:     char        type[256];

1832:     PetscCall(PetscObjectGetComm((PetscObject)pc, &comm));
1833:     PetscCallMPI(MPI_Comm_rank(comm, &rank));
1834:     if (rank == 0) {
1835:       PetscCall(PetscViewerBinaryWrite(viewer, &classid, 1, PETSC_INT));
1836:       PetscCall(PetscStrncpy(type, ((PetscObject)pc)->type_name, 256));
1837:       PetscCall(PetscViewerBinaryWrite(viewer, type, 256, PETSC_CHAR));
1838:     }
1839:     PetscTryTypeMethod(pc, view, viewer);
1840:   } else if (isdraw) {
1841:     PetscDraw draw;
1842:     char      str[25];
1843:     PetscReal x, y, bottom, h;
1844:     PetscInt  n;

1846:     PetscCall(PetscViewerDrawGetDraw(viewer, 0, &draw));
1847:     PetscCall(PetscDrawGetCurrentPoint(draw, &x, &y));
1848:     if (pc->mat) {
1849:       PetscCall(MatGetSize(pc->mat, &n, NULL));
1850:       PetscCall(PetscSNPrintf(str, 25, "PC: %s (%" PetscInt_FMT ")", ((PetscObject)pc)->type_name, n));
1851:     } else {
1852:       PetscCall(PetscSNPrintf(str, 25, "PC: %s", ((PetscObject)pc)->type_name));
1853:     }
1854:     PetscCall(PetscDrawStringBoxed(draw, x, y, PETSC_DRAW_RED, PETSC_DRAW_BLACK, str, NULL, &h));
1855:     bottom = y - h;
1856:     PetscCall(PetscDrawPushCurrentPoint(draw, x, bottom));
1857:     PetscTryTypeMethod(pc, view, viewer);
1858:     PetscCall(PetscDrawPopCurrentPoint(draw));
1859: #if defined(PETSC_HAVE_SAWS)
1860:   } else if (issaws) {
1861:     PetscMPIInt rank;

1863:     PetscCall(PetscObjectName((PetscObject)pc));
1864:     PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
1865:     if (!((PetscObject)pc)->amsmem && rank == 0) PetscCall(PetscObjectViewSAWs((PetscObject)pc, viewer));
1866:     if (pc->mat) PetscCall(MatView(pc->mat, viewer));
1867:     if (pc->pmat && pc->pmat != pc->mat) PetscCall(MatView(pc->pmat, viewer));
1868: #endif
1869:   }
1870:   PetscFunctionReturn(PETSC_SUCCESS);
1871: }

1873: /*@C
1874:   PCRegister -  Adds a method (`PCType`) to the preconditioner package.

1876:   Not collective. No Fortran Support

1878:   Input Parameters:
1879: + sname    - name of a new user-defined solver
1880: - function - routine to create method context

1882:   Example Usage:
1883: .vb
1884:    PCRegister("my_solver", MySolverCreate);
1885: .ve

1887:   Then, your solver can be chosen with the procedural interface via
1888: $     PCSetType(pc, "my_solver")
1889:   or at runtime via the option
1890: $     -pc_type my_solver

1892:   Level: advanced

1894:   Note:
1895:   `PCRegister()` may be called multiple times to add several user-defined preconditioners.

1897: .seealso: [](ch_ksp), `PC`, `PCType`, `PCRegisterAll()`
1898: @*/
1899: PetscErrorCode PCRegister(const char sname[], PetscErrorCode (*function)(PC))
1900: {
1901:   PetscFunctionBegin;
1902:   PetscCall(PCInitializePackage());
1903:   PetscCall(PetscFunctionListAdd(&PCList, sname, function));
1904:   PetscFunctionReturn(PETSC_SUCCESS);
1905: }

1907: static PetscErrorCode MatMult_PC(Mat A, Vec X, Vec Y)
1908: {
1909:   PC pc;

1911:   PetscFunctionBegin;
1912:   PetscCall(MatShellGetContext(A, &pc));
1913:   PetscCall(PCApply(pc, X, Y));
1914:   PetscFunctionReturn(PETSC_SUCCESS);
1915: }

1917: /*@C
1918:   PCComputeOperator - Computes the explicit preconditioned operator.

1920:   Collective

1922:   Input Parameters:
1923: + pc      - the preconditioner object
1924: - mattype - the `MatType` to be used for the operator

1926:   Output Parameter:
1927: . mat - the explicit preconditioned operator

1929:   Level: advanced

1931:   Note:
1932:   This computation is done by applying the operators to columns of the identity matrix.
1933:   This routine is costly in general, and is recommended for use only with relatively small systems.
1934:   Currently, this routine uses a dense matrix format when `mattype` == `NULL`

1936: .seealso: [](ch_ksp), `PC`, `KSPComputeOperator()`, `MatType`
1937: @*/
1938: PetscErrorCode PCComputeOperator(PC pc, MatType mattype, Mat *mat)
1939: {
1940:   PetscInt N, M, m, n;
1941:   Mat      A, Apc;

1943:   PetscFunctionBegin;
1945:   PetscAssertPointer(mat, 3);
1946:   PetscCall(PCGetOperators(pc, &A, NULL));
1947:   PetscCall(MatGetLocalSize(A, &m, &n));
1948:   PetscCall(MatGetSize(A, &M, &N));
1949:   PetscCall(MatCreateShell(PetscObjectComm((PetscObject)pc), m, n, M, N, pc, &Apc));
1950:   PetscCall(MatShellSetOperation(Apc, MATOP_MULT, (void (*)(void))MatMult_PC));
1951:   PetscCall(MatComputeOperator(Apc, mattype, mat));
1952:   PetscCall(MatDestroy(&Apc));
1953:   PetscFunctionReturn(PETSC_SUCCESS);
1954: }

1956: /*@
1957:   PCSetCoordinates - sets the coordinates of all the nodes on the local process

1959:   Collective

1961:   Input Parameters:
1962: + pc     - the solver context
1963: . dim    - the dimension of the coordinates 1, 2, or 3
1964: . nloc   - the blocked size of the coordinates array
1965: - coords - the coordinates array

1967:   Level: intermediate

1969:   Note:
1970:   `coords` is an array of the dim coordinates for the nodes on
1971:   the local processor, of size `dim`*`nloc`.
1972:   If there are 108 equation on a processor
1973:   for a displacement finite element discretization of elasticity (so
1974:   that there are nloc = 36 = 108/3 nodes) then the array must have 108
1975:   double precision values (ie, 3 * 36).  These x y z coordinates
1976:   should be ordered for nodes 0 to N-1 like so: [ 0.x, 0.y, 0.z, 1.x,
1977:   ... , N-1.z ].

1979: .seealso: [](ch_ksp), `PC`, `MatSetNearNullSpace()`
1980: @*/
1981: PetscErrorCode PCSetCoordinates(PC pc, PetscInt dim, PetscInt nloc, PetscReal coords[])
1982: {
1983:   PetscFunctionBegin;
1986:   PetscTryMethod(pc, "PCSetCoordinates_C", (PC, PetscInt, PetscInt, PetscReal[]), (pc, dim, nloc, coords));
1987:   PetscFunctionReturn(PETSC_SUCCESS);
1988: }

1990: /*@
1991:   PCGetInterpolations - Gets interpolation matrices for all levels (except level 0)

1993:   Logically Collective

1995:   Input Parameter:
1996: . pc - the precondition context

1998:   Output Parameters:
1999: + num_levels     - the number of levels
2000: - interpolations - the interpolation matrices (size of `num_levels`-1)

2002:   Level: advanced

2004:   Developer Note:
2005:   Why is this here instead of in `PCMG` etc?

2007: .seealso: [](ch_ksp), `PC`, `PCMG`, `PCMGGetRestriction()`, `PCMGSetInterpolation()`, `PCMGGetInterpolation()`, `PCGetCoarseOperators()`
2008: @*/
2009: PetscErrorCode PCGetInterpolations(PC pc, PetscInt *num_levels, Mat *interpolations[])
2010: {
2011:   PetscFunctionBegin;
2013:   PetscAssertPointer(num_levels, 2);
2014:   PetscAssertPointer(interpolations, 3);
2015:   PetscUseMethod(pc, "PCGetInterpolations_C", (PC, PetscInt *, Mat *[]), (pc, num_levels, interpolations));
2016:   PetscFunctionReturn(PETSC_SUCCESS);
2017: }

2019: /*@
2020:   PCGetCoarseOperators - Gets coarse operator matrices for all levels (except the finest level)

2022:   Logically Collective

2024:   Input Parameter:
2025: . pc - the precondition context

2027:   Output Parameters:
2028: + num_levels      - the number of levels
2029: - coarseOperators - the coarse operator matrices (size of `num_levels`-1)

2031:   Level: advanced

2033:   Developer Note:
2034:   Why is this here instead of in `PCMG` etc?

2036: .seealso: [](ch_ksp), `PC`, `PCMG`, `PCMGGetRestriction()`, `PCMGSetInterpolation()`, `PCMGGetRScale()`, `PCMGGetInterpolation()`, `PCGetInterpolations()`
2037: @*/
2038: PetscErrorCode PCGetCoarseOperators(PC pc, PetscInt *num_levels, Mat *coarseOperators[])
2039: {
2040:   PetscFunctionBegin;
2042:   PetscAssertPointer(num_levels, 2);
2043:   PetscAssertPointer(coarseOperators, 3);
2044:   PetscUseMethod(pc, "PCGetCoarseOperators_C", (PC, PetscInt *, Mat *[]), (pc, num_levels, coarseOperators));
2045:   PetscFunctionReturn(PETSC_SUCCESS);
2046: }