Actual source code: slepcnep.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 nonlinear eigenvalue solvers
12: */
14: #pragma once
16: #include <slepceps.h>
17: #include <slepcpep.h>
18: #include <slepcfn.h>
20: /* SUBMANSEC = NEP */
22: SLEPC_EXTERN PetscErrorCode NEPInitializePackage(void);
23: SLEPC_EXTERN PetscErrorCode NEPFinalizePackage(void);
25: /*S
26: NEP - SLEPc object that manages all solvers for nonlinear eigenvalue problems.
28: Level: beginner
30: .seealso: [](ch:nep), `NEPCreate()`
31: S*/
32: typedef struct _p_NEP* NEP;
34: /*J
35: NEPType - String with the name of a nonlinear eigensolver.
37: Level: beginner
39: .seealso: [](ch:nep), `NEPSetType()`, `NEP`
40: J*/
41: typedef const char *NEPType;
42: #define NEPRII "rii"
43: #define NEPSLP "slp"
44: #define NEPNARNOLDI "narnoldi"
45: #define NEPNLEIGS "nleigs"
46: #define NEPCISS "ciss"
47: #define NEPINTERPOL "interpol"
49: /* Logging support */
50: SLEPC_EXTERN PetscClassId NEP_CLASSID;
52: /*E
53: NEPProblemType - Determines the type of the nonlinear eigenproblem.
55: Values:
56: + `NEP_GENERAL` - no particular structure
57: - `NEP_RATIONAL` - problem defined in split form with all $f_i$ rational
59: Note:
60: Currently, the `NEP_RATIONAL` case is only used in `NEPNLEIGS` to
61: determine the singularities automatically.
63: Level: intermediate
65: .seealso: [](ch:nep), `NEPSetProblemType()`, `NEPGetProblemType()`
66: E*/
67: typedef enum { NEP_GENERAL = 1,
68: NEP_RATIONAL = 2
69: } NEPProblemType;
71: /*E
72: NEPWhich - Determines which part of the spectrum is requested.
74: Values:
75: + `NEP_LARGEST_MAGNITUDE` - largest $|\lambda|$
76: . `NEP_SMALLEST_MAGNITUDE` - smallest $|\lambda|$
77: . `NEP_LARGEST_REAL` - largest $\mathrm{Re}(\lambda)$
78: . `NEP_SMALLEST_REAL` - smallest $\mathrm{Re}(\lambda)$
79: . `NEP_LARGEST_IMAGINARY` - largest $\mathrm{Im}(\lambda)$
80: . `NEP_SMALLEST_IMAGINARY` - smallest $\mathrm{Im}(\lambda)$
81: . `NEP_TARGET_MAGNITUDE` - smallest $|\lambda-\tau|$
82: . `NEP_TARGET_REAL` - smallest $|\mathrm{Re}(\lambda-\tau)|$
83: . `NEP_TARGET_IMAGINARY` - smallest $|\mathrm{Im}(\lambda-\tau)|$
84: . `NEP_ALL` - all $\lambda\in\Omega$
85: - `NEP_WHICH_USER` - user-defined sorting criterion
87: Notes:
88: The target $\tau$ is a scalar value provided with `NEPSetTarget()`.
90: The case `NEP_ALL` needs a region $\Omega$ specified with an `RG` object.
92: Level: intermediate
94: .seealso: [](ch:nep), `NEPSetWhichEigenpairs()`, `NEPSetTarget()`
95: E*/
96: typedef enum { NEP_LARGEST_MAGNITUDE = 1,
97: NEP_SMALLEST_MAGNITUDE = 2,
98: NEP_LARGEST_REAL = 3,
99: NEP_SMALLEST_REAL = 4,
100: NEP_LARGEST_IMAGINARY = 5,
101: NEP_SMALLEST_IMAGINARY = 6,
102: NEP_TARGET_MAGNITUDE = 7,
103: NEP_TARGET_REAL = 8,
104: NEP_TARGET_IMAGINARY = 9,
105: NEP_ALL = 10,
106: NEP_WHICH_USER = 11 } NEPWhich;
108: /*E
109: NEPErrorType - The error type used to assess the accuracy of computed solutions.
111: Values:
112: + `NEP_ERROR_ABSOLUTE` - compute error bound as $\|r\|$
113: . `NEP_ERROR_RELATIVE` - compute error bound as $\|r\|/|\lambda|$
114: - `NEP_ERROR_BACKWARD` - compute error bound as $\|r\|/(\sum_i|f_i(\lambda)|\|A_i\|)$
116: Level: intermediate
118: .seealso: [](ch:nep), `NEPComputeError()`
119: E*/
120: typedef enum { NEP_ERROR_ABSOLUTE,
121: NEP_ERROR_RELATIVE,
122: NEP_ERROR_BACKWARD } NEPErrorType;
123: SLEPC_EXTERN const char *NEPErrorTypes[];
125: /*E
126: NEPRefine - The type of Newton iterative refinement.
128: Values:
129: + `NEP_REFINE_NONE` - no refinement
130: . `NEP_REFINE_SIMPLE` - refinement of each converged eigenpair individually
131: - `NEP_REFINE_MULTIPLE` - refinement of the invariant pair as a whole
133: Note:
134: See section [](#sec:refine) for a discussion of the different refinement strategies.
136: Level: intermediate
138: .seealso: [](ch:nep), [](#sec:refine), `NEPSetRefine()`
139: E*/
140: typedef enum { NEP_REFINE_NONE,
141: NEP_REFINE_SIMPLE,
142: NEP_REFINE_MULTIPLE } NEPRefine;
143: SLEPC_EXTERN const char *NEPRefineTypes[];
145: /*E
146: NEPRefineScheme - The scheme used for solving linear systems during iterative refinement.
148: Values:
149: + `NEP_REFINE_SCHEME_SCHUR` - use the Schur complement
150: . `NEP_REFINE_SCHEME_MBE` - use the mixed block elimination (MBE) scheme
151: - `NEP_REFINE_SCHEME_EXPLICIT` - build the full matrix explicitly
153: Note:
154: Iterative refinement may be very costly, due to the expensive linear
155: solves. These linear systems have a particular structure that can be
156: exploited in different ways, as described in {cite:p}`Cam16b`. See
157: `NEPSetRefine()` for additional details.
159: Level: intermediate
161: .seealso: [](ch:nep), [](#sec:refine), `NEPSetRefine()`
162: E*/
163: typedef enum { NEP_REFINE_SCHEME_SCHUR = 1,
164: NEP_REFINE_SCHEME_MBE = 2,
165: NEP_REFINE_SCHEME_EXPLICIT = 3 } NEPRefineScheme;
166: SLEPC_EXTERN const char *NEPRefineSchemes[];
168: /*E
169: NEPConv - The convergence criterion to be used by the solver.
171: Values:
172: + `NEP_CONV_ABS` - absolute convergence criterion, $\|r\|$
173: . `NEP_CONV_REL` - convergence criterion relative to eigenvalue, $\|r\|/|\lambda|$
174: . `NEP_CONV_NORM` - convergence criterion relative to matrix norms, $\|r\|/(\sum_j|f_j(\lambda)|\|A_j\|)$
175: - `NEP_CONV_USER` - convergence dictated by user-provided function
177: Level: intermediate
179: .seealso: [](ch:nep), `NEPSetConvergenceTest()`, `NEPSetConvergenceTestFunction()`
180: E*/
181: typedef enum { NEP_CONV_ABS,
182: NEP_CONV_REL,
183: NEP_CONV_NORM,
184: NEP_CONV_USER } NEPConv;
186: /*E
187: NEPStop - The stopping test to decide the termination of the outer loop
188: of the eigensolver.
190: Values:
191: + `NEP_STOP_BASIC` - default stopping test
192: - `NEP_STOP_USER` - user-provided stopping test
194: Level: advanced
196: .seealso: [](ch:nep), `NEPSetStoppingTest()`, `NEPSetStoppingTestFunction()`
197: E*/
198: typedef enum { NEP_STOP_BASIC,
199: NEP_STOP_USER } NEPStop;
201: /*MC
202: NEP_STOP_BASIC - The default stopping test.
204: Note:
205: By default, the termination of the outer loop is decided by calling
206: `NEPStoppingBasic()`, which will stop if all requested eigenvalues are converged,
207: or if the maximum number of iterations has been reached.
209: Level: advanced
211: .seealso: [](ch:nep), `NEPStop`, `NEPSetStoppingTest()`, `NEPStoppingBasic()`
212: M*/
214: /*MC
215: NEP_STOP_USER - The user-provided stopping test.
217: Note:
218: Customized stopping test using the user-provided function given with
219: `NEPSetStoppingTestFunction()`.
221: Level: advanced
223: .seealso: [](ch:nep), `NEPStop`, `NEPSetStoppingTest()`, `NEPSetStoppingTestFunction()`
224: M*/
226: /*E
227: NEPConvergedReason - Reason a nonlinear eigensolver was determined to have converged
228: or diverged.
230: Values:
231: + `NEP_CONVERGED_TOL` - converged up to tolerance
232: . `NEP_CONVERGED_USER` - converged due to a user-defined condition
233: . `NEP_DIVERGED_ITS` - exceeded the maximum number of allowed iterations
234: . `NEP_DIVERGED_BREAKDOWN` - generic breakdown in method
235: . `NEP_DIVERGED_LINEAR_SOLVE` - inner linear solve failed
236: . `NEP_DIVERGED_SUBSPACE_EXHAUSTED` - run out of space for the basis in an unrestarted solver
237: - `NEP_CONVERGED_ITERATING` - the solver is still running
239: Level: intermediate
241: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPSetTolerances()`
242: E*/
243: typedef enum {/* converged */
244: NEP_CONVERGED_TOL = 1,
245: NEP_CONVERGED_USER = 2,
246: /* diverged */
247: NEP_DIVERGED_ITS = -1,
248: NEP_DIVERGED_BREAKDOWN = -2,
249: /* unused = -3 */
250: NEP_DIVERGED_LINEAR_SOLVE = -4,
251: NEP_DIVERGED_SUBSPACE_EXHAUSTED = -5,
252: NEP_CONVERGED_ITERATING = 0} NEPConvergedReason;
253: SLEPC_EXTERN const char *const*NEPConvergedReasons;
255: /*MC
256: NEP_CONVERGED_TOL - The computed error estimates, based on residual norms,
257: for all requested eigenvalues are below the tolerance.
259: Level: intermediate
261: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPConvergedReason`
262: M*/
264: /*MC
265: NEP_CONVERGED_USER - The solver was declared converged due to a user-defined condition.
267: Note:
268: This happens only when a user-defined stopping test has been set with
269: `NEPSetStoppingTestFunction()`.
271: Level: intermediate
273: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPConvergedReason`, `NEPSetStoppingTestFunction()`
274: M*/
276: /*MC
277: NEP_DIVERGED_ITS - Exceeded the maximum number of allowed iterations
278: before the convergence criterion was satisfied.
280: Level: intermediate
282: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPConvergedReason`
283: M*/
285: /*MC
286: NEP_DIVERGED_BREAKDOWN - A breakdown in the solver was detected so the
287: method could not continue.
289: Level: intermediate
291: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPConvergedReason`
292: M*/
294: /*MC
295: NEP_DIVERGED_LINEAR_SOLVE - The inner linear solve failed so the nonlinear
296: eigensolver could not continue.
298: Level: intermediate
300: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPConvergedReason`
301: M*/
303: /*MC
304: NEP_DIVERGED_SUBSPACE_EXHAUSTED - The solver has run out of space for the
305: basis in the case of an unrestarted method.
307: Level: intermediate
309: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPConvergedReason`
310: M*/
312: /*MC
313: NEP_CONVERGED_ITERATING - This value is returned if `NEPGetConvergedReason()` is called
314: while `NEPSolve()` is still running.
316: Level: intermediate
318: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPConvergedReason`
319: M*/
321: SLEPC_EXTERN PetscErrorCode NEPCreate(MPI_Comm,NEP*);
322: SLEPC_EXTERN PetscErrorCode NEPDestroy(NEP*);
323: SLEPC_EXTERN PetscErrorCode NEPReset(NEP);
324: SLEPC_EXTERN PetscErrorCode NEPSetType(NEP,NEPType);
325: SLEPC_EXTERN PetscErrorCode NEPGetType(NEP,NEPType*);
326: SLEPC_EXTERN PetscErrorCode NEPSetProblemType(NEP,NEPProblemType);
327: SLEPC_EXTERN PetscErrorCode NEPGetProblemType(NEP,NEPProblemType*);
328: SLEPC_EXTERN PetscErrorCode NEPSetTarget(NEP,PetscScalar);
329: SLEPC_EXTERN PetscErrorCode NEPGetTarget(NEP,PetscScalar*);
330: SLEPC_EXTERN PetscErrorCode NEPSetFromOptions(NEP);
331: SLEPC_EXTERN PetscErrorCode NEPSetDSType(NEP);
332: SLEPC_EXTERN PetscErrorCode NEPSetUp(NEP);
333: SLEPC_EXTERN PetscErrorCode NEPSolve(NEP);
334: SLEPC_EXTERN PetscErrorCode NEPView(NEP,PetscViewer);
335: SLEPC_EXTERN PetscErrorCode NEPViewFromOptions(NEP,PetscObject,const char[]);
336: SLEPC_EXTERN PetscErrorCode NEPErrorView(NEP,NEPErrorType,PetscViewer);
337: SLEPC_EXTERN PetscErrorCode NEPErrorViewFromOptions(NEP);
338: SLEPC_EXTERN PetscErrorCode NEPConvergedReasonView(NEP,PetscViewer);
339: SLEPC_EXTERN PetscErrorCode NEPConvergedReasonViewFromOptions(NEP);
340: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "NEPConvergedReasonView()", ) static inline PetscErrorCode NEPReasonView(NEP nep,PetscViewer v) {return NEPConvergedReasonView(nep,v);}
341: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "NEPConvergedReasonViewFromOptions()", ) static inline PetscErrorCode NEPReasonViewFromOptions(NEP nep) {return NEPConvergedReasonViewFromOptions(nep);}
342: SLEPC_EXTERN PetscErrorCode NEPValuesView(NEP,PetscViewer);
343: SLEPC_EXTERN PetscErrorCode NEPValuesViewFromOptions(NEP);
344: SLEPC_EXTERN PetscErrorCode NEPVectorsView(NEP,PetscViewer);
345: SLEPC_EXTERN PetscErrorCode NEPVectorsViewFromOptions(NEP);
347: /*S
348: NEPFunctionFn - A prototype of a `NEP` function evaluation function that
349: would be passed to `NEPSetFunction()`.
351: Calling Sequence:
352: + nep - the nonlinear eigensolver context
353: . lambda - the scalar argument where $T(\cdot)$ must be evaluated
354: . T - matrix that will contain $T(\lambda)$
355: . P - [optional] different matrix to build the preconditioner
356: - ctx - [optional] user-defined context for private data for the
357: function evaluation routine (may be `NULL`)
359: Level: beginner
361: .seealso: [](ch:nep), `NEPSetFunction()`
362: S*/
363: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPFunctionFn(NEP nep,PetscScalar lambda,Mat T,Mat P,void *ctx);
365: /*S
366: NEPJacobianFn - A prototype of a `NEP` Jacobian evaluation function that
367: would be passed to `NEPSetJacobian()`.
369: Calling Sequence:
370: + nep - the nonlinear eigensolver context
371: . lambda - the scalar argument where $T'(\cdot)$ must be evaluated
372: . J - matrix that will contain $T'(\lambda)$
373: - ctx - [optional] user-defined context for private data for the
374: Jacobian evaluation routine (may be `NULL`)
376: Level: beginner
378: .seealso: [](ch:nep), `NEPSetJacobian()`
379: S*/
380: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPJacobianFn(NEP nep,PetscScalar lambda,Mat J,void *ctx);
382: SLEPC_EXTERN PetscErrorCode NEPSetFunction(NEP,Mat,Mat,NEPFunctionFn*,void*);
383: SLEPC_EXTERN PetscErrorCode NEPGetFunction(NEP,Mat*,Mat*,NEPFunctionFn**,void**);
384: SLEPC_EXTERN PetscErrorCode NEPSetJacobian(NEP,Mat,NEPJacobianFn*,void*);
385: SLEPC_EXTERN PetscErrorCode NEPGetJacobian(NEP,Mat*,NEPJacobianFn**,void**);
386: PETSC_DEPRECATED_FUNCTION(3, 12, 0, "NEPSetFunction() and NEPSetJacobian()", ) static inline PetscErrorCode NEPSetDerivatives(NEP nep,Mat A,PetscErrorCode (*fun)(NEP,PetscScalar,PetscInt,Mat,void*),void *ctx) {(void)A;(void)fun;(void)ctx;SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"Not implemented in this version");}
387: PETSC_DEPRECATED_FUNCTION(3, 12, 0, "NEPGetFunction() and NEPGetJacobian()", ) static inline PetscErrorCode NEPGetDerivatives(NEP nep,Mat *A,PetscErrorCode (**fun)(NEP,PetscScalar,PetscInt,Mat,void*),void **ctx) {(void)A;(void)fun;(void)ctx;SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"Not implemented in this version");}
388: SLEPC_EXTERN PetscErrorCode NEPSetSplitOperator(NEP,PetscInt,Mat[],FN[],MatStructure);
389: SLEPC_EXTERN PetscErrorCode NEPGetSplitOperatorTerm(NEP,PetscInt,Mat*,FN*);
390: SLEPC_EXTERN PetscErrorCode NEPGetSplitOperatorInfo(NEP,PetscInt*,MatStructure*);
391: SLEPC_EXTERN PetscErrorCode NEPSetSplitPreconditioner(NEP,PetscInt,Mat[],MatStructure);
392: SLEPC_EXTERN PetscErrorCode NEPGetSplitPreconditionerTerm(NEP,PetscInt,Mat*);
393: SLEPC_EXTERN PetscErrorCode NEPGetSplitPreconditionerInfo(NEP,PetscInt*,MatStructure*);
395: SLEPC_EXTERN PetscErrorCode NEPSetBV(NEP,BV);
396: SLEPC_EXTERN PetscErrorCode NEPGetBV(NEP,BV*);
397: SLEPC_EXTERN PetscErrorCode NEPSetRG(NEP,RG);
398: SLEPC_EXTERN PetscErrorCode NEPGetRG(NEP,RG*);
399: SLEPC_EXTERN PetscErrorCode NEPSetDS(NEP,DS);
400: SLEPC_EXTERN PetscErrorCode NEPGetDS(NEP,DS*);
401: SLEPC_EXTERN PetscErrorCode NEPRefineGetKSP(NEP,KSP*);
402: SLEPC_EXTERN PetscErrorCode NEPSetTolerances(NEP,PetscReal,PetscInt);
403: SLEPC_EXTERN PetscErrorCode NEPGetTolerances(NEP,PetscReal*,PetscInt*);
404: SLEPC_EXTERN PetscErrorCode NEPSetDimensions(NEP,PetscInt,PetscInt,PetscInt);
405: SLEPC_EXTERN PetscErrorCode NEPGetDimensions(NEP,PetscInt*,PetscInt*,PetscInt*);
406: SLEPC_EXTERN PetscErrorCode NEPSetRefine(NEP,NEPRefine,PetscInt,PetscReal,PetscInt,NEPRefineScheme);
407: SLEPC_EXTERN PetscErrorCode NEPGetRefine(NEP,NEPRefine*,PetscInt*,PetscReal*,PetscInt*,NEPRefineScheme*);
409: SLEPC_EXTERN PetscErrorCode NEPGetConverged(NEP,PetscInt*);
410: SLEPC_EXTERN PetscErrorCode NEPGetEigenpair(NEP,PetscInt,PetscScalar*,PetscScalar*,Vec,Vec);
411: SLEPC_EXTERN PetscErrorCode NEPGetLeftEigenvector(NEP,PetscInt,Vec,Vec);
413: SLEPC_EXTERN PetscErrorCode NEPComputeError(NEP,PetscInt,NEPErrorType,PetscReal*);
414: PETSC_DEPRECATED_FUNCTION(3, 6, 0, "NEPComputeError()", ) static inline PetscErrorCode NEPComputeRelativeError(NEP nep,PetscInt i,PetscReal *r) {return NEPComputeError(nep,i,NEP_ERROR_RELATIVE,r);}
415: PETSC_DEPRECATED_FUNCTION(3, 6, 0, "NEPComputeError() with NEP_ERROR_ABSOLUTE", ) static inline PetscErrorCode NEPComputeResidualNorm(NEP nep,PetscInt i,PetscReal *r) {return NEPComputeError(nep,i,NEP_ERROR_ABSOLUTE,r);}
416: SLEPC_EXTERN PetscErrorCode NEPGetErrorEstimate(NEP,PetscInt,PetscReal*);
418: SLEPC_EXTERN PetscErrorCode NEPComputeFunction(NEP,PetscScalar,Mat,Mat);
419: SLEPC_EXTERN PetscErrorCode NEPComputeJacobian(NEP,PetscScalar,Mat);
420: SLEPC_EXTERN PetscErrorCode NEPApplyFunction(NEP,PetscScalar,Vec,Vec,Vec,Mat,Mat);
421: SLEPC_EXTERN PetscErrorCode NEPApplyAdjoint(NEP,PetscScalar,Vec,Vec,Vec,Mat,Mat);
422: SLEPC_EXTERN PetscErrorCode NEPApplyJacobian(NEP,PetscScalar,Vec,Vec,Vec,Mat);
423: SLEPC_EXTERN PetscErrorCode NEPProjectOperator(NEP,PetscInt,PetscInt);
424: SLEPC_EXTERN PetscErrorCode NEPGetIterationNumber(NEP,PetscInt*);
426: SLEPC_EXTERN PetscErrorCode NEPSetInitialSpace(NEP,PetscInt,Vec[]);
427: SLEPC_EXTERN PetscErrorCode NEPSetWhichEigenpairs(NEP,NEPWhich);
428: SLEPC_EXTERN PetscErrorCode NEPGetWhichEigenpairs(NEP,NEPWhich*);
429: SLEPC_EXTERN PetscErrorCode NEPSetTwoSided(NEP,PetscBool);
430: SLEPC_EXTERN PetscErrorCode NEPGetTwoSided(NEP,PetscBool*);
432: SLEPC_EXTERN PetscErrorCode NEPApplyResolvent(NEP,RG,PetscScalar,Vec,Vec);
434: SLEPC_EXTERN PetscErrorCode NEPSetTrackAll(NEP,PetscBool);
435: SLEPC_EXTERN PetscErrorCode NEPGetTrackAll(NEP,PetscBool*);
437: SLEPC_EXTERN PetscErrorCode NEPGetConvergedReason(NEP,NEPConvergedReason*);
439: /*S
440: NEPMonitorFn - A function prototype for functions provided to `NEPMonitorSet()`.
442: Calling Sequence:
443: + nep - the nonlinear eigensolver context
444: . its - iteration number
445: . nconv - number of converged eigenpairs
446: . eigr - real part of the eigenvalues
447: . eigi - imaginary part of the eigenvalues
448: . errest - relative error estimates for each eigenpair
449: . nest - number of error estimates
450: - ctx - optional monitoring context, as provided with `NEPMonitorSet()`
452: Level: intermediate
454: .seealso: [](ch:nep), `NEPMonitorSet()`
455: S*/
456: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorFn(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,void *ctx);
458: /*S
459: NEPMonitorRegisterFn - A function prototype for functions provided to `NEPMonitorRegister()`.
461: Calling Sequence:
462: + nep - the nonlinear eigensolver context
463: . its - iteration number
464: . nconv - number of converged eigenpairs
465: . eigr - real part of the eigenvalues
466: . eigi - imaginary part of the eigenvalues
467: . errest - relative error estimates for each eigenpair
468: . nest - number of error estimates
469: - ctx - `PetscViewerAndFormat` object
471: Level: advanced
473: Note:
474: This is a `NEPMonitorFn` specialized for a context of `PetscViewerAndFormat`.
476: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorRegister()`, `NEPMonitorFn`, `NEPMonitorRegisterCreateFn`, `NEPMonitorRegisterDestroyFn`
477: S*/
478: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorRegisterFn(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *ctx);
480: /*S
481: NEPMonitorRegisterCreateFn - A function prototype for functions that do the
482: creation when provided to `NEPMonitorRegister()`.
484: Calling Sequence:
485: + viewer - the viewer to be used with the `NEPMonitorRegisterFn`
486: . format - the format of the viewer
487: . ctx - a context for the monitor
488: - result - a `PetscViewerAndFormat` object
490: Level: advanced
492: .seealso: [](ch:nep), `NEPMonitorRegisterFn`, `NEPMonitorSet()`, `NEPMonitorRegister()`, `NEPMonitorFn`, `NEPMonitorRegisterDestroyFn`
493: S*/
494: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorRegisterCreateFn(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **result);
496: /*S
497: NEPMonitorRegisterDestroyFn - A function prototype for functions that do the after
498: use destruction when provided to `NEPMonitorRegister()`.
500: Calling Sequence:
501: . vf - a `PetscViewerAndFormat` object to be destroyed, including any context
503: Level: advanced
505: .seealso: [](ch:nep), `NEPMonitorRegisterFn`, `NEPMonitorSet()`, `NEPMonitorRegister()`, `NEPMonitorFn`, `NEPMonitorRegisterCreateFn`
506: S*/
507: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorRegisterDestroyFn(PetscViewerAndFormat **result);
509: SLEPC_EXTERN PetscErrorCode NEPMonitor(NEP,PetscInt,PetscInt,PetscScalar[],PetscScalar[],PetscReal[],PetscInt);
510: SLEPC_EXTERN PetscErrorCode NEPMonitorSet(NEP,NEPMonitorFn,void*,PetscCtxDestroyFn*);
511: SLEPC_EXTERN PetscErrorCode NEPMonitorCancel(NEP);
512: SLEPC_EXTERN PetscErrorCode NEPGetMonitorContext(NEP,void*);
514: SLEPC_EXTERN PetscErrorCode NEPMonitorSetFromOptions(NEP,const char[],const char[],void*,PetscBool);
515: SLEPC_EXTERN NEPMonitorRegisterFn NEPMonitorFirst;
516: SLEPC_EXTERN NEPMonitorRegisterFn NEPMonitorFirstDrawLG;
517: SLEPC_EXTERN NEPMonitorRegisterCreateFn NEPMonitorFirstDrawLGCreate;
518: SLEPC_EXTERN NEPMonitorRegisterFn NEPMonitorAll;
519: SLEPC_EXTERN NEPMonitorRegisterFn NEPMonitorAllDrawLG;
520: SLEPC_EXTERN NEPMonitorRegisterCreateFn NEPMonitorAllDrawLGCreate;
521: SLEPC_EXTERN NEPMonitorRegisterFn NEPMonitorConverged;
522: SLEPC_EXTERN NEPMonitorRegisterCreateFn NEPMonitorConvergedCreate;
523: SLEPC_EXTERN NEPMonitorRegisterFn NEPMonitorConvergedDrawLG;
524: SLEPC_EXTERN NEPMonitorRegisterCreateFn NEPMonitorConvergedDrawLGCreate;
525: SLEPC_EXTERN NEPMonitorRegisterDestroyFn NEPMonitorConvergedDestroy;
527: SLEPC_EXTERN PetscErrorCode NEPSetOptionsPrefix(NEP,const char[]);
528: SLEPC_EXTERN PetscErrorCode NEPAppendOptionsPrefix(NEP,const char[]);
529: SLEPC_EXTERN PetscErrorCode NEPGetOptionsPrefix(NEP,const char*[]);
531: SLEPC_EXTERN PetscFunctionList NEPList;
532: SLEPC_EXTERN PetscFunctionList NEPMonitorList;
533: SLEPC_EXTERN PetscFunctionList NEPMonitorCreateList;
534: SLEPC_EXTERN PetscFunctionList NEPMonitorDestroyList;
535: SLEPC_EXTERN PetscErrorCode NEPRegister(const char[],PetscErrorCode(*)(NEP));
536: SLEPC_EXTERN PetscErrorCode NEPMonitorRegister(const char[],PetscViewerType,PetscViewerFormat,NEPMonitorRegisterFn*,NEPMonitorRegisterCreateFn*,NEPMonitorRegisterDestroyFn*);
538: SLEPC_EXTERN PetscErrorCode NEPSetWorkVecs(NEP,PetscInt);
539: SLEPC_EXTERN PetscErrorCode NEPAllocateSolution(NEP,PetscInt);
541: /*S
542: NEPConvergenceTestFn - A prototype of a `NEP` convergence test function that
543: would be passed to `NEPSetConvergenceTestFunction()`.
545: Calling Sequence:
546: + nep - the nonlinear eigensolver context
547: . eigr - real part of the eigenvalue
548: . eigi - imaginary part of the eigenvalue
549: . res - residual norm associated to the eigenpair
550: . errest - [output] computed error estimate
551: - ctx - optional convergence context, as set by `NEPSetConvergenceTestFunction()`
553: Level: advanced
555: .seealso: [](ch:nep), `NEPSetConvergenceTest()`, `NEPSetConvergenceTestFunction()`
556: S*/
557: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPConvergenceTestFn(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx);
559: SLEPC_EXTERN PetscErrorCode NEPSetConvergenceTest(NEP,NEPConv);
560: SLEPC_EXTERN PetscErrorCode NEPGetConvergenceTest(NEP,NEPConv*);
561: SLEPC_EXTERN NEPConvergenceTestFn NEPConvergedAbsolute;
562: SLEPC_EXTERN NEPConvergenceTestFn NEPConvergedRelative;
563: SLEPC_EXTERN NEPConvergenceTestFn NEPConvergedNorm;
564: SLEPC_EXTERN PetscErrorCode NEPSetConvergenceTestFunction(NEP,NEPConvergenceTestFn*,void*,PetscCtxDestroyFn*);
566: /*S
567: NEPStoppingTestFn - A prototype of a `NEP` stopping test function that would
568: be passed to `NEPSetStoppingTestFunction()`.
570: Calling Sequence:
571: + nep - the nonlinear eigensolver context
572: . its - current number of iterations
573: . max_it - maximum number of iterations
574: . nconv - number of currently converged eigenpairs
575: . nev - number of requested eigenpairs
576: . reason - [output] result of the stopping test
577: - ctx - optional stopping context, as set by `NEPSetStoppingTestFunction()`
579: Note:
580: A positive value of `reason` indicates that the iteration has finished successfully
581: (converged), and a negative value indicates an error condition (diverged). If
582: the iteration needs to be continued, `reason` must be set to `NEP_CONVERGED_ITERATING`
583: (zero).
585: Level: advanced
587: .seealso: [](ch:nep), `NEPSetStoppingTestFunction()`
588: S*/
589: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPStoppingTestFn(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx);
591: SLEPC_EXTERN PetscErrorCode NEPSetStoppingTest(NEP,NEPStop);
592: SLEPC_EXTERN PetscErrorCode NEPGetStoppingTest(NEP,NEPStop*);
593: SLEPC_EXTERN NEPStoppingTestFn NEPStoppingBasic;
594: SLEPC_EXTERN PetscErrorCode NEPSetStoppingTestFunction(NEP,NEPStoppingTestFn*,void*,PetscCtxDestroyFn*);
596: SLEPC_EXTERN PetscErrorCode NEPSetEigenvalueComparison(NEP,SlepcEigenvalueComparisonFn*,void*);
598: /* --------- options specific to particular eigensolvers -------- */
600: SLEPC_EXTERN PetscErrorCode NEPRIISetMaximumIterations(NEP,PetscInt);
601: SLEPC_EXTERN PetscErrorCode NEPRIIGetMaximumIterations(NEP,PetscInt*);
602: SLEPC_EXTERN PetscErrorCode NEPRIISetLagPreconditioner(NEP,PetscInt);
603: SLEPC_EXTERN PetscErrorCode NEPRIIGetLagPreconditioner(NEP,PetscInt*);
604: SLEPC_EXTERN PetscErrorCode NEPRIISetConstCorrectionTol(NEP,PetscBool);
605: SLEPC_EXTERN PetscErrorCode NEPRIIGetConstCorrectionTol(NEP,PetscBool*);
606: SLEPC_EXTERN PetscErrorCode NEPRIISetHermitian(NEP,PetscBool);
607: SLEPC_EXTERN PetscErrorCode NEPRIIGetHermitian(NEP,PetscBool*);
608: SLEPC_EXTERN PetscErrorCode NEPRIISetDeflationThreshold(NEP,PetscReal);
609: SLEPC_EXTERN PetscErrorCode NEPRIIGetDeflationThreshold(NEP,PetscReal*);
610: SLEPC_EXTERN PetscErrorCode NEPRIISetKSP(NEP,KSP);
611: SLEPC_EXTERN PetscErrorCode NEPRIIGetKSP(NEP,KSP*);
613: SLEPC_EXTERN PetscErrorCode NEPSLPSetDeflationThreshold(NEP,PetscReal);
614: SLEPC_EXTERN PetscErrorCode NEPSLPGetDeflationThreshold(NEP,PetscReal*);
615: SLEPC_EXTERN PetscErrorCode NEPSLPSetEPS(NEP,EPS);
616: SLEPC_EXTERN PetscErrorCode NEPSLPGetEPS(NEP,EPS*);
617: SLEPC_EXTERN PetscErrorCode NEPSLPSetEPSLeft(NEP,EPS);
618: SLEPC_EXTERN PetscErrorCode NEPSLPGetEPSLeft(NEP,EPS*);
619: SLEPC_EXTERN PetscErrorCode NEPSLPSetKSP(NEP,KSP);
620: SLEPC_EXTERN PetscErrorCode NEPSLPGetKSP(NEP,KSP*);
622: SLEPC_EXTERN PetscErrorCode NEPNArnoldiSetKSP(NEP,KSP);
623: SLEPC_EXTERN PetscErrorCode NEPNArnoldiGetKSP(NEP,KSP*);
624: SLEPC_EXTERN PetscErrorCode NEPNArnoldiSetLagPreconditioner(NEP,PetscInt);
625: SLEPC_EXTERN PetscErrorCode NEPNArnoldiGetLagPreconditioner(NEP,PetscInt*);
627: /*E
628: NEPCISSExtraction - The extraction technique used in the CISS solver.
630: Values:
631: + `NEP_CISS_EXTRACTION_RITZ` - Rayleigh-Ritz extraction
632: . `NEP_CISS_EXTRACTION_HANKEL` - block Hankel method
633: - `NEP_CISS_EXTRACTION_CAA` - communication-avoiding Arnoldi method
635: Level: advanced
637: .seealso: [](ch:nep), `NEPCISSSetExtraction()`, `NEPCISSGetExtraction()`
638: E*/
639: typedef enum { NEP_CISS_EXTRACTION_RITZ,
640: NEP_CISS_EXTRACTION_HANKEL,
641: NEP_CISS_EXTRACTION_CAA } NEPCISSExtraction;
642: SLEPC_EXTERN const char *NEPCISSExtractions[];
644: #if defined(PETSC_USE_COMPLEX) || defined(PETSC_CLANG_STATIC_ANALYZER)
645: SLEPC_EXTERN PetscErrorCode NEPCISSSetExtraction(NEP,NEPCISSExtraction);
646: SLEPC_EXTERN PetscErrorCode NEPCISSGetExtraction(NEP,NEPCISSExtraction*);
647: SLEPC_EXTERN PetscErrorCode NEPCISSSetSizes(NEP,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscBool);
648: SLEPC_EXTERN PetscErrorCode NEPCISSGetSizes(NEP,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscBool*);
649: SLEPC_EXTERN PetscErrorCode NEPCISSSetThreshold(NEP,PetscReal,PetscReal);
650: SLEPC_EXTERN PetscErrorCode NEPCISSGetThreshold(NEP,PetscReal*,PetscReal*);
651: SLEPC_EXTERN PetscErrorCode NEPCISSSetRefinement(NEP,PetscInt,PetscInt);
652: SLEPC_EXTERN PetscErrorCode NEPCISSGetRefinement(NEP,PetscInt*,PetscInt*);
653: SLEPC_EXTERN PetscErrorCode NEPCISSGetKSPs(NEP,PetscInt*,KSP*[]);
654: #else
655: #define SlepcNEPCISSUnavailable(nep) do { \
656: PetscFunctionBegin; \
657: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"%s() not available with real scalars",PETSC_FUNCTION_NAME); \
658: } while (0)
659: static inline PetscErrorCode NEPCISSSetExtraction(NEP nep,PETSC_UNUSED NEPCISSExtraction ex) {SlepcNEPCISSUnavailable(nep);}
660: static inline PetscErrorCode NEPCISSGetExtraction(NEP nep,PETSC_UNUSED NEPCISSExtraction *ex) {SlepcNEPCISSUnavailable(nep);}
661: static inline PetscErrorCode NEPCISSSetSizes(NEP nep,PETSC_UNUSED PetscInt ip,PETSC_UNUSED PetscInt bs,PETSC_UNUSED PetscInt ms,PETSC_UNUSED PetscInt npart,PETSC_UNUSED PetscInt bsmax,PETSC_UNUSED PetscBool realmats) {SlepcNEPCISSUnavailable(nep);}
662: static inline PetscErrorCode NEPCISSGetSizes(NEP nep,PETSC_UNUSED PetscInt *ip,PETSC_UNUSED PetscInt *bs,PETSC_UNUSED PetscInt *ms,PETSC_UNUSED PetscInt *npart,PETSC_UNUSED PetscInt *bsmak,PETSC_UNUSED PetscBool *realmats) {SlepcNEPCISSUnavailable(nep);}
663: static inline PetscErrorCode NEPCISSSetThreshold(NEP nep,PETSC_UNUSED PetscReal delta,PETSC_UNUSED PetscReal spur) {SlepcNEPCISSUnavailable(nep);}
664: static inline PetscErrorCode NEPCISSGetThreshold(NEP nep,PETSC_UNUSED PetscReal *delta,PETSC_UNUSED PetscReal *spur) {SlepcNEPCISSUnavailable(nep);}
665: static inline PetscErrorCode NEPCISSSetRefinement(NEP nep,PETSC_UNUSED PetscInt inner,PETSC_UNUSED PetscInt blsize) {SlepcNEPCISSUnavailable(nep);}
666: static inline PetscErrorCode NEPCISSGetRefinement(NEP nep,PETSC_UNUSED PetscInt *inner,PETSC_UNUSED PetscInt *blsize) {SlepcNEPCISSUnavailable(nep);}
667: static inline PetscErrorCode NEPCISSGetKSPs(NEP nep,PETSC_UNUSED PetscInt *nsolve,PETSC_UNUSED KSP *ksp[]) {SlepcNEPCISSUnavailable(nep);}
668: #undef SlepcNEPCISSUnavailable
669: #endif
671: SLEPC_EXTERN PetscErrorCode NEPInterpolSetPEP(NEP,PEP);
672: SLEPC_EXTERN PetscErrorCode NEPInterpolGetPEP(NEP,PEP*);
673: SLEPC_EXTERN PetscErrorCode NEPInterpolSetInterpolation(NEP,PetscReal,PetscInt);
674: SLEPC_EXTERN PetscErrorCode NEPInterpolGetInterpolation(NEP,PetscReal*,PetscInt*);
676: /*S
677: NEPNLEIGSSingularitiesFn - A prototype of a function that would be passed
678: to `NEPNLEIGSSetSingularitiesFunction()`.
680: Calling Sequence:
681: + nep - the nonlinear eigensolver context
682: . maxnp - on input the number of requested points in the discretization (can be modified)
683: . xi - computed values of the discretization
684: - ctx - optional context, as set by `NEPNLEIGSSetSingularitiesFunction()`
686: Note:
687: The user-defined function can set a smaller value of `maxnp` if necessary.
688: It is wrong to return a larger value.
690: Level: intermediate
692: .seealso: [](ch:nep), `NEPNLEIGSSetSingularitiesFunction()`
693: S*/
694: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPNLEIGSSingularitiesFn(NEP nep,PetscInt *maxnp,PetscScalar *xi,void *ctx);
696: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetSingularitiesFunction(NEP,NEPNLEIGSSingularitiesFn*,void*);
697: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetSingularitiesFunction(NEP,NEPNLEIGSSingularitiesFn**,void**);
698: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetRestart(NEP,PetscReal);
699: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetRestart(NEP,PetscReal*);
700: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetLocking(NEP,PetscBool);
701: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetLocking(NEP,PetscBool*);
702: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetInterpolation(NEP,PetscReal,PetscInt);
703: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetInterpolation(NEP,PetscReal*,PetscInt*);
704: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetRKShifts(NEP,PetscInt,PetscScalar[]);
705: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetRKShifts(NEP,PetscInt*,PetscScalar*[]);
706: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetKSPs(NEP,PetscInt*,KSP*[]);
707: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetFullBasis(NEP,PetscBool);
708: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetFullBasis(NEP,PetscBool*);
709: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetEPS(NEP,EPS);
710: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetEPS(NEP,EPS*);