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-2007, Universidad Politecnica de Valencia, Spain
10: This file is part of SLEPc. See the README file for conditions of use
11: and additional information.
12: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
13: */
15: #include src/st/stimpl.h
16: #include slepceps.h
19: typedef struct {
20: void *ctx; /* user provided context */
21: PetscErrorCode (*apply)(void *,Vec,Vec);
22: PetscErrorCode (*applytrans)(void *,Vec,Vec);
23: PetscErrorCode (*backtr)(void *,PetscScalar*,PetscScalar*);
24: char *name;
25: } ST_Shell;
30: /*@C
31: STShellGetContext - Returns the user-provided context associated with a shell ST
33: Not Collective
35: Input Parameter:
36: . st - spectral transformation context
38: Output Parameter:
39: . ctx - the user provided context
41: Level: advanced
43: Notes:
44: This routine is intended for use within various shell routines
45:
46: .seealso: STShellSetContext()
47: @*/
48: PetscErrorCode STShellGetContext(ST st,void **ctx)
49: {
51: PetscTruth flg;
56: PetscTypeCompare((PetscObject)st,STSHELL,&flg);
57: if (!flg) *ctx = 0;
58: else *ctx = ((ST_Shell*)(st->data))->ctx;
59: return(0);
60: }
64: /*@C
65: STShellSetContext - sets the context for a shell ST
67: Collective on ST
69: Input Parameters:
70: + st - the shell ST
71: - ctx - the context
73: Level: advanced
75: Fortran Notes: The context can only be an integer or a PetscObject;
76: unfortunately it cannot be a Fortran array or derived type.
78: .seealso: STShellGetContext()
79: @*/
80: PetscErrorCode STShellSetContext(ST st,void *ctx)
81: {
82: ST_Shell *shell = (ST_Shell*)st->data;
84: PetscTruth flg;
88: PetscTypeCompare((PetscObject)st,STSHELL,&flg);
89: if (flg) {
90: shell->ctx = ctx;
91: }
92: return(0);
93: }
97: PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
98: {
100: ST_Shell *shell = (ST_Shell*)st->data;
103: if (!shell->apply) SETERRQ(PETSC_ERR_USER,"No apply() routine provided to Shell ST");
104: PetscStackPush("PCSHELL user function");
105: CHKMEMQ;
106: (*shell->apply)(shell->ctx,x,y);
107: CHKMEMQ;
108: PetscStackPop;
109: return(0);
110: }
114: PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
115: {
117: ST_Shell *shell = (ST_Shell*)st->data;
120: if (!shell->applytrans) SETERRQ(PETSC_ERR_USER,"No applytranspose() routine provided to Shell ST");
121: (*shell->applytrans)(shell->ctx,x,y);
122: return(0);
123: }
127: PetscErrorCode STBackTransform_Shell(ST st,PetscScalar *eigr,PetscScalar *eigi)
128: {
130: ST_Shell *shell = (ST_Shell*)st->data;
133: if (shell->backtr) {
134: (*shell->backtr)(shell->ctx,eigr,eigi);
135: }
136: return(0);
137: }
141: PetscErrorCode STDestroy_Shell(ST st)
142: {
144: ST_Shell *shell = (ST_Shell*)st->data;
147: PetscFree(shell->name);
148: PetscFree(shell);
149: return(0);
150: }
154: PetscErrorCode STView_Shell(ST st,PetscViewer viewer)
155: {
157: ST_Shell *ctx = (ST_Shell*)st->data;
158: PetscTruth isascii;
161: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
162: if (isascii) {
163: if (ctx->name) {PetscViewerASCIIPrintf(viewer," ST Shell: %s\n",ctx->name);}
164: else {PetscViewerASCIIPrintf(viewer," ST Shell: no name\n");}
165: } else {
166: SETERRQ1(1,"Viewer type %s not supported for STShell",((PetscObject)viewer)->type_name);
167: }
168: return(0);
169: }
174: PetscErrorCode STShellSetApply_Shell(ST st,PetscErrorCode (*apply)(void*,Vec,Vec))
175: {
176: ST_Shell *shell = (ST_Shell*)st->data;
179: shell->apply = apply;
180: return(0);
181: }
187: PetscErrorCode STShellSetApplyTranspose_Shell(ST st,PetscErrorCode (*applytrans)(void*,Vec,Vec))
188: {
189: ST_Shell *shell = (ST_Shell*)st->data;
192: shell->applytrans = applytrans;
193: return(0);
194: }
200: PetscErrorCode STShellSetBackTransform_Shell(ST st,PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*))
201: {
202: ST_Shell *shell = (ST_Shell *) st->data;
205: shell->backtr = backtr;
206: return(0);
207: }
213: PetscErrorCode STShellSetName_Shell(ST st,const char name[])
214: {
215: ST_Shell *shell = (ST_Shell*)st->data;
219: PetscStrfree(shell->name);
220: PetscStrallocpy(name,&shell->name);
221: return(0);
222: }
228: PetscErrorCode STShellGetName_Shell(ST st,char *name[])
229: {
230: ST_Shell *shell = (ST_Shell*)st->data;
233: *name = shell->name;
234: return(0);
235: }
240: /*@C
241: STShellSetApply - Sets routine to use as the application of the
242: operator to a vector in the user-defined spectral transformation.
244: Collective on ST
246: Input Parameters:
247: + st - the spectral transformation context
248: - apply - the application-provided transformation routine
250: Calling sequence of apply:
251: .vb
252: PetscErrorCode apply (void *ptr,Vec xin,Vec xout)
253: .ve
255: + ptr - the application context
256: . xin - input vector
257: - xout - output vector
259: Level: developer
261: .seealso: STShellSetBackTransform(), STShellSetApplyTranspose()
262: @*/
263: PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(void*,Vec,Vec))
264: {
265: PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(void*,Vec,Vec));
269: PetscObjectQueryFunction((PetscObject)st,"STShellSetApply_C",(void (**)(void))&f);
270: if (f) {
271: (*f)(st,apply);
272: }
273: return(0);
274: }
278: /*@C
279: STShellSetApplyTranspose - Sets routine to use as the application of the
280: transposed operator to a vector in the user-defined spectral transformation.
282: Collective on ST
284: Input Parameters:
285: + st - the spectral transformation context
286: - applytrans - the application-provided transformation routine
288: Calling sequence of apply:
289: .vb
290: PetscErrorCode applytrans (void *ptr,Vec xin,Vec xout)
291: .ve
293: + ptr - the application context
294: . xin - input vector
295: - xout - output vector
297: Level: developer
299: .seealso: STShellSetApply(), STShellSetBackTransform()
300: @*/
301: PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(void*,Vec,Vec))
302: {
303: PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(void*,Vec,Vec));
307: PetscObjectQueryFunction((PetscObject)st,"STShellSetApplyTranspose_C",(void (**)(void))&f);
308: if (f) {
309: (*f)(st,applytrans);
310: }
311: return(0);
312: }
316: /*@C
317: STShellSetBackTransform - Sets the routine to be called after the
318: eigensolution process has finished in order to transform back the
319: computed eigenvalues.
321: Collective on ST
323: Input Parameters:
324: + st - the spectral transformation context
325: - backtr - the application-provided backtransform routine
327: Calling sequence of backtr:
328: .vb
329: PetscErrorCode backtr (void *ptr,PetscScalar *eigr,PetscScalar *eigi)
330: .ve
332: + ptr - the application context
333: . eigr - pointer ot the real part of the eigenvalue to transform back
334: - eigi - pointer ot the imaginary part
336: Level: developer
338: .seealso: STShellSetApply(), STShellSetApplyTranspose()
339: @*/
340: PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*))
341: {
342: PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(void*,PetscScalar*,PetscScalar*));
346: PetscObjectQueryFunction((PetscObject)st,"STShellSetBackTransform_C",(void (**)(void))&f);
347: if (f) {
348: (*f)(st,(int (*)(void*,PetscScalar*,PetscScalar*))backtr);
349: }
350: return(0);
351: }
355: /*@C
356: STShellSetName - Sets an optional name to associate with a shell
357: spectral transformation.
359: Not Collective
361: Input Parameters:
362: + st - the spectral transformation context
363: - name - character string describing the shell spectral transformation
365: Level: developer
367: .seealso: STShellGetName()
368: @*/
369: PetscErrorCode STShellSetName(ST st,const char name[])
370: {
371: PetscErrorCode ierr, (*f)(ST,const char []);
375: PetscObjectQueryFunction((PetscObject)st,"STShellSetName_C",(void (**)(void))&f);
376: if (f) {
377: (*f)(st,name);
378: }
379: return(0);
380: }
384: /*@C
385: STShellGetName - Gets an optional name that the user has set for a shell
386: spectral transformation.
388: Not Collective
390: Input Parameter:
391: . st - the spectral transformation context
393: Output Parameter:
394: . name - character string describing the shell spectral transformation
395: (you should not free this)
397: Level: developer
399: .seealso: STShellSetName()
400: @*/
401: PetscErrorCode STShellGetName(ST st,char *name[])
402: {
403: PetscErrorCode ierr, (*f)(ST,char *[]);
407: PetscObjectQueryFunction((PetscObject)st,"STShellGetName_C",(void (**)(void))&f);
408: if (f) {
409: (*f)(st,name);
410: } else {
411: SETERRQ(PETSC_ERR_ARG_WRONG,"Not shell spectral transformation, cannot get name");
412: }
413: return(0);
414: }
416: /*MC
417: STSHELL - Creates a new spectral transformation class.
418: This is intended to provide a simple class to use with EPS.
419: You should not use this if you plan to make a complete class.
421: Level: advanced
423: Usage:
424: $ PetscErrorCode (*apply)(void*,Vec,Vec);
425: $ PetscErrorCode (*applytrans)(void*,Vec,Vec);
426: $ PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*);
427: $ STCreate(comm,&st);
428: $ STSetType(st,STSHELL);
429: $ STShellSetApply(st,apply);
430: $ STShellSetApplyTranspose(st,applytrans);
431: $ STShellSetBackTransform(st,backtr); (optional)
433: M*/
438: PetscErrorCode STCreate_Shell(ST st)
439: {
441: ST_Shell *shell;
444: st->ops->destroy = STDestroy_Shell;
445: PetscNew(ST_Shell,&shell);
446: PetscLogObjectMemory(st,sizeof(ST_Shell));
448: st->data = (void *) shell;
449: st->name = 0;
451: st->ops->apply = STApply_Shell;
452: st->ops->applytrans= STApplyTranspose_Shell;
453: st->ops->backtr = STBackTransform_Shell;
454: st->ops->view = STView_Shell;
456: shell->apply = 0;
457: shell->applytrans = 0;
458: shell->backtr = 0;
459: shell->name = 0;
460: shell->ctx = 0;
462: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApply_C","STShellSetApply_Shell",
463: STShellSetApply_Shell);
464: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApplyTranspose_C","STShellSetApplyTranspose_Shell",
465: STShellSetApplyTranspose_Shell);
466: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetBackTransform_C","STShellSetBackTransform_Shell",
467: STShellSetBackTransform_Shell);
468: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetName_C","STShellSetName_Shell",
469: STShellSetName_Shell);
470: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellGetName_C","STShellGetName_Shell",
471: STShellGetName_Shell);
473: return(0);
474: }