Actual source code: epsimpl.h

  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2011, Universitat 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>

 27: extern PetscFList EPSList;
 28: extern PetscLogEvent EPS_SetUp, EPS_Solve, EPS_Dense;

 30: typedef struct _EPSOps *EPSOps;

 32: struct _EPSOps {
 33:   PetscErrorCode  (*solve)(EPS);
 34:   PetscErrorCode  (*setup)(EPS);
 35:   PetscErrorCode  (*setfromoptions)(EPS);
 36:   PetscErrorCode  (*publishoptions)(EPS);
 37:   PetscErrorCode  (*destroy)(EPS);
 38:   PetscErrorCode  (*reset)(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:                  nini, ninil,      /* number of initial vectors (negative means not copied yet) */
 60:                  nds;              /* number of basis vectors of deflation space */
 61:   PetscScalar    target;           /* target value */
 62:   PetscReal      tol;              /* tolerance */
 63:   EPSConv        conv;             /* convergence test */
 64:   PetscErrorCode (*conv_func)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
 65:   void           *conv_ctx;
 66:   EPSWhich       which;            /* which part of the spectrum to be sought */
 67:   PetscBool      leftvecs;         /* if left eigenvectors are requested */
 68:   PetscErrorCode (*which_func)(EPS,PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*);
 69:   void           *which_ctx;
 70:   PetscReal      inta, intb;       /* interval [a,b] for spectrum slicing */
 71:   EPSProblemType problem_type;     /* which kind of problem to be solved */
 72:   EPSExtraction  extraction;       /* which kind of extraction to be applied */
 73:   EPSBalance     balance;          /* the balancing method */
 74:   PetscInt       balance_its;      /* number of iterations of the balancing method */
 75:   PetscReal      balance_cutoff;   /* cutoff value for balancing */
 76:   PetscReal      nrma, nrmb;       /* matrix norms */
 77:   PetscBool      adaptive;         /* whether matrix norms are adaptively improved */
 78:   PetscBool      trueres;          /* whether the true residual norm must be computed */
 79:   PetscBool      trackall;         /* whether all the residuals must be computed */

 81:   /*------------------------- Working data --------------------------*/
 82:   Vec         D,                /* diagonal matrix for balancing */
 83:               *V,               /* set of basis vectors and computed eigenvectors */
 84:               *W,               /* set of left basis vectors and computed left eigenvectors */
 85:               *IS, *ISL,        /* placeholder for references to user-provided initial space */
 86:               *DS;              /* deflation space */
 87:   PetscScalar *eigr, *eigi,     /* real and imaginary parts of eigenvalues */
 88:               *T, *Tl;          /* projected matrices */
 89:   PetscReal   *errest,          /* error estimates */
 90:               *errest_left;     /* left error estimates */
 91:   ST          OP;               /* spectral transformation object */
 92:   IP          ip;               /* innerproduct object */
 93:   void        *data;            /* placeholder for misc stuff associated 
 94:                                    with a particular solver */
 95:   PetscInt    nconv,            /* number of converged eigenvalues */
 96:               its,              /* number of iterations so far computed */
 97:               *perm,            /* permutation for eigenvalue ordering */
 98:               nv,               /* size of current Schur decomposition */
 99:               n, nloc,          /* problem dimensions (global, local) */
100:               allocated_ncv;    /* number of basis vectors allocated */
101:   PetscBool   evecsavailable;   /* computed eigenvectors */
102:   PetscRandom rand;             /* random number generator */
103:   Vec         t;                /* template vector */

105:   /* ---------------- Default work-area and status vars -------------------- */
106:   PetscInt   nwork;
107:   Vec        *work;

109:   PetscBool  ds_ortho;         /* if DS vectors have been stored and orthonormalized */
110:   PetscInt   setupcalled;
111:   PetscBool  isgeneralized,
112:              ispositive,
113:              ishermitian;
114:   EPSConvergedReason reason;

116:   PetscErrorCode (*monitor[MAXEPSMONITORS])(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
117:   PetscErrorCode (*monitordestroy[MAXEPSMONITORS])(void**);
118:   void       *monitorcontext[MAXEPSMONITORS];
119:   PetscInt    numbermonitors;
120: };

122: extern PetscErrorCode EPSMonitor(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt);

124: extern PetscErrorCode EPSReset_Default(EPS);
125: extern PetscErrorCode EPSDefaultGetWork(EPS,PetscInt);
126: extern PetscErrorCode EPSDefaultFreeWork(EPS);
127: extern PetscErrorCode EPSDefaultSetWhich(EPS);
128: extern PetscErrorCode EPSAllocateSolution(EPS);
129: extern PetscErrorCode EPSFreeSolution(EPS);
130: extern PetscErrorCode EPSBackTransform_Default(EPS);
131: extern PetscErrorCode EPSComputeVectors_Default(EPS);
132: extern PetscErrorCode EPSComputeVectors_Hermitian(EPS);
133: extern PetscErrorCode EPSComputeVectors_Schur(EPS);
134: extern PetscErrorCode EPSComputeResidualNorm_Private(EPS,PetscScalar,PetscScalar,Vec,Vec,PetscReal*);
135: extern PetscErrorCode EPSComputeRelativeError_Private(EPS,PetscScalar,PetscScalar,Vec,Vec,PetscReal*);
136: extern PetscErrorCode EPSComputeTrueResidual(EPS,PetscScalar,PetscScalar,PetscScalar*,Vec*,PetscInt,PetscReal*);

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

140: extern PetscErrorCode EPSBasicArnoldi(EPS,PetscBool,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscBool*);
141: extern PetscErrorCode EPSDelayedArnoldi(EPS,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscBool*);
142: extern PetscErrorCode EPSDelayedArnoldi1(EPS,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscBool*);
143: extern PetscErrorCode EPSKrylovConvergence(EPS,PetscBool,PetscBool,PetscInt,PetscInt,PetscScalar*,PetscInt,PetscScalar*,Vec*,PetscInt,PetscReal,PetscReal,PetscInt*,PetscScalar*);
144: extern PetscErrorCode EPSFullLanczos(EPS,PetscReal*,PetscReal*,Vec*,PetscInt,PetscInt*,Vec,PetscBool*);
145: extern PetscErrorCode EPSTranslateHarmonic(PetscInt,PetscScalar*,PetscInt,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*);
146: extern PetscErrorCode EPSBuildBalance_Krylov(EPS);
147: extern PetscErrorCode EPSProjectedKSNonsym(EPS,PetscInt,PetscScalar*,PetscInt,PetscScalar*,PetscInt);

149: #endif