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  (*setfromoptions)(EPS);
 13:   int  (*publishoptions)(EPS);
 14:   int  (*destroy)(EPS);
 15:   int  (*view)(EPS,PetscViewer);
 16:   int  (*backtransform)(EPS);
 17:   int  (*computevectors)(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:              allocated_ncv,     /* number of basis vectors allocated */
 35:              nds;               /* number of basis vectors of deflation space */
 36:   PetscReal  tol;               /* tolerance */
 37:   EPSWhich   which;             /* which part of the spectrum to be sought */
 38:   PetscTruth evecsavailable;   /* computed eigenvectors */
 39:   EPSProblemType problem_type;  /* which kind of problem to be solved */

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

 57:   /* ---------------- Default work-area and status vars -------------------- */
 58:   int        nwork;
 59:   Vec        *work;

 61:   int        setupcalled;
 62:   PetscTruth isgeneralized,
 63:              ishermitian,
 64:              vec_initial_set;
 65:   EPSConvergedReason reason;

 67:   int        (*monitor[MAXEPSMONITORS])(EPS,int,int,PetscScalar*,PetscScalar*,PetscReal*,int,void*);
 68:   int        (*monitordestroy[MAXEPSMONITORS])(void*);
 69:   void       *monitorcontext[MAXEPSMONITORS];
 70:   int        numbermonitors;

 72:   /* --------------- Orthogonalization --------------------- */
 73:   EPSOrthogonalizationType           orthog_type; /* which orthogonalization to use */
 74:   EPSOrthogonalizationRefinementType orthog_ref;   /* refinement method */
 75:   PetscReal               orthog_eta;
 76:   PetscTruth              ds_ortho;    /* if vectors in DS have to be orthonormalized */
 77: };

 79: #define EPSMonitor(eps,it,nconv,eigr,eigi,errest,nest) \
 80:         { int _ierr,_i,_im = eps->numbermonitors; \
 81:           for ( _i=0; _i<_im; _i++ ) {\
 82:             _ierr=(*eps->monitor[_i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[_i]);\
 83:             CHKERRQ(_ierr); \
 84:           } \
 85:         }

 87: EXTERN PetscErrorCode EPSDestroy_Default(EPS);
 88: EXTERN PetscErrorCode EPSDefaultGetWork(EPS,int);
 89: EXTERN PetscErrorCode EPSDefaultFreeWork(EPS);
 90: EXTERN PetscErrorCode EPSAllocateSolution(EPS);
 91: EXTERN PetscErrorCode EPSFreeSolution(EPS);
 92: EXTERN PetscErrorCode EPSAllocateSolutionContiguous(EPS);
 93: EXTERN PetscErrorCode EPSFreeSolutionContiguous(EPS);
 94: EXTERN PetscErrorCode EPSBackTransform_Default(EPS);
 95: EXTERN PetscErrorCode EPSComputeVectors_Default(EPS);
 96: EXTERN PetscErrorCode EPSComputeVectors_Schur(EPS);

 98: #endif