Actual source code: shell.c
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: This provides a simple shell interface for programmers to create
12: their own spectral transformations without writing much interface code
13: */
15: #include <slepc/private/stimpl.h>
17: typedef struct {
18: void *ctx; /* user-provided context */
19: STShellDestroyFn *destroy; /* context destroy */
20: STShellApplyFn *apply;
21: STShellApplyTransposeFn *applytrans;
22: STShellApplyHermitianTransposeFn *applyhermtrans;
23: STShellBackTransformFn *backtransform;
24: } ST_SHELL;
26: static PetscErrorCode STShellGetContext_Shell(ST st,PetscCtxRt ctx)
27: {
28: ST_SHELL *shell = (ST_SHELL*)st->data;
30: PetscFunctionBegin;
31: *(void**)ctx = shell->ctx;
32: PetscFunctionReturn(PETSC_SUCCESS);
33: }
35: /*@C
36: STShellGetContext - Returns the user-provided context associated with an `STSHELL`.
38: Not Collective
40: Input Parameter:
41: . st - spectral transformation context
43: Output Parameter:
44: . ctx - the user provided context
46: Level: advanced
48: .seealso: [](ch:st), `STSHELL`, `STShellSetContext()`
49: @*/
50: PetscErrorCode STShellGetContext(ST st,PetscCtxRt ctx)
51: {
52: PetscFunctionBegin;
54: PetscAssertPointer(ctx,2);
55: PetscUseMethod(st,"STShellGetContext_C",(ST,PetscCtxRt),(st,ctx));
56: PetscFunctionReturn(PETSC_SUCCESS);
57: }
59: static PetscErrorCode STShellSetContext_Shell(ST st,PetscCtx ctx)
60: {
61: ST_SHELL *shell = (ST_SHELL*)st->data;
63: PetscFunctionBegin;
64: if (shell->destroy) PetscCall((*shell->destroy)(st));
65: shell->ctx = ctx;
66: PetscFunctionReturn(PETSC_SUCCESS);
67: }
69: /*@
70: STShellSetContext - Sets the user-defined context for an `STSHELL`.
72: Logically Collective
74: Input Parameters:
75: + st - the shell `ST`
76: - ctx - the context
78: Level: advanced
80: Fortran Notes:
81: To use this from Fortran you must write a Fortran interface definition
82: for this function that tells Fortran the Fortran derived data type that
83: you are passing in as the `ctx` argument.
85: .seealso: [](ch:st), `STSHELL`, `STShellGetContext()`, `STShellSetDestroy()`
86: @*/
87: PetscErrorCode STShellSetContext(ST st,PetscCtx ctx)
88: {
89: PetscFunctionBegin;
91: PetscTryMethod(st,"STShellSetContext_C",(ST,PetscCtx),(st,ctx));
92: PetscFunctionReturn(PETSC_SUCCESS);
93: }
95: static PetscErrorCode STShellSetDestroy_Shell(ST st,STShellDestroyFn *destroy)
96: {
97: ST_SHELL *shell = (ST_SHELL*)st->data;
99: PetscFunctionBegin;
100: shell->destroy = destroy;
101: PetscFunctionReturn(PETSC_SUCCESS);
102: }
104: /*@C
105: STShellSetDestroy - Set a context destroy function for the `STSHELL` context.
107: Logically Collective
109: Input Parameters:
110: + st - the shell `ST`
111: - destroy - context destroy function, see `STShellDestroyFn` for its calling sequence
113: Level: advanced
115: .seealso: [](ch:st), `STSHELL`, `STShellSetContext()`
116: @*/
117: PetscErrorCode STShellSetDestroy(ST st,STShellDestroyFn *destroy)
118: {
119: PetscFunctionBegin;
121: PetscTryMethod(st,"STShellSetDestroy_C",(ST,STShellDestroyFn*),(st,destroy));
122: PetscFunctionReturn(PETSC_SUCCESS);
123: }
125: static PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
126: {
127: ST_SHELL *shell = (ST_SHELL*)st->data;
128: PetscObjectState instate,outstate;
130: PetscFunctionBegin;
131: PetscCheck(shell->apply,PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No apply() routine provided to Shell ST");
132: PetscCall(VecGetState(y,&instate));
133: PetscCallBack("STSHELL user function apply()",(*shell->apply)(st,x,y));
134: PetscCall(VecGetState(y,&outstate));
135: if (instate == outstate) {
136: /* user forgot to increase the state of the output vector */
137: PetscCall(PetscObjectStateIncrease((PetscObject)y));
138: }
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: static PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
143: {
144: ST_SHELL *shell = (ST_SHELL*)st->data;
145: PetscObjectState instate,outstate;
147: PetscFunctionBegin;
148: PetscCheck(shell->applytrans,PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No applytrans() routine provided to Shell ST");
149: PetscCall(VecGetState(y,&instate));
150: PetscCallBack("STSHELL user function applytrans()",(*shell->applytrans)(st,x,y));
151: PetscCall(VecGetState(y,&outstate));
152: if (instate == outstate) {
153: /* user forgot to increase the state of the output vector */
154: PetscCall(PetscObjectStateIncrease((PetscObject)y));
155: }
156: PetscFunctionReturn(PETSC_SUCCESS);
157: }
159: #if defined(PETSC_USE_COMPLEX)
160: static PetscErrorCode STApplyHermitianTranspose_Shell(ST st,Vec x,Vec y)
161: {
162: ST_SHELL *shell = (ST_SHELL*)st->data;
163: PetscObjectState instate,outstate;
164: Vec w;
166: PetscFunctionBegin;
167: if (shell->applyhermtrans) {
168: PetscCall(VecGetState(y,&instate));
169: PetscCallBack("STSHELL user function applyhermtrans()",(*shell->applyhermtrans)(st,x,y));
170: PetscCall(VecGetState(y,&outstate));
171: if (instate == outstate) {
172: /* user forgot to increase the state of the output vector */
173: PetscCall(PetscObjectStateIncrease((PetscObject)y));
174: }
175: } else {
176: PetscCall(VecDuplicate(x,&w));
177: PetscCall(VecCopy(x,w));
178: PetscCall(VecConjugate(w));
179: PetscCall(STApplyTranspose_Shell(st,w,y));
180: PetscCall(VecDestroy(&w));
181: PetscCall(VecConjugate(y));
182: }
183: PetscFunctionReturn(PETSC_SUCCESS);
184: }
185: #endif
187: static PetscErrorCode STBackTransform_Shell(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi)
188: {
189: ST_SHELL *shell = (ST_SHELL*)st->data;
191: PetscFunctionBegin;
192: if (shell->backtransform) PetscCallBack("STSHELL user function backtransform()",(*shell->backtransform)(st,n,eigr,eigi));
193: PetscFunctionReturn(PETSC_SUCCESS);
194: }
196: /*
197: STIsInjective_Shell - Check if the user has provided the backtransform operation.
198: */
199: PetscErrorCode STIsInjective_Shell(ST st,PetscBool* is)
200: {
201: ST_SHELL *shell = (ST_SHELL*)st->data;
203: PetscFunctionBegin;
204: *is = shell->backtransform? PETSC_TRUE: PETSC_FALSE;
205: PetscFunctionReturn(PETSC_SUCCESS);
206: }
208: static PetscErrorCode STDestroy_Shell(ST st)
209: {
210: ST_SHELL *shell = (ST_SHELL*)st->data;
212: PetscFunctionBegin;
213: if (shell->destroy) PetscCall((*shell->destroy)(st));
214: PetscCall(PetscFree(st->data));
215: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellGetContext_C",NULL));
216: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetContext_C",NULL));
217: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetDestroy_C",NULL));
218: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",NULL));
219: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",NULL));
220: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyHermitianTranspose_C",NULL));
221: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",NULL));
222: PetscFunctionReturn(PETSC_SUCCESS);
223: }
225: static PetscErrorCode STShellSetApply_Shell(ST st,STShellApplyFn *apply)
226: {
227: ST_SHELL *shell = (ST_SHELL*)st->data;
229: PetscFunctionBegin;
230: shell->apply = apply;
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: /*@C
235: STShellSetApply - Sets routine to use as the application of the
236: operator to a vector in the user-defined spectral transformation.
238: Logically Collective
240: Input Parameters:
241: + st - the spectral transformation context
242: - apply - the application-provided transformation routine
244: Level: advanced
246: .seealso: [](ch:st), `STSHELL`, `STShellSetBackTransform()`, `STShellSetApplyTranspose()`, `STShellSetApplyHermitianTranspose()`
247: @*/
248: PetscErrorCode STShellSetApply(ST st,STShellApplyFn *apply)
249: {
250: PetscFunctionBegin;
252: PetscTryMethod(st,"STShellSetApply_C",(ST,STShellApplyFn*),(st,apply));
253: PetscFunctionReturn(PETSC_SUCCESS);
254: }
256: static PetscErrorCode STShellSetApplyTranspose_Shell(ST st,STShellApplyTransposeFn *applytrans)
257: {
258: ST_SHELL *shell = (ST_SHELL*)st->data;
260: PetscFunctionBegin;
261: shell->applytrans = applytrans;
262: PetscFunctionReturn(PETSC_SUCCESS);
263: }
265: /*@C
266: STShellSetApplyTranspose - Sets routine to use as the application of the
267: transposed operator to a vector in the user-defined spectral transformation.
269: Logically Collective
271: Input Parameters:
272: + st - the spectral transformation context
273: - applytrans - the application-provided transformation routine
275: Level: advanced
277: .seealso: [](ch:st), `STSHELL`, `STShellSetApply()`, `STShellSetBackTransform()`
278: @*/
279: PetscErrorCode STShellSetApplyTranspose(ST st,STShellApplyTransposeFn *applytrans)
280: {
281: PetscFunctionBegin;
283: PetscTryMethod(st,"STShellSetApplyTranspose_C",(ST,STShellApplyTransposeFn*),(st,applytrans));
284: PetscFunctionReturn(PETSC_SUCCESS);
285: }
287: #if defined(PETSC_USE_COMPLEX)
288: static PetscErrorCode STShellSetApplyHermitianTranspose_Shell(ST st,STShellApplyHermitianTransposeFn *applyhermtrans)
289: {
290: ST_SHELL *shell = (ST_SHELL*)st->data;
292: PetscFunctionBegin;
293: shell->applyhermtrans = applyhermtrans;
294: PetscFunctionReturn(PETSC_SUCCESS);
295: }
296: #endif
298: /*@C
299: STShellSetApplyHermitianTranspose - Sets routine to use as the application of the
300: conjugate-transposed operator to a vector in the user-defined spectral transformation.
302: Logically Collective
304: Input Parameters:
305: + st - the spectral transformation context
306: - applyhermtrans - the application-provided transformation routine
308: Note:
309: If configured with real scalars, this function has the same effect as `STShellSetApplyTranspose()`,
310: so no need to call both.
312: Level: advanced
314: .seealso: [](ch:st), `STSHELL`, `STShellSetApply()`, `STShellSetApplyTranspose()`, `STShellSetBackTransform()`
315: @*/
316: PetscErrorCode STShellSetApplyHermitianTranspose(ST st,STShellApplyHermitianTransposeFn *applyhermtrans)
317: {
318: PetscFunctionBegin;
320: PetscTryMethod(st,"STShellSetApplyHermitianTranspose_C",(ST,STShellApplyHermitianTransposeFn*),(st,applyhermtrans));
321: PetscFunctionReturn(PETSC_SUCCESS);
322: }
324: static PetscErrorCode STShellSetBackTransform_Shell(ST st,STShellBackTransformFn *backtr)
325: {
326: ST_SHELL *shell = (ST_SHELL*)st->data;
328: PetscFunctionBegin;
329: shell->backtransform = backtr;
330: PetscFunctionReturn(PETSC_SUCCESS);
331: }
333: /*@C
334: STShellSetBackTransform - Sets the routine to be called after the
335: eigensolution process has finished in order to transform back the
336: computed eigenvalues.
338: Logically Collective
340: Input Parameters:
341: + st - the spectral transformation context
342: - backtr - the application-provided backtransform routine
344: Level: advanced
346: .seealso: [](ch:st), `STSHELL`, `STShellSetApply()`, `STShellSetApplyTranspose()`
347: @*/
348: PetscErrorCode STShellSetBackTransform(ST st,STShellBackTransformFn *backtr)
349: {
350: PetscFunctionBegin;
352: PetscTryMethod(st,"STShellSetBackTransform_C",(ST,STShellBackTransformFn*),(st,backtr));
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: /*MC
357: STSHELL - STSHELL = "shell" - User-defined spectral transformation via callback
358: functions for the application of the operator to a vector and (optionally) the
359: backtransform operation.
361: Level: beginner
363: Note:
364: In order to define a `shell` spectral transformation, the user has to provide
365: the `apply` operation via `STShellSetApply()` and related functions, and
366: optionally a `backtransform` operation via `STShellSetBackTransform()`, and
367: in some cases a user-defined context containing relevant data via
368: `STShellSetContext()`.
370: .seealso: [](ch:st), `ST`, `STType`, `STSetType()`, `STShellSetApply()`, `STShellSetBackTransform()`, `STShellSetContext()`
371: M*/
373: SLEPC_EXTERN PetscErrorCode STCreate_Shell(ST st)
374: {
375: ST_SHELL *ctx;
377: PetscFunctionBegin;
378: PetscCall(PetscNew(&ctx));
379: st->data = (void*)ctx;
381: st->usesksp = PETSC_FALSE;
383: st->ops->apply = STApply_Shell;
384: st->ops->applytrans = STApplyTranspose_Shell;
385: #if defined(PETSC_USE_COMPLEX)
386: st->ops->applyhermtrans = STApplyHermitianTranspose_Shell;
387: #else
388: st->ops->applyhermtrans = STApplyTranspose_Shell;
389: #endif
390: st->ops->backtransform = STBackTransform_Shell;
391: st->ops->destroy = STDestroy_Shell;
393: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellGetContext_C",STShellGetContext_Shell));
394: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetContext_C",STShellSetContext_Shell));
395: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetDestroy_C",STShellSetDestroy_Shell));
396: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",STShellSetApply_Shell));
397: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",STShellSetApplyTranspose_Shell));
398: #if defined(PETSC_USE_COMPLEX)
399: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyHermitianTranspose_C",STShellSetApplyHermitianTranspose_Shell));
400: #else
401: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyHermitianTranspose_C",STShellSetApplyTranspose_Shell));
402: #endif
403: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",STShellSetBackTransform_Shell));
404: PetscFunctionReturn(PETSC_SUCCESS);
405: }