Actual source code: slepcpep.h

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

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    User interface for SLEPc's polynomial eigenvalue solvers
 12: */

 14: #pragma once

 16: #include <slepceps.h>

 18: /* SUBMANSEC = PEP */

 20: SLEPC_EXTERN PetscErrorCode PEPInitializePackage(void);
 21: SLEPC_EXTERN PetscErrorCode PEPFinalizePackage(void);

 23: /*S
 24:    PEP - SLEPc object that manages all the polynomial eigenvalue problem solvers.

 26:    Level: beginner

 28: .seealso: [](ch:pep), `PEPCreate()`
 29: S*/
 30: typedef struct _p_PEP* PEP;

 32: /*J
 33:    PEPType - String with the name of a polynomial eigensolver.

 35:    Level: beginner

 37: .seealso: [](ch:pep), `PEPSetType()`, `PEP`
 38: J*/
 39: typedef const char *PEPType;
 40: #define PEPTOAR      "toar"
 41: #define PEPSTOAR     "stoar"
 42: #define PEPQARNOLDI  "qarnoldi"
 43: #define PEPLINEAR    "linear"
 44: #define PEPJD        "jd"
 45: #define PEPCISS      "ciss"

 47: /* Logging support */
 48: SLEPC_EXTERN PetscClassId PEP_CLASSID;

 50: /*E
 51:    PEPProblemType - Determines the type of the polynomial eigenproblem

 53:    PEP_HERMITIAN is used when all A_i are Hermitian,
 54:    PEP_HYPERBOLIC is reserved for a QEP with Hermitian matrices, M>0, (x'Cx)^2 > 4(x'Mx)(x'Kx),
 55:    PEP_GYROSCOPIC is for aQEP with M, K  Hermitian, M>0, C skew-Hermitian.

 57:    Level: intermediate

 59: .seealso: [](ch:pep), `PEPSetProblemType()`, `PEPGetProblemType()`
 60: E*/
 61: typedef enum { PEP_GENERAL    = 1,
 62:                PEP_HERMITIAN  = 2,
 63:                PEP_HYPERBOLIC = 3,
 64:                PEP_GYROSCOPIC = 4
 65:              } PEPProblemType;

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

 70:    Level: intermediate

 72: .seealso: [](ch:pep), `PEPSetWhichEigenpairs()`, `PEPGetWhichEigenpairs()`
 73: E*/
 74: typedef enum { PEP_LARGEST_MAGNITUDE  = 1,
 75:                PEP_SMALLEST_MAGNITUDE = 2,
 76:                PEP_LARGEST_REAL       = 3,
 77:                PEP_SMALLEST_REAL      = 4,
 78:                PEP_LARGEST_IMAGINARY  = 5,
 79:                PEP_SMALLEST_IMAGINARY = 6,
 80:                PEP_TARGET_MAGNITUDE   = 7,
 81:                PEP_TARGET_REAL        = 8,
 82:                PEP_TARGET_IMAGINARY   = 9,
 83:                PEP_ALL                = 10,
 84:                PEP_WHICH_USER         = 11 } PEPWhich;

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

 90:    Level: intermediate

 92: .seealso: [](ch:pep), `PEPSetBasis()`
 93: E*/
 94: typedef enum { PEP_BASIS_MONOMIAL,
 95:                PEP_BASIS_CHEBYSHEV1,
 96:                PEP_BASIS_CHEBYSHEV2,
 97:                PEP_BASIS_LEGENDRE,
 98:                PEP_BASIS_LAGUERRE,
 99:                PEP_BASIS_HERMITE } PEPBasis;
100: SLEPC_EXTERN const char *PEPBasisTypes[];

102: /*E
103:    PEPScale - The scaling strategy

105:    Level: intermediate

107: .seealso: [](ch:pep), `PEPSetScale()`
108: E*/
109: typedef enum { PEP_SCALE_NONE,
110:                PEP_SCALE_SCALAR,
111:                PEP_SCALE_DIAGONAL,
112:                PEP_SCALE_BOTH } PEPScale;
113: SLEPC_EXTERN const char *PEPScaleTypes[];

115: /*E
116:    PEPRefine - The refinement type

118:    Level: intermediate

120: .seealso: [](ch:pep), `PEPSetRefine()`
121: E*/
122: typedef enum { PEP_REFINE_NONE,
123:                PEP_REFINE_SIMPLE,
124:                PEP_REFINE_MULTIPLE } PEPRefine;
125: SLEPC_EXTERN const char *PEPRefineTypes[];

127: /*E
128:    PEPRefineScheme - The scheme used for solving linear systems during iterative refinement

130:    Level: intermediate

132: .seealso: [](ch:pep), `PEPSetRefine()`
133: E*/
134: typedef enum { PEP_REFINE_SCHEME_SCHUR    = 1,
135:                PEP_REFINE_SCHEME_MBE      = 2,
136:                PEP_REFINE_SCHEME_EXPLICIT = 3 } PEPRefineScheme;
137: SLEPC_EXTERN const char *PEPRefineSchemes[];

139: /*E
140:    PEPExtract - The extraction type

142:    Level: intermediate

144: .seealso: [](ch:pep), `PEPSetExtract()`
145: E*/
146: typedef enum { PEP_EXTRACT_NONE       = 1,
147:                PEP_EXTRACT_NORM       = 2,
148:                PEP_EXTRACT_RESIDUAL   = 3,
149:                PEP_EXTRACT_STRUCTURED = 4 } PEPExtract;
150: SLEPC_EXTERN const char *PEPExtractTypes[];

152: /*E
153:    PEPErrorType - The error type used to assess accuracy of computed solutions

155:    Level: intermediate

157: .seealso: [](ch:pep), `PEPComputeError()`
158: E*/
159: typedef enum { PEP_ERROR_ABSOLUTE,
160:                PEP_ERROR_RELATIVE,
161:                PEP_ERROR_BACKWARD } PEPErrorType;
162: SLEPC_EXTERN const char *PEPErrorTypes[];

164: /*E
165:    PEPConv - Determines the convergence test

167:    Level: intermediate

169: .seealso: [](ch:pep), `PEPSetConvergenceTest()`, `PEPSetConvergenceTestFunction()`
170: E*/
171: typedef enum { PEP_CONV_ABS,
172:                PEP_CONV_REL,
173:                PEP_CONV_NORM,
174:                PEP_CONV_USER } PEPConv;

176: /*E
177:    PEPStop - Determines the stopping test

179:    Level: advanced

181: .seealso: [](ch:pep), `PEPSetStoppingTest()`, `PEPSetStoppingTestFunction()`
182: E*/
183: typedef enum { PEP_STOP_BASIC,
184:                PEP_STOP_USER } PEPStop;

186: /*E
187:    PEPConvergedReason - Reason an eigensolver was determined to have converged
188:    or diverged.

190:    Values:
191: +  `PEP_CONVERGED_TOL`          - converged up to tolerance
192: .  `PEP_CONVERGED_USER`         - converged due to a user-defined condition
193: .  `PEP_DIVERGED_ITS`           - exceeded the maximum number of allowed iterations
194: .  `PEP_DIVERGED_BREAKDOWN`     - generic breakdown in method
195: -  `PEP_DIVERGED_SYMMETRY_LOST` - pseudo-Lanczos was not able to keep symmetry

197:    Level: intermediate

199: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPSetTolerances()`
200: E*/
201: typedef enum {/* converged */
202:               PEP_CONVERGED_TOL                =  1,
203:               PEP_CONVERGED_USER               =  2,
204:               /* diverged */
205:               PEP_DIVERGED_ITS                 = -1,
206:               PEP_DIVERGED_BREAKDOWN           = -2,
207:               PEP_DIVERGED_SYMMETRY_LOST       = -3,
208:               PEP_CONVERGED_ITERATING          =  0} PEPConvergedReason;
209: SLEPC_EXTERN const char *const*PEPConvergedReasons;

211: SLEPC_EXTERN PetscErrorCode PEPCreate(MPI_Comm,PEP*);
212: SLEPC_EXTERN PetscErrorCode PEPDestroy(PEP*);
213: SLEPC_EXTERN PetscErrorCode PEPReset(PEP);
214: SLEPC_EXTERN PetscErrorCode PEPSetType(PEP,PEPType);
215: SLEPC_EXTERN PetscErrorCode PEPGetType(PEP,PEPType*);
216: SLEPC_EXTERN PetscErrorCode PEPSetProblemType(PEP,PEPProblemType);
217: SLEPC_EXTERN PetscErrorCode PEPGetProblemType(PEP,PEPProblemType*);
218: SLEPC_EXTERN PetscErrorCode PEPSetOperators(PEP,PetscInt,Mat[]);
219: SLEPC_EXTERN PetscErrorCode PEPGetOperators(PEP,PetscInt,Mat*);
220: SLEPC_EXTERN PetscErrorCode PEPGetNumMatrices(PEP,PetscInt*);
221: SLEPC_EXTERN PetscErrorCode PEPSetTarget(PEP,PetscScalar);
222: SLEPC_EXTERN PetscErrorCode PEPGetTarget(PEP,PetscScalar*);
223: SLEPC_EXTERN PetscErrorCode PEPSetInterval(PEP,PetscReal,PetscReal);
224: SLEPC_EXTERN PetscErrorCode PEPGetInterval(PEP,PetscReal*,PetscReal*);
225: SLEPC_EXTERN PetscErrorCode PEPSetFromOptions(PEP);
226: SLEPC_EXTERN PetscErrorCode PEPSetDSType(PEP);
227: SLEPC_EXTERN PetscErrorCode PEPSetUp(PEP);
228: SLEPC_EXTERN PetscErrorCode PEPSolve(PEP);
229: SLEPC_EXTERN PetscErrorCode PEPView(PEP,PetscViewer);
230: SLEPC_EXTERN PetscErrorCode PEPViewFromOptions(PEP,PetscObject,const char[]);
231: SLEPC_EXTERN PetscErrorCode PEPErrorView(PEP,PEPErrorType,PetscViewer);
232: PETSC_DEPRECATED_FUNCTION(3, 6, 0, "PEPErrorView()", ) static inline PetscErrorCode PEPPrintSolution(PEP pep,PetscViewer v) {return PEPErrorView(pep,PEP_ERROR_BACKWARD,v);}
233: SLEPC_EXTERN PetscErrorCode PEPErrorViewFromOptions(PEP);
234: SLEPC_EXTERN PetscErrorCode PEPConvergedReasonView(PEP,PetscViewer);
235: SLEPC_EXTERN PetscErrorCode PEPConvergedReasonViewFromOptions(PEP);
236: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "PEPConvergedReasonView()", ) static inline PetscErrorCode PEPReasonView(PEP pep,PetscViewer v) {return PEPConvergedReasonView(pep,v);}
237: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "PEPConvergedReasonViewFromOptions()", ) static inline PetscErrorCode PEPReasonViewFromOptions(PEP pep) {return PEPConvergedReasonViewFromOptions(pep);}
238: SLEPC_EXTERN PetscErrorCode PEPValuesView(PEP,PetscViewer);
239: SLEPC_EXTERN PetscErrorCode PEPValuesViewFromOptions(PEP);
240: SLEPC_EXTERN PetscErrorCode PEPVectorsView(PEP,PetscViewer);
241: SLEPC_EXTERN PetscErrorCode PEPVectorsViewFromOptions(PEP);
242: SLEPC_EXTERN PetscErrorCode PEPSetBV(PEP,BV);
243: SLEPC_EXTERN PetscErrorCode PEPGetBV(PEP,BV*);
244: SLEPC_EXTERN PetscErrorCode PEPSetRG(PEP,RG);
245: SLEPC_EXTERN PetscErrorCode PEPGetRG(PEP,RG*);
246: SLEPC_EXTERN PetscErrorCode PEPSetDS(PEP,DS);
247: SLEPC_EXTERN PetscErrorCode PEPGetDS(PEP,DS*);
248: SLEPC_EXTERN PetscErrorCode PEPSetST(PEP,ST);
249: SLEPC_EXTERN PetscErrorCode PEPGetST(PEP,ST*);
250: SLEPC_EXTERN PetscErrorCode PEPRefineGetKSP(PEP,KSP*);

252: SLEPC_EXTERN PetscErrorCode PEPSetTolerances(PEP,PetscReal,PetscInt);
253: SLEPC_EXTERN PetscErrorCode PEPGetTolerances(PEP,PetscReal*,PetscInt*);
254: SLEPC_EXTERN PetscErrorCode PEPGetConvergedReason(PEP,PEPConvergedReason*);

256: SLEPC_EXTERN PetscErrorCode PEPSetDimensions(PEP,PetscInt,PetscInt,PetscInt);
257: SLEPC_EXTERN PetscErrorCode PEPGetDimensions(PEP,PetscInt*,PetscInt*,PetscInt*);
258: SLEPC_EXTERN PetscErrorCode PEPSetScale(PEP,PEPScale,PetscReal,Vec,Vec,PetscInt,PetscReal);
259: SLEPC_EXTERN PetscErrorCode PEPGetScale(PEP,PEPScale*,PetscReal*,Vec*,Vec*,PetscInt*,PetscReal*);
260: SLEPC_EXTERN PetscErrorCode PEPSetRefine(PEP,PEPRefine,PetscInt,PetscReal,PetscInt,PEPRefineScheme);
261: SLEPC_EXTERN PetscErrorCode PEPGetRefine(PEP,PEPRefine*,PetscInt*,PetscReal*,PetscInt*,PEPRefineScheme*);
262: SLEPC_EXTERN PetscErrorCode PEPSetExtract(PEP,PEPExtract);
263: SLEPC_EXTERN PetscErrorCode PEPGetExtract(PEP,PEPExtract*);
264: SLEPC_EXTERN PetscErrorCode PEPSetBasis(PEP,PEPBasis);
265: SLEPC_EXTERN PetscErrorCode PEPGetBasis(PEP,PEPBasis*);

267: SLEPC_EXTERN PetscErrorCode PEPGetConverged(PEP,PetscInt*);
268: SLEPC_EXTERN PetscErrorCode PEPGetEigenpair(PEP,PetscInt,PetscScalar*,PetscScalar*,Vec,Vec);
269: SLEPC_EXTERN PetscErrorCode PEPComputeError(PEP,PetscInt,PEPErrorType,PetscReal*);
270: PETSC_DEPRECATED_FUNCTION(3, 6, 0, "PEPComputeError()", ) static inline PetscErrorCode PEPComputeRelativeError(PEP pep,PetscInt i,PetscReal *r) {return PEPComputeError(pep,i,PEP_ERROR_BACKWARD,r);}
271: PETSC_DEPRECATED_FUNCTION(3, 6, 0, "PEPComputeError() with PEP_ERROR_ABSOLUTE", ) static inline PetscErrorCode PEPComputeResidualNorm(PEP pep,PetscInt i,PetscReal *r) {return PEPComputeError(pep,i,PEP_ERROR_ABSOLUTE,r);}
272: SLEPC_EXTERN PetscErrorCode PEPGetErrorEstimate(PEP,PetscInt,PetscReal*);
273: SLEPC_EXTERN PetscErrorCode PEPGetIterationNumber(PEP,PetscInt*);

275: SLEPC_EXTERN PetscErrorCode PEPSetInitialSpace(PEP,PetscInt,Vec[]);
276: SLEPC_EXTERN PetscErrorCode PEPSetWhichEigenpairs(PEP,PEPWhich);
277: SLEPC_EXTERN PetscErrorCode PEPGetWhichEigenpairs(PEP,PEPWhich*);

279: SLEPC_EXTERN PetscErrorCode PEPSetTrackAll(PEP,PetscBool);
280: SLEPC_EXTERN PetscErrorCode PEPGetTrackAll(PEP,PetscBool*);

282: /*S
283:    PEPMonitorFn - A function prototype for functions provided to `PEPMonitorSet()`.

285:    Calling Sequence:
286: +  pep    - the polynomial eigensolver context
287: .  its    - iteration number
288: .  nconv  - number of converged eigenpairs
289: .  eigr   - real part of the eigenvalues
290: .  eigi   - imaginary part of the eigenvalues
291: .  errest - relative error estimates for each eigenpair
292: .  nest   - number of error estimates
293: -  ctx    - optional monitoring context, as provided with `PEPMonitorSet()`

295:    Level: intermediate

297: .seealso: [](ch:pep), `PEPMonitorSet()`
298: S*/
299: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorFn(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,void *ctx);

301: /*S
302:    PEPMonitorRegisterFn - A function prototype for functions provided to PEPMonitorRegister()

304:    Calling Sequence:
305: +  pep    - the polynomial eigensolver context
306: .  its    - iteration number
307: .  nconv  - number of converged eigenpairs
308: .  eigr   - real part of the eigenvalues
309: .  eigi   - imaginary part of the eigenvalues
310: .  errest - relative error estimates for each eigenpair
311: .  nest   - number of error estimates
312: -  ctx    - `PetscViewerAndFormat` object

314:    Level: advanced

316:    Note:
317:    This is a `PEPMonitorFn` specialized for a context of `PetscViewerAndFormat`.

319: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorRegister()`, `PEPMonitorFn`, `PEPMonitorRegisterCreateFn`, `PEPMonitorRegisterDestroyFn`
320: S*/
321: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorRegisterFn(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *ctx);

323: /*S
324:    PEPMonitorRegisterCreateFn - A function prototype for functions that do the
325:    creation when provided to `PEPMonitorRegister()`.

327:    Calling Sequence:
328: +  viewer - the viewer to be used with the `PEPMonitorRegisterFn`
329: .  format - the format of the viewer
330: .  ctx    - a context for the monitor
331: -  result - a `PetscViewerAndFormat` object

333:    Level: advanced

335: .seealso: [](ch:pep), `PEPMonitorRegisterFn`, `PEPMonitorSet()`, `PEPMonitorRegister()`, `PEPMonitorFn`, `PEPMonitorRegisterDestroyFn`
336: S*/
337: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorRegisterCreateFn(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **result);

339: /*S
340:    PEPMonitorRegisterDestroyFn - A function prototype for functions that do the after
341:    use destruction when provided to `PEPMonitorRegister()`.

343:    Calling Sequence:
344: .  vf - a `PetscViewerAndFormat` object to be destroyed, including any context

346:    Level: advanced

348: .seealso: [](ch:pep), `PEPMonitorRegisterFn`, `PEPMonitorSet()`, `PEPMonitorRegister()`, `PEPMonitorFn`, `PEPMonitorRegisterCreateFn`
349: S*/
350: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorRegisterDestroyFn(PetscViewerAndFormat **result);

352: SLEPC_EXTERN PetscErrorCode PEPMonitor(PEP,PetscInt,PetscInt,PetscScalar[],PetscScalar[],PetscReal[],PetscInt);
353: SLEPC_EXTERN PetscErrorCode PEPMonitorSet(PEP,PEPMonitorFn,void*,PetscCtxDestroyFn*);
354: SLEPC_EXTERN PetscErrorCode PEPMonitorCancel(PEP);
355: SLEPC_EXTERN PetscErrorCode PEPGetMonitorContext(PEP,void*);

357: SLEPC_EXTERN PetscErrorCode PEPMonitorSetFromOptions(PEP,const char[],const char[],void*,PetscBool);
358: SLEPC_EXTERN PEPMonitorRegisterFn        PEPMonitorFirst;
359: SLEPC_EXTERN PEPMonitorRegisterFn        PEPMonitorFirstDrawLG;
360: SLEPC_EXTERN PEPMonitorRegisterCreateFn  PEPMonitorFirstDrawLGCreate;
361: SLEPC_EXTERN PEPMonitorRegisterFn        PEPMonitorAll;
362: SLEPC_EXTERN PEPMonitorRegisterFn        PEPMonitorAllDrawLG;
363: SLEPC_EXTERN PEPMonitorRegisterCreateFn  PEPMonitorAllDrawLGCreate;
364: SLEPC_EXTERN PEPMonitorRegisterFn        PEPMonitorConverged;
365: SLEPC_EXTERN PEPMonitorRegisterCreateFn  PEPMonitorConvergedCreate;
366: SLEPC_EXTERN PEPMonitorRegisterFn        PEPMonitorConvergedDrawLG;
367: SLEPC_EXTERN PEPMonitorRegisterCreateFn  PEPMonitorConvergedDrawLGCreate;
368: SLEPC_EXTERN PEPMonitorRegisterDestroyFn PEPMonitorConvergedDestroy;

370: SLEPC_EXTERN PetscErrorCode PEPSetOptionsPrefix(PEP,const char[]);
371: SLEPC_EXTERN PetscErrorCode PEPAppendOptionsPrefix(PEP,const char[]);
372: SLEPC_EXTERN PetscErrorCode PEPGetOptionsPrefix(PEP,const char*[]);

374: SLEPC_EXTERN PetscFunctionList PEPList;
375: SLEPC_EXTERN PetscFunctionList PEPMonitorList;
376: SLEPC_EXTERN PetscFunctionList PEPMonitorCreateList;
377: SLEPC_EXTERN PetscFunctionList PEPMonitorDestroyList;
378: SLEPC_EXTERN PetscErrorCode PEPRegister(const char[],PetscErrorCode(*)(PEP));
379: SLEPC_EXTERN PetscErrorCode PEPMonitorRegister(const char[],PetscViewerType,PetscViewerFormat,PEPMonitorRegisterFn*,PEPMonitorRegisterCreateFn*,PEPMonitorRegisterDestroyFn*);

381: SLEPC_EXTERN PetscErrorCode PEPSetWorkVecs(PEP,PetscInt);
382: SLEPC_EXTERN PetscErrorCode PEPAllocateSolution(PEP,PetscInt);

384: /*S
385:    PEPConvergenceTestFn - A prototype of a PEP convergence test function that
386:    would be passed to PEPSetConvergenceTestFunction()

388:    Calling Sequence:
389: +  pep    - the polynomial eigensolver context
390: .  eigr   - real part of the eigenvalue
391: .  eigi   - imaginary part of the eigenvalue
392: .  res    - residual norm associated to the eigenpair
393: .  errest - [output] computed error estimate
394: -  ctx    - [optional] user-defined context for private data for the
395:             convergence test routine (may be NULL)

397:    Level: advanced

399: .seealso: [](ch:pep), `PEPSetConvergenceTestFunction()`
400: S*/
401: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPConvergenceTestFn(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx);

403: SLEPC_EXTERN PetscErrorCode PEPSetConvergenceTest(PEP,PEPConv);
404: SLEPC_EXTERN PetscErrorCode PEPGetConvergenceTest(PEP,PEPConv*);
405: SLEPC_EXTERN PEPConvergenceTestFn PEPConvergedAbsolute;
406: SLEPC_EXTERN PEPConvergenceTestFn PEPConvergedRelative;
407: SLEPC_EXTERN PEPConvergenceTestFn PEPConvergedNorm;
408: SLEPC_EXTERN PetscErrorCode PEPSetConvergenceTestFunction(PEP,PEPConvergenceTestFn*,void*,PetscCtxDestroyFn*);

410: /*S
411:    PEPStoppingTestFn - A prototype of a PEP stopping test function that would be
412:    passed to PEPSetStoppingTestFunction()

414:    Calling Sequence:
415: +  pep    - the polynomial eigensolver context
416: .  its    - current number of iterations
417: .  max_it - maximum number of iterations
418: .  nconv  - number of currently converged eigenpairs
419: .  nev    - number of requested eigenpairs
420: .  reason - [output] result of the stopping test
421: -  ctx    - [optional] user-defined context for private data for the
422:             stopping test routine (may be NULL)

424:    Level: advanced

426: .seealso: [](ch:pep), `PEPSetStoppingTestFunction()`
427: S*/
428: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPStoppingTestFn(PEP pep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,PEPConvergedReason *reason,void *ctx);

430: SLEPC_EXTERN PetscErrorCode PEPSetStoppingTest(PEP,PEPStop);
431: SLEPC_EXTERN PetscErrorCode PEPGetStoppingTest(PEP,PEPStop*);
432: SLEPC_EXTERN PEPStoppingTestFn PEPStoppingBasic;
433: SLEPC_EXTERN PetscErrorCode PEPSetStoppingTestFunction(PEP,PEPStoppingTestFn*,void*,PetscCtxDestroyFn*);

435: SLEPC_EXTERN PetscErrorCode PEPSetEigenvalueComparison(PEP,SlepcEigenvalueComparisonFn*,void*);

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

439: SLEPC_EXTERN PetscErrorCode PEPLinearSetLinearization(PEP,PetscReal,PetscReal);
440: SLEPC_EXTERN PetscErrorCode PEPLinearGetLinearization(PEP,PetscReal*,PetscReal*);
441: SLEPC_EXTERN PetscErrorCode PEPLinearSetExplicitMatrix(PEP,PetscBool);
442: SLEPC_EXTERN PetscErrorCode PEPLinearGetExplicitMatrix(PEP,PetscBool*);
443: SLEPC_EXTERN PetscErrorCode PEPLinearSetEPS(PEP,EPS);
444: SLEPC_EXTERN PetscErrorCode PEPLinearGetEPS(PEP,EPS*);
445: PETSC_DEPRECATED_FUNCTION(3, 10, 0, "PEPLinearSetLinearization()", ) static inline PetscErrorCode PEPLinearSetCompanionForm(PEP pep,PetscInt cform) {return (cform==1)?PEPLinearSetLinearization(pep,1.0,0.0):PEPLinearSetLinearization(pep,0.0,1.0);}
446: PETSC_DEPRECATED_FUNCTION(3, 10, 0, "PEPLinearGetLinearization()", ) static inline PetscErrorCode PEPLinearGetCompanionForm(PEP pep,PetscInt *cform) {(void)pep; if (cform) *cform=1; return PETSC_SUCCESS;}

448: SLEPC_EXTERN PetscErrorCode PEPQArnoldiSetRestart(PEP,PetscReal);
449: SLEPC_EXTERN PetscErrorCode PEPQArnoldiGetRestart(PEP,PetscReal*);
450: SLEPC_EXTERN PetscErrorCode PEPQArnoldiSetLocking(PEP,PetscBool);
451: SLEPC_EXTERN PetscErrorCode PEPQArnoldiGetLocking(PEP,PetscBool*);

453: SLEPC_EXTERN PetscErrorCode PEPTOARSetRestart(PEP,PetscReal);
454: SLEPC_EXTERN PetscErrorCode PEPTOARGetRestart(PEP,PetscReal*);
455: SLEPC_EXTERN PetscErrorCode PEPTOARSetLocking(PEP,PetscBool);
456: SLEPC_EXTERN PetscErrorCode PEPTOARGetLocking(PEP,PetscBool*);

458: SLEPC_EXTERN PetscErrorCode PEPSTOARSetLinearization(PEP,PetscReal,PetscReal);
459: SLEPC_EXTERN PetscErrorCode PEPSTOARGetLinearization(PEP,PetscReal*,PetscReal*);
460: SLEPC_EXTERN PetscErrorCode PEPSTOARSetLocking(PEP,PetscBool);
461: SLEPC_EXTERN PetscErrorCode PEPSTOARGetLocking(PEP,PetscBool*);
462: SLEPC_EXTERN PetscErrorCode PEPSTOARSetDetectZeros(PEP,PetscBool);
463: SLEPC_EXTERN PetscErrorCode PEPSTOARGetDetectZeros(PEP,PetscBool*);
464: SLEPC_EXTERN PetscErrorCode PEPSTOARGetInertias(PEP,PetscInt*,PetscReal*[],PetscInt*[]);
465: SLEPC_EXTERN PetscErrorCode PEPSTOARSetDimensions(PEP,PetscInt,PetscInt,PetscInt);
466: SLEPC_EXTERN PetscErrorCode PEPSTOARGetDimensions(PEP,PetscInt*,PetscInt*,PetscInt*);
467: SLEPC_EXTERN PetscErrorCode PEPSTOARSetCheckEigenvalueType(PEP,PetscBool);
468: SLEPC_EXTERN PetscErrorCode PEPSTOARGetCheckEigenvalueType(PEP,PetscBool*);
469: SLEPC_EXTERN PetscErrorCode PEPCheckDefiniteQEP(PEP,PetscReal*,PetscReal*,PetscInt*,PetscInt*);

471: /*E
472:    PEPJDProjection - The type of projection to be used in the Jacobi-Davidson solver

474:    Level: intermediate

476: .seealso: [](ch:pep), `PEPJDSetProjection()`
477: E*/
478: typedef enum { PEP_JD_PROJECTION_HARMONIC,
479:                PEP_JD_PROJECTION_ORTHOGONAL } PEPJDProjection;
480: SLEPC_EXTERN const char *PEPJDProjectionTypes[];

482: SLEPC_EXTERN PetscErrorCode PEPJDSetRestart(PEP,PetscReal);
483: SLEPC_EXTERN PetscErrorCode PEPJDGetRestart(PEP,PetscReal*);
484: SLEPC_EXTERN PetscErrorCode PEPJDSetFix(PEP,PetscReal);
485: SLEPC_EXTERN PetscErrorCode PEPJDGetFix(PEP,PetscReal*);
486: SLEPC_EXTERN PetscErrorCode PEPJDSetReusePreconditioner(PEP,PetscBool);
487: SLEPC_EXTERN PetscErrorCode PEPJDGetReusePreconditioner(PEP,PetscBool*);
488: SLEPC_EXTERN PetscErrorCode PEPJDSetMinimalityIndex(PEP,PetscInt);
489: SLEPC_EXTERN PetscErrorCode PEPJDGetMinimalityIndex(PEP,PetscInt*);
490: SLEPC_EXTERN PetscErrorCode PEPJDSetProjection(PEP,PEPJDProjection);
491: SLEPC_EXTERN PetscErrorCode PEPJDGetProjection(PEP,PEPJDProjection*);

493: /*E
494:    PEPCISSExtraction - determines the extraction technique in the CISS solver

496:    Level: advanced

498: .seealso: [](ch:pep), `PEPCISSSetExtraction()`, `PEPCISSGetExtraction()`
499: E*/
500: typedef enum { PEP_CISS_EXTRACTION_RITZ,
501:                PEP_CISS_EXTRACTION_HANKEL,
502:                PEP_CISS_EXTRACTION_CAA    } PEPCISSExtraction;
503: SLEPC_EXTERN const char *PEPCISSExtractions[];

505: #if defined(PETSC_USE_COMPLEX) || defined(PETSC_CLANG_STATIC_ANALYZER)
506: SLEPC_EXTERN PetscErrorCode PEPCISSSetExtraction(PEP,PEPCISSExtraction);
507: SLEPC_EXTERN PetscErrorCode PEPCISSGetExtraction(PEP,PEPCISSExtraction*);
508: SLEPC_EXTERN PetscErrorCode PEPCISSSetSizes(PEP,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscBool);
509: SLEPC_EXTERN PetscErrorCode PEPCISSGetSizes(PEP,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscBool*);
510: SLEPC_EXTERN PetscErrorCode PEPCISSSetThreshold(PEP,PetscReal,PetscReal);
511: SLEPC_EXTERN PetscErrorCode PEPCISSGetThreshold(PEP,PetscReal*,PetscReal*);
512: SLEPC_EXTERN PetscErrorCode PEPCISSSetRefinement(PEP,PetscInt,PetscInt);
513: SLEPC_EXTERN PetscErrorCode PEPCISSGetRefinement(PEP,PetscInt*,PetscInt*);
514: SLEPC_EXTERN PetscErrorCode PEPCISSGetKSPs(PEP,PetscInt*,KSP*[]);
515: #else
516: #define SlepcPEPCISSUnavailable(pep) do { \
517:     PetscFunctionBegin; \
518:     SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"%s() not available with real scalars",PETSC_FUNCTION_NAME); \
519:     } while (0)
520: static inline PetscErrorCode PEPCISSSetExtraction(PEP pep,PETSC_UNUSED PEPCISSExtraction ex) {SlepcPEPCISSUnavailable(pep);}
521: static inline PetscErrorCode PEPCISSGetExtraction(PEP pep,PETSC_UNUSED PEPCISSExtraction *ex) {SlepcPEPCISSUnavailable(pep);}
522: static inline PetscErrorCode PEPCISSSetSizes(PEP pep,PETSC_UNUSED PetscInt ip,PETSC_UNUSED PetscInt bs,PETSC_UNUSED PetscInt ms,PETSC_UNUSED PetscInt npart,PETSC_UNUSED PetscInt bsmax,PETSC_UNUSED PetscBool realmats) {SlepcPEPCISSUnavailable(pep);}
523: static inline PetscErrorCode PEPCISSGetSizes(PEP pep,PETSC_UNUSED PetscInt *ip,PETSC_UNUSED PetscInt *bs,PETSC_UNUSED PetscInt *ms,PETSC_UNUSED PetscInt *npart,PETSC_UNUSED PetscInt *bsmak,PETSC_UNUSED PetscBool *realmats) {SlepcPEPCISSUnavailable(pep);}
524: static inline PetscErrorCode PEPCISSSetThreshold(PEP pep,PETSC_UNUSED PetscReal delta,PETSC_UNUSED PetscReal spur) {SlepcPEPCISSUnavailable(pep);}
525: static inline PetscErrorCode PEPCISSGetThreshold(PEP pep,PETSC_UNUSED PetscReal *delta,PETSC_UNUSED PetscReal *spur) {SlepcPEPCISSUnavailable(pep);}
526: static inline PetscErrorCode PEPCISSSetRefinement(PEP pep,PETSC_UNUSED PetscInt inner,PETSC_UNUSED PetscInt blsize) {SlepcPEPCISSUnavailable(pep);}
527: static inline PetscErrorCode PEPCISSGetRefinement(PEP pep,PETSC_UNUSED PetscInt *inner,PETSC_UNUSED PetscInt *blsize) {SlepcPEPCISSUnavailable(pep);}
528: static inline PetscErrorCode PEPCISSGetKSPs(PEP pep,PETSC_UNUSED PetscInt *nsolve,PETSC_UNUSED KSP *ksp[]) {SlepcPEPCISSUnavailable(pep);}
529: #undef SlepcPEPCISSUnavailable
530: #endif