Actual source code: epsimpl.h

  2: #ifndef _EPSIMPL
  3: #define _EPSIMPL

 5:  #include slepceps.h

  7: typedef struct _EPSOps *EPSOps;

  9: struct _EPSOps {
 10:   int  (*solve)(EPS);                   /* actual solver */
 11:   int  (*setup)(EPS);
 12:   int  (*setdefaults)(EPS);
 13:   int  (*setfromoptions)(EPS);
 14:   int  (*publishoptions)(EPS);
 15:   int  (*destroy)(EPS);
 16:   int  (*view)(EPS,PetscViewer);
 17:   int  (*backtransform)(EPS);
 18: };

 20: /*
 21:      Maximum number of monitors you can run with a single EPS
 22: */
 23: #define MAXEPSMONITORS 5 

 25: /*
 26:    Defines the EPS data structure.
 27: */
 28: struct _p_EPS {
 29:   PETSCHEADER(struct _EPSOps)
 30:   /*------------------------- User parameters --------------------------*/
 31:   int        max_it,            /* maximum number of iterations */
 32:              nev,               /* number of eigenvalues to compute */
 33:              ncv;               /* number of basis vectors */
 34:   PetscReal  tol;               /* tolerance */
 35:   EPSWhich   which;             /* which part of the spectrum to be sought */
 36:   PetscTruth dropvectors;       /* do not compute eigenvectors */
 37:   EPSProblemType problem_type;  /* which kind of problem to be solved */

 39:   /*------------------------- Working data --------------------------*/
 40:   Vec         vec_initial;      /* initial vector for iterative methods */
 41:   Vec         *V;               /* set of basis vectors */
 42:   PetscScalar *eigr, *eigi;     /* real and imaginary parts of eigenvalues */
 43:   PetscReal  *errest;           /* error estimates */
 44:   ST         OP;                /* spectral transformation object */
 45:   void       *data;             /* holder for misc stuff associated 
 46:                                    with a particular solver */
 47:   int        nconv,             /* number of converged eigenvalues */
 48:              its;               /* number of iterations so far computed */
 49:   int        *perm;             /* permutation for eigenvalue ordering */

 51:   /* ---------------- Default work-area and status vars -------------------- */
 52:   int        nwork;
 53:   Vec        *work;

 55:   int        setupcalled;
 56:   PetscTruth isgeneralized,
 57:              ishermitian;
 58:   EPSConvergedReason reason;

 60:   int        (*monitor[MAXEPSMONITORS])(EPS,int,int,PetscReal*,int,void*);
 61:   int        (*monitordestroy[MAXEPSMONITORS])(void*);
 62:   void       *monitorcontext[MAXEPSMONITORS];
 63:   int        numbermonitors;
 64:   int        (*vmonitor[MAXEPSMONITORS])(EPS,int,int,PetscScalar*,PetscScalar*,int,void*);
 65:   int        (*vmonitordestroy[MAXEPSMONITORS])(void*);
 66:   void       *vmonitorcontext[MAXEPSMONITORS];
 67:   int        numbervmonitors;

 69:   /* --------------- Orthogonalization --------------------- */
 70:   int        (*orthog)(EPS,int,PetscScalar*,PetscReal*);
 71:   PetscReal  orth_eta;
 72:   EPSOrthogonalizationType orth_type;   /* which orthogonalization to use */
 73: 
 74: };

 76: #define EPSMonitorEstimates(eps,it,nconv,errest,nest) \
 77:         { int _ierr,_i,_im = eps->numbermonitors; \
 78:           for ( _i=0; _i<_im; _i++ ) {\
 79:             _ierr=(*eps->monitor[_i])(eps,it,nconv,errest,nest,eps->monitorcontext[_i]);\
 80:             CHKERRQ(_ierr); \
 81:           } \
 82:         }

 84: #define EPSMonitorValues(eps,it,nconv,eigr,eigi,neig) \
 85:         { int _ierr,_i,_im = eps->numbervmonitors; \
 86:           for ( _i=0; _i<_im; _i++ ) {\
 87:             _ierr=(*eps->vmonitor[_i])(eps,it,nconv,eigr,eigi,neig,eps->monitorcontext[_i]);\
 88:             CHKERRQ(_ierr); \
 89:           } \
 90:         }

 92: extern int EPSDefaultDestroy(EPS);
 93: extern int EPSDefaultGetWork(EPS,int);
 94: extern int EPSDefaultFreeWork(EPS);
 95: extern int EPSModifiedGramSchmidtOrthogonalization(EPS,int,PetscScalar*,PetscReal*);
 96: extern int EPSClassicalGramSchmidtOrthogonalization(EPS,int,PetscScalar*,PetscReal*);
 97: extern int EPSBackTransform_Default(EPS);

 99: #endif