Actual source code: shellpc.c

petsc-master 2016-12-06
Report Typos and Errors

  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 (*applysymmetricleft)(PC,Vec,Vec);
 16:   PetscErrorCode (*applysymmetricright)(PC,Vec,Vec);
 17:   PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec);
 18:   PetscErrorCode (*presolve)(PC,KSP,Vec,Vec);
 19:   PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec);
 20:   PetscErrorCode (*view)(PC,PetscViewer);
 21:   PetscErrorCode (*applytranspose)(PC,Vec,Vec);
 22:   PetscErrorCode (*applyrich)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*);

 24:   char *name;
 25: } PC_Shell;

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

 32:     Not Collective

 34:     Input Parameter:
 35: .   pc - should have been created with PCSetType(pc,shell)

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

 40:     Level: advanced

 42:     Notes:
 43:     This routine is intended for use within various shell routines

 45:    Fortran Notes: 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: .keywords: PC, shell, get, context

 50: .seealso: PCShellSetContext()
 51: @*/
 52: PetscErrorCode  PCShellGetContext(PC pc,void **ctx)
 53: {
 55:   PetscBool      flg;

 60:   PetscObjectTypeCompare((PetscObject)pc,PCSHELL,&flg);
 61:   if (!flg) *ctx = 0;
 62:   else      *ctx = ((PC_Shell*)(pc->data))->ctx;
 63:   return(0);
 64: }

 68: /*@
 69:     PCShellSetContext - sets the context for a shell PC

 71:    Logically Collective on PC

 73:     Input Parameters:
 74: +   pc - the shell PC
 75: -   ctx - the context

 77:    Level: advanced

 79:    Fortran Notes: To use this from Fortran you must write a Fortran interface definition for this
 80:     function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.



 84: .seealso: PCShellGetContext(), PCSHELL
 85: @*/
 86: PetscErrorCode  PCShellSetContext(PC pc,void *ctx)
 87: {
 88:   PC_Shell       *shell = (PC_Shell*)pc->data;
 90:   PetscBool      flg;

 94:   PetscObjectTypeCompare((PetscObject)pc,PCSHELL,&flg);
 95:   if (flg) shell->ctx = ctx;
 96:   return(0);
 97: }

101: static PetscErrorCode PCSetUp_Shell(PC pc)
102: {
103:   PC_Shell       *shell = (PC_Shell*)pc->data;

107:   if (!shell->setup) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No setup() routine provided to Shell PC");
108:   PetscStackCall("PCSHELL user function setup()",(*shell->setup)(pc);CHKERRQ(ierr));
109:   return(0);
110: }

114: static PetscErrorCode PCApply_Shell(PC pc,Vec x,Vec y)
115: {
116:   PC_Shell         *shell = (PC_Shell*)pc->data;
117:   PetscErrorCode   ierr;
118:   PetscObjectState instate,outstate;

121:   if (!shell->apply) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No apply() routine provided to Shell PC");
122:   PetscObjectStateGet((PetscObject)y, &instate);
123:   PetscStackCall("PCSHELL user function apply()",(*shell->apply)(pc,x,y);CHKERRQ(ierr));
124:   PetscObjectStateGet((PetscObject)y, &outstate);
125:   if (instate == outstate) {
126:     /* increase the state of the output vector since the user did not update its state themselve as should have been done */
127:     PetscObjectStateIncrease((PetscObject)y);
128:   }
129:   return(0);
130: }

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

140:   if (!shell->applysymmetricleft) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No apply() routine provided to Shell PC");
141:   PetscStackCall("PCSHELL user function apply()",(*shell->applysymmetricleft)(pc,x,y);CHKERRQ(ierr));
142:   return(0);
143: }

147: static PetscErrorCode PCApplySymmetricRight_Shell(PC pc,Vec x,Vec y)
148: {
149:   PC_Shell       *shell = (PC_Shell*)pc->data;

153:   if (!shell->applysymmetricright) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No apply() routine provided to Shell PC");
154:   PetscStackCall("PCSHELL user function apply()",(*shell->applysymmetricright)(pc,x,y);CHKERRQ(ierr));
155:   return(0);
156: }

160: static PetscErrorCode PCApplyBA_Shell(PC pc,PCSide side,Vec x,Vec y,Vec w)
161: {
162:   PC_Shell         *shell = (PC_Shell*)pc->data;
163:   PetscErrorCode   ierr;
164:   PetscObjectState instate,outstate;

167:   if (!shell->applyBA) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applyBA() routine provided to Shell PC");
168:   PetscObjectStateGet((PetscObject)w, &instate);
169:   PetscStackCall("PCSHELL user function applyBA()",(*shell->applyBA)(pc,side,x,y,w);CHKERRQ(ierr));
170:   PetscObjectStateGet((PetscObject)w, &outstate);
171:   if (instate == outstate) {
172:     /* increase the state of the output vector since the user did not update its state themselve as should have been done */
173:     PetscObjectStateIncrease((PetscObject)w);
174:   }
175:   return(0);
176: }

180: static PetscErrorCode PCPreSolveChangeRHS_Shell(PC pc,PetscBool* change)
181: {
183:   *change = PETSC_TRUE;
184:   return(0);
185: }

189: static PetscErrorCode PCPreSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
190: {
191:   PC_Shell       *shell = (PC_Shell*)pc->data;

195:   if (!shell->presolve) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No presolve() routine provided to Shell PC");
196:   PetscStackCall("PCSHELL user function presolve()",(*shell->presolve)(pc,ksp,b,x);CHKERRQ(ierr));
197:   return(0);
198: }

202: static PetscErrorCode PCPostSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
203: {
204:   PC_Shell       *shell = (PC_Shell*)pc->data;

208:   if (!shell->postsolve) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No postsolve() routine provided to Shell PC");
209:   PetscStackCall("PCSHELL user function postsolve()",(*shell->postsolve)(pc,ksp,b,x);CHKERRQ(ierr));
210:   return(0);
211: }

215: static PetscErrorCode PCApplyTranspose_Shell(PC pc,Vec x,Vec y)
216: {
217:   PC_Shell         *shell = (PC_Shell*)pc->data;
218:   PetscErrorCode   ierr;
219:   PetscObjectState instate,outstate;

222:   if (!shell->applytranspose) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applytranspose() routine provided to Shell PC");
223:   PetscObjectStateGet((PetscObject)y, &instate);
224:   PetscStackCall("PCSHELL user function applytranspose()",(*shell->applytranspose)(pc,x,y);CHKERRQ(ierr));
225:   PetscObjectStateGet((PetscObject)y, &outstate);
226:   if (instate == outstate) {
227:     /* increase the state of the output vector since the user did not update its state themself as should have been done */
228:     PetscObjectStateIncrease((PetscObject)y);
229:   }
230:   return(0);
231: }

235: 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)
236: {
237:   PetscErrorCode   ierr;
238:   PC_Shell         *shell = (PC_Shell*)pc->data;
239:   PetscObjectState instate,outstate;

242:   if (!shell->applyrich) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applyrichardson() routine provided to Shell PC");
243:   PetscObjectStateGet((PetscObject)y, &instate);
244:   PetscStackCall("PCSHELL user function applyrichardson()",(*shell->applyrich)(pc,x,y,w,rtol,abstol,dtol,it,guesszero,outits,reason);CHKERRQ(ierr));
245:   PetscObjectStateGet((PetscObject)y, &outstate);
246:   if (instate == outstate) {
247:     /* increase the state of the output vector since the user did not update its state themself as should have been done */
248:     PetscObjectStateIncrease((PetscObject)y);
249:   }
250:   return(0);
251: }

255: static PetscErrorCode PCDestroy_Shell(PC pc)
256: {
257:   PC_Shell       *shell = (PC_Shell*)pc->data;

261:   PetscFree(shell->name);
262:   if (shell->destroy) PetscStackCall("PCSHELL user function destroy()",(*shell->destroy)(pc);CHKERRQ(ierr));
263:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetDestroy_C",NULL);
264:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetSetUp_C",NULL);
265:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApply_C",NULL);
266:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricLeft_C",NULL);
267:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricRight_C",NULL);
268:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyBA_C",NULL);
269:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPreSolve_C",NULL);
270:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPostSolve_C",NULL);
271:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetView_C",NULL);
272:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",NULL);
273:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetName_C",NULL);
274:   PetscObjectComposeFunction((PetscObject)pc,"PCShellGetName_C",NULL);
275:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",NULL);
276:   PetscObjectComposeFunction((PetscObject)pc,"PCPreSolveChangeRHS_C",NULL);
277:   PetscFree(pc->data);
278:   return(0);
279: }

283: static PetscErrorCode PCView_Shell(PC pc,PetscViewer viewer)
284: {
285:   PC_Shell       *shell = (PC_Shell*)pc->data;
287:   PetscBool      iascii;

290:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
291:   if (iascii) {
292:     if (shell->name) {
293:       PetscViewerASCIIPrintf(viewer,"  Shell: %s\n",shell->name);
294:     } else {
295:       PetscViewerASCIIPrintf(viewer,"  Shell: no name\n");
296:     }
297:   }
298:   if (shell->view) {
299:     PetscViewerASCIIPushTab(viewer);
300:     (*shell->view)(pc,viewer);
301:     PetscViewerASCIIPopTab(viewer);
302:   }
303:   return(0);
304: }

306: /* ------------------------------------------------------------------------------*/
309: static PetscErrorCode  PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
310: {
311:   PC_Shell *shell= (PC_Shell*)pc->data;

314:   shell->destroy = destroy;
315:   return(0);
316: }

320: static PetscErrorCode  PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
321: {
322:   PC_Shell *shell = (PC_Shell*)pc->data;;

325:   shell->setup = setup;
326:   if (setup) pc->ops->setup = PCSetUp_Shell;
327:   else       pc->ops->setup = 0;
328:   return(0);
329: }

333: static PetscErrorCode  PCShellSetApply_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
334: {
335:   PC_Shell *shell = (PC_Shell*)pc->data;

338:   shell->apply = apply;
339:   return(0);
340: }

344: static PetscErrorCode  PCShellSetApplySymmetricLeft_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
345: {
346:   PC_Shell *shell = (PC_Shell*)pc->data;

349:   shell->applysymmetricleft = apply;
350:   return(0);
351: }

355: static PetscErrorCode  PCShellSetApplySymmetricRight_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
356: {
357:   PC_Shell *shell = (PC_Shell*)pc->data;

360:   shell->applysymmetricright = apply;
361:   return(0);
362: }

366: static PetscErrorCode  PCShellSetApplyBA_Shell(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
367: {
368:   PC_Shell *shell = (PC_Shell*)pc->data;

371:   shell->applyBA = applyBA;
372:   if (applyBA) pc->ops->applyBA  = PCApplyBA_Shell;
373:   else         pc->ops->applyBA  = 0;
374:   return(0);
375: }

379: static PetscErrorCode  PCShellSetPreSolve_Shell(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
380: {
381:   PC_Shell       *shell = (PC_Shell*)pc->data;

385:   shell->presolve = presolve;
386:   if (presolve) {
387:     pc->ops->presolve = PCPreSolve_Shell;
388:     PetscObjectComposeFunction((PetscObject)pc,"PCPreSolveChangeRHS_C",PCPreSolveChangeRHS_Shell);
389:   } else {
390:     pc->ops->presolve = 0;
391:     PetscObjectComposeFunction((PetscObject)pc,"PCPreSolveChangeRHS_C",NULL);
392:   }
393:   return(0);
394: }

398: static PetscErrorCode  PCShellSetPostSolve_Shell(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
399: {
400:   PC_Shell *shell = (PC_Shell*)pc->data;

403:   shell->postsolve = postsolve;
404:   if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
405:   else           pc->ops->postsolve = 0;
406:   return(0);
407: }

411: static PetscErrorCode  PCShellSetView_Shell(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
412: {
413:   PC_Shell *shell = (PC_Shell*)pc->data;

416:   shell->view = view;
417:   return(0);
418: }

422: static PetscErrorCode  PCShellSetApplyTranspose_Shell(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
423: {
424:   PC_Shell *shell = (PC_Shell*)pc->data;

427:   shell->applytranspose = applytranspose;
428:   if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
429:   else                pc->ops->applytranspose = 0;
430:   return(0);
431: }

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

440:   shell->applyrich = applyrich;
441:   if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
442:   else           pc->ops->applyrichardson = 0;
443:   return(0);
444: }

448: static PetscErrorCode  PCShellSetName_Shell(PC pc,const char name[])
449: {
450:   PC_Shell       *shell = (PC_Shell*)pc->data;

454:   PetscFree(shell->name);
455:   PetscStrallocpy(name,&shell->name);
456:   return(0);
457: }

461: static PetscErrorCode  PCShellGetName_Shell(PC pc,const char *name[])
462: {
463:   PC_Shell *shell = (PC_Shell*)pc->data;

466:   *name = shell->name;
467:   return(0);
468: }

470: /* -------------------------------------------------------------------------------*/

474: /*@C
475:    PCShellSetDestroy - Sets routine to use to destroy the user-provided
476:    application context.

478:    Logically Collective on PC

480:    Input Parameters:
481: +  pc - the preconditioner context
482: .  destroy - the application-provided destroy routine

484:    Calling sequence of destroy:
485: .vb
486:    PetscErrorCode destroy (PC)
487: .ve

489: .  ptr - the application context

491:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

493:    Level: developer

495: .keywords: PC, shell, set, destroy, user-provided

497: .seealso: PCShellSetApply(), PCShellSetContext()
498: @*/
499: PetscErrorCode  PCShellSetDestroy(PC pc,PetscErrorCode (*destroy)(PC))
500: {

505:   PetscTryMethod(pc,"PCShellSetDestroy_C",(PC,PetscErrorCode (*)(PC)),(pc,destroy));
506:   return(0);
507: }


512: /*@C
513:    PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
514:    matrix operator is changed.

516:    Logically Collective on PC

518:    Input Parameters:
519: +  pc - the preconditioner context
520: .  setup - the application-provided setup routine

522:    Calling sequence of setup:
523: .vb
524:    PetscErrorCode setup (PC pc)
525: .ve

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

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

531:    Level: developer

533: .keywords: PC, shell, set, setup, user-provided

535: .seealso: PCShellSetApplyRichardson(), PCShellSetApply(), PCShellSetContext()
536: @*/
537: PetscErrorCode  PCShellSetSetUp(PC pc,PetscErrorCode (*setup)(PC))
538: {

543:   PetscTryMethod(pc,"PCShellSetSetUp_C",(PC,PetscErrorCode (*)(PC)),(pc,setup));
544:   return(0);
545: }


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

553:    Logically Collective on PC

555:    Input Parameters:
556: +  pc - the preconditioner context
557: -  view - the application-provided view routine

559:    Calling sequence of apply:
560: .vb
561:    PetscErrorCode view(PC pc,PetscViewer v)
562: .ve

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

567:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

569:    Level: developer

571: .keywords: PC, shell, set, apply, user-provided

573: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose()
574: @*/
575: PetscErrorCode  PCShellSetView(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
576: {

581:   PetscTryMethod(pc,"PCShellSetView_C",(PC,PetscErrorCode (*)(PC,PetscViewer)),(pc,view));
582:   return(0);
583: }

587: /*@C
588:    PCShellSetApply - Sets routine to use as preconditioner.

590:    Logically Collective on PC

592:    Input Parameters:
593: +  pc - the preconditioner context
594: -  apply - the application-provided preconditioning routine

596:    Calling sequence of apply:
597: .vb
598:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
599: .ve

601: +  pc - the preconditioner, get the application context with PCShellGetContext()
602: .  xin - input vector
603: -  xout - output vector

605:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

607:    Level: developer

609: .keywords: PC, shell, set, apply, user-provided

611: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApplyBA(), PCShellSetApplySymmetricRight(),PCShellSetApplySymmetricLeft()
612: @*/
613: PetscErrorCode  PCShellSetApply(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
614: {

619:   PetscTryMethod(pc,"PCShellSetApply_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
620:   return(0);
621: }

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

628:    Logically Collective on PC

630:    Input Parameters:
631: +  pc - the preconditioner context
632: -  apply - the application-provided left preconditioning routine

634:    Calling sequence of apply:
635: .vb
636:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
637: .ve

639: +  pc - the preconditioner, get the application context with PCShellGetContext()
640: .  xin - input vector
641: -  xout - output vector

643:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

645:    Level: developer

647: .keywords: PC, shell, set, apply, user-provided

649: .seealso: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
650: @*/
651: PetscErrorCode  PCShellSetApplySymmetricLeft(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
652: {

657:   PetscTryMethod(pc,"PCShellSetApplySymmetricLeft_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
658:   return(0);
659: }

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

666:    Logically Collective on PC

668:    Input Parameters:
669: +  pc - the preconditioner context
670: -  apply - the application-provided right preconditioning routine

672:    Calling sequence of apply:
673: .vb
674:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
675: .ve

677: +  pc - the preconditioner, get the application context with PCShellGetContext()
678: .  xin - input vector
679: -  xout - output vector

681:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

683:    Level: developer

685: .keywords: PC, shell, set, apply, user-provided

687: .seealso: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
688: @*/
689: PetscErrorCode  PCShellSetApplySymmetricRight(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
690: {

695:   PetscTryMethod(pc,"PCShellSetApplySymmetricRight_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
696:   return(0);
697: }

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

704:    Logically Collective on PC

706:    Input Parameters:
707: +  pc - the preconditioner context
708: -  applyBA - the application-provided BA routine

710:    Calling sequence of apply:
711: .vb
712:    PetscErrorCode applyBA (PC pc,Vec xin,Vec xout)
713: .ve

715: +  pc - the preconditioner, get the application context with PCShellGetContext()
716: .  xin - input vector
717: -  xout - output vector

719:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

721:    Level: developer

723: .keywords: PC, shell, set, apply, user-provided

725: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApply()
726: @*/
727: PetscErrorCode  PCShellSetApplyBA(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
728: {

733:   PetscTryMethod(pc,"PCShellSetApplyBA_C",(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec)),(pc,applyBA));
734:   return(0);
735: }

739: /*@C
740:    PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.

742:    Logically Collective on PC

744:    Input Parameters:
745: +  pc - the preconditioner context
746: -  apply - the application-provided preconditioning transpose routine

748:    Calling sequence of apply:
749: .vb
750:    PetscErrorCode applytranspose (PC pc,Vec xin,Vec xout)
751: .ve

753: +  pc - the preconditioner, get the application context with PCShellGetContext()
754: .  xin - input vector
755: -  xout - output vector

757:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

759:    Level: developer

761:    Notes:
762:    Uses the same context variable as PCShellSetApply().

764: .keywords: PC, shell, set, apply, user-provided

766: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApply(), PCSetContext(), PCShellSetApplyBA()
767: @*/
768: PetscErrorCode  PCShellSetApplyTranspose(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
769: {

774:   PetscTryMethod(pc,"PCShellSetApplyTranspose_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,applytranspose));
775:   return(0);
776: }

780: /*@C
781:    PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
782:       applied. This usually does something like scale the linear system in some application
783:       specific way.

785:    Logically Collective on PC

787:    Input Parameters:
788: +  pc - the preconditioner context
789: -  presolve - the application-provided presolve routine

791:    Calling sequence of presolve:
792: .vb
793:    PetscErrorCode presolve (PC,KSP ksp,Vec b,Vec x)
794: .ve

796: +  pc - the preconditioner, get the application context with PCShellGetContext()
797: .  xin - input vector
798: -  xout - output vector

800:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

802:    Level: developer

804: .keywords: PC, shell, set, apply, user-provided

806: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPostSolve(), PCShellSetContext()
807: @*/
808: PetscErrorCode  PCShellSetPreSolve(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
809: {

814:   PetscTryMethod(pc,"PCShellSetPreSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,presolve));
815:   return(0);
816: }

820: /*@C
821:    PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
822:       applied. This usually does something like scale the linear system in some application
823:       specific way.

825:    Logically Collective on PC

827:    Input Parameters:
828: +  pc - the preconditioner context
829: -  postsolve - the application-provided presolve routine

831:    Calling sequence of postsolve:
832: .vb
833:    PetscErrorCode postsolve(PC,KSP ksp,Vec b,Vec x)
834: .ve

836: +  pc - the preconditioner, get the application context with PCShellGetContext()
837: .  xin - input vector
838: -  xout - output vector

840:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

842:    Level: developer

844: .keywords: PC, shell, set, apply, user-provided

846: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPreSolve(), PCShellSetContext()
847: @*/
848: PetscErrorCode  PCShellSetPostSolve(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
849: {

854:   PetscTryMethod(pc,"PCShellSetPostSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,postsolve));
855:   return(0);
856: }

860: /*@C
861:    PCShellSetName - Sets an optional name to associate with a shell
862:    preconditioner.

864:    Not Collective

866:    Input Parameters:
867: +  pc - the preconditioner context
868: -  name - character string describing shell preconditioner

870:    Level: developer

872: .keywords: PC, shell, set, name, user-provided

874: .seealso: PCShellGetName()
875: @*/
876: PetscErrorCode  PCShellSetName(PC pc,const char name[])
877: {

882:   PetscTryMethod(pc,"PCShellSetName_C",(PC,const char []),(pc,name));
883:   return(0);
884: }

888: /*@C
889:    PCShellGetName - Gets an optional name that the user has set for a shell
890:    preconditioner.

892:    Not Collective

894:    Input Parameter:
895: .  pc - the preconditioner context

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

900:    Level: developer

902: .keywords: PC, shell, get, name, user-provided

904: .seealso: PCShellSetName()
905: @*/
906: PetscErrorCode  PCShellGetName(PC pc,const char *name[])
907: {

913:   PetscUseMethod(pc,"PCShellGetName_C",(PC,const char*[]),(pc,name));
914:   return(0);
915: }

919: /*@C
920:    PCShellSetApplyRichardson - Sets routine to use as preconditioner
921:    in Richardson iteration.

923:    Logically Collective on PC

925:    Input Parameters:
926: +  pc - the preconditioner context
927: -  apply - the application-provided preconditioning routine

929:    Calling sequence of apply:
930: .vb
931:    PetscErrorCode apply (PC pc,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
932: .ve

934: +  pc - the preconditioner, get the application context with PCShellGetContext()
935: .  b - right-hand-side
936: .  x - current iterate
937: .  r - work space
938: .  rtol - relative tolerance of residual norm to stop at
939: .  abstol - absolute tolerance of residual norm to stop at
940: .  dtol - if residual norm increases by this factor than return
941: -  maxits - number of iterations to run

943:    Notes: the function MUST return an error code of 0 on success and nonzero on failure.

945:    Level: developer

947: .keywords: PC, shell, set, apply, Richardson, user-provided

949: .seealso: PCShellSetApply(), PCShellSetContext()
950: @*/
951: PetscErrorCode  PCShellSetApplyRichardson(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*))
952: {

957:   PetscTryMethod(pc,"PCShellSetApplyRichardson_C",(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*)),(pc,apply));
958:   return(0);
959: }

961: /*MC
962:    PCSHELL - Creates a new preconditioner class for use with your
963:               own private data storage format.

965:    Level: advanced
966: >
967:    Concepts: providing your own preconditioner

969:   Usage:
970: $             extern PetscErrorCode apply(PC,Vec,Vec);
971: $             extern PetscErrorCode applyba(PC,PCSide,Vec,Vec,Vec);
972: $             extern PetscErrorCode applytranspose(PC,Vec,Vec);
973: $             extern PetscErrorCode setup(PC);
974: $             extern PetscErrorCode destroy(PC);
975: $
976: $             PCCreate(comm,&pc);
977: $             PCSetType(pc,PCSHELL);
978: $             PCShellSetContext(pc,ctx)
979: $             PCShellSetApply(pc,apply);
980: $             PCShellSetApplyBA(pc,applyba);               (optional)
981: $             PCShellSetApplyTranspose(pc,applytranspose); (optional)
982: $             PCShellSetSetUp(pc,setup);                   (optional)
983: $             PCShellSetDestroy(pc,destroy);               (optional)

985: .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC,
986:            MATSHELL, PCShellSetSetUp(), PCShellSetApply(), PCShellSetView(),
987:            PCShellSetApplyTranspose(), PCShellSetName(), PCShellSetApplyRichardson(),
988:            PCShellGetName(), PCShellSetContext(), PCShellGetContext(), PCShellSetApplyBA()
989: M*/

993: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
994: {
996:   PC_Shell       *shell;

999:   PetscNewLog(pc,&shell);
1000:   pc->data = (void*)shell;

1002:   pc->ops->destroy         = PCDestroy_Shell;
1003:   pc->ops->view            = PCView_Shell;
1004:   pc->ops->apply           = PCApply_Shell;
1005:   pc->ops->applysymmetricleft  = PCApplySymmetricLeft_Shell;
1006:   pc->ops->applysymmetricright = PCApplySymmetricRight_Shell;
1007:   pc->ops->applytranspose  = 0;
1008:   pc->ops->applyrichardson = 0;
1009:   pc->ops->setup           = 0;
1010:   pc->ops->presolve        = 0;
1011:   pc->ops->postsolve       = 0;

1013:   shell->apply          = 0;
1014:   shell->applytranspose = 0;
1015:   shell->name           = 0;
1016:   shell->applyrich      = 0;
1017:   shell->presolve       = 0;
1018:   shell->postsolve      = 0;
1019:   shell->ctx            = 0;
1020:   shell->setup          = 0;
1021:   shell->view           = 0;
1022:   shell->destroy        = 0;
1023:   shell->applysymmetricleft  = 0;
1024:   shell->applysymmetricright = 0;

1026:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetDestroy_C",PCShellSetDestroy_Shell);
1027:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetSetUp_C",PCShellSetSetUp_Shell);
1028:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApply_C",PCShellSetApply_Shell);
1029:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricLeft_C",PCShellSetApplySymmetricLeft_Shell);
1030:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricRight_C",PCShellSetApplySymmetricRight_Shell);
1031:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyBA_C",PCShellSetApplyBA_Shell);
1032:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPreSolve_C",PCShellSetPreSolve_Shell);
1033:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPostSolve_C",PCShellSetPostSolve_Shell);
1034:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetView_C",PCShellSetView_Shell);
1035:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",PCShellSetApplyTranspose_Shell);
1036:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetName_C",PCShellSetName_Shell);
1037:   PetscObjectComposeFunction((PetscObject)pc,"PCShellGetName_C",PCShellGetName_Shell);
1038:   PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",PCShellSetApplyRichardson_Shell);
1039:   return(0);
1040: }