Actual source code: shellpc.c
petsc-3.7.1 2016-05-15
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 PCPreSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
182: {
183: PC_Shell *shell = (PC_Shell*)pc->data;
187: if (!shell->presolve) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No presolve() routine provided to Shell PC");
188: PetscStackCall("PCSHELL user function presolve()",(*shell->presolve)(pc,ksp,b,x);CHKERRQ(ierr));
189: return(0);
190: }
194: static PetscErrorCode PCPostSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
195: {
196: PC_Shell *shell = (PC_Shell*)pc->data;
200: if (!shell->postsolve) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No postsolve() routine provided to Shell PC");
201: PetscStackCall("PCSHELL user function postsolve()",(*shell->postsolve)(pc,ksp,b,x);CHKERRQ(ierr));
202: return(0);
203: }
207: static PetscErrorCode PCApplyTranspose_Shell(PC pc,Vec x,Vec y)
208: {
209: PC_Shell *shell = (PC_Shell*)pc->data;
210: PetscErrorCode ierr;
211: PetscObjectState instate,outstate;
214: if (!shell->applytranspose) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applytranspose() routine provided to Shell PC");
215: PetscObjectStateGet((PetscObject)y, &instate);
216: PetscStackCall("PCSHELL user function applytranspose()",(*shell->applytranspose)(pc,x,y);CHKERRQ(ierr));
217: PetscObjectStateGet((PetscObject)y, &outstate);
218: if (instate == outstate) {
219: /* increase the state of the output vector since the user did not update its state themself as should have been done */
220: PetscObjectStateIncrease((PetscObject)y);
221: }
222: return(0);
223: }
227: 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)
228: {
229: PetscErrorCode ierr;
230: PC_Shell *shell = (PC_Shell*)pc->data;
231: PetscObjectState instate,outstate;
234: if (!shell->applyrich) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applyrichardson() routine provided to Shell PC");
235: PetscObjectStateGet((PetscObject)y, &instate);
236: PetscStackCall("PCSHELL user function applyrichardson()",(*shell->applyrich)(pc,x,y,w,rtol,abstol,dtol,it,guesszero,outits,reason);CHKERRQ(ierr));
237: PetscObjectStateGet((PetscObject)y, &outstate);
238: if (instate == outstate) {
239: /* increase the state of the output vector since the user did not update its state themself as should have been done */
240: PetscObjectStateIncrease((PetscObject)y);
241: }
242: return(0);
243: }
247: static PetscErrorCode PCDestroy_Shell(PC pc)
248: {
249: PC_Shell *shell = (PC_Shell*)pc->data;
253: PetscFree(shell->name);
254: if (shell->destroy) PetscStackCall("PCSHELL user function destroy()",(*shell->destroy)(pc);CHKERRQ(ierr));
255: PetscFree(pc->data);
256: return(0);
257: }
261: static PetscErrorCode PCView_Shell(PC pc,PetscViewer viewer)
262: {
263: PC_Shell *shell = (PC_Shell*)pc->data;
265: PetscBool iascii;
268: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
269: if (iascii) {
270: if (shell->name) {
271: PetscViewerASCIIPrintf(viewer," Shell: %s\n",shell->name);
272: } else {
273: PetscViewerASCIIPrintf(viewer," Shell: no name\n");
274: }
275: }
276: if (shell->view) {
277: PetscViewerASCIIPushTab(viewer);
278: (*shell->view)(pc,viewer);
279: PetscViewerASCIIPopTab(viewer);
280: }
281: return(0);
282: }
284: /* ------------------------------------------------------------------------------*/
287: static PetscErrorCode PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
288: {
289: PC_Shell *shell= (PC_Shell*)pc->data;
292: shell->destroy = destroy;
293: return(0);
294: }
298: static PetscErrorCode PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
299: {
300: PC_Shell *shell = (PC_Shell*)pc->data;;
303: shell->setup = setup;
304: if (setup) pc->ops->setup = PCSetUp_Shell;
305: else pc->ops->setup = 0;
306: return(0);
307: }
311: static PetscErrorCode PCShellSetApply_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
312: {
313: PC_Shell *shell = (PC_Shell*)pc->data;
316: shell->apply = apply;
317: return(0);
318: }
322: static PetscErrorCode PCShellSetApplySymmetricLeft_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
323: {
324: PC_Shell *shell = (PC_Shell*)pc->data;
327: shell->applysymmetricleft = apply;
328: return(0);
329: }
333: static PetscErrorCode PCShellSetApplySymmetricRight_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
334: {
335: PC_Shell *shell = (PC_Shell*)pc->data;
338: shell->applysymmetricright = apply;
339: return(0);
340: }
344: static PetscErrorCode PCShellSetApplyBA_Shell(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
345: {
346: PC_Shell *shell = (PC_Shell*)pc->data;
349: shell->applyBA = applyBA;
350: if (applyBA) pc->ops->applyBA = PCApplyBA_Shell;
351: else pc->ops->applyBA = 0;
352: return(0);
353: }
357: static PetscErrorCode PCShellSetPreSolve_Shell(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
358: {
359: PC_Shell *shell = (PC_Shell*)pc->data;
362: shell->presolve = presolve;
363: if (presolve) pc->ops->presolve = PCPreSolve_Shell;
364: else pc->ops->presolve = 0;
365: return(0);
366: }
370: static PetscErrorCode PCShellSetPostSolve_Shell(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
371: {
372: PC_Shell *shell = (PC_Shell*)pc->data;
375: shell->postsolve = postsolve;
376: if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
377: else pc->ops->postsolve = 0;
378: return(0);
379: }
383: static PetscErrorCode PCShellSetView_Shell(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
384: {
385: PC_Shell *shell = (PC_Shell*)pc->data;
388: shell->view = view;
389: return(0);
390: }
394: static PetscErrorCode PCShellSetApplyTranspose_Shell(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
395: {
396: PC_Shell *shell = (PC_Shell*)pc->data;
399: shell->applytranspose = applytranspose;
400: if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
401: else pc->ops->applytranspose = 0;
402: return(0);
403: }
407: static PetscErrorCode PCShellSetApplyRichardson_Shell(PC pc,PetscErrorCode (*applyrich)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*))
408: {
409: PC_Shell *shell = (PC_Shell*)pc->data;
412: shell->applyrich = applyrich;
413: if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
414: else pc->ops->applyrichardson = 0;
415: return(0);
416: }
420: static PetscErrorCode PCShellSetName_Shell(PC pc,const char name[])
421: {
422: PC_Shell *shell = (PC_Shell*)pc->data;
426: PetscFree(shell->name);
427: PetscStrallocpy(name,&shell->name);
428: return(0);
429: }
433: static PetscErrorCode PCShellGetName_Shell(PC pc,const char *name[])
434: {
435: PC_Shell *shell = (PC_Shell*)pc->data;
438: *name = shell->name;
439: return(0);
440: }
442: /* -------------------------------------------------------------------------------*/
446: /*@C
447: PCShellSetDestroy - Sets routine to use to destroy the user-provided
448: application context.
450: Logically Collective on PC
452: Input Parameters:
453: + pc - the preconditioner context
454: . destroy - the application-provided destroy routine
456: Calling sequence of destroy:
457: .vb
458: PetscErrorCode destroy (PC)
459: .ve
461: . ptr - the application context
463: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
465: Level: developer
467: .keywords: PC, shell, set, destroy, user-provided
469: .seealso: PCShellSetApply(), PCShellSetContext()
470: @*/
471: PetscErrorCode PCShellSetDestroy(PC pc,PetscErrorCode (*destroy)(PC))
472: {
477: PetscTryMethod(pc,"PCShellSetDestroy_C",(PC,PetscErrorCode (*)(PC)),(pc,destroy));
478: return(0);
479: }
484: /*@C
485: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
486: matrix operator is changed.
488: Logically Collective on PC
490: Input Parameters:
491: + pc - the preconditioner context
492: . setup - the application-provided setup routine
494: Calling sequence of setup:
495: .vb
496: PetscErrorCode setup (PC pc)
497: .ve
499: . pc - the preconditioner, get the application context with PCShellGetContext()
501: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
503: Level: developer
505: .keywords: PC, shell, set, setup, user-provided
507: .seealso: PCShellSetApplyRichardson(), PCShellSetApply(), PCShellSetContext()
508: @*/
509: PetscErrorCode PCShellSetSetUp(PC pc,PetscErrorCode (*setup)(PC))
510: {
515: PetscTryMethod(pc,"PCShellSetSetUp_C",(PC,PetscErrorCode (*)(PC)),(pc,setup));
516: return(0);
517: }
522: /*@C
523: PCShellSetView - Sets routine to use as viewer of shell preconditioner
525: Logically Collective on PC
527: Input Parameters:
528: + pc - the preconditioner context
529: - view - the application-provided view routine
531: Calling sequence of apply:
532: .vb
533: PetscErrorCode view(PC pc,PetscViewer v)
534: .ve
536: + pc - the preconditioner, get the application context with PCShellGetContext()
537: - v - viewer
539: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
541: Level: developer
543: .keywords: PC, shell, set, apply, user-provided
545: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose()
546: @*/
547: PetscErrorCode PCShellSetView(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
548: {
553: PetscTryMethod(pc,"PCShellSetView_C",(PC,PetscErrorCode (*)(PC,PetscViewer)),(pc,view));
554: return(0);
555: }
559: /*@C
560: PCShellSetApply - Sets routine to use as preconditioner.
562: Logically Collective on PC
564: Input Parameters:
565: + pc - the preconditioner context
566: - apply - the application-provided preconditioning routine
568: Calling sequence of apply:
569: .vb
570: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
571: .ve
573: + pc - the preconditioner, get the application context with PCShellGetContext()
574: . xin - input vector
575: - xout - output vector
577: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
579: Developer Notes: There should also be a PCShellSetApplySymmetricRight() and PCShellSetApplySymmetricLeft().
581: Level: developer
583: .keywords: PC, shell, set, apply, user-provided
585: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApplyBA()
586: @*/
587: PetscErrorCode PCShellSetApply(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
588: {
593: PetscTryMethod(pc,"PCShellSetApply_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
594: return(0);
595: }
599: /*@C
600: PCShellSetApplySymmetricLeft - Sets routine to use as left preconditioner (when the PC_SYMMETRIC is used).
602: Logically Collective on PC
604: Input Parameters:
605: + pc - the preconditioner context
606: - apply - the application-provided left preconditioning routine
608: Calling sequence of apply:
609: .vb
610: PetscErrorCode apply (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: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
624: @*/
625: PetscErrorCode PCShellSetApplySymmetricLeft(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
626: {
631: PetscTryMethod(pc,"PCShellSetApplySymmetricLeft_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
632: return(0);
633: }
637: /*@C
638: PCShellSetApply - Sets routine to use as right preconditioner (when the PC_SYMMETRIC is used).
640: Logically Collective on PC
642: Input Parameters:
643: + pc - the preconditioner context
644: - apply - the application-provided right preconditioning routine
646: Calling sequence of apply:
647: .vb
648: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
649: .ve
651: + pc - the preconditioner, get the application context with PCShellGetContext()
652: . xin - input vector
653: - xout - output vector
655: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
657: Level: developer
659: .keywords: PC, shell, set, apply, user-provided
661: .seealso: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
662: @*/
663: PetscErrorCode PCShellSetApplySymmetricRight(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
664: {
669: PetscTryMethod(pc,"PCShellSetApplySymmetricRight_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
670: return(0);
671: }
675: /*@C
676: PCShellSetApplyBA - Sets routine to use as preconditioner times operator.
678: Logically Collective on PC
680: Input Parameters:
681: + pc - the preconditioner context
682: - applyBA - the application-provided BA routine
684: Calling sequence of apply:
685: .vb
686: PetscErrorCode applyBA (PC pc,Vec xin,Vec xout)
687: .ve
689: + pc - the preconditioner, get the application context with PCShellGetContext()
690: . xin - input vector
691: - xout - output vector
693: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
695: Level: developer
697: .keywords: PC, shell, set, apply, user-provided
699: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApply()
700: @*/
701: PetscErrorCode PCShellSetApplyBA(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
702: {
707: PetscTryMethod(pc,"PCShellSetApplyBA_C",(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec)),(pc,applyBA));
708: return(0);
709: }
713: /*@C
714: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
716: Logically Collective on PC
718: Input Parameters:
719: + pc - the preconditioner context
720: - apply - the application-provided preconditioning transpose routine
722: Calling sequence of apply:
723: .vb
724: PetscErrorCode applytranspose (PC pc,Vec xin,Vec xout)
725: .ve
727: + pc - the preconditioner, get the application context with PCShellGetContext()
728: . xin - input vector
729: - xout - output vector
731: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
733: Level: developer
735: Notes:
736: Uses the same context variable as PCShellSetApply().
738: .keywords: PC, shell, set, apply, user-provided
740: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApply(), PCSetContext(), PCShellSetApplyBA()
741: @*/
742: PetscErrorCode PCShellSetApplyTranspose(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
743: {
748: PetscTryMethod(pc,"PCShellSetApplyTranspose_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,applytranspose));
749: return(0);
750: }
754: /*@C
755: PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
756: applied. This usually does something like scale the linear system in some application
757: specific way.
759: Logically Collective on PC
761: Input Parameters:
762: + pc - the preconditioner context
763: - presolve - the application-provided presolve routine
765: Calling sequence of presolve:
766: .vb
767: PetscErrorCode presolve (PC,KSP ksp,Vec b,Vec x)
768: .ve
770: + pc - the preconditioner, get the application context with PCShellGetContext()
771: . xin - input vector
772: - xout - output vector
774: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
776: Level: developer
778: .keywords: PC, shell, set, apply, user-provided
780: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPostSolve(), PCShellSetContext()
781: @*/
782: PetscErrorCode PCShellSetPreSolve(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
783: {
788: PetscTryMethod(pc,"PCShellSetPreSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,presolve));
789: return(0);
790: }
794: /*@C
795: PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
796: applied. This usually does something like scale the linear system in some application
797: specific way.
799: Logically Collective on PC
801: Input Parameters:
802: + pc - the preconditioner context
803: - postsolve - the application-provided presolve routine
805: Calling sequence of postsolve:
806: .vb
807: PetscErrorCode postsolve(PC,KSP ksp,Vec b,Vec x)
808: .ve
810: + pc - the preconditioner, get the application context with PCShellGetContext()
811: . xin - input vector
812: - xout - output vector
814: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
816: Level: developer
818: .keywords: PC, shell, set, apply, user-provided
820: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPreSolve(), PCShellSetContext()
821: @*/
822: PetscErrorCode PCShellSetPostSolve(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
823: {
828: PetscTryMethod(pc,"PCShellSetPostSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,postsolve));
829: return(0);
830: }
834: /*@C
835: PCShellSetName - Sets an optional name to associate with a shell
836: preconditioner.
838: Not Collective
840: Input Parameters:
841: + pc - the preconditioner context
842: - name - character string describing shell preconditioner
844: Level: developer
846: .keywords: PC, shell, set, name, user-provided
848: .seealso: PCShellGetName()
849: @*/
850: PetscErrorCode PCShellSetName(PC pc,const char name[])
851: {
856: PetscTryMethod(pc,"PCShellSetName_C",(PC,const char []),(pc,name));
857: return(0);
858: }
862: /*@C
863: PCShellGetName - Gets an optional name that the user has set for a shell
864: preconditioner.
866: Not Collective
868: Input Parameter:
869: . pc - the preconditioner context
871: Output Parameter:
872: . name - character string describing shell preconditioner (you should not free this)
874: Level: developer
876: .keywords: PC, shell, get, name, user-provided
878: .seealso: PCShellSetName()
879: @*/
880: PetscErrorCode PCShellGetName(PC pc,const char *name[])
881: {
887: PetscUseMethod(pc,"PCShellGetName_C",(PC,const char*[]),(pc,name));
888: return(0);
889: }
893: /*@C
894: PCShellSetApplyRichardson - Sets routine to use as preconditioner
895: in Richardson iteration.
897: Logically Collective on PC
899: Input Parameters:
900: + pc - the preconditioner context
901: - apply - the application-provided preconditioning routine
903: Calling sequence of apply:
904: .vb
905: PetscErrorCode apply (PC pc,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
906: .ve
908: + pc - the preconditioner, get the application context with PCShellGetContext()
909: . b - right-hand-side
910: . x - current iterate
911: . r - work space
912: . rtol - relative tolerance of residual norm to stop at
913: . abstol - absolute tolerance of residual norm to stop at
914: . dtol - if residual norm increases by this factor than return
915: - maxits - number of iterations to run
917: Notes: the function MUST return an error code of 0 on success and nonzero on failure.
919: Level: developer
921: .keywords: PC, shell, set, apply, Richardson, user-provided
923: .seealso: PCShellSetApply(), PCShellSetContext()
924: @*/
925: PetscErrorCode PCShellSetApplyRichardson(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*))
926: {
931: PetscTryMethod(pc,"PCShellSetApplyRichardson_C",(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*)),(pc,apply));
932: return(0);
933: }
935: /*MC
936: PCSHELL - Creates a new preconditioner class for use with your
937: own private data storage format.
939: Level: advanced
940: >
941: Concepts: providing your own preconditioner
943: Usage:
944: $ extern PetscErrorCode apply(PC,Vec,Vec);
945: $ extern PetscErrorCode applyba(PC,PCSide,Vec,Vec,Vec);
946: $ extern PetscErrorCode applytranspose(PC,Vec,Vec);
947: $ extern PetscErrorCode setup(PC);
948: $ extern PetscErrorCode destroy(PC);
949: $
950: $ PCCreate(comm,&pc);
951: $ PCSetType(pc,PCSHELL);
952: $ PCShellSetContext(pc,ctx)
953: $ PCShellSetApply(pc,apply);
954: $ PCShellSetApplyBA(pc,applyba); (optional)
955: $ PCShellSetApplyTranspose(pc,applytranspose); (optional)
956: $ PCShellSetSetUp(pc,setup); (optional)
957: $ PCShellSetDestroy(pc,destroy); (optional)
959: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
960: MATSHELL, PCShellSetSetUp(), PCShellSetApply(), PCShellSetView(),
961: PCShellSetApplyTranspose(), PCShellSetName(), PCShellSetApplyRichardson(),
962: PCShellGetName(), PCShellSetContext(), PCShellGetContext(), PCShellSetApplyBA()
963: M*/
967: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
968: {
970: PC_Shell *shell;
973: PetscNewLog(pc,&shell);
974: pc->data = (void*)shell;
976: pc->ops->destroy = PCDestroy_Shell;
977: pc->ops->view = PCView_Shell;
978: pc->ops->apply = PCApply_Shell;
979: pc->ops->applysymmetricleft = PCApplySymmetricLeft_Shell;
980: pc->ops->applysymmetricright = PCApplySymmetricRight_Shell;
981: pc->ops->applytranspose = 0;
982: pc->ops->applyrichardson = 0;
983: pc->ops->setup = 0;
984: pc->ops->presolve = 0;
985: pc->ops->postsolve = 0;
987: shell->apply = 0;
988: shell->applytranspose = 0;
989: shell->name = 0;
990: shell->applyrich = 0;
991: shell->presolve = 0;
992: shell->postsolve = 0;
993: shell->ctx = 0;
994: shell->setup = 0;
995: shell->view = 0;
996: shell->destroy = 0;
997: shell->applysymmetricleft = 0;
998: shell->applysymmetricright = 0;
1000: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetDestroy_C",PCShellSetDestroy_Shell);
1001: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetSetUp_C",PCShellSetSetUp_Shell);
1002: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApply_C",PCShellSetApply_Shell);
1003: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricLeft_C",PCShellSetApplySymmetricLeft_Shell);
1004: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricRight_C",PCShellSetApplySymmetricRight_Shell);
1005: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyBA_C",PCShellSetApplyBA_Shell);
1006: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPreSolve_C",PCShellSetPreSolve_Shell);
1007: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPostSolve_C",PCShellSetPostSolve_Shell);
1008: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetView_C",PCShellSetView_Shell);
1009: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",PCShellSetApplyTranspose_Shell);
1010: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetName_C",PCShellSetName_Shell);
1011: PetscObjectComposeFunction((PetscObject)pc,"PCShellGetName_C",PCShellGetName_Shell);
1012: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",PCShellSetApplyRichardson_Shell);
1013: return(0);
1014: }