Actual source code: shellpc.c
petsc-3.11.4 2019-09-28
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>
9: typedef struct {
10: void *ctx; /* user provided contexts for preconditioner */
12: PetscErrorCode (*destroy)(PC);
13: PetscErrorCode (*setup)(PC);
14: PetscErrorCode (*apply)(PC,Vec,Vec);
15: PetscErrorCode (*applysymmetricleft)(PC,Vec,Vec);
16: PetscErrorCode (*applysymmetricright)(PC,Vec,Vec);
17: PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec);
18: PetscErrorCode (*presolve)(PC,KSP,Vec,Vec);
19: PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec);
20: PetscErrorCode (*view)(PC,PetscViewer);
21: PetscErrorCode (*applytranspose)(PC,Vec,Vec);
22: PetscErrorCode (*applyrich)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*);
24: char *name;
25: } PC_Shell;
27: /*@C
28: PCShellGetContext - Returns the user-provided context associated with a shell PC
30: Not Collective
32: Input Parameter:
33: . pc - should have been created with PCSetType(pc,shell)
35: Output Parameter:
36: . ctx - the user provided context
38: Level: advanced
40: Notes:
41: This routine is intended for use within various shell routines
43: Fortran Notes:
44: To use this from Fortran you must write a Fortran interface definition for this
45: function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.
47: .keywords: PC, shell, get, context
49: .seealso: PCShellSetContext()
50: @*/
51: PetscErrorCode PCShellGetContext(PC pc,void **ctx)
52: {
54: PetscBool flg;
59: PetscObjectTypeCompare((PetscObject)pc,PCSHELL,&flg);
60: if (!flg) *ctx = 0;
61: else *ctx = ((PC_Shell*)(pc->data))->ctx;
62: return(0);
63: }
65: /*@
66: PCShellSetContext - sets the context for a shell PC
68: Logically Collective on PC
70: Input Parameters:
71: + pc - the shell PC
72: - ctx - the context
74: Level: advanced
76: Fortran Notes:
77: To use this from Fortran you must write a Fortran interface definition for this
78: function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.
82: .seealso: PCShellGetContext(), PCSHELL
83: @*/
84: PetscErrorCode PCShellSetContext(PC pc,void *ctx)
85: {
86: PC_Shell *shell = (PC_Shell*)pc->data;
88: PetscBool flg;
92: PetscObjectTypeCompare((PetscObject)pc,PCSHELL,&flg);
93: if (flg) shell->ctx = ctx;
94: return(0);
95: }
97: static PetscErrorCode PCSetUp_Shell(PC pc)
98: {
99: PC_Shell *shell = (PC_Shell*)pc->data;
103: if (!shell->setup) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No setup() routine provided to Shell PC");
104: PetscStackCall("PCSHELL user function setup()",(*shell->setup)(pc);CHKERRQ(ierr));
105: return(0);
106: }
108: static PetscErrorCode PCApply_Shell(PC pc,Vec x,Vec y)
109: {
110: PC_Shell *shell = (PC_Shell*)pc->data;
111: PetscErrorCode ierr;
112: PetscObjectState instate,outstate;
115: if (!shell->apply) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No apply() routine provided to Shell PC");
116: PetscObjectStateGet((PetscObject)y, &instate);
117: PetscStackCall("PCSHELL user function apply()",(*shell->apply)(pc,x,y);CHKERRQ(ierr));
118: PetscObjectStateGet((PetscObject)y, &outstate);
119: if (instate == outstate) {
120: /* increase the state of the output vector since the user did not update its state themselve as should have been done */
121: PetscObjectStateIncrease((PetscObject)y);
122: }
123: return(0);
124: }
126: static PetscErrorCode PCApplySymmetricLeft_Shell(PC pc,Vec x,Vec y)
127: {
128: PC_Shell *shell = (PC_Shell*)pc->data;
132: if (!shell->applysymmetricleft) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No apply() routine provided to Shell PC");
133: PetscStackCall("PCSHELL user function apply()",(*shell->applysymmetricleft)(pc,x,y);CHKERRQ(ierr));
134: return(0);
135: }
137: static PetscErrorCode PCApplySymmetricRight_Shell(PC pc,Vec x,Vec y)
138: {
139: PC_Shell *shell = (PC_Shell*)pc->data;
143: if (!shell->applysymmetricright) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No apply() routine provided to Shell PC");
144: PetscStackCall("PCSHELL user function apply()",(*shell->applysymmetricright)(pc,x,y);CHKERRQ(ierr));
145: return(0);
146: }
148: static PetscErrorCode PCApplyBA_Shell(PC pc,PCSide side,Vec x,Vec y,Vec w)
149: {
150: PC_Shell *shell = (PC_Shell*)pc->data;
151: PetscErrorCode ierr;
152: PetscObjectState instate,outstate;
155: if (!shell->applyBA) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applyBA() routine provided to Shell PC");
156: PetscObjectStateGet((PetscObject)w, &instate);
157: PetscStackCall("PCSHELL user function applyBA()",(*shell->applyBA)(pc,side,x,y,w);CHKERRQ(ierr));
158: PetscObjectStateGet((PetscObject)w, &outstate);
159: if (instate == outstate) {
160: /* increase the state of the output vector since the user did not update its state themselve as should have been done */
161: PetscObjectStateIncrease((PetscObject)w);
162: }
163: return(0);
164: }
166: static PetscErrorCode PCPreSolveChangeRHS_Shell(PC pc,PetscBool* change)
167: {
169: *change = PETSC_TRUE;
170: return(0);
171: }
173: static PetscErrorCode PCPreSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
174: {
175: PC_Shell *shell = (PC_Shell*)pc->data;
179: if (!shell->presolve) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No presolve() routine provided to Shell PC");
180: PetscStackCall("PCSHELL user function presolve()",(*shell->presolve)(pc,ksp,b,x);CHKERRQ(ierr));
181: return(0);
182: }
184: static PetscErrorCode PCPostSolve_Shell(PC pc,KSP ksp,Vec b,Vec x)
185: {
186: PC_Shell *shell = (PC_Shell*)pc->data;
190: if (!shell->postsolve) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No postsolve() routine provided to Shell PC");
191: PetscStackCall("PCSHELL user function postsolve()",(*shell->postsolve)(pc,ksp,b,x);CHKERRQ(ierr));
192: return(0);
193: }
195: static PetscErrorCode PCApplyTranspose_Shell(PC pc,Vec x,Vec y)
196: {
197: PC_Shell *shell = (PC_Shell*)pc->data;
198: PetscErrorCode ierr;
199: PetscObjectState instate,outstate;
202: if (!shell->applytranspose) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applytranspose() routine provided to Shell PC");
203: PetscObjectStateGet((PetscObject)y, &instate);
204: PetscStackCall("PCSHELL user function applytranspose()",(*shell->applytranspose)(pc,x,y);CHKERRQ(ierr));
205: PetscObjectStateGet((PetscObject)y, &outstate);
206: if (instate == outstate) {
207: /* increase the state of the output vector since the user did not update its state themself as should have been done */
208: PetscObjectStateIncrease((PetscObject)y);
209: }
210: return(0);
211: }
213: 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)
214: {
215: PetscErrorCode ierr;
216: PC_Shell *shell = (PC_Shell*)pc->data;
217: PetscObjectState instate,outstate;
220: if (!shell->applyrich) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"No applyrichardson() routine provided to Shell PC");
221: PetscObjectStateGet((PetscObject)y, &instate);
222: PetscStackCall("PCSHELL user function applyrichardson()",(*shell->applyrich)(pc,x,y,w,rtol,abstol,dtol,it,guesszero,outits,reason);CHKERRQ(ierr));
223: PetscObjectStateGet((PetscObject)y, &outstate);
224: if (instate == outstate) {
225: /* increase the state of the output vector since the user did not update its state themself as should have been done */
226: PetscObjectStateIncrease((PetscObject)y);
227: }
228: return(0);
229: }
231: static PetscErrorCode PCDestroy_Shell(PC pc)
232: {
233: PC_Shell *shell = (PC_Shell*)pc->data;
237: PetscFree(shell->name);
238: if (shell->destroy) PetscStackCall("PCSHELL user function destroy()",(*shell->destroy)(pc);CHKERRQ(ierr));
239: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetDestroy_C",NULL);
240: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetSetUp_C",NULL);
241: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApply_C",NULL);
242: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricLeft_C",NULL);
243: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricRight_C",NULL);
244: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyBA_C",NULL);
245: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPreSolve_C",NULL);
246: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPostSolve_C",NULL);
247: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetView_C",NULL);
248: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",NULL);
249: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetName_C",NULL);
250: PetscObjectComposeFunction((PetscObject)pc,"PCShellGetName_C",NULL);
251: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",NULL);
252: PetscObjectComposeFunction((PetscObject)pc,"PCPreSolveChangeRHS_C",NULL);
253: PetscFree(pc->data);
254: return(0);
255: }
257: static PetscErrorCode PCView_Shell(PC pc,PetscViewer viewer)
258: {
259: PC_Shell *shell = (PC_Shell*)pc->data;
261: PetscBool iascii;
264: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
265: if (iascii) {
266: if (shell->name) {
267: PetscViewerASCIIPrintf(viewer," %s\n",shell->name);
268: } else {
269: PetscViewerASCIIPrintf(viewer," no name\n");
270: }
271: }
272: if (shell->view) {
273: PetscViewerASCIIPushTab(viewer);
274: (*shell->view)(pc,viewer);
275: PetscViewerASCIIPopTab(viewer);
276: }
277: return(0);
278: }
280: /* ------------------------------------------------------------------------------*/
281: static PetscErrorCode PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
282: {
283: PC_Shell *shell= (PC_Shell*)pc->data;
286: shell->destroy = destroy;
287: return(0);
288: }
290: static PetscErrorCode PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
291: {
292: PC_Shell *shell = (PC_Shell*)pc->data;;
295: shell->setup = setup;
296: if (setup) pc->ops->setup = PCSetUp_Shell;
297: else pc->ops->setup = 0;
298: return(0);
299: }
301: static PetscErrorCode PCShellSetApply_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
302: {
303: PC_Shell *shell = (PC_Shell*)pc->data;
306: shell->apply = apply;
307: return(0);
308: }
310: static PetscErrorCode PCShellSetApplySymmetricLeft_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
311: {
312: PC_Shell *shell = (PC_Shell*)pc->data;
315: shell->applysymmetricleft = apply;
316: return(0);
317: }
319: static PetscErrorCode PCShellSetApplySymmetricRight_Shell(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
320: {
321: PC_Shell *shell = (PC_Shell*)pc->data;
324: shell->applysymmetricright = apply;
325: return(0);
326: }
328: static PetscErrorCode PCShellSetApplyBA_Shell(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
329: {
330: PC_Shell *shell = (PC_Shell*)pc->data;
333: shell->applyBA = applyBA;
334: if (applyBA) pc->ops->applyBA = PCApplyBA_Shell;
335: else pc->ops->applyBA = 0;
336: return(0);
337: }
339: static PetscErrorCode PCShellSetPreSolve_Shell(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
340: {
341: PC_Shell *shell = (PC_Shell*)pc->data;
345: shell->presolve = presolve;
346: if (presolve) {
347: pc->ops->presolve = PCPreSolve_Shell;
348: PetscObjectComposeFunction((PetscObject)pc,"PCPreSolveChangeRHS_C",PCPreSolveChangeRHS_Shell);
349: } else {
350: pc->ops->presolve = 0;
351: PetscObjectComposeFunction((PetscObject)pc,"PCPreSolveChangeRHS_C",NULL);
352: }
353: return(0);
354: }
356: static PetscErrorCode PCShellSetPostSolve_Shell(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
357: {
358: PC_Shell *shell = (PC_Shell*)pc->data;
361: shell->postsolve = postsolve;
362: if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
363: else pc->ops->postsolve = 0;
364: return(0);
365: }
367: static PetscErrorCode PCShellSetView_Shell(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
368: {
369: PC_Shell *shell = (PC_Shell*)pc->data;
372: shell->view = view;
373: return(0);
374: }
376: static PetscErrorCode PCShellSetApplyTranspose_Shell(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
377: {
378: PC_Shell *shell = (PC_Shell*)pc->data;
381: shell->applytranspose = applytranspose;
382: if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
383: else pc->ops->applytranspose = 0;
384: return(0);
385: }
387: static PetscErrorCode PCShellSetApplyRichardson_Shell(PC pc,PetscErrorCode (*applyrich)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*))
388: {
389: PC_Shell *shell = (PC_Shell*)pc->data;
392: shell->applyrich = applyrich;
393: if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
394: else pc->ops->applyrichardson = 0;
395: return(0);
396: }
398: static PetscErrorCode PCShellSetName_Shell(PC pc,const char name[])
399: {
400: PC_Shell *shell = (PC_Shell*)pc->data;
404: PetscFree(shell->name);
405: PetscStrallocpy(name,&shell->name);
406: return(0);
407: }
409: static PetscErrorCode PCShellGetName_Shell(PC pc,const char *name[])
410: {
411: PC_Shell *shell = (PC_Shell*)pc->data;
414: *name = shell->name;
415: return(0);
416: }
418: /* -------------------------------------------------------------------------------*/
420: /*@C
421: PCShellSetDestroy - Sets routine to use to destroy the user-provided
422: application context.
424: Logically Collective on PC
426: Input Parameters:
427: + pc - the preconditioner context
428: . destroy - the application-provided destroy routine
430: Calling sequence of destroy:
431: .vb
432: PetscErrorCode destroy (PC)
433: .ve
435: . ptr - the application context
437: Notes:
438: the function MUST return an error code of 0 on success and nonzero on failure.
440: Level: developer
442: .keywords: PC, shell, set, destroy, user-provided
444: .seealso: PCShellSetApply(), PCShellSetContext()
445: @*/
446: PetscErrorCode PCShellSetDestroy(PC pc,PetscErrorCode (*destroy)(PC))
447: {
452: PetscTryMethod(pc,"PCShellSetDestroy_C",(PC,PetscErrorCode (*)(PC)),(pc,destroy));
453: return(0);
454: }
457: /*@C
458: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
459: matrix operator is changed.
461: Logically Collective on PC
463: Input Parameters:
464: + pc - the preconditioner context
465: . setup - the application-provided setup routine
467: Calling sequence of setup:
468: .vb
469: PetscErrorCode setup (PC pc)
470: .ve
472: . pc - the preconditioner, get the application context with PCShellGetContext()
474: Notes:
475: the function MUST return an error code of 0 on success and nonzero on failure.
477: Level: developer
479: .keywords: PC, shell, set, setup, user-provided
481: .seealso: PCShellSetApplyRichardson(), PCShellSetApply(), PCShellSetContext()
482: @*/
483: PetscErrorCode PCShellSetSetUp(PC pc,PetscErrorCode (*setup)(PC))
484: {
489: PetscTryMethod(pc,"PCShellSetSetUp_C",(PC,PetscErrorCode (*)(PC)),(pc,setup));
490: return(0);
491: }
494: /*@C
495: PCShellSetView - Sets routine to use as viewer of shell preconditioner
497: Logically Collective on PC
499: Input Parameters:
500: + pc - the preconditioner context
501: - view - the application-provided view routine
503: Calling sequence of apply:
504: .vb
505: PetscErrorCode view(PC pc,PetscViewer v)
506: .ve
508: + pc - the preconditioner, get the application context with PCShellGetContext()
509: - v - viewer
511: Notes:
512: the function MUST return an error code of 0 on success and nonzero on failure.
514: Level: developer
516: .keywords: PC, shell, set, apply, user-provided
518: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose()
519: @*/
520: PetscErrorCode PCShellSetView(PC pc,PetscErrorCode (*view)(PC,PetscViewer))
521: {
526: PetscTryMethod(pc,"PCShellSetView_C",(PC,PetscErrorCode (*)(PC,PetscViewer)),(pc,view));
527: return(0);
528: }
530: /*@C
531: PCShellSetApply - Sets routine to use as preconditioner.
533: Logically Collective on PC
535: Input Parameters:
536: + pc - the preconditioner context
537: - apply - the application-provided preconditioning routine
539: Calling sequence of apply:
540: .vb
541: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
542: .ve
544: + pc - the preconditioner, get the application context with PCShellGetContext()
545: . xin - input vector
546: - xout - output vector
548: Notes:
549: the function MUST return an error code of 0 on success and nonzero on failure.
551: Level: developer
553: .keywords: PC, shell, set, apply, user-provided
555: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApplyBA(), PCShellSetApplySymmetricRight(),PCShellSetApplySymmetricLeft()
556: @*/
557: PetscErrorCode PCShellSetApply(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
558: {
563: PetscTryMethod(pc,"PCShellSetApply_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
564: return(0);
565: }
567: /*@C
568: PCShellSetApplySymmetricLeft - Sets routine to use as left preconditioner (when the PC_SYMMETRIC is used).
570: Logically Collective on PC
572: Input Parameters:
573: + pc - the preconditioner context
574: - apply - the application-provided left preconditioning routine
576: Calling sequence of apply:
577: .vb
578: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
579: .ve
581: + pc - the preconditioner, get the application context with PCShellGetContext()
582: . xin - input vector
583: - xout - output vector
585: Notes:
586: the function MUST return an error code of 0 on success and nonzero on failure.
588: Level: developer
590: .keywords: PC, shell, set, apply, user-provided
592: .seealso: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
593: @*/
594: PetscErrorCode PCShellSetApplySymmetricLeft(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
595: {
600: PetscTryMethod(pc,"PCShellSetApplySymmetricLeft_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
601: return(0);
602: }
604: /*@C
605: PCShellSetApplySymmetricRight - Sets routine to use as right preconditioner (when the PC_SYMMETRIC is used).
607: Logically Collective on PC
609: Input Parameters:
610: + pc - the preconditioner context
611: - apply - the application-provided right preconditioning routine
613: Calling sequence of apply:
614: .vb
615: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
616: .ve
618: + pc - the preconditioner, get the application context with PCShellGetContext()
619: . xin - input vector
620: - xout - output vector
622: Notes:
623: the function MUST return an error code of 0 on success and nonzero on failure.
625: Level: developer
627: .keywords: PC, shell, set, apply, user-provided
629: .seealso: PCShellSetApply(), PCShellSetApplySymmetricLeft(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext()
630: @*/
631: PetscErrorCode PCShellSetApplySymmetricRight(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec))
632: {
637: PetscTryMethod(pc,"PCShellSetApplySymmetricRight_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,apply));
638: return(0);
639: }
641: /*@C
642: PCShellSetApplyBA - Sets routine to use as preconditioner times operator.
644: Logically Collective on PC
646: Input Parameters:
647: + pc - the preconditioner context
648: - applyBA - the application-provided BA routine
650: Calling sequence of apply:
651: .vb
652: PetscErrorCode applyBA (PC pc,Vec xin,Vec xout)
653: .ve
655: + pc - the preconditioner, get the application context with PCShellGetContext()
656: . xin - input vector
657: - xout - output vector
659: Notes:
660: the function MUST return an error code of 0 on success and nonzero on failure.
662: Level: developer
664: .keywords: PC, shell, set, apply, user-provided
666: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetContext(), PCShellSetApply()
667: @*/
668: PetscErrorCode PCShellSetApplyBA(PC pc,PetscErrorCode (*applyBA)(PC,PCSide,Vec,Vec,Vec))
669: {
674: PetscTryMethod(pc,"PCShellSetApplyBA_C",(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec)),(pc,applyBA));
675: return(0);
676: }
678: /*@C
679: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
681: Logically Collective on PC
683: Input Parameters:
684: + pc - the preconditioner context
685: - apply - the application-provided preconditioning transpose routine
687: Calling sequence of apply:
688: .vb
689: PetscErrorCode applytranspose (PC pc,Vec xin,Vec xout)
690: .ve
692: + pc - the preconditioner, get the application context with PCShellGetContext()
693: . xin - input vector
694: - xout - output vector
696: Notes:
697: the function MUST return an error code of 0 on success and nonzero on failure.
699: Level: developer
701: Notes:
702: Uses the same context variable as PCShellSetApply().
704: .keywords: PC, shell, set, apply, user-provided
706: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApply(), PCSetContext(), PCShellSetApplyBA()
707: @*/
708: PetscErrorCode PCShellSetApplyTranspose(PC pc,PetscErrorCode (*applytranspose)(PC,Vec,Vec))
709: {
714: PetscTryMethod(pc,"PCShellSetApplyTranspose_C",(PC,PetscErrorCode (*)(PC,Vec,Vec)),(pc,applytranspose));
715: return(0);
716: }
718: /*@C
719: PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
720: applied. This usually does something like scale the linear system in some application
721: specific way.
723: Logically Collective on PC
725: Input Parameters:
726: + pc - the preconditioner context
727: - presolve - the application-provided presolve routine
729: Calling sequence of presolve:
730: .vb
731: PetscErrorCode presolve (PC,KSP ksp,Vec b,Vec x)
732: .ve
734: + pc - the preconditioner, get the application context with PCShellGetContext()
735: . xin - input vector
736: - xout - output vector
738: Notes:
739: the function MUST return an error code of 0 on success and nonzero on failure.
741: Level: developer
743: .keywords: PC, shell, set, apply, user-provided
745: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPostSolve(), PCShellSetContext()
746: @*/
747: PetscErrorCode PCShellSetPreSolve(PC pc,PetscErrorCode (*presolve)(PC,KSP,Vec,Vec))
748: {
753: PetscTryMethod(pc,"PCShellSetPreSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,presolve));
754: return(0);
755: }
757: /*@C
758: PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a KSPSolve() is
759: applied. This usually does something like scale the linear system in some application
760: specific way.
762: Logically Collective on PC
764: Input Parameters:
765: + pc - the preconditioner context
766: - postsolve - the application-provided presolve routine
768: Calling sequence of postsolve:
769: .vb
770: PetscErrorCode postsolve(PC,KSP ksp,Vec b,Vec x)
771: .ve
773: + pc - the preconditioner, get the application context with PCShellGetContext()
774: . xin - input vector
775: - xout - output vector
777: Notes:
778: the function MUST return an error code of 0 on success and nonzero on failure.
780: Level: developer
782: .keywords: PC, shell, set, apply, user-provided
784: .seealso: PCShellSetApplyRichardson(), PCShellSetSetUp(), PCShellSetApplyTranspose(), PCShellSetPreSolve(), PCShellSetContext()
785: @*/
786: PetscErrorCode PCShellSetPostSolve(PC pc,PetscErrorCode (*postsolve)(PC,KSP,Vec,Vec))
787: {
792: PetscTryMethod(pc,"PCShellSetPostSolve_C",(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec)),(pc,postsolve));
793: return(0);
794: }
796: /*@C
797: PCShellSetName - Sets an optional name to associate with a shell
798: preconditioner.
800: Not Collective
802: Input Parameters:
803: + pc - the preconditioner context
804: - name - character string describing shell preconditioner
806: Level: developer
808: .keywords: PC, shell, set, name, user-provided
810: .seealso: PCShellGetName()
811: @*/
812: PetscErrorCode PCShellSetName(PC pc,const char name[])
813: {
818: PetscTryMethod(pc,"PCShellSetName_C",(PC,const char []),(pc,name));
819: return(0);
820: }
822: /*@C
823: PCShellGetName - Gets an optional name that the user has set for a shell
824: preconditioner.
826: Not Collective
828: Input Parameter:
829: . pc - the preconditioner context
831: Output Parameter:
832: . name - character string describing shell preconditioner (you should not free this)
834: Level: developer
836: .keywords: PC, shell, get, name, user-provided
838: .seealso: PCShellSetName()
839: @*/
840: PetscErrorCode PCShellGetName(PC pc,const char *name[])
841: {
847: PetscUseMethod(pc,"PCShellGetName_C",(PC,const char*[]),(pc,name));
848: return(0);
849: }
851: /*@C
852: PCShellSetApplyRichardson - Sets routine to use as preconditioner
853: in Richardson iteration.
855: Logically Collective on PC
857: Input Parameters:
858: + pc - the preconditioner context
859: - apply - the application-provided preconditioning routine
861: Calling sequence of apply:
862: .vb
863: PetscErrorCode apply (PC pc,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
864: .ve
866: + pc - the preconditioner, get the application context with PCShellGetContext()
867: . b - right-hand-side
868: . x - current iterate
869: . r - work space
870: . rtol - relative tolerance of residual norm to stop at
871: . abstol - absolute tolerance of residual norm to stop at
872: . dtol - if residual norm increases by this factor than return
873: - maxits - number of iterations to run
875: Notes:
876: the function MUST return an error code of 0 on success and nonzero on failure.
878: Level: developer
880: .keywords: PC, shell, set, apply, Richardson, user-provided
882: .seealso: PCShellSetApply(), PCShellSetContext()
883: @*/
884: PetscErrorCode PCShellSetApplyRichardson(PC pc,PetscErrorCode (*apply)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*))
885: {
890: PetscTryMethod(pc,"PCShellSetApplyRichardson_C",(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool,PetscInt*,PCRichardsonConvergedReason*)),(pc,apply));
891: return(0);
892: }
894: /*MC
895: PCSHELL - Creates a new preconditioner class for use with your
896: own private data storage format.
898: Level: advanced
899: >
900: Concepts: providing your own preconditioner
902: Usage:
903: $ extern PetscErrorCode apply(PC,Vec,Vec);
904: $ extern PetscErrorCode applyba(PC,PCSide,Vec,Vec,Vec);
905: $ extern PetscErrorCode applytranspose(PC,Vec,Vec);
906: $ extern PetscErrorCode setup(PC);
907: $ extern PetscErrorCode destroy(PC);
908: $
909: $ PCCreate(comm,&pc);
910: $ PCSetType(pc,PCSHELL);
911: $ PCShellSetContext(pc,ctx)
912: $ PCShellSetApply(pc,apply);
913: $ PCShellSetApplyBA(pc,applyba); (optional)
914: $ PCShellSetApplyTranspose(pc,applytranspose); (optional)
915: $ PCShellSetSetUp(pc,setup); (optional)
916: $ PCShellSetDestroy(pc,destroy); (optional)
918: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
919: MATSHELL, PCShellSetSetUp(), PCShellSetApply(), PCShellSetView(),
920: PCShellSetApplyTranspose(), PCShellSetName(), PCShellSetApplyRichardson(),
921: PCShellGetName(), PCShellSetContext(), PCShellGetContext(), PCShellSetApplyBA()
922: M*/
924: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
925: {
927: PC_Shell *shell;
930: PetscNewLog(pc,&shell);
931: pc->data = (void*)shell;
933: pc->ops->destroy = PCDestroy_Shell;
934: pc->ops->view = PCView_Shell;
935: pc->ops->apply = PCApply_Shell;
936: pc->ops->applysymmetricleft = PCApplySymmetricLeft_Shell;
937: pc->ops->applysymmetricright = PCApplySymmetricRight_Shell;
938: pc->ops->applytranspose = 0;
939: pc->ops->applyrichardson = 0;
940: pc->ops->setup = 0;
941: pc->ops->presolve = 0;
942: pc->ops->postsolve = 0;
944: shell->apply = 0;
945: shell->applytranspose = 0;
946: shell->name = 0;
947: shell->applyrich = 0;
948: shell->presolve = 0;
949: shell->postsolve = 0;
950: shell->ctx = 0;
951: shell->setup = 0;
952: shell->view = 0;
953: shell->destroy = 0;
954: shell->applysymmetricleft = 0;
955: shell->applysymmetricright = 0;
957: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetDestroy_C",PCShellSetDestroy_Shell);
958: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetSetUp_C",PCShellSetSetUp_Shell);
959: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApply_C",PCShellSetApply_Shell);
960: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricLeft_C",PCShellSetApplySymmetricLeft_Shell);
961: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplySymmetricRight_C",PCShellSetApplySymmetricRight_Shell);
962: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyBA_C",PCShellSetApplyBA_Shell);
963: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPreSolve_C",PCShellSetPreSolve_Shell);
964: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetPostSolve_C",PCShellSetPostSolve_Shell);
965: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetView_C",PCShellSetView_Shell);
966: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyTranspose_C",PCShellSetApplyTranspose_Shell);
967: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetName_C",PCShellSetName_Shell);
968: PetscObjectComposeFunction((PetscObject)pc,"PCShellGetName_C",PCShellGetName_Shell);
969: PetscObjectComposeFunction((PetscObject)pc,"PCShellSetApplyRichardson_C",PCShellSetApplyRichardson_Shell);
970: return(0);
971: }