Actual source code: slepcpep.h

slepc-3.5.2 2014-10-10
Report Typos and Errors
  1: /*
  2:    User interface for SLEPc's polynomial eigenvalue solvers.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.

 10:    SLEPc is free software: you can redistribute it and/or modify it under  the
 11:    terms of version 3 of the GNU Lesser General Public License as published by
 12:    the Free Software Foundation.

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

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

 26: #include <slepceps.h>

 28: PETSC_EXTERN PetscErrorCode PEPInitializePackage(void);

 30: /*S
 31:      PEP - Abstract SLEPc object that manages all the polynomial eigenvalue
 32:      problem solvers.

 34:    Level: beginner

 36: .seealso:  PEPCreate()
 37: S*/
 38: typedef struct _p_PEP* PEP;

 40: /*J
 41:     PEPType - String with the name of a polynomial eigensolver

 43:    Level: beginner

 45: .seealso: PEPSetType(), PEP
 46: J*/
 47: typedef const char* PEPType;
 48: #define PEPLINEAR    "linear"
 49: #define PEPQARNOLDI  "qarnoldi"
 50: #define PEPTOAR      "toar"

 52: /* Logging support */
 53: PETSC_EXTERN PetscClassId PEP_CLASSID;

 55: /*E
 56:     PEPProblemType - Determines the type of the polynomial eigenproblem

 58:     Level: intermediate

 60: .seealso: PEPSetProblemType(), PEPGetProblemType()
 61: E*/
 62: typedef enum { PEP_GENERAL=1,
 63:                PEP_HERMITIAN,   /* All A_i  Hermitian */
 64:                PEP_GYROSCOPIC   /* QEP with M, K  Hermitian, M>0, C skew-Hermitian */
 65:              } PEPProblemType;

 67: /*E
 68:     PEPWhich - Determines which part of the spectrum is requested

 70:     Level: intermediate

 72: .seealso: PEPSetWhichEigenpairs(), PEPGetWhichEigenpairs()
 73: E*/
 74: typedef enum { PEP_LARGEST_MAGNITUDE=1,
 75:                PEP_SMALLEST_MAGNITUDE,
 76:                PEP_LARGEST_REAL,
 77:                PEP_SMALLEST_REAL,
 78:                PEP_LARGEST_IMAGINARY,
 79:                PEP_SMALLEST_IMAGINARY,
 80:                PEP_TARGET_MAGNITUDE,
 81:                PEP_TARGET_REAL,
 82:                PEP_TARGET_IMAGINARY} PEPWhich;

 84: /*E
 85:     PEPBasis - The type of polynomial basis used to represent the polynomial
 86:     eigenproblem

 88:     Level: intermediate

 90: .seealso: PEPSetBasis()
 91: E*/
 92: typedef enum { PEP_BASIS_MONOMIAL,
 93:                PEP_BASIS_CHEBYSHEV1,
 94:                PEP_BASIS_CHEBYSHEV2,
 95:                PEP_BASIS_LEGENDRE,
 96:                PEP_BASIS_LAGUERRE,
 97:                PEP_BASIS_HERMITE } PEPBasis;
 98: PETSC_EXTERN const char *PEPBasisTypes[];

100: /*E
101:     PEPScale - The scaling strategy

103:     Level: intermediate

105: .seealso: PEPSetScale()
106: E*/
107: typedef enum { PEP_SCALE_NONE,
108:                PEP_SCALE_SCALAR,
109:                PEP_SCALE_DIAGONAL,
110:                PEP_SCALE_BOTH } PEPScale;
111: PETSC_EXTERN const char *PEPScaleTypes[];

113: /*E
114:     PEPRefine - The refinement type

116:     Level: intermediate

118: .seealso: PEPSetRefine()
119: E*/
120: typedef enum { PEP_REFINE_NONE,
121:                PEP_REFINE_SIMPLE,
122:                PEP_REFINE_MULTIPLE } PEPRefine;
123: PETSC_EXTERN const char *PEPRefineTypes[];

125: /*E
126:     PEPExtract - The extraction type

128:     Level: intermediate

130: .seealso: PEPSetExtract()
131: E*/
132: typedef enum { PEP_EXTRACT_NORM,
133:                PEP_EXTRACT_RESIDUAL,
134:                PEP_EXTRACT_STRUCTURED } PEPExtract;
135: PETSC_EXTERN const char *PEPExtractTypes[];

137: /*E
138:     PEPConv - Determines the convergence test

140:     Level: intermediate

142: .seealso: PEPSetConvergenceTest(), PEPSetConvergenceTestFunction()
143: E*/
144: typedef enum { PEP_CONV_ABS,
145:                PEP_CONV_EIG,
146:                PEP_CONV_NORM,
147:                PEP_CONV_USER } PEPConv;

149: PETSC_EXTERN PetscErrorCode PEPCreate(MPI_Comm,PEP*);
150: PETSC_EXTERN PetscErrorCode PEPDestroy(PEP*);
151: PETSC_EXTERN PetscErrorCode PEPReset(PEP);
152: PETSC_EXTERN PetscErrorCode PEPSetType(PEP,PEPType);
153: PETSC_EXTERN PetscErrorCode PEPGetType(PEP,PEPType*);
154: PETSC_EXTERN PetscErrorCode PEPSetProblemType(PEP,PEPProblemType);
155: PETSC_EXTERN PetscErrorCode PEPGetProblemType(PEP,PEPProblemType*);
156: PETSC_EXTERN PetscErrorCode PEPSetOperators(PEP,PetscInt,Mat[]);
157: PETSC_EXTERN PetscErrorCode PEPGetOperators(PEP,PetscInt,Mat*);
158: PETSC_EXTERN PetscErrorCode PEPGetNumMatrices(PEP,PetscInt*);
159: PETSC_EXTERN PetscErrorCode PEPSetTarget(PEP,PetscScalar);
160: PETSC_EXTERN PetscErrorCode PEPGetTarget(PEP,PetscScalar*);
161: PETSC_EXTERN PetscErrorCode PEPSetFromOptions(PEP);
162: PETSC_EXTERN PetscErrorCode PEPSetUp(PEP);
163: PETSC_EXTERN PetscErrorCode PEPSolve(PEP);
164: PETSC_EXTERN PetscErrorCode PEPView(PEP,PetscViewer);
165: PETSC_EXTERN PetscErrorCode PEPPrintSolution(PEP,PetscViewer);
166: PETSC_EXTERN PetscErrorCode PEPSetBV(PEP,BV);
167: PETSC_EXTERN PetscErrorCode PEPGetBV(PEP,BV*);
168: PETSC_EXTERN PetscErrorCode PEPSetRG(PEP,RG);
169: PETSC_EXTERN PetscErrorCode PEPGetRG(PEP,RG*);
170: PETSC_EXTERN PetscErrorCode PEPSetDS(PEP,DS);
171: PETSC_EXTERN PetscErrorCode PEPGetDS(PEP,DS*);
172: PETSC_EXTERN PetscErrorCode PEPSetST(PEP,ST);
173: PETSC_EXTERN PetscErrorCode PEPGetST(PEP,ST*);

175: PETSC_EXTERN PetscErrorCode PEPSetTolerances(PEP,PetscReal,PetscInt);
176: PETSC_EXTERN PetscErrorCode PEPGetTolerances(PEP,PetscReal*,PetscInt*);
177: PETSC_EXTERN PetscErrorCode PEPSetConvergenceTestFunction(PEP,PetscErrorCode (*)(PEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void*,PetscErrorCode (*)(void*));
178: PETSC_EXTERN PetscErrorCode PEPSetConvergenceTest(PEP,PEPConv);
179: PETSC_EXTERN PetscErrorCode PEPGetConvergenceTest(PEP,PEPConv*);
180: PETSC_EXTERN PetscErrorCode PEPConvergedEigRelative(PEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
181: PETSC_EXTERN PetscErrorCode PEPConvergedNormRelative(PEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
182: PETSC_EXTERN PetscErrorCode PEPConvergedAbsolute(PEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
183: PETSC_EXTERN PetscErrorCode PEPSetDimensions(PEP,PetscInt,PetscInt,PetscInt);
184: PETSC_EXTERN PetscErrorCode PEPGetDimensions(PEP,PetscInt*,PetscInt*,PetscInt*);
185: PETSC_EXTERN PetscErrorCode PEPSetScale(PEP,PEPScale,PetscReal,PetscInt,PetscReal);
186: PETSC_EXTERN PetscErrorCode PEPGetScale(PEP,PEPScale*,PetscReal*,PetscInt*,PetscReal*);
187: PETSC_EXTERN PetscErrorCode PEPSetRefine(PEP,PEPRefine,PetscInt,PetscReal,PetscInt,PetscBool);
188: PETSC_EXTERN PetscErrorCode PEPGetRefine(PEP,PEPRefine*,PetscInt*,PetscReal*,PetscInt*,PetscBool*);
189: PETSC_EXTERN PetscErrorCode PEPSetExtract(PEP,PEPExtract);
190: PETSC_EXTERN PetscErrorCode PEPGetExtract(PEP,PEPExtract*);
191: PETSC_EXTERN PetscErrorCode PEPSetBasis(PEP,PEPBasis);
192: PETSC_EXTERN PetscErrorCode PEPGetBasis(PEP,PEPBasis*);

194: PETSC_EXTERN PetscErrorCode PEPGetConverged(PEP,PetscInt*);
195: PETSC_EXTERN PetscErrorCode PEPGetEigenpair(PEP,PetscInt,PetscScalar*,PetscScalar*,Vec,Vec);
196: PETSC_EXTERN PetscErrorCode PEPComputeRelativeError(PEP,PetscInt,PetscReal*);
197: PETSC_EXTERN PetscErrorCode PEPComputeResidualNorm(PEP,PetscInt,PetscReal*);
198: PETSC_EXTERN PetscErrorCode PEPGetErrorEstimate(PEP,PetscInt,PetscReal*);

200: PETSC_EXTERN PetscErrorCode PEPMonitor(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt);
201: PETSC_EXTERN PetscErrorCode PEPMonitorSet(PEP,PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void*,PetscErrorCode (*)(void**));
202: PETSC_EXTERN PetscErrorCode PEPMonitorCancel(PEP);
203: PETSC_EXTERN PetscErrorCode PEPGetMonitorContext(PEP,void **);
204: PETSC_EXTERN PetscErrorCode PEPGetIterationNumber(PEP,PetscInt*);

206: PETSC_EXTERN PetscErrorCode PEPSetInitialSpace(PEP,PetscInt,Vec*);
207: PETSC_EXTERN PetscErrorCode PEPSetWhichEigenpairs(PEP,PEPWhich);
208: PETSC_EXTERN PetscErrorCode PEPGetWhichEigenpairs(PEP,PEPWhich*);
209: PETSC_EXTERN PetscErrorCode PEPSetEigenvalueComparison(PEP,PetscErrorCode (*func)(PEP,PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void*);

211: PETSC_EXTERN PetscErrorCode PEPMonitorAll(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
212: PETSC_EXTERN PetscErrorCode PEPMonitorFirst(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
213: PETSC_EXTERN PetscErrorCode PEPMonitorConverged(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
214: PETSC_EXTERN PetscErrorCode PEPMonitorLG(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
215: PETSC_EXTERN PetscErrorCode PEPMonitorLGAll(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);

217: PETSC_EXTERN PetscErrorCode PEPSetTrackAll(PEP,PetscBool);
218: PETSC_EXTERN PetscErrorCode PEPGetTrackAll(PEP,PetscBool*);

220: PETSC_EXTERN PetscErrorCode PEPSetOptionsPrefix(PEP,const char*);
221: PETSC_EXTERN PetscErrorCode PEPAppendOptionsPrefix(PEP,const char*);
222: PETSC_EXTERN PetscErrorCode PEPGetOptionsPrefix(PEP,const char*[]);

224: /*E
225:     PEPConvergedReason - Reason an eigensolver was said to
226:          have converged or diverged

228:     Level: beginner

230: .seealso: PEPSolve(), PEPGetConvergedReason(), PEPSetTolerances()
231: E*/
232: typedef enum {/* converged */
233:               PEP_CONVERGED_TOL                =  2,
234:               /* diverged */
235:               PEP_DIVERGED_ITS                 = -3,
236:               PEP_DIVERGED_BREAKDOWN           = -4,
237:               PEP_CONVERGED_ITERATING          =  0} PEPConvergedReason;

239: PETSC_EXTERN PetscErrorCode PEPGetConvergedReason(PEP,PEPConvergedReason *);

241: PETSC_EXTERN PetscFunctionList PEPList;
242: PETSC_EXTERN PetscBool         PEPRegisterAllCalled;
243: PETSC_EXTERN PetscErrorCode PEPRegisterAll(void);
244: PETSC_EXTERN PetscErrorCode PEPRegister(const char[],PetscErrorCode(*)(PEP));

246: PETSC_EXTERN PetscErrorCode PEPSetWorkVecs(PEP,PetscInt);
247: PETSC_EXTERN PetscErrorCode PEPAllocateSolution(PEP,PetscInt);

249: /* --------- options specific to particular eigensolvers -------- */

251: PETSC_EXTERN PetscErrorCode PEPLinearSetCompanionForm(PEP,PetscInt);
252: PETSC_EXTERN PetscErrorCode PEPLinearGetCompanionForm(PEP,PetscInt*);
253: PETSC_EXTERN PetscErrorCode PEPLinearSetExplicitMatrix(PEP,PetscBool);
254: PETSC_EXTERN PetscErrorCode PEPLinearGetExplicitMatrix(PEP,PetscBool*);
255: PETSC_EXTERN PetscErrorCode PEPLinearSetEPS(PEP,EPS);
256: PETSC_EXTERN PetscErrorCode PEPLinearGetEPS(PEP,EPS*);

258: PETSC_EXTERN PetscErrorCode PEPQArnoldiSetRestart(PEP,PetscReal);
259: PETSC_EXTERN PetscErrorCode PEPQArnoldiGetRestart(PEP,PetscReal*);

261: PETSC_EXTERN PetscErrorCode PEPTOARSetRestart(PEP,PetscReal);
262: PETSC_EXTERN PetscErrorCode PEPTOARGetRestart(PEP,PetscReal*);

264: #endif