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(_PEPIMPL)
23: #define _PEPIMPL 25: #include <slepcpep.h>
26: #include <slepc/private/slepcimpl.h>
28: PETSC_EXTERN PetscBool PEPRegisterAllCalled;
29: PETSC_EXTERN PetscErrorCode PEPRegisterAll(void);
30: PETSC_EXTERN PetscLogEvent PEP_SetUp,PEP_Solve,PEP_Refine;
32: typedef struct _PEPOps *PEPOps;
34: struct _PEPOps {
35: PetscErrorCode (*solve)(PEP);
36: PetscErrorCode (*setup)(PEP);
37: PetscErrorCode (*setfromoptions)(PetscOptionItems*,PEP);
38: PetscErrorCode (*publishoptions)(PEP);
39: PetscErrorCode (*destroy)(PEP);
40: PetscErrorCode (*reset)(PEP);
41: PetscErrorCode (*view)(PEP,PetscViewer);
42: PetscErrorCode (*backtransform)(PEP);
43: PetscErrorCode (*computevectors)(PEP);
44: PetscErrorCode (*extractvectors)(PEP);
45: };
47: /*
48: Maximum number of monitors you can run with a single PEP 49: */
50: #define MAXPEPMONITORS 5 52: typedef enum { PEP_STATE_INITIAL,
53: PEP_STATE_SETUP,
54: PEP_STATE_SOLVED,
55: PEP_STATE_EIGENVECTORS } PEPStateType;
57: /*
58: Defines the PEP data structure.
59: */
60: struct _p_PEP {
61: PETSCHEADER(struct _PEPOps);
62: /*------------------------- User parameters ---------------------------*/
63: PetscInt max_it; /* maximum number of iterations */
64: PetscInt nev; /* number of eigenvalues to compute */
65: PetscInt ncv; /* number of basis vectors */
66: PetscInt mpd; /* maximum dimension of projected problem */
67: PetscInt nini; /* number of initial vectors (negative means not copied yet) */
68: PetscScalar target; /* target value */
69: PetscReal tol; /* tolerance */
70: PEPConv conv; /* convergence test */
71: PEPStop stop; /* stopping test */
72: PEPWhich which; /* which part of the spectrum to be sought */
73: PEPBasis basis; /* polynomial basis used to represent the problem */
74: PEPProblemType problem_type; /* which kind of problem to be solved */
75: PEPScale scale; /* scaling strategy to be used */
76: PetscReal sfactor,dsfactor; /* scaling factors */
77: PetscInt sits; /* number of iterations of the scaling method */
78: PetscReal slambda; /* norm eigenvalue approximation for scaling */
79: PEPRefine refine; /* type of refinement to be applied after solve */
80: PetscInt npart; /* number of partitions of the communicator */
81: PetscReal rtol; /* tolerance for refinement */
82: PetscInt rits; /* number of iterations of the refinement method */
83: PEPRefineScheme scheme; /* scheme for solving linear systems within refinement */
84: PEPExtract extract; /* type of extraction used */
85: PetscBool trackall; /* whether all the residuals must be computed */
87: /*-------------- User-provided functions and contexts -----------------*/
88: PetscErrorCode (*converged)(PEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
89: PetscErrorCode (*convergeddestroy)(void*);
90: PetscErrorCode (*stopping)(PEP,PetscInt,PetscInt,PetscInt,PetscInt,PEPConvergedReason*,void*);
91: PetscErrorCode (*stoppingdestroy)(void*);
92: void *convergedctx;
93: void *stoppingctx;
94: PetscErrorCode (*monitor[MAXPEPMONITORS])(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
95: PetscErrorCode (*monitordestroy[MAXPEPMONITORS])(void**);
96: void *monitorcontext[MAXPEPMONITORS];
97: PetscInt numbermonitors;
99: /*----------------- Child objects and working data -------------------*/
100: ST st; /* spectral transformation object */
101: DS ds; /* direct solver object */
102: BV V; /* set of basis vectors and computed eigenvectors */
103: RG rg; /* optional region for filtering */
104: SlepcSC sc; /* sorting criterion data */
105: Mat *A; /* coefficient matrices of the polynomial */
106: PetscInt nmat; /* number of matrices */
107: Vec Dl,Dr; /* diagonal matrices for balancing */
108: Vec *IS; /* references to user-provided initial space */
109: PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */
110: PetscReal *errest; /* error estimates */
111: PetscInt *perm; /* permutation for eigenvalue ordering */
112: PetscReal *pbc; /* coefficients defining the polynomial basis */
113: PetscScalar *solvematcoeffs; /* coefficients to compute the matrix to be inverted */
114: PetscInt nwork; /* number of work vectors */
115: Vec *work; /* work vectors */
116: KSP refineksp; /* ksp used in refinement */
117: PetscSubcomm refinesubc; /* context for sub-communicators */
118: void *data; /* placeholder for solver-specific stuff */
120: /* ----------------------- Status variables --------------------------*/
121: PEPStateType state; /* initial -> setup -> solved -> eigenvectors */
122: PetscInt nconv; /* number of converged eigenvalues */
123: PetscInt its; /* number of iterations so far computed */
124: PetscInt n,nloc; /* problem dimensions (global, local) */
125: PetscReal *nrma; /* computed matrix norms */
126: PetscReal nrml[2]; /* computed matrix norms for the linearization */
127: PetscBool sfactor_set; /* flag to indicate the user gave sfactor */
128: PetscBool lineariz; /* current solver is based on linearization */
129: PEPConvergedReason reason;
130: };
132: /*
133: Macros to test valid PEP arguments
134: */
135: #if !defined(PETSC_USE_DEBUG)
137: #define PEPCheckSolved(h,arg) do {} while (0)139: #else
141: #define PEPCheckSolved(h,arg) \142: do { \143: if (h->state<PEP_STATE_SOLVED) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"Must call PEPSolve() first: Parameter #%d",arg); \144: } while (0)146: #endif
148: PETSC_INTERN PetscErrorCode PEPSetDimensions_Default(PEP,PetscInt,PetscInt*,PetscInt*);
149: PETSC_INTERN PetscErrorCode PEPExtractVectors(PEP);
150: PETSC_INTERN PetscErrorCode PEPBackTransform_Default(PEP);
151: PETSC_INTERN PetscErrorCode PEPComputeVectors(PEP);
152: PETSC_INTERN PetscErrorCode PEPComputeVectors_Default(PEP);
153: PETSC_INTERN PetscErrorCode PEPComputeVectors_Indefinite(PEP);
154: PETSC_INTERN PetscErrorCode PEPComputeResidualNorm_Private(PEP,PetscScalar,PetscScalar,Vec,Vec,Vec*,PetscReal*);
155: PETSC_INTERN PetscErrorCode PEPKrylovConvergence(PEP,PetscBool,PetscInt,PetscInt,PetscReal,PetscInt*);
156: PETSC_INTERN PetscErrorCode PEPComputeScaleFactor(PEP);
157: PETSC_INTERN PetscErrorCode PEPBuildDiagonalScaling(PEP);
158: PETSC_INTERN PetscErrorCode PEPBasisCoefficients(PEP,PetscReal*);
159: PETSC_INTERN PetscErrorCode PEPEvaluateBasis(PEP,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*);
160: PETSC_INTERN PetscErrorCode PEPNewtonRefinement_TOAR(PEP,PetscScalar,PetscInt*,PetscReal*,PetscInt,PetscScalar*,PetscInt,PetscInt*);
161: PETSC_INTERN PetscErrorCode PEPNewtonRefinementSimple(PEP,PetscInt*,PetscReal,PetscInt);
163: #endif