Actual source code: shell.c
slepc-3.6.1 2015-09-03
1: /*
2: This provides a simple shell interface for programmers to
3: create their own spectral transformations without writing much
4: interface code.
6: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
7: SLEPc - Scalable Library for Eigenvalue Problem Computations
8: Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
10: This file is part of SLEPc.
12: SLEPc is free software: you can redistribute it and/or modify it under the
13: terms of version 3 of the GNU Lesser General Public License as published by
14: the Free Software Foundation.
16: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
17: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
19: more details.
21: You should have received a copy of the GNU Lesser General Public License
22: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
23: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
24: */
26: #include <slepc/private/stimpl.h> /*I "slepcst.h" I*/
28: typedef struct {
29: void *ctx; /* user provided context */
30: PetscErrorCode (*apply)(ST,Vec,Vec);
31: PetscErrorCode (*applytrans)(ST,Vec,Vec);
32: PetscErrorCode (*backtransform)(ST,PetscInt n,PetscScalar*,PetscScalar*);
33: } ST_SHELL;
37: /*@C
38: STShellGetContext - Returns the user-provided context associated with a shell ST
40: Not Collective
42: Input Parameter:
43: . st - spectral transformation context
45: Output Parameter:
46: . ctx - the user provided context
48: Level: advanced
50: Notes:
51: This routine is intended for use within various shell routines
53: .seealso: STShellSetContext()
54: @*/
55: PetscErrorCode STShellGetContext(ST st,void **ctx)
56: {
58: PetscBool flg;
63: PetscObjectTypeCompare((PetscObject)st,STSHELL,&flg);
64: if (!flg) *ctx = 0;
65: else *ctx = ((ST_SHELL*)(st->data))->ctx;
66: return(0);
67: }
71: /*@
72: STShellSetContext - Sets the context for a shell ST
74: Logically Collective on ST
76: Input Parameters:
77: + st - the shell ST
78: - ctx - the context
80: Level: advanced
82: Fortran Notes: The context can only be an integer or a PetscObject;
83: unfortunately it cannot be a Fortran array or derived type.
85: .seealso: STShellGetContext()
86: @*/
87: PetscErrorCode STShellSetContext(ST st,void *ctx)
88: {
89: ST_SHELL *shell = (ST_SHELL*)st->data;
91: PetscBool flg;
95: PetscObjectTypeCompare((PetscObject)st,STSHELL,&flg);
96: if (flg) {
97: shell->ctx = ctx;
98: }
99: return(0);
100: }
104: PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
105: {
107: ST_SHELL *shell = (ST_SHELL*)st->data;
110: if (!shell->apply) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No apply() routine provided to Shell ST");
111: PetscStackPush("STSHELL apply() user function");
112: CHKMEMQ;
113: (*shell->apply)(st,x,y);
114: CHKMEMQ;
115: PetscStackPop;
116: return(0);
117: }
121: PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
122: {
124: ST_SHELL *shell = (ST_SHELL*)st->data;
127: if (!shell->applytrans) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No applytranspose() routine provided to Shell ST");
128: PetscStackPush("STSHELL applytranspose() user function");
129: CHKMEMQ;
130: (*shell->applytrans)(st,x,y);
131: CHKMEMQ;
132: PetscStackPop;
133: return(0);
134: }
138: PetscErrorCode STBackTransform_Shell(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi)
139: {
141: ST_SHELL *shell = (ST_SHELL*)st->data;
144: if (shell->backtransform) {
145: PetscStackPush("STSHELL backtransform() user function");
146: CHKMEMQ;
147: (*shell->backtransform)(st,n,eigr,eigi);
148: CHKMEMQ;
149: PetscStackPop;
150: }
151: return(0);
152: }
156: PetscErrorCode STDestroy_Shell(ST st)
157: {
161: PetscFree(st->data);
162: PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",NULL);
163: PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",NULL);
164: PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",NULL);
165: return(0);
166: }
170: static PetscErrorCode STShellSetApply_Shell(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
171: {
172: ST_SHELL *shell = (ST_SHELL*)st->data;
175: shell->apply = apply;
176: return(0);
177: }
181: static PetscErrorCode STShellSetApplyTranspose_Shell(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
182: {
183: ST_SHELL *shell = (ST_SHELL*)st->data;
186: shell->applytrans = applytrans;
187: return(0);
188: }
192: static PetscErrorCode STShellSetBackTransform_Shell(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
193: {
194: ST_SHELL *shell = (ST_SHELL*)st->data;
197: shell->backtransform = backtr;
198: return(0);
199: }
203: /*@C
204: STShellSetApply - Sets routine to use as the application of the
205: operator to a vector in the user-defined spectral transformation.
207: Logically Collective on ST
209: Input Parameters:
210: + st - the spectral transformation context
211: - apply - the application-provided transformation routine
213: Calling sequence of apply:
214: .vb
215: PetscErrorCode apply (ST st,Vec xin,Vec xout)
216: .ve
218: + st - the spectral transformation context
219: . xin - input vector
220: - xout - output vector
222: Level: advanced
224: .seealso: STShellSetBackTransform(), STShellSetApplyTranspose()
225: @*/
226: PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
227: {
232: PetscTryMethod(st,"STShellSetApply_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,apply));
233: return(0);
234: }
238: /*@C
239: STShellSetApplyTranspose - Sets routine to use as the application of the
240: transposed operator to a vector in the user-defined spectral transformation.
242: Logically Collective on ST
244: Input Parameters:
245: + st - the spectral transformation context
246: - applytrans - the application-provided transformation routine
248: Calling sequence of apply:
249: .vb
250: PetscErrorCode applytrans (ST st,Vec xin,Vec xout)
251: .ve
253: + st - the spectral transformation context
254: . xin - input vector
255: - xout - output vector
257: Level: advanced
259: .seealso: STShellSetApply(), STShellSetBackTransform()
260: @*/
261: PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
262: {
267: PetscTryMethod(st,"STShellSetApplyTranspose_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,applytrans));
268: return(0);
269: }
273: /*@C
274: STShellSetBackTransform - Sets the routine to be called after the
275: eigensolution process has finished in order to transform back the
276: computed eigenvalues.
278: Logically Collective on ST
280: Input Parameters:
281: + st - the spectral transformation context
282: - backtr - the application-provided backtransform routine
284: Calling sequence of backtr:
285: .vb
286: PetscErrorCode backtr(ST st,PetscScalar *eigr,PetscScalar *eigi)
287: .ve
289: + st - the spectral transformation context
290: . eigr - pointer ot the real part of the eigenvalue to transform back
291: - eigi - pointer ot the imaginary part
293: Level: advanced
295: .seealso: STShellSetApply(), STShellSetApplyTranspose()
296: @*/
297: PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
298: {
303: PetscTryMethod(st,"STShellSetBackTransform_C",(ST,PetscErrorCode (*)(ST,PetscInt,PetscScalar*,PetscScalar*)),(st,backtr));
304: return(0);
305: }
309: PetscErrorCode STSetFromOptions_Shell(PetscOptions *PetscOptionsObject,ST st)
310: {
312: PC pc;
313: PCType pctype;
314: KSPType ksptype;
317: if (!st->ksp) { STGetKSP(st,&st->ksp); }
318: KSPGetPC(st->ksp,&pc);
319: KSPGetType(st->ksp,&ksptype);
320: PCGetType(pc,&pctype);
321: if (!pctype && !ksptype) {
322: if (st->shift_matrix == ST_MATMODE_SHELL) {
323: /* in shell mode use GMRES with Jacobi as the default */
324: KSPSetType(st->ksp,KSPGMRES);
325: PCSetType(pc,PCJACOBI);
326: } else {
327: /* use direct solver as default */
328: KSPSetType(st->ksp,KSPPREONLY);
329: PCSetType(pc,PCLU);
330: }
331: }
332: return(0);
333: }
335: /*MC
336: STSHELL - Creates a new spectral transformation class.
337: This is intended to provide a simple class to use with EPS.
338: You should not use this if you plan to make a complete class.
340: Level: advanced
342: Usage:
343: $ PetscErrorCode (*apply)(void*,Vec,Vec);
344: $ PetscErrorCode (*applytrans)(void*,Vec,Vec);
345: $ PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*);
346: $ STCreate(comm,&st);
347: $ STSetType(st,STSHELL);
348: $ STShellSetApply(st,apply);
349: $ STShellSetApplyTranspose(st,applytrans);
350: $ STShellSetBackTransform(st,backtr); (optional)
352: M*/
356: PETSC_EXTERN PetscErrorCode STCreate_Shell(ST st)
357: {
359: ST_SHELL *ctx;
362: PetscNewLog(st,&ctx);
363: st->data = (void*)ctx;
365: st->ops->apply = STApply_Shell;
366: st->ops->applytrans = STApplyTranspose_Shell;
367: st->ops->backtransform = STBackTransform_Shell;
368: st->ops->setfromoptions = STSetFromOptions_Shell;
369: st->ops->destroy = STDestroy_Shell;
370: PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",STShellSetApply_Shell);
371: PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",STShellSetApplyTranspose_Shell);
372: PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",STShellSetBackTransform_Shell);
373: return(0);
374: }