Actual source code: shellpc.c
petsc-3.7.2 2016-06-05
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: }