Actual source code: epsimpl.h

  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2009, Universidad Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:       
  8:    SLEPc is free software: you can redistribute it and/or modify it under  the
  9:    terms of version 3 of the GNU Lesser General Public License as published by
 10:    the Free Software Foundation.

 12:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY 
 13:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS 
 14:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for 
 15:    more details.

 17:    You  should have received a copy of the GNU Lesser General  Public  License
 18:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 19:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 20: */

 22: #ifndef _EPSIMPL
 23: #define _EPSIMPL

 25:  #include slepceps.h


 30: typedef struct _EPSOps *EPSOps;

 32: struct _EPSOps {
 33:   PetscErrorCode  (*solve)(EPS);            /* one-sided solver */
 34:   PetscErrorCode  (*solvets)(EPS);          /* two-sided solver */
 35:   PetscErrorCode  (*setup)(EPS);
 36:   PetscErrorCode  (*setfromoptions)(EPS);
 37:   PetscErrorCode  (*publishoptions)(EPS);
 38:   PetscErrorCode  (*destroy)(EPS);
 39:   PetscErrorCode  (*view)(EPS,PetscViewer);
 40:   PetscErrorCode  (*backtransform)(EPS);
 41:   PetscErrorCode  (*computevectors)(EPS);
 42: };

 44: /*
 45:      Maximum number of monitors you can run with a single EPS
 46: */
 47: #define MAXEPSMONITORS 5 

 49: /*
 50:    Defines the EPS data structure.
 51: */
 52: struct _p_EPS {
 53:   PETSCHEADER(struct _EPSOps);
 54:   /*------------------------- User parameters --------------------------*/
 55:   PetscInt       max_it,           /* maximum number of iterations */
 56:                  nev,              /* number of eigenvalues to compute */
 57:                  ncv,              /* number of basis vectors */
 58:                  mpd,              /* maximum dimension of projected problem */
 59:                  allocated_ncv,    /* number of basis vectors allocated */
 60:                  nds;              /* number of basis vectors of deflation space */
 61:   PetscScalar    target;           /* target value */
 62:   PetscTruth     target_set;       /* flag indicating if target was specified */
 63:   PetscReal      tol;              /* tolerance */
 64:   EPSWhich       which;            /* which part of the spectrum to be sought */
 65:   PetscTruth     evecsavailable;   /* computed eigenvectors */
 66:   EPSProblemType problem_type;     /* which kind of problem to be solved */
 67:   EPSExtraction  extraction;       /* which kind of extraction to be applied */
 68:   EPSClass       solverclass;      /* whether the selected solver is one- or two-sided */

 70:   /*------------------------- Working data --------------------------*/
 71:   Vec         vec_initial,      /* initial vector */
 72:               vec_initial_left, /* left initial vector for two-sided solvers */
 73:               *V,               /* set of basis vectors and computed eigenvectors */
 74:               *AV,              /* auxiliar set of basis vectors */
 75:               *W,               /* set of left basis vectors and computed left eigenvectors */
 76:               *DS,              /* deflation space */
 77:               *DSV;             /* deflation space and basis vectors*/
 78:   PetscScalar *eigr, *eigi,     /* real and imaginary parts of eigenvalues */
 79:               *T, *Tl;          /* projected matrices */
 80:   PetscReal   *errest,          /* error estimates */
 81:               *errest_left;     /* left error estimates */
 82:   ST          OP;               /* spectral transformation object */
 83:   IP          ip;               /* innerproduct object */
 84:   void        *data;            /* placeholder for misc stuff associated 
 85:                                    with a particular solver */
 86:   PetscInt    nconv,            /* number of converged eigenvalues */
 87:               its,              /* number of iterations so far computed */
 88:               *perm;            /* permutation for eigenvalue ordering */

 90:   /* ---------------- Default work-area and status vars -------------------- */
 91:   PetscInt   nwork;
 92:   Vec        *work;

 94:   PetscInt   setupcalled;
 95:   PetscTruth isgeneralized,
 96:              ispositive,
 97:              ishermitian;
 98:   EPSConvergedReason reason;

100:   PetscErrorCode (*monitor[MAXEPSMONITORS])(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
101:   PetscErrorCode (*monitordestroy[MAXEPSMONITORS])(void*);
102:   void       *monitorcontext[MAXEPSMONITORS];
103:   PetscInt    numbermonitors;

105:   PetscTruth ds_ortho;    /* if vectors in DS have to be orthonormalized */
106: };

108: #define EPSMonitor(eps,it,nconv,eigr,eigi,errest,nest) \
109:         { PetscErrorCode _ierr; PetscInt _i,_im = eps->numbermonitors; \
110:           for ( _i=0; _i<_im; _i++ ) {\
111:             _ierr=(*eps->monitor[_i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[_i]);\
112:             CHKERRQ(_ierr); \
113:           } \
114:         }

116: EXTERN PetscErrorCode EPSRegisterAll(char *);

118: EXTERN PetscErrorCode EPSDestroy_Default(EPS);
119: EXTERN PetscErrorCode EPSDefaultGetWork(EPS,PetscInt);
120: EXTERN PetscErrorCode EPSDefaultFreeWork(EPS);
121: EXTERN PetscErrorCode EPSAllocateSolution(EPS);
122: EXTERN PetscErrorCode EPSFreeSolution(EPS);
123: EXTERN PetscErrorCode EPSBackTransform_Default(EPS);
124: EXTERN PetscErrorCode EPSComputeVectors_Default(EPS);
125: EXTERN PetscErrorCode EPSComputeVectors_Hermitian(EPS);
126: EXTERN PetscErrorCode EPSComputeVectors_Schur(EPS);

128: /* Private functions of the solver implementations */

130: EXTERN PetscErrorCode EPSBasicArnoldi(EPS,PetscTruth,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscTruth*);
131: EXTERN PetscErrorCode EPSDelayedArnoldi(EPS,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscTruth*);
132: EXTERN PetscErrorCode EPSDelayedArnoldi1(EPS,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscTruth*);
133: EXTERN PetscErrorCode ArnoldiResiduals(PetscScalar*,PetscInt,PetscScalar*,PetscReal,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscScalar*);
134: EXTERN PetscErrorCode EPSFullLanczos(EPS,PetscReal*,PetscReal*,Vec*,PetscInt,PetscInt*,Vec,PetscTruth*);
135: EXTERN PetscErrorCode EPSTranslateHarmonic(PetscInt,PetscScalar*,PetscInt,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*);

137: #endif