Actual source code: redundant.c


  2: /*
  3:   This file defines a "solve the problem redundantly on each subgroup of processor" preconditioner.
  4: */
  5: #include <petsc/private/pcimpl.h>
  6: #include <petscksp.h>

  8: typedef struct {
  9:   KSP                ksp;
 10:   PC                 pc;                    /* actual preconditioner used on each processor */
 11:   Vec                xsub, ysub;            /* vectors of a subcommunicator to hold parallel vectors of PetscObjectComm((PetscObject)pc) */
 12:   Vec                xdup, ydup;            /* parallel vector that congregates xsub or ysub facilitating vector scattering */
 13:   Mat                pmats;                 /* matrix and optional preconditioner matrix belong to a subcommunicator */
 14:   VecScatter         scatterin, scatterout; /* scatter used to move all values to each processor group (subcommunicator) */
 15:   PetscBool          useparallelmat;
 16:   PetscSubcomm       psubcomm;
 17:   PetscInt           nsubcomm; /* num of data structure PetscSubcomm */
 18:   PetscBool          shifttypeset;
 19:   MatFactorShiftType shifttype;
 20: } PC_Redundant;

 22: PetscErrorCode PCFactorSetShiftType_Redundant(PC pc, MatFactorShiftType shifttype)
 23: {
 24:   PC_Redundant *red = (PC_Redundant *)pc->data;

 26:   PetscFunctionBegin;
 27:   if (red->ksp) {
 28:     PC pc;
 29:     PetscCall(KSPGetPC(red->ksp, &pc));
 30:     PetscCall(PCFactorSetShiftType(pc, shifttype));
 31:   } else {
 32:     red->shifttypeset = PETSC_TRUE;
 33:     red->shifttype    = shifttype;
 34:   }
 35:   PetscFunctionReturn(PETSC_SUCCESS);
 36: }

 38: static PetscErrorCode PCView_Redundant(PC pc, PetscViewer viewer)
 39: {
 40:   PC_Redundant *red = (PC_Redundant *)pc->data;
 41:   PetscBool     iascii, isstring;
 42:   PetscViewer   subviewer;

 44:   PetscFunctionBegin;
 45:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
 46:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring));
 47:   if (iascii) {
 48:     if (!red->psubcomm) {
 49:       PetscCall(PetscViewerASCIIPrintf(viewer, "  Not yet setup\n"));
 50:     } else {
 51:       PetscCall(PetscViewerASCIIPrintf(viewer, "  First (color=0) of %" PetscInt_FMT " PCs follows\n", red->nsubcomm));
 52:       PetscCall(PetscViewerGetSubViewer(viewer, ((PetscObject)red->pc)->comm, &subviewer));
 53:       if (!red->psubcomm->color) { /* only view first redundant pc */
 54:         PetscCall(PetscViewerASCIIPushTab(subviewer));
 55:         PetscCall(KSPView(red->ksp, subviewer));
 56:         PetscCall(PetscViewerASCIIPopTab(subviewer));
 57:       }
 58:       PetscCall(PetscViewerRestoreSubViewer(viewer, ((PetscObject)red->pc)->comm, &subviewer));
 59:     }
 60:   } else if (isstring) {
 61:     PetscCall(PetscViewerStringSPrintf(viewer, " Redundant solver preconditioner"));
 62:   }
 63:   PetscFunctionReturn(PETSC_SUCCESS);
 64: }

 66: #include <../src/mat/impls/aij/mpi/mpiaij.h>
 67: static PetscErrorCode PCSetUp_Redundant(PC pc)
 68: {
 69:   PC_Redundant *red = (PC_Redundant *)pc->data;
 70:   PetscInt      mstart, mend, mlocal, M;
 71:   PetscMPIInt   size;
 72:   MPI_Comm      comm, subcomm;
 73:   Vec           x;

 75:   PetscFunctionBegin;
 76:   PetscCall(PetscObjectGetComm((PetscObject)pc, &comm));

 78:   /* if pmatrix set by user is sequential then we do not need to gather the parallel matrix */
 79:   PetscCallMPI(MPI_Comm_size(comm, &size));
 80:   if (size == 1) red->useparallelmat = PETSC_FALSE;

 82:   if (!pc->setupcalled) {
 83:     PetscInt mloc_sub;
 84:     if (!red->psubcomm) { /* create red->psubcomm, new ksp and pc over subcomm */
 85:       KSP ksp;
 86:       PetscCall(PCRedundantGetKSP(pc, &ksp));
 87:     }
 88:     subcomm = PetscSubcommChild(red->psubcomm);

 90:     if (red->useparallelmat) {
 91:       /* grab the parallel matrix and put it into processors of a subcomminicator */
 92:       PetscCall(MatCreateRedundantMatrix(pc->pmat, red->psubcomm->n, subcomm, MAT_INITIAL_MATRIX, &red->pmats));

 94:       PetscCallMPI(MPI_Comm_size(subcomm, &size));
 95:       if (size > 1) {
 96:         PetscBool foundpack, issbaij;
 97:         PetscCall(PetscObjectTypeCompare((PetscObject)red->pmats, MATMPISBAIJ, &issbaij));
 98:         if (!issbaij) {
 99:           PetscCall(MatGetFactorAvailable(red->pmats, NULL, MAT_FACTOR_LU, &foundpack));
100:         } else {
101:           PetscCall(MatGetFactorAvailable(red->pmats, NULL, MAT_FACTOR_CHOLESKY, &foundpack));
102:         }
103:         if (!foundpack) { /* reset default ksp and pc */
104:           PetscCall(KSPSetType(red->ksp, KSPGMRES));
105:           PetscCall(PCSetType(red->pc, PCBJACOBI));
106:         } else {
107:           PetscCall(PCFactorSetMatSolverType(red->pc, NULL));
108:         }
109:       }

111:       PetscCall(KSPSetOperators(red->ksp, red->pmats, red->pmats));

113:       /* get working vectors xsub and ysub */
114:       PetscCall(MatCreateVecs(red->pmats, &red->xsub, &red->ysub));

116:       /* create working vectors xdup and ydup.
117:        xdup concatenates all xsub's contigously to form a mpi vector over dupcomm  (see PetscSubcommCreate_interlaced())
118:        ydup concatenates all ysub and has empty local arrays because ysub's arrays will be place into it.
119:        Note: we use communicator dupcomm, not PetscObjectComm((PetscObject)pc)! */
120:       PetscCall(MatGetLocalSize(red->pmats, &mloc_sub, NULL));
121:       PetscCall(VecCreateMPI(PetscSubcommContiguousParent(red->psubcomm), mloc_sub, PETSC_DECIDE, &red->xdup));
122:       PetscCall(VecCreateMPIWithArray(PetscSubcommContiguousParent(red->psubcomm), 1, mloc_sub, PETSC_DECIDE, NULL, &red->ydup));

124:       /* create vecscatters */
125:       if (!red->scatterin) { /* efficiency of scatterin is independent from psubcomm_type! */
126:         IS        is1, is2;
127:         PetscInt *idx1, *idx2, i, j, k;

129:         PetscCall(MatCreateVecs(pc->pmat, &x, NULL));
130:         PetscCall(VecGetSize(x, &M));
131:         PetscCall(VecGetOwnershipRange(x, &mstart, &mend));
132:         mlocal = mend - mstart;
133:         PetscCall(PetscMalloc2(red->psubcomm->n * mlocal, &idx1, red->psubcomm->n * mlocal, &idx2));
134:         j = 0;
135:         for (k = 0; k < red->psubcomm->n; k++) {
136:           for (i = mstart; i < mend; i++) {
137:             idx1[j]   = i;
138:             idx2[j++] = i + M * k;
139:           }
140:         }
141:         PetscCall(ISCreateGeneral(comm, red->psubcomm->n * mlocal, idx1, PETSC_COPY_VALUES, &is1));
142:         PetscCall(ISCreateGeneral(comm, red->psubcomm->n * mlocal, idx2, PETSC_COPY_VALUES, &is2));
143:         PetscCall(VecScatterCreate(x, is1, red->xdup, is2, &red->scatterin));
144:         PetscCall(ISDestroy(&is1));
145:         PetscCall(ISDestroy(&is2));

147:         /* Impl below is good for PETSC_SUBCOMM_INTERLACED (no inter-process communication) and PETSC_SUBCOMM_CONTIGUOUS (communication within subcomm) */
148:         PetscCall(ISCreateStride(comm, mlocal, mstart + red->psubcomm->color * M, 1, &is1));
149:         PetscCall(ISCreateStride(comm, mlocal, mstart, 1, &is2));
150:         PetscCall(VecScatterCreate(red->xdup, is1, x, is2, &red->scatterout));
151:         PetscCall(ISDestroy(&is1));
152:         PetscCall(ISDestroy(&is2));
153:         PetscCall(PetscFree2(idx1, idx2));
154:         PetscCall(VecDestroy(&x));
155:       }
156:     } else { /* !red->useparallelmat */
157:       PetscCall(KSPSetOperators(red->ksp, pc->mat, pc->pmat));
158:     }
159:   } else { /* pc->setupcalled */
160:     if (red->useparallelmat) {
161:       MatReuse reuse;
162:       /* grab the parallel matrix and put it into processors of a subcomminicator */
163:       /*--------------------------------------------------------------------------*/
164:       if (pc->flag == DIFFERENT_NONZERO_PATTERN) {
165:         /* destroy old matrices */
166:         PetscCall(MatDestroy(&red->pmats));
167:         reuse = MAT_INITIAL_MATRIX;
168:       } else {
169:         reuse = MAT_REUSE_MATRIX;
170:       }
171:       PetscCall(MatCreateRedundantMatrix(pc->pmat, red->psubcomm->n, PetscSubcommChild(red->psubcomm), reuse, &red->pmats));
172:       PetscCall(KSPSetOperators(red->ksp, red->pmats, red->pmats));
173:     } else { /* !red->useparallelmat */
174:       PetscCall(KSPSetOperators(red->ksp, pc->mat, pc->pmat));
175:     }
176:   }

178:   if (pc->setfromoptionscalled) PetscCall(KSPSetFromOptions(red->ksp));
179:   PetscCall(KSPSetUp(red->ksp));
180:   PetscFunctionReturn(PETSC_SUCCESS);
181: }

183: static PetscErrorCode PCApply_Redundant(PC pc, Vec x, Vec y)
184: {
185:   PC_Redundant *red = (PC_Redundant *)pc->data;
186:   PetscScalar  *array;

188:   PetscFunctionBegin;
189:   if (!red->useparallelmat) {
190:     PetscCall(KSPSolve(red->ksp, x, y));
191:     PetscCall(KSPCheckSolve(red->ksp, pc, y));
192:     PetscFunctionReturn(PETSC_SUCCESS);
193:   }

195:   /* scatter x to xdup */
196:   PetscCall(VecScatterBegin(red->scatterin, x, red->xdup, INSERT_VALUES, SCATTER_FORWARD));
197:   PetscCall(VecScatterEnd(red->scatterin, x, red->xdup, INSERT_VALUES, SCATTER_FORWARD));

199:   /* place xdup's local array into xsub */
200:   PetscCall(VecGetArray(red->xdup, &array));
201:   PetscCall(VecPlaceArray(red->xsub, (const PetscScalar *)array));

203:   /* apply preconditioner on each processor */
204:   PetscCall(KSPSolve(red->ksp, red->xsub, red->ysub));
205:   PetscCall(KSPCheckSolve(red->ksp, pc, red->ysub));
206:   PetscCall(VecResetArray(red->xsub));
207:   PetscCall(VecRestoreArray(red->xdup, &array));

209:   /* place ysub's local array into ydup */
210:   PetscCall(VecGetArray(red->ysub, &array));
211:   PetscCall(VecPlaceArray(red->ydup, (const PetscScalar *)array));

213:   /* scatter ydup to y */
214:   PetscCall(VecScatterBegin(red->scatterout, red->ydup, y, INSERT_VALUES, SCATTER_FORWARD));
215:   PetscCall(VecScatterEnd(red->scatterout, red->ydup, y, INSERT_VALUES, SCATTER_FORWARD));
216:   PetscCall(VecResetArray(red->ydup));
217:   PetscCall(VecRestoreArray(red->ysub, &array));
218:   PetscFunctionReturn(PETSC_SUCCESS);
219: }

221: static PetscErrorCode PCApplyTranspose_Redundant(PC pc, Vec x, Vec y)
222: {
223:   PC_Redundant *red = (PC_Redundant *)pc->data;
224:   PetscScalar  *array;

226:   PetscFunctionBegin;
227:   if (!red->useparallelmat) {
228:     PetscCall(KSPSolveTranspose(red->ksp, x, y));
229:     PetscCall(KSPCheckSolve(red->ksp, pc, y));
230:     PetscFunctionReturn(PETSC_SUCCESS);
231:   }

233:   /* scatter x to xdup */
234:   PetscCall(VecScatterBegin(red->scatterin, x, red->xdup, INSERT_VALUES, SCATTER_FORWARD));
235:   PetscCall(VecScatterEnd(red->scatterin, x, red->xdup, INSERT_VALUES, SCATTER_FORWARD));

237:   /* place xdup's local array into xsub */
238:   PetscCall(VecGetArray(red->xdup, &array));
239:   PetscCall(VecPlaceArray(red->xsub, (const PetscScalar *)array));

241:   /* apply preconditioner on each processor */
242:   PetscCall(KSPSolveTranspose(red->ksp, red->xsub, red->ysub));
243:   PetscCall(KSPCheckSolve(red->ksp, pc, red->ysub));
244:   PetscCall(VecResetArray(red->xsub));
245:   PetscCall(VecRestoreArray(red->xdup, &array));

247:   /* place ysub's local array into ydup */
248:   PetscCall(VecGetArray(red->ysub, &array));
249:   PetscCall(VecPlaceArray(red->ydup, (const PetscScalar *)array));

251:   /* scatter ydup to y */
252:   PetscCall(VecScatterBegin(red->scatterout, red->ydup, y, INSERT_VALUES, SCATTER_FORWARD));
253:   PetscCall(VecScatterEnd(red->scatterout, red->ydup, y, INSERT_VALUES, SCATTER_FORWARD));
254:   PetscCall(VecResetArray(red->ydup));
255:   PetscCall(VecRestoreArray(red->ysub, &array));
256:   PetscFunctionReturn(PETSC_SUCCESS);
257: }

259: static PetscErrorCode PCReset_Redundant(PC pc)
260: {
261:   PC_Redundant *red = (PC_Redundant *)pc->data;

263:   PetscFunctionBegin;
264:   if (red->useparallelmat) {
265:     PetscCall(VecScatterDestroy(&red->scatterin));
266:     PetscCall(VecScatterDestroy(&red->scatterout));
267:     PetscCall(VecDestroy(&red->ysub));
268:     PetscCall(VecDestroy(&red->xsub));
269:     PetscCall(VecDestroy(&red->xdup));
270:     PetscCall(VecDestroy(&red->ydup));
271:   }
272:   PetscCall(MatDestroy(&red->pmats));
273:   PetscCall(KSPReset(red->ksp));
274:   PetscFunctionReturn(PETSC_SUCCESS);
275: }

277: static PetscErrorCode PCDestroy_Redundant(PC pc)
278: {
279:   PC_Redundant *red = (PC_Redundant *)pc->data;

281:   PetscFunctionBegin;
282:   PetscCall(PCReset_Redundant(pc));
283:   PetscCall(KSPDestroy(&red->ksp));
284:   PetscCall(PetscSubcommDestroy(&red->psubcomm));
285:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantSetScatter_C", NULL));
286:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantSetNumber_C", NULL));
287:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantGetKSP_C", NULL));
288:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantGetOperators_C", NULL));
289:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetShiftType_C", NULL));
290:   PetscCall(PetscFree(pc->data));
291:   PetscFunctionReturn(PETSC_SUCCESS);
292: }

294: static PetscErrorCode PCSetFromOptions_Redundant(PC pc, PetscOptionItems *PetscOptionsObject)
295: {
296:   PC_Redundant *red = (PC_Redundant *)pc->data;

298:   PetscFunctionBegin;
299:   PetscOptionsHeadBegin(PetscOptionsObject, "Redundant options");
300:   PetscCall(PetscOptionsInt("-pc_redundant_number", "Number of redundant pc", "PCRedundantSetNumber", red->nsubcomm, &red->nsubcomm, NULL));
301:   PetscOptionsHeadEnd();
302:   PetscFunctionReturn(PETSC_SUCCESS);
303: }

305: static PetscErrorCode PCRedundantSetNumber_Redundant(PC pc, PetscInt nreds)
306: {
307:   PC_Redundant *red = (PC_Redundant *)pc->data;

309:   PetscFunctionBegin;
310:   red->nsubcomm = nreds;
311:   PetscFunctionReturn(PETSC_SUCCESS);
312: }

314: /*@
315:    PCRedundantSetNumber - Sets the number of redundant preconditioner contexts.

317:    Logically Collective

319:    Input Parameters:
320: +  pc - the preconditioner context
321: -  nredundant - number of redundant preconditioner contexts; for example if you are using 64 MPI processes and
322:                               use an nredundant of 4 there will be 4 parallel solves each on 16 = 64/4 processes.

324:    Level: advanced

326: .seealso: `PCREDUNDANT`
327: @*/
328: PetscErrorCode PCRedundantSetNumber(PC pc, PetscInt nredundant)
329: {
330:   PetscFunctionBegin;
332:   PetscCheck(nredundant > 0, PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_WRONG, "num of redundant pc %" PetscInt_FMT " must be positive", nredundant);
333:   PetscTryMethod(pc, "PCRedundantSetNumber_C", (PC, PetscInt), (pc, nredundant));
334:   PetscFunctionReturn(PETSC_SUCCESS);
335: }

337: static PetscErrorCode PCRedundantSetScatter_Redundant(PC pc, VecScatter in, VecScatter out)
338: {
339:   PC_Redundant *red = (PC_Redundant *)pc->data;

341:   PetscFunctionBegin;
342:   PetscCall(PetscObjectReference((PetscObject)in));
343:   PetscCall(VecScatterDestroy(&red->scatterin));

345:   red->scatterin = in;

347:   PetscCall(PetscObjectReference((PetscObject)out));
348:   PetscCall(VecScatterDestroy(&red->scatterout));
349:   red->scatterout = out;
350:   PetscFunctionReturn(PETSC_SUCCESS);
351: }

353: /*@
354:    PCRedundantSetScatter - Sets the scatter used to copy values into the
355:      redundant local solve and the scatter to move them back into the global
356:      vector.

358:    Logically Collective

360:    Input Parameters:
361: +  pc - the preconditioner context
362: .  in - the scatter to move the values in
363: -  out - the scatter to move them out

365:    Level: advanced

367: .seealso: `PCREDUNDANT`
368: @*/
369: PetscErrorCode PCRedundantSetScatter(PC pc, VecScatter in, VecScatter out)
370: {
371:   PetscFunctionBegin;
375:   PetscTryMethod(pc, "PCRedundantSetScatter_C", (PC, VecScatter, VecScatter), (pc, in, out));
376:   PetscFunctionReturn(PETSC_SUCCESS);
377: }

379: static PetscErrorCode PCRedundantGetKSP_Redundant(PC pc, KSP *innerksp)
380: {
381:   PC_Redundant *red = (PC_Redundant *)pc->data;
382:   MPI_Comm      comm, subcomm;
383:   const char   *prefix;
384:   PetscBool     issbaij;

386:   PetscFunctionBegin;
387:   if (!red->psubcomm) {
388:     PetscCall(PCGetOptionsPrefix(pc, &prefix));

390:     PetscCall(PetscObjectGetComm((PetscObject)pc, &comm));
391:     PetscCall(PetscSubcommCreate(comm, &red->psubcomm));
392:     PetscCall(PetscSubcommSetNumber(red->psubcomm, red->nsubcomm));
393:     PetscCall(PetscSubcommSetType(red->psubcomm, PETSC_SUBCOMM_CONTIGUOUS));

395:     PetscCall(PetscSubcommSetOptionsPrefix(red->psubcomm, prefix));
396:     PetscCall(PetscSubcommSetFromOptions(red->psubcomm));

398:     /* create a new PC that processors in each subcomm have copy of */
399:     subcomm = PetscSubcommChild(red->psubcomm);

401:     PetscCall(KSPCreate(subcomm, &red->ksp));
402:     PetscCall(KSPSetErrorIfNotConverged(red->ksp, pc->erroriffailure));
403:     PetscCall(PetscObjectIncrementTabLevel((PetscObject)red->ksp, (PetscObject)pc, 1));
404:     PetscCall(KSPSetType(red->ksp, KSPPREONLY));
405:     PetscCall(KSPGetPC(red->ksp, &red->pc));
406:     PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATSEQSBAIJ, &issbaij));
407:     if (!issbaij) PetscCall(PetscObjectTypeCompare((PetscObject)pc->pmat, MATMPISBAIJ, &issbaij));
408:     if (!issbaij) {
409:       PetscCall(PCSetType(red->pc, PCLU));
410:     } else {
411:       PetscCall(PCSetType(red->pc, PCCHOLESKY));
412:     }
413:     if (red->shifttypeset) {
414:       PetscCall(PCFactorSetShiftType(red->pc, red->shifttype));
415:       red->shifttypeset = PETSC_FALSE;
416:     }
417:     PetscCall(KSPSetOptionsPrefix(red->ksp, prefix));
418:     PetscCall(KSPAppendOptionsPrefix(red->ksp, "redundant_"));
419:   }
420:   *innerksp = red->ksp;
421:   PetscFunctionReturn(PETSC_SUCCESS);
422: }

424: /*@
425:    PCRedundantGetKSP - Gets the less parallel `KSP` created by the redundant `PC`.

427:    Not Collective

429:    Input Parameter:
430: .  pc - the preconditioner context

432:    Output Parameter:
433: .  innerksp - the `KSP` on the smaller set of processes

435:    Level: advanced

437: .seealso: `PCREDUNDANT`
438: @*/
439: PetscErrorCode PCRedundantGetKSP(PC pc, KSP *innerksp)
440: {
441:   PetscFunctionBegin;
444:   PetscUseMethod(pc, "PCRedundantGetKSP_C", (PC, KSP *), (pc, innerksp));
445:   PetscFunctionReturn(PETSC_SUCCESS);
446: }

448: static PetscErrorCode PCRedundantGetOperators_Redundant(PC pc, Mat *mat, Mat *pmat)
449: {
450:   PC_Redundant *red = (PC_Redundant *)pc->data;

452:   PetscFunctionBegin;
453:   if (mat) *mat = red->pmats;
454:   if (pmat) *pmat = red->pmats;
455:   PetscFunctionReturn(PETSC_SUCCESS);
456: }

458: /*@
459:    PCRedundantGetOperators - gets the sequential matrix and preconditioner matrix

461:    Not Collective

463:    Input Parameter:
464: .  pc - the preconditioner context

466:    Output Parameters:
467: +  mat - the matrix
468: -  pmat - the (possibly different) preconditioner matrix

470:    Level: advanced

472: .seealso: `PCREDUNDANT`
473: @*/
474: PetscErrorCode PCRedundantGetOperators(PC pc, Mat *mat, Mat *pmat)
475: {
476:   PetscFunctionBegin;
480:   PetscUseMethod(pc, "PCRedundantGetOperators_C", (PC, Mat *, Mat *), (pc, mat, pmat));
481:   PetscFunctionReturn(PETSC_SUCCESS);
482: }

484: /*MC
485:      PCREDUNDANT - Runs a `KSP` solver with preconditioner for the entire problem on subgroups of processors

487:   Options Database Key:
488: .  -pc_redundant_number <n> - number of redundant solves, for example if you are using 64 MPI processes and
489:                               use an n of 4 there will be 4 parallel solves each on 16 = 64/4 processes.

491:    Level: intermediate

493:    Notes:
494:    Options for the redundant preconditioners can be set using the options database prefix -redundant_

496:    The default `KSP` is preonly and the default `PC` is `PCLU` or `PCCHOLESKY` if Pmat is of type `MATSBAIJ`.

498:    `PCFactorSetShiftType()` applied to this `PC` will convey they shift type into the inner `PC` if it is factorization based.

500:    Developer Note:
501:    `PCSetInitialGuessNonzero()` is not used by this class but likely should be.

503: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PCRedundantSetScatter()`,
504:           `PCRedundantGetKSP()`, `PCRedundantGetOperators()`, `PCRedundantSetNumber()`, `PCREDISTRIBUTE`
505: M*/

507: PETSC_EXTERN PetscErrorCode PCCreate_Redundant(PC pc)
508: {
509:   PC_Redundant *red;
510:   PetscMPIInt   size;

512:   PetscFunctionBegin;
513:   PetscCall(PetscNew(&red));
514:   PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)pc), &size));

516:   red->nsubcomm       = size;
517:   red->useparallelmat = PETSC_TRUE;
518:   pc->data            = (void *)red;

520:   pc->ops->apply          = PCApply_Redundant;
521:   pc->ops->applytranspose = PCApplyTranspose_Redundant;
522:   pc->ops->setup          = PCSetUp_Redundant;
523:   pc->ops->destroy        = PCDestroy_Redundant;
524:   pc->ops->reset          = PCReset_Redundant;
525:   pc->ops->setfromoptions = PCSetFromOptions_Redundant;
526:   pc->ops->view           = PCView_Redundant;

528:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantSetScatter_C", PCRedundantSetScatter_Redundant));
529:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantSetNumber_C", PCRedundantSetNumber_Redundant));
530:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantGetKSP_C", PCRedundantGetKSP_Redundant));
531:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCRedundantGetOperators_C", PCRedundantGetOperators_Redundant));
532:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCFactorSetShiftType_C", PCFactorSetShiftType_Redundant));
533:   PetscFunctionReturn(PETSC_SUCCESS);
534: }