Actual source code: shellpc.c

petsc-3.7.2 2016-06-05
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>        /*I "petscpc.h" I*/
  8: #include <petsc/private/vecimpl.h>

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

 13:   PetscErrorCode (*destroy)(PC);
 14:   PetscErrorCode (*setup)(PC);
 15:   PetscErrorCode (*apply)(PC,Vec,Vec);
 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;

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

 33:     Not Collective

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

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

 41:     Level: advanced

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

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

 49: .keywords: PC, shell, get, context

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

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

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

 72:    Logically Collective on PC

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

 78:    Level: advanced

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



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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

471: /* -------------------------------------------------------------------------------*/

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

479:    Logically Collective on PC

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

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

490: .  ptr - the application context

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

494:    Level: developer

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

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

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


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

517:    Logically Collective on PC

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

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

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

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

532:    Level: developer

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

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

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


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

554:    Logically Collective on PC

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

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

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

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

570:    Level: developer

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

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

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

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

591:    Logically Collective on PC

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

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

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

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

608:    Developer Notes: There should also be a PCShellSetApplySymmetricRight() and PCShellSetApplySymmetricLeft().

610:    Level: developer

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

614: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApplyBA()
615: @*/
616: PetscErrorCode  PCShellSetApply(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
617: {

622:   PetscTryMethod(pc,"PCShellSetApply_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
623:   return(0);
624: }

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

631:    Logically Collective on PC

633:    Input Parameters:
634: +  pc - the preconditioner context
635: -  apply - the application-provided left preconditioning routine

637:    Calling sequence of apply:
638: .vb
639:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
640: .ve

642: +  pc - the preconditioner, get the application context with PCShellGetContext()
643: .  xin - input vector
644: -  xout - output vector

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

648:    Level: developer

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

652: .seealso: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
653: @*/
654: PetscErrorCode  PCShellSetApplySymmetricLeft(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
655: {

660:   PetscTryMethod(pc,"PCShellSetApplySymmetricLeft_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
661:   return(0);
662: }

666: /*@C
667:    PCShellSetApply - Sets routine to use as right preconditioner (when the PC_SYMMETRIC is used).

669:    Logically Collective on PC

671:    Input Parameters:
672: +  pc - the preconditioner context
673: -  apply - the application-provided right preconditioning routine

675:    Calling sequence of apply:
676: .vb
677:    PetscErrorCode apply (PC pc,Vec xin,Vec xout)
678: .ve

680: +  pc - the preconditioner, get the application context with PCShellGetContext()
681: .  xin - input vector
682: -  xout - output vector

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

686:    Level: developer

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

690: .seealso: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
691: @*/
692: PetscErrorCode  PCShellSetApplySymmetricRight(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
693: {

698:   PetscTryMethod(pc,"PCShellSetApplySymmetricRight_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
699:   return(0);
700: }

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

707:    Logically Collective on PC

709:    Input Parameters:
710: +  pc - the preconditioner context
711: -  applyBA - the application-provided BA routine

713:    Calling sequence of apply:
714: .vb
715:    PetscErrorCode applyBA (PC pc,Vec xin,Vec xout)
716: .ve

718: +  pc - the preconditioner, get the application context with PCShellGetContext()
719: .  xin - input vector
720: -  xout - output vector

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

724:    Level: developer

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

728: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApply()
729: @*/
730: PetscErrorCode  PCShellSetApplyBA(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
731: {

736:   PetscTryMethod(pc,"PCShellSetApplyBA_C",(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec)),(pc,applyBA));
737:   return(0);
738: }

742: /*@C
743:    PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.

745:    Logically Collective on PC

747:    Input Parameters:
748: +  pc - the preconditioner context
749: -  apply - the application-provided preconditioning transpose routine

751:    Calling sequence of apply:
752: .vb
753:    PetscErrorCode applytranspose (PC pc,Vec xin,Vec xout)
754: .ve

756: +  pc - the preconditioner, get the application context with PCShellGetContext()
757: .  xin - input vector
758: -  xout - output vector

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

762:    Level: developer

764:    Notes:
765:    Uses the same context variable as PCShellSetApply().

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

769: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApply(), PCSetContext(), PCShellSetApplyBA()
770: @*/
771: PetscErrorCode  PCShellSetApplyTranspose(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
772: {

777:   PetscTryMethod(pc,"PCShellSetApplyTranspose_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,applytranspose));
778:   return(0);
779: }

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

788:    Logically Collective on PC

790:    Input Parameters:
791: +  pc - the preconditioner context
792: -  presolve - the application-provided presolve routine

794:    Calling sequence of presolve:
795: .vb
796:    PetscErrorCode presolve (PC,KSP ksp,Vec b,Vec x)
797: .ve

799: +  pc - the preconditioner, get the application context with PCShellGetContext()
800: .  xin - input vector
801: -  xout - output vector

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

805:    Level: developer

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

809: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPostSolve(), PCShellSetContext()
810: @*/
811: PetscErrorCode  PCShellSetPreSolve(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
812: {

817:   PetscTryMethod(pc,"PCShellSetPreSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,presolve));
818:   return(0);
819: }

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

828:    Logically Collective on PC

830:    Input Parameters:
831: +  pc - the preconditioner context
832: -  postsolve - the application-provided presolve routine

834:    Calling sequence of postsolve:
835: .vb
836:    PetscErrorCode postsolve(PC,KSP ksp,Vec b,Vec x)
837: .ve

839: +  pc - the preconditioner, get the application context with PCShellGetContext()
840: .  xin - input vector
841: -  xout - output vector

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

845:    Level: developer

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

849: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPreSolve(), PCShellSetContext()
850: @*/
851: PetscErrorCode  PCShellSetPostSolve(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
852: {

857:   PetscTryMethod(pc,"PCShellSetPostSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,postsolve));
858:   return(0);
859: }

863: /*@C
864:    PCShellSetName - Sets an optional name to associate with a shell
865:    preconditioner.

867:    Not Collective

869:    Input Parameters:
870: +  pc - the preconditioner context
871: -  name - character string describing shell preconditioner

873:    Level: developer

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

877: .seealso: PCShellGetName()
878: @*/
879: PetscErrorCode  PCShellSetName(PC pc,const char name[])
880: {

885:   PetscTryMethod(pc,"PCShellSetName_C",(PC,const char []),(pc,name));
886:   return(0);
887: }

891: /*@C
892:    PCShellGetName - Gets an optional name that the user has set for a shell
893:    preconditioner.

895:    Not Collective

897:    Input Parameter:
898: .  pc - the preconditioner context

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

903:    Level: developer

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

907: .seealso: PCShellSetName()
908: @*/
909: PetscErrorCode  PCShellGetName(PC pc,const char *name[])
910: {

916:   PetscUseMethod(pc,"PCShellGetName_C",(PC,const char*[]),(pc,name));
917:   return(0);
918: }

922: /*@C
923:    PCShellSetApplyRichardson - Sets routine to use as preconditioner
924:    in Richardson iteration.

926:    Logically Collective on PC

928:    Input Parameters:
929: +  pc - the preconditioner context
930: -  apply - the application-provided preconditioning routine

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

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

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

948:    Level: developer

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

952: .seealso: PCShellSetApply(), PCShellSetContext()
953: @*/
954: PetscErrorCode  PCShellSetApplyRichardson(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*))
955: {

960:   PetscTryMethod(pc,"PCShellSetApplyRichardson_C",(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*)),(pc,apply));
961:   return(0);
962: }

964: /*MC
965:    PCSHELL - Creates a new preconditioner class for use with your
966:               own private data storage format.

968:    Level: advanced
969: >
970:    Concepts: providing your own preconditioner

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

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

996: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
997: {
999:   PC_Shell       *shell;

1002:   PetscNewLog(pc,&shell);
1003:   pc->data = (void*)shell;

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

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

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