Actual source code: shellpc.c
1: /*
2: This provides a simple shell for Fortran (and C programmers) to
3: create their own preconditioner without writing much interface code.
4: */
6: #include <petsc/private/pcimpl.h>
8: typedef struct {
9: void *ctx; /* user provided contexts for preconditioner */
11: PetscErrorCode (*destroy)(PC);
12: PetscErrorCode (*setup)(PC);
13: PetscErrorCode (*apply)(PC, Vec, Vec);
14: PetscErrorCode (*matapply)(PC, Mat, Mat);
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` that was provided with `PCShellSetContext()`
30: Not Collective
32: Input Parameter:
33: . pc - of type `PCSHELL`
35: Output Parameter:
36: . ctx - the user provided context
38: Level: advanced
40: Note:
41: This routine is intended for use within the various user-provided routines set with, for example, `PCShellSetApply()`
43: Fortran Note:
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: .seealso: [](ch_ksp), `PC`, `PCSHELL`, `PCShellSetContext()`, `PCShellSetApply()`, `PCShellSetDestroy()`
48: @*/
49: PetscErrorCode PCShellGetContext(PC pc, void *ctx)
50: {
51: PetscBool flg;
53: PetscFunctionBegin;
55: PetscAssertPointer(ctx, 2);
56: PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg));
57: if (!flg) *(void **)ctx = NULL;
58: else *(void **)ctx = ((PC_Shell *)(pc->data))->ctx;
59: PetscFunctionReturn(PETSC_SUCCESS);
60: }
62: /*@
63: PCShellSetContext - sets the context for a shell `PC` that can be accessed with `PCShellGetContext()`
65: Logically Collective
67: Input Parameters:
68: + pc - the `PC` of type `PCSHELL`
69: - ctx - the context
71: Level: advanced
73: Notes:
74: This routine is intended for use within the various user-provided routines set with, for example, `PCShellSetApply()`
76: One should also provide a routine to destroy the context when `pc` is destroyed with `PCShellSetDestroy()`
78: Fortran Notes:
79: To use this from Fortran you must write a Fortran interface definition for this
80: function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.
82: .seealso: [](ch_ksp), `PC`, `PCShellGetContext()`, `PCSHELL`, `PCShellSetApply()`, `PCShellSetDestroy()`
83: @*/
84: PetscErrorCode PCShellSetContext(PC pc, void *ctx)
85: {
86: PC_Shell *shell = (PC_Shell *)pc->data;
87: PetscBool flg;
89: PetscFunctionBegin;
91: PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg));
92: if (flg) shell->ctx = ctx;
93: PetscFunctionReturn(PETSC_SUCCESS);
94: }
96: static PetscErrorCode PCSetUp_Shell(PC pc)
97: {
98: PC_Shell *shell = (PC_Shell *)pc->data;
100: PetscFunctionBegin;
101: PetscCheck(shell->setup, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No setup() routine provided to Shell PC");
102: PetscCallBack("PCSHELL callback setup", (*shell->setup)(pc));
103: PetscFunctionReturn(PETSC_SUCCESS);
104: }
106: static PetscErrorCode PCApply_Shell(PC pc, Vec x, Vec y)
107: {
108: PC_Shell *shell = (PC_Shell *)pc->data;
109: PetscObjectState instate, outstate;
111: PetscFunctionBegin;
112: PetscCheck(shell->apply, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
113: PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
114: PetscCallBack("PCSHELL callback apply", (*shell->apply)(pc, x, y));
115: PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
116: /* increase the state of the output vector if the user did not update its state themself as should have been done */
117: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
118: PetscFunctionReturn(PETSC_SUCCESS);
119: }
121: static PetscErrorCode PCMatApply_Shell(PC pc, Mat X, Mat Y)
122: {
123: PC_Shell *shell = (PC_Shell *)pc->data;
124: PetscObjectState instate, outstate;
126: PetscFunctionBegin;
127: PetscCheck(shell->matapply, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
128: PetscCall(PetscObjectStateGet((PetscObject)Y, &instate));
129: PetscCallBack("PCSHELL callback apply", (*shell->matapply)(pc, X, Y));
130: PetscCall(PetscObjectStateGet((PetscObject)Y, &outstate));
131: /* increase the state of the output vector if the user did not update its state themself as should have been done */
132: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)Y));
133: PetscFunctionReturn(PETSC_SUCCESS);
134: }
136: static PetscErrorCode PCApplySymmetricLeft_Shell(PC pc, Vec x, Vec y)
137: {
138: PC_Shell *shell = (PC_Shell *)pc->data;
140: PetscFunctionBegin;
141: PetscCheck(shell->applysymmetricleft, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
142: PetscCallBack("PCSHELL callback apply symmetric left", (*shell->applysymmetricleft)(pc, x, y));
143: PetscFunctionReturn(PETSC_SUCCESS);
144: }
146: static PetscErrorCode PCApplySymmetricRight_Shell(PC pc, Vec x, Vec y)
147: {
148: PC_Shell *shell = (PC_Shell *)pc->data;
150: PetscFunctionBegin;
151: PetscCheck(shell->applysymmetricright, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
152: PetscCallBack("PCSHELL callback apply symmetric right", (*shell->applysymmetricright)(pc, x, y));
153: PetscFunctionReturn(PETSC_SUCCESS);
154: }
156: static PetscErrorCode PCApplyBA_Shell(PC pc, PCSide side, Vec x, Vec y, Vec w)
157: {
158: PC_Shell *shell = (PC_Shell *)pc->data;
159: PetscObjectState instate, outstate;
161: PetscFunctionBegin;
162: PetscCheck(shell->applyBA, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applyBA() routine provided to Shell PC");
163: PetscCall(PetscObjectStateGet((PetscObject)w, &instate));
164: PetscCallBack("PCSHELL callback applyBA", (*shell->applyBA)(pc, side, x, y, w));
165: PetscCall(PetscObjectStateGet((PetscObject)w, &outstate));
166: /* increase the state of the output vector if the user did not update its state themself as should have been done */
167: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)w));
168: PetscFunctionReturn(PETSC_SUCCESS);
169: }
171: static PetscErrorCode PCPreSolveChangeRHS_Shell(PC pc, PetscBool *change)
172: {
173: PetscFunctionBegin;
174: *change = PETSC_TRUE;
175: PetscFunctionReturn(PETSC_SUCCESS);
176: }
178: static PetscErrorCode PCPreSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
179: {
180: PC_Shell *shell = (PC_Shell *)pc->data;
182: PetscFunctionBegin;
183: PetscCheck(shell->presolve, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No presolve() routine provided to Shell PC");
184: PetscCallBack("PCSHELL callback presolve", (*shell->presolve)(pc, ksp, b, x));
185: PetscFunctionReturn(PETSC_SUCCESS);
186: }
188: static PetscErrorCode PCPostSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
189: {
190: PC_Shell *shell = (PC_Shell *)pc->data;
192: PetscFunctionBegin;
193: PetscCheck(shell->postsolve, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No postsolve() routine provided to Shell PC");
194: PetscCallBack("PCSHELL callback postsolve()", (*shell->postsolve)(pc, ksp, b, x));
195: PetscFunctionReturn(PETSC_SUCCESS);
196: }
198: static PetscErrorCode PCApplyTranspose_Shell(PC pc, Vec x, Vec y)
199: {
200: PC_Shell *shell = (PC_Shell *)pc->data;
201: PetscObjectState instate, outstate;
203: PetscFunctionBegin;
204: PetscCheck(shell->applytranspose, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applytranspose() routine provided to Shell PC");
205: PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
206: PetscCallBack("PCSHELL callback applytranspose", (*shell->applytranspose)(pc, x, y));
207: PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
208: /* increase the state of the output vector if the user did not update its state themself as should have been done */
209: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
210: PetscFunctionReturn(PETSC_SUCCESS);
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: PC_Shell *shell = (PC_Shell *)pc->data;
216: PetscObjectState instate, outstate;
218: PetscFunctionBegin;
219: PetscCheck(shell->applyrich, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applyrichardson() routine provided to Shell PC");
220: PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
221: PetscCallBack("PCSHELL callback applyrichardson", (*shell->applyrich)(pc, x, y, w, rtol, abstol, dtol, it, guesszero, outits, reason));
222: PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
223: /* increase the state of the output vector since the user did not update its state themself as should have been done */
224: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
225: PetscFunctionReturn(PETSC_SUCCESS);
226: }
228: static PetscErrorCode PCDestroy_Shell(PC pc)
229: {
230: PC_Shell *shell = (PC_Shell *)pc->data;
232: PetscFunctionBegin;
233: PetscCall(PetscFree(shell->name));
234: if (shell->destroy) PetscCallBack("PCSHELL callback destroy", (*shell->destroy)(pc));
235: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", NULL));
236: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", NULL));
237: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", NULL));
238: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", NULL));
239: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", NULL));
240: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", NULL));
241: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", NULL));
242: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", NULL));
243: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", NULL));
244: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", NULL));
245: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", NULL));
246: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", NULL));
247: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", NULL));
248: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", NULL));
249: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL));
250: PetscCall(PetscFree(pc->data));
251: PetscFunctionReturn(PETSC_SUCCESS);
252: }
254: static PetscErrorCode PCView_Shell(PC pc, PetscViewer viewer)
255: {
256: PC_Shell *shell = (PC_Shell *)pc->data;
257: PetscBool iascii;
259: PetscFunctionBegin;
260: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
261: if (iascii) {
262: if (shell->name) PetscCall(PetscViewerASCIIPrintf(viewer, " %s\n", shell->name));
263: else PetscCall(PetscViewerASCIIPrintf(viewer, " no name\n"));
264: }
265: if (shell->view) {
266: PetscCall(PetscViewerASCIIPushTab(viewer));
267: PetscCall((*shell->view)(pc, viewer));
268: PetscCall(PetscViewerASCIIPopTab(viewer));
269: }
270: PetscFunctionReturn(PETSC_SUCCESS);
271: }
273: static PetscErrorCode PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
274: {
275: PC_Shell *shell = (PC_Shell *)pc->data;
277: PetscFunctionBegin;
278: shell->destroy = destroy;
279: PetscFunctionReturn(PETSC_SUCCESS);
280: }
282: static PetscErrorCode PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
283: {
284: PC_Shell *shell = (PC_Shell *)pc->data;
286: PetscFunctionBegin;
287: shell->setup = setup;
288: if (setup) pc->ops->setup = PCSetUp_Shell;
289: else pc->ops->setup = NULL;
290: PetscFunctionReturn(PETSC_SUCCESS);
291: }
293: static PetscErrorCode PCShellSetApply_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
294: {
295: PC_Shell *shell = (PC_Shell *)pc->data;
297: PetscFunctionBegin;
298: shell->apply = apply;
299: PetscFunctionReturn(PETSC_SUCCESS);
300: }
302: static PetscErrorCode PCShellSetMatApply_Shell(PC pc, PetscErrorCode (*matapply)(PC, Mat, Mat))
303: {
304: PC_Shell *shell = (PC_Shell *)pc->data;
306: PetscFunctionBegin;
307: shell->matapply = matapply;
308: if (matapply) pc->ops->matapply = PCMatApply_Shell;
309: else pc->ops->matapply = NULL;
310: PetscFunctionReturn(PETSC_SUCCESS);
311: }
313: static PetscErrorCode PCShellSetApplySymmetricLeft_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
314: {
315: PC_Shell *shell = (PC_Shell *)pc->data;
317: PetscFunctionBegin;
318: shell->applysymmetricleft = apply;
319: PetscFunctionReturn(PETSC_SUCCESS);
320: }
322: static PetscErrorCode PCShellSetApplySymmetricRight_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
323: {
324: PC_Shell *shell = (PC_Shell *)pc->data;
326: PetscFunctionBegin;
327: shell->applysymmetricright = apply;
328: PetscFunctionReturn(PETSC_SUCCESS);
329: }
331: static PetscErrorCode PCShellSetApplyBA_Shell(PC pc, PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec))
332: {
333: PC_Shell *shell = (PC_Shell *)pc->data;
335: PetscFunctionBegin;
336: shell->applyBA = applyBA;
337: if (applyBA) pc->ops->applyBA = PCApplyBA_Shell;
338: else pc->ops->applyBA = NULL;
339: PetscFunctionReturn(PETSC_SUCCESS);
340: }
342: static PetscErrorCode PCShellSetPreSolve_Shell(PC pc, PetscErrorCode (*presolve)(PC, KSP, Vec, Vec))
343: {
344: PC_Shell *shell = (PC_Shell *)pc->data;
346: PetscFunctionBegin;
347: shell->presolve = presolve;
348: if (presolve) {
349: pc->ops->presolve = PCPreSolve_Shell;
350: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", PCPreSolveChangeRHS_Shell));
351: } else {
352: pc->ops->presolve = NULL;
353: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL));
354: }
355: PetscFunctionReturn(PETSC_SUCCESS);
356: }
358: static PetscErrorCode PCShellSetPostSolve_Shell(PC pc, PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec))
359: {
360: PC_Shell *shell = (PC_Shell *)pc->data;
362: PetscFunctionBegin;
363: shell->postsolve = postsolve;
364: if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
365: else pc->ops->postsolve = NULL;
366: PetscFunctionReturn(PETSC_SUCCESS);
367: }
369: static PetscErrorCode PCShellSetView_Shell(PC pc, PetscErrorCode (*view)(PC, PetscViewer))
370: {
371: PC_Shell *shell = (PC_Shell *)pc->data;
373: PetscFunctionBegin;
374: shell->view = view;
375: PetscFunctionReturn(PETSC_SUCCESS);
376: }
378: static PetscErrorCode PCShellSetApplyTranspose_Shell(PC pc, PetscErrorCode (*applytranspose)(PC, Vec, Vec))
379: {
380: PC_Shell *shell = (PC_Shell *)pc->data;
382: PetscFunctionBegin;
383: shell->applytranspose = applytranspose;
384: if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
385: else pc->ops->applytranspose = NULL;
386: PetscFunctionReturn(PETSC_SUCCESS);
387: }
389: static PetscErrorCode PCShellSetApplyRichardson_Shell(PC pc, PetscErrorCode (*applyrich)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *))
390: {
391: PC_Shell *shell = (PC_Shell *)pc->data;
393: PetscFunctionBegin;
394: shell->applyrich = applyrich;
395: if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
396: else pc->ops->applyrichardson = NULL;
397: PetscFunctionReturn(PETSC_SUCCESS);
398: }
400: static PetscErrorCode PCShellSetName_Shell(PC pc, const char name[])
401: {
402: PC_Shell *shell = (PC_Shell *)pc->data;
404: PetscFunctionBegin;
405: PetscCall(PetscFree(shell->name));
406: PetscCall(PetscStrallocpy(name, &shell->name));
407: PetscFunctionReturn(PETSC_SUCCESS);
408: }
410: static PetscErrorCode PCShellGetName_Shell(PC pc, const char *name[])
411: {
412: PC_Shell *shell = (PC_Shell *)pc->data;
414: PetscFunctionBegin;
415: *name = shell->name;
416: PetscFunctionReturn(PETSC_SUCCESS);
417: }
419: /*@C
420: PCShellSetDestroy - Sets routine to use to destroy the user-provided application context that was provided with `PCShellSetContext()`
422: Logically Collective
424: Input Parameters:
425: + pc - the preconditioner context
426: - destroy - the application-provided destroy routine
428: Calling sequence of `destroy`:
429: . pc - the preconditioner
431: Level: intermediate
433: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApply()`, `PCShellSetContext()`, `PCShellGetContext()`
434: @*/
435: PetscErrorCode PCShellSetDestroy(PC pc, PetscErrorCode (*destroy)(PC pc))
436: {
437: PetscFunctionBegin;
439: PetscTryMethod(pc, "PCShellSetDestroy_C", (PC, PetscErrorCode(*)(PC)), (pc, destroy));
440: PetscFunctionReturn(PETSC_SUCCESS);
441: }
443: /*@C
444: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
445: matrix operator is changed.
447: Logically Collective
449: Input Parameters:
450: + pc - the preconditioner context
451: - setup - the application-provided setup routine
453: Calling sequence of `setup`:
454: . pc - the preconditioner
456: Level: intermediate
458: Note:
459: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `setup`.
461: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetApply()`, `PCShellSetContext()`, , `PCShellGetContext()`
462: @*/
463: PetscErrorCode PCShellSetSetUp(PC pc, PetscErrorCode (*setup)(PC pc))
464: {
465: PetscFunctionBegin;
467: PetscTryMethod(pc, "PCShellSetSetUp_C", (PC, PetscErrorCode(*)(PC)), (pc, setup));
468: PetscFunctionReturn(PETSC_SUCCESS);
469: }
471: /*@C
472: PCShellSetView - Sets routine to use as viewer of a `PCSHELL` shell preconditioner
474: Logically Collective
476: Input Parameters:
477: + pc - the preconditioner context
478: - view - the application-provided view routine
480: Calling sequence of `view`:
481: + pc - the preconditioner
482: - v - viewer
484: Level: advanced
486: Note:
487: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `view`.
489: .seealso: [](ch_ksp), `PC`, `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellGetContext()`
490: @*/
491: PetscErrorCode PCShellSetView(PC pc, PetscErrorCode (*view)(PC pc, PetscViewer v))
492: {
493: PetscFunctionBegin;
495: PetscTryMethod(pc, "PCShellSetView_C", (PC, PetscErrorCode(*)(PC, PetscViewer)), (pc, view));
496: PetscFunctionReturn(PETSC_SUCCESS);
497: }
499: /*@C
500: PCShellSetApply - Sets routine to use as preconditioner.
502: Logically Collective
504: Input Parameters:
505: + pc - the preconditioner context
506: - apply - the application-provided preconditioning routine
508: Calling sequence of `apply`:
509: + pc - the preconditioner, get the application context with `PCShellGetContext()`
510: . xin - input vector
511: - xout - output vector
513: Level: intermediate
515: Note:
516: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `apply`.
518: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApplyBA()`, `PCShellSetApplySymmetricRight()`, `PCShellSetApplySymmetricLeft()`, `PCShellGetContext()`
519: @*/
520: PetscErrorCode PCShellSetApply(PC pc, PetscErrorCode (*apply)(PC pc, Vec xin, Vec xout))
521: {
522: PetscFunctionBegin;
524: PetscTryMethod(pc, "PCShellSetApply_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
525: PetscFunctionReturn(PETSC_SUCCESS);
526: }
528: /*@C
529: PCShellSetMatApply - Sets routine to use as preconditioner on a block of vectors.
531: Logically Collective
533: Input Parameters:
534: + pc - the preconditioner context
535: - matapply - the application-provided preconditioning routine
537: Calling sequence of `matapply`:
538: + pc - the preconditioner
539: . Xin - input block of vectors represented as a dense `Mat`
540: - Xout - output block of vectors represented as a dense `Mat`
542: Level: advanced
544: Note:
545: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `matapply`.
547: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApply()`, `PCShellSetContext()`, `PCShellGetContext()`
548: @*/
549: PetscErrorCode PCShellSetMatApply(PC pc, PetscErrorCode (*matapply)(PC pc, Mat Xin, Mat Xout))
550: {
551: PetscFunctionBegin;
553: PetscTryMethod(pc, "PCShellSetMatApply_C", (PC, PetscErrorCode(*)(PC, Mat, Mat)), (pc, matapply));
554: PetscFunctionReturn(PETSC_SUCCESS);
555: }
557: /*@C
558: PCShellSetApplySymmetricLeft - Sets routine to use as left preconditioner (when the `PC_SYMMETRIC` is used).
560: Logically Collective
562: Input Parameters:
563: + pc - the preconditioner context
564: - apply - the application-provided left preconditioning routine
566: Calling sequence of `apply`:
567: + pc - the preconditioner
568: . xin - input vector
569: - xout - output vector
571: Level: advanced
573: Note:
574: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `apply`.
576: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApply()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`
577: @*/
578: PetscErrorCode PCShellSetApplySymmetricLeft(PC pc, PetscErrorCode (*apply)(PC pc, Vec xin, Vec xout))
579: {
580: PetscFunctionBegin;
582: PetscTryMethod(pc, "PCShellSetApplySymmetricLeft_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
583: PetscFunctionReturn(PETSC_SUCCESS);
584: }
586: /*@C
587: PCShellSetApplySymmetricRight - Sets routine to use as right preconditioner (when the `PC_SYMMETRIC` is used).
589: Logically Collective
591: Input Parameters:
592: + pc - the preconditioner context
593: - apply - the application-provided right preconditioning routine
595: Calling sequence of `apply`:
596: + pc - the preconditioner
597: . xin - input vector
598: - xout - output vector
600: Level: advanced
602: Note:
603: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `apply`.
605: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApply()`, `PCShellSetApplySymmetricLeft()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellGetContext()`
606: @*/
607: PetscErrorCode PCShellSetApplySymmetricRight(PC pc, PetscErrorCode (*apply)(PC pc, Vec xin, Vec xout))
608: {
609: PetscFunctionBegin;
611: PetscTryMethod(pc, "PCShellSetApplySymmetricRight_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
612: PetscFunctionReturn(PETSC_SUCCESS);
613: }
615: /*@C
616: PCShellSetApplyBA - Sets routine to use as the preconditioner times the operator.
618: Logically Collective
620: Input Parameters:
621: + pc - the preconditioner context
622: - applyBA - the application-provided BA routine
624: Calling sequence of `applyBA`:
625: + pc - the preconditioner
626: . side - `PC_LEFT`, `PC_RIGHT`, or `PC_SYMMETRIC`
627: . xin - input vector
628: . xout - output vector
629: - w - work vector
631: Level: intermediate
633: Note:
634: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `applyBA`.
636: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApply()`, `PCShellGetContext()`, `PCSide`
637: @*/
638: PetscErrorCode PCShellSetApplyBA(PC pc, PetscErrorCode (*applyBA)(PC pc, PCSide side, Vec xin, Vec xout, Vec w))
639: {
640: PetscFunctionBegin;
642: PetscTryMethod(pc, "PCShellSetApplyBA_C", (PC, PetscErrorCode(*)(PC, PCSide, Vec, Vec, Vec)), (pc, applyBA));
643: PetscFunctionReturn(PETSC_SUCCESS);
644: }
646: /*@C
647: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
649: Logically Collective
651: Input Parameters:
652: + pc - the preconditioner context
653: - applytranspose - the application-provided preconditioning transpose routine
655: Calling sequence of `applytranspose`:
656: + pc - the preconditioner
657: . xin - input vector
658: - xout - output vector
660: Level: intermediate
662: Note:
663: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `applytranspose`.
665: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApply()`, `PCSetContext()`, `PCShellSetApplyBA()`, `PCGetContext()`
666: @*/
667: PetscErrorCode PCShellSetApplyTranspose(PC pc, PetscErrorCode (*applytranspose)(PC pc, Vec xin, Vec xout))
668: {
669: PetscFunctionBegin;
671: PetscTryMethod(pc, "PCShellSetApplyTranspose_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, applytranspose));
672: PetscFunctionReturn(PETSC_SUCCESS);
673: }
675: /*@C
676: PCShellSetPreSolve - Sets routine to apply to the operators/vectors before a `KSPSolve()` is
677: applied. This usually does something like scale the linear system in some application
678: specific way.
680: Logically Collective
682: Input Parameters:
683: + pc - the preconditioner context
684: - presolve - the application-provided presolve routine
686: Calling sequence of `presolve`:
687: + pc - the preconditioner
688: . ksp - the `KSP` that contains `pc`
689: . xin - input vector
690: - xout - output vector
692: Level: advanced
694: Note:
695: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `presolve`.
697: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPostSolve()`, `PCShellSetContext()`, `PCGetContext()`
698: @*/
699: PetscErrorCode PCShellSetPreSolve(PC pc, PetscErrorCode (*presolve)(PC pc, KSP ksp, Vec xin, Vec xout))
700: {
701: PetscFunctionBegin;
703: PetscTryMethod(pc, "PCShellSetPreSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, presolve));
704: PetscFunctionReturn(PETSC_SUCCESS);
705: }
707: /*@C
708: PCShellSetPostSolve - Sets routine to apply to the operators/vectors after a `KSPSolve()` is
709: applied. This usually does something like scale the linear system in some application
710: specific way.
712: Logically Collective
714: Input Parameters:
715: + pc - the preconditioner context
716: - postsolve - the application-provided presolve routine
718: Calling sequence of `postsolve`:
719: + pc - the preconditioner
720: . ksp - the `KSP` that contains `pc`
721: . xin - input vector
722: - xout - output vector
724: Level: advanced
726: Note:
727: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `postsolve`.
729: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPreSolve()`, `PCShellSetContext()`, `PCGetContext()`
730: @*/
731: PetscErrorCode PCShellSetPostSolve(PC pc, PetscErrorCode (*postsolve)(PC pc, KSP ksp, Vec xin, Vec xout))
732: {
733: PetscFunctionBegin;
735: PetscTryMethod(pc, "PCShellSetPostSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, postsolve));
736: PetscFunctionReturn(PETSC_SUCCESS);
737: }
739: /*@C
740: PCShellSetName - Sets an optional name to associate with a `PCSHELL`
741: preconditioner.
743: Not Collective
745: Input Parameters:
746: + pc - the preconditioner context
747: - name - character string describing shell preconditioner
749: Level: intermediate
751: Note:
752: This is separate from the name you can provide with `PetscObjectSetName()`
754: .seealso: [](ch_ksp), `PCSHELL`, `PCShellGetName()`, `PetscObjectSetName()`, `PetscObjectGetName()`
755: @*/
756: PetscErrorCode PCShellSetName(PC pc, const char name[])
757: {
758: PetscFunctionBegin;
760: PetscTryMethod(pc, "PCShellSetName_C", (PC, const char[]), (pc, name));
761: PetscFunctionReturn(PETSC_SUCCESS);
762: }
764: /*@C
765: PCShellGetName - Gets an optional name that the user has set for a `PCSHELL` with `PCShellSetName()`
766: preconditioner.
768: Not Collective
770: Input Parameter:
771: . pc - the preconditioner context
773: Output Parameter:
774: . name - character string describing shell preconditioner (you should not free this)
776: Level: intermediate
778: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetName()`, `PetscObjectSetName()`, `PetscObjectGetName()`
779: @*/
780: PetscErrorCode PCShellGetName(PC pc, const char *name[])
781: {
782: PetscFunctionBegin;
784: PetscAssertPointer(name, 2);
785: PetscUseMethod(pc, "PCShellGetName_C", (PC, const char *[]), (pc, name));
786: PetscFunctionReturn(PETSC_SUCCESS);
787: }
789: /*@C
790: PCShellSetApplyRichardson - Sets routine to use as preconditioner
791: in Richardson iteration.
793: Logically Collective
795: Input Parameters:
796: + pc - the preconditioner context
797: - apply - the application-provided preconditioning routine
799: Calling sequence of `apply`:
800: + pc - the preconditioner
801: . b - right-hand-side
802: . x - current iterate
803: . r - work space
804: . rtol - relative tolerance of residual norm to stop at
805: . abstol - absolute tolerance of residual norm to stop at
806: . dtol - if residual norm increases by this factor than return
807: . maxits - number of iterations to run
808: . zeroinitialguess - `PETSC_TRUE` if `x` is known to be initially zero
809: . its - returns the number of iterations used
810: - reason - returns the reason the iteration has converged
812: Level: advanced
814: Note:
815: You can get the `PCSHELL` context set with `PCShellSetContext()` using `PCShellGetContext()` if needed by `apply`.
817: .seealso: [](ch_ksp), `PCSHELL`, `PCShellSetApply()`, `PCShellSetContext()`, `PCRichardsonConvergedReason()`, `PCShellGetContext()`
818: @*/
819: PetscErrorCode PCShellSetApplyRichardson(PC pc, PetscErrorCode (*apply)(PC pc, Vec b, Vec x, Vec r, PetscReal rtol, PetscReal abstol, PetscReal dtol, PetscInt maxits, PetscBool zeroinitialguess, PetscInt *its, PCRichardsonConvergedReason *reason))
820: {
821: PetscFunctionBegin;
823: PetscTryMethod(pc, "PCShellSetApplyRichardson_C", (PC, PetscErrorCode(*)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *)), (pc, apply));
824: PetscFunctionReturn(PETSC_SUCCESS);
825: }
827: /*MC
828: PCSHELL - Creates a new preconditioner class for use with a users
829: own private data storage format and preconditioner application code
831: Level: advanced
833: Usage:
834: .vb
835: extern PetscErrorCode apply(PC,Vec,Vec);
836: extern PetscErrorCode applyba(PC,PCSide,Vec,Vec,Vec);
837: extern PetscErrorCode applytranspose(PC,Vec,Vec);
838: extern PetscErrorCode setup(PC);
839: extern PetscErrorCode destroy(PC);
841: PCCreate(comm,&pc);
842: PCSetType(pc,PCSHELL);
843: PCShellSetContext(pc,ctx)
844: PCShellSetApply(pc,apply);
845: PCShellSetApplyBA(pc,applyba); (optional)
846: PCShellSetApplyTranspose(pc,applytranspose); (optional)
847: PCShellSetSetUp(pc,setup); (optional)
848: PCShellSetDestroy(pc,destroy); (optional)
849: .ve
851: Note:
852: Information required for the preconditioner and its internal datastructures can be set with `PCShellSetContext()` and then accessed
853: with `PCShellGetContext()` inside the routines provided above
855: .seealso: [](ch_ksp), `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
856: `MATSHELL`, `PCShellSetSetUp()`, `PCShellSetApply()`, `PCShellSetView()`, `PCShellSetDestroy()`, `PCShellSetPostSolve()`,
857: `PCShellSetApplyTranspose()`, `PCShellSetName()`, `PCShellSetApplyRichardson()`, `PCShellSetPreSolve()`, `PCShellSetView()`,
858: `PCShellGetName()`, `PCShellSetContext()`, `PCShellGetContext()`, `PCShellSetApplyBA()`, `MATSHELL`, `PCShellSetMatApply()`,
859: M*/
861: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
862: {
863: PC_Shell *shell;
865: PetscFunctionBegin;
866: PetscCall(PetscNew(&shell));
867: pc->data = (void *)shell;
869: pc->ops->destroy = PCDestroy_Shell;
870: pc->ops->view = PCView_Shell;
871: pc->ops->apply = PCApply_Shell;
872: pc->ops->applysymmetricleft = PCApplySymmetricLeft_Shell;
873: pc->ops->applysymmetricright = PCApplySymmetricRight_Shell;
874: pc->ops->matapply = NULL;
875: pc->ops->applytranspose = NULL;
876: pc->ops->applyrichardson = NULL;
877: pc->ops->setup = NULL;
878: pc->ops->presolve = NULL;
879: pc->ops->postsolve = NULL;
881: shell->apply = NULL;
882: shell->applytranspose = NULL;
883: shell->name = NULL;
884: shell->applyrich = NULL;
885: shell->presolve = NULL;
886: shell->postsolve = NULL;
887: shell->ctx = NULL;
888: shell->setup = NULL;
889: shell->view = NULL;
890: shell->destroy = NULL;
891: shell->applysymmetricleft = NULL;
892: shell->applysymmetricright = NULL;
894: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", PCShellSetDestroy_Shell));
895: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", PCShellSetSetUp_Shell));
896: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", PCShellSetApply_Shell));
897: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", PCShellSetMatApply_Shell));
898: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", PCShellSetApplySymmetricLeft_Shell));
899: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", PCShellSetApplySymmetricRight_Shell));
900: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", PCShellSetApplyBA_Shell));
901: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", PCShellSetPreSolve_Shell));
902: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", PCShellSetPostSolve_Shell));
903: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", PCShellSetView_Shell));
904: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", PCShellSetApplyTranspose_Shell));
905: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", PCShellSetName_Shell));
906: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", PCShellGetName_Shell));
907: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", PCShellSetApplyRichardson_Shell));
908: PetscFunctionReturn(PETSC_SUCCESS);
909: }