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-2016, 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:
83: To use this from Fortran you must write a Fortran interface definition
84: for this function that tells Fortran the Fortran derived data type that
85: you are passing in as the ctx argument.
87: .seealso: STShellGetContext()
88: @*/
89: PetscErrorCode STShellSetContext(ST st,void *ctx) 90: {
91: ST_SHELL *shell = (ST_SHELL*)st->data;
93: PetscBool flg;
97: PetscObjectTypeCompare((PetscObject)st,STSHELL,&flg);
98: if (flg) shell->ctx = ctx;
99: return(0);
100: }
104: PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)105: {
106: PetscErrorCode ierr;
107: ST_SHELL *shell = (ST_SHELL*)st->data;
108: PetscObjectState instate,outstate;
111: if (!shell->apply) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No apply() routine provided to Shell ST");
112: PetscObjectStateGet((PetscObject)y,&instate);
113: PetscStackCall("STSHELL user function apply()",(*shell->apply)(st,x,y);CHKERRQ(ierr));
114: PetscObjectStateGet((PetscObject)y,&outstate);
115: if (instate == outstate) {
116: /* user forgot to increase the state of the output vector */
117: PetscObjectStateIncrease((PetscObject)y);
118: }
119: return(0);
120: }
124: PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)125: {
127: ST_SHELL *shell = (ST_SHELL*)st->data;
128: PetscObjectState instate,outstate;
131: if (!shell->applytrans) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_USER,"No applytranspose() routine provided to Shell ST");
132: PetscObjectStateGet((PetscObject)y,&instate);
133: PetscStackCall("STSHELL user function applytrans()",(*shell->applytrans)(st,x,y);CHKERRQ(ierr));
134: PetscObjectStateGet((PetscObject)y,&outstate);
135: if (instate == outstate) {
136: /* user forgot to increase the state of the output vector */
137: PetscObjectStateIncrease((PetscObject)y);
138: }
139: return(0);
140: }
144: PetscErrorCode STBackTransform_Shell(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi)145: {
147: ST_SHELL *shell = (ST_SHELL*)st->data;
150: if (shell->backtransform) PetscStackCall("STSHELL user function backtransform()",(*shell->backtransform)(st,n,eigr,eigi);CHKERRQ(ierr));
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 ST209: Input Parameters:
210: + st - the spectral transformation context
211: - apply - the application-provided transformation routine
213: Calling sequence of apply:
214: $ apply(ST st,Vec xin,Vec xout)
216: + st - the spectral transformation context
217: . xin - input vector
218: - xout - output vector
220: Level: advanced
222: .seealso: STShellSetBackTransform(), STShellSetApplyTranspose()
223: @*/
224: PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))225: {
230: PetscTryMethod(st,"STShellSetApply_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,apply));
231: return(0);
232: }
236: /*@C
237: STShellSetApplyTranspose - Sets routine to use as the application of the
238: transposed operator to a vector in the user-defined spectral transformation.
240: Logically Collective on ST242: Input Parameters:
243: + st - the spectral transformation context
244: - applytrans - the application-provided transformation routine
246: Calling sequence of applytrans:
247: $ applytrans(ST st,Vec xin,Vec xout)
249: + st - the spectral transformation context
250: . xin - input vector
251: - xout - output vector
253: Level: advanced
255: .seealso: STShellSetApply(), STShellSetBackTransform()
256: @*/
257: PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))258: {
263: PetscTryMethod(st,"STShellSetApplyTranspose_C",(ST,PetscErrorCode (*)(ST,Vec,Vec)),(st,applytrans));
264: return(0);
265: }
269: /*@C
270: STShellSetBackTransform - Sets the routine to be called after the
271: eigensolution process has finished in order to transform back the
272: computed eigenvalues.
274: Logically Collective on ST276: Input Parameters:
277: + st - the spectral transformation context
278: - backtr - the application-provided backtransform routine
280: Calling sequence of backtr:
281: $ backtr(ST st,PetscScalar *eigr,PetscScalar *eigi)
283: + st - the spectral transformation context
284: . eigr - pointer ot the real part of the eigenvalue to transform back
285: - eigi - pointer ot the imaginary part
287: Level: advanced
289: .seealso: STShellSetApply(), STShellSetApplyTranspose()
290: @*/
291: PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))292: {
297: PetscTryMethod(st,"STShellSetBackTransform_C",(ST,PetscErrorCode (*)(ST,PetscInt,PetscScalar*,PetscScalar*)),(st,backtr));
298: return(0);
299: }
303: PetscErrorCode STSetFromOptions_Shell(PetscOptionItems *PetscOptionsObject,ST st)304: {
306: PC pc;
307: PCType pctype;
308: KSPType ksptype;
311: if (!st->ksp) { STGetKSP(st,&st->ksp); }
312: KSPGetPC(st->ksp,&pc);
313: KSPGetType(st->ksp,&ksptype);
314: PCGetType(pc,&pctype);
315: if (!pctype && !ksptype) {
316: if (st->shift_matrix == ST_MATMODE_SHELL) {
317: /* in shell mode use GMRES with Jacobi as the default */
318: KSPSetType(st->ksp,KSPGMRES);
319: PCSetType(pc,PCJACOBI);
320: } else {
321: /* use direct solver as default */
322: KSPSetType(st->ksp,KSPPREONLY);
323: PCSetType(pc,PCLU);
324: }
325: }
326: return(0);
327: }
329: /*MC
330: STSHELL - Creates a new spectral transformation class.
331: This is intended to provide a simple class to use with EPS.
332: You should not use this if you plan to make a complete class.
334: Level: advanced
336: Usage:
337: $ PetscErrorCode (*apply)(void*,Vec,Vec);
338: $ PetscErrorCode (*applytrans)(void*,Vec,Vec);
339: $ PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*);
340: $ STCreate(comm,&st);
341: $ STSetType(st,STSHELL);
342: $ STShellSetApply(st,apply);
343: $ STShellSetApplyTranspose(st,applytrans);
344: $ STShellSetBackTransform(st,backtr); (optional)
346: M*/
350: PETSC_EXTERN PetscErrorCode STCreate_Shell(ST st)351: {
353: ST_SHELL *ctx;
356: PetscNewLog(st,&ctx);
357: st->data = (void*)ctx;
359: st->ops->apply = STApply_Shell;
360: st->ops->applytrans = STApplyTranspose_Shell;
361: st->ops->backtransform = STBackTransform_Shell;
362: st->ops->setfromoptions = STSetFromOptions_Shell;
363: st->ops->destroy = STDestroy_Shell;
364: PetscObjectComposeFunction((PetscObject)st,"STShellSetApply_C",STShellSetApply_Shell);
365: PetscObjectComposeFunction((PetscObject)st,"STShellSetApplyTranspose_C",STShellSetApplyTranspose_Shell);
366: PetscObjectComposeFunction((PetscObject)st,"STShellSetBackTransform_C",STShellSetBackTransform_Shell);
367: return(0);
368: }