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: Values:
54: + `PEP_GENERAL` - polynomial eigenproblem with no particular structure
55: . `PEP_HERMITIAN` - polynomial eigenproblem with all coefficient matrices Hermitian
56: . `PEP_HYPERBOLIC` - quadratic eigenproblem with hyperbolic structure
57: - `PEP_GYROSCOPIC` - quadratic eigenproblem with gyroscopic structure
59: Note:
60: By default, no particular structure is assumed (`PEP_GENERAL`).
62: Level: intermediate
64: .seealso: [](ch:pep), `PEPSetProblemType()`, `PEPGetProblemType()`
65: E*/
66: typedef enum { PEP_GENERAL = 1,
67: PEP_HERMITIAN = 2,
68: PEP_HYPERBOLIC = 3,
69: PEP_GYROSCOPIC = 4
70: } PEPProblemType;
72: /*MC
73: PEP_GENERAL - A polynomial eigenproblem with no particular structure.
75: Note:
76: This is the default problem type.
78: Level: intermediate
80: .seealso: [](ch:pep), `PEPProblemType`, `PEPSetProblemType()`, `PEP_HERMITIAN`, `PEP_HYPERBOLIC`, `PEP_GYROSCOPIC`
81: M*/
83: /*MC
84: PEP_HERMITIAN - A polynomial eigenvalue problem with all coefficient matrices Hermitian.
86: Notes:
87: This is used when all $A_i$ matrices passed in `PEPSetOperators()`
88: are Hermitian.
90: Currently there is support only for quadratic eigenvalue problems,
91: $(K+\lambda C+\lambda^2M)x=0$ with $K$, $C$, $M$ Hermitian.
93: Level: intermediate
95: .seealso: [](ch:pep), `PEPProblemType`, `PEPSetProblemType()`, `PEP_GENERAL`, `PEP_HYPERBOLIC`, `PEP_GYROSCOPIC`
96: M*/
98: /*MC
99: PEP_HYPERBOLIC - A quadratic eigenvalue problem with hyperbolic structure.
101: Note:
102: This is reserved for the case of a quadratic eigenvalue problem
103: $(K+\lambda C+\lambda^2M)x=0$ with Hermitian coefficient matrices, and in
104: addition $M$ is positive definite and $(x^*Cx)^2>4(x^*Mx)(x^*Kx)$ for all
105: nonzero $x\in\mathbb{C}^n$. All eigenvalues are real, and form two separate
106: groups of $n$ eigenvalues, each of them having linearly independent eigenvectors.
108: Level: intermediate
110: .seealso: [](ch:pep), `PEPProblemType`, `PEPSetProblemType()`, `PEP_GENERAL`, `PEP_HERMITIAN`, `PEP_GYROSCOPIC`
111: M*/
113: /*MC
114: PEP_GYROSCOPIC - A quadratic eigenvalue problem with gyroscopic structure.
116: Notes:
117: This is reserved for the case of a quadratic eigenvalue problem
118: $(K+\lambda C+\lambda^2M)x=0$ with $M$, $K$ Hermitian, $M>0$, and
119: $C$ skew-Hermitian.
121: Currently there is support for this problem type only in `PEPLINEAR`, using
122: a general eigensolver, without exploiting the structure.
124: Level: intermediate
126: .seealso: [](ch:pep), `PEPProblemType`, `PEPSetProblemType()`, `PEP_GENERAL`, `PEP_HERMITIAN`, `PEP_HYPERBOLIC`, `PEPLINEAR`
127: M*/
129: /*E
130: PEPWhich - Determines which part of the spectrum is requested.
132: Values:
133: + `PEP_LARGEST_MAGNITUDE` - largest $|\lambda|$
134: . `PEP_SMALLEST_MAGNITUDE` - smallest $|\lambda|$
135: . `PEP_LARGEST_REAL` - largest $\mathrm{Re}(\lambda)$
136: . `PEP_SMALLEST_REAL` - smallest $\mathrm{Re}(\lambda)$
137: . `PEP_LARGEST_IMAGINARY` - largest $\mathrm{Im}(\lambda)$
138: . `PEP_SMALLEST_IMAGINARY` - smallest $\mathrm{Im}(\lambda)$
139: . `PEP_TARGET_MAGNITUDE` - smallest $|\lambda-\tau|$
140: . `PEP_TARGET_REAL` - smallest $|\mathrm{Re}(\lambda-\tau)|$
141: . `PEP_TARGET_IMAGINARY` - smallest $|\mathrm{Im}(\lambda-\tau)|$
142: . `PEP_ALL` - all $\lambda\in[a,b]$ or $\lambda\in\Omega$
143: - `PEP_WHICH_USER` - user-defined sorting criterion
145: Notes:
146: If SLEPc is compiled for real scalars `PEP_LARGEST_IMAGINARY` and
147: `PEP_SMALLEST_IMAGINARY` use the absolute value of the imaginary part
148: for eigenvalue selection.
150: The target $\tau$ is a scalar value provided with `PEPSetTarget()`.
152: The case `PEP_ALL` needs an interval $[a,b]$ given with `PEPSetInterval()`
153: or a region $\Omega$ specified with an `RG` object.
155: Level: intermediate
157: .seealso: [](ch:pep), `PEPSetWhichEigenpairs()`, `PEPSetTarget()`, `PEPSetInterval()`
158: E*/
159: typedef enum { PEP_LARGEST_MAGNITUDE = 1,
160: PEP_SMALLEST_MAGNITUDE = 2,
161: PEP_LARGEST_REAL = 3,
162: PEP_SMALLEST_REAL = 4,
163: PEP_LARGEST_IMAGINARY = 5,
164: PEP_SMALLEST_IMAGINARY = 6,
165: PEP_TARGET_MAGNITUDE = 7,
166: PEP_TARGET_REAL = 8,
167: PEP_TARGET_IMAGINARY = 9,
168: PEP_ALL = 10,
169: PEP_WHICH_USER = 11 } PEPWhich;
171: /*E
172: PEPBasis - The type of polynomial basis used to represent the polynomial
173: eigenproblem.
175: Values:
176: + `PEP_BASIS_MONOMIAL` - monomial basis
177: . `PEP_BASIS_CHEBYSHEV1` - Chebyshev polynomials of the 1st kind
178: . `PEP_BASIS_CHEBYSHEV2` - Chebyshev polynomials of the 2nd kind
179: . `PEP_BASIS_LEGENDRE` - Legendre polynomials
180: . `PEP_BASIS_LAGUERRE` - Laguerre polynomials
181: - `PEP_BASIS_HERMITE` - Hermite polynomials
183: Notes:
184: The default is to work with the monomial basis to represent the polynomial,
185: i.e., $1, x, x^2, \dots, x^d$. For large degree $d$, numerical
186: difficulties may arise. In that case, a different basis is recommended.
187: The user is responsible for providing the coefficient matrices to
188: `PEPSetOperators()` represented in the selected basis.
190: Level: intermediate
192: .seealso: [](ch:pep), `PEPSetBasis()`
193: E*/
194: typedef enum { PEP_BASIS_MONOMIAL,
195: PEP_BASIS_CHEBYSHEV1,
196: PEP_BASIS_CHEBYSHEV2,
197: PEP_BASIS_LEGENDRE,
198: PEP_BASIS_LAGUERRE,
199: PEP_BASIS_HERMITE } PEPBasis;
200: SLEPC_EXTERN const char *PEPBasisTypes[];
202: /*E
203: PEPScale - The scaling strategy.
205: Values:
206: + `PEP_SCALE_NONE` - no scaling
207: . `PEP_SCALE_SCALAR` - multiply by a scalar value
208: . `PEP_SCALE_DIAGONAL` - multiply by two diagonal matrices
209: - `PEP_SCALE_BOTH` - both scalar and diagonal scaling
211: Note:
212: See section [](#sec:scaling) for a discussion of the different scaling strategies.
214: Level: intermediate
216: .seealso: [](ch:pep), [](#sec:scaling), `PEPSetScale()`
217: E*/
218: typedef enum { PEP_SCALE_NONE,
219: PEP_SCALE_SCALAR,
220: PEP_SCALE_DIAGONAL,
221: PEP_SCALE_BOTH } PEPScale;
222: SLEPC_EXTERN const char *PEPScaleTypes[];
224: /*E
225: PEPRefine - The type of Newton iterative refinement.
227: Values:
228: + `PEP_REFINE_NONE` - no refinement
229: . `PEP_REFINE_SIMPLE` - refinement of each converged eigenpair individually
230: - `PEP_REFINE_MULTIPLE` - refinement of the invariant pair as a whole
232: Note:
233: See section [](#sec:refine) for a discussion of the different refinement strategies.
235: Level: intermediate
237: .seealso: [](ch:pep), [](#sec:refine), `PEPSetRefine()`
238: E*/
239: typedef enum { PEP_REFINE_NONE,
240: PEP_REFINE_SIMPLE,
241: PEP_REFINE_MULTIPLE } PEPRefine;
242: SLEPC_EXTERN const char *PEPRefineTypes[];
244: /*E
245: PEPRefineScheme - The scheme used for solving linear systems during iterative refinement.
247: Values:
248: + `PEP_REFINE_SCHEME_SCHUR` - use the Schur complement
249: . `PEP_REFINE_SCHEME_MBE` - use the mixed block elimination (MBE) scheme
250: - `PEP_REFINE_SCHEME_EXPLICIT` - build the full matrix explicitly
252: Note:
253: Iterative refinement may be very costly, due to the expensive linear
254: solves. These linear systems have a particular structure that can be
255: exploited in different ways, as described in {cite:p}`Cam16b`. See
256: `PEPSetRefine()` for additional details.
258: Level: intermediate
260: .seealso: [](ch:pep), [](#sec:refine), `PEPSetRefine()`
261: E*/
262: typedef enum { PEP_REFINE_SCHEME_SCHUR = 1,
263: PEP_REFINE_SCHEME_MBE = 2,
264: PEP_REFINE_SCHEME_EXPLICIT = 3 } PEPRefineScheme;
265: SLEPC_EXTERN const char *PEPRefineSchemes[];
267: /*E
268: PEPExtract - The eigenvector extraction strategy.
270: Values:
271: + `PEP_EXTRACT_NONE` - trivial extraction
272: . `PEP_EXTRACT_NORM` - extraction based on the norm
273: . `PEP_EXTRACT_RESIDUAL` - extraction based on the residual
274: - `PEP_EXTRACT_STRUCTURED` - extraction using a linear combination of all the blocks
276: Note:
277: This is relevant for solvers based on linearization. Once the solver has
278: converged, the polynomial eigenvectors can be extracted from the
279: eigenvectors of the linearized problem in different ways. See the
280: discussion in section [](#sec:pepextr).
282: Level: intermediate
284: .seealso: [](ch:pep), [](#sec:pepextr), `PEPSetExtract()`
285: E*/
286: typedef enum { PEP_EXTRACT_NONE = 1,
287: PEP_EXTRACT_NORM = 2,
288: PEP_EXTRACT_RESIDUAL = 3,
289: PEP_EXTRACT_STRUCTURED = 4 } PEPExtract;
290: SLEPC_EXTERN const char *PEPExtractTypes[];
292: /*MC
293: PEP_EXTRACT_NONE - Trivial eigenvector extraction.
295: Note:
296: Given the eigenvector of the linearization, $y$, the eigenvector of the
297: polynomial eigenproblem $x$ is taken from the first block of $y$.
299: Level: intermediate
301: .seealso: [](ch:pep), [](#sec:pepextr), `PEPExtract`, `PEPSetExtract()`, `PEP_EXTRACT_NORM`, `PEP_EXTRACT_RESIDUAL`, `PEP_EXTRACT_STRUCTURED`
302: M*/
304: /*MC
305: PEP_EXTRACT_NORM - Eigenvector extraction based on the norm.
307: Notes:
308: Given the eigenvector of the linearization, $y$, the eigenvector of the
309: polynomial eigenproblem $x$ is obtained from the $i$th block for which
310: $|\phi_i(\lambda)|$ is maximum, where $\phi_i$ is the $i$th element
311: of the polynomial basis. In the case of the monomial basis, it will
312: select the first or last block depending on $|\lambda|\geq 1$ or $|\lambda|<1$.
314: This is the default extraction.
316: Level: intermediate
318: .seealso: [](ch:pep), [](#sec:pepextr), `PEPExtract`, `PEPSetExtract()`, `PEP_EXTRACT_NONE`, `PEP_EXTRACT_RESIDUAL`, `PEP_EXTRACT_STRUCTURED`
319: M*/
321: /*MC
322: PEP_EXTRACT_RESIDUAL - Eigenvector extraction based on the residual.
324: Note:
325: Given the eigenvector of the linearization, $y$, the eigenvector of the
326: polynomial eigenproblem $x$ is the block of $y$ that minimizes the
327: residual norm.
329: Level: intermediate
331: .seealso: [](ch:pep), [](#sec:pepextr), `PEPExtract`, `PEPSetExtract()`, `PEP_EXTRACT_NONE`, `PEP_EXTRACT_NORM`, `PEP_EXTRACT_STRUCTURED`
332: M*/
334: /*MC
335: PEP_EXTRACT_STRUCTURED - Eigenvector extraction using a linear combination of
336: all the blocks.
338: Note:
339: Given the eigenvector of the linearization, $y$, the eigenvector of the
340: polynomial eigenproblem $x$ is obtained as a linear combination of all
341: the blocks, such that it minimizes a certain norm.
343: Level: intermediate
345: .seealso: [](ch:pep), [](#sec:pepextr), `PEPExtract`, `PEPSetExtract()`, `PEP_EXTRACT_NONE`, `PEP_EXTRACT_NORM`, `PEP_EXTRACT_RESIDUAL`
346: M*/
348: /*E
349: PEPErrorType - The error type used to assess the accuracy of computed solutions.
351: Values:
352: + `PEP_ERROR_ABSOLUTE` - compute error bound as $\|r\|$
353: . `PEP_ERROR_RELATIVE` - compute error bound as $\|r\|/|\lambda|$
354: - `PEP_ERROR_BACKWARD` - compute error bound as $\|r\|/(\sum_j\|A_j\||\lambda_i|^j)$
356: Level: intermediate
358: .seealso: [](ch:pep), `PEPComputeError()`
359: E*/
360: typedef enum { PEP_ERROR_ABSOLUTE,
361: PEP_ERROR_RELATIVE,
362: PEP_ERROR_BACKWARD } PEPErrorType;
363: SLEPC_EXTERN const char *PEPErrorTypes[];
365: /*E
366: PEPConv - The convergence criterion to be used by the solver.
368: Values:
369: + `PEP_CONV_ABS` - absolute convergence criterion, $\|r\|$
370: . `PEP_CONV_REL` - convergence criterion relative to eigenvalue, $\|r\|/|\lambda|$
371: . `PEP_CONV_NORM` - convergence criterion relative to matrix norms, $\|r\|/(\sum_j\|A_j\||\lambda|^j)$
372: - `PEP_CONV_USER` - convergence dictated by user-provided function
374: Level: intermediate
376: .seealso: [](ch:pep), `PEPSetConvergenceTest()`, `PEPSetConvergenceTestFunction()`
377: E*/
378: typedef enum { PEP_CONV_ABS,
379: PEP_CONV_REL,
380: PEP_CONV_NORM,
381: PEP_CONV_USER } PEPConv;
383: /*E
384: PEPStop - The stopping test to decide the termination of the outer loop
385: of the eigensolver.
387: Values:
388: + `PEP_STOP_BASIC` - default stopping test
389: - `PEP_STOP_USER` - user-provided stopping test
391: Level: advanced
393: .seealso: [](ch:pep), `PEPSetStoppingTest()`, `PEPSetStoppingTestFunction()`
394: E*/
395: typedef enum { PEP_STOP_BASIC,
396: PEP_STOP_USER } PEPStop;
398: /*MC
399: PEP_STOP_BASIC - The default stopping test.
401: Note:
402: By default, the termination of the outer loop is decided by calling
403: `PEPStoppingBasic()`, which will stop if all requested eigenvalues are converged,
404: or if the maximum number of iterations has been reached.
406: Level: advanced
408: .seealso: [](ch:pep), `PEPStop`, `PEPSetStoppingTest()`, `PEPStoppingBasic()`
409: M*/
411: /*MC
412: PEP_STOP_USER - The user-provided stopping test.
414: Note:
415: Customized stopping test using the user-provided function given with
416: `PEPSetStoppingTestFunction()`.
418: Level: advanced
420: .seealso: [](ch:pep), `PEPStop`, `PEPSetStoppingTest()`, `PEPSetStoppingTestFunction()`
421: M*/
423: /*E
424: PEPConvergedReason - Reason a polynomial eigensolver was determined to have converged
425: or diverged.
427: Values:
428: + `PEP_CONVERGED_TOL` - converged up to tolerance
429: . `PEP_CONVERGED_USER` - converged due to a user-defined condition
430: . `PEP_DIVERGED_ITS` - exceeded the maximum number of allowed iterations
431: . `PEP_DIVERGED_BREAKDOWN` - generic breakdown in method
432: . `PEP_DIVERGED_SYMMETRY_LOST` - pseudo-Lanczos was not able to keep symmetry
433: - `PEP_CONVERGED_ITERATING` - the solver is still running
435: Level: intermediate
437: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPSetTolerances()`
438: E*/
439: typedef enum {/* converged */
440: PEP_CONVERGED_TOL = 1,
441: PEP_CONVERGED_USER = 2,
442: /* diverged */
443: PEP_DIVERGED_ITS = -1,
444: PEP_DIVERGED_BREAKDOWN = -2,
445: PEP_DIVERGED_SYMMETRY_LOST = -3,
446: PEP_CONVERGED_ITERATING = 0} PEPConvergedReason;
447: SLEPC_EXTERN const char *const*PEPConvergedReasons;
449: /*MC
450: PEP_CONVERGED_TOL - The computed error estimates, based on residual norms,
451: for all requested eigenvalues are below the tolerance.
453: Level: intermediate
455: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPConvergedReason`
456: M*/
458: /*MC
459: PEP_CONVERGED_USER - The solver was declared converged due to a user-defined condition.
461: Note:
462: This happens only when a user-defined stopping test has been set with
463: `PEPSetStoppingTestFunction()`.
465: Level: intermediate
467: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPConvergedReason`, `PEPSetStoppingTestFunction()`
468: M*/
470: /*MC
471: PEP_DIVERGED_ITS - Exceeded the maximum number of allowed iterations
472: before the convergence criterion was satisfied.
474: Level: intermediate
476: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPConvergedReason`
477: M*/
479: /*MC
480: PEP_DIVERGED_BREAKDOWN - A breakdown in the solver was detected so the
481: method could not continue.
483: Level: intermediate
485: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPConvergedReason`
486: M*/
488: /*MC
489: PEP_DIVERGED_SYMMETRY_LOST - The selected solver uses a pseudo-Lanczos recurrence,
490: which is numerically unstable, and a symmetry test revealed that instability
491: had appeared so the solver could not continue.
493: Level: intermediate
495: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPConvergedReason`
496: M*/
498: /*MC
499: PEP_CONVERGED_ITERATING - This value is returned if `PEPGetConvergedReason()` is called
500: while `PEPSolve()` is still running.
502: Level: intermediate
504: .seealso: [](ch:pep), `PEPSolve()`, `PEPGetConvergedReason()`, `PEPConvergedReason`
505: M*/
507: SLEPC_EXTERN PetscErrorCode PEPCreate(MPI_Comm,PEP*);
508: SLEPC_EXTERN PetscErrorCode PEPDestroy(PEP*);
509: SLEPC_EXTERN PetscErrorCode PEPReset(PEP);
510: SLEPC_EXTERN PetscErrorCode PEPSetType(PEP,PEPType);
511: SLEPC_EXTERN PetscErrorCode PEPGetType(PEP,PEPType*);
512: SLEPC_EXTERN PetscErrorCode PEPSetProblemType(PEP,PEPProblemType);
513: SLEPC_EXTERN PetscErrorCode PEPGetProblemType(PEP,PEPProblemType*);
514: SLEPC_EXTERN PetscErrorCode PEPSetOperators(PEP,PetscInt,Mat[]);
515: SLEPC_EXTERN PetscErrorCode PEPGetOperators(PEP,PetscInt,Mat*);
516: SLEPC_EXTERN PetscErrorCode PEPGetNumMatrices(PEP,PetscInt*);
517: SLEPC_EXTERN PetscErrorCode PEPSetTarget(PEP,PetscScalar);
518: SLEPC_EXTERN PetscErrorCode PEPGetTarget(PEP,PetscScalar*);
519: SLEPC_EXTERN PetscErrorCode PEPSetInterval(PEP,PetscReal,PetscReal);
520: SLEPC_EXTERN PetscErrorCode PEPGetInterval(PEP,PetscReal*,PetscReal*);
521: SLEPC_EXTERN PetscErrorCode PEPSetFromOptions(PEP);
522: SLEPC_EXTERN PetscErrorCode PEPSetDSType(PEP);
523: SLEPC_EXTERN PetscErrorCode PEPSetUp(PEP);
524: SLEPC_EXTERN PetscErrorCode PEPSolve(PEP);
525: SLEPC_EXTERN PetscErrorCode PEPView(PEP,PetscViewer);
526: SLEPC_EXTERN PetscErrorCode PEPViewFromOptions(PEP,PetscObject,const char[]);
527: SLEPC_EXTERN PetscErrorCode PEPErrorView(PEP,PEPErrorType,PetscViewer);
528: PETSC_DEPRECATED_FUNCTION(3, 6, 0, "PEPErrorView()", ) static inline PetscErrorCode PEPPrintSolution(PEP pep,PetscViewer v) {return PEPErrorView(pep,PEP_ERROR_BACKWARD,v);}
529: SLEPC_EXTERN PetscErrorCode PEPErrorViewFromOptions(PEP);
530: SLEPC_EXTERN PetscErrorCode PEPConvergedReasonView(PEP,PetscViewer);
531: SLEPC_EXTERN PetscErrorCode PEPConvergedReasonViewFromOptions(PEP);
532: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "PEPConvergedReasonView()", ) static inline PetscErrorCode PEPReasonView(PEP pep,PetscViewer v) {return PEPConvergedReasonView(pep,v);}
533: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "PEPConvergedReasonViewFromOptions()", ) static inline PetscErrorCode PEPReasonViewFromOptions(PEP pep) {return PEPConvergedReasonViewFromOptions(pep);}
534: SLEPC_EXTERN PetscErrorCode PEPValuesView(PEP,PetscViewer);
535: SLEPC_EXTERN PetscErrorCode PEPValuesViewFromOptions(PEP);
536: SLEPC_EXTERN PetscErrorCode PEPVectorsView(PEP,PetscViewer);
537: SLEPC_EXTERN PetscErrorCode PEPVectorsViewFromOptions(PEP);
538: SLEPC_EXTERN PetscErrorCode PEPSetBV(PEP,BV);
539: SLEPC_EXTERN PetscErrorCode PEPGetBV(PEP,BV*);
540: SLEPC_EXTERN PetscErrorCode PEPSetRG(PEP,RG);
541: SLEPC_EXTERN PetscErrorCode PEPGetRG(PEP,RG*);
542: SLEPC_EXTERN PetscErrorCode PEPSetDS(PEP,DS);
543: SLEPC_EXTERN PetscErrorCode PEPGetDS(PEP,DS*);
544: SLEPC_EXTERN PetscErrorCode PEPSetST(PEP,ST);
545: SLEPC_EXTERN PetscErrorCode PEPGetST(PEP,ST*);
546: SLEPC_EXTERN PetscErrorCode PEPRefineGetKSP(PEP,KSP*);
548: SLEPC_EXTERN PetscErrorCode PEPSetTolerances(PEP,PetscReal,PetscInt);
549: SLEPC_EXTERN PetscErrorCode PEPGetTolerances(PEP,PetscReal*,PetscInt*);
550: SLEPC_EXTERN PetscErrorCode PEPGetConvergedReason(PEP,PEPConvergedReason*);
552: SLEPC_EXTERN PetscErrorCode PEPSetDimensions(PEP,PetscInt,PetscInt,PetscInt);
553: SLEPC_EXTERN PetscErrorCode PEPGetDimensions(PEP,PetscInt*,PetscInt*,PetscInt*);
554: SLEPC_EXTERN PetscErrorCode PEPSetScale(PEP,PEPScale,PetscReal,Vec,Vec,PetscInt,PetscReal);
555: SLEPC_EXTERN PetscErrorCode PEPGetScale(PEP,PEPScale*,PetscReal*,Vec*,Vec*,PetscInt*,PetscReal*);
556: SLEPC_EXTERN PetscErrorCode PEPSetRefine(PEP,PEPRefine,PetscInt,PetscReal,PetscInt,PEPRefineScheme);
557: SLEPC_EXTERN PetscErrorCode PEPGetRefine(PEP,PEPRefine*,PetscInt*,PetscReal*,PetscInt*,PEPRefineScheme*);
558: SLEPC_EXTERN PetscErrorCode PEPSetExtract(PEP,PEPExtract);
559: SLEPC_EXTERN PetscErrorCode PEPGetExtract(PEP,PEPExtract*);
560: SLEPC_EXTERN PetscErrorCode PEPSetBasis(PEP,PEPBasis);
561: SLEPC_EXTERN PetscErrorCode PEPGetBasis(PEP,PEPBasis*);
563: SLEPC_EXTERN PetscErrorCode PEPGetConverged(PEP,PetscInt*);
564: SLEPC_EXTERN PetscErrorCode PEPGetEigenpair(PEP,PetscInt,PetscScalar*,PetscScalar*,Vec,Vec);
565: SLEPC_EXTERN PetscErrorCode PEPComputeError(PEP,PetscInt,PEPErrorType,PetscReal*);
566: 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);}
567: 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);}
568: SLEPC_EXTERN PetscErrorCode PEPGetErrorEstimate(PEP,PetscInt,PetscReal*);
569: SLEPC_EXTERN PetscErrorCode PEPGetIterationNumber(PEP,PetscInt*);
571: SLEPC_EXTERN PetscErrorCode PEPSetInitialSpace(PEP,PetscInt,Vec[]);
572: SLEPC_EXTERN PetscErrorCode PEPSetWhichEigenpairs(PEP,PEPWhich);
573: SLEPC_EXTERN PetscErrorCode PEPGetWhichEigenpairs(PEP,PEPWhich*);
575: SLEPC_EXTERN PetscErrorCode PEPSetTrackAll(PEP,PetscBool);
576: SLEPC_EXTERN PetscErrorCode PEPGetTrackAll(PEP,PetscBool*);
578: /*S
579: PEPMonitorFn - A function prototype for functions provided to `PEPMonitorSet()`.
581: Calling Sequence:
582: + pep - the polynomial eigensolver context
583: . its - iteration number
584: . nconv - number of converged eigenpairs
585: . eigr - real part of the eigenvalues
586: . eigi - imaginary part of the eigenvalues
587: . errest - relative error estimates for each eigenpair
588: . nest - number of error estimates
589: - ctx - optional monitoring context, as provided with `PEPMonitorSet()`
591: Level: intermediate
593: .seealso: [](ch:pep), `PEPMonitorSet()`
594: S*/
595: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorFn(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,void *ctx);
597: /*S
598: PEPMonitorRegisterFn - A function prototype for functions provided to `PEPMonitorRegister()`.
600: Calling Sequence:
601: + pep - the polynomial eigensolver context
602: . its - iteration number
603: . nconv - number of converged eigenpairs
604: . eigr - real part of the eigenvalues
605: . eigi - imaginary part of the eigenvalues
606: . errest - relative error estimates for each eigenpair
607: . nest - number of error estimates
608: - ctx - `PetscViewerAndFormat` object
610: Level: advanced
612: Note:
613: This is a `PEPMonitorFn` specialized for a context of `PetscViewerAndFormat`.
615: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorRegister()`, `PEPMonitorFn`, `PEPMonitorRegisterCreateFn`, `PEPMonitorRegisterDestroyFn`
616: S*/
617: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorRegisterFn(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *ctx);
619: /*S
620: PEPMonitorRegisterCreateFn - A function prototype for functions that do the
621: creation when provided to `PEPMonitorRegister()`.
623: Calling Sequence:
624: + viewer - the viewer to be used with the `PEPMonitorRegisterFn`
625: . format - the format of the viewer
626: . ctx - a context for the monitor
627: - result - a `PetscViewerAndFormat` object
629: Level: advanced
631: .seealso: [](ch:pep), `PEPMonitorRegisterFn`, `PEPMonitorSet()`, `PEPMonitorRegister()`, `PEPMonitorFn`, `PEPMonitorRegisterDestroyFn`
632: S*/
633: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorRegisterCreateFn(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **result);
635: /*S
636: PEPMonitorRegisterDestroyFn - A function prototype for functions that do the after
637: use destruction when provided to `PEPMonitorRegister()`.
639: Calling Sequence:
640: . vf - a `PetscViewerAndFormat` object to be destroyed, including any context
642: Level: advanced
644: .seealso: [](ch:pep), `PEPMonitorRegisterFn`, `PEPMonitorSet()`, `PEPMonitorRegister()`, `PEPMonitorFn`, `PEPMonitorRegisterCreateFn`
645: S*/
646: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPMonitorRegisterDestroyFn(PetscViewerAndFormat **result);
648: SLEPC_EXTERN PetscErrorCode PEPMonitor(PEP,PetscInt,PetscInt,PetscScalar[],PetscScalar[],PetscReal[],PetscInt);
649: SLEPC_EXTERN PetscErrorCode PEPMonitorSet(PEP,PEPMonitorFn,void*,PetscCtxDestroyFn*);
650: SLEPC_EXTERN PetscErrorCode PEPMonitorCancel(PEP);
651: SLEPC_EXTERN PetscErrorCode PEPGetMonitorContext(PEP,void*);
653: SLEPC_EXTERN PetscErrorCode PEPMonitorSetFromOptions(PEP,const char[],const char[],void*,PetscBool);
654: SLEPC_EXTERN PEPMonitorRegisterFn PEPMonitorFirst;
655: SLEPC_EXTERN PEPMonitorRegisterFn PEPMonitorFirstDrawLG;
656: SLEPC_EXTERN PEPMonitorRegisterCreateFn PEPMonitorFirstDrawLGCreate;
657: SLEPC_EXTERN PEPMonitorRegisterFn PEPMonitorAll;
658: SLEPC_EXTERN PEPMonitorRegisterFn PEPMonitorAllDrawLG;
659: SLEPC_EXTERN PEPMonitorRegisterCreateFn PEPMonitorAllDrawLGCreate;
660: SLEPC_EXTERN PEPMonitorRegisterFn PEPMonitorConverged;
661: SLEPC_EXTERN PEPMonitorRegisterCreateFn PEPMonitorConvergedCreate;
662: SLEPC_EXTERN PEPMonitorRegisterFn PEPMonitorConvergedDrawLG;
663: SLEPC_EXTERN PEPMonitorRegisterCreateFn PEPMonitorConvergedDrawLGCreate;
664: SLEPC_EXTERN PEPMonitorRegisterDestroyFn PEPMonitorConvergedDestroy;
666: SLEPC_EXTERN PetscErrorCode PEPSetOptionsPrefix(PEP,const char[]);
667: SLEPC_EXTERN PetscErrorCode PEPAppendOptionsPrefix(PEP,const char[]);
668: SLEPC_EXTERN PetscErrorCode PEPGetOptionsPrefix(PEP,const char*[]);
670: SLEPC_EXTERN PetscFunctionList PEPList;
671: SLEPC_EXTERN PetscFunctionList PEPMonitorList;
672: SLEPC_EXTERN PetscFunctionList PEPMonitorCreateList;
673: SLEPC_EXTERN PetscFunctionList PEPMonitorDestroyList;
674: SLEPC_EXTERN PetscErrorCode PEPRegister(const char[],PetscErrorCode(*)(PEP));
675: SLEPC_EXTERN PetscErrorCode PEPMonitorRegister(const char[],PetscViewerType,PetscViewerFormat,PEPMonitorRegisterFn*,PEPMonitorRegisterCreateFn*,PEPMonitorRegisterDestroyFn*);
677: SLEPC_EXTERN PetscErrorCode PEPSetWorkVecs(PEP,PetscInt);
678: SLEPC_EXTERN PetscErrorCode PEPAllocateSolution(PEP,PetscInt);
680: /*S
681: PEPConvergenceTestFn - A prototype of a `PEP` convergence test function that
682: would be passed to `PEPSetConvergenceTestFunction()`.
684: Calling Sequence:
685: + pep - the polynomial eigensolver context
686: . eigr - real part of the eigenvalue
687: . eigi - imaginary part of the eigenvalue
688: . res - residual norm associated to the eigenpair
689: . errest - [output] computed error estimate
690: - ctx - optional convergence context, as set by `PEPSetConvergenceTestFunction()`
692: Level: advanced
694: .seealso: [](ch:pep), `PEPSetConvergenceTest()`, `PEPSetConvergenceTestFunction()`
695: S*/
696: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPConvergenceTestFn(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx);
698: SLEPC_EXTERN PetscErrorCode PEPSetConvergenceTest(PEP,PEPConv);
699: SLEPC_EXTERN PetscErrorCode PEPGetConvergenceTest(PEP,PEPConv*);
700: SLEPC_EXTERN PEPConvergenceTestFn PEPConvergedAbsolute;
701: SLEPC_EXTERN PEPConvergenceTestFn PEPConvergedRelative;
702: SLEPC_EXTERN PEPConvergenceTestFn PEPConvergedNorm;
703: SLEPC_EXTERN PetscErrorCode PEPSetConvergenceTestFunction(PEP,PEPConvergenceTestFn*,void*,PetscCtxDestroyFn*);
705: /*S
706: PEPStoppingTestFn - A prototype of a `PEP` stopping test function that would be
707: passed to `PEPSetStoppingTestFunction()`.
709: Calling Sequence:
710: + pep - the polynomial eigensolver context
711: . its - current number of iterations
712: . max_it - maximum number of iterations
713: . nconv - number of currently converged eigenpairs
714: . nev - number of requested eigenpairs
715: . reason - [output] result of the stopping test
716: - ctx - optional stopping context, as set by `PEPSetStoppingTestFunction()`
718: Note:
719: A positive value of `reason` indicates that the iteration has finished successfully
720: (converged), and a negative value indicates an error condition (diverged). If
721: the iteration needs to be continued, `reason` must be set to `PEP_CONVERGED_ITERATING`
722: (zero).
724: Level: advanced
726: .seealso: [](ch:pep), `PEPSetStoppingTestFunction()`
727: S*/
728: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode PEPStoppingTestFn(PEP pep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,PEPConvergedReason *reason,void *ctx);
730: SLEPC_EXTERN PetscErrorCode PEPSetStoppingTest(PEP,PEPStop);
731: SLEPC_EXTERN PetscErrorCode PEPGetStoppingTest(PEP,PEPStop*);
732: SLEPC_EXTERN PEPStoppingTestFn PEPStoppingBasic;
733: SLEPC_EXTERN PetscErrorCode PEPSetStoppingTestFunction(PEP,PEPStoppingTestFn*,void*,PetscCtxDestroyFn*);
735: SLEPC_EXTERN PetscErrorCode PEPSetEigenvalueComparison(PEP,SlepcEigenvalueComparisonFn*,void*);
737: /* --------- options specific to particular eigensolvers -------- */
739: SLEPC_EXTERN PetscErrorCode PEPLinearSetLinearization(PEP,PetscReal,PetscReal);
740: SLEPC_EXTERN PetscErrorCode PEPLinearGetLinearization(PEP,PetscReal*,PetscReal*);
741: SLEPC_EXTERN PetscErrorCode PEPLinearSetExplicitMatrix(PEP,PetscBool);
742: SLEPC_EXTERN PetscErrorCode PEPLinearGetExplicitMatrix(PEP,PetscBool*);
743: SLEPC_EXTERN PetscErrorCode PEPLinearSetEPS(PEP,EPS);
744: SLEPC_EXTERN PetscErrorCode PEPLinearGetEPS(PEP,EPS*);
745: 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);}
746: PETSC_DEPRECATED_FUNCTION(3, 10, 0, "PEPLinearGetLinearization()", ) static inline PetscErrorCode PEPLinearGetCompanionForm(PEP pep,PetscInt *cform) {(void)pep; if (cform) *cform=1; return PETSC_SUCCESS;}
748: SLEPC_EXTERN PetscErrorCode PEPQArnoldiSetRestart(PEP,PetscReal);
749: SLEPC_EXTERN PetscErrorCode PEPQArnoldiGetRestart(PEP,PetscReal*);
750: SLEPC_EXTERN PetscErrorCode PEPQArnoldiSetLocking(PEP,PetscBool);
751: SLEPC_EXTERN PetscErrorCode PEPQArnoldiGetLocking(PEP,PetscBool*);
753: SLEPC_EXTERN PetscErrorCode PEPTOARSetRestart(PEP,PetscReal);
754: SLEPC_EXTERN PetscErrorCode PEPTOARGetRestart(PEP,PetscReal*);
755: SLEPC_EXTERN PetscErrorCode PEPTOARSetLocking(PEP,PetscBool);
756: SLEPC_EXTERN PetscErrorCode PEPTOARGetLocking(PEP,PetscBool*);
758: SLEPC_EXTERN PetscErrorCode PEPSTOARSetLinearization(PEP,PetscReal,PetscReal);
759: SLEPC_EXTERN PetscErrorCode PEPSTOARGetLinearization(PEP,PetscReal*,PetscReal*);
760: SLEPC_EXTERN PetscErrorCode PEPSTOARSetLocking(PEP,PetscBool);
761: SLEPC_EXTERN PetscErrorCode PEPSTOARGetLocking(PEP,PetscBool*);
762: SLEPC_EXTERN PetscErrorCode PEPSTOARSetDetectZeros(PEP,PetscBool);
763: SLEPC_EXTERN PetscErrorCode PEPSTOARGetDetectZeros(PEP,PetscBool*);
764: SLEPC_EXTERN PetscErrorCode PEPSTOARGetInertias(PEP,PetscInt*,PetscReal*[],PetscInt*[]);
765: SLEPC_EXTERN PetscErrorCode PEPSTOARSetDimensions(PEP,PetscInt,PetscInt,PetscInt);
766: SLEPC_EXTERN PetscErrorCode PEPSTOARGetDimensions(PEP,PetscInt*,PetscInt*,PetscInt*);
767: SLEPC_EXTERN PetscErrorCode PEPSTOARSetCheckEigenvalueType(PEP,PetscBool);
768: SLEPC_EXTERN PetscErrorCode PEPSTOARGetCheckEigenvalueType(PEP,PetscBool*);
769: SLEPC_EXTERN PetscErrorCode PEPCheckDefiniteQEP(PEP,PetscReal*,PetscReal*,PetscInt*,PetscInt*);
771: /*E
772: PEPJDProjection - The type of projection to be used in the Jacobi-Davidson solver.
774: Values:
775: + `PEP_JD_PROJECTION_HARMONIC` - oblique projection
776: - `PEP_JD_PROJECTION_ORTHOGONAL` - orthogonal projection
778: Level: advanced
780: .seealso: [](ch:pep), `PEPJDSetProjection()`
781: E*/
782: typedef enum { PEP_JD_PROJECTION_HARMONIC,
783: PEP_JD_PROJECTION_ORTHOGONAL } PEPJDProjection;
784: SLEPC_EXTERN const char *PEPJDProjectionTypes[];
786: SLEPC_EXTERN PetscErrorCode PEPJDSetRestart(PEP,PetscReal);
787: SLEPC_EXTERN PetscErrorCode PEPJDGetRestart(PEP,PetscReal*);
788: SLEPC_EXTERN PetscErrorCode PEPJDSetFix(PEP,PetscReal);
789: SLEPC_EXTERN PetscErrorCode PEPJDGetFix(PEP,PetscReal*);
790: SLEPC_EXTERN PetscErrorCode PEPJDSetReusePreconditioner(PEP,PetscBool);
791: SLEPC_EXTERN PetscErrorCode PEPJDGetReusePreconditioner(PEP,PetscBool*);
792: SLEPC_EXTERN PetscErrorCode PEPJDSetMinimalityIndex(PEP,PetscInt);
793: SLEPC_EXTERN PetscErrorCode PEPJDGetMinimalityIndex(PEP,PetscInt*);
794: SLEPC_EXTERN PetscErrorCode PEPJDSetProjection(PEP,PEPJDProjection);
795: SLEPC_EXTERN PetscErrorCode PEPJDGetProjection(PEP,PEPJDProjection*);
797: /*E
798: PEPCISSExtraction - The extraction technique used in the CISS solver.
800: Values:
801: + `PEP_CISS_EXTRACTION_RITZ` - Rayleigh-Ritz extraction
802: . `PEP_CISS_EXTRACTION_HANKEL` - block Hankel method
803: - `PEP_CISS_EXTRACTION_CAA` - communication-avoiding Arnoldi method
805: Level: advanced
807: .seealso: [](ch:pep), `PEPCISSSetExtraction()`, `PEPCISSGetExtraction()`
808: E*/
809: typedef enum { PEP_CISS_EXTRACTION_RITZ,
810: PEP_CISS_EXTRACTION_HANKEL,
811: PEP_CISS_EXTRACTION_CAA } PEPCISSExtraction;
812: SLEPC_EXTERN const char *PEPCISSExtractions[];
814: #if defined(PETSC_USE_COMPLEX) || defined(PETSC_CLANG_STATIC_ANALYZER)
815: SLEPC_EXTERN PetscErrorCode PEPCISSSetExtraction(PEP,PEPCISSExtraction);
816: SLEPC_EXTERN PetscErrorCode PEPCISSGetExtraction(PEP,PEPCISSExtraction*);
817: SLEPC_EXTERN PetscErrorCode PEPCISSSetSizes(PEP,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscBool);
818: SLEPC_EXTERN PetscErrorCode PEPCISSGetSizes(PEP,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscBool*);
819: SLEPC_EXTERN PetscErrorCode PEPCISSSetThreshold(PEP,PetscReal,PetscReal);
820: SLEPC_EXTERN PetscErrorCode PEPCISSGetThreshold(PEP,PetscReal*,PetscReal*);
821: SLEPC_EXTERN PetscErrorCode PEPCISSSetRefinement(PEP,PetscInt,PetscInt);
822: SLEPC_EXTERN PetscErrorCode PEPCISSGetRefinement(PEP,PetscInt*,PetscInt*);
823: SLEPC_EXTERN PetscErrorCode PEPCISSGetKSPs(PEP,PetscInt*,KSP*[]);
824: #else
825: #define SlepcPEPCISSUnavailable(pep) do { \
826: PetscFunctionBegin; \
827: SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"%s() not available with real scalars",PETSC_FUNCTION_NAME); \
828: } while (0)
829: static inline PetscErrorCode PEPCISSSetExtraction(PEP pep,PETSC_UNUSED PEPCISSExtraction ex) {SlepcPEPCISSUnavailable(pep);}
830: static inline PetscErrorCode PEPCISSGetExtraction(PEP pep,PETSC_UNUSED PEPCISSExtraction *ex) {SlepcPEPCISSUnavailable(pep);}
831: 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);}
832: 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);}
833: static inline PetscErrorCode PEPCISSSetThreshold(PEP pep,PETSC_UNUSED PetscReal delta,PETSC_UNUSED PetscReal spur) {SlepcPEPCISSUnavailable(pep);}
834: static inline PetscErrorCode PEPCISSGetThreshold(PEP pep,PETSC_UNUSED PetscReal *delta,PETSC_UNUSED PetscReal *spur) {SlepcPEPCISSUnavailable(pep);}
835: static inline PetscErrorCode PEPCISSSetRefinement(PEP pep,PETSC_UNUSED PetscInt inner,PETSC_UNUSED PetscInt blsize) {SlepcPEPCISSUnavailable(pep);}
836: static inline PetscErrorCode PEPCISSGetRefinement(PEP pep,PETSC_UNUSED PetscInt *inner,PETSC_UNUSED PetscInt *blsize) {SlepcPEPCISSUnavailable(pep);}
837: static inline PetscErrorCode PEPCISSGetKSPs(PEP pep,PETSC_UNUSED PetscInt *nsolve,PETSC_UNUSED KSP *ksp[]) {SlepcPEPCISSUnavailable(pep);}
838: #undef SlepcPEPCISSUnavailable
839: #endif