Actual source code: composite.c
1: /*
2: Defines a preconditioner that can consist of a collection of PCs
3: */
4: #include <petsc/private/pcimpl.h>
5: #include <petscksp.h>
7: typedef struct _PC_CompositeLink *PC_CompositeLink;
8: struct _PC_CompositeLink {
9: PC pc;
10: PC_CompositeLink next;
11: PC_CompositeLink previous;
12: };
14: typedef struct {
15: PC_CompositeLink head;
16: PCCompositeType type;
17: Vec work1;
18: Vec work2;
19: PetscScalar alpha;
20: } PC_Composite;
22: static PetscErrorCode PCApply_Composite_Multiplicative(PC pc, Vec x, Vec y)
23: {
24: PC_Composite *jac = (PC_Composite *)pc->data;
25: PC_CompositeLink next = jac->head;
26: Mat mat = pc->pmat;
28: PetscFunctionBegin;
29: PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
31: /* Set the reuse flag on children PCs */
32: while (next) {
33: PetscCall(PCSetReusePreconditioner(next->pc, pc->reusepreconditioner));
34: next = next->next;
35: }
36: next = jac->head;
38: if (next->next && !jac->work2) { /* allocate second work vector */
39: PetscCall(VecDuplicate(jac->work1, &jac->work2));
40: }
41: if (pc->useAmat) mat = pc->mat;
42: PetscCall(PCApply(next->pc, x, y)); /* y <- B x */
43: while (next->next) {
44: next = next->next;
45: PetscCall(MatMult(mat, y, jac->work1)); /* work1 <- A y */
46: PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x)); /* work2 <- x - work1 */
47: PetscCall(PCApply(next->pc, jac->work2, jac->work1)); /* work1 <- C work2 */
48: PetscCall(VecAXPY(y, 1.0, jac->work1)); /* y <- y + work1 = B x + C (x - A B x) = (B + C (1 - A B)) x */
49: }
50: if (jac->type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
51: while (next->previous) {
52: next = next->previous;
53: PetscCall(MatMult(mat, y, jac->work1));
54: PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x));
55: PetscCall(PCApply(next->pc, jac->work2, jac->work1));
56: PetscCall(VecAXPY(y, 1.0, jac->work1));
57: }
58: }
59: PetscFunctionReturn(PETSC_SUCCESS);
60: }
62: static PetscErrorCode PCApplyTranspose_Composite_Multiplicative(PC pc, Vec x, Vec y)
63: {
64: PC_Composite *jac = (PC_Composite *)pc->data;
65: PC_CompositeLink next = jac->head;
66: Mat mat = pc->pmat;
68: PetscFunctionBegin;
69: PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
70: if (next->next && !jac->work2) { /* allocate second work vector */
71: PetscCall(VecDuplicate(jac->work1, &jac->work2));
72: }
73: if (pc->useAmat) mat = pc->mat;
74: /* locate last PC */
75: while (next->next) next = next->next;
76: PetscCall(PCApplyTranspose(next->pc, x, y));
77: while (next->previous) {
78: next = next->previous;
79: PetscCall(MatMultTranspose(mat, y, jac->work1));
80: PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x));
81: PetscCall(PCApplyTranspose(next->pc, jac->work2, jac->work1));
82: PetscCall(VecAXPY(y, 1.0, jac->work1));
83: }
84: if (jac->type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
85: next = jac->head;
86: while (next->next) {
87: next = next->next;
88: PetscCall(MatMultTranspose(mat, y, jac->work1));
89: PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x));
90: PetscCall(PCApplyTranspose(next->pc, jac->work2, jac->work1));
91: PetscCall(VecAXPY(y, 1.0, jac->work1));
92: }
93: }
94: PetscFunctionReturn(PETSC_SUCCESS);
95: }
97: /*
98: This is very special for a matrix of the form alpha I + R + S
99: where first preconditioner is built from alpha I + S and second from
100: alpha I + R
101: */
102: static PetscErrorCode PCApply_Composite_Special(PC pc, Vec x, Vec y)
103: {
104: PC_Composite *jac = (PC_Composite *)pc->data;
105: PC_CompositeLink next = jac->head;
107: PetscFunctionBegin;
108: PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
109: PetscCheck(next->next && !next->next->next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "Special composite preconditioners requires exactly two PCs");
111: /* Set the reuse flag on children PCs */
112: PetscCall(PCSetReusePreconditioner(next->pc, pc->reusepreconditioner));
113: PetscCall(PCSetReusePreconditioner(next->next->pc, pc->reusepreconditioner));
115: PetscCall(PCApply(next->pc, x, jac->work1));
116: PetscCall(PCApply(next->next->pc, jac->work1, y));
117: PetscFunctionReturn(PETSC_SUCCESS);
118: }
120: static PetscErrorCode PCApply_Composite_Additive(PC pc, Vec x, Vec y)
121: {
122: PC_Composite *jac = (PC_Composite *)pc->data;
123: PC_CompositeLink next = jac->head;
125: PetscFunctionBegin;
126: PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
128: /* Set the reuse flag on children PCs */
129: while (next) {
130: PetscCall(PCSetReusePreconditioner(next->pc, pc->reusepreconditioner));
131: next = next->next;
132: }
133: next = jac->head;
135: PetscCall(PCApply(next->pc, x, y));
136: while (next->next) {
137: next = next->next;
138: PetscCall(PCApply(next->pc, x, jac->work1));
139: PetscCall(VecAXPY(y, 1.0, jac->work1));
140: }
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
144: static PetscErrorCode PCApplyTranspose_Composite_Additive(PC pc, Vec x, Vec y)
145: {
146: PC_Composite *jac = (PC_Composite *)pc->data;
147: PC_CompositeLink next = jac->head;
149: PetscFunctionBegin;
150: PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
151: PetscCall(PCApplyTranspose(next->pc, x, y));
152: while (next->next) {
153: next = next->next;
154: PetscCall(PCApplyTranspose(next->pc, x, jac->work1));
155: PetscCall(VecAXPY(y, 1.0, jac->work1));
156: }
157: PetscFunctionReturn(PETSC_SUCCESS);
158: }
160: static PetscErrorCode PCSetUp_Composite(PC pc)
161: {
162: PC_Composite *jac = (PC_Composite *)pc->data;
163: PC_CompositeLink next = jac->head;
164: DM dm;
166: PetscFunctionBegin;
167: if (!jac->work1) PetscCall(MatCreateVecs(pc->pmat, &jac->work1, NULL));
168: PetscCall(PCGetDM(pc, &dm));
169: while (next) {
170: if (!next->pc->dm) PetscCall(PCSetDM(next->pc, dm));
171: if (!next->pc->mat) PetscCall(PCSetOperators(next->pc, pc->mat, pc->pmat));
172: next = next->next;
173: }
174: PetscFunctionReturn(PETSC_SUCCESS);
175: }
177: static PetscErrorCode PCSetUpOnBlocks_Composite(PC pc)
178: {
179: PC_Composite *jac = (PC_Composite *)pc->data;
180: PC_CompositeLink next = jac->head;
181: PCFailedReason reason;
183: PetscFunctionBegin;
184: while (next) {
185: PetscCall(PCSetUp(next->pc));
186: PetscCall(PCGetFailedReasonRank(next->pc, &reason));
187: if (reason) pc->failedreason = reason;
188: next = next->next;
189: }
190: PetscFunctionReturn(PETSC_SUCCESS);
191: }
193: static PetscErrorCode PCReset_Composite(PC pc)
194: {
195: PC_Composite *jac = (PC_Composite *)pc->data;
196: PC_CompositeLink next = jac->head;
198: PetscFunctionBegin;
199: while (next) {
200: PetscCall(PCReset(next->pc));
201: next = next->next;
202: }
203: PetscCall(VecDestroy(&jac->work1));
204: PetscCall(VecDestroy(&jac->work2));
205: PetscFunctionReturn(PETSC_SUCCESS);
206: }
208: static PetscErrorCode PCDestroy_Composite(PC pc)
209: {
210: PC_Composite *jac = (PC_Composite *)pc->data;
211: PC_CompositeLink next = jac->head, next_tmp;
213: PetscFunctionBegin;
214: PetscCall(PCReset_Composite(pc));
215: while (next) {
216: PetscCall(PCDestroy(&next->pc));
217: next_tmp = next;
218: next = next->next;
219: PetscCall(PetscFree(next_tmp));
220: }
221: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSetType_C", NULL));
222: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetType_C", NULL));
223: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPCType_C", NULL));
224: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPC_C", NULL));
225: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetNumberPC_C", NULL));
226: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetPC_C", NULL));
227: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSpecialSetAlpha_C", NULL));
228: PetscCall(PetscFree(pc->data));
229: PetscFunctionReturn(PETSC_SUCCESS);
230: }
232: static PetscErrorCode PCSetFromOptions_Composite(PC pc, PetscOptionItems *PetscOptionsObject)
233: {
234: PC_Composite *jac = (PC_Composite *)pc->data;
235: PetscInt nmax = 8, i;
236: PC_CompositeLink next;
237: char *pcs[8];
238: PetscBool flg;
240: PetscFunctionBegin;
241: PetscOptionsHeadBegin(PetscOptionsObject, "Composite preconditioner options");
242: PetscCall(PetscOptionsEnum("-pc_composite_type", "Type of composition", "PCCompositeSetType", PCCompositeTypes, (PetscEnum)jac->type, (PetscEnum *)&jac->type, &flg));
243: if (flg) PetscCall(PCCompositeSetType(pc, jac->type));
244: PetscCall(PetscOptionsStringArray("-pc_composite_pcs", "List of composite solvers", "PCCompositeAddPCType", pcs, &nmax, &flg));
245: if (flg) {
246: for (i = 0; i < nmax; i++) {
247: PetscCall(PCCompositeAddPCType(pc, pcs[i]));
248: PetscCall(PetscFree(pcs[i])); /* deallocate string pcs[i], which is allocated in PetscOptionsStringArray() */
249: }
250: }
251: PetscOptionsHeadEnd();
253: next = jac->head;
254: while (next) {
255: PetscCall(PCSetFromOptions(next->pc));
256: next = next->next;
257: }
258: PetscFunctionReturn(PETSC_SUCCESS);
259: }
261: static PetscErrorCode PCView_Composite(PC pc, PetscViewer viewer)
262: {
263: PC_Composite *jac = (PC_Composite *)pc->data;
264: PC_CompositeLink next = jac->head;
265: PetscBool iascii;
267: PetscFunctionBegin;
268: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
269: if (iascii) {
270: PetscCall(PetscViewerASCIIPrintf(viewer, "Composite PC type - %s\n", PCCompositeTypes[jac->type]));
271: PetscCall(PetscViewerASCIIPrintf(viewer, "PCs on composite preconditioner follow\n"));
272: PetscCall(PetscViewerASCIIPrintf(viewer, "---------------------------------\n"));
273: }
274: if (iascii) PetscCall(PetscViewerASCIIPushTab(viewer));
275: while (next) {
276: PetscCall(PCView(next->pc, viewer));
277: next = next->next;
278: }
279: if (iascii) {
280: PetscCall(PetscViewerASCIIPopTab(viewer));
281: PetscCall(PetscViewerASCIIPrintf(viewer, "---------------------------------\n"));
282: }
283: PetscFunctionReturn(PETSC_SUCCESS);
284: }
286: static PetscErrorCode PCCompositeSpecialSetAlpha_Composite(PC pc, PetscScalar alpha)
287: {
288: PC_Composite *jac = (PC_Composite *)pc->data;
290: PetscFunctionBegin;
291: jac->alpha = alpha;
292: PetscFunctionReturn(PETSC_SUCCESS);
293: }
295: static PetscErrorCode PCCompositeSetType_Composite(PC pc, PCCompositeType type)
296: {
297: PC_Composite *jac = (PC_Composite *)pc->data;
299: PetscFunctionBegin;
300: if (type == PC_COMPOSITE_ADDITIVE) {
301: pc->ops->apply = PCApply_Composite_Additive;
302: pc->ops->applytranspose = PCApplyTranspose_Composite_Additive;
303: } else if (type == PC_COMPOSITE_MULTIPLICATIVE || type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
304: pc->ops->apply = PCApply_Composite_Multiplicative;
305: pc->ops->applytranspose = PCApplyTranspose_Composite_Multiplicative;
306: } else if (type == PC_COMPOSITE_SPECIAL) {
307: pc->ops->apply = PCApply_Composite_Special;
308: pc->ops->applytranspose = NULL;
309: } else SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONG, "Unknown composite preconditioner type");
310: jac->type = type;
311: PetscFunctionReturn(PETSC_SUCCESS);
312: }
314: static PetscErrorCode PCCompositeGetType_Composite(PC pc, PCCompositeType *type)
315: {
316: PC_Composite *jac = (PC_Composite *)pc->data;
318: PetscFunctionBegin;
319: *type = jac->type;
320: PetscFunctionReturn(PETSC_SUCCESS);
321: }
323: static PetscErrorCode PCCompositeAddPC_Composite(PC pc, PC subpc)
324: {
325: PC_Composite *jac;
326: PC_CompositeLink next, ilink;
327: PetscInt cnt = 0;
328: const char *prefix;
329: char newprefix[20];
331: PetscFunctionBegin;
332: PetscCall(PetscNew(&ilink));
333: ilink->next = NULL;
334: ilink->pc = subpc;
336: jac = (PC_Composite *)pc->data;
337: next = jac->head;
338: if (!next) {
339: jac->head = ilink;
340: ilink->previous = NULL;
341: } else {
342: cnt++;
343: while (next->next) {
344: next = next->next;
345: cnt++;
346: }
347: next->next = ilink;
348: ilink->previous = next;
349: }
350: PetscCall(PCGetOptionsPrefix(pc, &prefix));
351: PetscCall(PCSetOptionsPrefix(subpc, prefix));
352: PetscCall(PetscSNPrintf(newprefix, 20, "sub_%d_", (int)cnt));
353: PetscCall(PCAppendOptionsPrefix(subpc, newprefix));
354: PetscCall(PetscObjectReference((PetscObject)subpc));
355: PetscFunctionReturn(PETSC_SUCCESS);
356: }
358: static PetscErrorCode PCCompositeAddPCType_Composite(PC pc, PCType type)
359: {
360: PC subpc;
362: PetscFunctionBegin;
363: PetscCall(PCCreate(PetscObjectComm((PetscObject)pc), &subpc));
364: PetscCall(PetscObjectIncrementTabLevel((PetscObject)subpc, (PetscObject)pc, 1));
365: PetscCall(PCCompositeAddPC_Composite(pc, subpc));
366: /* type is set after prefix, because some methods may modify prefix, e.g. pcksp */
367: PetscCall(PCSetType(subpc, type));
368: PetscCall(PCDestroy(&subpc));
369: PetscFunctionReturn(PETSC_SUCCESS);
370: }
372: static PetscErrorCode PCCompositeGetNumberPC_Composite(PC pc, PetscInt *n)
373: {
374: PC_Composite *jac;
375: PC_CompositeLink next;
377: PetscFunctionBegin;
378: jac = (PC_Composite *)pc->data;
379: next = jac->head;
380: *n = 0;
381: while (next) {
382: next = next->next;
383: (*n)++;
384: }
385: PetscFunctionReturn(PETSC_SUCCESS);
386: }
388: static PetscErrorCode PCCompositeGetPC_Composite(PC pc, PetscInt n, PC *subpc)
389: {
390: PC_Composite *jac;
391: PC_CompositeLink next;
392: PetscInt i;
394: PetscFunctionBegin;
395: jac = (PC_Composite *)pc->data;
396: next = jac->head;
397: for (i = 0; i < n; i++) {
398: PetscCheck(next->next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Not enough PCs in composite preconditioner");
399: next = next->next;
400: }
401: *subpc = next->pc;
402: PetscFunctionReturn(PETSC_SUCCESS);
403: }
405: /*@
406: PCCompositeSetType - Sets the type of composite preconditioner.
408: Logically Collective
410: Input Parameters:
411: + pc - the preconditioner context
412: - type - `PC_COMPOSITE_ADDITIVE` (default), `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`
414: Options Database Key:
415: . -pc_composite_type <type: one of multiplicative, additive, special> - Sets composite preconditioner type
417: Level: advanced
419: .seealso: [](ch_ksp), `PCCOMPOSITE`, `PC_COMPOSITE_ADDITIVE`, `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`, `PCCompositeType`,
420: `PCCompositeGetType()`
421: @*/
422: PetscErrorCode PCCompositeSetType(PC pc, PCCompositeType type)
423: {
424: PetscFunctionBegin;
427: PetscTryMethod(pc, "PCCompositeSetType_C", (PC, PCCompositeType), (pc, type));
428: PetscFunctionReturn(PETSC_SUCCESS);
429: }
431: /*@
432: PCCompositeGetType - Gets the type of composite preconditioner.
434: Logically Collective
436: Input Parameter:
437: . pc - the preconditioner context
439: Output Parameter:
440: . type - `PC_COMPOSITE_ADDITIVE` (default), `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`
442: Level: advanced
444: .seealso: [](ch_ksp), `PCCOMPOSITE`, `PC_COMPOSITE_ADDITIVE`, `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`, `PCCompositeType`,
445: `PCCompositeSetType()`
446: @*/
447: PetscErrorCode PCCompositeGetType(PC pc, PCCompositeType *type)
448: {
449: PetscFunctionBegin;
451: PetscUseMethod(pc, "PCCompositeGetType_C", (PC, PCCompositeType *), (pc, type));
452: PetscFunctionReturn(PETSC_SUCCESS);
453: }
455: /*@
456: PCCompositeSpecialSetAlpha - Sets alpha for the special composite preconditioner, `PC_COMPOSITE_SPECIAL`,
457: for $\alpha I + R + S$
459: Logically Collective
461: Input Parameters:
462: + pc - the preconditioner context
463: - alpha - scale on identity
465: Level: developer
467: .seealso: [](ch_ksp), `PCCOMPOSITE`, `PC_COMPOSITE_ADDITIVE`, `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`, `PCCompositeType`,
468: `PCCompositeSetType()`, `PCCompositeGetType()`
469: @*/
470: PetscErrorCode PCCompositeSpecialSetAlpha(PC pc, PetscScalar alpha)
471: {
472: PetscFunctionBegin;
475: PetscTryMethod(pc, "PCCompositeSpecialSetAlpha_C", (PC, PetscScalar), (pc, alpha));
476: PetscFunctionReturn(PETSC_SUCCESS);
477: }
479: /*@C
480: PCCompositeAddPCType - Adds another `PC` of the given type to the composite `PC`.
482: Collective
484: Input Parameters:
485: + pc - the preconditioner context
486: - type - the type of the new preconditioner
488: Level: intermediate
490: .seealso: [](ch_ksp), `PCCOMPOSITE`, `PCCompositeAddPC()`, `PCCompositeGetNumberPC()`
491: @*/
492: PetscErrorCode PCCompositeAddPCType(PC pc, PCType type)
493: {
494: PetscFunctionBegin;
496: PetscTryMethod(pc, "PCCompositeAddPCType_C", (PC, PCType), (pc, type));
497: PetscFunctionReturn(PETSC_SUCCESS);
498: }
500: /*@
501: PCCompositeAddPC - Adds another `PC` to the composite `PC`.
503: Collective
505: Input Parameters:
506: + pc - the preconditioner context
507: - subpc - the new preconditioner
509: Level: intermediate
511: .seealso: [](ch_ksp), `PCCOMPOSITE`, `PCCompositeAddPCType()`, `PCCompositeGetNumberPC()`
512: @*/
513: PetscErrorCode PCCompositeAddPC(PC pc, PC subpc)
514: {
515: PetscFunctionBegin;
518: PetscTryMethod(pc, "PCCompositeAddPC_C", (PC, PC), (pc, subpc));
519: PetscFunctionReturn(PETSC_SUCCESS);
520: }
522: /*@
523: PCCompositeGetNumberPC - Gets the number of `PC` objects in the composite `PC`.
525: Not Collective
527: Input Parameter:
528: . pc - the preconditioner context
530: Output Parameter:
531: . num - the number of sub pcs
533: Level: developer
535: .seealso: [](ch_ksp), `PCCOMPOSITE`, `PCCompositeGetPC()`, `PCCompositeAddPC()`, `PCCompositeAddPCType()`
536: @*/
537: PetscErrorCode PCCompositeGetNumberPC(PC pc, PetscInt *num)
538: {
539: PetscFunctionBegin;
541: PetscAssertPointer(num, 2);
542: PetscUseMethod(pc, "PCCompositeGetNumberPC_C", (PC, PetscInt *), (pc, num));
543: PetscFunctionReturn(PETSC_SUCCESS);
544: }
546: /*@
547: PCCompositeGetPC - Gets one of the `PC` objects in the composite `PC`.
549: Not Collective
551: Input Parameters:
552: + pc - the preconditioner context
553: - n - the number of the pc requested
555: Output Parameter:
556: . subpc - the PC requested
558: Level: intermediate
560: Note:
561: To use a different operator to construct one of the inner preconditioners first call `PCCompositeGetPC()`, then
562: call `PCSetOperators()` on that `PC`.
564: .seealso: [](ch_ksp), `PCCOMPOSITE`, `PCCompositeAddPCType()`, `PCCompositeGetNumberPC()`, `PCSetOperators()`
565: @*/
566: PetscErrorCode PCCompositeGetPC(PC pc, PetscInt n, PC *subpc)
567: {
568: PetscFunctionBegin;
570: PetscAssertPointer(subpc, 3);
571: PetscUseMethod(pc, "PCCompositeGetPC_C", (PC, PetscInt, PC *), (pc, n, subpc));
572: PetscFunctionReturn(PETSC_SUCCESS);
573: }
575: /*MC
576: PCCOMPOSITE - Build a preconditioner by composing together several preconditioners
578: Options Database Keys:
579: + -pc_composite_type <type: one of multiplicative, additive, symmetric_multiplicative, special> - Sets composite preconditioner type
580: . -pc_use_amat - activates `PCSetUseAmat()`
581: - -pc_composite_pcs - <pc0,pc1,...> list of PCs to compose
583: Level: intermediate
585: Notes:
586: To use a Krylov method inside the composite preconditioner, set the `PCType` of one or more
587: inner `PC`s to be `PCKSP`. Using a Krylov method inside another Krylov method can be dangerous (you get divergence or
588: the incorrect answer) unless you use `KSPFGMRES` as the outer Krylov method
590: To use a different operator to construct one of the inner preconditioners first call `PCCompositeGetPC()`, then
591: call `PCSetOperators()` on that `PC`.
593: .seealso: [](ch_ksp), `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
594: `PCSHELL`, `PCKSP`, `PCCompositeSetType()`, `PCCompositeSpecialSetAlpha()`, `PCCompositeAddPCType()`,
595: `PCCompositeGetPC()`, `PCSetUseAmat()`, `PCCompositeAddPC()`, `PCCompositeGetNumberPC()`
596: M*/
598: PETSC_EXTERN PetscErrorCode PCCreate_Composite(PC pc)
599: {
600: PC_Composite *jac;
602: PetscFunctionBegin;
603: PetscCall(PetscNew(&jac));
605: pc->ops->apply = PCApply_Composite_Additive;
606: pc->ops->applytranspose = PCApplyTranspose_Composite_Additive;
607: pc->ops->setup = PCSetUp_Composite;
608: pc->ops->setuponblocks = PCSetUpOnBlocks_Composite;
609: pc->ops->reset = PCReset_Composite;
610: pc->ops->destroy = PCDestroy_Composite;
611: pc->ops->setfromoptions = PCSetFromOptions_Composite;
612: pc->ops->view = PCView_Composite;
613: pc->ops->applyrichardson = NULL;
615: pc->data = (void *)jac;
616: jac->type = PC_COMPOSITE_ADDITIVE;
617: jac->work1 = NULL;
618: jac->work2 = NULL;
619: jac->head = NULL;
621: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSetType_C", PCCompositeSetType_Composite));
622: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetType_C", PCCompositeGetType_Composite));
623: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPCType_C", PCCompositeAddPCType_Composite));
624: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPC_C", PCCompositeAddPC_Composite));
625: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetNumberPC_C", PCCompositeGetNumberPC_Composite));
626: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetPC_C", PCCompositeGetPC_Composite));
627: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSpecialSetAlpha_C", PCCompositeSpecialSetAlpha_Composite));
628: PetscFunctionReturn(PETSC_SUCCESS);
629: }