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: /*@C
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()`, `MatDiagonalScale()`
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()`, `MatDiagonalScale()`
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`, `PCMG`, `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`, `PCMG`, `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: /*@C
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: /*@C
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: /*@C
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: /*@C
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: /*@C
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: /*@C
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
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: }