Actual source code: shell.c
2: /*
3: This provides a simple shell interface for programmers to
4: create their own spectral transformations without writing much
5: interface code.
6: */
8: #include src/st/stimpl.h
9: #include slepceps.h
12: typedef struct {
13: void *ctx; /* user provided context */
14: PetscErrorCode (*apply)(void *,Vec,Vec);
15: PetscErrorCode (*applytrans)(void *,Vec,Vec);
16: PetscErrorCode (*backtr)(void *,PetscScalar*,PetscScalar*);
17: char *name;
18: } ST_Shell;
23: /*@C
24: STShellGetContext - Returns the user-provided context associated with a shell ST
26: Not Collective
28: Input Parameter:
29: . st - spectral transformation context
31: Output Parameter:
32: . ctx - the user provided context
34: Level: advanced
36: Notes:
37: This routine is intended for use within various shell routines
38:
39: .seealso: STShellSetContext()
40: @*/
41: PetscErrorCode STShellGetContext(ST st,void **ctx)
42: {
44: PetscTruth flg;
49: PetscTypeCompare((PetscObject)st,STSHELL,&flg);
50: if (!flg) *ctx = 0;
51: else *ctx = ((ST_Shell*)(st->data))->ctx;
52: return(0);
53: }
57: /*@C
58: STShellSetContext - sets the context for a shell ST
60: Collective on ST
62: Input Parameters:
63: + st - the shell ST
64: - ctx - the context
66: Level: advanced
68: Fortran Notes: The context can only be an integer or a PetscObject;
69: unfortunately it cannot be a Fortran array or derived type.
71: .seealso: STShellGetContext()
72: @*/
73: PetscErrorCode STShellSetContext(ST st,void *ctx)
74: {
75: ST_Shell *shell = (ST_Shell*)st->data;
77: PetscTruth flg;
81: PetscTypeCompare((PetscObject)st,STSHELL,&flg);
82: if (flg) {
83: shell->ctx = ctx;
84: }
85: return(0);
86: }
90: PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
91: {
93: ST_Shell *shell = (ST_Shell*)st->data;
96: if (!shell->apply) SETERRQ(PETSC_ERR_USER,"No apply() routine provided to Shell ST");
97: PetscStackPush("PCSHELL user function");
98: CHKMEMQ;
99: (*shell->apply)(shell->ctx,x,y);
100: CHKMEMQ;
101: PetscStackPop;
102: return(0);
103: }
107: PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
108: {
110: ST_Shell *shell = (ST_Shell*)st->data;
113: if (!shell->applytrans) SETERRQ(PETSC_ERR_USER,"No applytranspose() routine provided to Shell ST");
114: (*shell->applytrans)(shell->ctx,x,y);
115: return(0);
116: }
120: PetscErrorCode STBackTransform_Shell(ST st,PetscScalar *eigr,PetscScalar *eigi)
121: {
123: ST_Shell *shell = (ST_Shell*)st->data;
126: if (shell->backtr) {
127: (*shell->backtr)(shell->ctx,eigr,eigi);
128: }
129: return(0);
130: }
134: PetscErrorCode STDestroy_Shell(ST st)
135: {
137: ST_Shell *shell = (ST_Shell*)st->data;
140: PetscFree(shell->name);
141: PetscFree(shell);
142: return(0);
143: }
147: PetscErrorCode STView_Shell(ST st,PetscViewer viewer)
148: {
150: ST_Shell *ctx = (ST_Shell*)st->data;
151: PetscTruth isascii;
154: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
155: if (isascii) {
156: if (ctx->name) {PetscViewerASCIIPrintf(viewer," ST Shell: %s\n",ctx->name);}
157: else {PetscViewerASCIIPrintf(viewer," ST Shell: no name\n");}
158: } else {
159: SETERRQ1(1,"Viewer type %s not supported for STShell",((PetscObject)viewer)->type_name);
160: }
161: return(0);
162: }
167: PetscErrorCode STShellSetApply_Shell(ST st,PetscErrorCode (*apply)(void*,Vec,Vec))
168: {
169: ST_Shell *shell = (ST_Shell*)st->data;
172: shell->apply = apply;
173: return(0);
174: }
180: PetscErrorCode STShellSetApplyTranspose_Shell(ST st,PetscErrorCode (*applytrans)(void*,Vec,Vec))
181: {
182: ST_Shell *shell = (ST_Shell*)st->data;
185: shell->applytrans = applytrans;
186: return(0);
187: }
193: PetscErrorCode STShellSetBackTransform_Shell(ST st,PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*))
194: {
195: ST_Shell *shell = (ST_Shell *) st->data;
198: shell->backtr = backtr;
199: return(0);
200: }
206: PetscErrorCode STShellSetName_Shell(ST st,const char name[])
207: {
208: ST_Shell *shell = (ST_Shell*)st->data;
212: PetscStrfree(shell->name);
213: PetscStrallocpy(name,&shell->name);
214: return(0);
215: }
221: PetscErrorCode STShellGetName_Shell(ST st,char *name[])
222: {
223: ST_Shell *shell = (ST_Shell*)st->data;
226: *name = shell->name;
227: return(0);
228: }
233: /*@C
234: STShellSetApply - Sets routine to use as the application of the
235: operator to a vector in the user-defined spectral transformation.
237: Collective on ST
239: Input Parameters:
240: + st - the spectral transformation context
241: - apply - the application-provided transformation routine
243: Calling sequence of apply:
244: .vb
245: PetscErrorCode apply (void *ptr,Vec xin,Vec xout)
246: .ve
248: + ptr - the application context
249: . xin - input vector
250: - xout - output vector
252: Level: developer
254: .seealso: STShellSetBackTransform(), STShellSetApplyTranspose()
255: @*/
256: PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(void*,Vec,Vec))
257: {
258: PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(void*,Vec,Vec));
262: PetscObjectQueryFunction((PetscObject)st,"STShellSetApply_C",(void (**)(void))&f);
263: if (f) {
264: (*f)(st,apply);
265: }
266: return(0);
267: }
271: /*@C
272: STShellSetApplyTranspose - Sets routine to use as the application of the
273: transposed operator to a vector in the user-defined spectral transformation.
275: Collective on ST
277: Input Parameters:
278: + st - the spectral transformation context
279: - applytrans - the application-provided transformation routine
281: Calling sequence of apply:
282: .vb
283: PetscErrorCode applytrans (void *ptr,Vec xin,Vec xout)
284: .ve
286: + ptr - the application context
287: . xin - input vector
288: - xout - output vector
290: Level: developer
292: .seealso: STShellSetApply(), STShellSetBackTransform()
293: @*/
294: PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(void*,Vec,Vec))
295: {
296: PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(void*,Vec,Vec));
300: PetscObjectQueryFunction((PetscObject)st,"STShellSetApplyTranspose_C",(void (**)(void))&f);
301: if (f) {
302: (*f)(st,applytrans);
303: }
304: return(0);
305: }
309: /*@C
310: STShellSetBackTransform - Sets the routine to be called after the
311: eigensolution process has finished in order to transform back the
312: computed eigenvalues.
314: Collective on ST
316: Input Parameters:
317: + st - the spectral transformation context
318: - backtr - the application-provided backtransform routine
320: Calling sequence of backtr:
321: .vb
322: PetscErrorCode backtr (void *ptr,PetscScalar *eigr,PetscScalar *eigi)
323: .ve
325: + ptr - the application context
326: . eigr - pointer ot the real part of the eigenvalue to transform back
327: - eigi - pointer ot the imaginary part
329: Level: developer
331: .seealso: STShellSetApply(), STShellSetApplyTranspose()
332: @*/
333: PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*))
334: {
335: PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(void*,PetscScalar*,PetscScalar*));
339: PetscObjectQueryFunction((PetscObject)st,"STShellSetBackTransform_C",(void (**)(void))&f);
340: if (f) {
341: (*f)(st,(int (*)(void*,PetscScalar*,PetscScalar*))backtr);
342: }
343: return(0);
344: }
348: /*@C
349: STShellSetName - Sets an optional name to associate with a shell
350: spectral transformation.
352: Not Collective
354: Input Parameters:
355: + st - the spectral transformation context
356: - name - character string describing the shell spectral transformation
358: Level: developer
360: .seealso: STShellGetName()
361: @*/
362: PetscErrorCode STShellSetName(ST st,const char name[])
363: {
364: PetscErrorCode ierr, (*f)(ST,const char []);
368: PetscObjectQueryFunction((PetscObject)st,"STShellSetName_C",(void (**)(void))&f);
369: if (f) {
370: (*f)(st,name);
371: }
372: return(0);
373: }
377: /*@C
378: STShellGetName - Gets an optional name that the user has set for a shell
379: spectral transformation.
381: Not Collective
383: Input Parameter:
384: . st - the spectral transformation context
386: Output Parameter:
387: . name - character string describing the shell spectral transformation
388: (you should not free this)
390: Level: developer
392: .seealso: STShellSetName()
393: @*/
394: PetscErrorCode STShellGetName(ST st,char *name[])
395: {
396: PetscErrorCode ierr, (*f)(ST,char *[]);
400: PetscObjectQueryFunction((PetscObject)st,"STShellGetName_C",(void (**)(void))&f);
401: if (f) {
402: (*f)(st,name);
403: } else {
404: SETERRQ(PETSC_ERR_ARG_WRONG,"Not shell spectral transformation, cannot get name");
405: }
406: return(0);
407: }
409: /*MC
410: STSHELL - Creates a new spectral transformation class.
411: This is intended to provide a simple class to use with EPS.
412: You should not use this if you plan to make a complete class.
414: Level: advanced
416: Usage:
417: $ PetscErrorCode (*apply)(void*,Vec,Vec);
418: $ PetscErrorCode (*applytrans)(void*,Vec,Vec);
419: $ PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*);
420: $ STCreate(comm,&st);
421: $ STSetType(st,STSHELL);
422: $ STShellSetApply(st,apply);
423: $ STShellSetApplyTranspose(st,applytrans);
424: $ STShellSetBackTransform(st,backtr); (optional)
426: M*/
431: PetscErrorCode STCreate_Shell(ST st)
432: {
434: ST_Shell *shell;
437: st->ops->destroy = STDestroy_Shell;
438: PetscNew(ST_Shell,&shell);
439: PetscLogObjectMemory(st,sizeof(ST_Shell));
441: st->data = (void *) shell;
442: st->name = 0;
444: st->ops->apply = STApply_Shell;
445: st->ops->applytrans= STApplyTranspose_Shell;
446: st->ops->backtr = STBackTransform_Shell;
447: st->ops->view = STView_Shell;
449: shell->apply = 0;
450: shell->applytrans = 0;
451: shell->backtr = 0;
452: shell->name = 0;
453: shell->ctx = 0;
455: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApply_C","STShellSetApply_Shell",
456: STShellSetApply_Shell);
457: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApplyTranspose_C","STShellSetApplyTranspose_Shell",
458: STShellSetApplyTranspose_Shell);
459: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetBackTransform_C","STShellSetBackTransform_Shell",
460: STShellSetBackTransform_Shell);
461: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetName_C","STShellSetName_Shell",
462: STShellSetName_Shell);
463: PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellGetName_C","STShellGetName_Shell",
464: STShellGetName_Shell);
466: return(0);
467: }