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-2010, Universidad 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
 27:  #include slepceps.h

 30: typedef struct {
 31:   void           *ctx;                       /* user provided context */
 32:   PetscErrorCode (*apply)(ST,Vec,Vec);
 33:   PetscErrorCode (*applytrans)(ST,Vec,Vec);
 34:   PetscErrorCode (*backtr)(ST,PetscInt n,PetscScalar*,PetscScalar*);
 35:   char           *name;
 36: } ST_Shell;

 41: /*@C
 42:     STShellGetContext - Returns the user-provided context associated with a shell ST

 44:     Not Collective

 46:     Input Parameter:
 47: .   st - spectral transformation context

 49:     Output Parameter:
 50: .   ctx - the user provided context

 52:     Level: advanced

 54:     Notes:
 55:     This routine is intended for use within various shell routines
 56:     
 57: .seealso: STShellSetContext()
 58: @*/
 59: PetscErrorCode STShellGetContext(ST st,void **ctx)
 60: {
 62:   PetscTruth     flg;

 67:   PetscTypeCompare((PetscObject)st,STSHELL,&flg);
 68:   if (!flg) *ctx = 0;
 69:   else      *ctx = ((ST_Shell*)(st->data))->ctx;
 70:   return(0);
 71: }

 75: /*@
 76:     STShellSetContext - sets the context for a shell ST

 78:    Collective on ST

 80:     Input Parameters:
 81: +   st - the shell ST
 82: -   ctx - the context

 84:    Level: advanced

 86:    Fortran Notes: The context can only be an integer or a PetscObject;
 87:       unfortunately it cannot be a Fortran array or derived type.

 89: .seealso: STShellGetContext()
 90: @*/
 91: PetscErrorCode STShellSetContext(ST st,void *ctx)
 92: {
 93:   ST_Shell      *shell = (ST_Shell*)st->data;
 95:   PetscTruth     flg;

 99:   PetscTypeCompare((PetscObject)st,STSHELL,&flg);
100:   if (flg) {
101:     shell->ctx = ctx;
102:   }
103:   return(0);
104: }

108: PetscErrorCode STApply_Shell(ST st,Vec x,Vec y)
109: {
111:   ST_Shell       *shell = (ST_Shell*)st->data;

114:   if (!shell->apply) SETERRQ(PETSC_ERR_USER,"No apply() routine provided to Shell ST");
115:   PetscStackPush("STSHELL apply() user function");
116:   CHKMEMQ;
117:   (*shell->apply)(st,x,y);
118:   CHKMEMQ;
119:   PetscStackPop;
120:   return(0);
121: }

125: PetscErrorCode STApplyTranspose_Shell(ST st,Vec x,Vec y)
126: {
128:   ST_Shell       *shell = (ST_Shell*)st->data;

131:   if (!shell->applytrans) SETERRQ(PETSC_ERR_USER,"No applytranspose() routine provided to Shell ST");
132:   PetscStackPush("STSHELL applytranspose() user function");
133:   CHKMEMQ;
134:   (*shell->applytrans)(st,x,y);
135:   CHKMEMQ;
136:   PetscStackPop;
137:   return(0);
138: }

142: PetscErrorCode STBackTransform_Shell(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi)
143: {
145:   ST_Shell       *shell = (ST_Shell*)st->data;

148:   if (shell->backtr) {
149:     PetscStackPush("STSHELL backtransform() user function");
150:     CHKMEMQ;
151:     (*shell->backtr)(st,n,eigr,eigi);
152:     CHKMEMQ;
153:     PetscStackPop;
154:   }
155:   return(0);
156: }

160: PetscErrorCode STDestroy_Shell(ST st)
161: {
163:   ST_Shell       *shell = (ST_Shell*)st->data;

166:   PetscFree(shell->name);
167:   PetscFree(shell);
168:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApply_C","",PETSC_NULL);
169:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApplyTranspose_C","",PETSC_NULL);
170:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetBackTransform_C","",PETSC_NULL);
171:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetName_C","",PETSC_NULL);
172:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellGetName_C","",PETSC_NULL);
173:   return(0);
174: }

178: PetscErrorCode STView_Shell(ST st,PetscViewer viewer)
179: {
181:   ST_Shell       *ctx = (ST_Shell*)st->data;
182:   PetscTruth     isascii;

185:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
186:   if (isascii) {
187:     if (ctx->name) {PetscViewerASCIIPrintf(viewer,"  ST Shell: %s\n",ctx->name);}
188:     else           {PetscViewerASCIIPrintf(viewer,"  ST Shell: no name\n");}
189:   } else {
190:     SETERRQ1(1,"Viewer type %s not supported for STShell",((PetscObject)viewer)->type_name);
191:   }
192:   return(0);
193: }

198: PetscErrorCode STShellSetApply_Shell(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
199: {
200:   ST_Shell *shell = (ST_Shell*)st->data;

203:   shell->apply = apply;
204:   return(0);
205: }

211: PetscErrorCode STShellSetApplyTranspose_Shell(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
212: {
213:   ST_Shell *shell = (ST_Shell*)st->data;

216:   shell->applytrans = applytrans;
217:   return(0);
218: }

224: PetscErrorCode STShellSetBackTransform_Shell(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
225: {
226:   ST_Shell *shell = (ST_Shell *) st->data;

229:   shell->backtr = backtr;
230:   return(0);
231: }

237: PetscErrorCode STShellSetName_Shell(ST st,const char name[])
238: {
239:   ST_Shell *shell = (ST_Shell*)st->data;

243:   PetscFree(shell->name);
244:   PetscStrallocpy(name,&shell->name);
245:   return(0);
246: }

252: PetscErrorCode STShellGetName_Shell(ST st,char *name[])
253: {
254:   ST_Shell *shell = (ST_Shell*)st->data;

257:   *name  = shell->name;
258:   return(0);
259: }

264: /*@C
265:    STShellSetApply - Sets routine to use as the application of the 
266:    operator to a vector in the user-defined spectral transformation.

268:    Collective on ST

270:    Input Parameters:
271: +  st    - the spectral transformation context
272: -  apply - the application-provided transformation routine

274:    Calling sequence of apply:
275: .vb
276:    PetscErrorCode apply (ST st,Vec xin,Vec xout)
277: .ve

279: +  st   - the spectral transformation context
280: .  xin  - input vector
281: -  xout - output vector

283:    Level: developer

285: .seealso: STShellSetBackTransform(), STShellSetApplyTranspose()
286: @*/
287: PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(ST,Vec,Vec))
288: {
289:   PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(ST,Vec,Vec));

293:   PetscObjectQueryFunction((PetscObject)st,"STShellSetApply_C",(void (**)(void))&f);
294:   if (f) {
295:     (*f)(st,apply);
296:   }
297:   return(0);
298: }

302: /*@C
303:    STShellSetApplyTranspose - Sets routine to use as the application of the 
304:    transposed operator to a vector in the user-defined spectral transformation.

306:    Collective on ST

308:    Input Parameters:
309: +  st    - the spectral transformation context
310: -  applytrans - the application-provided transformation routine

312:    Calling sequence of apply:
313: .vb
314:    PetscErrorCode applytrans (ST st,Vec xin,Vec xout)
315: .ve

317: +  st   - the spectral transformation context
318: .  xin  - input vector
319: -  xout - output vector

321:    Level: developer

323: .seealso: STShellSetApply(), STShellSetBackTransform()
324: @*/
325: PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec))
326: {
327:   PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(ST,Vec,Vec));

331:   PetscObjectQueryFunction((PetscObject)st,"STShellSetApplyTranspose_C",(void (**)(void))&f);
332:   if (f) {
333:     (*f)(st,applytrans);
334:   }
335:   return(0);
336: }

340: /*@C
341:    STShellSetBackTransform - Sets the routine to be called after the 
342:    eigensolution process has finished in order to transform back the
343:    computed eigenvalues.

345:    Collective on ST

347:    Input Parameters:
348: +  st     - the spectral transformation context
349: -  backtr - the application-provided backtransform routine

351:    Calling sequence of backtr:
352: .vb
353:    PetscErrorCode backtr (ST st,PetscScalar *eigr,PetscScalar *eigi)
354: .ve

356: +  st   - the spectral transformation context
357: .  eigr - pointer ot the real part of the eigenvalue to transform back
358: -  eigi - pointer ot the imaginary part 

360:    Level: developer

362: .seealso: STShellSetApply(), STShellSetApplyTranspose()
363: @*/
364: PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*))
365: {
366:   PetscErrorCode ierr, (*f)(ST,PetscErrorCode (*)(ST,PetscInt,PetscScalar*,PetscScalar*));

370:   PetscObjectQueryFunction((PetscObject)st,"STShellSetBackTransform_C",(void (**)(void))&f);
371:   if (f) {
372:     (*f)(st,backtr);
373:   }
374:   return(0);
375: }

379: /*@C
380:    STShellSetName - Sets an optional name to associate with a shell
381:    spectral transformation.

383:    Not Collective

385:    Input Parameters:
386: +  st   - the spectral transformation context
387: -  name - character string describing the shell spectral transformation

389:    Level: developer

391: .seealso: STShellGetName()
392: @*/
393: PetscErrorCode STShellSetName(ST st,const char name[])
394: {
395:   PetscErrorCode ierr, (*f)(ST,const char []);

399:   PetscObjectQueryFunction((PetscObject)st,"STShellSetName_C",(void (**)(void))&f);
400:   if (f) {
401:     (*f)(st,name);
402:   }
403:   return(0);
404: }

408: /*@C
409:    STShellGetName - Gets an optional name that the user has set for a shell
410:    spectral transformation.

412:    Not Collective

414:    Input Parameter:
415: .  st - the spectral transformation context

417:    Output Parameter:
418: .  name - character string describing the shell spectral transformation 
419:           (you should not free this)

421:    Level: developer

423: .seealso: STShellSetName()
424: @*/
425: PetscErrorCode STShellGetName(ST st,char *name[])
426: {
427:   PetscErrorCode ierr, (*f)(ST,char *[]);

431:   PetscObjectQueryFunction((PetscObject)st,"STShellGetName_C",(void (**)(void))&f);
432:   if (f) {
433:     (*f)(st,name);
434:   } else {
435:     SETERRQ(PETSC_ERR_ARG_WRONG,"Not shell spectral transformation, cannot get name");
436:   }
437:   return(0);
438: }

442: PetscErrorCode STSetFromOptions_Shell(ST st)
443: {
445:   PC             pc;
446:   const PCType   pctype;
447:   const KSPType  ksptype;


451:   KSPGetPC(st->ksp,&pc);
452:   KSPGetType(st->ksp,&ksptype);
453:   PCGetType(pc,&pctype);
454:   if (!pctype && !ksptype) {
455:     if (st->shift_matrix == ST_MATMODE_SHELL) {
456:       /* in shell mode use GMRES with Jacobi as the default */
457:       KSPSetType(st->ksp,KSPGMRES);
458:       PCSetType(pc,PCJACOBI);
459:     } else {
460:       /* use direct solver as default */
461:       KSPSetType(st->ksp,KSPPREONLY);
462:       PCSetType(pc,PCREDUNDANT);
463:     }
464:   }

466:   return(0);
467: }

469: /*MC
470:    STSHELL - Creates a new spectral transformation class.
471:           This is intended to provide a simple class to use with EPS.
472:           You should not use this if you plan to make a complete class.

474:   Level: advanced

476:   Usage:
477: $             PetscErrorCode (*apply)(void*,Vec,Vec);
478: $             PetscErrorCode (*applytrans)(void*,Vec,Vec);
479: $             PetscErrorCode (*backtr)(void*,PetscScalar*,PetscScalar*);
480: $             STCreate(comm,&st);
481: $             STSetType(st,STSHELL);
482: $             STShellSetApply(st,apply);
483: $             STShellSetApplyTranspose(st,applytrans);
484: $             STShellSetBackTransform(st,backtr);    (optional)

486: M*/

491: PetscErrorCode STCreate_Shell(ST st)
492: {
494:   ST_Shell       *shell;

497:   st->ops->destroy = STDestroy_Shell;
498:   PetscNew(ST_Shell,&shell);
499:   PetscLogObjectMemory(st,sizeof(ST_Shell));

501:   st->data           = (void *) shell;
502:   ((PetscObject)st)->name           = 0;

504:   st->ops->apply          = STApply_Shell;
505:   st->ops->applytrans     = STApplyTranspose_Shell;
506:   st->ops->backtr         = STBackTransform_Shell;
507:   st->ops->view           = STView_Shell;
508:   st->ops->setfromoptions = STSetFromOptions_Shell;

510:   shell->apply       = 0;
511:   shell->applytrans  = 0;
512:   shell->backtr      = 0;
513:   shell->name        = 0;
514:   shell->ctx         = 0;

516:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApply_C","STShellSetApply_Shell",STShellSetApply_Shell);
517:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetApplyTranspose_C","STShellSetApplyTranspose_Shell",STShellSetApplyTranspose_Shell);
518:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetBackTransform_C","STShellSetBackTransform_Shell",STShellSetBackTransform_Shell);
519:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellSetName_C","STShellSetName_Shell",STShellSetName_Shell);
520:   PetscObjectComposeFunctionDynamic((PetscObject)st,"STShellGetName_C","STShellGetName_Shell",STShellGetName_Shell);

522:   return(0);
523: }