Actual source code: bvimpl.h
slepc-3.5.2 2014-10-10
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) , 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(_BVIMPL)
23: #define _BVIMPL
25: #include <slepcbv.h>
26: #include <slepc-private/slepcimpl.h>
28: PETSC_EXTERN PetscLogEvent BV_Create,BV_Copy,BV_Mult,BV_Dot,BV_Orthogonalize,BV_Scale,BV_Norm,BV_SetRandom,BV_MatMult,BV_MatProject,BV_AXPY;
30: typedef struct _BVOps *BVOps;
32: struct _BVOps {
33: PetscErrorCode (*mult)(BV,PetscScalar,PetscScalar,BV,Mat);
34: PetscErrorCode (*multvec)(BV,PetscScalar,PetscScalar,Vec,PetscScalar*);
35: PetscErrorCode (*multinplace)(BV,Mat,PetscInt,PetscInt);
36: PetscErrorCode (*multinplacetrans)(BV,Mat,PetscInt,PetscInt);
37: PetscErrorCode (*axpy)(BV,PetscScalar,BV);
38: PetscErrorCode (*dot)(BV,BV,Mat);
39: PetscErrorCode (*dotvec)(BV,Vec,PetscScalar*);
40: PetscErrorCode (*scale)(BV,PetscInt,PetscScalar);
41: PetscErrorCode (*norm)(BV,PetscInt,NormType,PetscReal*);
42: PetscErrorCode (*orthogonalize)(BV,Mat);
43: PetscErrorCode (*matmult)(BV,Mat,BV);
44: PetscErrorCode (*copy)(BV,BV);
45: PetscErrorCode (*resize)(BV,PetscInt,PetscBool);
46: PetscErrorCode (*getcolumn)(BV,PetscInt,Vec*);
47: PetscErrorCode (*restorecolumn)(BV,PetscInt,Vec*);
48: PetscErrorCode (*getarray)(BV,PetscScalar**);
49: PetscErrorCode (*restorearray)(BV,PetscScalar**);
50: PetscErrorCode (*setfromoptions)(BV);
51: PetscErrorCode (*create)(BV);
52: PetscErrorCode (*view)(BV,PetscViewer);
53: PetscErrorCode (*destroy)(BV);
54: };
56: struct _p_BV {
57: PETSCHEADER(struct _BVOps);
58: /*------------------------- User parameters --------------------------*/
59: Vec t; /* template vector */
60: PetscInt n,N; /* dimensions of vectors (local, global) */
61: PetscInt m; /* number of vectors */
62: PetscInt l; /* number of leading columns */
63: PetscInt k; /* number of active columns */
64: PetscInt nc; /* number of constraints */
65: BVOrthogType orthog_type; /* which orthogonalization to use */
66: BVOrthogRefineType orthog_ref; /* refinement method */
67: PetscReal orthog_eta; /* refinement threshold */
68: Mat matrix; /* inner product matrix */
69: PetscBool indef; /* matrix is indefinite */
71: /*---------------------- Cached data and workspace -------------------*/
72: Vec Bx; /* result of matrix times a vector x */
73: PetscInt xid; /* object id of vector x */
74: PetscInt xstate; /* state of vector x */
75: Vec cv[2]; /* column vectors obtained with BVGetColumn() */
76: PetscInt ci[2]; /* column indices of obtained vectors */
77: PetscObjectState st[2]; /* state of obtained vectors */
78: PetscObjectId id[2]; /* object id of obtained vectors */
79: PetscScalar *h,*c; /* orthogonalization coefficients */
80: PetscReal *omega; /* signature matrix values for indefinite case */
81: PetscScalar *work;
82: PetscInt lwork;
83: void *data;
84: };
88: /*
89: BV_IPMatMult - Multiply a vector x by the inner-product matrix, cache the
90: result in Bx.
91: */
92: PETSC_STATIC_INLINE PetscErrorCode BV_IPMatMult(BV bv,Vec x)
93: {
97: if (((PetscObject)x)->id != bv->xid || ((PetscObject)x)->state != bv->xstate) {
98: MatMult(bv->matrix,x,bv->Bx);
99: bv->xid = ((PetscObject)x)->id;
100: bv->xstate = ((PetscObject)x)->state;
101: }
102: return(0);
103: }
107: /*
108: BV_AllocateCoeffs - Allocate orthogonalization coefficients if not done already.
109: */
110: PETSC_STATIC_INLINE PetscErrorCode BV_AllocateCoeffs(BV bv)
111: {
115: if (!bv->h) {
116: PetscMalloc2(bv->nc+bv->m,&bv->h,bv->nc+bv->m,&bv->c);
117: PetscLogObjectMemory((PetscObject)bv,2*bv->m*sizeof(PetscScalar));
118: }
119: return(0);
120: }
124: /*
125: BV_AllocateSignature - Allocate signature coefficients if not done already.
126: */
127: PETSC_STATIC_INLINE PetscErrorCode BV_AllocateSignature(BV bv)
128: {
130: PetscInt i;
133: if (bv->indef && !bv->omega) {
134: PetscMalloc1(bv->nc+bv->m,&bv->omega);
135: PetscLogObjectMemory((PetscObject)bv,bv->m*sizeof(PetscReal));
136: for (i=-bv->nc;i<bv->m;i++) bv->omega[i] = 1.0;
137: }
138: return(0);
139: }
141: /*
142: BVAvailableVec: First (0) or second (1) vector available for
143: getcolumn operation (or -1 if both vectors already fetched).
144: */
145: #define BVAvailableVec (((bv->ci[0]==-bv->nc-1)? 0: (bv->ci[1]==-bv->nc-1)? 1: -1))
147: /*
148: Macros to test valid BV arguments
149: */
150: #if !defined(PETSC_USE_DEBUG)
152: #define BVCheckSizes(h,arg) do {} while (0)
154: #else
156: #define BVCheckSizes(h,arg) \
157: do { \
158: if (!h->m) SETERRQ1(PetscObjectComm((PetscObject)h),PETSC_ERR_ARG_WRONGSTATE,"BV sizes have not been defined: Parameter #%d",arg); \
159: } while (0)
161: #endif
163: PETSC_INTERN PetscErrorCode BVView_Vecs(BV,PetscViewer);
165: PETSC_INTERN PetscErrorCode BVAllocateWork_Private(BV,PetscInt);
167: PETSC_INTERN PetscErrorCode BVMult_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar,PetscScalar*,PetscScalar*,PetscScalar,PetscScalar*);
168: PETSC_INTERN PetscErrorCode BVMultVec_BLAS_Private(BV,PetscInt,PetscInt,PetscScalar,PetscScalar*,PetscScalar*,PetscScalar,PetscScalar*);
169: PETSC_INTERN PetscErrorCode BVMultInPlace_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscBool);
170: PETSC_INTERN PetscErrorCode BVMultInPlace_Vecs_Private(BV,PetscInt,PetscInt,PetscInt,Vec*,PetscScalar*,PetscBool);
171: PETSC_INTERN PetscErrorCode BVAXPY_BLAS_Private(BV,PetscInt,PetscInt,PetscScalar,PetscScalar*,PetscScalar*);
172: PETSC_INTERN PetscErrorCode BVDot_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscScalar*,PetscBool);
173: PETSC_INTERN PetscErrorCode BVDotVec_BLAS_Private(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscScalar*,PetscBool);
174: PETSC_INTERN PetscErrorCode BVScale_BLAS_Private(BV,PetscInt,PetscScalar*,PetscScalar);
175: PETSC_INTERN PetscErrorCode BVNorm_LAPACK_Private(BV,PetscInt,PetscInt,PetscScalar*,NormType,PetscReal*,PetscBool);
176: PETSC_INTERN PetscErrorCode BVOrthogonalize_LAPACK_Private(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscBool);
178: #endif