Actual source code: shellpc.c
1: #define PETSCKSP_DLL
3: /*
4: This provides a simple shell for Fortran (and C programmers) to
5: create their own preconditioner without writing much interface code.
6: */
8: #include private/pcimpl.h
9: #include private/vecimpl.h
12: typedef struct {
13: void *ctx; /* user provided contexts for preconditioner */
14: PetscErrorCode (*destroy)(void*);
15: PetscErrorCode (*setup)(void*);
16: PetscErrorCode (*apply)(void*,Vec,Vec);
17: PetscErrorCode (*presolve)(void*,KSP,Vec,Vec);
18: PetscErrorCode (*postsolve)(void*,KSP,Vec,Vec);
19: PetscErrorCode (*view)(void*,PetscViewer);
20: PetscErrorCode (*applytranspose)(void*,Vec,Vec);
21: PetscErrorCode (*applyrich)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt);
22: char *name;
23: } PC_Shell;
28: /*@
29: PCShellGetContext - Returns the user-provided context associated with a shell PC
31: Not Collective
33: Input Parameter:
34: . pc - should have been created with PCCreateShell()
36: Output Parameter:
37: . ctx - the user provided context
39: Level: advanced
41: Notes:
42: This routine is intended for use within various shell routines
43:
44: .keywords: PC, shell, get, context
46: .seealso: PCCreateShell(), PCShellSetContext()
47: @*/
48: PetscErrorCode PETSCKSP_DLLEXPORT PCShellGetContext(PC pc,void **ctx)
49: {
51: PetscTruth flg;
56: PetscTypeCompare((PetscObject)pc,PCSHELL,&flg);
57: if (!flg) *ctx = 0;
58: else *ctx = ((PC_Shell*)(pc->data))->ctx;
59: return(0);
60: }
64: /*@C
65: PCShellSetContext - sets the context for a shell PC
67: Collective on PC
69: Input Parameters:
70: + pc - the shell PC
71: - ctx - the context
73: Level: advanced
75: Fortran Notes: The context can only be an integer or a PetscObject
76: unfortunately it cannot be a Fortran array or derived type.
78: .seealso: PCCreateShell(), PCShellGetContext()
79: @*/
80: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetContext(PC pc,void *ctx)
81: {
82: PC_Shell *shell = (PC_Shell*)pc->data;
84: PetscTruth flg;
88: PetscTypeCompare((PetscObject)pc,PCSHELL,&flg);
89: if (flg) {
90: shell->ctx = ctx;
91: }
92: return(0);
93: }
97: static PetscErrorCode PCSetUp_Shell(PC pc)
98: {
99: PC_Shell *shell;
103: shell = (PC_Shell*)pc->data;
104: if (shell->setup) {
105: CHKMEMQ;
106: (*shell->setup)(shell->ctx);
107: CHKMEMQ;
108: }
109: return(0);
110: }
114: static PetscErrorCode PCApply_Shell(PC pc,Vec x,Vec y)
115: {
116: PC_Shell *shell;
120: shell = (PC_Shell*)pc->data;
121: if (!shell->apply) SETERRQ(PETSC_ERR_USER,"No apply() routine provided to Shell PC");
122: PetscStackPush("PCSHELL user function");
123: CHKMEMQ;
124: (*shell->apply)(shell->ctx,x,y);
125: CHKMEMQ;
126: PetscStackPop;
127: return(0);
128: }
132: static PetscErrorCode PCPreSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
133: {
134: PC_Shell *shell;
138: shell = (PC_Shell*)pc->data;
139: if (!shell->presolve) SETERRQ(PETSC_ERR_USER,"No presolve() routine provided to Shell PC");
140: (*shell->presolve)(shell->ctx,ksp,b,x);
141: return(0);
142: }
146: static PetscErrorCode PCPostSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
147: {
148: PC_Shell *shell;
152: shell = (PC_Shell*)pc->data;
153: if (!shell->postsolve) SETERRQ(PETSC_ERR_USER,"No postsolve() routine provided to Shell PC");
154: (*shell->postsolve)(shell->ctx,ksp,b,x);
155: return(0);
156: }
160: static PetscErrorCode PCApplyTranspose_Shell(PC pc,Vec x,Vec y)
161: {
162: PC_Shell *shell;
166: shell = (PC_Shell*)pc->data;
167: if (!shell->applytranspose) SETERRQ(PETSC_ERR_USER,"No applytranspose() routine provided to Shell PC");
168: (*shell->applytranspose)(shell->ctx,x,y);
169: return(0);
170: }
174: static PetscErrorCode PCApplyRichardson_Shell(PC pc,Vec x,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt it)
175: {
177: PC_Shell *shell;
180: shell = (PC_Shell*)pc->data;
181: (*shell->applyrich)(shell->ctx,x,y,w,rtol,abstol,dtol,it);
182: return(0);
183: }
187: static PetscErrorCode PCDestroy_Shell(PC pc)
188: {
189: PC_Shell *shell = (PC_Shell*)pc->data;
193: PetscStrfree(shell->name);
194: if (shell->destroy) {
195: (*shell->destroy)(shell->ctx);
196: }
197: PetscFree(shell);
198: return(0);
199: }
203: static PetscErrorCode PCView_Shell(PC pc,PetscViewer viewer)
204: {
205: PC_Shell *shell = (PC_Shell*)pc->data;
207: PetscTruth iascii;
210: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
211: if (iascii) {
212: if (shell->name) {PetscViewerASCIIPrintf(viewer," Shell: %s\n",shell->name);}
213: else {PetscViewerASCIIPrintf(viewer," Shell: no name\n");}
214: }
215: if (shell->view) {
216: PetscViewerASCIIPushTab(viewer);
217: (*shell->view)(shell->ctx,viewer);
218: PetscViewerASCIIPopTab(viewer);
219: }
220: return(0);
221: }
223: /* ------------------------------------------------------------------------------*/
227: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(void*))
228: {
229: PC_Shell *shell;
232: shell = (PC_Shell*)pc->data;
233: shell->destroy = destroy;
234: return(0);
235: }
241: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(void*))
242: {
243: PC_Shell *shell;
246: shell = (PC_Shell*)pc->data;
247: shell->setup = setup;
248: return(0);
249: }
255: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApply_Shell(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec))
256: {
257: PC_Shell *shell;
260: shell = (PC_Shell*)pc->data;
261: shell->apply = apply;
262: return(0);
263: }
269: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPreSolve_Shell(PC pc,PetscErrorCode (*presolve)(void*,KSP,Vec,Vec))
270: {
271: PC_Shell *shell;
274: shell = (PC_Shell*)pc->data;
275: shell->presolve = presolve;
276: if (presolve) {
277: pc->ops->presolve = PCPreSolve_Shell;
278: } else {
279: pc->ops->presolve = 0;
280: }
281: return(0);
282: }
288: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPostSolve_Shell(PC pc,PetscErrorCode (*postsolve)(void*,KSP,Vec,Vec))
289: {
290: PC_Shell *shell;
293: shell = (PC_Shell*)pc->data;
294: shell->postsolve = postsolve;
295: if (postsolve) {
296: pc->ops->postsolve = PCPostSolve_Shell;
297: } else {
298: pc->ops->postsolve = 0;
299: }
300: return(0);
301: }
307: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetView_Shell(PC pc,PetscErrorCode (*view)(void*,PetscViewer))
308: {
309: PC_Shell *shell;
312: shell = (PC_Shell*)pc->data;
313: shell->view = view;
314: return(0);
315: }
321: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyTranspose_Shell(PC pc,PetscErrorCode (*applytranspose)(void*,Vec,Vec))
322: {
323: PC_Shell *shell;
326: shell = (PC_Shell*)pc->data;
327: shell->applytranspose = applytranspose;
328: return(0);
329: }
335: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetName_Shell(PC pc,const char name[])
336: {
337: PC_Shell *shell;
341: shell = (PC_Shell*)pc->data;
342: PetscStrfree(shell->name);
343: PetscStrallocpy(name,&shell->name);
344: return(0);
345: }
351: PetscErrorCode PETSCKSP_DLLEXPORT PCShellGetName_Shell(PC pc,char *name[])
352: {
353: PC_Shell *shell;
356: shell = (PC_Shell*)pc->data;
357: *name = shell->name;
358: return(0);
359: }
365: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyRichardson_Shell(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt))
366: {
367: PC_Shell *shell;
370: shell = (PC_Shell*)pc->data;
371: pc->ops->applyrichardson = PCApplyRichardson_Shell;
372: shell->applyrich = apply;
373: return(0);
374: }
377: /* -------------------------------------------------------------------------------*/
381: /*@C
382: PCShellSetDestroy - Sets routine to use to destroy the user-provided
383: application context.
385: Collective on PC
387: Input Parameters:
388: + pc - the preconditioner context
389: . destroy - the application-provided destroy routine
391: Calling sequence of destroy:
392: .vb
393: PetscErrorCode destroy (void *ptr)
394: .ve
396: . ptr - the application context
398: Level: developer
400: .keywords: PC, shell, set, destroy, user-provided
402: .seealso: PCShellSetApply(), PCShellSetContext()
403: @*/
404: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetDestroy(PC pc,PetscErrorCode (*destroy)(void*))
405: {
406: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*));
410: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetDestroy_C",(void (**)(void))&f);
411: if (f) {
412: (*f)(pc,destroy);
413: }
414: return(0);
415: }
420: /*@C
421: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
422: matrix operator is changed.
424: Collective on PC
426: Input Parameters:
427: + pc - the preconditioner context
428: . setup - the application-provided setup routine
430: Calling sequence of setup:
431: .vb
432: PetscErrorCode setup (void *ptr)
433: .ve
435: . ptr - the application context
437: Level: developer
439: .keywords: PC, shell, set, setup, user-provided
441: .seealso: PCShellSetApplyRichardson(), PCShellSetApply(), PCShellSetContext()
442: @*/
443: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetSetUp(PC pc,PetscErrorCode (*setup)(void*))
444: {
445: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*));
449: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetSetUp_C",(void (**)(void))&f);
450: if (f) {
451: (*f)(pc,setup);
452: }
453: return(0);
454: }
459: /*@C
460: PCShellSetView - Sets routine to use as viewer of shell preconditioner
462: Collective on PC
464: Input Parameters:
465: + pc - the preconditioner context
466: - view - the application-provided view routine
468: Calling sequence of apply:
469: .vb
470: PetscErrorCode view(void *ptr,PetscViewer v)
471: .ve
473: + ptr - the application context
474: - v - viewer
476: Level: developer
478: .keywords: PC, shell, set, apply, user-provided
480: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose()
481: @*/
482: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetView(PC pc,PetscErrorCode (*view)(void*,PetscViewer))
483: {
484: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,PetscViewer));
488: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetView_C",(void (**)(void))&f);
489: if (f) {
490: (*f)(pc,view);
491: }
492: return(0);
493: }
497: /*@C
498: PCShellSetApply - Sets routine to use as preconditioner.
500: Collective on PC
502: Input Parameters:
503: + pc - the preconditioner context
504: - apply - the application-provided preconditioning routine
506: Calling sequence of apply:
507: .vb
508: PetscErrorCode apply (void *ptr,Vec xin,Vec xout)
509: .ve
511: + ptr - the application context
512: . xin - input vector
513: - xout - output vector
515: Level: developer
517: .keywords: PC, shell, set, apply, user-provided
519: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
520: @*/
521: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApply(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec))
522: {
523: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,Vec,Vec));
527: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApply_C",(void (**)(void))&f);
528: if (f) {
529: (*f)(pc,apply);
530: }
531: return(0);
532: }
536: /*@C
537: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
539: Collective on PC
541: Input Parameters:
542: + pc - the preconditioner context
543: - apply - the application-provided preconditioning transpose routine
545: Calling sequence of apply:
546: .vb
547: PetscErrorCode applytranspose (void *ptr,Vec xin,Vec xout)
548: .ve
550: + ptr - the application context
551: . xin - input vector
552: - xout - output vector
554: Level: developer
556: Notes:
557: Uses the same context variable as PCShellSetApply().
559: .keywords: PC, shell, set, apply, user-provided
561: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApply(), PCSetContext()
562: @*/
563: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyTranspose(PC pc,PetscErrorCode (*applytranspose)(void*,Vec,Vec))
564: {
565: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,Vec,Vec));
569: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",(void (**)(void))&f);
570: if (f) {
571: (*f)(pc,applytranspose);
572: }
573: return(0);
574: }
578: /*@C
579: PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
580: applied. This usually does something like scale the linear system in some application
581: specific way.
583: Collective on PC
585: Input Parameters:
586: + pc - the preconditioner context
587: - presolve - the application-provided presolve routine
589: Calling sequence of presolve:
590: .vb
591: PetscErrorCode presolve (void *ptr,KSP ksp,Vec b,Vec x)
592: .ve
594: + ptr - the application context
595: . xin - input vector
596: - xout - output vector
598: Level: developer
600: .keywords: PC, shell, set, apply, user-provided
602: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPostSolve(), PCShellSetContext()
603: @*/
604: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPreSolve(PC pc,PetscErrorCode (*presolve)(void*,KSP,Vec,Vec))
605: {
606: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,KSP,Vec,Vec));
610: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetPreSolve_C",(void (**)(void))&f);
611: if (f) {
612: (*f)(pc,presolve);
613: }
614: return(0);
615: }
619: /*@C
620: PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
621: applied. This usually does something like scale the linear system in some application
622: specific way.
624: Collective on PC
626: Input Parameters:
627: + pc - the preconditioner context
628: - postsolve - the application-provided presolve routine
630: Calling sequence of postsolve:
631: .vb
632: PetscErrorCode postsolve(void *ptr,KSP ksp,Vec b,Vec x)
633: .ve
635: + ptr - the application context
636: . xin - input vector
637: - xout - output vector
639: Level: developer
641: .keywords: PC, shell, set, apply, user-provided
643: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPreSolve(), PCShellSetContext()
644: @*/
645: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPostSolve(PC pc,PetscErrorCode (*postsolve)(void*,KSP,Vec,Vec))
646: {
647: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,KSP,Vec,Vec));
651: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetPostSolve_C",(void (**)(void))&f);
652: if (f) {
653: (*f)(pc,postsolve);
654: }
655: return(0);
656: }
660: /*@C
661: PCShellSetName - Sets an optional name to associate with a shell
662: preconditioner.
664: Not Collective
666: Input Parameters:
667: + pc - the preconditioner context
668: - name - character string describing shell preconditioner
670: Level: developer
672: .keywords: PC, shell, set, name, user-provided
674: .seealso: PCShellGetName()
675: @*/
676: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetName(PC pc,const char name[])
677: {
678: PetscErrorCode ierr,(*f)(PC,const char []);
682: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetName_C",(void (**)(void))&f);
683: if (f) {
684: (*f)(pc,name);
685: }
686: return(0);
687: }
691: /*@C
692: PCShellGetName - Gets an optional name that the user has set for a shell
693: preconditioner.
695: Not Collective
697: Input Parameter:
698: . pc - the preconditioner context
700: Output Parameter:
701: . name - character string describing shell preconditioner (you should not free this)
703: Level: developer
705: .keywords: PC, shell, get, name, user-provided
707: .seealso: PCShellSetName()
708: @*/
709: PetscErrorCode PETSCKSP_DLLEXPORT PCShellGetName(PC pc,char *name[])
710: {
711: PetscErrorCode ierr,(*f)(PC,char *[]);
716: PetscObjectQueryFunction((PetscObject)pc,"PCShellGetName_C",(void (**)(void))&f);
717: if (f) {
718: (*f)(pc,name);
719: } else {
720: SETERRQ(PETSC_ERR_ARG_WRONG,"Not shell preconditioner, cannot get name");
721: }
722: return(0);
723: }
727: /*@C
728: PCShellSetApplyRichardson - Sets routine to use as preconditioner
729: in Richardson iteration.
731: Collective on PC
733: Input Parameters:
734: + pc - the preconditioner context
735: - apply - the application-provided preconditioning routine
737: Calling sequence of apply:
738: .vb
739: PetscErrorCode apply (void *ptr,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
740: .ve
742: + ptr - the application context
743: . b - right-hand-side
744: . x - current iterate
745: . r - work space
746: . rtol - relative tolerance of residual norm to stop at
747: . abstol - absolute tolerance of residual norm to stop at
748: . dtol - if residual norm increases by this factor than return
749: - maxits - number of iterations to run
751: Level: developer
753: .keywords: PC, shell, set, apply, Richardson, user-provided
755: .seealso: PCShellSetApply(), PCShellSetContext()
756: @*/
757: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyRichardson(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt))
758: {
759: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt));
763: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",(void (**)(void))&f);
764: if (f) {
765: (*f)(pc,apply);
766: }
767: return(0);
768: }
770: /*MC
771: PCSHELL - Creates a new preconditioner class for use with your
772: own private data storage format.
774: Level: advanced
776: Concepts: providing your own preconditioner
778: Usage:
779: $ PetscErrorCode (*mult)(void*,Vec,Vec);
780: $ PetscErrorCode (*setup)(void*);
781: $ PCCreate(comm,&pc);
782: $ PCSetType(pc,PCSHELL);
783: $ PCShellSetApply(pc,mult);
784: $ PCShellSetContext(pc,ctx)
785: $ PCShellSetSetUp(pc,setup); (optional)
787: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
788: MATSHELL, PCShellSetSetUp(), PCShellSetApply(), PCShellSetView(),
789: PCShellSetApplyTranspose(), PCShellSetName(), PCShellSetApplyRichardson(),
790: PCShellGetName(), PCShellSetContext(), PCShellGetContext()
791: M*/
796: PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_Shell(PC pc)
797: {
799: PC_Shell *shell;
802: pc->ops->destroy = PCDestroy_Shell;
803: PetscNew(PC_Shell,&shell);
804: PetscLogObjectMemory(pc,sizeof(PC_Shell));
805: pc->data = (void*)shell;
806: pc->name = 0;
808: pc->ops->apply = PCApply_Shell;
809: pc->ops->view = PCView_Shell;
810: pc->ops->applytranspose = PCApplyTranspose_Shell;
811: pc->ops->applyrichardson = 0;
812: pc->ops->setup = PCSetUp_Shell;
813: pc->ops->presolve = 0;
814: pc->ops->postsolve = 0;
815: pc->ops->view = PCView_Shell;
817: shell->apply = 0;
818: shell->applytranspose = 0;
819: shell->name = 0;
820: shell->applyrich = 0;
821: shell->presolve = 0;
822: shell->postsolve = 0;
823: shell->ctx = 0;
824: shell->setup = 0;
825: shell->view = 0;
826: shell->destroy = 0;
828: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetDestroy_C","PCShellSetDestroy_Shell",
829: PCShellSetDestroy_Shell);
830: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetSetUp_C","PCShellSetSetUp_Shell",
831: PCShellSetSetUp_Shell);
832: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApply_C","PCShellSetApply_Shell",
833: PCShellSetApply_Shell);
834: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetPreSolve_C","PCShellSetPreSolve_Shell",
835: PCShellSetPreSolve_Shell);
836: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetPostSolve_C","PCShellSetPostSolve_Shell",
837: PCShellSetPostSolve_Shell);
838: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetView_C","PCShellSetView_Shell",
839: PCShellSetView_Shell);
840: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApplyTranspose_C","PCShellSetApplyTranspose_Shell",
841: PCShellSetApplyTranspose_Shell);
842: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetName_C","PCShellSetName_Shell",
843: PCShellSetName_Shell);
844: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellGetName_C","PCShellGetName_Shell",
845: PCShellGetName_Shell);
846: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApplyRichardson_C","PCShellSetApplyRichardson_Shell",
847: PCShellSetApplyRichardson_Shell);
848: return(0);
849: }