Actual source code: nepimpl.h
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2013, 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 PetscLogEvent NEP_SetUp,NEP_Solve,NEP_Dense,NEP_FunctionEval,NEP_JacobianEval;
30: typedef struct _NEPOps *NEPOps;
32: struct _NEPOps {
33: PetscErrorCode (*solve)(NEP);
34: PetscErrorCode (*setup)(NEP);
35: PetscErrorCode (*setfromoptions)(NEP);
36: PetscErrorCode (*publishoptions)(NEP);
37: PetscErrorCode (*destroy)(NEP);
38: PetscErrorCode (*reset)(NEP);
39: PetscErrorCode (*view)(NEP,PetscViewer);
40: };
42: /*
43: Maximum number of monitors you can run with a single NEP
44: */
45: #define MAXNEPMONITORS 5
47: /*
48: Defines the NEP data structure.
49: */
50: struct _p_NEP {
51: PETSCHEADER(struct _NEPOps);
52: /*------------------------- User parameters --------------------------*/
53: PetscInt max_it; /* maximum number of iterations */
54: PetscInt max_funcs; /* maximum number of function evaluations */
55: PetscInt nev; /* number of eigenvalues to compute */
56: PetscInt ncv; /* number of basis vectors */
57: PetscInt mpd; /* maximum dimension of projected problem */
58: PetscInt lag; /* interval to rebuild preconditioner */
59: PetscInt nini; /* number of initial vectors (negative means not copied yet) */
60: PetscScalar target; /* target value */
61: PetscReal abstol,rtol,stol; /* user tolerances */
62: PetscReal ktol; /* tolerance for linear solver */
63: PetscBool cctol; /* constant correction tolerance */
64: PetscReal ttol; /* tolerance used in the convergence criterion */
65: PetscBool trackall; /* whether all the residuals must be computed */
66: NEPWhich which; /* which part of the spectrum to be sought */
67: PetscBool split; /* the nonlinear operator has been set in
68: split form, otherwise user callbacks are used */
70: /*-------------- User-provided functions and contexts -----------------*/
71: PetscErrorCode (*computefunction)(NEP,PetscScalar,Mat*,Mat*,MatStructure*,void*);
72: PetscErrorCode (*computejacobian)(NEP,PetscScalar,Mat*,MatStructure*,void*);
73: PetscErrorCode (*comparison)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*);
74: PetscErrorCode (*converged)(NEP,PetscInt,PetscReal,PetscReal,PetscReal,NEPConvergedReason*,void*);
75: PetscErrorCode (*convergeddestroy)(void*);
76: void *comparisonctx;
77: void *convergedctx;
78: Mat function,function_pre;
79: void *functionctx;
80: Mat jacobian;
81: void *jacobianctx;
82: PetscInt nt; /* number of terms in split form */
83: MatStructure mstr; /* pattern of split matrices */
84: Mat *A; /* matrix coefficients of split form */
85: FN *f; /* matrix functions of split form */
87: /*------------------------- Working data --------------------------*/
88: Vec *V; /* set of basis vectors and computed eigenvectors */
89: Vec *IS; /* placeholder for references to user-provided initial space */
90: PetscScalar *eig; /* computed eigenvalues */
91: PetscReal *errest; /* error estimates */
92: IP ip; /* innerproduct object */
93: DS ds; /* direct solver object */
94: KSP ksp; /* linear solver object */
95: void *data; /* placeholder for misc stuff associated
96: with a particular solver */
97: PetscInt nconv; /* number of converged eigenvalues */
98: PetscInt its; /* number of iterations so far computed */
99: PetscInt *perm; /* permutation for eigenvalue ordering */
100: PetscInt nfuncs,linits; /* operation counters */
101: PetscInt n,nloc; /* problem dimensions (global, local) */
102: PetscRandom rand; /* random number generator */
103: Vec t; /* template vector */
104: PetscInt allocated_ncv; /* number of basis vectors allocated */
106: /* ---------------- Default work-area and status vars -------------------- */
107: PetscInt nwork;
108: Vec *work;
110: PetscInt setupcalled;
111: NEPConvergedReason reason;
113: PetscErrorCode (*monitor[MAXNEPMONITORS])(NEP,PetscInt,PetscInt,PetscScalar*,PetscReal*,PetscInt,void*);
114: PetscErrorCode (*monitordestroy[MAXNEPMONITORS])(void**);
115: void *monitorcontext[MAXNEPMONITORS];
116: PetscInt numbermonitors;
117: };
119: PETSC_INTERN PetscErrorCode NEPReset_Default(NEP);
120: PETSC_INTERN PetscErrorCode NEPGetDefaultShift(NEP,PetscScalar*);
121: PETSC_INTERN PetscErrorCode NEPAllocateSolution(NEP);
122: PETSC_INTERN PetscErrorCode NEPFreeSolution(NEP);
123: PETSC_INTERN PetscErrorCode NEP_KSPSolve(NEP,Vec,Vec);
124: PETSC_INTERN PetscErrorCode NEPComputeResidualNorm_Private(NEP,PetscScalar,Vec,PetscReal*);
125: PETSC_INTERN PetscErrorCode NEPComputeRelativeError_Private(NEP,PetscScalar,Vec,PetscReal*);
127: #endif