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:    Level: intermediate

 57: .seealso: [](ch:nep), `NEPSetProblemType()`, `NEPGetProblemType()`
 58: E*/
 59: typedef enum { NEP_GENERAL  = 1,
 60:                NEP_RATIONAL = 2     /* NEP defined in split form with all f_i rational */
 61:              } NEPProblemType;

 63: /*E
 64:    NEPWhich - Determines which part of the spectrum is requested

 66:    Level: intermediate

 68: .seealso: [](ch:nep), `NEPSetWhichEigenpairs()`, `NEPGetWhichEigenpairs()`
 69: E*/
 70: typedef enum { NEP_LARGEST_MAGNITUDE  = 1,
 71:                NEP_SMALLEST_MAGNITUDE = 2,
 72:                NEP_LARGEST_REAL       = 3,
 73:                NEP_SMALLEST_REAL      = 4,
 74:                NEP_LARGEST_IMAGINARY  = 5,
 75:                NEP_SMALLEST_IMAGINARY = 6,
 76:                NEP_TARGET_MAGNITUDE   = 7,
 77:                NEP_TARGET_REAL        = 8,
 78:                NEP_TARGET_IMAGINARY   = 9,
 79:                NEP_ALL                = 10,
 80:                NEP_WHICH_USER         = 11 } NEPWhich;

 82: /*E
 83:    NEPErrorType - The error type used to assess accuracy of computed solutions

 85:    Level: intermediate

 87: .seealso: [](ch:nep), `NEPComputeError()`
 88: E*/
 89: typedef enum { NEP_ERROR_ABSOLUTE,
 90:                NEP_ERROR_RELATIVE,
 91:                NEP_ERROR_BACKWARD } NEPErrorType;
 92: SLEPC_EXTERN const char *NEPErrorTypes[];

 94: /*E
 95:    NEPRefine - The refinement type

 97:    Level: intermediate

 99: .seealso: [](ch:nep), `NEPSetRefine()`
100: E*/
101: typedef enum { NEP_REFINE_NONE,
102:                NEP_REFINE_SIMPLE,
103:                NEP_REFINE_MULTIPLE } NEPRefine;
104: SLEPC_EXTERN const char *NEPRefineTypes[];

106: /*E
107:    NEPRefineScheme - The scheme used for solving linear systems during iterative refinement

109:    Level: intermediate

111: .seealso: [](ch:nep), `NEPSetRefine()`
112: E*/
113: typedef enum { NEP_REFINE_SCHEME_SCHUR    = 1,
114:                NEP_REFINE_SCHEME_MBE      = 2,
115:                NEP_REFINE_SCHEME_EXPLICIT = 3 } NEPRefineScheme;
116: SLEPC_EXTERN const char *NEPRefineSchemes[];

118: /*E
119:    NEPConv - Determines the convergence test

121:    Level: intermediate

123: .seealso: [](ch:nep), `NEPSetConvergenceTest()`, `NEPSetConvergenceTestFunction()`
124: E*/
125: typedef enum { NEP_CONV_ABS,
126:                NEP_CONV_REL,
127:                NEP_CONV_NORM,
128:                NEP_CONV_USER } NEPConv;

130: /*E
131:    NEPStop - Determines the stopping test

133:    Level: advanced

135: .seealso: [](ch:nep), `NEPSetStoppingTest()`, `NEPSetStoppingTestFunction()`
136: E*/
137: typedef enum { NEP_STOP_BASIC,
138:                NEP_STOP_USER } NEPStop;

140: /*E
141:    NEPConvergedReason - Reason a nonlinear eigensolver was determined to have converged
142:    or diverged.

144:    Values:
145: +  `NEP_CONVERGED_TOL`               - converged up to tolerance
146: .  `NEP_CONVERGED_USER`              - converged due to a user-defined condition
147: .  `NEP_DIVERGED_ITS`                - exceeded the maximum number of allowed iterations
148: .  `NEP_DIVERGED_BREAKDOWN`          - generic breakdown in method
149: .  `NEP_DIVERGED_LINEAR_SOLVE`       - inner linear solve failed
150: -  `NEP_DIVERGED_SUBSPACE_EXHAUSTED` - run out of space for the basis in an unrestarted solver

152:    Level: intermediate

154: .seealso: [](ch:nep), `NEPSolve()`, `NEPGetConvergedReason()`, `NEPSetTolerances()`
155: E*/
156: typedef enum {/* converged */
157:               NEP_CONVERGED_TOL                =  1,
158:               NEP_CONVERGED_USER               =  2,
159:               /* diverged */
160:               NEP_DIVERGED_ITS                 = -1,
161:               NEP_DIVERGED_BREAKDOWN           = -2,
162:                     /* unused                  = -3 */
163:               NEP_DIVERGED_LINEAR_SOLVE        = -4,
164:               NEP_DIVERGED_SUBSPACE_EXHAUSTED  = -5,
165:               NEP_CONVERGED_ITERATING          =  0} NEPConvergedReason;
166: SLEPC_EXTERN const char *const*NEPConvergedReasons;

168: SLEPC_EXTERN PetscErrorCode NEPCreate(MPI_Comm,NEP*);
169: SLEPC_EXTERN PetscErrorCode NEPDestroy(NEP*);
170: SLEPC_EXTERN PetscErrorCode NEPReset(NEP);
171: SLEPC_EXTERN PetscErrorCode NEPSetType(NEP,NEPType);
172: SLEPC_EXTERN PetscErrorCode NEPGetType(NEP,NEPType*);
173: SLEPC_EXTERN PetscErrorCode NEPSetProblemType(NEP,NEPProblemType);
174: SLEPC_EXTERN PetscErrorCode NEPGetProblemType(NEP,NEPProblemType*);
175: SLEPC_EXTERN PetscErrorCode NEPSetTarget(NEP,PetscScalar);
176: SLEPC_EXTERN PetscErrorCode NEPGetTarget(NEP,PetscScalar*);
177: SLEPC_EXTERN PetscErrorCode NEPSetFromOptions(NEP);
178: SLEPC_EXTERN PetscErrorCode NEPSetDSType(NEP);
179: SLEPC_EXTERN PetscErrorCode NEPSetUp(NEP);
180: SLEPC_EXTERN PetscErrorCode NEPSolve(NEP);
181: SLEPC_EXTERN PetscErrorCode NEPView(NEP,PetscViewer);
182: SLEPC_EXTERN PetscErrorCode NEPViewFromOptions(NEP,PetscObject,const char[]);
183: SLEPC_EXTERN PetscErrorCode NEPErrorView(NEP,NEPErrorType,PetscViewer);
184: SLEPC_EXTERN PetscErrorCode NEPErrorViewFromOptions(NEP);
185: SLEPC_EXTERN PetscErrorCode NEPConvergedReasonView(NEP,PetscViewer);
186: SLEPC_EXTERN PetscErrorCode NEPConvergedReasonViewFromOptions(NEP);
187: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "NEPConvergedReasonView()", ) static inline PetscErrorCode NEPReasonView(NEP nep,PetscViewer v) {return NEPConvergedReasonView(nep,v);}
188: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "NEPConvergedReasonViewFromOptions()", ) static inline PetscErrorCode NEPReasonViewFromOptions(NEP nep) {return NEPConvergedReasonViewFromOptions(nep);}
189: SLEPC_EXTERN PetscErrorCode NEPValuesView(NEP,PetscViewer);
190: SLEPC_EXTERN PetscErrorCode NEPValuesViewFromOptions(NEP);
191: SLEPC_EXTERN PetscErrorCode NEPVectorsView(NEP,PetscViewer);
192: SLEPC_EXTERN PetscErrorCode NEPVectorsViewFromOptions(NEP);

194: /*S
195:    NEPFunctionFn - A prototype of a NEP function evaluation function that
196:    would be passed to NEPSetFunction()

198:    Calling Sequence:
199: +  nep    - the nonlinear eigensolver context
200: .  lambda - the scalar argument where T(.) must be evaluated
201: .  T      - matrix that will contain T(lambda)
202: .  P      - [optional] different matrix to build the preconditioner
203: -  ctx    - [optional] user-defined context for private data for the
204:             function evaluation routine (may be NULL)

206:    Level: beginner

208: .seealso: [](ch:nep), `NEPSetFunction()`
209: S*/
210: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPFunctionFn(NEP nep,PetscScalar lambda,Mat T,Mat P,void *ctx);

212: /*S
213:    NEPJacobianFn - A prototype of a NEP Jacobian evaluation function that would be passed to NEPSetJacobian()

215:    Calling Sequence:
216: +  nep    - the nonlinear eigensolver context
217: .  lambda - the scalar argument where T'(.) must be evaluated
218: .  J      - matrix that will contain T'(lambda)
219: -  ctx    - [optional] user-defined context for private data for the
220:             Jacobian evaluation routine (may be NULL)

222:    Level: beginner

224: .seealso: [](ch:nep), `NEPSetJacobian()`
225: S*/
226: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPJacobianFn(NEP nep,PetscScalar lambda,Mat J,void *ctx);

228: SLEPC_EXTERN PetscErrorCode NEPSetFunction(NEP,Mat,Mat,NEPFunctionFn*,void*);
229: SLEPC_EXTERN PetscErrorCode NEPGetFunction(NEP,Mat*,Mat*,NEPFunctionFn**,void**);
230: SLEPC_EXTERN PetscErrorCode NEPSetJacobian(NEP,Mat,NEPJacobianFn*,void*);
231: SLEPC_EXTERN PetscErrorCode NEPGetJacobian(NEP,Mat*,NEPJacobianFn**,void**);
232: 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");}
233: 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");}
234: SLEPC_EXTERN PetscErrorCode NEPSetSplitOperator(NEP,PetscInt,Mat[],FN[],MatStructure);
235: SLEPC_EXTERN PetscErrorCode NEPGetSplitOperatorTerm(NEP,PetscInt,Mat*,FN*);
236: SLEPC_EXTERN PetscErrorCode NEPGetSplitOperatorInfo(NEP,PetscInt*,MatStructure*);
237: SLEPC_EXTERN PetscErrorCode NEPSetSplitPreconditioner(NEP,PetscInt,Mat[],MatStructure);
238: SLEPC_EXTERN PetscErrorCode NEPGetSplitPreconditionerTerm(NEP,PetscInt,Mat*);
239: SLEPC_EXTERN PetscErrorCode NEPGetSplitPreconditionerInfo(NEP,PetscInt*,MatStructure*);

241: SLEPC_EXTERN PetscErrorCode NEPSetBV(NEP,BV);
242: SLEPC_EXTERN PetscErrorCode NEPGetBV(NEP,BV*);
243: SLEPC_EXTERN PetscErrorCode NEPSetRG(NEP,RG);
244: SLEPC_EXTERN PetscErrorCode NEPGetRG(NEP,RG*);
245: SLEPC_EXTERN PetscErrorCode NEPSetDS(NEP,DS);
246: SLEPC_EXTERN PetscErrorCode NEPGetDS(NEP,DS*);
247: SLEPC_EXTERN PetscErrorCode NEPRefineGetKSP(NEP,KSP*);
248: SLEPC_EXTERN PetscErrorCode NEPSetTolerances(NEP,PetscReal,PetscInt);
249: SLEPC_EXTERN PetscErrorCode NEPGetTolerances(NEP,PetscReal*,PetscInt*);
250: SLEPC_EXTERN PetscErrorCode NEPSetDimensions(NEP,PetscInt,PetscInt,PetscInt);
251: SLEPC_EXTERN PetscErrorCode NEPGetDimensions(NEP,PetscInt*,PetscInt*,PetscInt*);
252: SLEPC_EXTERN PetscErrorCode NEPSetRefine(NEP,NEPRefine,PetscInt,PetscReal,PetscInt,NEPRefineScheme);
253: SLEPC_EXTERN PetscErrorCode NEPGetRefine(NEP,NEPRefine*,PetscInt*,PetscReal*,PetscInt*,NEPRefineScheme*);

255: SLEPC_EXTERN PetscErrorCode NEPGetConverged(NEP,PetscInt*);
256: SLEPC_EXTERN PetscErrorCode NEPGetEigenpair(NEP,PetscInt,PetscScalar*,PetscScalar*,Vec,Vec);
257: SLEPC_EXTERN PetscErrorCode NEPGetLeftEigenvector(NEP,PetscInt,Vec,Vec);

259: SLEPC_EXTERN PetscErrorCode NEPComputeError(NEP,PetscInt,NEPErrorType,PetscReal*);
260: 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);}
261: 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);}
262: SLEPC_EXTERN PetscErrorCode NEPGetErrorEstimate(NEP,PetscInt,PetscReal*);

264: SLEPC_EXTERN PetscErrorCode NEPComputeFunction(NEP,PetscScalar,Mat,Mat);
265: SLEPC_EXTERN PetscErrorCode NEPComputeJacobian(NEP,PetscScalar,Mat);
266: SLEPC_EXTERN PetscErrorCode NEPApplyFunction(NEP,PetscScalar,Vec,Vec,Vec,Mat,Mat);
267: SLEPC_EXTERN PetscErrorCode NEPApplyAdjoint(NEP,PetscScalar,Vec,Vec,Vec,Mat,Mat);
268: SLEPC_EXTERN PetscErrorCode NEPApplyJacobian(NEP,PetscScalar,Vec,Vec,Vec,Mat);
269: SLEPC_EXTERN PetscErrorCode NEPProjectOperator(NEP,PetscInt,PetscInt);
270: SLEPC_EXTERN PetscErrorCode NEPGetIterationNumber(NEP,PetscInt*);

272: SLEPC_EXTERN PetscErrorCode NEPSetInitialSpace(NEP,PetscInt,Vec[]);
273: SLEPC_EXTERN PetscErrorCode NEPSetWhichEigenpairs(NEP,NEPWhich);
274: SLEPC_EXTERN PetscErrorCode NEPGetWhichEigenpairs(NEP,NEPWhich*);
275: SLEPC_EXTERN PetscErrorCode NEPSetTwoSided(NEP,PetscBool);
276: SLEPC_EXTERN PetscErrorCode NEPGetTwoSided(NEP,PetscBool*);

278: SLEPC_EXTERN PetscErrorCode NEPApplyResolvent(NEP,RG,PetscScalar,Vec,Vec);

280: SLEPC_EXTERN PetscErrorCode NEPSetTrackAll(NEP,PetscBool);
281: SLEPC_EXTERN PetscErrorCode NEPGetTrackAll(NEP,PetscBool*);

283: SLEPC_EXTERN PetscErrorCode NEPGetConvergedReason(NEP,NEPConvergedReason*);

285: /*S
286:    NEPMonitorFn - A function prototype for functions provided to `NEPMonitorSet()`.

288:    Calling Sequence:
289: +  nep    - the nonlinear eigensolver context
290: .  its    - iteration number
291: .  nconv  - number of converged eigenpairs
292: .  eigr   - real part of the eigenvalues
293: .  eigi   - imaginary part of the eigenvalues
294: .  errest - relative error estimates for each eigenpair
295: .  nest   - number of error estimates
296: -  ctx    - optional monitoring context, as provided with `NEPMonitorSet()`

298:    Level: intermediate

300: .seealso: [](ch:nep), `NEPMonitorSet()`
301: S*/
302: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorFn(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,void *ctx);

304: /*S
305:    NEPMonitorRegisterFn - A function prototype for functions provided to `NEPMonitorRegister()`.

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

317:    Level: advanced

319:    Note:
320:    This is a `NEPMonitorFn` specialized for a context of `PetscViewerAndFormat`.

322: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorRegister()`, `NEPMonitorFn`, `NEPMonitorRegisterCreateFn`, `NEPMonitorRegisterDestroyFn`
323: S*/
324: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorRegisterFn(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *ctx);

326: /*S
327:    NEPMonitorRegisterCreateFn - A function prototype for functions that do the
328:    creation when provided to `NEPMonitorRegister()`.

330:    Calling Sequence:
331: +  viewer - the viewer to be used with the `NEPMonitorRegisterFn`
332: .  format - the format of the viewer
333: .  ctx    - a context for the monitor
334: -  result - a `PetscViewerAndFormat` object

336:    Level: advanced

338: .seealso: [](ch:nep), `NEPMonitorRegisterFn`, `NEPMonitorSet()`, `NEPMonitorRegister()`, `NEPMonitorFn`, `NEPMonitorRegisterDestroyFn`
339: S*/
340: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorRegisterCreateFn(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **result);

342: /*S
343:    NEPMonitorRegisterDestroyFn - A function prototype for functions that do the after
344:    use destruction when provided to `NEPMonitorRegister()`.

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

349:    Level: advanced

351: .seealso: [](ch:nep), `NEPMonitorRegisterFn`, `NEPMonitorSet()`, `NEPMonitorRegister()`, `NEPMonitorFn`, `NEPMonitorRegisterCreateFn`
352: S*/
353: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPMonitorRegisterDestroyFn(PetscViewerAndFormat **result);

355: SLEPC_EXTERN PetscErrorCode NEPMonitor(NEP,PetscInt,PetscInt,PetscScalar[],PetscScalar[],PetscReal[],PetscInt);
356: SLEPC_EXTERN PetscErrorCode NEPMonitorSet(NEP,NEPMonitorFn,void*,PetscCtxDestroyFn*);
357: SLEPC_EXTERN PetscErrorCode NEPMonitorCancel(NEP);
358: SLEPC_EXTERN PetscErrorCode NEPGetMonitorContext(NEP,void*);

360: SLEPC_EXTERN PetscErrorCode NEPMonitorSetFromOptions(NEP,const char[],const char[],void*,PetscBool);
361: SLEPC_EXTERN NEPMonitorRegisterFn        NEPMonitorFirst;
362: SLEPC_EXTERN NEPMonitorRegisterFn        NEPMonitorFirstDrawLG;
363: SLEPC_EXTERN NEPMonitorRegisterCreateFn  NEPMonitorFirstDrawLGCreate;
364: SLEPC_EXTERN NEPMonitorRegisterFn        NEPMonitorAll;
365: SLEPC_EXTERN NEPMonitorRegisterFn        NEPMonitorAllDrawLG;
366: SLEPC_EXTERN NEPMonitorRegisterCreateFn  NEPMonitorAllDrawLGCreate;
367: SLEPC_EXTERN NEPMonitorRegisterFn        NEPMonitorConverged;
368: SLEPC_EXTERN NEPMonitorRegisterCreateFn  NEPMonitorConvergedCreate;
369: SLEPC_EXTERN NEPMonitorRegisterFn        NEPMonitorConvergedDrawLG;
370: SLEPC_EXTERN NEPMonitorRegisterCreateFn  NEPMonitorConvergedDrawLGCreate;
371: SLEPC_EXTERN NEPMonitorRegisterDestroyFn NEPMonitorConvergedDestroy;

373: SLEPC_EXTERN PetscErrorCode NEPSetOptionsPrefix(NEP,const char[]);
374: SLEPC_EXTERN PetscErrorCode NEPAppendOptionsPrefix(NEP,const char[]);
375: SLEPC_EXTERN PetscErrorCode NEPGetOptionsPrefix(NEP,const char*[]);

377: SLEPC_EXTERN PetscFunctionList NEPList;
378: SLEPC_EXTERN PetscFunctionList NEPMonitorList;
379: SLEPC_EXTERN PetscFunctionList NEPMonitorCreateList;
380: SLEPC_EXTERN PetscFunctionList NEPMonitorDestroyList;
381: SLEPC_EXTERN PetscErrorCode NEPRegister(const char[],PetscErrorCode(*)(NEP));
382: SLEPC_EXTERN PetscErrorCode NEPMonitorRegister(const char[],PetscViewerType,PetscViewerFormat,NEPMonitorRegisterFn*,NEPMonitorRegisterCreateFn*,NEPMonitorRegisterDestroyFn*);

384: SLEPC_EXTERN PetscErrorCode NEPSetWorkVecs(NEP,PetscInt);
385: SLEPC_EXTERN PetscErrorCode NEPAllocateSolution(NEP,PetscInt);

387: /*S
388:    NEPConvergenceTestFn - A prototype of a NEP convergence test function that
389:    would be passed to NEPSetConvergenceTestFunction()

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

400:    Level: advanced

402: .seealso: [](ch:nep), `NEPSetConvergenceTestFunction()`
403: S*/
404: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPConvergenceTestFn(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx);

406: SLEPC_EXTERN PetscErrorCode NEPSetConvergenceTest(NEP,NEPConv);
407: SLEPC_EXTERN PetscErrorCode NEPGetConvergenceTest(NEP,NEPConv*);
408: SLEPC_EXTERN NEPConvergenceTestFn NEPConvergedAbsolute;
409: SLEPC_EXTERN NEPConvergenceTestFn NEPConvergedRelative;
410: SLEPC_EXTERN NEPConvergenceTestFn NEPConvergedNorm;
411: SLEPC_EXTERN PetscErrorCode NEPSetConvergenceTestFunction(NEP,NEPConvergenceTestFn*,void*,PetscCtxDestroyFn*);

413: /*S
414:    NEPStoppingTestFn - A prototype of a NEP stopping test function that would
415:    be passed to NEPSetStoppingTestFunction()

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

427:    Level: advanced

429: .seealso: [](ch:nep), `NEPSetStoppingTestFunction()`
430: S*/
431: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPStoppingTestFn(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx);

433: SLEPC_EXTERN PetscErrorCode NEPSetStoppingTest(NEP,NEPStop);
434: SLEPC_EXTERN PetscErrorCode NEPGetStoppingTest(NEP,NEPStop*);
435: SLEPC_EXTERN NEPStoppingTestFn NEPStoppingBasic;
436: SLEPC_EXTERN PetscErrorCode NEPSetStoppingTestFunction(NEP,NEPStoppingTestFn*,void*,PetscCtxDestroyFn*);

438: SLEPC_EXTERN PetscErrorCode NEPSetEigenvalueComparison(NEP,SlepcEigenvalueComparisonFn*,void*);

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

442: SLEPC_EXTERN PetscErrorCode NEPRIISetMaximumIterations(NEP,PetscInt);
443: SLEPC_EXTERN PetscErrorCode NEPRIIGetMaximumIterations(NEP,PetscInt*);
444: SLEPC_EXTERN PetscErrorCode NEPRIISetLagPreconditioner(NEP,PetscInt);
445: SLEPC_EXTERN PetscErrorCode NEPRIIGetLagPreconditioner(NEP,PetscInt*);
446: SLEPC_EXTERN PetscErrorCode NEPRIISetConstCorrectionTol(NEP,PetscBool);
447: SLEPC_EXTERN PetscErrorCode NEPRIIGetConstCorrectionTol(NEP,PetscBool*);
448: SLEPC_EXTERN PetscErrorCode NEPRIISetHermitian(NEP,PetscBool);
449: SLEPC_EXTERN PetscErrorCode NEPRIIGetHermitian(NEP,PetscBool*);
450: SLEPC_EXTERN PetscErrorCode NEPRIISetDeflationThreshold(NEP,PetscReal);
451: SLEPC_EXTERN PetscErrorCode NEPRIIGetDeflationThreshold(NEP,PetscReal*);
452: SLEPC_EXTERN PetscErrorCode NEPRIISetKSP(NEP,KSP);
453: SLEPC_EXTERN PetscErrorCode NEPRIIGetKSP(NEP,KSP*);

455: SLEPC_EXTERN PetscErrorCode NEPSLPSetDeflationThreshold(NEP,PetscReal);
456: SLEPC_EXTERN PetscErrorCode NEPSLPGetDeflationThreshold(NEP,PetscReal*);
457: SLEPC_EXTERN PetscErrorCode NEPSLPSetEPS(NEP,EPS);
458: SLEPC_EXTERN PetscErrorCode NEPSLPGetEPS(NEP,EPS*);
459: SLEPC_EXTERN PetscErrorCode NEPSLPSetEPSLeft(NEP,EPS);
460: SLEPC_EXTERN PetscErrorCode NEPSLPGetEPSLeft(NEP,EPS*);
461: SLEPC_EXTERN PetscErrorCode NEPSLPSetKSP(NEP,KSP);
462: SLEPC_EXTERN PetscErrorCode NEPSLPGetKSP(NEP,KSP*);

464: SLEPC_EXTERN PetscErrorCode NEPNArnoldiSetKSP(NEP,KSP);
465: SLEPC_EXTERN PetscErrorCode NEPNArnoldiGetKSP(NEP,KSP*);
466: SLEPC_EXTERN PetscErrorCode NEPNArnoldiSetLagPreconditioner(NEP,PetscInt);
467: SLEPC_EXTERN PetscErrorCode NEPNArnoldiGetLagPreconditioner(NEP,PetscInt*);

469: /*E
470:    NEPCISSExtraction - determines the extraction technique in the CISS solver

472:    Level: advanced

474: .seealso: [](ch:nep), `NEPCISSSetExtraction()`, `NEPCISSGetExtraction()`
475: E*/
476: typedef enum { NEP_CISS_EXTRACTION_RITZ,
477:                NEP_CISS_EXTRACTION_HANKEL,
478:                NEP_CISS_EXTRACTION_CAA    } NEPCISSExtraction;
479: SLEPC_EXTERN const char *NEPCISSExtractions[];

481: #if defined(PETSC_USE_COMPLEX) || defined(PETSC_CLANG_STATIC_ANALYZER)
482: SLEPC_EXTERN PetscErrorCode NEPCISSSetExtraction(NEP,NEPCISSExtraction);
483: SLEPC_EXTERN PetscErrorCode NEPCISSGetExtraction(NEP,NEPCISSExtraction*);
484: SLEPC_EXTERN PetscErrorCode NEPCISSSetSizes(NEP,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscBool);
485: SLEPC_EXTERN PetscErrorCode NEPCISSGetSizes(NEP,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscInt*,PetscBool*);
486: SLEPC_EXTERN PetscErrorCode NEPCISSSetThreshold(NEP,PetscReal,PetscReal);
487: SLEPC_EXTERN PetscErrorCode NEPCISSGetThreshold(NEP,PetscReal*,PetscReal*);
488: SLEPC_EXTERN PetscErrorCode NEPCISSSetRefinement(NEP,PetscInt,PetscInt);
489: SLEPC_EXTERN PetscErrorCode NEPCISSGetRefinement(NEP,PetscInt*,PetscInt*);
490: SLEPC_EXTERN PetscErrorCode NEPCISSGetKSPs(NEP,PetscInt*,KSP*[]);
491: #else
492: #define SlepcNEPCISSUnavailable(nep) do { \
493:     PetscFunctionBegin; \
494:     SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"%s() not available with real scalars",PETSC_FUNCTION_NAME); \
495:     } while (0)
496: static inline PetscErrorCode NEPCISSSetExtraction(NEP nep,PETSC_UNUSED NEPCISSExtraction ex) {SlepcNEPCISSUnavailable(nep);}
497: static inline PetscErrorCode NEPCISSGetExtraction(NEP nep,PETSC_UNUSED NEPCISSExtraction *ex) {SlepcNEPCISSUnavailable(nep);}
498: 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);}
499: 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);}
500: static inline PetscErrorCode NEPCISSSetThreshold(NEP nep,PETSC_UNUSED PetscReal delta,PETSC_UNUSED PetscReal spur) {SlepcNEPCISSUnavailable(nep);}
501: static inline PetscErrorCode NEPCISSGetThreshold(NEP nep,PETSC_UNUSED PetscReal *delta,PETSC_UNUSED PetscReal *spur) {SlepcNEPCISSUnavailable(nep);}
502: static inline PetscErrorCode NEPCISSSetRefinement(NEP nep,PETSC_UNUSED PetscInt inner,PETSC_UNUSED PetscInt blsize) {SlepcNEPCISSUnavailable(nep);}
503: static inline PetscErrorCode NEPCISSGetRefinement(NEP nep,PETSC_UNUSED PetscInt *inner,PETSC_UNUSED PetscInt *blsize) {SlepcNEPCISSUnavailable(nep);}
504: static inline PetscErrorCode NEPCISSGetKSPs(NEP nep,PETSC_UNUSED PetscInt *nsolve,PETSC_UNUSED KSP *ksp[]) {SlepcNEPCISSUnavailable(nep);}
505: #undef SlepcNEPCISSUnavailable
506: #endif

508: SLEPC_EXTERN PetscErrorCode NEPInterpolSetPEP(NEP,PEP);
509: SLEPC_EXTERN PetscErrorCode NEPInterpolGetPEP(NEP,PEP*);
510: SLEPC_EXTERN PetscErrorCode NEPInterpolSetInterpolation(NEP,PetscReal,PetscInt);
511: SLEPC_EXTERN PetscErrorCode NEPInterpolGetInterpolation(NEP,PetscReal*,PetscInt*);

513: /*S
514:    NEPNLEIGSSingularitiesFn - A prototype of a function that would be passed
515:    to NEPNLEIGSSetSingularitiesFunction()

517:    Calling Sequence:
518: +  nep   - the nonlinear eigensolver context
519: .  maxnp - on input number of requested points in the discretization (can be set)
520: .  xi    - computed values of the discretization
521: -  ctx   - optional context, as set by NEPNLEIGSSetSingularitiesFunction()

523:    Level: advanced

525: .seealso: [](ch:nep), `NEPNLEIGSSetSingularitiesFunction()`
526: S*/
527: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode NEPNLEIGSSingularitiesFn(NEP nep,PetscInt *maxnp,PetscScalar *xi,void *ctx);

529: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetSingularitiesFunction(NEP,NEPNLEIGSSingularitiesFn*,void*);
530: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetSingularitiesFunction(NEP,NEPNLEIGSSingularitiesFn**,void**);
531: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetRestart(NEP,PetscReal);
532: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetRestart(NEP,PetscReal*);
533: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetLocking(NEP,PetscBool);
534: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetLocking(NEP,PetscBool*);
535: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetInterpolation(NEP,PetscReal,PetscInt);
536: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetInterpolation(NEP,PetscReal*,PetscInt*);
537: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetRKShifts(NEP,PetscInt,PetscScalar[]);
538: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetRKShifts(NEP,PetscInt*,PetscScalar*[]);
539: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetKSPs(NEP,PetscInt*,KSP*[]);
540: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetFullBasis(NEP,PetscBool);
541: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetFullBasis(NEP,PetscBool*);
542: SLEPC_EXTERN PetscErrorCode NEPNLEIGSSetEPS(NEP,EPS);
543: SLEPC_EXTERN PetscErrorCode NEPNLEIGSGetEPS(NEP,EPS*);