1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
8: SLEPc is free software: you can redistribute it and/or modify it under the
9: terms of version 3 of the GNU Lesser General Public License as published by
10: the Free Software Foundation.
12: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
13: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15: more details.
17: You should have received a copy of the GNU Lesser General Public License
18: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
19: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: */
22: #if !defined(_NEPIMPL)
23: #define _NEPIMPL 25: #include <slepcnep.h>
26: #include <slepc/private/slepcimpl.h>
28: PETSC_EXTERN PetscBool NEPRegisterAllCalled;
29: PETSC_EXTERN PetscErrorCode NEPRegisterAll(void);
30: PETSC_EXTERN PetscLogEvent NEP_SetUp,NEP_Solve,NEP_Refine,NEP_FunctionEval,NEP_JacobianEval,NEP_DerivativesEval;
32: typedef struct _NEPOps *NEPOps;
34: struct _NEPOps {
35: PetscErrorCode (*solve)(NEP);
36: PetscErrorCode (*setup)(NEP);
37: PetscErrorCode (*setfromoptions)(PetscOptionItems*,NEP);
38: PetscErrorCode (*publishoptions)(NEP);
39: PetscErrorCode (*destroy)(NEP);
40: PetscErrorCode (*reset)(NEP);
41: PetscErrorCode (*view)(NEP,PetscViewer);
42: PetscErrorCode (*computevectors)(NEP);
43: };
45: /*
46: Maximum number of monitors you can run with a single NEP 47: */
48: #define MAXNEPMONITORS 5 50: typedef enum { NEP_STATE_INITIAL,
51: NEP_STATE_SETUP,
52: NEP_STATE_SOLVED,
53: NEP_STATE_EIGENVECTORS } NEPStateType;
55: /*
56: How the problem function T(lambda) has been defined by the user
57: - Callback: one callback to build the function matrix, another one for the Jacobian
58: - Split: in split form sum_j(A_j*f_j(lambda))
59: - Derivatives: a single callback for all the derivatives (including the 0th derivative)
60: */
61: typedef enum { NEP_USER_INTERFACE_CALLBACK=1,
62: NEP_USER_INTERFACE_SPLIT,
63: NEP_USER_INTERFACE_DERIVATIVES } NEPUserInterface;
65: /*
66: Defines the NEP data structure.
67: */
68: struct _p_NEP {
69: PETSCHEADER(struct _NEPOps);
70: /*------------------------- User parameters ---------------------------*/
71: PetscInt max_it; /* maximum number of iterations */
72: PetscInt nev; /* number of eigenvalues to compute */
73: PetscInt ncv; /* number of basis vectors */
74: PetscInt mpd; /* maximum dimension of projected problem */
75: PetscInt nini; /* number of initial vectors (negative means not copied yet) */
76: PetscScalar target; /* target value */
77: PetscReal tol; /* tolerance */
78: NEPConv conv; /* convergence test */
79: NEPStop stop; /* stopping test */
80: NEPWhich which; /* which part of the spectrum to be sought */
81: NEPRefine refine; /* type of refinement to be applied after solve */
82: PetscInt npart; /* number of partitions of the communicator */
83: PetscReal rtol; /* tolerance for refinement */
84: PetscInt rits; /* number of iterations of the refinement method */
85: NEPRefineScheme scheme; /* scheme for solving linear systems within refinement */
86: PetscBool trackall; /* whether all the residuals must be computed */
88: /*-------------- User-provided functions and contexts -----------------*/
89: PetscErrorCode (*computefunction)(NEP,PetscScalar,Mat,Mat,void*);
90: PetscErrorCode (*computejacobian)(NEP,PetscScalar,Mat,void*);
91: void *functionctx;
92: void *jacobianctx;
93: PetscErrorCode (*computederivatives)(NEP,PetscScalar,PetscInt,Mat,void*);
94: void *derivativesctx;
95: PetscErrorCode (*converged)(NEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
96: PetscErrorCode (*convergeddestroy)(void*);
97: PetscErrorCode (*stopping)(NEP,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*);
98: PetscErrorCode (*stoppingdestroy)(void*);
99: void *convergedctx;
100: void *stoppingctx;
101: PetscErrorCode (*monitor[MAXNEPMONITORS])(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
102: PetscErrorCode (*monitordestroy[MAXNEPMONITORS])(void**);
103: void *monitorcontext[MAXNEPMONITORS];
104: PetscInt numbermonitors;
106: /*----------------- Child objects and working data -------------------*/
107: DS ds; /* direct solver object */
108: BV V; /* set of basis vectors and computed eigenvectors */
109: RG rg; /* optional region for filtering */
110: SlepcSC sc; /* sorting criterion data */
111: Mat function; /* function matrix */
112: Mat function_pre; /* function matrix (preconditioner) */
113: Mat jacobian; /* Jacobian matrix */
114: Mat derivatives; /* derivatives matrix */
115: Mat *A; /* matrix coefficients of split form */
116: FN *f; /* matrix functions of split form */
117: PetscInt nt; /* number of terms in split form */
118: MatStructure mstr; /* pattern of split matrices */
119: Vec *IS; /* references to user-provided initial space */
120: PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */
121: PetscReal *errest; /* error estimates */
122: PetscInt *perm; /* permutation for eigenvalue ordering */
123: PetscInt nwork; /* number of work vectors */
124: Vec *work; /* work vectors */
125: KSP refineksp; /* ksp used in refinement */
126: PetscSubcomm refinesubc; /* context for sub-communicators */
127: void *data; /* placeholder for solver-specific stuff */
129: /* ----------------------- Status variables --------------------------*/
130: NEPStateType state; /* initial -> setup -> solved -> eigenvectors */
131: PetscInt nconv; /* number of converged eigenvalues */
132: PetscInt its; /* number of iterations so far computed */
133: PetscInt n,nloc; /* problem dimensions (global, local) */
134: PetscReal *nrma; /* computed matrix norms */
135: NEPUserInterface fui; /* how the user has defined the nonlinear operator */
136: NEPConvergedReason reason;
137: };
139: /*
140: Macros to test valid NEP arguments
141: */
142: #if !defined(PETSC_USE_DEBUG)
144: #define NEPCheckProblem(h,arg) do {} while (0)145: #define NEPCheckCallback(h,arg) do {} while (0)146: #define NEPCheckSplit(h,arg) do {} while (0)147: #define NEPCheckDerivatives(h,arg) do {} while (0)148: #define NEPCheckSolved(h,arg) do {} while (0)150: #else
152: #define NEPCheckProblem(h,arg) \153: do { \154: if (!(h->fui)) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"The nonlinear eigenproblem has not been specified yet. Parameter #%d",arg); \155: } while (0)157: #define NEPCheckCallback(h,arg) \158: do { \159: 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); \160: } while (0)162: #define NEPCheckSplit(h,arg) \163: do { \164: 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); \165: } while (0)167: #define NEPCheckDerivatives(h,arg) \168: do { \169: 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); \170: } while (0)172: #define NEPCheckSolved(h,arg) \173: do { \174: if (h->state<NEP_STATE_SOLVED) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"Must call NEPSolve() first: Parameter #%d",arg); \175: } while (0)177: #endif
179: PETSC_INTERN PetscErrorCode NEPSetDimensions_Default(NEP,PetscInt,PetscInt*,PetscInt*);
180: PETSC_INTERN PetscErrorCode NEPComputeVectors(NEP);
181: PETSC_INTERN PetscErrorCode NEPReset_Problem(NEP);
182: PETSC_INTERN PetscErrorCode NEPGetDefaultShift(NEP,PetscScalar*);
183: PETSC_INTERN PetscErrorCode NEPComputeVectors_Schur(NEP);
184: PETSC_INTERN PetscErrorCode NEPComputeResidualNorm_Private(NEP,PetscScalar,Vec,Vec*,PetscReal*);
185: PETSC_INTERN PetscErrorCode NEPNewtonRefinementSimple(NEP,PetscInt*,PetscReal,PetscInt);
187: #endif