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;
54: PetscFunctionBegin;
57: PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg));
58: if (!flg) *(void **)ctx = NULL;
59: else *(void **)ctx = ((PC_Shell *)(pc->data))->ctx;
60: PetscFunctionReturn(PETSC_SUCCESS);
61: }
63: /*@
64: PCShellSetContext - sets the context for a shell `PC`
66: Logically Collective
68: Input Parameters:
69: + pc - the `PC` of type `PCSHELL`
70: - ctx - the context
72: Level: advanced
74: Fortran Note:
75: To use this from Fortran you must write a Fortran interface definition for this
76: function that tells Fortran the Fortran derived data type that you are passing in as the ctx argument.
78: .seealso: `PCShellGetContext()`, `PCSHELL`
79: @*/
80: PetscErrorCode PCShellSetContext(PC pc, void *ctx)
81: {
82: PC_Shell *shell = (PC_Shell *)pc->data;
83: PetscBool flg;
85: PetscFunctionBegin;
87: PetscCall(PetscObjectTypeCompare((PetscObject)pc, PCSHELL, &flg));
88: if (flg) shell->ctx = ctx;
89: PetscFunctionReturn(PETSC_SUCCESS);
90: }
92: static PetscErrorCode PCSetUp_Shell(PC pc)
93: {
94: PC_Shell *shell = (PC_Shell *)pc->data;
96: PetscFunctionBegin;
97: PetscCheck(shell->setup, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No setup() routine provided to Shell PC");
98: PetscCallBack("PCSHELL callback setup", (*shell->setup)(pc));
99: PetscFunctionReturn(PETSC_SUCCESS);
100: }
102: static PetscErrorCode PCApply_Shell(PC pc, Vec x, Vec y)
103: {
104: PC_Shell *shell = (PC_Shell *)pc->data;
105: PetscObjectState instate, outstate;
107: PetscFunctionBegin;
108: PetscCheck(shell->apply, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
109: PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
110: PetscCallBack("PCSHELL callback apply", (*shell->apply)(pc, x, y));
111: PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
112: /* increase the state of the output vector if the user did not update its state themself as should have been done */
113: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
114: PetscFunctionReturn(PETSC_SUCCESS);
115: }
117: static PetscErrorCode PCMatApply_Shell(PC pc, Mat X, Mat Y)
118: {
119: PC_Shell *shell = (PC_Shell *)pc->data;
120: PetscObjectState instate, outstate;
122: PetscFunctionBegin;
123: PetscCheck(shell->matapply, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
124: PetscCall(PetscObjectStateGet((PetscObject)Y, &instate));
125: PetscCallBack("PCSHELL callback apply", (*shell->matapply)(pc, X, Y));
126: PetscCall(PetscObjectStateGet((PetscObject)Y, &outstate));
127: /* increase the state of the output vector if the user did not update its state themself as should have been done */
128: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)Y));
129: PetscFunctionReturn(PETSC_SUCCESS);
130: }
132: static PetscErrorCode PCApplySymmetricLeft_Shell(PC pc, Vec x, Vec y)
133: {
134: PC_Shell *shell = (PC_Shell *)pc->data;
136: PetscFunctionBegin;
137: PetscCheck(shell->applysymmetricleft, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
138: PetscCallBack("PCSHELL callback apply symmetric left", (*shell->applysymmetricleft)(pc, x, y));
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: static PetscErrorCode PCApplySymmetricRight_Shell(PC pc, Vec x, Vec y)
143: {
144: PC_Shell *shell = (PC_Shell *)pc->data;
146: PetscFunctionBegin;
147: PetscCheck(shell->applysymmetricright, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No apply() routine provided to Shell PC");
148: PetscCallBack("PCSHELL callback apply symmetric right", (*shell->applysymmetricright)(pc, x, y));
149: PetscFunctionReturn(PETSC_SUCCESS);
150: }
152: static PetscErrorCode PCApplyBA_Shell(PC pc, PCSide side, Vec x, Vec y, Vec w)
153: {
154: PC_Shell *shell = (PC_Shell *)pc->data;
155: PetscObjectState instate, outstate;
157: PetscFunctionBegin;
158: PetscCheck(shell->applyBA, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applyBA() routine provided to Shell PC");
159: PetscCall(PetscObjectStateGet((PetscObject)w, &instate));
160: PetscCallBack("PCSHELL callback applyBA", (*shell->applyBA)(pc, side, x, y, w));
161: PetscCall(PetscObjectStateGet((PetscObject)w, &outstate));
162: /* increase the state of the output vector if the user did not update its state themself as should have been done */
163: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)w));
164: PetscFunctionReturn(PETSC_SUCCESS);
165: }
167: static PetscErrorCode PCPreSolveChangeRHS_Shell(PC pc, PetscBool *change)
168: {
169: PetscFunctionBegin;
170: *change = PETSC_TRUE;
171: PetscFunctionReturn(PETSC_SUCCESS);
172: }
174: static PetscErrorCode PCPreSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
175: {
176: PC_Shell *shell = (PC_Shell *)pc->data;
178: PetscFunctionBegin;
179: PetscCheck(shell->presolve, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No presolve() routine provided to Shell PC");
180: PetscCallBack("PCSHELL callback presolve", (*shell->presolve)(pc, ksp, b, x));
181: PetscFunctionReturn(PETSC_SUCCESS);
182: }
184: static PetscErrorCode PCPostSolve_Shell(PC pc, KSP ksp, Vec b, Vec x)
185: {
186: PC_Shell *shell = (PC_Shell *)pc->data;
188: PetscFunctionBegin;
189: PetscCheck(shell->postsolve, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No postsolve() routine provided to Shell PC");
190: PetscCallBack("PCSHELL callback postsolve()", (*shell->postsolve)(pc, ksp, b, x));
191: PetscFunctionReturn(PETSC_SUCCESS);
192: }
194: static PetscErrorCode PCApplyTranspose_Shell(PC pc, Vec x, Vec y)
195: {
196: PC_Shell *shell = (PC_Shell *)pc->data;
197: PetscObjectState instate, outstate;
199: PetscFunctionBegin;
200: PetscCheck(shell->applytranspose, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applytranspose() routine provided to Shell PC");
201: PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
202: PetscCallBack("PCSHELL callback applytranspose", (*shell->applytranspose)(pc, x, y));
203: PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
204: /* increase the state of the output vector if the user did not update its state themself as should have been done */
205: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
206: PetscFunctionReturn(PETSC_SUCCESS);
207: }
209: 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)
210: {
211: PC_Shell *shell = (PC_Shell *)pc->data;
212: PetscObjectState instate, outstate;
214: PetscFunctionBegin;
215: PetscCheck(shell->applyrich, PetscObjectComm((PetscObject)pc), PETSC_ERR_USER, "No applyrichardson() routine provided to Shell PC");
216: PetscCall(PetscObjectStateGet((PetscObject)y, &instate));
217: PetscCallBack("PCSHELL callback applyrichardson", (*shell->applyrich)(pc, x, y, w, rtol, abstol, dtol, it, guesszero, outits, reason));
218: PetscCall(PetscObjectStateGet((PetscObject)y, &outstate));
219: /* increase the state of the output vector since the user did not update its state themself as should have been done */
220: if (instate == outstate) PetscCall(PetscObjectStateIncrease((PetscObject)y));
221: PetscFunctionReturn(PETSC_SUCCESS);
222: }
224: static PetscErrorCode PCDestroy_Shell(PC pc)
225: {
226: PC_Shell *shell = (PC_Shell *)pc->data;
228: PetscFunctionBegin;
229: PetscCall(PetscFree(shell->name));
230: if (shell->destroy) PetscCallBack("PCSHELL callback destroy", (*shell->destroy)(pc));
231: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", NULL));
232: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", NULL));
233: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", NULL));
234: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", NULL));
235: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", NULL));
236: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", NULL));
237: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", NULL));
238: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", NULL));
239: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", NULL));
240: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", NULL));
241: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", NULL));
242: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", NULL));
243: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", NULL));
244: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", NULL));
245: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL));
246: PetscCall(PetscFree(pc->data));
247: PetscFunctionReturn(PETSC_SUCCESS);
248: }
250: static PetscErrorCode PCView_Shell(PC pc, PetscViewer viewer)
251: {
252: PC_Shell *shell = (PC_Shell *)pc->data;
253: PetscBool iascii;
255: PetscFunctionBegin;
256: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
257: if (iascii) {
258: if (shell->name) PetscCall(PetscViewerASCIIPrintf(viewer, " %s\n", shell->name));
259: else PetscCall(PetscViewerASCIIPrintf(viewer, " no name\n"));
260: }
261: if (shell->view) {
262: PetscCall(PetscViewerASCIIPushTab(viewer));
263: PetscCall((*shell->view)(pc, viewer));
264: PetscCall(PetscViewerASCIIPopTab(viewer));
265: }
266: PetscFunctionReturn(PETSC_SUCCESS);
267: }
269: static PetscErrorCode PCShellSetDestroy_Shell(PC pc, PetscErrorCode (*destroy)(PC))
270: {
271: PC_Shell *shell = (PC_Shell *)pc->data;
273: PetscFunctionBegin;
274: shell->destroy = destroy;
275: PetscFunctionReturn(PETSC_SUCCESS);
276: }
278: static PetscErrorCode PCShellSetSetUp_Shell(PC pc, PetscErrorCode (*setup)(PC))
279: {
280: PC_Shell *shell = (PC_Shell *)pc->data;
282: PetscFunctionBegin;
283: shell->setup = setup;
284: if (setup) pc->ops->setup = PCSetUp_Shell;
285: else pc->ops->setup = NULL;
286: PetscFunctionReturn(PETSC_SUCCESS);
287: }
289: static PetscErrorCode PCShellSetApply_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
290: {
291: PC_Shell *shell = (PC_Shell *)pc->data;
293: PetscFunctionBegin;
294: shell->apply = apply;
295: PetscFunctionReturn(PETSC_SUCCESS);
296: }
298: static PetscErrorCode PCShellSetMatApply_Shell(PC pc, PetscErrorCode (*matapply)(PC, Mat, Mat))
299: {
300: PC_Shell *shell = (PC_Shell *)pc->data;
302: PetscFunctionBegin;
303: shell->matapply = matapply;
304: if (matapply) pc->ops->matapply = PCMatApply_Shell;
305: else pc->ops->matapply = NULL;
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }
309: static PetscErrorCode PCShellSetApplySymmetricLeft_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
310: {
311: PC_Shell *shell = (PC_Shell *)pc->data;
313: PetscFunctionBegin;
314: shell->applysymmetricleft = apply;
315: PetscFunctionReturn(PETSC_SUCCESS);
316: }
318: static PetscErrorCode PCShellSetApplySymmetricRight_Shell(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
319: {
320: PC_Shell *shell = (PC_Shell *)pc->data;
322: PetscFunctionBegin;
323: shell->applysymmetricright = apply;
324: PetscFunctionReturn(PETSC_SUCCESS);
325: }
327: static PetscErrorCode PCShellSetApplyBA_Shell(PC pc, PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec))
328: {
329: PC_Shell *shell = (PC_Shell *)pc->data;
331: PetscFunctionBegin;
332: shell->applyBA = applyBA;
333: if (applyBA) pc->ops->applyBA = PCApplyBA_Shell;
334: else pc->ops->applyBA = NULL;
335: PetscFunctionReturn(PETSC_SUCCESS);
336: }
338: static PetscErrorCode PCShellSetPreSolve_Shell(PC pc, PetscErrorCode (*presolve)(PC, KSP, Vec, Vec))
339: {
340: PC_Shell *shell = (PC_Shell *)pc->data;
342: PetscFunctionBegin;
343: shell->presolve = presolve;
344: if (presolve) {
345: pc->ops->presolve = PCPreSolve_Shell;
346: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", PCPreSolveChangeRHS_Shell));
347: } else {
348: pc->ops->presolve = NULL;
349: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCPreSolveChangeRHS_C", NULL));
350: }
351: PetscFunctionReturn(PETSC_SUCCESS);
352: }
354: static PetscErrorCode PCShellSetPostSolve_Shell(PC pc, PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec))
355: {
356: PC_Shell *shell = (PC_Shell *)pc->data;
358: PetscFunctionBegin;
359: shell->postsolve = postsolve;
360: if (postsolve) pc->ops->postsolve = PCPostSolve_Shell;
361: else pc->ops->postsolve = NULL;
362: PetscFunctionReturn(PETSC_SUCCESS);
363: }
365: static PetscErrorCode PCShellSetView_Shell(PC pc, PetscErrorCode (*view)(PC, PetscViewer))
366: {
367: PC_Shell *shell = (PC_Shell *)pc->data;
369: PetscFunctionBegin;
370: shell->view = view;
371: PetscFunctionReturn(PETSC_SUCCESS);
372: }
374: static PetscErrorCode PCShellSetApplyTranspose_Shell(PC pc, PetscErrorCode (*applytranspose)(PC, Vec, Vec))
375: {
376: PC_Shell *shell = (PC_Shell *)pc->data;
378: PetscFunctionBegin;
379: shell->applytranspose = applytranspose;
380: if (applytranspose) pc->ops->applytranspose = PCApplyTranspose_Shell;
381: else pc->ops->applytranspose = NULL;
382: PetscFunctionReturn(PETSC_SUCCESS);
383: }
385: static PetscErrorCode PCShellSetApplyRichardson_Shell(PC pc, PetscErrorCode (*applyrich)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *))
386: {
387: PC_Shell *shell = (PC_Shell *)pc->data;
389: PetscFunctionBegin;
390: shell->applyrich = applyrich;
391: if (applyrich) pc->ops->applyrichardson = PCApplyRichardson_Shell;
392: else pc->ops->applyrichardson = NULL;
393: PetscFunctionReturn(PETSC_SUCCESS);
394: }
396: static PetscErrorCode PCShellSetName_Shell(PC pc, const char name[])
397: {
398: PC_Shell *shell = (PC_Shell *)pc->data;
400: PetscFunctionBegin;
401: PetscCall(PetscFree(shell->name));
402: PetscCall(PetscStrallocpy(name, &shell->name));
403: PetscFunctionReturn(PETSC_SUCCESS);
404: }
406: static PetscErrorCode PCShellGetName_Shell(PC pc, const char *name[])
407: {
408: PC_Shell *shell = (PC_Shell *)pc->data;
410: PetscFunctionBegin;
411: *name = shell->name;
412: PetscFunctionReturn(PETSC_SUCCESS);
413: }
415: /*@C
416: PCShellSetDestroy - Sets routine to use to destroy the user-provided
417: application context.
419: Logically Collective
421: Input Parameters:
422: + pc - the preconditioner context
423: - destroy - the application-provided destroy routine
425: Calling sequence of `destroy`:
426: $ PetscErrorCode destroy(PC pc)
427: . pc - the preconditioner, get the application context with `PCShellGetContext()`
429: Level: intermediate
431: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetContext()`
432: @*/
433: PetscErrorCode PCShellSetDestroy(PC pc, PetscErrorCode (*destroy)(PC))
434: {
435: PetscFunctionBegin;
437: PetscTryMethod(pc, "PCShellSetDestroy_C", (PC, PetscErrorCode(*)(PC)), (pc, destroy));
438: PetscFunctionReturn(PETSC_SUCCESS);
439: }
441: /*@C
442: PCShellSetSetUp - Sets routine to use to "setup" the preconditioner whenever the
443: matrix operator is changed.
445: Logically Collective
447: Input Parameters:
448: + pc - the preconditioner context
449: - setup - the application-provided setup routine
451: Calling sequence of `setup`:
452: $ PetscErrorCode setup(PC pc)
453: . pc - the preconditioner, get the application context with `PCShellGetContext()`
455: Note:
456: the function MUST return an error code of 0 on success and nonzero on failure.
458: Level: intermediate
460: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetApply()`, `PCShellSetContext()`
461: @*/
462: PetscErrorCode PCShellSetSetUp(PC pc, PetscErrorCode (*setup)(PC))
463: {
464: PetscFunctionBegin;
466: PetscTryMethod(pc, "PCShellSetSetUp_C", (PC, PetscErrorCode(*)(PC)), (pc, setup));
467: PetscFunctionReturn(PETSC_SUCCESS);
468: }
470: /*@C
471: PCShellSetView - Sets routine to use as viewer of shell preconditioner
473: Logically Collective
475: Input Parameters:
476: + pc - the preconditioner context
477: - view - the application-provided view routine
479: Calling sequence of `view`:
480: .vb
481: PetscErrorCode view(PC pc, PetscViewer v)
482: .ve
483: + pc - the preconditioner, get the application context with `PCShellGetContext()`
484: - v - viewer
486: Level: advanced
488: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`
489: @*/
490: PetscErrorCode PCShellSetView(PC pc, PetscErrorCode (*view)(PC, PetscViewer))
491: {
492: PetscFunctionBegin;
494: PetscTryMethod(pc, "PCShellSetView_C", (PC, PetscErrorCode(*)(PC, PetscViewer)), (pc, view));
495: PetscFunctionReturn(PETSC_SUCCESS);
496: }
498: /*@C
499: PCShellSetApply - Sets routine to use as preconditioner.
501: Logically Collective
503: Input Parameters:
504: + pc - the preconditioner context
505: - apply - the application-provided preconditioning routine
507: Calling sequence of `apply`:
508: .vb
509: PetscErrorCode apply(PC pc, Vec xin, Vec xout)
510: .ve
511: + pc - the preconditioner, get the application context with `PCShellGetContext()`
512: . xin - input vector
513: - xout - output vector
515: Level: intermediate
517: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApplyBA()`, `PCShellSetApplySymmetricRight()`, `PCShellSetApplySymmetricLeft()`
518: @*/
519: PetscErrorCode PCShellSetApply(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
520: {
521: PetscFunctionBegin;
523: PetscTryMethod(pc, "PCShellSetApply_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
524: PetscFunctionReturn(PETSC_SUCCESS);
525: }
527: /*@C
528: PCShellSetMatApply - Sets routine to use as preconditioner on a block of vectors.
530: Logically Collective
532: Input Parameters:
533: + pc - the preconditioner context
534: - apply - the application-provided preconditioning routine
536: Calling sequence of `apply`:
537: .vb
538: PetscErrorCode apply(PC pc, Mat Xin, Mat Xout)
539: .ve
540: + pc - the preconditioner, get the application context with `PCShellGetContext()`
541: . Xin - input block of vectors
542: - Xout - output block of vectors
544: Level: advanced
546: .seealso: `PCSHELL`, `PCShellSetApply()`
547: @*/
548: PetscErrorCode PCShellSetMatApply(PC pc, PetscErrorCode (*matapply)(PC, Mat, Mat))
549: {
550: PetscFunctionBegin;
552: PetscTryMethod(pc, "PCShellSetMatApply_C", (PC, PetscErrorCode(*)(PC, Mat, Mat)), (pc, matapply));
553: PetscFunctionReturn(PETSC_SUCCESS);
554: }
556: /*@C
557: PCShellSetApplySymmetricLeft - Sets routine to use as left preconditioner (when the `PC_SYMMETRIC` is used).
559: Logically Collective
561: Input Parameters:
562: + pc - the preconditioner context
563: - apply - the application-provided left preconditioning routine
565: Calling sequence of `apply`:
566: .vb
567: PetscErrorCode apply(PC pc, Vec xin, Vec xout)
568: .ve
569: + pc - the preconditioner, get the application context with `PCShellGetContext()`
570: . xin - input vector
571: - xout - output vector
573: Level: advanced
575: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetApplySymmetricLeft()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`
576: @*/
577: PetscErrorCode PCShellSetApplySymmetricLeft(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
578: {
579: PetscFunctionBegin;
581: PetscTryMethod(pc, "PCShellSetApplySymmetricLeft_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
582: PetscFunctionReturn(PETSC_SUCCESS);
583: }
585: /*@C
586: PCShellSetApplySymmetricRight - Sets routine to use as right preconditioner (when the `PC_SYMMETRIC` is used).
588: Logically Collective
590: Input Parameters:
591: + pc - the preconditioner context
592: - apply - the application-provided right preconditioning routine
594: Calling sequence of `apply`:
595: .vb
596: PetscErrorCode apply(PC pc, Vec xin, Vec xout)
597: .ve
598: + pc - the preconditioner, get the application context with PCShellGetContext()
599: . xin - input vector
600: - xout - output vector
602: Level: advanced
604: .seealso: `PCSHELL`, `PCShellSetApply()`, `PCShellSetApplySymmetricLeft()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`
605: @*/
606: PetscErrorCode PCShellSetApplySymmetricRight(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec))
607: {
608: PetscFunctionBegin;
610: PetscTryMethod(pc, "PCShellSetApplySymmetricRight_C", (PC, PetscErrorCode(*)(PC, Vec, Vec)), (pc, apply));
611: PetscFunctionReturn(PETSC_SUCCESS);
612: }
614: /*@C
615: PCShellSetApplyBA - Sets routine to use as preconditioner times operator.
617: Logically Collective
619: Input Parameters:
620: + pc - the preconditioner context
621: - applyBA - the application-provided BA routine
623: Calling sequence of `applyBA`:
624: .vb
625: PetscErrorCode applyBA(PC pc, Vec xin, Vec xout)
626: .ve
627: + pc - the preconditioner, get the application context with `PCShellGetContext()`
628: . xin - input vector
629: - xout - output vector
631: Level: intermediate
633: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetContext()`, `PCShellSetApply()`
634: @*/
635: PetscErrorCode PCShellSetApplyBA(PC pc, PetscErrorCode (*applyBA)(PC, PCSide, Vec, Vec, Vec))
636: {
637: PetscFunctionBegin;
639: PetscTryMethod(pc, "PCShellSetApplyBA_C", (PC, PetscErrorCode(*)(PC, PCSide, Vec, Vec, Vec)), (pc, applyBA));
640: PetscFunctionReturn(PETSC_SUCCESS);
641: }
643: /*@C
644: PCShellSetApplyTranspose - Sets routine to use as preconditioner transpose.
646: Logically Collective
648: Input Parameters:
649: + pc - the preconditioner context
650: - applytranspose - the application-provided preconditioning transpose routine
652: Calling sequence of `applytranspose`:
653: .vb
654: PetscErrorCode applytranspose(PC pc, Vec xin, Vec xout)
655: .ve
656: + pc - the preconditioner, get the application context with `PCShellGetContext()`
657: . xin - input vector
658: - xout - output vector
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: {
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: .vb
688: PetscErrorCode presolve(PC pc, KSP ksp, Vec b, Vec x)
689: .ve
690: + pc - the preconditioner, get the application context with `PCShellGetContext()`
691: . xin - input vector
692: - xout - output vector
694: Level: advanced
696: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPostSolve()`, `PCShellSetContext()`
697: @*/
698: PetscErrorCode PCShellSetPreSolve(PC pc, PetscErrorCode (*presolve)(PC, KSP, Vec, Vec))
699: {
700: PetscFunctionBegin;
702: PetscTryMethod(pc, "PCShellSetPreSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, presolve));
703: PetscFunctionReturn(PETSC_SUCCESS);
704: }
706: /*@C
707: PCShellSetPostSolve - Sets routine to apply to the operators/vectors before a `KSPSolve()` is
708: applied. This usually does something like scale the linear system in some application
709: specific way.
711: Logically Collective
713: Input Parameters:
714: + pc - the preconditioner context
715: - postsolve - the application-provided presolve routine
717: Calling sequence of `postsolve`:
718: .vb
719: PetscErrorCode postsolve(PC pc, KSP ksp, Vec b, Vec x)
720: .ve
721: + pc - the preconditioner, get the application context with `PCShellGetContext()`
722: . xin - input vector
723: - xout - output vector
725: Level: advanced
727: .seealso: `PCSHELL`, `PCShellSetApplyRichardson()`, `PCShellSetSetUp()`, `PCShellSetApplyTranspose()`, `PCShellSetPreSolve()`, `PCShellSetContext()`
728: @*/
729: PetscErrorCode PCShellSetPostSolve(PC pc, PetscErrorCode (*postsolve)(PC, KSP, Vec, Vec))
730: {
731: PetscFunctionBegin;
733: PetscTryMethod(pc, "PCShellSetPostSolve_C", (PC, PetscErrorCode(*)(PC, KSP, Vec, Vec)), (pc, postsolve));
734: PetscFunctionReturn(PETSC_SUCCESS);
735: }
737: /*@C
738: PCShellSetName - Sets an optional name to associate with a `PCSHELL`
739: preconditioner.
741: Not Collective
743: Input Parameters:
744: + pc - the preconditioner context
745: - name - character string describing shell preconditioner
747: Level: intermediate
749: .seealso: `PCSHELL`, `PCShellGetName()`
750: @*/
751: PetscErrorCode PCShellSetName(PC pc, const char name[])
752: {
753: PetscFunctionBegin;
755: PetscTryMethod(pc, "PCShellSetName_C", (PC, const char[]), (pc, name));
756: PetscFunctionReturn(PETSC_SUCCESS);
757: }
759: /*@C
760: PCShellGetName - Gets an optional name that the user has set for a `PCSHELL`
761: preconditioner.
763: Not Collective
765: Input Parameter:
766: . pc - the preconditioner context
768: Output Parameter:
769: . name - character string describing shell preconditioner (you should not free this)
771: Level: intermediate
773: .seealso: `PCSHELL`, `PCShellSetName()`
774: @*/
775: PetscErrorCode PCShellGetName(PC pc, const char *name[])
776: {
777: PetscFunctionBegin;
780: PetscUseMethod(pc, "PCShellGetName_C", (PC, const char *[]), (pc, name));
781: PetscFunctionReturn(PETSC_SUCCESS);
782: }
784: /*@C
785: PCShellSetApplyRichardson - Sets routine to use as preconditioner
786: in Richardson iteration.
788: Logically Collective
790: Input Parameters:
791: + pc - the preconditioner context
792: - apply - the application-provided preconditioning routine
794: Calling sequence of `apply`:
795: .vb
796: PetscErrorCode apply(PC pc, Vec b, Vec x, Vec r, PetscReal rtol, PetscReal abstol, PetscReal dtol, PetscInt maxits)
797: .ve
798: + pc - the preconditioner, get the application context with `PCShellGetContext()`
799: . b - right-hand-side
800: . x - current iterate
801: . r - work space
802: . rtol - relative tolerance of residual norm to stop at
803: . abstol - absolute tolerance of residual norm to stop at
804: . dtol - if residual norm increases by this factor than return
805: - maxits - number of iterations to run
807: Level: advanced
809: .seealso: `PCShellSetApply()`, `PCShellSetContext()`
810: @*/
811: PetscErrorCode PCShellSetApplyRichardson(PC pc, PetscErrorCode (*apply)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *))
812: {
813: PetscFunctionBegin;
815: PetscTryMethod(pc, "PCShellSetApplyRichardson_C", (PC, PetscErrorCode(*)(PC, Vec, Vec, Vec, PetscReal, PetscReal, PetscReal, PetscInt, PetscBool, PetscInt *, PCRichardsonConvergedReason *)), (pc, apply));
816: PetscFunctionReturn(PETSC_SUCCESS);
817: }
819: /*MC
820: PCSHELL - Creates a new preconditioner class for use with a users
821: own private data storage format and preconditioner application code
823: Level: advanced
825: Usage:
826: .vb
827: extern PetscErrorCode apply(PC,Vec,Vec);
828: extern PetscErrorCode applyba(PC,PCSide,Vec,Vec,Vec);
829: extern PetscErrorCode applytranspose(PC,Vec,Vec);
830: extern PetscErrorCode setup(PC);
831: extern PetscErrorCode destroy(PC);
833: PCCreate(comm,&pc);
834: PCSetType(pc,PCSHELL);
835: PCShellSetContext(pc,ctx)
836: PCShellSetApply(pc,apply);
837: PCShellSetApplyBA(pc,applyba); (optional)
838: PCShellSetApplyTranspose(pc,applytranspose); (optional)
839: PCShellSetSetUp(pc,setup); (optional)
840: PCShellSetDestroy(pc,destroy); (optional)
841: .ve
843: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
844: `MATSHELL`, `PCShellSetSetUp()`, `PCShellSetApply()`, `PCShellSetView()`, `PCShellSetDestroy()`, `PCShellSetPostSolve()`,
845: `PCShellSetApplyTranspose()`, `PCShellSetName()`, `PCShellSetApplyRichardson()`, `PCShellSetPreSolve()`, `PCShellSetView()`,
846: `PCShellGetName()`, `PCShellSetContext()`, `PCShellGetContext()`, `PCShellSetApplyBA()`, `MATSHELL`, `PCShellSetMatApply()`,
847: M*/
849: PETSC_EXTERN PetscErrorCode PCCreate_Shell(PC pc)
850: {
851: PC_Shell *shell;
853: PetscFunctionBegin;
854: PetscCall(PetscNew(&shell));
855: pc->data = (void *)shell;
857: pc->ops->destroy = PCDestroy_Shell;
858: pc->ops->view = PCView_Shell;
859: pc->ops->apply = PCApply_Shell;
860: pc->ops->applysymmetricleft = PCApplySymmetricLeft_Shell;
861: pc->ops->applysymmetricright = PCApplySymmetricRight_Shell;
862: pc->ops->matapply = NULL;
863: pc->ops->applytranspose = NULL;
864: pc->ops->applyrichardson = NULL;
865: pc->ops->setup = NULL;
866: pc->ops->presolve = NULL;
867: pc->ops->postsolve = NULL;
869: shell->apply = NULL;
870: shell->applytranspose = NULL;
871: shell->name = NULL;
872: shell->applyrich = NULL;
873: shell->presolve = NULL;
874: shell->postsolve = NULL;
875: shell->ctx = NULL;
876: shell->setup = NULL;
877: shell->view = NULL;
878: shell->destroy = NULL;
879: shell->applysymmetricleft = NULL;
880: shell->applysymmetricright = NULL;
882: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetDestroy_C", PCShellSetDestroy_Shell));
883: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetSetUp_C", PCShellSetSetUp_Shell));
884: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApply_C", PCShellSetApply_Shell));
885: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetMatApply_C", PCShellSetMatApply_Shell));
886: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricLeft_C", PCShellSetApplySymmetricLeft_Shell));
887: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplySymmetricRight_C", PCShellSetApplySymmetricRight_Shell));
888: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyBA_C", PCShellSetApplyBA_Shell));
889: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPreSolve_C", PCShellSetPreSolve_Shell));
890: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetPostSolve_C", PCShellSetPostSolve_Shell));
891: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetView_C", PCShellSetView_Shell));
892: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyTranspose_C", PCShellSetApplyTranspose_Shell));
893: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetName_C", PCShellSetName_Shell));
894: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellGetName_C", PCShellGetName_Shell));
895: PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCShellSetApplyRichardson_C", PCShellSetApplyRichardson_Shell));
896: PetscFunctionReturn(PETSC_SUCCESS);
897: }