Actual source code: nepimpl.h
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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: */
11: #if !defined(SLEPCNEPIMPL_H)
12: #define SLEPCNEPIMPL_H
14: #include <slepcnep.h>
15: #include <slepc/private/slepcimpl.h>
17: SLEPC_EXTERN PetscBool NEPRegisterAllCalled;
18: SLEPC_EXTERN PetscErrorCode NEPRegisterAll(void);
19: SLEPC_EXTERN PetscLogEvent NEP_SetUp,NEP_Solve,NEP_Refine,NEP_FunctionEval,NEP_JacobianEval,NEP_DerivativesEval,NEP_Resolvent;
21: typedef struct _NEPOps *NEPOps;
23: struct _NEPOps {
24: PetscErrorCode (*solve)(NEP);
25: PetscErrorCode (*setup)(NEP);
26: PetscErrorCode (*setfromoptions)(PetscOptionItems*,NEP);
27: PetscErrorCode (*publishoptions)(NEP);
28: PetscErrorCode (*destroy)(NEP);
29: PetscErrorCode (*reset)(NEP);
30: PetscErrorCode (*view)(NEP,PetscViewer);
31: PetscErrorCode (*computevectors)(NEP);
32: };
34: /*
35: Maximum number of monitors you can run with a single NEP
36: */
37: #define MAXNEPMONITORS 5
39: typedef enum { NEP_STATE_INITIAL,
40: NEP_STATE_SETUP,
41: NEP_STATE_SOLVED,
42: NEP_STATE_EIGENVECTORS } NEPStateType;
44: /*
45: How the problem function T(lambda) has been defined by the user
46: - Callback: one callback to build the function matrix, another one for the Jacobian
47: - Split: in split form sum_j(A_j*f_j(lambda))
48: - Derivatives: a single callback for all the derivatives (including the 0th derivative)
49: */
50: typedef enum { NEP_USER_INTERFACE_CALLBACK=1,
51: NEP_USER_INTERFACE_SPLIT,
52: NEP_USER_INTERFACE_DERIVATIVES } NEPUserInterface;
54: /*
55: Defines the NEP data structure.
56: */
57: struct _p_NEP {
58: PETSCHEADER(struct _NEPOps);
59: /*------------------------- User parameters ---------------------------*/
60: PetscInt max_it; /* maximum number of iterations */
61: PetscInt nev; /* number of eigenvalues to compute */
62: PetscInt ncv; /* number of basis vectors */
63: PetscInt mpd; /* maximum dimension of projected problem */
64: PetscInt nini; /* number of initial vectors (negative means not copied yet) */
65: PetscScalar target; /* target value */
66: PetscReal tol; /* tolerance */
67: NEPConv conv; /* convergence test */
68: NEPStop stop; /* stopping test */
69: NEPWhich which; /* which part of the spectrum to be sought */
70: NEPProblemType problem_type; /* which kind of problem to be solved */
71: NEPRefine refine; /* type of refinement to be applied after solve */
72: PetscInt npart; /* number of partitions of the communicator */
73: PetscReal rtol; /* tolerance for refinement */
74: PetscInt rits; /* number of iterations of the refinement method */
75: NEPRefineScheme scheme; /* scheme for solving linear systems within refinement */
76: PetscBool trackall; /* whether all the residuals must be computed */
77: PetscBool twosided; /* whether to compute left eigenvectors (two-sided solver) */
79: /*-------------- User-provided functions and contexts -----------------*/
80: PetscErrorCode (*computefunction)(NEP,PetscScalar,Mat,Mat,void*);
81: PetscErrorCode (*computejacobian)(NEP,PetscScalar,Mat,void*);
82: void *functionctx;
83: void *jacobianctx;
84: PetscErrorCode (*computederivatives)(NEP,PetscScalar,PetscInt,Mat,void*);
85: void *derivativesctx;
86: PetscErrorCode (*converged)(NEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
87: PetscErrorCode (*convergeduser)(NEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
88: PetscErrorCode (*convergeddestroy)(void*);
89: PetscErrorCode (*stopping)(NEP,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*);
90: PetscErrorCode (*stoppinguser)(NEP,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*);
91: PetscErrorCode (*stoppingdestroy)(void*);
92: void *convergedctx;
93: void *stoppingctx;
94: PetscErrorCode (*monitor[MAXNEPMONITORS])(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
95: PetscErrorCode (*monitordestroy[MAXNEPMONITORS])(void**);
96: void *monitorcontext[MAXNEPMONITORS];
97: PetscInt numbermonitors;
99: /*----------------- Child objects and working data -------------------*/
100: DS ds; /* direct solver object */
101: BV V; /* set of basis vectors and computed eigenvectors */
102: BV W; /* left basis vectors (if left eigenvectors requested) */
103: RG rg; /* optional region for filtering */
104: SlepcSC sc; /* sorting criterion data */
105: Mat function; /* function matrix */
106: Mat function_pre; /* function matrix (preconditioner) */
107: Mat jacobian; /* Jacobian matrix */
108: Mat derivatives; /* derivatives matrix */
109: Mat *A; /* matrix coefficients of split form */
110: FN *f; /* matrix functions of split form */
111: PetscInt nt; /* number of terms in split form */
112: MatStructure mstr; /* pattern of split matrices */
113: Vec *IS; /* references to user-provided initial space */
114: PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */
115: PetscReal *errest; /* error estimates */
116: PetscInt *perm; /* permutation for eigenvalue ordering */
117: PetscInt nwork; /* number of work vectors */
118: Vec *work; /* work vectors */
119: KSP refineksp; /* ksp used in refinement */
120: PetscSubcomm refinesubc; /* context for sub-communicators */
121: void *data; /* placeholder for solver-specific stuff */
123: /* ----------------------- Status variables --------------------------*/
124: NEPStateType state; /* initial -> setup -> solved -> eigenvectors */
125: PetscInt nconv; /* number of converged eigenvalues */
126: PetscInt its; /* number of iterations so far computed */
127: PetscInt n,nloc; /* problem dimensions (global, local) */
128: PetscReal *nrma; /* computed matrix norms */
129: NEPUserInterface fui; /* how the user has defined the nonlinear operator */
130: PetscBool useds; /* whether the solver uses the DS object or not */
131: PetscBool hasts; /* whether the solver has two-sided variant */
132: Mat resolvent; /* shell matrix to be used in NEPApplyResolvent */
133: NEPConvergedReason reason;
134: };
136: /*
137: Macros to test valid NEP arguments
138: */
139: #if !defined(PETSC_USE_DEBUG)
141: #define NEPCheckProblem(h,arg) do {} while (0)
142: #define NEPCheckCallback(h,arg) do {} while (0)
143: #define NEPCheckSplit(h,arg) do {} while (0)
144: #define NEPCheckDerivatives(h,arg) do {} while (0)
145: #define NEPCheckSolved(h,arg) do {} while (0)
147: #else
149: #define NEPCheckProblem(h,arg) \
150: do { \
151: if (!(h->fui)) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"The nonlinear eigenproblem has not been specified yet. Parameter #%d",arg); \
152: } while (0)
154: #define NEPCheckCallback(h,arg) \
155: do { \
156: if (h->fui!=NEP_USER_INTERFACE_CALLBACK) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"This operation requires the nonlinear eigenproblem specified with callbacks. Parameter #%d",arg); \
157: } while (0)
159: #define NEPCheckSplit(h,arg) \
160: do { \
161: if (h->fui!=NEP_USER_INTERFACE_SPLIT) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"This operation requires the nonlinear eigenproblem in split form. Parameter #%d",arg); \
162: } while (0)
164: #define NEPCheckDerivatives(h,arg) \
165: do { \
166: if (h->fui!=NEP_USER_INTERFACE_DERIVATIVES) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"This operation requires the nonlinear eigenproblem specified with derivatives callback. Parameter #%d",arg); \
167: } while (0)
169: #define NEPCheckSolved(h,arg) \
170: do { \
171: if (h->state<NEP_STATE_SOLVED) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"Must call NEPSolve() first: Parameter #%d",arg); \
172: } while (0)
174: #endif
176: SLEPC_INTERN PetscErrorCode NEPSetDimensions_Default(NEP,PetscInt,PetscInt*,PetscInt*);
177: SLEPC_INTERN PetscErrorCode NEPComputeVectors(NEP);
178: SLEPC_INTERN PetscErrorCode NEPReset_Problem(NEP);
179: SLEPC_INTERN PetscErrorCode NEPGetDefaultShift(NEP,PetscScalar*);
180: SLEPC_INTERN PetscErrorCode NEPComputeVectors_Schur(NEP);
181: SLEPC_INTERN PetscErrorCode NEPComputeResidualNorm_Private(NEP,PetscBool,PetscScalar,Vec,Vec*,PetscReal*);
182: SLEPC_INTERN PetscErrorCode NEPNewtonRefinementSimple(NEP,PetscInt*,PetscReal,PetscInt);
184: #endif