Actual source code: shellpc.c
1: #define PETSCKSP_DLL
4: /*
5: This provides a simple shell for Fortran (and C programmers) to
6: create their own preconditioner without writing much interface code.
7: */
9: #include private/pcimpl.h
10: #include private/vecimpl.h
13: typedef struct {
14: void *ctx; /* user provided contexts for preconditioner */
15: PetscErrorCode (*destroy)(PC);
16: PetscErrorCode (*setup)(PC);
17: PetscErrorCode (*apply)(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,PetscTruth,PetscInt*,PCRichardsonConvergedReason*);
24: char *name;
25: } 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 PCCreateShell()
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
45:
46: .keywords: PC, shell, get, context
48: .seealso: PCCreateShell(), PCShellSetContext()
49: @*/
50: PetscErrorCode PCShellGetContext(PC pc,void **ctx)
51: {
53: PetscTruth flg;
58: PetscTypeCompare((PetscObject)pc,PCSHELL,&flg);
59: if (!flg) *ctx = 0;
60: else *ctx = ((PC_Shell*)(pc->data))->ctx;
61: return(0);
62: }
66: /*@
67: PCShellSetContext - sets the context for a shell PC
69: Collective on PC
71: Input Parameters:
72: + pc - the shell PC
73: - ctx - the context
75: Level: advanced
77: Fortran Notes: The context can only be an integer or a PetscObject
78: unfortunately it cannot be a Fortran array or derived type.
81: .seealso: PCCreateShell(), PCShellGetContext()
82: @*/
83: PetscErrorCode PCShellSetContext(PC pc,void *ctx)
84: {
85: PC_Shell *shell;
87: PetscTruth flg;
91: shell = (PC_Shell*)pc->data;
92: PetscTypeCompare((PetscObject)pc,PCSHELL,&flg);
93: if (flg) {
94: shell->ctx = ctx;
95: }
96: return(0);
97: }
101: static PetscErrorCode PCSetUp_Shell(PC pc)
102: {
103: PC_Shell *shell;
107: shell = (PC_Shell*)pc->data;
108: if (!shell->setup) SETERRQ(PETSC_ERR_USER,"No setup() routine provided to Shell PC");
109: PetscStackPush("PCSHELL user function setup()");
110: CHKMEMQ;
111: (*shell->setup)(pc);
112: CHKMEMQ;
113: PetscStackPop;
114: return(0);
115: }
119: static PetscErrorCode PCApply_Shell(PC pc,Vec x,Vec y)
120: {
121: PC_Shell *shell;
125: shell = (PC_Shell*)pc->data;
126: if (!shell->apply) SETERRQ(PETSC_ERR_USER,"No apply() routine provided to Shell PC");
127: PetscStackPush("PCSHELL user function apply()");
128: CHKMEMQ;
129: (*shell->apply)(pc,x,y);
130: CHKMEMQ;
131: PetscStackPop;
132: return(0);
133: }
137: static PetscErrorCode PCApplyBA_Shell(PC pc,PCSide side,Vec x,Vec y,Vec w)
138: {
139: PC_Shell *shell;
143: shell = (PC_Shell*)pc->data;
144: if (!shell->applyBA) SETERRQ(PETSC_ERR_USER,"No applyBA() routine provided to Shell PC");
145: PetscStackPush("PCSHELL user function applyBA()");
146: CHKMEMQ;
147: (*shell->applyBA)(pc,side,x,y,w);
148: CHKMEMQ;
149: PetscStackPop;
150: return(0);
151: }
155: static PetscErrorCode PCPreSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
156: {
157: PC_Shell *shell;
161: shell = (PC_Shell*)pc->data;
162: if (!shell->presolve) SETERRQ(PETSC_ERR_USER,"No presolve() routine provided to Shell PC");
163: PetscStackPush("PCSHELL user function presolve()");
164: CHKMEMQ;
165: (*shell->presolve)(pc,ksp,b,x);
166: CHKMEMQ;
167: PetscStackPop;
168: return(0);
169: }
173: static PetscErrorCode PCPostSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
174: {
175: PC_Shell *shell;
179: shell = (PC_Shell*)pc->data;
180: if (!shell->postsolve) SETERRQ(PETSC_ERR_USER,"No postsolve() routine provided to Shell PC");
181: PetscStackPush("PCSHELL user function postsolve()");
182: CHKMEMQ;
183: (*shell->postsolve)(pc,ksp,b,x);
184: CHKMEMQ;
185: PetscStackPop;
186: return(0);
187: }
191: static PetscErrorCode PCApplyTranspose_Shell(PC pc,Vec x,Vec y)
192: {
193: PC_Shell *shell;
197: shell = (PC_Shell*)pc->data;
198: if (!shell->applytranspose) SETERRQ(PETSC_ERR_USER,"No applytranspose() routine provided to Shell PC");
199: PetscStackPush("PCSHELL user function applytranspose()");
200: CHKMEMQ;
201: (*shell->applytranspose)(pc,x,y);
202: CHKMEMQ;
203: PetscStackPop;
204: return(0);
205: }
209: static PetscErrorCode PCApplyRichardson_Shell(PC pc,Vec x,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt it,PetscTruth guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason)
210: {
212: PC_Shell *shell;
215: shell = (PC_Shell*)pc->data;
216: if (!shell->applyrich) SETERRQ(PETSC_ERR_USER,"No applyrichardson() routine provided to Shell PC");
217: PetscStackPush("PCSHELL user function applyrichardson()");
218: CHKMEMQ;
219: (*shell->applyrich)(pc,x,y,w,rtol,abstol,dtol,it,guesszero,outits,reason);
220: CHKMEMQ;
221: PetscStackPop;
222: return(0);
223: }
227: static PetscErrorCode PCDestroy_Shell(PC pc)
228: {
229: PC_Shell *shell = (PC_Shell*)pc->data;
233: PetscStrfree(shell->name);
234: if (shell->destroy) {
235: (*shell->destroy)(pc);
236: }
237: PetscFree(shell);
238: return(0);
239: }
243: static PetscErrorCode PCView_Shell(PC pc,PetscViewer viewer)
244: {
245: PC_Shell *shell = (PC_Shell*)pc->data;
247: PetscTruth iascii;
250: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
251: if (iascii) {
252: if (shell->name) {PetscViewerASCIIPrintf(viewer," Shell: %s\n",shell->name);}
253: else {PetscViewerASCIIPrintf(viewer," Shell: no name\n");}
254: }
255: if (shell->view) {
256: PetscViewerASCIIPushTab(viewer);
257: (*shell->view)(pc,viewer);
258: PetscViewerASCIIPopTab(viewer);
259: }
260: return(0);
261: }
263: /* ------------------------------------------------------------------------------*/
267: PetscErrorCode PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
268: {
269: PC_Shell *shell;
272: shell = (PC_Shell*)pc->data;
273: shell->destroy = destroy;
274: return(0);
275: }
281: PetscErrorCode PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
282: {
283: PC_Shell *shell;
286: shell = (PC_Shell*)pc->data;
287: shell->setup = setup;
288: if (setup) pc->ops->setup = PCSetUp_Shell;
289: else pc->ops->setup = 0;
290: return(0);
291: }
297: PetscErrorCode PCShellSetApply_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
298: {
299: PC_Shell *shell;
302: shell = (PC_Shell*)pc->data;
303: shell->apply = apply;
304: return(0);
305: }
311: PetscErrorCode PCShellSetApplyBA_Shell(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
312: {
313: PC_Shell *shell;
316: shell = (PC_Shell*)pc->data;
317: shell->applyBA = applyBA;
318: if (applyBA) pc->ops->applyBA = PCApplyBA_Shell;
319: else pc->ops->applyBA = 0;
320: return(0);
321: }
327: PetscErrorCode PCShellSetPreSolve_Shell(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
328: {
329: PC_Shell *shell;
332: shell = (PC_Shell*)pc->data;
333: shell->presolve = presolve;
334: if (presolve) pc->ops->presolve = PCPreSolve_Shell;
335: else pc->ops->presolve = 0;
336: return(0);
337: }
343: PetscErrorCode PCShellSetPostSolve_Shell(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
344: {
345: PC_Shell *shell;
348: shell = (PC_Shell*)pc->data;
349: shell->postsolve = postsolve;
350: if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
351: else pc->ops->postsolve = 0;
352: return(0);
353: }
359: PetscErrorCode PCShellSetView_Shell(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
360: {
361: PC_Shell *shell;
364: shell = (PC_Shell*)pc->data;
365: shell->view = view;
366: return(0);
367: }
373: PetscErrorCode PCShellSetApplyTranspose_Shell(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
374: {
375: PC_Shell *shell;
378: shell = (PC_Shell*)pc->data;
379: shell->applytranspose = applytranspose;
380: if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
381: else pc->ops->applytranspose = 0;
382: return(0);
383: }
389: PetscErrorCode PCShellSetApplyRichardson_Shell(PC pc,PetscErrorCode (*applyrich)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscTruth,PetscInt*,PCRichardsonConvergedReason*))
390: {
391: PC_Shell *shell;
394: shell = (PC_Shell*)pc->data;
395: shell->applyrich = applyrich;
396: if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
397: else pc->ops->applyrichardson = 0;
398: return(0);
399: }
405: PetscErrorCode PCShellSetName_Shell(PC pc,const char name[])
406: {
407: PC_Shell *shell;
411: shell = (PC_Shell*)pc->data;
412: PetscStrfree(shell->name);
413: PetscStrallocpy(name,&shell->name);
414: return(0);
415: }
421: PetscErrorCode PCShellGetName_Shell(PC pc,char *name[])
422: {
423: PC_Shell *shell;
426: shell = (PC_Shell*)pc->data;
427: *name = shell->name;
428: return(0);
429: }
432: /* -------------------------------------------------------------------------------*/
436: /*@C
437: PCShellSetDestroy - Sets routine to use to destroy the user-provided
438: application context.
440: Collective on PC
442: Input Parameters:
443: + pc - the preconditioner context
444: . destroy - the application-provided destroy routine
446: Calling sequence of destroy:
447: .vb
448: PetscErrorCode destroy (PC)
449: .ve
451: . ptr - the application context
453: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
455: Level: developer
457: .keywords: PC, shell, set, destroy, user-provided
459: .seealso: PCShellSetApply(), PCShellSetContext()
460: @*/
461: PetscErrorCode PCShellSetDestroy(PC pc,PetscErrorCode (*destroy)(PC))
462: {
463: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC));
467: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetDestroy_C",(void (**)(void))&f);
468: if (f) {
469: (*f)(pc,destroy);
470: }
471: return(0);
472: }
477: /*@C
478: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
479: matrix operator is changed.
481: Collective on PC
483: Input Parameters:
484: + pc - the preconditioner context
485: . setup - the application-provided setup routine
487: Calling sequence of setup:
488: .vb
489: PetscErrorCode setup (PC pc)
490: .ve
492: . pc - the preconditioner, get the application context with PCShellGetContext()
494: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
496: Level: developer
498: .keywords: PC, shell, set, setup, user-provided
500: .seealso: PCShellSetApplyRichardson(), PCShellSetApply(), PCShellSetContext()
501: @*/
502: PetscErrorCode PCShellSetSetUp(PC pc,PetscErrorCode (*setup)(PC))
503: {
504: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC));
508: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetSetUp_C",(void (**)(void))&f);
509: if (f) {
510: (*f)(pc,setup);
511: }
512: return(0);
513: }
518: /*@C
519: PCShellSetView - Sets routine to use as viewer of shell preconditioner
521: Collective on PC
523: Input Parameters:
524: + pc - the preconditioner context
525: - view - the application-provided view routine
527: Calling sequence of apply:
528: .vb
529: PetscErrorCode view(PC pc,PetscViewer v)
530: .ve
532: + pc - the preconditioner, get the application context with PCShellGetContext()
533: - v - viewer
535: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
537: Level: developer
539: .keywords: PC, shell, set, apply, user-provided
541: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose()
542: @*/
543: PetscErrorCode PCShellSetView(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
544: {
545: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC,PetscViewer));
549: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetView_C",(void (**)(void))&f);
550: if (f) {
551: (*f)(pc,view);
552: }
553: return(0);
554: }
558: /*@C
559: PCShellSetApply - Sets routine to use as preconditioner.
561: Collective on PC
563: Input Parameters:
564: + pc - the preconditioner context
565: - apply - the application-provided preconditioning routine
567: Calling sequence of apply:
568: .vb
569: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
570: .ve
572: + pc - the preconditioner, get the application context with PCShellGetContext()
573: . xin - input vector
574: - xout - output vector
576: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
578: Level: developer
580: .keywords: PC, shell, set, apply, user-provided
582: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApplyBA()
583: @*/
584: PetscErrorCode PCShellSetApply(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
585: {
586: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC,Vec,Vec));
590: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApply_C",(void (**)(void))&f);
591: if (f) {
592: (*f)(pc,apply);
593: }
594: return(0);
595: }
599: /*@C
600: PCShellSetApplyBA - Sets routine to use as preconditioner times operator.
602: Collective on PC
604: Input Parameters:
605: + pc - the preconditioner context
606: - applyBA - the application-provided BA routine
608: Calling sequence of apply:
609: .vb
610: PetscErrorCode applyBA (PC pc,Vec xin,Vec xout)
611: .ve
613: + pc - the preconditioner, get the application context with PCShellGetContext()
614: . xin - input vector
615: - xout - output vector
617: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
619: Level: developer
621: .keywords: PC, shell, set, apply, user-provided
623: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApply()
624: @*/
625: PetscErrorCode PCShellSetApplyBA(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
626: {
627: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec));
631: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApplyBA_C",(void (**)(void))&f);
632: if (f) {
633: (*f)(pc,applyBA);
634: }
635: return(0);
636: }
640: /*@C
641: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
643: Collective on PC
645: Input Parameters:
646: + pc - the preconditioner context
647: - apply - the application-provided preconditioning transpose routine
649: Calling sequence of apply:
650: .vb
651: PetscErrorCode applytranspose (PC pc,Vec xin,Vec xout)
652: .ve
654: + pc - the preconditioner, get the application context with PCShellGetContext()
655: . xin - input vector
656: - xout - output vector
658: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
660: Level: developer
662: Notes:
663: Uses the same context variable as PCShellSetApply().
665: .keywords: PC, shell, set, apply, user-provided
667: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApply(), PCSetContext(), PCShellSetApplyBA()
668: @*/
669: PetscErrorCode PCShellSetApplyTranspose(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
670: {
671: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC,Vec,Vec));
675: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",(void (**)(void))&f);
676: if (f) {
677: (*f)(pc,applytranspose);
678: }
679: return(0);
680: }
684: /*@C
685: PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
686: applied. This usually does something like scale the linear system in some application
687: specific way.
689: Collective on PC
691: Input Parameters:
692: + pc - the preconditioner context
693: - presolve - the application-provided presolve routine
695: Calling sequence of presolve:
696: .vb
697: PetscErrorCode presolve (PC,KSP ksp,Vec b,Vec x)
698: .ve
700: + pc - the preconditioner, get the application context with PCShellGetContext()
701: . xin - input vector
702: - xout - output vector
704: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
706: Level: developer
708: .keywords: PC, shell, set, apply, user-provided
710: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPostSolve(), PCShellSetContext()
711: @*/
712: PetscErrorCode PCShellSetPreSolve(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
713: {
714: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec));
718: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetPreSolve_C",(void (**)(void))&f);
719: if (f) {
720: (*f)(pc,presolve);
721: }
722: return(0);
723: }
727: /*@C
728: PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
729: applied. This usually does something like scale the linear system in some application
730: specific way.
732: Collective on PC
734: Input Parameters:
735: + pc - the preconditioner context
736: - postsolve - the application-provided presolve routine
738: Calling sequence of postsolve:
739: .vb
740: PetscErrorCode postsolve(PC,KSP ksp,Vec b,Vec x)
741: .ve
743: + pc - the preconditioner, get the application context with PCShellGetContext()
744: . xin - input vector
745: - xout - output vector
747: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
749: Level: developer
751: .keywords: PC, shell, set, apply, user-provided
753: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPreSolve(), PCShellSetContext()
754: @*/
755: PetscErrorCode PCShellSetPostSolve(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
756: {
757: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec));
761: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetPostSolve_C",(void (**)(void))&f);
762: if (f) {
763: (*f)(pc,postsolve);
764: }
765: return(0);
766: }
770: /*@C
771: PCShellSetName - Sets an optional name to associate with a shell
772: preconditioner.
774: Not Collective
776: Input Parameters:
777: + pc - the preconditioner context
778: - name - character string describing shell preconditioner
780: Level: developer
782: .keywords: PC, shell, set, name, user-provided
784: .seealso: PCShellGetName()
785: @*/
786: PetscErrorCode PCShellSetName(PC pc,const char name[])
787: {
788: PetscErrorCode ierr,(*f)(PC,const char []);
792: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetName_C",(void (**)(void))&f);
793: if (f) {
794: (*f)(pc,name);
795: }
796: return(0);
797: }
801: /*@C
802: PCShellGetName - Gets an optional name that the user has set for a shell
803: preconditioner.
805: Not Collective
807: Input Parameter:
808: . pc - the preconditioner context
810: Output Parameter:
811: . name - character string describing shell preconditioner (you should not free this)
813: Level: developer
815: .keywords: PC, shell, get, name, user-provided
817: .seealso: PCShellSetName()
818: @*/
819: PetscErrorCode PCShellGetName(PC pc,char *name[])
820: {
821: PetscErrorCode ierr,(*f)(PC,char *[]);
826: PetscObjectQueryFunction((PetscObject)pc,"PCShellGetName_C",(void (**)(void))&f);
827: if (f) {
828: (*f)(pc,name);
829: } else {
830: SETERRQ(PETSC_ERR_ARG_WRONG,"Not shell preconditioner, cannot get name");
831: }
832: return(0);
833: }
837: /*@C
838: PCShellSetApplyRichardson - Sets routine to use as preconditioner
839: in Richardson iteration.
841: Collective on PC
843: Input Parameters:
844: + pc - the preconditioner context
845: - apply - the application-provided preconditioning routine
847: Calling sequence of apply:
848: .vb
849: PetscErrorCode apply (PC pc,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
850: .ve
852: + pc - the preconditioner, get the application context with PCShellGetContext()
853: . b - right-hand-side
854: . x - current iterate
855: . r - work space
856: . rtol - relative tolerance of residual norm to stop at
857: . abstol - absolute tolerance of residual norm to stop at
858: . dtol - if residual norm increases by this factor than return
859: - maxits - number of iterations to run
861: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
863: Level: developer
865: .keywords: PC, shell, set, apply, Richardson, user-provided
867: .seealso: PCShellSetApply(), PCShellSetContext()
868: @*/
869: PetscErrorCode PCShellSetApplyRichardson(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscTruth,PetscInt*,PCRichardsonConvergedReason*))
870: {
871: PetscErrorCode ierr,(*f)(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscTruth,PetscInt*,PCRichardsonConvergedReason*));
875: PetscObjectQueryFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",(void (**)(void))&f);
876: if (f) {
877: (*f)(pc,apply);
878: }
879: return(0);
880: }
882: /*MC
883: PCSHELL - Creates a new preconditioner class for use with your
884: own private data storage format.
886: Level: advanced
887: >
888: Concepts: providing your own preconditioner
890: Usage:
896: $
897: $ PCCreate(comm,&pc);
898: $ PCSetType(pc,PCSHELL);
899: $ PCShellSetContext(pc,ctx)
900: $ PCShellSetApply(pc,apply);
901: $ PCShellSetApplyBA(pc,applyba); (optional)
902: $ PCShellSetApplyTranspose(pc,applytranspose); (optional)
903: $ PCShellSetSetUp(pc,setup); (optional)
904: $ PCShellSetDestroy(pc,destroy); (optional)
906: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
907: MATSHELL, PCShellSetSetUp(), PCShellSetApply(), PCShellSetView(),
908: PCShellSetApplyTranspose(), PCShellSetName(), PCShellSetApplyRichardson(),
909: PCShellGetName(), PCShellSetContext(), PCShellGetContext(), PCShellSetApplyBA()
910: M*/
915: PetscErrorCode PCCreate_Shell(PC pc)
916: {
918: PC_Shell *shell;
921: PetscNewLog(pc,PC_Shell,&shell);
922: pc->data = (void*)shell;
924: pc->ops->destroy = PCDestroy_Shell;
925: pc->ops->view = PCView_Shell;
926: pc->ops->apply = PCApply_Shell;
927: pc->ops->applytranspose = 0;
928: pc->ops->applyrichardson = 0;
929: pc->ops->setup = 0;
930: pc->ops->presolve = 0;
931: pc->ops->postsolve = 0;
933: shell->apply = 0;
934: shell->applytranspose = 0;
935: shell->name = 0;
936: shell->applyrich = 0;
937: shell->presolve = 0;
938: shell->postsolve = 0;
939: shell->ctx = 0;
940: shell->setup = 0;
941: shell->view = 0;
942: shell->destroy = 0;
944: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetDestroy_C","PCShellSetDestroy_Shell",
945: PCShellSetDestroy_Shell);
946: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetSetUp_C","PCShellSetSetUp_Shell",
947: PCShellSetSetUp_Shell);
948: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApply_C","PCShellSetApply_Shell",
949: PCShellSetApply_Shell);
950: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApplyBA_C","PCShellSetApplyBA_Shell",
951: PCShellSetApplyBA_Shell);
952: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetPreSolve_C","PCShellSetPreSolve_Shell",
953: PCShellSetPreSolve_Shell);
954: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetPostSolve_C","PCShellSetPostSolve_Shell",
955: PCShellSetPostSolve_Shell);
956: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetView_C","PCShellSetView_Shell",
957: PCShellSetView_Shell);
958: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApplyTranspose_C","PCShellSetApplyTranspose_Shell",
959: PCShellSetApplyTranspose_Shell);
960: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetName_C","PCShellSetName_Shell",
961: PCShellSetName_Shell);
962: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellGetName_C","PCShellGetName_Shell",
963: PCShellGetName_Shell);
964: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCShellSetApplyRichardson_C","PCShellSetApplyRichardson_Shell",
965: PCShellSetApplyRichardson_Shell);
966: return(0);
967: }