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 src/ksp/pc/pcimpl.h
9: #include 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: (*shell->setup)(shell->ctx);
106: }
107: return(0);
108: }
112: static PetscErrorCode PCApply_Shell(PC pc,Vec x,Vec y)
113: {
114: PC_Shell *shell;
118: shell = (PC_Shell*)pc->data;
119: if (!shell->apply) SETERRQ(PETSC_ERR_USER,"No apply() routine provided to Shell PC");
120: (*shell->apply)(shell->ctx,x,y);
121: return(0);
122: }
126: static PetscErrorCode PCPreSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
127: {
128: PC_Shell *shell;
132: shell = (PC_Shell*)pc->data;
133: if (!shell->presolve) SETERRQ(PETSC_ERR_USER,"No presolve() routine provided to Shell PC");
134: (*shell->presolve)(shell->ctx,ksp,b,x);
135: return(0);
136: }
140: static PetscErrorCode PCPostSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
141: {
142: PC_Shell *shell;
146: shell = (PC_Shell*)pc->data;
147: if (!shell->postsolve) SETERRQ(PETSC_ERR_USER,"No postsolve() routine provided to Shell PC");
148: (*shell->postsolve)(shell->ctx,ksp,b,x);
149: return(0);
150: }
154: static PetscErrorCode PCApplyTranspose_Shell(PC pc,Vec x,Vec y)
155: {
156: PC_Shell *shell;
160: shell = (PC_Shell*)pc->data;
161: if (!shell->applytranspose) SETERRQ(PETSC_ERR_USER,"No applytranspose() routine provided to Shell PC");
162: (*shell->applytranspose)(shell->ctx,x,y);
163: return(0);
164: }
168: static PetscErrorCode PCApplyRichardson_Shell(PC pc,Vec x,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt it)
169: {
171: PC_Shell *shell;
174: shell = (PC_Shell*)pc->data;
175: (*shell->applyrich)(shell->ctx,x,y,w,rtol,abstol,dtol,it);
176: return(0);
177: }
181: static PetscErrorCode PCDestroy_Shell(PC pc)
182: {
183: PC_Shell *shell = (PC_Shell*)pc->data;
187: if (shell->name) {PetscFree(shell->name);}
188: if (shell->destroy) {
189: (*shell->destroy)(shell->ctx);
190: }
191: PetscFree(shell);
192: return(0);
193: }
197: static PetscErrorCode PCView_Shell(PC pc,PetscViewer viewer)
198: {
199: PC_Shell *shell = (PC_Shell*)pc->data;
201: PetscTruth iascii;
204: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
205: if (iascii) {
206: if (shell->name) {PetscViewerASCIIPrintf(viewer," Shell: %s\n",shell->name);}
207: else {PetscViewerASCIIPrintf(viewer," Shell: no name\n");}
208: }
209: if (shell->view) {
210: PetscViewerASCIIPushTab(viewer);
211: (*shell->view)(shell->ctx,viewer);
212: PetscViewerASCIIPopTab(viewer);
213: }
214: return(0);
215: }
217: /* ------------------------------------------------------------------------------*/
221: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(void*))
222: {
223: PC_Shell *shell;
226: shell = (PC_Shell*)pc->data;
227: shell->destroy = destroy;
228: return(0);
229: }
235: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(void*))
236: {
237: PC_Shell *shell;
240: shell = (PC_Shell*)pc->data;
241: shell->setup = setup;
242: return(0);
243: }
249: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApply_Shell(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec))
250: {
251: PC_Shell *shell;
254: shell = (PC_Shell*)pc->data;
255: shell->apply = apply;
256: return(0);
257: }
263: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPreSolve_Shell(PC pc,PetscErrorCode (*presolve)(void*,KSP,Vec,Vec))
264: {
265: PC_Shell *shell;
268: shell = (PC_Shell*)pc->data;
269: shell->presolve = presolve;
270: if (presolve) {
271: pc->ops->presolve = PCPreSolve_Shell;
272: } else {
273: pc->ops->presolve = 0;
274: }
275: return(0);
276: }
282: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPostSolve_Shell(PC pc,PetscErrorCode (*postsolve)(void*,KSP,Vec,Vec))
283: {
284: PC_Shell *shell;
287: shell = (PC_Shell*)pc->data;
288: shell->postsolve = postsolve;
289: if (postsolve) {
290: pc->ops->postsolve = PCPostSolve_Shell;
291: } else {
292: pc->ops->postsolve = 0;
293: }
294: return(0);
295: }
301: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetView_Shell(PC pc,PetscErrorCode (*view)(void*,PetscViewer))
302: {
303: PC_Shell *shell;
306: shell = (PC_Shell*)pc->data;
307: shell->view = view;
308: return(0);
309: }
315: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyTranspose_Shell(PC pc,PetscErrorCode (*applytranspose)(void*,Vec,Vec))
316: {
317: PC_Shell *shell;
320: shell = (PC_Shell*)pc->data;
321: shell->applytranspose = applytranspose;
322: return(0);
323: }
329: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetName_Shell(PC pc,const char name[])
330: {
331: PC_Shell *shell;
335: shell = (PC_Shell*)pc->data;
336: PetscStrfree(shell->name);
337: PetscStrallocpy(name,&shell->name);
338: return(0);
339: }
345: PetscErrorCode PETSCKSP_DLLEXPORT PCShellGetName_Shell(PC pc,char *name[])
346: {
347: PC_Shell *shell;
350: shell = (PC_Shell*)pc->data;
351: *name = shell->name;
352: return(0);
353: }
359: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyRichardson_Shell(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt))
360: {
361: PC_Shell *shell;
364: shell = (PC_Shell*)pc->data;
365: pc->ops->applyrichardson = PCApplyRichardson_Shell;
366: shell->applyrich = apply;
367: return(0);
368: }
371: /* -------------------------------------------------------------------------------*/
375: /*@C
376: PCShellSetDestroy - Sets routine to use to destroy the user-provided
377: application context.
379: Collective on PC
381: Input Parameters:
382: + pc - the preconditioner context
383: . destroy - the application-provided destroy routine
385: Calling sequence of destroy:
386: .vb
387: PetscErrorCode destroy (void *ptr)
388: .ve
390: . ptr - the application context
392: Level: developer
394: .keywords: PC, shell, set, destroy, user-provided
396: .seealso: PCShellSetApply(), PCShellSetContext()
397: @*/
398: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetDestroy(PC pc,PetscErrorCode (*destroy)(void*))
399: {
400: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*));
404: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetDestroy_C",(void (**)(void))&f);
405: if (f) {
406: (*f)(pc,destroy);
407: }
408: return(0);
409: }
414: /*@C
415: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
416: matrix operator is changed.
418: Collective on PC
420: Input Parameters:
421: + pc - the preconditioner context
422: . setup - the application-provided setup routine
424: Calling sequence of setup:
425: .vb
426: PetscErrorCode setup (void *ptr)
427: .ve
429: . ptr - the application context
431: Level: developer
433: .keywords: PC, shell, set, setup, user-provided
435: .seealso: PCShellSetApplyRichardson(), PCShellSetApply(), PCShellSetContext()
436: @*/
437: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetSetUp(PC pc,PetscErrorCode (*setup)(void*))
438: {
439: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*));
443: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetSetUp_C",(void (**)(void))&f);
444: if (f) {
445: (*f)(pc,setup);
446: }
447: return(0);
448: }
453: /*@C
454: PCShellSetView - Sets routine to use as viewer of shell preconditioner
456: Collective on PC
458: Input Parameters:
459: + pc - the preconditioner context
460: - view - the application-provided view routine
462: Calling sequence of apply:
463: .vb
464: PetscErrorCode view(void *ptr,PetscViewer v)
465: .ve
467: + ptr - the application context
468: - v - viewer
470: Level: developer
472: .keywords: PC, shell, set, apply, user-provided
474: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose()
475: @*/
476: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetView(PC pc,PetscErrorCode (*view)(void*,PetscViewer))
477: {
478: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,PetscViewer));
482: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetView_C",(void (**)(void))&f);
483: if (f) {
484: (*f)(pc,view);
485: }
486: return(0);
487: }
491: /*@C
492: PCShellSetApply - Sets routine to use as preconditioner.
494: Collective on PC
496: Input Parameters:
497: + pc - the preconditioner context
498: - apply - the application-provided preconditioning routine
500: Calling sequence of apply:
501: .vb
502: PetscErrorCode apply (void *ptr,Vec xin,Vec xout)
503: .ve
505: + ptr - the application context
506: . xin - input vector
507: - xout - output vector
509: Level: developer
511: .keywords: PC, shell, set, apply, user-provided
513: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
514: @*/
515: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApply(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec))
516: {
517: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,Vec,Vec));
521: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApply_C",(void (**)(void))&f);
522: if (f) {
523: (*f)(pc,apply);
524: }
525: return(0);
526: }
530: /*@C
531: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
533: Collective on PC
535: Input Parameters:
536: + pc - the preconditioner context
537: - apply - the application-provided preconditioning transpose routine
539: Calling sequence of apply:
540: .vb
541: PetscErrorCode applytranspose (void *ptr,Vec xin,Vec xout)
542: .ve
544: + ptr - the application context
545: . xin - input vector
546: - xout - output vector
548: Level: developer
550: Notes:
551: Uses the same context variable as PCShellSetApply().
553: .keywords: PC, shell, set, apply, user-provided
555: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApply(), PCSetContext()
556: @*/
557: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyTranspose(PC pc,PetscErrorCode (*applytranspose)(void*,Vec,Vec))
558: {
559: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,Vec,Vec));
563: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",(void (**)(void))&f);
564: if (f) {
565: (*f)(pc,applytranspose);
566: }
567: return(0);
568: }
572: /*@C
573: PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
574: applied. This usually does something like scale the linear system in some application
575: specific way.
577: Collective on PC
579: Input Parameters:
580: + pc - the preconditioner context
581: - presolve - the application-provided presolve routine
583: Calling sequence of presolve:
584: .vb
585: PetscErrorCode presolve (void *ptr,KSP ksp,Vec b,Vec x)
586: .ve
588: + ptr - the application context
589: . xin - input vector
590: - xout - output vector
592: Level: developer
594: .keywords: PC, shell, set, apply, user-provided
596: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPostSolve(), PCShellSetContext()
597: @*/
598: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPreSolve(PC pc,PetscErrorCode (*presolve)(void*,KSP,Vec,Vec))
599: {
600: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,KSP,Vec,Vec));
604: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetPreSolve_C",(void (**)(void))&f);
605: if (f) {
606: (*f)(pc,presolve);
607: }
608: return(0);
609: }
613: /*@C
614: PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
615: applied. This usually does something like scale the linear system in some application
616: specific way.
618: Collective on PC
620: Input Parameters:
621: + pc - the preconditioner context
622: - postsolve - the application-provided presolve routine
624: Calling sequence of postsolve:
625: .vb
626: PetscErrorCode postsolve(void *ptr,KSP ksp,Vec b,Vec x)
627: .ve
629: + ptr - the application context
630: . xin - input vector
631: - xout - output vector
633: Level: developer
635: .keywords: PC, shell, set, apply, user-provided
637: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPreSolve(), PCShellSetContext()
638: @*/
639: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetPostSolve(PC pc,PetscErrorCode (*postsolve)(void*,KSP,Vec,Vec))
640: {
641: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,KSP,Vec,Vec));
645: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetPostSolve_C",(void (**)(void))&f);
646: if (f) {
647: (*f)(pc,postsolve);
648: }
649: return(0);
650: }
654: /*@C
655: PCShellSetName - Sets an optional name to associate with a shell
656: preconditioner.
658: Not Collective
660: Input Parameters:
661: + pc - the preconditioner context
662: - name - character string describing shell preconditioner
664: Level: developer
666: .keywords: PC, shell, set, name, user-provided
668: .seealso: PCShellGetName()
669: @*/
670: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetName(PC pc,const char name[])
671: {
672: PetscErrorCode ierr,(*f)(PC,const char []);
676: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetName_C",(void (**)(void))&f);
677: if (f) {
678: (*f)(pc,name);
679: }
680: return(0);
681: }
685: /*@C
686: PCShellGetName - Gets an optional name that the user has set for a shell
687: preconditioner.
689: Not Collective
691: Input Parameter:
692: . pc - the preconditioner context
694: Output Parameter:
695: . name - character string describing shell preconditioner (you should not free this)
697: Level: developer
699: .keywords: PC, shell, get, name, user-provided
701: .seealso: PCShellSetName()
702: @*/
703: PetscErrorCode PETSCKSP_DLLEXPORT PCShellGetName(PC pc,char *name[])
704: {
705: PetscErrorCode ierr,(*f)(PC,char *[]);
710: PetscObjectQueryFunction((PetscObject)pc,"PCShellGetName_C",(void (**)(void))&f);
711: if (f) {
712: (*f)(pc,name);
713: } else {
714: SETERRQ(PETSC_ERR_ARG_WRONG,"Not shell preconditioner, cannot get name");
715: }
716: return(0);
717: }
721: /*@C
722: PCShellSetApplyRichardson - Sets routine to use as preconditioner
723: in Richardson iteration.
725: Collective on PC
727: Input Parameters:
728: + pc - the preconditioner context
729: - apply - the application-provided preconditioning routine
731: Calling sequence of apply:
732: .vb
733: PetscErrorCode apply (void *ptr,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
734: .ve
736: + ptr - the application context
737: . b - right-hand-side
738: . x - current iterate
739: . r - work space
740: . rtol - relative tolerance of residual norm to stop at
741: . abstol - absolute tolerance of residual norm to stop at
742: . dtol - if residual norm increases by this factor than return
743: - maxits - number of iterations to run
745: Level: developer
747: .keywords: PC, shell, set, apply, Richardson, user-provided
749: .seealso: PCShellSetApply(), PCShellSetContext()
750: @*/
751: PetscErrorCode PETSCKSP_DLLEXPORT PCShellSetApplyRichardson(PC pc,PetscErrorCode (*apply)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt))
752: {
753: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt));
757: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",(void (**)(void))&f);
758: if (f) {
759: (*f)(pc,apply);
760: }
761: return(0);
762: }
764: /*MC
765: PCSHELL - Creates a new preconditioner class for use with your
766: own private data storage format.
768: Level: advanced
770: Concepts: providing your own preconditioner
772: Usage:
773: $ PetscErrorCode (*mult)(void*,Vec,Vec);
774: $ PetscErrorCode (*setup)(void*);
775: $ PCCreate(comm,&pc);
776: $ PCSetType(pc,PCSHELL);
777: $ PCShellSetApply(pc,mult);
778: $ PCShellSetContext(pc,ctx)
779: $ PCShellSetSetUp(pc,setup); (optional)
781: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
782: MATSHELL, PCShellSetUp(), PCShellSetApply(), PCShellSetView(),
783: PCShellSetApplyTranpose(), PCShellSetName(), PCShellSetApplyRichardson(),
784: PCShellGetName(), PCShellSetContext(), PCShellGetContext()
785: M*/
790: PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_Shell(PC pc)
791: {
793: PC_Shell *shell;
796: pc->ops->destroy = PCDestroy_Shell;
797: PetscNew(PC_Shell,&shell);
798: PetscLogObjectMemory(pc,sizeof(PC_Shell));
799: pc->data = (void*)shell;
800: pc->name = 0;
802: pc->ops->apply = PCApply_Shell;
803: pc->ops->view = PCView_Shell;
804: pc->ops->applytranspose = PCApplyTranspose_Shell;
805: pc->ops->applyrichardson = 0;
806: pc->ops->setup = PCSetUp_Shell;
807: pc->ops->presolve = 0;
808: pc->ops->postsolve = 0;
809: pc->ops->view = PCView_Shell;
811: shell->apply = 0;
812: shell->applytranspose = 0;
813: shell->name = 0;
814: shell->applyrich = 0;
815: shell->presolve = 0;
816: shell->postsolve = 0;
817: shell->ctx = 0;
818: shell->setup = 0;
819: shell->view = 0;
820: shell->destroy = 0;
822: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetDestroy_C","PCShellSetDestroy_Shell",
823: PCShellSetDestroy_Shell);
824: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetSetUp_C","PCShellSetSetUp_Shell",
825: PCShellSetSetUp_Shell);
826: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApply_C","PCShellSetApply_Shell",
827: PCShellSetApply_Shell);
828: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetPreSolve_C","PCShellSetPreSolve_Shell",
829: PCShellSetPreSolve_Shell);
830: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetPostSolve_C","PCShellSetPostSolve_Shell",
831: PCShellSetPostSolve_Shell);
832: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetView_C","PCShellSetView_Shell",
833: PCShellSetView_Shell);
834: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApplyTranspose_C","PCShellSetApplyTranspose_Shell",
835: PCShellSetApplyTranspose_Shell);
836: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetName_C","PCShellSetName_Shell",
837: PCShellSetName_Shell);
838: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellGetName_C","PCShellGetName_Shell",
839: PCShellGetName_Shell);
840: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApplyRichardson_C","PCShellSetApplyRichardson_Shell",
841: PCShellSetApplyRichardson_Shell);
842: return(0);
843: }