Actual source code: shellpc.c
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 (*matapply)(PC, Mat, Mat);
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;
28: /*@C
29: PCShellGetContext - Returns the user-provided context associated with a shell `PC`
31: Not Collective
33: Input Parameter:
34: . pc - of type `PCSHELL` created with `PCSetType`(pc,shell)
36: Output Parameter:
37: . ctx - the user provided context
39: Level: advanced
41: Note:
42: This routine is intended for use within various shell routines
44: Fortran Note:
45: To use this from Fortran you must write a Fortran interface definition for this
46: function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.
48: .seealso: `PCSHELL`, `PCShellSetContext()`
49: @*/
50: PetscErrorCode PCShellGetContext(PC pc, void *ctx)
51: {
52: PetscBool flg;
56: PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg);
57: if (!flg) *(void **)ctx = NULL;
58: else *(void **)ctx = ((PC_Shell *)(pc->data))->ctx;
59: return 0;
60: }
62: /*@
63: PCShellSetContext - sets the context for a shell `PC`
65: Logically Collective
67: Input Parameters:
68: + pc - the `PC` of type `PCSHELL`
69: - ctx - the context
71: Level: advanced
73: Fortran Note:
74: To use this from Fortran you must write a Fortran interface definition for this
75: function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.
77: .seealso: `PCShellGetContext()`, `PCSHELL`
78: @*/
79: PetscErrorCode PCShellSetContext(PC pc, void *ctx)
80: {
81: PC_Shell *shell = (PC_Shell *)pc->data;
82: PetscBool flg;
85: PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg);
86: if (flg) shell->ctx = ctx;
87: return 0;
88: }
90: static PetscErrorCode PCSetUp_Shell(PC pc)
91: {
92: PC_Shell *shell = (PC_Shell *)pc->data;
95: PetscCallBack("PCSHELL callback setup", (*shell->setup)(pc));
96: return 0;
97: }
99: static PetscErrorCode PCApply_Shell(PC pc, Vec x, Vec y)
100: {
101: PC_Shell *shell = (PC_Shell *)pc->data;
102: PetscObjectState instate, outstate;
105: PetscObjectStateGet((PetscObject)y, &instate);
106: PetscCallBack("PCSHELL callback apply", (*shell->apply)(pc, x, y));
107: PetscObjectStateGet((PetscObject)y, &outstate);
108: /* increase the state of the output vector if the user did not update its state themself as should have been done */
109: if (instate == outstate) PetscObjectStateIncrease((PetscObject)y);
110: return 0;
111: }
113: static PetscErrorCode PCMatApply_Shell(PC pc, Mat X, Mat Y)
114: {
115: PC_Shell *shell = (PC_Shell *)pc->data;
116: PetscObjectState instate, outstate;
119: PetscObjectStateGet((PetscObject)Y, &instate);
120: PetscCallBack("PCSHELL callback apply", (*shell->matapply)(pc, X, Y));
121: PetscObjectStateGet((PetscObject)Y, &outstate);
122: /* increase the state of the output vector if the user did not update its state themself as should have been done */
123: if (instate == outstate) PetscObjectStateIncrease((PetscObject)Y);
124: return 0;
125: }
127: static PetscErrorCode PCApplySymmetricLeft_Shell(PC pc, Vec x, Vec y)
128: {
129: PC_Shell *shell = (PC_Shell *)pc->data;
132: PetscCallBack("PCSHELL callback apply symmetric left", (*shell->applysymmetricleft)(pc, x, y));
133: return 0;
134: }
136: static PetscErrorCode PCApplySymmetricRight_Shell(PC pc, Vec x, Vec y)
137: {
138: PC_Shell *shell = (PC_Shell *)pc->data;
141: PetscCallBack("PCSHELL callback apply symmetric right", (*shell->applysymmetricright)(pc, x, y));
142: return 0;
143: }
145: static PetscErrorCode PCApplyBA_Shell(PC pc, PCSide side, Vec x, Vec y, Vec w)
146: {
147: PC_Shell *shell = (PC_Shell *)pc->data;
148: PetscObjectState instate, outstate;
151: PetscObjectStateGet((PetscObject)w, &instate);
152: PetscCallBack("PCSHELL callback applyBA", (*shell->applyBA)(pc, side, x, y, w));
153: PetscObjectStateGet((PetscObject)w, &outstate);
154: /* increase the state of the output vector if the user did not update its state themself as should have been done */
155: if (instate == outstate) PetscObjectStateIncrease((PetscObject)w);
156: return 0;
157: }
159: static PetscErrorCode PCPreSolveChangeRHS_Shell(PC pc, PetscBool *change)
160: {
161: *change = PETSC_TRUE;
162: return 0;
163: }
165: static PetscErrorCode PCPreSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
166: {
167: PC_Shell *shell = (PC_Shell *)pc->data;
170: PetscCallBack("PCSHELL callback presolve", (*shell->presolve)(pc, ksp, b, x));
171: return 0;
172: }
174: static PetscErrorCode PCPostSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
175: {
176: PC_Shell *shell = (PC_Shell *)pc->data;
179: PetscCallBack("PCSHELL callback postsolve()", (*shell->postsolve)(pc, ksp, b, x));
180: return 0;
181: }
183: static PetscErrorCode PCApplyTranspose_Shell(PC pc, Vec x, Vec y)
184: {
185: PC_Shell *shell = (PC_Shell *)pc->data;
186: PetscObjectState instate, outstate;
189: PetscObjectStateGet((PetscObject)y, &instate);
190: PetscCallBack("PCSHELL callback applytranspose", (*shell->applytranspose)(pc, x, y));
191: PetscObjectStateGet((PetscObject)y, &outstate);
192: /* increase the state of the output vector if the user did not update its state themself as should have been done */
193: if (instate == outstate) PetscObjectStateIncrease((PetscObject)y);
194: return 0;
195: }
197: 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)
198: {
199: PC_Shell *shell = (PC_Shell *)pc->data;
200: PetscObjectState instate, outstate;
203: PetscObjectStateGet((PetscObject)y, &instate);
204: PetscCallBack("PCSHELL callback applyrichardson", (*shell->applyrich)(pc, x, y, w, rtol, abstol, dtol, it, guesszero, outits, reason));
205: PetscObjectStateGet((PetscObject)y, &outstate);
206: /* increase the state of the output vector since the user did not update its state themself as should have been done */
207: if (instate == outstate) PetscObjectStateIncrease((PetscObject)y);
208: return 0;
209: }
211: static PetscErrorCode PCDestroy_Shell(PC pc)
212: {
213: PC_Shell *shell = (PC_Shell *)pc->data;
215: PetscFree(shell->name);
216: if (shell->destroy) PetscCallBack("PCSHELL callback destroy", (*shell->destroy)(pc));
217: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", NULL);
218: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", NULL);
219: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", NULL);
220: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", NULL);
221: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", NULL);
222: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", NULL);
223: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", NULL);
224: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", NULL);
225: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", NULL);
226: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", NULL);
227: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", NULL);
228: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", NULL);
229: PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", NULL);
230: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", NULL);
231: PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL);
232: PetscFree(pc->data);
233: return 0;
234: }
236: static PetscErrorCode PCView_Shell(PC pc, PetscViewer viewer)
237: {
238: PC_Shell *shell = (PC_Shell *)pc->data;
239: PetscBool iascii;
241: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
242: if (iascii) {
243: if (shell->name) PetscViewerASCIIPrintf(viewer, " %s\n", shell->name);
244: else PetscViewerASCIIPrintf(viewer, " no name\n");
245: }
246: if (shell->view) {
247: PetscViewerASCIIPushTab(viewer);
248: (*shell->view)(pc, viewer);
249: PetscViewerASCIIPopTab(viewer);
250: }
251: return 0;
252: }
254: static PetscErrorCode PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
255: {
256: PC_Shell *shell = (PC_Shell *)pc->data;
258: shell->destroy = destroy;
259: return 0;
260: }
262: static PetscErrorCode PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
263: {
264: PC_Shell *shell = (PC_Shell *)pc->data;
266: shell->setup = setup;
267: if (setup) pc->ops->setup = PCSetUp_Shell;
268: else pc->ops->setup = NULL;
269: return 0;
270: }
272: static PetscErrorCode PCShellSetApply_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
273: {
274: PC_Shell *shell = (PC_Shell *)pc->data;
276: shell->apply = apply;
277: return 0;
278: }
280: static PetscErrorCode PCShellSetMatApply_Shell(PC pc, PetscErrorCode (*matapply)(PC, Mat, Mat))
281: {
282: PC_Shell *shell = (PC_Shell *)pc->data;
284: shell->matapply = matapply;
285: if (matapply) pc->ops->matapply = PCMatApply_Shell;
286: else pc->ops->matapply = NULL;
287: return 0;
288: }
290: static PetscErrorCode PCShellSetApplySymmetricLeft_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
291: {
292: PC_Shell *shell = (PC_Shell *)pc->data;
294: shell->applysymmetricleft = apply;
295: return 0;
296: }
298: static PetscErrorCode PCShellSetApplySymmetricRight_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
299: {
300: PC_Shell *shell = (PC_Shell *)pc->data;
302: shell->applysymmetricright = apply;
303: return 0;
304: }
306: static PetscErrorCode PCShellSetApplyBA_Shell(PC pc, PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec))
307: {
308: PC_Shell *shell = (PC_Shell *)pc->data;
310: shell->applyBA = applyBA;
311: if (applyBA) pc->ops->applyBA = PCApplyBA_Shell;
312: else pc->ops->applyBA = NULL;
313: return 0;
314: }
316: static PetscErrorCode PCShellSetPreSolve_Shell(PC pc, PetscErrorCode (*presolve)(PC, KSP, Vec, Vec))
317: {
318: PC_Shell *shell = (PC_Shell *)pc->data;
320: shell->presolve = presolve;
321: if (presolve) {
322: pc->ops->presolve = PCPreSolve_Shell;
323: PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", PCPreSolveChangeRHS_Shell);
324: } else {
325: pc->ops->presolve = NULL;
326: PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL);
327: }
328: return 0;
329: }
331: static PetscErrorCode PCShellSetPostSolve_Shell(PC pc, PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec))
332: {
333: PC_Shell *shell = (PC_Shell *)pc->data;
335: shell->postsolve = postsolve;
336: if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
337: else pc->ops->postsolve = NULL;
338: return 0;
339: }
341: static PetscErrorCode PCShellSetView_Shell(PC pc, PetscErrorCode (*view)(PC, PetscViewer))
342: {
343: PC_Shell *shell = (PC_Shell *)pc->data;
345: shell->view = view;
346: return 0;
347: }
349: static PetscErrorCode PCShellSetApplyTranspose_Shell(PC pc, PetscErrorCode (*applytranspose)(PC, Vec, Vec))
350: {
351: PC_Shell *shell = (PC_Shell *)pc->data;
353: shell->applytranspose = applytranspose;
354: if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
355: else pc->ops->applytranspose = NULL;
356: return 0;
357: }
359: static PetscErrorCode PCShellSetApplyRichardson_Shell(PC pc, PetscErrorCode (*applyrich)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *))
360: {
361: PC_Shell *shell = (PC_Shell *)pc->data;
363: shell->applyrich = applyrich;
364: if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
365: else pc->ops->applyrichardson = NULL;
366: return 0;
367: }
369: static PetscErrorCode PCShellSetName_Shell(PC pc, const char name[])
370: {
371: PC_Shell *shell = (PC_Shell *)pc->data;
373: PetscFree(shell->name);
374: PetscStrallocpy(name, &shell->name);
375: return 0;
376: }
378: static PetscErrorCode PCShellGetName_Shell(PC pc, const char *name[])
379: {
380: PC_Shell *shell = (PC_Shell *)pc->data;
382: *name = shell->name;
383: return 0;
384: }
386: /*@C
387: PCShellSetDestroy - Sets routine to use to destroy the user-provided
388: application context.
390: Logically Collective
392: Input Parameters:
393: + pc - the preconditioner context
394: - destroy - the application-provided destroy routine
396: Calling sequence of destroy:
397: .vb
398: PetscErrorCode destroy (PC)
399: .ve
401: . ptr - the application context
403: Note:
404: the function MUST return an error code of 0 on success and nonzero on failure.
406: Level: intermediate
408: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetContext()`
409: @*/
410: PetscErrorCode PCShellSetDestroy(PC pc, PetscErrorCode (*destroy)(PC))
411: {
413: PetscTryMethod(pc, "PCShellSetDestroy_C", (PC, PetscErrorCode(*)(PC)), (pc, destroy));
414: return 0;
415: }
417: /*@C
418: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
419: matrix operator is changed.
421: Logically Collective
423: Input Parameters:
424: + pc - the preconditioner context
425: - setup - the application-provided setup routine
427: Calling sequence of setup:
428: .vb
429: PetscErrorCode setup (PC pc)
430: .ve
432: . pc - the preconditioner, get the application context with PCShellGetContext()
434: Note:
435: the function MUST return an error code of 0 on success and nonzero on failure.
437: Level: intermediate
439: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetApply()`, `PCShellSetContext()`
440: @*/
441: PetscErrorCode PCShellSetSetUp(PC pc, PetscErrorCode (*setup)(PC))
442: {
444: PetscTryMethod(pc, "PCShellSetSetUp_C", (PC, PetscErrorCode(*)(PC)), (pc, setup));
445: return 0;
446: }
448: /*@C
449: PCShellSetView - Sets routine to use as viewer of shell preconditioner
451: Logically Collective
453: Input Parameters:
454: + pc - the preconditioner context
455: - view - the application-provided view routine
457: Calling sequence of view:
458: .vb
459: PetscErrorCode view(PC pc,PetscViewer v)
460: .ve
462: + pc - the preconditioner, get the application context with PCShellGetContext()
463: - v - viewer
465: Note:
466: the function MUST return an error code of 0 on success and nonzero on failure.
468: Level: advanced
470: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`
471: @*/
472: PetscErrorCode PCShellSetView(PC pc, PetscErrorCode (*view)(PC, PetscViewer))
473: {
475: PetscTryMethod(pc, "PCShellSetView_C", (PC, PetscErrorCode(*)(PC, PetscViewer)), (pc, view));
476: return 0;
477: }
479: /*@C
480: PCShellSetApply - Sets routine to use as preconditioner.
482: Logically Collective
484: Input Parameters:
485: + pc - the preconditioner context
486: - apply - the application-provided preconditioning routine
488: Calling sequence of apply:
489: .vb
490: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
491: .ve
493: + pc - the preconditioner, get the application context with PCShellGetContext()
494: . xin - input vector
495: - xout - output vector
497: Note:
498: the function MUST return an error code of 0 on success and nonzero on failure.
500: Level: intermediate
502: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApplyBA()`, `PCShellSetApplySymmetricRight()`, `PCShellSetApplySymmetricLeft()`
503: @*/
504: PetscErrorCode PCShellSetApply(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
505: {
507: PetscTryMethod(pc, "PCShellSetApply_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
508: return 0;
509: }
511: /*@C
512: PCShellSetMatApply - Sets routine to use as preconditioner on a block of vectors.
514: Logically Collective
516: Input Parameters:
517: + pc - the preconditioner context
518: - apply - the application-provided preconditioning routine
520: Calling sequence of apply:
521: .vb
522: PetscErrorCode apply (PC pc,Mat Xin,Mat Xout)
523: .ve
525: + pc - the preconditioner, get the application context with PCShellGetContext()
526: . Xin - input block of vectors
527: - Xout - output block of vectors
529: Note:
530: The function MUST return an error code of 0 on success and nonzero on failure.
532: Level: advanced
534: .seealso: `PCSHELL`, `PCShellSetApply()`
535: @*/
536: PetscErrorCode PCShellSetMatApply(PC pc, PetscErrorCode (*matapply)(PC, Mat, Mat))
537: {
539: PetscTryMethod(pc, "PCShellSetMatApply_C", (PC, PetscErrorCode(*)(PC, Mat, Mat)), (pc, matapply));
540: return 0;
541: }
543: /*@C
544: PCShellSetApplySymmetricLeft - Sets routine to use as left preconditioner (when the `PC_SYMMETRIC` is used).
546: Logically Collective
548: Input Parameters:
549: + pc - the preconditioner context
550: - apply - the application-provided left preconditioning routine
552: Calling sequence of apply:
553: .vb
554: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
555: .ve
557: + pc - the preconditioner, get the application context with `PCShellGetContext()`
558: . xin - input vector
559: - xout - output vector
561: Note:
562: The function MUST return an error code of 0 on success and nonzero on failure.
564: Level: advanced
566: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetApplySymmetricLeft()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`
567: @*/
568: PetscErrorCode PCShellSetApplySymmetricLeft(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
569: {
571: PetscTryMethod(pc, "PCShellSetApplySymmetricLeft_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
572: return 0;
573: }
575: /*@C
576: PCShellSetApplySymmetricRight - Sets routine to use as right preconditioner (when the `PC_SYMMETRIC` is used).
578: Logically Collective
580: Input Parameters:
581: + pc - the preconditioner context
582: - apply - the application-provided right preconditioning routine
584: Calling sequence of apply:
585: .vb
586: PetscErrorCode apply (PC pc,Vec xin,Vec xout)
587: .ve
589: + pc - the preconditioner, get the application context with PCShellGetContext()
590: . xin - input vector
591: - xout - output vector
593: Note:
594: The function MUST return an error code of 0 on success and nonzero on failure.
596: Level: advanced
598: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetApplySymmetricLeft()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`
599: @*/
600: PetscErrorCode PCShellSetApplySymmetricRight(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
601: {
603: PetscTryMethod(pc, "PCShellSetApplySymmetricRight_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
604: return 0;
605: }
607: /*@C
608: PCShellSetApplyBA - Sets routine to use as preconditioner times operator.
610: Logically Collective
612: Input Parameters:
613: + pc - the preconditioner context
614: - applyBA - the application-provided BA routine
616: Calling sequence of applyBA:
617: .vb
618: PetscErrorCode applyBA (PC pc,Vec xin,Vec xout)
619: .ve
621: + pc - the preconditioner, get the application context with PCShellGetContext()
622: . xin - input vector
623: - xout - output vector
625: Note:
626: The function MUST return an error code of 0 on success and nonzero on failure.
628: Level: intermediate
630: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApply()`
631: @*/
632: PetscErrorCode PCShellSetApplyBA(PC pc, PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec))
633: {
635: PetscTryMethod(pc, "PCShellSetApplyBA_C", (PC, PetscErrorCode(*)(PC, PCSide, Vec, Vec, Vec)), (pc, applyBA));
636: return 0;
637: }
639: /*@C
640: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
642: Logically Collective
644: Input Parameters:
645: + pc - the preconditioner context
646: - apply - the application-provided preconditioning transpose routine
648: Calling sequence of apply:
649: .vb
650: PetscErrorCode applytranspose (PC pc,Vec xin,Vec xout)
651: .ve
653: + pc - the preconditioner, get the application context with PCShellGetContext()
654: . xin - input vector
655: - xout - output vector
657: Note:
658: the function MUST return an error code of 0 on success and nonzero on failure.
660: Level: intermediate
662: Note:
663: Uses the same context variable as `PCShellSetApply()`.
665: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApply()`, `PCSetContext()`, `PCShellSetApplyBA()`
666: @*/
667: PetscErrorCode PCShellSetApplyTranspose(PC pc, PetscErrorCode (*applytranspose)(PC, Vec, Vec))
668: {
670: PetscTryMethod(pc, "PCShellSetApplyTranspose_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, applytranspose));
671: return 0;
672: }
674: /*@C
675: PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a `KSPSolve()` is
676: applied. This usually does something like scale the linear system in some application
677: specific way.
679: Logically Collective
681: Input Parameters:
682: + pc - the preconditioner context
683: - presolve - the application-provided presolve routine
685: Calling sequence of presolve:
686: .vb
687: PetscErrorCode presolve (PC,KSP ksp,Vec b,Vec x)
688: .ve
690: + pc - the preconditioner, get the application context with PCShellGetContext()
691: . xin - input vector
692: - xout - output vector
694: Note:
695: The function MUST return an error code of 0 on success and nonzero on failure.
697: Level: advanced
699: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPostSolve()`, `PCShellSetContext()`
700: @*/
701: PetscErrorCode PCShellSetPreSolve(PC pc, PetscErrorCode (*presolve)(PC, KSP, Vec, Vec))
702: {
704: PetscTryMethod(pc, "PCShellSetPreSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, presolve));
705: return 0;
706: }
708: /*@C
709: PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a `KSPSolve()` is
710: applied. This usually does something like scale the linear system in some application
711: specific way.
713: Logically Collective
715: Input Parameters:
716: + pc - the preconditioner context
717: - postsolve - the application-provided presolve routine
719: Calling sequence of postsolve:
720: .vb
721: PetscErrorCode postsolve(PC,KSP ksp,Vec b,Vec x)
722: .ve
724: + pc - the preconditioner, get the application context with `PCShellGetContext()`
725: . xin - input vector
726: - xout - output vector
728: Note:
729: the function MUST return an error code of 0 on success and nonzero on failure.
731: Level: advanced
733: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPreSolve()`, `PCShellSetContext()`
734: @*/
735: PetscErrorCode PCShellSetPostSolve(PC pc, PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec))
736: {
738: PetscTryMethod(pc, "PCShellSetPostSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, postsolve));
739: return 0;
740: }
742: /*@C
743: PCShellSetName - Sets an optional name to associate with a `PCSHELL`
744: preconditioner.
746: Not Collective
748: Input Parameters:
749: + pc - the preconditioner context
750: - name - character string describing shell preconditioner
752: Level: intermediate
754: .seealso: `PCSHELL`, `PCShellGetName()`
755: @*/
756: PetscErrorCode PCShellSetName(PC pc, const char name[])
757: {
759: PetscTryMethod(pc, "PCShellSetName_C", (PC, const char[]), (pc, name));
760: return 0;
761: }
763: /*@C
764: PCShellGetName - Gets an optional name that the user has set for a `PCSHELL`
765: preconditioner.
767: Not Collective
769: Input Parameter:
770: . pc - the preconditioner context
772: Output Parameter:
773: . name - character string describing shell preconditioner (you should not free this)
775: Level: intermediate
777: .seealso: `PCSHELL`, `PCShellSetName()`
778: @*/
779: PetscErrorCode PCShellGetName(PC pc, const char *name[])
780: {
783: PetscUseMethod(pc, "PCShellGetName_C", (PC, const char *[]), (pc, name));
784: return 0;
785: }
787: /*@C
788: PCShellSetApplyRichardson - Sets routine to use as preconditioner
789: in Richardson iteration.
791: Logically Collective
793: Input Parameters:
794: + pc - the preconditioner context
795: - apply - the application-provided preconditioning routine
797: Calling sequence of apply:
798: .vb
799: PetscErrorCode apply (PC pc,Vec b,Vec x,Vec r,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
800: .ve
802: + pc - the preconditioner, get the application context with PCShellGetContext()
803: . b - right-hand-side
804: . x - current iterate
805: . r - work space
806: . rtol - relative tolerance of residual norm to stop at
807: . abstol - absolute tolerance of residual norm to stop at
808: . dtol - if residual norm increases by this factor than return
809: - maxits - number of iterations to run
811: Note:
812: the function MUST return an error code of 0 on success and nonzero on failure.
814: Level: advanced
816: .seealso: `PCShellSetApply()`, `PCShellSetContext()`
817: @*/
818: PetscErrorCode PCShellSetApplyRichardson(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *))
819: {
821: PetscTryMethod(pc, "PCShellSetApplyRichardson_C", (PC, PetscErrorCode(*)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *)), (pc, apply));
822: return 0;
823: }
825: /*MC
826: PCSHELL - Creates a new preconditioner class for use with a users
827: own private data storage format and preconditioner application code
829: Level: advanced
831: Usage:
832: .vb
833: extern PetscErrorCode apply(PC,Vec,Vec);
834: extern PetscErrorCode applyba(PC,PCSide,Vec,Vec,Vec);
835: extern PetscErrorCode applytranspose(PC,Vec,Vec);
836: extern PetscErrorCode setup(PC);
837: extern PetscErrorCode destroy(PC);
839: PCCreate(comm,&pc);
840: PCSetType(pc,PCSHELL);
841: PCShellSetContext(pc,ctx)
842: PCShellSetApply(pc,apply);
843: PCShellSetApplyBA(pc,applyba); (optional)
844: PCShellSetApplyTranspose(pc,applytranspose); (optional)
845: PCShellSetSetUp(pc,setup); (optional)
846: PCShellSetDestroy(pc,destroy); (optional)
847: .ve
849: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
850: `MATSHELL`, `PCShellSetSetUp()`, `PCShellSetApply()`, `PCShellSetView()`, `PCShellSetDestroy()`, `PCShellSetPostSolve()`,
851: `PCShellSetApplyTranspose()`, `PCShellSetName()`, `PCShellSetApplyRichardson()`, `PCShellSetPreSolve()`, `PCShellSetView()`,
852: `PCShellGetName()`, `PCShellSetContext()`, `PCShellGetContext()`, `PCShellSetApplyBA()`, `MATSHELL`, `PCShellSetMatApply()`,
853: M*/
855: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
856: {
857: PC_Shell *shell;
859: PetscNew(&shell);
860: pc->data = (void *)shell;
862: pc->ops->destroy = PCDestroy_Shell;
863: pc->ops->view = PCView_Shell;
864: pc->ops->apply = PCApply_Shell;
865: pc->ops->applysymmetricleft = PCApplySymmetricLeft_Shell;
866: pc->ops->applysymmetricright = PCApplySymmetricRight_Shell;
867: pc->ops->matapply = NULL;
868: pc->ops->applytranspose = NULL;
869: pc->ops->applyrichardson = NULL;
870: pc->ops->setup = NULL;
871: pc->ops->presolve = NULL;
872: pc->ops->postsolve = NULL;
874: shell->apply = NULL;
875: shell->applytranspose = NULL;
876: shell->name = NULL;
877: shell->applyrich = NULL;
878: shell->presolve = NULL;
879: shell->postsolve = NULL;
880: shell->ctx = NULL;
881: shell->setup = NULL;
882: shell->view = NULL;
883: shell->destroy = NULL;
884: shell->applysymmetricleft = NULL;
885: shell->applysymmetricright = NULL;
887: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", PCShellSetDestroy_Shell);
888: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", PCShellSetSetUp_Shell);
889: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", PCShellSetApply_Shell);
890: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", PCShellSetMatApply_Shell);
891: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", PCShellSetApplySymmetricLeft_Shell);
892: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", PCShellSetApplySymmetricRight_Shell);
893: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", PCShellSetApplyBA_Shell);
894: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", PCShellSetPreSolve_Shell);
895: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", PCShellSetPostSolve_Shell);
896: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", PCShellSetView_Shell);
897: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", PCShellSetApplyTranspose_Shell);
898: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", PCShellSetName_Shell);
899: PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", PCShellGetName_Shell);
900: PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", PCShellSetApplyRichardson_Shell);
901: return 0;
902: }