1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2015, 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;
32: typedef struct _NEPOps *NEPOps;
34: struct _NEPOps {
35: PetscErrorCode (*solve)(NEP);
36: PetscErrorCode (*setup)(NEP);
37: PetscErrorCode (*setfromoptions)(PetscOptions*,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: Defines the NEP data structure.
57: */
58: struct _p_NEP {
59: PETSCHEADER(struct _NEPOps);
60: /*------------------------- User parameters ---------------------------*/
61: PetscInt max_it; /* maximum number of iterations */
62: PetscInt max_funcs; /* maximum number of function evaluations */
63: PetscInt nev; /* number of eigenvalues to compute */
64: PetscInt ncv; /* number of basis vectors */
65: PetscInt mpd; /* maximum dimension of projected problem */
66: PetscInt lag; /* interval to rebuild preconditioner */
67: PetscInt nini; /* number of initial vectors (negative means not copied yet) */
68: PetscScalar target; /* target value */
69: PetscReal abstol,rtol,stol; /* user tolerances */
70: PetscReal ktol; /* tolerance for linear solver */
71: PetscBool cctol; /* constant correction tolerance */
72: PetscReal ttol; /* tolerance used in the convergence criterion */
73: NEPWhich which; /* which part of the spectrum to be sought */
74: NEPRefine refine; /* type of refinement to be applied after solve */
75: PetscInt npart; /* number of partitions of the communicator */
76: PetscReal reftol; /* tolerance for refinement */
77: PetscInt rits; /* number of iterations of the refinement method */
78: PetscBool trackall; /* whether all the residuals must be computed */
80: /*-------------- User-provided functions and contexts -----------------*/
81: PetscErrorCode (*computefunction)(NEP,PetscScalar,Mat,Mat,void*);
82: PetscErrorCode (*computejacobian)(NEP,PetscScalar,Mat,void*);
83: void *functionctx;
84: void *jacobianctx;
85: PetscErrorCode (*converged)(NEP,PetscInt,PetscReal,PetscReal,PetscReal,NEPConvergedReason*,void*);
86: PetscErrorCode (*convergeddestroy)(void*);
87: void *convergedctx;
88: PetscErrorCode (*monitor[MAXNEPMONITORS])(NEP,PetscInt,PetscInt,PetscScalar*,PetscReal*,PetscInt,void*);
89: PetscErrorCode (*monitordestroy[MAXNEPMONITORS])(void**);
90: void *monitorcontext[MAXNEPMONITORS];
91: PetscInt numbermonitors;
93: /*----------------- Child objects and working data -------------------*/
94: DS ds; /* direct solver object */
95: BV V; /* set of basis vectors and computed eigenvectors */
96: RG rg; /* optional region for filtering */
97: PetscRandom rand; /* random number generator */
98: SlepcSC sc; /* sorting criterion data */
99: KSP ksp; /* linear solver object */
100: Mat function; /* function matrix */
101: Mat function_pre; /* function matrix (preconditioner) */
102: Mat jacobian; /* Jacobian matrix */
103: Mat *A; /* matrix coefficients of split form */
104: FN *f; /* matrix functions of split form */
105: PetscInt nt; /* number of terms in split form */
106: MatStructure mstr; /* pattern of split matrices */
107: Vec *IS; /* references to user-provided initial space */
108: PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */
109: PetscReal *errest; /* error estimates */
110: PetscInt *perm; /* permutation for eigenvalue ordering */
111: PetscInt nwork; /* number of work vectors */
112: Vec *work; /* work vectors */
113: void *data; /* placeholder for solver-specific stuff */
115: /* ----------------------- Status variables --------------------------*/
116: NEPStateType state; /* initial -> setup -> solved -> eigenvectors */
117: PetscInt nconv; /* number of converged eigenvalues */
118: PetscInt its; /* number of iterations so far computed */
119: PetscInt n,nloc; /* problem dimensions (global, local) */
120: PetscInt nfuncs; /* number of function evaluations */
121: PetscBool split; /* the nonlinear operator has been set in
122: split form, otherwise user callbacks are used */
123: NEPConvergedReason reason;
124: };
126: /*
127: Macros to test valid NEP arguments
128: */
129: #if !defined(PETSC_USE_DEBUG)
131: #define NEPCheckSolved(h,arg) do {} while (0)133: #else
135: #define NEPCheckSolved(h,arg) \136: do { \137: if (h->state<NEP_STATE_SOLVED) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"Must call NEPSolve() first: Parameter #%d",arg); \138: } while (0)140: #endif
144: PETSC_STATIC_INLINE PetscErrorCode NEP_KSPSolve(NEP nep,Vec b,Vec x)145: {
147: PetscInt lits;
150: KSPSolve(nep->ksp,b,x);
151: KSPGetIterationNumber(nep->ksp,&lits);
152: PetscInfo2(nep,"iter=%D, linear solve iterations=%D\n",nep->its,lits);
153: return(0);
154: }
156: PETSC_INTERN PetscErrorCode NEPComputeVectors(NEP);
157: PETSC_INTERN PetscErrorCode NEPGetDefaultShift(NEP,PetscScalar*);
158: PETSC_INTERN PetscErrorCode NEPComputeResidualNorm_Private(NEP,PetscScalar,Vec,Vec*,PetscReal*);
159: PETSC_INTERN PetscErrorCode NEPNewtonRefinementSimple(NEP,PetscInt*,PetscReal*,PetscInt);
161: #endif