Actual source code: composite.c


  2: /*
  3:       Defines a preconditioner that can consist of a collection of PCs
  4: */
  5: #include <petsc/private/pcimpl.h>
  6: #include <petscksp.h>

  8: typedef struct _PC_CompositeLink *PC_CompositeLink;
  9: struct _PC_CompositeLink {
 10:   PC               pc;
 11:   PC_CompositeLink next;
 12:   PC_CompositeLink previous;
 13: };

 15: typedef struct {
 16:   PC_CompositeLink head;
 17:   PCCompositeType  type;
 18:   Vec              work1;
 19:   Vec              work2;
 20:   PetscScalar      alpha;
 21: } PC_Composite;

 23: static PetscErrorCode PCApply_Composite_Multiplicative(PC pc, Vec x, Vec y)
 24: {
 25:   PC_Composite    *jac  = (PC_Composite *)pc->data;
 26:   PC_CompositeLink next = jac->head;
 27:   Mat              mat  = pc->pmat;

 29:   PetscFunctionBegin;

 31:   PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");

 33:   /* Set the reuse flag on children PCs */
 34:   while (next) {
 35:     PetscCall(PCSetReusePreconditioner(next->pc, pc->reusepreconditioner));
 36:     next = next->next;
 37:   }
 38:   next = jac->head;

 40:   if (next->next && !jac->work2) { /* allocate second work vector */
 41:     PetscCall(VecDuplicate(jac->work1, &jac->work2));
 42:   }
 43:   if (pc->useAmat) mat = pc->mat;
 44:   PetscCall(PCApply(next->pc, x, y)); /* y <- B x */
 45:   while (next->next) {
 46:     next = next->next;
 47:     PetscCall(MatMult(mat, y, jac->work1));               /* work1 <- A y */
 48:     PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x)); /* work2 <- x - work1 */
 49:     PetscCall(PCApply(next->pc, jac->work2, jac->work1)); /* work1 <- C work2 */
 50:     PetscCall(VecAXPY(y, 1.0, jac->work1));               /* y <- y + work1 = B x + C (x - A B x) = (B + C (1 - A B)) x */
 51:   }
 52:   if (jac->type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
 53:     while (next->previous) {
 54:       next = next->previous;
 55:       PetscCall(MatMult(mat, y, jac->work1));
 56:       PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x));
 57:       PetscCall(PCApply(next->pc, jac->work2, jac->work1));
 58:       PetscCall(VecAXPY(y, 1.0, jac->work1));
 59:     }
 60:   }
 61:   PetscFunctionReturn(PETSC_SUCCESS);
 62: }

 64: static PetscErrorCode PCApplyTranspose_Composite_Multiplicative(PC pc, Vec x, Vec y)
 65: {
 66:   PC_Composite    *jac  = (PC_Composite *)pc->data;
 67:   PC_CompositeLink next = jac->head;
 68:   Mat              mat  = pc->pmat;

 70:   PetscFunctionBegin;
 71:   PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
 72:   if (next->next && !jac->work2) { /* allocate second work vector */
 73:     PetscCall(VecDuplicate(jac->work1, &jac->work2));
 74:   }
 75:   if (pc->useAmat) mat = pc->mat;
 76:   /* locate last PC */
 77:   while (next->next) next = next->next;
 78:   PetscCall(PCApplyTranspose(next->pc, x, y));
 79:   while (next->previous) {
 80:     next = next->previous;
 81:     PetscCall(MatMultTranspose(mat, y, jac->work1));
 82:     PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x));
 83:     PetscCall(PCApplyTranspose(next->pc, jac->work2, jac->work1));
 84:     PetscCall(VecAXPY(y, 1.0, jac->work1));
 85:   }
 86:   if (jac->type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
 87:     next = jac->head;
 88:     while (next->next) {
 89:       next = next->next;
 90:       PetscCall(MatMultTranspose(mat, y, jac->work1));
 91:       PetscCall(VecWAXPY(jac->work2, -1.0, jac->work1, x));
 92:       PetscCall(PCApplyTranspose(next->pc, jac->work2, jac->work1));
 93:       PetscCall(VecAXPY(y, 1.0, jac->work1));
 94:     }
 95:   }
 96:   PetscFunctionReturn(PETSC_SUCCESS);
 97: }

 99: /*
100:     This is very special for a matrix of the form alpha I + R + S
101: where first preconditioner is built from alpha I + S and second from
102: alpha I + R
103: */
104: static PetscErrorCode PCApply_Composite_Special(PC pc, Vec x, Vec y)
105: {
106:   PC_Composite    *jac  = (PC_Composite *)pc->data;
107:   PC_CompositeLink next = jac->head;

109:   PetscFunctionBegin;
110:   PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
111:   PetscCheck(next->next && !next->next->next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "Special composite preconditioners requires exactly two PCs");

113:   /* Set the reuse flag on children PCs */
114:   PetscCall(PCSetReusePreconditioner(next->pc, pc->reusepreconditioner));
115:   PetscCall(PCSetReusePreconditioner(next->next->pc, pc->reusepreconditioner));

117:   PetscCall(PCApply(next->pc, x, jac->work1));
118:   PetscCall(PCApply(next->next->pc, jac->work1, y));
119:   PetscFunctionReturn(PETSC_SUCCESS);
120: }

122: static PetscErrorCode PCApply_Composite_Additive(PC pc, Vec x, Vec y)
123: {
124:   PC_Composite    *jac  = (PC_Composite *)pc->data;
125:   PC_CompositeLink next = jac->head;

127:   PetscFunctionBegin;
128:   PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");

130:   /* Set the reuse flag on children PCs */
131:   while (next) {
132:     PetscCall(PCSetReusePreconditioner(next->pc, pc->reusepreconditioner));
133:     next = next->next;
134:   }
135:   next = jac->head;

137:   PetscCall(PCApply(next->pc, x, y));
138:   while (next->next) {
139:     next = next->next;
140:     PetscCall(PCApply(next->pc, x, jac->work1));
141:     PetscCall(VecAXPY(y, 1.0, jac->work1));
142:   }
143:   PetscFunctionReturn(PETSC_SUCCESS);
144: }

146: static PetscErrorCode PCApplyTranspose_Composite_Additive(PC pc, Vec x, Vec y)
147: {
148:   PC_Composite    *jac  = (PC_Composite *)pc->data;
149:   PC_CompositeLink next = jac->head;

151:   PetscFunctionBegin;
152:   PetscCheck(next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONGSTATE, "No composite preconditioners supplied via PCCompositeAddPCType() or -pc_composite_pcs");
153:   PetscCall(PCApplyTranspose(next->pc, x, y));
154:   while (next->next) {
155:     next = next->next;
156:     PetscCall(PCApplyTranspose(next->pc, x, jac->work1));
157:     PetscCall(VecAXPY(y, 1.0, jac->work1));
158:   }
159:   PetscFunctionReturn(PETSC_SUCCESS);
160: }

162: static PetscErrorCode PCSetUp_Composite(PC pc)
163: {
164:   PC_Composite    *jac  = (PC_Composite *)pc->data;
165:   PC_CompositeLink next = jac->head;
166:   DM               dm;

168:   PetscFunctionBegin;
169:   if (!jac->work1) PetscCall(MatCreateVecs(pc->pmat, &jac->work1, NULL));
170:   PetscCall(PCGetDM(pc, &dm));
171:   while (next) {
172:     if (!next->pc->dm) PetscCall(PCSetDM(next->pc, dm));
173:     if (!next->pc->mat) PetscCall(PCSetOperators(next->pc, pc->mat, pc->pmat));
174:     next = next->next;
175:   }
176:   PetscFunctionReturn(PETSC_SUCCESS);
177: }

179: static PetscErrorCode PCReset_Composite(PC pc)
180: {
181:   PC_Composite    *jac  = (PC_Composite *)pc->data;
182:   PC_CompositeLink next = jac->head;

184:   PetscFunctionBegin;
185:   while (next) {
186:     PetscCall(PCReset(next->pc));
187:     next = next->next;
188:   }
189:   PetscCall(VecDestroy(&jac->work1));
190:   PetscCall(VecDestroy(&jac->work2));
191:   PetscFunctionReturn(PETSC_SUCCESS);
192: }

194: static PetscErrorCode PCDestroy_Composite(PC pc)
195: {
196:   PC_Composite    *jac  = (PC_Composite *)pc->data;
197:   PC_CompositeLink next = jac->head, next_tmp;

199:   PetscFunctionBegin;
200:   PetscCall(PCReset_Composite(pc));
201:   while (next) {
202:     PetscCall(PCDestroy(&next->pc));
203:     next_tmp = next;
204:     next     = next->next;
205:     PetscCall(PetscFree(next_tmp));
206:   }
207:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSetType_C", NULL));
208:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetType_C", NULL));
209:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPCType_C", NULL));
210:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPC_C", NULL));
211:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetNumberPC_C", NULL));
212:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetPC_C", NULL));
213:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSpecialSetAlpha_C", NULL));
214:   PetscCall(PetscFree(pc->data));
215:   PetscFunctionReturn(PETSC_SUCCESS);
216: }

218: static PetscErrorCode PCSetFromOptions_Composite(PC pc, PetscOptionItems *PetscOptionsObject)
219: {
220:   PC_Composite    *jac  = (PC_Composite *)pc->data;
221:   PetscInt         nmax = 8, i;
222:   PC_CompositeLink next;
223:   char            *pcs[8];
224:   PetscBool        flg;

226:   PetscFunctionBegin;
227:   PetscOptionsHeadBegin(PetscOptionsObject, "Composite preconditioner options");
228:   PetscCall(PetscOptionsEnum("-pc_composite_type", "Type of composition", "PCCompositeSetType", PCCompositeTypes, (PetscEnum)jac->type, (PetscEnum *)&jac->type, &flg));
229:   if (flg) PetscCall(PCCompositeSetType(pc, jac->type));
230:   PetscCall(PetscOptionsStringArray("-pc_composite_pcs", "List of composite solvers", "PCCompositeAddPCType", pcs, &nmax, &flg));
231:   if (flg) {
232:     for (i = 0; i < nmax; i++) {
233:       PetscCall(PCCompositeAddPCType(pc, pcs[i]));
234:       PetscCall(PetscFree(pcs[i])); /* deallocate string pcs[i], which is allocated in PetscOptionsStringArray() */
235:     }
236:   }
237:   PetscOptionsHeadEnd();

239:   next = jac->head;
240:   while (next) {
241:     PetscCall(PCSetFromOptions(next->pc));
242:     next = next->next;
243:   }
244:   PetscFunctionReturn(PETSC_SUCCESS);
245: }

247: static PetscErrorCode PCView_Composite(PC pc, PetscViewer viewer)
248: {
249:   PC_Composite    *jac  = (PC_Composite *)pc->data;
250:   PC_CompositeLink next = jac->head;
251:   PetscBool        iascii;

253:   PetscFunctionBegin;
254:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
255:   if (iascii) {
256:     PetscCall(PetscViewerASCIIPrintf(viewer, "Composite PC type - %s\n", PCCompositeTypes[jac->type]));
257:     PetscCall(PetscViewerASCIIPrintf(viewer, "PCs on composite preconditioner follow\n"));
258:     PetscCall(PetscViewerASCIIPrintf(viewer, "---------------------------------\n"));
259:   }
260:   if (iascii) PetscCall(PetscViewerASCIIPushTab(viewer));
261:   while (next) {
262:     PetscCall(PCView(next->pc, viewer));
263:     next = next->next;
264:   }
265:   if (iascii) {
266:     PetscCall(PetscViewerASCIIPopTab(viewer));
267:     PetscCall(PetscViewerASCIIPrintf(viewer, "---------------------------------\n"));
268:   }
269:   PetscFunctionReturn(PETSC_SUCCESS);
270: }

272: static PetscErrorCode PCCompositeSpecialSetAlpha_Composite(PC pc, PetscScalar alpha)
273: {
274:   PC_Composite *jac = (PC_Composite *)pc->data;

276:   PetscFunctionBegin;
277:   jac->alpha = alpha;
278:   PetscFunctionReturn(PETSC_SUCCESS);
279: }

281: static PetscErrorCode PCCompositeSetType_Composite(PC pc, PCCompositeType type)
282: {
283:   PC_Composite *jac = (PC_Composite *)pc->data;

285:   PetscFunctionBegin;
286:   if (type == PC_COMPOSITE_ADDITIVE) {
287:     pc->ops->apply          = PCApply_Composite_Additive;
288:     pc->ops->applytranspose = PCApplyTranspose_Composite_Additive;
289:   } else if (type == PC_COMPOSITE_MULTIPLICATIVE || type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
290:     pc->ops->apply          = PCApply_Composite_Multiplicative;
291:     pc->ops->applytranspose = PCApplyTranspose_Composite_Multiplicative;
292:   } else if (type == PC_COMPOSITE_SPECIAL) {
293:     pc->ops->apply          = PCApply_Composite_Special;
294:     pc->ops->applytranspose = NULL;
295:   } else SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONG, "Unknown composite preconditioner type");
296:   jac->type = type;
297:   PetscFunctionReturn(PETSC_SUCCESS);
298: }

300: static PetscErrorCode PCCompositeGetType_Composite(PC pc, PCCompositeType *type)
301: {
302:   PC_Composite *jac = (PC_Composite *)pc->data;

304:   PetscFunctionBegin;
305:   *type = jac->type;
306:   PetscFunctionReturn(PETSC_SUCCESS);
307: }

309: static PetscErrorCode PCCompositeAddPC_Composite(PC pc, PC subpc)
310: {
311:   PC_Composite    *jac;
312:   PC_CompositeLink next, ilink;
313:   PetscInt         cnt = 0;
314:   const char      *prefix;
315:   char             newprefix[20];

317:   PetscFunctionBegin;
318:   PetscCall(PetscNew(&ilink));
319:   ilink->next = NULL;
320:   ilink->pc   = subpc;

322:   jac  = (PC_Composite *)pc->data;
323:   next = jac->head;
324:   if (!next) {
325:     jac->head       = ilink;
326:     ilink->previous = NULL;
327:   } else {
328:     cnt++;
329:     while (next->next) {
330:       next = next->next;
331:       cnt++;
332:     }
333:     next->next      = ilink;
334:     ilink->previous = next;
335:   }
336:   PetscCall(PCGetOptionsPrefix(pc, &prefix));
337:   PetscCall(PCSetOptionsPrefix(subpc, prefix));
338:   PetscCall(PetscSNPrintf(newprefix, 20, "sub_%d_", (int)cnt));
339:   PetscCall(PCAppendOptionsPrefix(subpc, newprefix));
340:   PetscCall(PetscObjectReference((PetscObject)subpc));
341:   PetscFunctionReturn(PETSC_SUCCESS);
342: }

344: static PetscErrorCode PCCompositeAddPCType_Composite(PC pc, PCType type)
345: {
346:   PC subpc;

348:   PetscFunctionBegin;
349:   PetscCall(PCCreate(PetscObjectComm((PetscObject)pc), &subpc));
350:   PetscCall(PetscObjectIncrementTabLevel((PetscObject)subpc, (PetscObject)pc, 1));
351:   PetscCall(PCCompositeAddPC_Composite(pc, subpc));
352:   /* type is set after prefix, because some methods may modify prefix, e.g. pcksp */
353:   PetscCall(PCSetType(subpc, type));
354:   PetscCall(PCDestroy(&subpc));
355:   PetscFunctionReturn(PETSC_SUCCESS);
356: }

358: static PetscErrorCode PCCompositeGetNumberPC_Composite(PC pc, PetscInt *n)
359: {
360:   PC_Composite    *jac;
361:   PC_CompositeLink next;

363:   PetscFunctionBegin;
364:   jac  = (PC_Composite *)pc->data;
365:   next = jac->head;
366:   *n   = 0;
367:   while (next) {
368:     next = next->next;
369:     (*n)++;
370:   }
371:   PetscFunctionReturn(PETSC_SUCCESS);
372: }

374: static PetscErrorCode PCCompositeGetPC_Composite(PC pc, PetscInt n, PC *subpc)
375: {
376:   PC_Composite    *jac;
377:   PC_CompositeLink next;
378:   PetscInt         i;

380:   PetscFunctionBegin;
381:   jac  = (PC_Composite *)pc->data;
382:   next = jac->head;
383:   for (i = 0; i < n; i++) {
384:     PetscCheck(next->next, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_INCOMP, "Not enough PCs in composite preconditioner");
385:     next = next->next;
386:   }
387:   *subpc = next->pc;
388:   PetscFunctionReturn(PETSC_SUCCESS);
389: }

391: /*@
392:    PCCompositeSetType - Sets the type of composite preconditioner.

394:    Logically Collective

396:    Input Parameters:
397: +  pc - the preconditioner context
398: -  type - `PC_COMPOSITE_ADDITIVE` (default), `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`

400:    Options Database Key:
401: .  -pc_composite_type <type: one of multiplicative, additive, special> - Sets composite preconditioner type

403:    Level: advanced

405: .seealso: `PCCOMPOSITE`, `PC_COMPOSITE_ADDITIVE`, `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`, `PCCompositeType`,
406:           `PCCompositeGetType()`
407: @*/
408: PetscErrorCode PCCompositeSetType(PC pc, PCCompositeType type)
409: {
410:   PetscFunctionBegin;
413:   PetscTryMethod(pc, "PCCompositeSetType_C", (PC, PCCompositeType), (pc, type));
414:   PetscFunctionReturn(PETSC_SUCCESS);
415: }

417: /*@
418:    PCCompositeGetType - Gets the type of composite preconditioner.

420:    Logically Collective

422:    Input Parameter:
423: .  pc - the preconditioner context

425:    Output Parameter:
426: .  type - `PC_COMPOSITE_ADDITIVE` (default), `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`

428:    Level: advanced

430: .seealso: `PCCOMPOSITE`, `PC_COMPOSITE_ADDITIVE`, `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`, `PCCompositeType`,
431:           `PCCompositeSetType()`
432: @*/
433: PetscErrorCode PCCompositeGetType(PC pc, PCCompositeType *type)
434: {
435:   PetscFunctionBegin;
437:   PetscUseMethod(pc, "PCCompositeGetType_C", (PC, PCCompositeType *), (pc, type));
438:   PetscFunctionReturn(PETSC_SUCCESS);
439: }

441: /*@
442:    PCCompositeSpecialSetAlpha - Sets alpha for the special composite preconditioner, `PC_COMPOSITE_SPECIAL`,
443:      for alphaI + R + S

445:    Logically Collective

447:    Input Parameters:
448: +  pc - the preconditioner context
449: -  alpha - scale on identity

451:    Level: Developer

453: .seealso: `PCCOMPOSITE`, `PC_COMPOSITE_ADDITIVE`, `PC_COMPOSITE_MULTIPLICATIVE`, `PC_COMPOSITE_SPECIAL`, `PCCompositeType`,
454:           `PCCompositeSetType()`, `PCCompositeGetType()`
455: @*/
456: PetscErrorCode PCCompositeSpecialSetAlpha(PC pc, PetscScalar alpha)
457: {
458:   PetscFunctionBegin;
461:   PetscTryMethod(pc, "PCCompositeSpecialSetAlpha_C", (PC, PetscScalar), (pc, alpha));
462:   PetscFunctionReturn(PETSC_SUCCESS);
463: }

465: /*@C
466:   PCCompositeAddPCType - Adds another `PC` of the given type to the composite `PC`.

468:   Collective

470:   Input Parameters:
471: + pc - the preconditioner context
472: - type - the type of the new preconditioner

474:   Level: intermediate

476: .seealso: `PCCOMPOSITE`, `PCCompositeAddPC()`, `PCCompositeGetNumberPC()`
477: @*/
478: PetscErrorCode PCCompositeAddPCType(PC pc, PCType type)
479: {
480:   PetscFunctionBegin;
482:   PetscTryMethod(pc, "PCCompositeAddPCType_C", (PC, PCType), (pc, type));
483:   PetscFunctionReturn(PETSC_SUCCESS);
484: }

486: /*@
487:   PCCompositeAddPC - Adds another `PC` to the composite `PC`.

489:   Collective

491:   Input Parameters:
492: + pc    - the preconditioner context
493: - subpc - the new preconditioner

495:    Level: intermediate

497: .seealso: `PCCOMPOSITE`, `PCCompositeAddPCType()`, `PCCompositeGetNumberPC()`
498: @*/
499: PetscErrorCode PCCompositeAddPC(PC pc, PC subpc)
500: {
501:   PetscFunctionBegin;
504:   PetscTryMethod(pc, "PCCompositeAddPC_C", (PC, PC), (pc, subpc));
505:   PetscFunctionReturn(PETSC_SUCCESS);
506: }

508: /*@
509:    PCCompositeGetNumberPC - Gets the number of `PC` objects in the composite `PC`.

511:    Not Collective

513:    Input Parameter:
514: .  pc - the preconditioner context

516:    Output Parameter:
517: .  num - the number of sub pcs

519:    Level: Developer

521: .seealso: `PCCOMPOSITE`, `PCCompositeGetPC()`, `PCCompositeAddPC()`, `PCCompositeAddPCType()`
522: @*/
523: PetscErrorCode PCCompositeGetNumberPC(PC pc, PetscInt *num)
524: {
525:   PetscFunctionBegin;
528:   PetscUseMethod(pc, "PCCompositeGetNumberPC_C", (PC, PetscInt *), (pc, num));
529:   PetscFunctionReturn(PETSC_SUCCESS);
530: }

532: /*@
533:    PCCompositeGetPC - Gets one of the `PC` objects in the composite `PC`.

535:    Not Collective

537:    Input Parameters:
538: +  pc - the preconditioner context
539: -  n - the number of the pc requested

541:    Output Parameter:
542: .  subpc - the PC requested

544:    Level: intermediate

546:     Note:
547:     To use a different operator to construct one of the inner preconditioners first call `PCCompositeGetPC()`, then
548:     call `PCSetOperators()` on that `PC`.

550: .seealso: `PCCOMPOSITE`, `PCCompositeAddPCType()`, `PCCompositeGetNumberPC()`, `PCSetOperators()`
551: @*/
552: PetscErrorCode PCCompositeGetPC(PC pc, PetscInt n, PC *subpc)
553: {
554:   PetscFunctionBegin;
557:   PetscUseMethod(pc, "PCCompositeGetPC_C", (PC, PetscInt, PC *), (pc, n, subpc));
558:   PetscFunctionReturn(PETSC_SUCCESS);
559: }

561: /*MC
562:      PCCOMPOSITE - Build a preconditioner by composing together several preconditioners

564:    Options Database Keys:
565: +  -pc_composite_type <type: one of multiplicative, additive, symmetric_multiplicative, special> - Sets composite preconditioner type
566: .  -pc_use_amat - activates `PCSetUseAmat()`
567: -  -pc_composite_pcs - <pc0,pc1,...> list of PCs to compose

569:    Level: intermediate

571:    Notes:
572:    To use a Krylov method inside the composite preconditioner, set the `PCType` of one or more
573:    inner `PC`s to be `PCKSP`. Using a Krylov method inside another Krylov method can be dangerous (you get divergence or
574:    the incorrect answer) unless you use `KSPFGMRES` as the outer Krylov method

576:    To use a different operator to construct one of the inner preconditioners first call `PCCompositeGetPC()`, then
577:    call `PCSetOperators()` on that `PC`.

579: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
580:           `PCSHELL`, `PCKSP`, `PCCompositeSetType()`, `PCCompositeSpecialSetAlpha()`, `PCCompositeAddPCType()`,
581:           `PCCompositeGetPC()`, `PCSetUseAmat()`, `PCCompositeAddPC()`, `PCCompositeGetNumberPC()`
582: M*/

584: PETSC_EXTERN PetscErrorCode PCCreate_Composite(PC pc)
585: {
586:   PC_Composite *jac;

588:   PetscFunctionBegin;
589:   PetscCall(PetscNew(&jac));

591:   pc->ops->apply           = PCApply_Composite_Additive;
592:   pc->ops->applytranspose  = PCApplyTranspose_Composite_Additive;
593:   pc->ops->setup           = PCSetUp_Composite;
594:   pc->ops->reset           = PCReset_Composite;
595:   pc->ops->destroy         = PCDestroy_Composite;
596:   pc->ops->setfromoptions  = PCSetFromOptions_Composite;
597:   pc->ops->view            = PCView_Composite;
598:   pc->ops->applyrichardson = NULL;

600:   pc->data   = (void *)jac;
601:   jac->type  = PC_COMPOSITE_ADDITIVE;
602:   jac->work1 = NULL;
603:   jac->work2 = NULL;
604:   jac->head  = NULL;

606:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSetType_C", PCCompositeSetType_Composite));
607:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetType_C", PCCompositeGetType_Composite));
608:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPCType_C", PCCompositeAddPCType_Composite));
609:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeAddPC_C", PCCompositeAddPC_Composite));
610:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetNumberPC_C", PCCompositeGetNumberPC_Composite));
611:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeGetPC_C", PCCompositeGetPC_Composite));
612:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCCompositeSpecialSetAlpha_C", PCCompositeSpecialSetAlpha_Composite));
613:   PetscFunctionReturn(PETSC_SUCCESS);
614: }