Actual source code: shell.c
slepc-3.17.2 2022-08-09
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 (*backtransform)(ST,PetscInt n,PetscScalar*,PetscScalar*);
22: } ST_SHELL;
24: /*@C
25: STShellGetContext - Returns the user-provided context associated with a shell ST
27: Not Collective
29: Input Parameter:
30: . st - spectral transformation context
32: Output Parameter:
33: . ctx - the user provided context
35: Level: advanced
37: Notes:
38: This routine is intended for use within various shell routines
40: .seealso: STShellSetContext()
41: @*/
42: PetscErrorCode STShellGetContext(ST st,void *ctx)
43: {
44: PetscBool flg;
48: PetscObjectTypeCompare((PetscObject)st,STSHELL,&flg);
49: if (!flg) *(void**)ctx = NULL;
50: else *(void**)ctx = ((ST_SHELL*)(st->data))->ctx;
51: PetscFunctionReturn(0);
52: }
54: /*@
55: STShellSetContext - Sets the context for a shell ST
57: Logically Collective on st
59: Input Parameters:
60: + st - the shell ST
61: - ctx - the context
63: Level: advanced
65: Fortran Notes:
66: To use this from Fortran you must write a Fortran interface definition
67: for this function that tells Fortran the Fortran derived data type that
68: you are passing in as the ctx argument.
70: .seealso: STShellGetContext()
71: @*/
72: PetscErrorCode STShellSetContext(ST st,void *ctx)
73: {
74: ST_SHELL *shell = (ST_SHELL*)st->data;
75: PetscBool flg;
78: PetscObjectTypeCompare((PetscObject)st,STSHELL,&flg);
79: if (flg) shell->ctx = ctx;
80: PetscFunctionReturn(0);
81: }
83: PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
84: {
85: ST_SHELL *shell = (ST_SHELL*)st->data;
86: PetscObjectState instate,outstate;
89: PetscObjectStateGet((PetscObject)y,&instate);
90: PetscStackCall("STSHELL user function apply()",(*shell->apply)(st,x,y));
91: PetscObjectStateGet((PetscObject)y,&outstate);
92: if (instate == outstate) {
93: /* user forgot to increase the state of the output vector */
94: PetscObjectStateIncrease((PetscObject)y);
95: }
96: PetscFunctionReturn(0);
97: }
99: PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
100: {
101: ST_SHELL *shell = (ST_SHELL*)st->data;
102: PetscObjectState instate,outstate;
105: PetscObjectStateGet((PetscObject)y,&instate);
106: PetscStackCall("STSHELL user function applytrans()",(*shell->applytrans)(st,x,y));
107: PetscObjectStateGet((PetscObject)y,&outstate);
108: if (instate == outstate) {
109: /* user forgot to increase the state of the output vector */
110: PetscObjectStateIncrease((PetscObject)y);
111: }
112: PetscFunctionReturn(0);
113: }
115: PetscErrorCode STBackTransform_Shell(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi)
116: {
117: ST_SHELL *shell = (ST_SHELL*)st->data;
119: if (shell->backtransform) PetscStackCall("STSHELL user function backtransform()",(*shell->backtransform)(st,n,eigr,eigi));
120: PetscFunctionReturn(0);
121: }
123: /*
124: STIsInjective_Shell - Check if the user has provided the backtransform operation.
125: */
126: PetscErrorCode STIsInjective_Shell(ST st,PetscBool* is)
127: {
128: ST_SHELL *shell = (ST_SHELL*)st->data;
130: *is = shell->backtransform? PETSC_TRUE: PETSC_FALSE;
131: PetscFunctionReturn(0);
132: }
134: PetscErrorCode STDestroy_Shell(ST st)
135: {
136: PetscFree(st->data);
137: PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",NULL);
138: PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",NULL);
139: PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",NULL);
140: PetscFunctionReturn(0);
141: }
143: static PetscErrorCode STShellSetApply_Shell(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
144: {
145: ST_SHELL *shell = (ST_SHELL*)st->data;
147: shell->apply = apply;
148: PetscFunctionReturn(0);
149: }
151: /*@C
152: STShellSetApply - Sets routine to use as the application of the
153: operator to a vector in the user-defined spectral transformation.
155: Logically Collective on st
157: Input Parameters:
158: + st - the spectral transformation context
159: - apply - the application-provided transformation routine
161: Calling sequence of apply:
162: $ PetscErrorCode apply(ST st,Vec xin,Vec xout)
164: + st - the spectral transformation context
165: . xin - input vector
166: - xout - output vector
168: Level: advanced
170: .seealso: STShellSetBackTransform(), STShellSetApplyTranspose()
171: @*/
172: PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
173: {
175: PetscTryMethod(st,"STShellSetApply_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,apply));
176: PetscFunctionReturn(0);
177: }
179: static PetscErrorCode STShellSetApplyTranspose_Shell(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
180: {
181: ST_SHELL *shell = (ST_SHELL*)st->data;
183: shell->applytrans = applytrans;
184: PetscFunctionReturn(0);
185: }
187: /*@C
188: STShellSetApplyTranspose - Sets routine to use as the application of the
189: transposed operator to a vector in the user-defined spectral transformation.
191: Logically Collective on st
193: Input Parameters:
194: + st - the spectral transformation context
195: - applytrans - the application-provided transformation routine
197: Calling sequence of applytrans:
198: $ PetscErrorCode applytrans(ST st,Vec xin,Vec xout)
200: + st - the spectral transformation context
201: . xin - input vector
202: - xout - output vector
204: Level: advanced
206: .seealso: STShellSetApply(), STShellSetBackTransform()
207: @*/
208: PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
209: {
211: PetscTryMethod(st,"STShellSetApplyTranspose_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,applytrans));
212: PetscFunctionReturn(0);
213: }
215: static PetscErrorCode STShellSetBackTransform_Shell(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
216: {
217: ST_SHELL *shell = (ST_SHELL*)st->data;
219: shell->backtransform = backtr;
220: PetscFunctionReturn(0);
221: }
223: /*@C
224: STShellSetBackTransform - Sets the routine to be called after the
225: eigensolution process has finished in order to transform back the
226: computed eigenvalues.
228: Logically Collective on st
230: Input Parameters:
231: + st - the spectral transformation context
232: - backtr - the application-provided backtransform routine
234: Calling sequence of backtr:
235: $ PetscErrorCode backtr(ST st,PetscScalar *eigr,PetscScalar *eigi)
237: + st - the spectral transformation context
238: . eigr - pointer ot the real part of the eigenvalue to transform back
239: - eigi - pointer ot the imaginary part
241: Level: advanced
243: .seealso: STShellSetApply(), STShellSetApplyTranspose()
244: @*/
245: PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
246: {
248: PetscTryMethod(st,"STShellSetBackTransform_C",(ST,PetscErrorCode (*)(ST,PetscInt,PetscScalar*,PetscScalar*)),(st,backtr));
249: PetscFunctionReturn(0);
250: }
252: /*MC
253: STSHELL - User-defined spectral transformation via callback functions
254: for the application of the operator to a vector and (optionally) the
255: backtransform operation.
257: Level: advanced
259: Usage:
260: $ extern PetscErrorCode (*apply)(void*,Vec,Vec);
261: $ extern PetscErrorCode (*applytrans)(void*,Vec,Vec);
262: $ extern PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*);
263: $
264: $ STCreate(comm,&st);
265: $ STSetType(st,STSHELL);
266: $ STShellSetContext(st,ctx);
267: $ STShellSetApply(st,apply);
268: $ STShellSetApplyTranspose(st,applytrans); (optional)
269: $ STShellSetBackTransform(st,backtr); (optional)
271: M*/
273: SLEPC_EXTERN PetscErrorCode STCreate_Shell(ST st)
274: {
275: ST_SHELL *ctx;
277: PetscNewLog(st,&ctx);
278: st->data = (void*)ctx;
280: st->usesksp = PETSC_FALSE;
282: st->ops->apply = STApply_Shell;
283: st->ops->applytrans = STApplyTranspose_Shell;
284: st->ops->backtransform = STBackTransform_Shell;
285: st->ops->destroy = STDestroy_Shell;
287: PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",STShellSetApply_Shell);
288: PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",STShellSetApplyTranspose_Shell);
289: PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",STShellSetBackTransform_Shell);
290: PetscFunctionReturn(0);
291: }