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: PetscErrorCode (*apply)(ST,Vec,Vec);
20: PetscErrorCode (*applytrans)(ST,Vec,Vec);
21: PetscErrorCode (*applyhermtrans)(ST,Vec,Vec);
22: PetscErrorCode (*backtransform)(ST,PetscInt n,PetscScalar*,PetscScalar*);
23: } ST_SHELL;
25: /*@C
26: STShellGetContext - Returns the user-provided context associated with an `STSHELL`.
28: Not Collective
30: Input Parameter:
31: . st - spectral transformation context
33: Output Parameter:
34: . ctx - the user provided context
36: Level: advanced
38: .seealso: [](ch:st), `STShellSetContext()`
39: @*/
40: PetscErrorCode STShellGetContext(ST st,void *ctx)
41: {
42: PetscBool flg;
44: PetscFunctionBegin;
46: PetscAssertPointer(ctx,2);
47: PetscCall(PetscObjectTypeCompare((PetscObject)st,STSHELL,&flg));
48: if (!flg) *(void**)ctx = NULL;
49: else *(void**)ctx = ((ST_SHELL*)st->data)->ctx;
50: PetscFunctionReturn(PETSC_SUCCESS);
51: }
53: /*@
54: STShellSetContext - Sets the user-defined context for an `STSHELL`.
56: Logically Collective
58: Input Parameters:
59: + st - the shell ST
60: - ctx - the context
62: Level: advanced
64: Fortran Notes:
65: To use this from Fortran you must write a Fortran interface definition
66: for this function that tells Fortran the Fortran derived data type that
67: you are passing in as the `ctx` argument.
69: .seealso: [](ch:st), `STShellGetContext()`
70: @*/
71: PetscErrorCode STShellSetContext(ST st,void *ctx)
72: {
73: ST_SHELL *shell = (ST_SHELL*)st->data;
74: PetscBool flg;
76: PetscFunctionBegin;
78: PetscCall(PetscObjectTypeCompare((PetscObject)st,STSHELL,&flg));
79: if (flg) shell->ctx = ctx;
80: PetscFunctionReturn(PETSC_SUCCESS);
81: }
83: static PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
84: {
85: ST_SHELL *shell = (ST_SHELL*)st->data;
86: PetscObjectState instate,outstate;
88: PetscFunctionBegin;
89: PetscCheck(shell->apply,PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No apply() routine provided to Shell ST");
90: PetscCall(VecGetState(y,&instate));
91: PetscCallBack("STSHELL user function apply()",(*shell->apply)(st,x,y));
92: PetscCall(VecGetState(y,&outstate));
93: if (instate == outstate) {
94: /* user forgot to increase the state of the output vector */
95: PetscCall(PetscObjectStateIncrease((PetscObject)y));
96: }
97: PetscFunctionReturn(PETSC_SUCCESS);
98: }
100: static PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
101: {
102: ST_SHELL *shell = (ST_SHELL*)st->data;
103: PetscObjectState instate,outstate;
105: PetscFunctionBegin;
106: PetscCheck(shell->applytrans,PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No applytrans() routine provided to Shell ST");
107: PetscCall(VecGetState(y,&instate));
108: PetscCallBack("STSHELL user function applytrans()",(*shell->applytrans)(st,x,y));
109: PetscCall(VecGetState(y,&outstate));
110: if (instate == outstate) {
111: /* user forgot to increase the state of the output vector */
112: PetscCall(PetscObjectStateIncrease((PetscObject)y));
113: }
114: PetscFunctionReturn(PETSC_SUCCESS);
115: }
117: #if defined(PETSC_USE_COMPLEX)
118: static PetscErrorCode STApplyHermitianTranspose_Shell(ST st,Vec x,Vec y)
119: {
120: ST_SHELL *shell = (ST_SHELL*)st->data;
121: PetscObjectState instate,outstate;
122: Vec w;
124: PetscFunctionBegin;
125: if (shell->applyhermtrans) {
126: PetscCall(VecGetState(y,&instate));
127: PetscCallBack("STSHELL user function applyhermtrans()",(*shell->applyhermtrans)(st,x,y));
128: PetscCall(VecGetState(y,&outstate));
129: if (instate == outstate) {
130: /* user forgot to increase the state of the output vector */
131: PetscCall(PetscObjectStateIncrease((PetscObject)y));
132: }
133: } else {
134: PetscCall(VecDuplicate(x,&w));
135: PetscCall(VecCopy(x,w));
136: PetscCall(VecConjugate(w));
137: PetscCall(STApplyTranspose_Shell(st,w,y));
138: PetscCall(VecDestroy(&w));
139: PetscCall(VecConjugate(y));
140: }
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
143: #endif
145: static PetscErrorCode STBackTransform_Shell(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi)
146: {
147: ST_SHELL *shell = (ST_SHELL*)st->data;
149: PetscFunctionBegin;
150: if (shell->backtransform) PetscCallBack("STSHELL user function backtransform()",(*shell->backtransform)(st,n,eigr,eigi));
151: PetscFunctionReturn(PETSC_SUCCESS);
152: }
154: /*
155: STIsInjective_Shell - Check if the user has provided the backtransform operation.
156: */
157: PetscErrorCode STIsInjective_Shell(ST st,PetscBool* is)
158: {
159: ST_SHELL *shell = (ST_SHELL*)st->data;
161: PetscFunctionBegin;
162: *is = shell->backtransform? PETSC_TRUE: PETSC_FALSE;
163: PetscFunctionReturn(PETSC_SUCCESS);
164: }
166: static PetscErrorCode STDestroy_Shell(ST st)
167: {
168: PetscFunctionBegin;
169: PetscCall(PetscFree(st->data));
170: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",NULL));
171: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",NULL));
172: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyHermitianTranspose_C",NULL));
173: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",NULL));
174: PetscFunctionReturn(PETSC_SUCCESS);
175: }
177: static PetscErrorCode STShellSetApply_Shell(ST st,STShellApplyFn *apply)
178: {
179: ST_SHELL *shell = (ST_SHELL*)st->data;
181: PetscFunctionBegin;
182: shell->apply = apply;
183: PetscFunctionReturn(PETSC_SUCCESS);
184: }
186: /*@C
187: STShellSetApply - Sets routine to use as the application of the
188: operator to a vector in the user-defined spectral transformation.
190: Logically Collective
192: Input Parameters:
193: + st - the spectral transformation context
194: - apply - the application-provided transformation routine
196: Level: advanced
198: .seealso: [](ch:st), `STShellSetBackTransform()`, `STShellSetApplyTranspose()`, `STShellSetApplyHermitianTranspose()`
199: @*/
200: PetscErrorCode STShellSetApply(ST st,STShellApplyFn *apply)
201: {
202: PetscFunctionBegin;
204: PetscTryMethod(st,"STShellSetApply_C",(ST,STShellApplyFn*),(st,apply));
205: PetscFunctionReturn(PETSC_SUCCESS);
206: }
208: static PetscErrorCode STShellSetApplyTranspose_Shell(ST st,STShellApplyTransposeFn *applytrans)
209: {
210: ST_SHELL *shell = (ST_SHELL*)st->data;
212: PetscFunctionBegin;
213: shell->applytrans = applytrans;
214: PetscFunctionReturn(PETSC_SUCCESS);
215: }
217: /*@C
218: STShellSetApplyTranspose - Sets routine to use as the application of the
219: transposed operator to a vector in the user-defined spectral transformation.
221: Logically Collective
223: Input Parameters:
224: + st - the spectral transformation context
225: - applytrans - the application-provided transformation routine
227: Level: advanced
229: .seealso: [](ch:st), `STShellSetApply()`, `STShellSetBackTransform()`
230: @*/
231: PetscErrorCode STShellSetApplyTranspose(ST st,STShellApplyTransposeFn *applytrans)
232: {
233: PetscFunctionBegin;
235: PetscTryMethod(st,"STShellSetApplyTranspose_C",(ST,STShellApplyTransposeFn*),(st,applytrans));
236: PetscFunctionReturn(PETSC_SUCCESS);
237: }
239: #if defined(PETSC_USE_COMPLEX)
240: static PetscErrorCode STShellSetApplyHermitianTranspose_Shell(ST st,STShellApplyHermitianTransposeFn *applyhermtrans)
241: {
242: ST_SHELL *shell = (ST_SHELL*)st->data;
244: PetscFunctionBegin;
245: shell->applyhermtrans = applyhermtrans;
246: PetscFunctionReturn(PETSC_SUCCESS);
247: }
248: #endif
250: /*@C
251: STShellSetApplyHermitianTranspose - Sets routine to use as the application of the
252: conjugate-transposed operator to a vector in the user-defined spectral transformation.
254: Logically Collective
256: Input Parameters:
257: + st - the spectral transformation context
258: - applyhermtrans - the application-provided transformation routine
260: Note:
261: If configured with real scalars, this function has the same effect as `STShellSetApplyTranspose()`,
262: so no need to call both.
264: Level: advanced
266: .seealso: [](ch:st), `STShellSetApply()`, `STShellSetApplyTranspose()`, `STShellSetBackTransform()`
267: @*/
268: PetscErrorCode STShellSetApplyHermitianTranspose(ST st,STShellApplyHermitianTransposeFn *applyhermtrans)
269: {
270: PetscFunctionBegin;
272: PetscTryMethod(st,"STShellSetApplyHermitianTranspose_C",(ST,STShellApplyHermitianTransposeFn*),(st,applyhermtrans));
273: PetscFunctionReturn(PETSC_SUCCESS);
274: }
276: static PetscErrorCode STShellSetBackTransform_Shell(ST st,STShellBackTransformFn *backtr)
277: {
278: ST_SHELL *shell = (ST_SHELL*)st->data;
280: PetscFunctionBegin;
281: shell->backtransform = backtr;
282: PetscFunctionReturn(PETSC_SUCCESS);
283: }
285: /*@C
286: STShellSetBackTransform - Sets the routine to be called after the
287: eigensolution process has finished in order to transform back the
288: computed eigenvalues.
290: Logically Collective
292: Input Parameters:
293: + st - the spectral transformation context
294: - backtr - the application-provided backtransform routine
296: Level: advanced
298: .seealso: [](ch:st), `STShellSetApply()`, `STShellSetApplyTranspose()`
299: @*/
300: PetscErrorCode STShellSetBackTransform(ST st,STShellBackTransformFn *backtr)
301: {
302: PetscFunctionBegin;
304: PetscTryMethod(st,"STShellSetBackTransform_C",(ST,STShellBackTransformFn*),(st,backtr));
305: PetscFunctionReturn(PETSC_SUCCESS);
306: }
308: /*MC
309: STSHELL - STSHELL = "shell" - User-defined spectral transformation via callback
310: functions for the application of the operator to a vector and (optionally) the
311: backtransform operation.
313: Level: beginner
315: Note:
316: In order to define a `shell` spectral transformation, the user has to provide
317: the `apply` operation via `STShellSetApply()` and related functions, and
318: optionally a `backtransform` operation via `STShellSetBackTransform()`, and
319: in some cases a user-defined context containing relevant data via
320: `STShellSetContext()`.
322: .seealso: [](ch:st), `ST`, `STType`, `STSetType()`, `STShellSetApply()`, `STShellSetBackTransform()`, `STShellSetContext()`
323: M*/
325: SLEPC_EXTERN PetscErrorCode STCreate_Shell(ST st)
326: {
327: ST_SHELL *ctx;
329: PetscFunctionBegin;
330: PetscCall(PetscNew(&ctx));
331: st->data = (void*)ctx;
333: st->usesksp = PETSC_FALSE;
335: st->ops->apply = STApply_Shell;
336: st->ops->applytrans = STApplyTranspose_Shell;
337: #if defined(PETSC_USE_COMPLEX)
338: st->ops->applyhermtrans = STApplyHermitianTranspose_Shell;
339: #else
340: st->ops->applyhermtrans = STApplyTranspose_Shell;
341: #endif
342: st->ops->backtransform = STBackTransform_Shell;
343: st->ops->destroy = STDestroy_Shell;
345: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",STShellSetApply_Shell));
346: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",STShellSetApplyTranspose_Shell));
347: #if defined(PETSC_USE_COMPLEX)
348: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyHermitianTranspose_C",STShellSetApplyHermitianTranspose_Shell));
349: #else
350: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyHermitianTranspose_C",STShellSetApplyTranspose_Shell));
351: #endif
352: PetscCall(PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",STShellSetBackTransform_Shell));
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }