Actual source code: qepimpl.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(_QEPIMPL)
23: #define _QEPIMPL
25: #include <slepcqep.h>
26: #include <slepc-private/slepcimpl.h>
28: PETSC_EXTERN PetscLogEvent QEP_SetUp, QEP_Solve, QEP_Dense;
30: typedef struct _QEPOps *QEPOps;
32: struct _QEPOps {
33: PetscErrorCode (*solve)(QEP);
34: PetscErrorCode (*setup)(QEP);
35: PetscErrorCode (*setfromoptions)(QEP);
36: PetscErrorCode (*publishoptions)(QEP);
37: PetscErrorCode (*destroy)(QEP);
38: PetscErrorCode (*reset)(QEP);
39: PetscErrorCode (*view)(QEP,PetscViewer);
40: };
42: /*
43: Maximum number of monitors you can run with a single QEP
44: */
45: #define MAXQEPMONITORS 5
47: /*
48: Defines the QEP data structure.
49: */
50: struct _p_QEP {
51: PETSCHEADER(struct _QEPOps);
52: /*------------------------- User parameters --------------------------*/
53: PetscInt max_it; /* maximum number of iterations */
54: PetscInt nev; /* number of eigenvalues to compute */
55: PetscInt ncv; /* number of basis vectors */
56: PetscInt mpd; /* maximum dimension of projected problem */
57: PetscInt nini,ninil; /* number of initial vectors (negative means not copied yet) */
58: PetscScalar target; /* target value */
59: PetscReal tol; /* tolerance */
60: PetscReal sfactor; /* scaling factor of the quadratic problem */
61: PetscBool sfactor_set; /* flag to indicate the user gave sfactor */
62: QEPWhich which; /* which part of the spectrum to be sought */
63: PetscBool leftvecs; /* if left eigenvectors are requested */
64: QEPProblemType problem_type; /* which kind of problem to be solved */
65: PetscBool trackall; /* whether all the residuals must be computed */
67: /*-------------- User-provided functions and contexts -----------------*/
68: PetscErrorCode (*comparison)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*);
69: PetscErrorCode (*converged)(QEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
70: void *comparisonctx;
71: void *convergedctx;
73: /*------------------------- Working data --------------------------*/
74: Mat M,C,K; /* problem matrices */
75: Vec *V; /* set of basis vectors and computed eigenvectors */
76: Vec *W; /* set of left basis vectors and computed left eigenvectors */
77: Vec *IS,*ISL; /* placeholder for references to user-provided initial space */
78: PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */
79: PetscReal *errest; /* error estimates */
80: IP ip; /* innerproduct object */
81: DS ds; /* direct solver object */
82: ST st; /* spectral transformation object */
83: void *data; /* placeholder for misc stuff associated
84: with a particular solver */
85: PetscInt allocated_ncv; /* number of basis vectors allocated */
86: PetscInt nconv; /* number of converged eigenvalues */
87: PetscInt its; /* number of iterations so far computed */
88: PetscInt *perm; /* permutation for eigenvalue ordering */
89: PetscInt matvecs,linits; /* operation counters */
90: PetscInt n,nloc; /* problem dimensions (global, local) */
91: PetscRandom rand; /* random number generator */
92: Vec t; /* template vector */
94: /* ---------------- Default work-area and status vars -------------------- */
95: PetscInt nwork;
96: Vec *work;
98: PetscInt setupcalled;
99: QEPConvergedReason reason;
101: PetscErrorCode (*monitor[MAXQEPMONITORS])(QEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
102: PetscErrorCode (*monitordestroy[MAXQEPMONITORS])(void**);
103: void *monitorcontext[MAXQEPMONITORS];
104: PetscInt numbermonitors;
105: };
107: PETSC_INTERN PetscErrorCode QEPReset_Default(QEP);
108: PETSC_INTERN PetscErrorCode QEPAllocateSolution(QEP);
109: PETSC_INTERN PetscErrorCode QEPFreeSolution(QEP);
110: PETSC_INTERN PetscErrorCode QEPComputeVectors_Schur(QEP);
111: PETSC_INTERN PetscErrorCode QEPComputeResidualNorm_Private(QEP,PetscScalar,PetscScalar,Vec,Vec,PetscReal*);
112: PETSC_INTERN PetscErrorCode QEPComputeRelativeError_Private(QEP,PetscScalar,PetscScalar,Vec,Vec,PetscReal*);
113: PETSC_INTERN PetscErrorCode QEPKrylovConvergence(QEP,PetscBool,PetscInt,PetscInt,PetscInt,PetscReal,PetscInt*);
115: #endif