Actual source code: shell.c
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-2011, Universitat Politecnica de Valencia, Spain
10: This file is part of SLEPc.
11:
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 private/stimpl.h
29: typedef struct {
30: void *ctx; /* user provided context */
31: PetscErrorCode (*apply)(ST,Vec,Vec);
32: PetscErrorCode (*applytrans)(ST,Vec,Vec);
33: PetscErrorCode (*backtr)(ST,PetscInt n,PetscScalar*,PetscScalar*);
34: } ST_Shell;
39: /*@C
40: STShellGetContext - Returns the user-provided context associated with a shell ST
42: Not Collective
44: Input Parameter:
45: . st - spectral transformation context
47: Output Parameter:
48: . ctx - the user provided context
50: Level: advanced
52: Notes:
53: This routine is intended for use within various shell routines
54:
55: .seealso: STShellSetContext()
56: @*/
57: PetscErrorCode STShellGetContext(ST st,void **ctx)
58: {
60: PetscBool flg;
65: PetscTypeCompare((PetscObject)st,STSHELL,&flg);
66: if (!flg) *ctx = 0;
67: else *ctx = ((ST_Shell*)(st->data))->ctx;
68: return(0);
69: }
73: /*@
74: STShellSetContext - sets the context for a shell ST
76: Logically Collective on ST
78: Input Parameters:
79: + st - the shell ST
80: - ctx - the context
82: Level: advanced
84: Fortran Notes: The context can only be an integer or a PetscObject;
85: unfortunately it cannot be a Fortran array or derived type.
87: .seealso: STShellGetContext()
88: @*/
89: PetscErrorCode STShellSetContext(ST st,void *ctx)
90: {
91: ST_Shell *shell = (ST_Shell*)st->data;
93: PetscBool flg;
97: PetscTypeCompare((PetscObject)st,STSHELL,&flg);
98: if (flg) {
99: shell->ctx = ctx;
100: }
101: return(0);
102: }
106: PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
107: {
109: ST_Shell *shell = (ST_Shell*)st->data;
112: if (!shell->apply) SETERRQ(((PetscObject)st)->comm,PETSC_ERR_USER,"No apply() routine provided to Shell ST");
113: PetscStackPush("STSHELL apply() user function");
114: CHKMEMQ;
115: (*shell->apply)(st,x,y);
116: CHKMEMQ;
117: PetscStackPop;
118: return(0);
119: }
123: PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
124: {
126: ST_Shell *shell = (ST_Shell*)st->data;
129: if (!shell->applytrans) SETERRQ(((PetscObject)st)->comm,PETSC_ERR_USER,"No applytranspose() routine provided to Shell ST");
130: PetscStackPush("STSHELL applytranspose() user function");
131: CHKMEMQ;
132: (*shell->applytrans)(st,x,y);
133: CHKMEMQ;
134: PetscStackPop;
135: return(0);
136: }
140: PetscErrorCode STBackTransform_Shell(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi)
141: {
143: ST_Shell *shell = (ST_Shell*)st->data;
146: if (shell->backtr) {
147: PetscStackPush("STSHELL backtransform() user function");
148: CHKMEMQ;
149: (*shell->backtr)(st,n,eigr,eigi);
150: CHKMEMQ;
151: PetscStackPop;
152: }
153: return(0);
154: }
158: PetscErrorCode STDestroy_Shell(ST st)
159: {
163: PetscFree(st->data);
164: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApply_C","",PETSC_NULL);
165: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApplyTranspose_C","",PETSC_NULL);
166: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetBackTransform_C","",PETSC_NULL);
167: return(0);
168: }
173: PetscErrorCode STShellSetApply_Shell(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
174: {
175: ST_Shell *shell = (ST_Shell*)st->data;
178: shell->apply = apply;
179: return(0);
180: }
186: PetscErrorCode STShellSetApplyTranspose_Shell(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
187: {
188: ST_Shell *shell = (ST_Shell*)st->data;
191: shell->applytrans = applytrans;
192: return(0);
193: }
199: PetscErrorCode STShellSetBackTransform_Shell(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
200: {
201: ST_Shell *shell = (ST_Shell*)st->data;
204: shell->backtr = backtr;
205: return(0);
206: }
211: /*@C
212: STShellSetApply - Sets routine to use as the application of the
213: operator to a vector in the user-defined spectral transformation.
215: Logically Collective on ST
217: Input Parameters:
218: + st - the spectral transformation context
219: - apply - the application-provided transformation routine
221: Calling sequence of apply:
222: .vb
223: PetscErrorCode apply (ST st,Vec xin,Vec xout)
224: .ve
226: + st - the spectral transformation context
227: . xin - input vector
228: - xout - output vector
230: Level: developer
232: .seealso: STShellSetBackTransform(), STShellSetApplyTranspose()
233: @*/
234: PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
235: {
240: PetscTryMethod(st,"STShellSetApply_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,apply));
241: return(0);
242: }
246: /*@C
247: STShellSetApplyTranspose - Sets routine to use as the application of the
248: transposed operator to a vector in the user-defined spectral transformation.
250: Logically Collective on ST
252: Input Parameters:
253: + st - the spectral transformation context
254: - applytrans - the application-provided transformation routine
256: Calling sequence of apply:
257: .vb
258: PetscErrorCode applytrans (ST st,Vec xin,Vec xout)
259: .ve
261: + st - the spectral transformation context
262: . xin - input vector
263: - xout - output vector
265: Level: developer
267: .seealso: STShellSetApply(), STShellSetBackTransform()
268: @*/
269: PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
270: {
275: PetscTryMethod(st,"STShellSetApplyTranspose_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,applytrans));
276: return(0);
277: }
281: /*@C
282: STShellSetBackTransform - Sets the routine to be called after the
283: eigensolution process has finished in order to transform back the
284: computed eigenvalues.
286: Logically Collective on ST
288: Input Parameters:
289: + st - the spectral transformation context
290: - backtr - the application-provided backtransform routine
292: Calling sequence of backtr:
293: .vb
294: PetscErrorCode backtr (ST st,PetscScalar *eigr,PetscScalar *eigi)
295: .ve
297: + st - the spectral transformation context
298: . eigr - pointer ot the real part of the eigenvalue to transform back
299: - eigi - pointer ot the imaginary part
301: Level: developer
303: .seealso: STShellSetApply(), STShellSetApplyTranspose()
304: @*/
305: PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
306: {
311: PetscTryMethod(st,"STShellSetBackTransform_C",(ST,PetscErrorCode (*)(ST,PetscInt,PetscScalar*,PetscScalar*)),(st,backtr));
312: return(0);
313: }
317: PetscErrorCode STSetFromOptions_Shell(ST st)
318: {
320: PC pc;
321: const PCType pctype;
322: const KSPType ksptype;
325: if (!st->ksp) { STGetKSP(st,&st->ksp); }
326: KSPGetPC(st->ksp,&pc);
327: KSPGetType(st->ksp,&ksptype);
328: PCGetType(pc,&pctype);
329: if (!pctype && !ksptype) {
330: if (st->shift_matrix == ST_MATMODE_SHELL) {
331: /* in shell mode use GMRES with Jacobi as the default */
332: KSPSetType(st->ksp,KSPGMRES);
333: PCSetType(pc,PCJACOBI);
334: } else {
335: /* use direct solver as default */
336: KSPSetType(st->ksp,KSPPREONLY);
337: PCSetType(pc,PCREDUNDANT);
338: }
339: }
340: return(0);
341: }
343: /*MC
344: STSHELL - Creates a new spectral transformation class.
345: This is intended to provide a simple class to use with EPS.
346: You should not use this if you plan to make a complete class.
348: Level: advanced
350: Usage:
351: $ PetscErrorCode (*apply)(void*,Vec,Vec);
352: $ PetscErrorCode (*applytrans)(void*,Vec,Vec);
353: $ PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*);
354: $ STCreate(comm,&st);
355: $ STSetType(st,STSHELL);
356: $ STShellSetApply(st,apply);
357: $ STShellSetApplyTranspose(st,applytrans);
358: $ STShellSetBackTransform(st,backtr); (optional)
360: M*/
365: PetscErrorCode STCreate_Shell(ST st)
366: {
370: PetscNewLog(st,ST_Shell,&st->data);
371: st->ops->apply = STApply_Shell;
372: st->ops->applytrans = STApplyTranspose_Shell;
373: st->ops->backtr = STBackTransform_Shell;
374: st->ops->setfromoptions = STSetFromOptions_Shell;
375: st->ops->destroy = STDestroy_Shell;
376: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApply_C","STShellSetApply_Shell",STShellSetApply_Shell);
377: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApplyTranspose_C","STShellSetApplyTranspose_Shell",STShellSetApplyTranspose_Shell);
378: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetBackTransform_C","STShellSetBackTransform_Shell",STShellSetBackTransform_Shell);
379: return(0);
380: }