Actual source code: shellpc.c


  2: /*
  3:    This provides a simple shell for Fortran (and C programmers) to
  4:   create their own preconditioner without writing much interface code.
  5: */

  7: #include <petsc/private/pcimpl.h>

  9: typedef struct {
 10:   void *ctx; /* user provided contexts for preconditioner */

 12:   PetscErrorCode (*destroy)(PC);
 13:   PetscErrorCode (*setup)(PC);
 14:   PetscErrorCode (*apply)(PC, Vec, Vec);
 15:   PetscErrorCode (*matapply)(PC, Mat, Mat);
 16:   PetscErrorCode (*applysymmetricleft)(PC, Vec, Vec);
 17:   PetscErrorCode (*applysymmetricright)(PC, Vec, Vec);
 18:   PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec);
 19:   PetscErrorCode (*presolve)(PC, KSP, Vec, Vec);
 20:   PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec);
 21:   PetscErrorCode (*view)(PC, PetscViewer);
 22:   PetscErrorCode (*applytranspose)(PC, Vec, Vec);
 23:   PetscErrorCode (*applyrich)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *);

 25:   char *name;
 26: } PC_Shell;

 28: /*@C
 29:     PCShellGetContext - Returns the user-provided context associated with a shell `PC`

 31:     Not Collective

 33:     Input Parameter:
 34: .   pc - of type `PCSHELL` created with `PCSetType`(pc,shell)

 36:     Output Parameter:
 37: .   ctx - the user provided context

 39:     Level: advanced

 41:     Note:
 42:     This routine is intended for use within various shell routines

 44:    Fortran Note:
 45:     To use this from Fortran you must write a Fortran interface definition for this
 46:     function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.

 48: .seealso: `PCSHELL`, `PCShellSetContext()`
 49: @*/
 50: PetscErrorCode PCShellGetContext(PC pc, void *ctx)
 51: {
 52:   PetscBool flg;

 54:   PetscFunctionBegin;
 57:   PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg));
 58:   if (!flg) *(void **)ctx = NULL;
 59:   else *(void **)ctx = ((PC_Shell *)(pc->data))->ctx;
 60:   PetscFunctionReturn(PETSC_SUCCESS);
 61: }

 63: /*@
 64:     PCShellSetContext - sets the context for a shell `PC`

 66:    Logically Collective

 68:     Input Parameters:
 69: +   pc - the `PC` of type `PCSHELL`
 70: -   ctx - the context

 72:    Level: advanced

 74:    Fortran Note:
 75:     To use this from Fortran you must write a Fortran interface definition for this
 76:     function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.

 78: .seealso: `PCShellGetContext()`, `PCSHELL`
 79: @*/
 80: PetscErrorCode PCShellSetContext(PC pc, void *ctx)
 81: {
 82:   PC_Shell *shell = (PC_Shell *)pc->data;
 83:   PetscBool flg;

 85:   PetscFunctionBegin;
 87:   PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg));
 88:   if (flg) shell->ctx = ctx;
 89:   PetscFunctionReturn(PETSC_SUCCESS);
 90: }

 92: static PetscErrorCode PCSetUp_Shell(PC pc)
 93: {
 94:   PC_Shell *shell = (PC_Shell *)pc->data;

 96:   PetscFunctionBegin;
 97:   PetscCheck(shell->setup, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No setup() routine provided to Shell PC");
 98:   PetscCallBack("PCSHELL callback setup", (*shell->setup)(pc));
 99:   PetscFunctionReturn(PETSC_SUCCESS);
100: }

102: static PetscErrorCode PCApply_Shell(PC pc, Vec x, Vec y)
103: {
104:   PC_Shell        *shell = (PC_Shell *)pc->data;
105:   PetscObjectState instate, outstate;

107:   PetscFunctionBegin;
108:   PetscCheck(shell->apply, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
109:   PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
110:   PetscCallBack("PCSHELL callback apply", (*shell->apply)(pc, x, y));
111:   PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
112:   /* increase the state of the output vector if the user did not update its state themself as should have been done */
113:   if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
114:   PetscFunctionReturn(PETSC_SUCCESS);
115: }

117: static PetscErrorCode PCMatApply_Shell(PC pc, Mat X, Mat Y)
118: {
119:   PC_Shell        *shell = (PC_Shell *)pc->data;
120:   PetscObjectState instate, outstate;

122:   PetscFunctionBegin;
123:   PetscCheck(shell->matapply, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
124:   PetscCall(PetscObjectStateGet((PetscObject)Y, &instate));
125:   PetscCallBack("PCSHELL callback apply", (*shell->matapply)(pc, X, Y));
126:   PetscCall(PetscObjectStateGet((PetscObject)Y, &outstate));
127:   /* increase the state of the output vector if the user did not update its state themself as should have been done */
128:   if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)Y));
129:   PetscFunctionReturn(PETSC_SUCCESS);
130: }

132: static PetscErrorCode PCApplySymmetricLeft_Shell(PC pc, Vec x, Vec y)
133: {
134:   PC_Shell *shell = (PC_Shell *)pc->data;

136:   PetscFunctionBegin;
137:   PetscCheck(shell->applysymmetricleft, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
138:   PetscCallBack("PCSHELL callback apply symmetric left", (*shell->applysymmetricleft)(pc, x, y));
139:   PetscFunctionReturn(PETSC_SUCCESS);
140: }

142: static PetscErrorCode PCApplySymmetricRight_Shell(PC pc, Vec x, Vec y)
143: {
144:   PC_Shell *shell = (PC_Shell *)pc->data;

146:   PetscFunctionBegin;
147:   PetscCheck(shell->applysymmetricright, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
148:   PetscCallBack("PCSHELL callback apply symmetric right", (*shell->applysymmetricright)(pc, x, y));
149:   PetscFunctionReturn(PETSC_SUCCESS);
150: }

152: static PetscErrorCode PCApplyBA_Shell(PC pc, PCSide side, Vec x, Vec y, Vec w)
153: {
154:   PC_Shell        *shell = (PC_Shell *)pc->data;
155:   PetscObjectState instate, outstate;

157:   PetscFunctionBegin;
158:   PetscCheck(shell->applyBA, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applyBA() routine provided to Shell PC");
159:   PetscCall(PetscObjectStateGet((PetscObject)w, &instate));
160:   PetscCallBack("PCSHELL callback applyBA", (*shell->applyBA)(pc, side, x, y, w));
161:   PetscCall(PetscObjectStateGet((PetscObject)w, &outstate));
162:   /* increase the state of the output vector if the user did not update its state themself as should have been done */
163:   if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)w));
164:   PetscFunctionReturn(PETSC_SUCCESS);
165: }

167: static PetscErrorCode PCPreSolveChangeRHS_Shell(PC pc, PetscBool *change)
168: {
169:   PetscFunctionBegin;
170:   *change = PETSC_TRUE;
171:   PetscFunctionReturn(PETSC_SUCCESS);
172: }

174: static PetscErrorCode PCPreSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
175: {
176:   PC_Shell *shell = (PC_Shell *)pc->data;

178:   PetscFunctionBegin;
179:   PetscCheck(shell->presolve, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No presolve() routine provided to Shell PC");
180:   PetscCallBack("PCSHELL callback presolve", (*shell->presolve)(pc, ksp, b, x));
181:   PetscFunctionReturn(PETSC_SUCCESS);
182: }

184: static PetscErrorCode PCPostSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
185: {
186:   PC_Shell *shell = (PC_Shell *)pc->data;

188:   PetscFunctionBegin;
189:   PetscCheck(shell->postsolve, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No postsolve() routine provided to Shell PC");
190:   PetscCallBack("PCSHELL callback postsolve()", (*shell->postsolve)(pc, ksp, b, x));
191:   PetscFunctionReturn(PETSC_SUCCESS);
192: }

194: static PetscErrorCode PCApplyTranspose_Shell(PC pc, Vec x, Vec y)
195: {
196:   PC_Shell        *shell = (PC_Shell *)pc->data;
197:   PetscObjectState instate, outstate;

199:   PetscFunctionBegin;
200:   PetscCheck(shell->applytranspose, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applytranspose() routine provided to Shell PC");
201:   PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
202:   PetscCallBack("PCSHELL callback applytranspose", (*shell->applytranspose)(pc, x, y));
203:   PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
204:   /* increase the state of the output vector if the user did not update its state themself as should have been done */
205:   if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
206:   PetscFunctionReturn(PETSC_SUCCESS);
207: }

209: static PetscErrorCode PCApplyRichardson_Shell(PC pc, Vec x, Vec y, Vec w, PetscReal rtol, PetscReal abstol, PetscReal dtol, PetscInt it, PetscBool guesszero, PetscInt *outits, PCRichardsonConvergedReason *reason)
210: {
211:   PC_Shell        *shell = (PC_Shell *)pc->data;
212:   PetscObjectState instate, outstate;

214:   PetscFunctionBegin;
215:   PetscCheck(shell->applyrich, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applyrichardson() routine provided to Shell PC");
216:   PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
217:   PetscCallBack("PCSHELL callback applyrichardson", (*shell->applyrich)(pc, x, y, w, rtol, abstol, dtol, it, guesszero, outits, reason));
218:   PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
219:   /* increase the state of the output vector since the user did not update its state themself as should have been done */
220:   if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
221:   PetscFunctionReturn(PETSC_SUCCESS);
222: }

224: static PetscErrorCode PCDestroy_Shell(PC pc)
225: {
226:   PC_Shell *shell = (PC_Shell *)pc->data;

228:   PetscFunctionBegin;
229:   PetscCall(PetscFree(shell->name));
230:   if (shell->destroy) PetscCallBack("PCSHELL callback destroy", (*shell->destroy)(pc));
231:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", NULL));
232:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", NULL));
233:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", NULL));
234:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", NULL));
235:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", NULL));
236:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", NULL));
237:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", NULL));
238:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", NULL));
239:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", NULL));
240:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", NULL));
241:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", NULL));
242:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", NULL));
243:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", NULL));
244:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", NULL));
245:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL));
246:   PetscCall(PetscFree(pc->data));
247:   PetscFunctionReturn(PETSC_SUCCESS);
248: }

250: static PetscErrorCode PCView_Shell(PC pc, PetscViewer viewer)
251: {
252:   PC_Shell *shell = (PC_Shell *)pc->data;
253:   PetscBool iascii;

255:   PetscFunctionBegin;
256:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
257:   if (iascii) {
258:     if (shell->name) PetscCall(PetscViewerASCIIPrintf(viewer, "  %s\n", shell->name));
259:     else PetscCall(PetscViewerASCIIPrintf(viewer, "  no name\n"));
260:   }
261:   if (shell->view) {
262:     PetscCall(PetscViewerASCIIPushTab(viewer));
263:     PetscCall((*shell->view)(pc, viewer));
264:     PetscCall(PetscViewerASCIIPopTab(viewer));
265:   }
266:   PetscFunctionReturn(PETSC_SUCCESS);
267: }

269: static PetscErrorCode PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
270: {
271:   PC_Shell *shell = (PC_Shell *)pc->data;

273:   PetscFunctionBegin;
274:   shell->destroy = destroy;
275:   PetscFunctionReturn(PETSC_SUCCESS);
276: }

278: static PetscErrorCode PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
279: {
280:   PC_Shell *shell = (PC_Shell *)pc->data;

282:   PetscFunctionBegin;
283:   shell->setup = setup;
284:   if (setup) pc->ops->setup = PCSetUp_Shell;
285:   else pc->ops->setup = NULL;
286:   PetscFunctionReturn(PETSC_SUCCESS);
287: }

289: static PetscErrorCode PCShellSetApply_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
290: {
291:   PC_Shell *shell = (PC_Shell *)pc->data;

293:   PetscFunctionBegin;
294:   shell->apply = apply;
295:   PetscFunctionReturn(PETSC_SUCCESS);
296: }

298: static PetscErrorCode PCShellSetMatApply_Shell(PC pc, PetscErrorCode (*matapply)(PC, Mat, Mat))
299: {
300:   PC_Shell *shell = (PC_Shell *)pc->data;

302:   PetscFunctionBegin;
303:   shell->matapply = matapply;
304:   if (matapply) pc->ops->matapply = PCMatApply_Shell;
305:   else pc->ops->matapply = NULL;
306:   PetscFunctionReturn(PETSC_SUCCESS);
307: }

309: static PetscErrorCode PCShellSetApplySymmetricLeft_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
310: {
311:   PC_Shell *shell = (PC_Shell *)pc->data;

313:   PetscFunctionBegin;
314:   shell->applysymmetricleft = apply;
315:   PetscFunctionReturn(PETSC_SUCCESS);
316: }

318: static PetscErrorCode PCShellSetApplySymmetricRight_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
319: {
320:   PC_Shell *shell = (PC_Shell *)pc->data;

322:   PetscFunctionBegin;
323:   shell->applysymmetricright = apply;
324:   PetscFunctionReturn(PETSC_SUCCESS);
325: }

327: static PetscErrorCode PCShellSetApplyBA_Shell(PC pc, PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec))
328: {
329:   PC_Shell *shell = (PC_Shell *)pc->data;

331:   PetscFunctionBegin;
332:   shell->applyBA = applyBA;
333:   if (applyBA) pc->ops->applyBA = PCApplyBA_Shell;
334:   else pc->ops->applyBA = NULL;
335:   PetscFunctionReturn(PETSC_SUCCESS);
336: }

338: static PetscErrorCode PCShellSetPreSolve_Shell(PC pc, PetscErrorCode (*presolve)(PC, KSP, Vec, Vec))
339: {
340:   PC_Shell *shell = (PC_Shell *)pc->data;

342:   PetscFunctionBegin;
343:   shell->presolve = presolve;
344:   if (presolve) {
345:     pc->ops->presolve = PCPreSolve_Shell;
346:     PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", PCPreSolveChangeRHS_Shell));
347:   } else {
348:     pc->ops->presolve = NULL;
349:     PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL));
350:   }
351:   PetscFunctionReturn(PETSC_SUCCESS);
352: }

354: static PetscErrorCode PCShellSetPostSolve_Shell(PC pc, PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec))
355: {
356:   PC_Shell *shell = (PC_Shell *)pc->data;

358:   PetscFunctionBegin;
359:   shell->postsolve = postsolve;
360:   if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
361:   else pc->ops->postsolve = NULL;
362:   PetscFunctionReturn(PETSC_SUCCESS);
363: }

365: static PetscErrorCode PCShellSetView_Shell(PC pc, PetscErrorCode (*view)(PC, PetscViewer))
366: {
367:   PC_Shell *shell = (PC_Shell *)pc->data;

369:   PetscFunctionBegin;
370:   shell->view = view;
371:   PetscFunctionReturn(PETSC_SUCCESS);
372: }

374: static PetscErrorCode PCShellSetApplyTranspose_Shell(PC pc, PetscErrorCode (*applytranspose)(PC, Vec, Vec))
375: {
376:   PC_Shell *shell = (PC_Shell *)pc->data;

378:   PetscFunctionBegin;
379:   shell->applytranspose = applytranspose;
380:   if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
381:   else pc->ops->applytranspose = NULL;
382:   PetscFunctionReturn(PETSC_SUCCESS);
383: }

385: static PetscErrorCode PCShellSetApplyRichardson_Shell(PC pc, PetscErrorCode (*applyrich)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *))
386: {
387:   PC_Shell *shell = (PC_Shell *)pc->data;

389:   PetscFunctionBegin;
390:   shell->applyrich = applyrich;
391:   if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
392:   else pc->ops->applyrichardson = NULL;
393:   PetscFunctionReturn(PETSC_SUCCESS);
394: }

396: static PetscErrorCode PCShellSetName_Shell(PC pc, const char name[])
397: {
398:   PC_Shell *shell = (PC_Shell *)pc->data;

400:   PetscFunctionBegin;
401:   PetscCall(PetscFree(shell->name));
402:   PetscCall(PetscStrallocpy(name, &shell->name));
403:   PetscFunctionReturn(PETSC_SUCCESS);
404: }

406: static PetscErrorCode PCShellGetName_Shell(PC pc, const char *name[])
407: {
408:   PC_Shell *shell = (PC_Shell *)pc->data;

410:   PetscFunctionBegin;
411:   *name = shell->name;
412:   PetscFunctionReturn(PETSC_SUCCESS);
413: }

415: /*@C
416:    PCShellSetDestroy - Sets routine to use to destroy the user-provided
417:    application context.

419:    Logically Collective

421:    Input Parameters:
422: +  pc - the preconditioner context
423: -  destroy - the application-provided destroy routine

425:    Calling sequence of destroy:
426: .vb
427:    PetscErrorCode destroy (PC)
428: .ve

430: .  ptr - the application context

432:    Note:
433:     the function MUST return an error code of 0 on success and nonzero on failure.

435:    Level: intermediate

437: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetContext()`
438: @*/
439: PetscErrorCode PCShellSetDestroy(PC pc, PetscErrorCode (*destroy)(PC))
440: {
441:   PetscFunctionBegin;
443:   PetscTryMethod(pc, "PCShellSetDestroy_C", (PC, PetscErrorCode(*)(PC)), (pc, destroy));
444:   PetscFunctionReturn(PETSC_SUCCESS);
445: }

447: /*@C
448:    PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
449:    matrix operator is changed.

451:    Logically Collective

453:    Input Parameters:
454: +  pc - the preconditioner context
455: -  setup - the application-provided setup routine

457:    Calling sequence of setup:
458: .vb
459:    PetscErrorCode setup (PC pc)
460: .ve

462: .  pc - the preconditioner, get the application context with PCShellGetContext()

464:    Note:
465:     the function MUST return an error code of 0 on success and nonzero on failure.

467:    Level: intermediate

469: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetApply()`, `PCShellSetContext()`
470: @*/
471: PetscErrorCode PCShellSetSetUp(PC pc, PetscErrorCode (*setup)(PC))
472: {
473:   PetscFunctionBegin;
475:   PetscTryMethod(pc, "PCShellSetSetUp_C", (PC, PetscErrorCode(*)(PC)), (pc, setup));
476:   PetscFunctionReturn(PETSC_SUCCESS);
477: }

479: /*@C
480:    PCShellSetView - Sets routine to use as viewer of shell preconditioner

482:    Logically Collective

484:    Input Parameters:
485: +  pc - the preconditioner context
486: -  view - the application-provided view routine

488:    Calling sequence of view:
489: .vb
490:    PetscErrorCode view(PC pc,PetscViewer v)
491: .ve

493: +  pc - the preconditioner, get the application context with PCShellGetContext()
494: -  v   - viewer

496:    Note:
497:     the function MUST return an error code of 0 on success and nonzero on failure.

499:    Level: advanced

501: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`
502: @*/
503: PetscErrorCode PCShellSetView(PC pc, PetscErrorCode (*view)(PC, PetscViewer))
504: {
505:   PetscFunctionBegin;
507:   PetscTryMethod(pc, "PCShellSetView_C", (PC, PetscErrorCode(*)(PC, PetscViewer)), (pc, view));
508:   PetscFunctionReturn(PETSC_SUCCESS);
509: }

511: /*@C
512:    PCShellSetApply - Sets routine to use as preconditioner.

514:    Logically Collective

516:    Input Parameters:
517: +  pc - the preconditioner context
518: -  apply - the application-provided preconditioning routine

520:    Calling sequence of apply:
521: .vb
522:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
523: .ve

525: +  pc - the preconditioner, get the application context with PCShellGetContext()
526: .  xin - input vector
527: -  xout - output vector

529:    Note:
530:     the function MUST return an error code of 0 on success and nonzero on failure.

532:    Level: intermediate

534: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApplyBA()`, `PCShellSetApplySymmetricRight()`, `PCShellSetApplySymmetricLeft()`
535: @*/
536: PetscErrorCode PCShellSetApply(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
537: {
538:   PetscFunctionBegin;
540:   PetscTryMethod(pc, "PCShellSetApply_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
541:   PetscFunctionReturn(PETSC_SUCCESS);
542: }

544: /*@C
545:    PCShellSetMatApply - Sets routine to use as preconditioner on a block of vectors.

547:    Logically Collective

549:    Input Parameters:
550: +  pc - the preconditioner context
551: -  apply - the application-provided preconditioning routine

553:    Calling sequence of apply:
554: .vb
555:    PetscErrorCode apply (PC pc,Mat Xin,Mat Xout)
556: .ve

558: +  pc - the preconditioner, get the application context with PCShellGetContext()
559: .  Xin - input block of vectors
560: -  Xout - output block of vectors

562:    Note:
563:    The function MUST return an error code of 0 on success and nonzero on failure.

565:    Level: advanced

567: .seealso: `PCSHELL`, `PCShellSetApply()`
568: @*/
569: PetscErrorCode PCShellSetMatApply(PC pc, PetscErrorCode (*matapply)(PC, Mat, Mat))
570: {
571:   PetscFunctionBegin;
573:   PetscTryMethod(pc, "PCShellSetMatApply_C", (PC, PetscErrorCode(*)(PC, Mat, Mat)), (pc, matapply));
574:   PetscFunctionReturn(PETSC_SUCCESS);
575: }

577: /*@C
578:    PCShellSetApplySymmetricLeft - Sets routine to use as left preconditioner (when the `PC_SYMMETRIC` is used).

580:    Logically Collective

582:    Input Parameters:
583: +  pc - the preconditioner context
584: -  apply - the application-provided left preconditioning routine

586:    Calling sequence of apply:
587: .vb
588:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
589: .ve

591: +  pc - the preconditioner, get the application context with `PCShellGetContext()`
592: .  xin - input vector
593: -  xout - output vector

595:    Note:
596:   The function MUST return an error code of 0 on success and nonzero on failure.

598:    Level: advanced

600: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetApplySymmetricLeft()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`
601: @*/
602: PetscErrorCode PCShellSetApplySymmetricLeft(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
603: {
604:   PetscFunctionBegin;
606:   PetscTryMethod(pc, "PCShellSetApplySymmetricLeft_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
607:   PetscFunctionReturn(PETSC_SUCCESS);
608: }

610: /*@C
611:    PCShellSetApplySymmetricRight - Sets routine to use as right preconditioner (when the `PC_SYMMETRIC` is used).

613:    Logically Collective

615:    Input Parameters:
616: +  pc - the preconditioner context
617: -  apply - the application-provided right preconditioning routine

619:    Calling sequence of apply:
620: .vb
621:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
622: .ve

624: +  pc - the preconditioner, get the application context with PCShellGetContext()
625: .  xin - input vector
626: -  xout - output vector

628:    Note:
629:    The function MUST return an error code of 0 on success and nonzero on failure.

631:    Level: advanced

633: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetApplySymmetricLeft()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`
634: @*/
635: PetscErrorCode PCShellSetApplySymmetricRight(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
636: {
637:   PetscFunctionBegin;
639:   PetscTryMethod(pc, "PCShellSetApplySymmetricRight_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
640:   PetscFunctionReturn(PETSC_SUCCESS);
641: }

643: /*@C
644:    PCShellSetApplyBA - Sets routine to use as preconditioner times operator.

646:    Logically Collective

648:    Input Parameters:
649: +  pc - the preconditioner context
650: -  applyBA - the application-provided BA routine

652:    Calling sequence of applyBA:
653: .vb
654:    PetscErrorCode applyBA (PC pc,Vec xin,Vec xout)
655: .ve

657: +  pc - the preconditioner, get the application context with PCShellGetContext()
658: .  xin - input vector
659: -  xout - output vector

661:    Note:
662:    The function MUST return an error code of 0 on success and nonzero on failure.

664:    Level: intermediate

666: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApply()`
667: @*/
668: PetscErrorCode PCShellSetApplyBA(PC pc, PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec))
669: {
670:   PetscFunctionBegin;
672:   PetscTryMethod(pc, "PCShellSetApplyBA_C", (PC, PetscErrorCode(*)(PC, PCSide, Vec, Vec, Vec)), (pc, applyBA));
673:   PetscFunctionReturn(PETSC_SUCCESS);
674: }

676: /*@C
677:    PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.

679:    Logically Collective

681:    Input Parameters:
682: +  pc - the preconditioner context
683: -  apply - the application-provided preconditioning transpose routine

685:    Calling sequence of apply:
686: .vb
687:    PetscErrorCode applytranspose (PC pc,Vec xin,Vec xout)
688: .ve

690: +  pc - the preconditioner, get the application context with PCShellGetContext()
691: .  xin - input vector
692: -  xout - output vector

694:    Note:
695:     the function MUST return an error code of 0 on success and nonzero on failure.

697:    Level: intermediate

699:    Note:
700:    Uses the same context variable as `PCShellSetApply()`.

702: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApply()`, `PCSetContext()`, `PCShellSetApplyBA()`
703: @*/
704: PetscErrorCode PCShellSetApplyTranspose(PC pc, PetscErrorCode (*applytranspose)(PC, Vec, Vec))
705: {
706:   PetscFunctionBegin;
708:   PetscTryMethod(pc, "PCShellSetApplyTranspose_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, applytranspose));
709:   PetscFunctionReturn(PETSC_SUCCESS);
710: }

712: /*@C
713:    PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a `KSPSolve()` is
714:       applied. This usually does something like scale the linear system in some application
715:       specific way.

717:    Logically Collective

719:    Input Parameters:
720: +  pc - the preconditioner context
721: -  presolve - the application-provided presolve routine

723:    Calling sequence of presolve:
724: .vb
725:    PetscErrorCode presolve (PC,KSP ksp,Vec b,Vec x)
726: .ve

728: +  pc - the preconditioner, get the application context with PCShellGetContext()
729: .  xin - input vector
730: -  xout - output vector

732:    Note:
733:    The function MUST return an error code of 0 on success and nonzero on failure.

735:    Level: advanced

737: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPostSolve()`, `PCShellSetContext()`
738: @*/
739: PetscErrorCode PCShellSetPreSolve(PC pc, PetscErrorCode (*presolve)(PC, KSP, Vec, Vec))
740: {
741:   PetscFunctionBegin;
743:   PetscTryMethod(pc, "PCShellSetPreSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, presolve));
744:   PetscFunctionReturn(PETSC_SUCCESS);
745: }

747: /*@C
748:    PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a `KSPSolve()` is
749:       applied. This usually does something like scale the linear system in some application
750:       specific way.

752:    Logically Collective

754:    Input Parameters:
755: +  pc - the preconditioner context
756: -  postsolve - the application-provided presolve routine

758:    Calling sequence of postsolve:
759: .vb
760:    PetscErrorCode postsolve(PC,KSP ksp,Vec b,Vec x)
761: .ve

763: +  pc - the preconditioner, get the application context with `PCShellGetContext()`
764: .  xin - input vector
765: -  xout - output vector

767:    Note:
768:     the function MUST return an error code of 0 on success and nonzero on failure.

770:    Level: advanced

772: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPreSolve()`, `PCShellSetContext()`
773: @*/
774: PetscErrorCode PCShellSetPostSolve(PC pc, PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec))
775: {
776:   PetscFunctionBegin;
778:   PetscTryMethod(pc, "PCShellSetPostSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, postsolve));
779:   PetscFunctionReturn(PETSC_SUCCESS);
780: }

782: /*@C
783:    PCShellSetName - Sets an optional name to associate with a `PCSHELL`
784:    preconditioner.

786:    Not Collective

788:    Input Parameters:
789: +  pc - the preconditioner context
790: -  name - character string describing shell preconditioner

792:    Level: intermediate

794: .seealso: `PCSHELL`, `PCShellGetName()`
795: @*/
796: PetscErrorCode PCShellSetName(PC pc, const char name[])
797: {
798:   PetscFunctionBegin;
800:   PetscTryMethod(pc, "PCShellSetName_C", (PC, const char[]), (pc, name));
801:   PetscFunctionReturn(PETSC_SUCCESS);
802: }

804: /*@C
805:    PCShellGetName - Gets an optional name that the user has set for a `PCSHELL`
806:    preconditioner.

808:    Not Collective

810:    Input Parameter:
811: .  pc - the preconditioner context

813:    Output Parameter:
814: .  name - character string describing shell preconditioner (you should not free this)

816:    Level: intermediate

818: .seealso: `PCSHELL`, `PCShellSetName()`
819: @*/
820: PetscErrorCode PCShellGetName(PC pc, const char *name[])
821: {
822:   PetscFunctionBegin;
825:   PetscUseMethod(pc, "PCShellGetName_C", (PC, const char *[]), (pc, name));
826:   PetscFunctionReturn(PETSC_SUCCESS);
827: }

829: /*@C
830:    PCShellSetApplyRichardson - Sets routine to use as preconditioner
831:    in Richardson iteration.

833:    Logically Collective

835:    Input Parameters:
836: +  pc - the preconditioner context
837: -  apply - the application-provided preconditioning routine

839:    Calling sequence of apply:
840: .vb
841:    PetscErrorCode apply (PC pc,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
842: .ve

844: +  pc - the preconditioner, get the application context with PCShellGetContext()
845: .  b - right-hand-side
846: .  x - current iterate
847: .  r - work space
848: .  rtol - relative tolerance of residual norm to stop at
849: .  abstol - absolute tolerance of residual norm to stop at
850: .  dtol - if residual norm increases by this factor than return
851: -  maxits - number of iterations to run

853:    Note:
854:     the function MUST return an error code of 0 on success and nonzero on failure.

856:    Level: advanced

858: .seealso: `PCShellSetApply()`, `PCShellSetContext()`
859: @*/
860: PetscErrorCode PCShellSetApplyRichardson(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *))
861: {
862:   PetscFunctionBegin;
864:   PetscTryMethod(pc, "PCShellSetApplyRichardson_C", (PC, PetscErrorCode(*)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *)), (pc, apply));
865:   PetscFunctionReturn(PETSC_SUCCESS);
866: }

868: /*MC
869:    PCSHELL - Creates a new preconditioner class for use with a users
870:               own private data storage format and preconditioner application code

872:    Level: advanced

874:   Usage:
875: .vb
876:        extern PetscErrorCode apply(PC,Vec,Vec);
877:        extern PetscErrorCode applyba(PC,PCSide,Vec,Vec,Vec);
878:        extern PetscErrorCode applytranspose(PC,Vec,Vec);
879:        extern PetscErrorCode setup(PC);
880:        extern PetscErrorCode destroy(PC);

882:        PCCreate(comm,&pc);
883:        PCSetType(pc,PCSHELL);
884:        PCShellSetContext(pc,ctx)
885:        PCShellSetApply(pc,apply);
886:        PCShellSetApplyBA(pc,applyba);               (optional)
887:        PCShellSetApplyTranspose(pc,applytranspose); (optional)
888:        PCShellSetSetUp(pc,setup);                   (optional)
889:        PCShellSetDestroy(pc,destroy);               (optional)
890: .ve

892: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
893:           `MATSHELL`, `PCShellSetSetUp()`, `PCShellSetApply()`, `PCShellSetView()`, `PCShellSetDestroy()`, `PCShellSetPostSolve()`,
894:           `PCShellSetApplyTranspose()`, `PCShellSetName()`, `PCShellSetApplyRichardson()`, `PCShellSetPreSolve()`, `PCShellSetView()`,
895:           `PCShellGetName()`, `PCShellSetContext()`, `PCShellGetContext()`, `PCShellSetApplyBA()`, `MATSHELL`, `PCShellSetMatApply()`,
896: M*/

898: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
899: {
900:   PC_Shell *shell;

902:   PetscFunctionBegin;
903:   PetscCall(PetscNew(&shell));
904:   pc->data = (void *)shell;

906:   pc->ops->destroy             = PCDestroy_Shell;
907:   pc->ops->view                = PCView_Shell;
908:   pc->ops->apply               = PCApply_Shell;
909:   pc->ops->applysymmetricleft  = PCApplySymmetricLeft_Shell;
910:   pc->ops->applysymmetricright = PCApplySymmetricRight_Shell;
911:   pc->ops->matapply            = NULL;
912:   pc->ops->applytranspose      = NULL;
913:   pc->ops->applyrichardson     = NULL;
914:   pc->ops->setup               = NULL;
915:   pc->ops->presolve            = NULL;
916:   pc->ops->postsolve           = NULL;

918:   shell->apply               = NULL;
919:   shell->applytranspose      = NULL;
920:   shell->name                = NULL;
921:   shell->applyrich           = NULL;
922:   shell->presolve            = NULL;
923:   shell->postsolve           = NULL;
924:   shell->ctx                 = NULL;
925:   shell->setup               = NULL;
926:   shell->view                = NULL;
927:   shell->destroy             = NULL;
928:   shell->applysymmetricleft  = NULL;
929:   shell->applysymmetricright = NULL;

931:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", PCShellSetDestroy_Shell));
932:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", PCShellSetSetUp_Shell));
933:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", PCShellSetApply_Shell));
934:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", PCShellSetMatApply_Shell));
935:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", PCShellSetApplySymmetricLeft_Shell));
936:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", PCShellSetApplySymmetricRight_Shell));
937:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", PCShellSetApplyBA_Shell));
938:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", PCShellSetPreSolve_Shell));
939:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", PCShellSetPostSolve_Shell));
940:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", PCShellSetView_Shell));
941:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", PCShellSetApplyTranspose_Shell));
942:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", PCShellSetName_Shell));
943:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", PCShellGetName_Shell));
944:   PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", PCShellSetApplyRichardson_Shell));
945:   PetscFunctionReturn(PETSC_SUCCESS);
946: }