Actual source code: bvimpl.h
slepc-3.15.1 2021-05-28
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2021, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #if !defined(SLEPCBVIMPL_H)
12: #define SLEPCBVIMPL_H
14: #include <slepcbv.h>
15: #include <slepc/private/slepcimpl.h>
17: SLEPC_EXTERN PetscBool BVRegisterAllCalled;
18: SLEPC_EXTERN PetscErrorCode BVRegisterAll(void);
20: SLEPC_EXTERN PetscLogEvent BV_Create,BV_Copy,BV_Mult,BV_MultVec,BV_MultInPlace,BV_Dot,BV_DotVec,BV_Orthogonalize,BV_OrthogonalizeVec,BV_Scale,BV_Norm,BV_NormVec,BV_Normalize,BV_SetRandom,BV_MatMult,BV_MatMultVec,BV_MatProject;
22: typedef struct _BVOps *BVOps;
24: struct _BVOps {
25: PetscErrorCode (*mult)(BV,PetscScalar,PetscScalar,BV,Mat);
26: PetscErrorCode (*multvec)(BV,PetscScalar,PetscScalar,Vec,PetscScalar*);
27: PetscErrorCode (*multinplace)(BV,Mat,PetscInt,PetscInt);
28: PetscErrorCode (*multinplacetrans)(BV,Mat,PetscInt,PetscInt);
29: PetscErrorCode (*dot)(BV,BV,Mat);
30: PetscErrorCode (*dotvec)(BV,Vec,PetscScalar*);
31: PetscErrorCode (*dotvec_local)(BV,Vec,PetscScalar*);
32: PetscErrorCode (*dotvec_begin)(BV,Vec,PetscScalar*);
33: PetscErrorCode (*dotvec_end)(BV,Vec,PetscScalar*);
34: PetscErrorCode (*scale)(BV,PetscInt,PetscScalar);
35: PetscErrorCode (*norm)(BV,PetscInt,NormType,PetscReal*);
36: PetscErrorCode (*norm_local)(BV,PetscInt,NormType,PetscReal*);
37: PetscErrorCode (*norm_begin)(BV,PetscInt,NormType,PetscReal*);
38: PetscErrorCode (*norm_end)(BV,PetscInt,NormType,PetscReal*);
39: PetscErrorCode (*normalize)(BV,PetscScalar*);
40: PetscErrorCode (*matmult)(BV,Mat,BV);
41: PetscErrorCode (*copy)(BV,BV);
42: PetscErrorCode (*copycolumn)(BV,PetscInt,PetscInt);
43: PetscErrorCode (*resize)(BV,PetscInt,PetscBool);
44: PetscErrorCode (*getcolumn)(BV,PetscInt,Vec*);
45: PetscErrorCode (*restorecolumn)(BV,PetscInt,Vec*);
46: PetscErrorCode (*getarray)(BV,PetscScalar**);
47: PetscErrorCode (*restorearray)(BV,PetscScalar**);
48: PetscErrorCode (*getarrayread)(BV,const PetscScalar**);
49: PetscErrorCode (*restorearrayread)(BV,const PetscScalar**);
50: PetscErrorCode (*restoresplit)(BV,BV*,BV*);
51: PetscErrorCode (*gramschmidt)(BV,PetscInt,Vec,PetscBool*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*);
52: PetscErrorCode (*getmat)(BV,Mat*);
53: PetscErrorCode (*restoremat)(BV,Mat*);
54: PetscErrorCode (*duplicate)(BV,BV);
55: PetscErrorCode (*create)(BV);
56: PetscErrorCode (*setfromoptions)(PetscOptionItems*,BV);
57: PetscErrorCode (*view)(BV,PetscViewer);
58: PetscErrorCode (*destroy)(BV);
59: };
61: struct _p_BV {
62: PETSCHEADER(struct _BVOps);
63: /*------------------------- User parameters --------------------------*/
64: Vec t; /* template vector */
65: PetscInt n,N; /* dimensions of vectors (local, global) */
66: PetscInt m; /* number of vectors */
67: PetscInt l; /* number of leading columns */
68: PetscInt k; /* number of active columns */
69: PetscInt nc; /* number of constraints */
70: BVOrthogType orthog_type; /* the method of vector orthogonalization */
71: BVOrthogRefineType orthog_ref; /* refinement method */
72: PetscReal orthog_eta; /* refinement threshold */
73: BVOrthogBlockType orthog_block; /* the method of block orthogonalization */
74: Mat matrix; /* inner product matrix */
75: PetscBool indef; /* matrix is indefinite */
76: BVMatMultType vmm; /* version of matmult operation */
77: PetscBool rrandom; /* reproducible random vectors */
78: PetscReal deftol; /* tolerance for BV_SafeSqrt */
80: /*---------------------- Cached data and workspace -------------------*/
81: Vec buffer; /* buffer vector used in orthogonalization */
82: Mat Abuffer; /* auxiliary seqdense matrix that wraps the buffer */
83: Vec Bx; /* result of matrix times a vector x */
84: PetscObjectId xid; /* object id of vector x */
85: PetscObjectState xstate; /* state of vector x */
86: Vec cv[2]; /* column vectors obtained with BVGetColumn() */
87: PetscInt ci[2]; /* column indices of obtained vectors */
88: PetscObjectState st[2]; /* state of obtained vectors */
89: PetscObjectId id[2]; /* object id of obtained vectors */
90: PetscScalar *h,*c; /* orthogonalization coefficients */
91: Vec omega; /* signature matrix values for indefinite case */
92: PetscBool defersfo; /* deferred call to setfromoptions */
93: BV cached; /* cached BV to store result of matrix times BV */
94: PetscObjectState bvstate; /* state of BV when BVApplyMatrixBV() was called */
95: BV L,R; /* BV objects obtained with BVGetSplit() */
96: PetscObjectState lstate,rstate;/* state of L and R when BVGetSplit() was called */
97: PetscInt lsplit; /* the value of l when BVGetSplit() was called */
98: PetscInt issplit; /* >0 if this BV has been created by splitting (1=left, 2=right) */
99: BV splitparent; /* my parent if I am a split BV */
100: PetscRandom rand; /* random number generator */
101: Mat Acreate; /* matrix given at BVCreateFromMat() */
102: Mat Aget; /* matrix returned for BVGetMat() */
103: PetscBool cuda; /* true if GPU must be used in SVEC */
104: PetscBool sfocalled; /* setfromoptions has been called */
105: PetscScalar *work;
106: PetscInt lwork;
107: void *data;
108: };
110: /*
111: BV_SafeSqrt - Computes the square root of a scalar value alpha, which is
112: assumed to be z'*B*z. The result is
113: if definite inner product: res = sqrt(alpha)
114: if indefinite inner product: res = sgn(alpha)*sqrt(abs(alpha))
115: */
116: PETSC_STATIC_INLINE PetscErrorCode BV_SafeSqrt(BV bv,PetscScalar alpha,PetscReal *res)
117: {
119: PetscReal absal,realp;
122: absal = PetscAbsScalar(alpha);
123: realp = PetscRealPart(alpha);
124: if (absal<PETSC_MACHINE_EPSILON) {
125: PetscInfo(bv,"Zero norm, either the vector is zero or a semi-inner product is being used\n");
126: }
127: #if defined(PETSC_USE_COMPLEX)
128: if (PetscAbsReal(PetscImaginaryPart(alpha))>bv->deftol && PetscAbsReal(PetscImaginaryPart(alpha))/absal>10*bv->deftol) SETERRQ1(PetscObjectComm((PetscObject)bv),1,"The inner product is not well defined: nonzero imaginary part %g",PetscImaginaryPart(alpha));
129: #endif
130: if (bv->indef) {
131: *res = (realp<0.0)? -PetscSqrtReal(-realp): PetscSqrtReal(realp);
132: } else {
133: if (realp<-bv->deftol) SETERRQ(PetscObjectComm((PetscObject)bv),1,"The inner product is not well defined: indefinite matrix");
134: *res = (realp<0.0)? 0.0: PetscSqrtReal(realp);
135: }
136: return(0);
137: }
139: /*
140: BV_IPMatMult - Multiply a vector x by the inner-product matrix, cache the
141: result in Bx.
142: */
143: PETSC_STATIC_INLINE PetscErrorCode BV_IPMatMult(BV bv,Vec x)
144: {
148: if (((PetscObject)x)->id != bv->xid || ((PetscObject)x)->state != bv->xstate) {
149: if (!bv->Bx) {
150: MatCreateVecs(bv->matrix,&bv->Bx,NULL);
151: PetscLogObjectParent((PetscObject)bv,(PetscObject)bv->Bx);
152: }
153: MatMult(bv->matrix,x,bv->Bx);
154: PetscObjectGetId((PetscObject)x,&bv->xid);
155: PetscObjectStateGet((PetscObject)x,&bv->xstate);
156: }
157: return(0);
158: }
160: /*
161: BV_IPMatMultBV - Multiply BV by the inner-product matrix, cache the
162: result internally in bv->cached.
163: */
164: PETSC_STATIC_INLINE PetscErrorCode BV_IPMatMultBV(BV bv)
165: {
169: BVGetCachedBV(bv,&bv->cached);
170: if (((PetscObject)bv)->state != bv->bvstate || bv->l != bv->cached->l || bv->k != bv->cached->k) {
171: BVSetActiveColumns(bv->cached,bv->l,bv->k);
172: if (bv->matrix) {
173: BVMatMult(bv,bv->matrix,bv->cached);
174: } else {
175: BVCopy(bv,bv->cached);
176: }
177: bv->bvstate = ((PetscObject)bv)->state;
178: }
179: return(0);
180: }
182: /*
183: BV_AllocateCoeffs - Allocate orthogonalization coefficients if not done already.
184: */
185: PETSC_STATIC_INLINE PetscErrorCode BV_AllocateCoeffs(BV bv)
186: {
190: if (!bv->h) {
191: PetscMalloc2(bv->nc+bv->m,&bv->h,bv->nc+bv->m,&bv->c);
192: PetscLogObjectMemory((PetscObject)bv,2*bv->m*sizeof(PetscScalar));
193: }
194: return(0);
195: }
197: /*
198: BV_AllocateSignature - Allocate signature coefficients if not done already.
199: */
200: PETSC_STATIC_INLINE PetscErrorCode BV_AllocateSignature(BV bv)
201: {
205: if (bv->indef && !bv->omega) {
206: if (bv->cuda) {
207: #if defined(PETSC_HAVE_CUDA)
208: VecCreateSeqCUDA(PETSC_COMM_SELF,bv->nc+bv->m,&bv->omega);
209: #else
210: SETERRQ(PetscObjectComm((PetscObject)bv),1,"Something wrong happened");
211: #endif
212: } else {
213: VecCreateSeq(PETSC_COMM_SELF,bv->nc+bv->m,&bv->omega);
214: }
215: PetscLogObjectParent((PetscObject)bv,(PetscObject)bv->omega);
216: VecSet(bv->omega,1.0);
217: }
218: return(0);
219: }
221: /*
222: BVAvailableVec: First (0) or second (1) vector available for
223: getcolumn operation (or -1 if both vectors already fetched).
224: */
225: #define BVAvailableVec (((bv->ci[0]==-bv->nc-1)? 0: (bv->ci[1]==-bv->nc-1)? 1: -1))
227: /*
228: Macros to test valid BV arguments
229: */
230: #if !defined(PETSC_USE_DEBUG)
232: #define BVCheckSizes(h,arg) do {(void)(h);} while (0)
233: #define BVCheckOp(h,arg,op) do {(void)(h);} while (0)
235: #else
237: #define BVCheckSizes(h,arg) \
238: do { \
239: if (!(h)->m) SETERRQ1(PetscObjectComm((PetscObject)(h)),PETSC_ERR_ARG_WRONGSTATE,"BV sizes have not been defined: Parameter #%d",arg); \
240: } while (0)
242: #define BVCheckOp(h,arg,op) \
243: do { \
244: if (!(h)->ops->op) SETERRQ1(PetscObjectComm((PetscObject)(h)),PETSC_ERR_SUP,"Operation not implemented in this BV type: Parameter #%d",arg); \
245: } while (0)
247: #endif
249: SLEPC_INTERN PetscErrorCode BVView_Vecs(BV,PetscViewer);
251: SLEPC_INTERN PetscErrorCode BVAllocateWork_Private(BV,PetscInt);
253: SLEPC_INTERN PetscErrorCode BVMult_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar,const PetscScalar*,const PetscScalar*,PetscScalar,PetscScalar*);
254: SLEPC_INTERN PetscErrorCode BVMultVec_BLAS_Private(BV,PetscInt,PetscInt,PetscScalar,const PetscScalar*,const PetscScalar*,PetscScalar,PetscScalar*);
255: SLEPC_INTERN PetscErrorCode BVMultInPlace_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar*,const PetscScalar*,PetscBool);
256: SLEPC_INTERN PetscErrorCode BVMultInPlace_Vecs_Private(BV,PetscInt,PetscInt,PetscInt,Vec*,const PetscScalar*,PetscBool);
257: SLEPC_INTERN PetscErrorCode BVAXPY_BLAS_Private(BV,PetscInt,PetscInt,PetscScalar,const PetscScalar*,PetscScalar,PetscScalar*);
258: SLEPC_INTERN PetscErrorCode BVDot_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,const PetscScalar*,const PetscScalar*,PetscScalar*,PetscBool);
259: SLEPC_INTERN PetscErrorCode BVDotVec_BLAS_Private(BV,PetscInt,PetscInt,const PetscScalar*,const PetscScalar*,PetscScalar*,PetscBool);
260: SLEPC_INTERN PetscErrorCode BVScale_BLAS_Private(BV,PetscInt,PetscScalar*,PetscScalar);
261: SLEPC_INTERN PetscErrorCode BVNorm_LAPACK_Private(BV,PetscInt,PetscInt,const PetscScalar*,NormType,PetscReal*,PetscBool);
262: SLEPC_INTERN PetscErrorCode BVNormalize_LAPACK_Private(BV,PetscInt,PetscInt,const PetscScalar*,PetscScalar*,PetscBool);
263: SLEPC_INTERN PetscErrorCode BVMatCholInv_LAPACK_Private(BV,Mat,Mat);
264: SLEPC_INTERN PetscErrorCode BVMatTriInv_LAPACK_Private(BV,Mat,Mat);
265: SLEPC_INTERN PetscErrorCode BVMatSVQB_LAPACK_Private(BV,Mat,Mat);
266: SLEPC_INTERN PetscErrorCode BVOrthogonalize_LAPACK_TSQR(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscInt);
267: SLEPC_INTERN PetscErrorCode BVOrthogonalize_LAPACK_TSQR_OnlyR(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscInt);
269: /* reduction operations used in BVOrthogonalize and BVNormalize */
270: SLEPC_EXTERN MPI_Op MPIU_TSQR, MPIU_LAPY2;
271: SLEPC_EXTERN void MPIAPI SlepcGivensPacked(void*,void*,PetscMPIInt*,MPI_Datatype*);
272: SLEPC_EXTERN void MPIAPI SlepcPythag(void*,void*,PetscMPIInt*,MPI_Datatype*);
274: /*
275: BV_CleanCoefficients_Default - Sets to zero all entries of column j of the bv buffer
276: */
277: PETSC_STATIC_INLINE PetscErrorCode BV_CleanCoefficients_Default(BV bv,PetscInt j,PetscScalar *h)
278: {
280: PetscScalar *hh=h,*a;
281: PetscInt i;
284: if (!h) {
285: VecGetArray(bv->buffer,&a);
286: hh = a + j*(bv->nc+bv->m);
287: }
288: for (i=0;i<bv->nc+j;i++) hh[i] = 0.0;
289: if (!h) { VecRestoreArray(bv->buffer,&a); }
290: return(0);
291: }
293: /*
294: BV_AddCoefficients_Default - Add the contents of the scratch (0-th column) of the bv buffer
295: into column j of the bv buffer
296: */
297: PETSC_STATIC_INLINE PetscErrorCode BV_AddCoefficients_Default(BV bv,PetscInt j,PetscScalar *h,PetscScalar *c)
298: {
300: PetscScalar *hh=h,*cc=c;
301: PetscInt i;
304: if (!h) {
305: VecGetArray(bv->buffer,&cc);
306: hh = cc + j*(bv->nc+bv->m);
307: }
308: for (i=0;i<bv->nc+j;i++) hh[i] += cc[i];
309: if (!h) { VecRestoreArray(bv->buffer,&cc); }
310: return(0);
311: }
313: /*
314: BV_SetValue_Default - Sets value in row j (counted after the constraints) of column k
315: of the coefficients array
316: */
317: PETSC_STATIC_INLINE PetscErrorCode BV_SetValue_Default(BV bv,PetscInt j,PetscInt k,PetscScalar *h,PetscScalar value)
318: {
320: PetscScalar *hh=h,*a;
323: if (!h) {
324: VecGetArray(bv->buffer,&a);
325: hh = a + k*(bv->nc+bv->m);
326: }
327: hh[bv->nc+j] = value;
328: if (!h) { VecRestoreArray(bv->buffer,&a); }
329: return(0);
330: }
332: /*
333: BV_SquareSum_Default - Returns the value h'*h, where h represents the contents of the
334: coefficients array (up to position j)
335: */
336: PETSC_STATIC_INLINE PetscErrorCode BV_SquareSum_Default(BV bv,PetscInt j,PetscScalar *h,PetscReal *sum)
337: {
339: PetscScalar *hh=h;
340: PetscInt i;
343: *sum = 0.0;
344: if (!h) { VecGetArray(bv->buffer,&hh); }
345: for (i=0;i<bv->nc+j;i++) *sum += PetscRealPart(hh[i]*PetscConj(hh[i]));
346: if (!h) { VecRestoreArray(bv->buffer,&hh); }
347: return(0);
348: }
350: /*
351: BV_ApplySignature_Default - Computes the pointwise product h*omega, where h represents
352: the contents of the coefficients array (up to position j) and omega is the signature;
353: if inverse=TRUE then the operation is h/omega
354: */
355: PETSC_STATIC_INLINE PetscErrorCode BV_ApplySignature_Default(BV bv,PetscInt j,PetscScalar *h,PetscBool inverse)
356: {
357: PetscErrorCode ierr;
358: PetscScalar *hh=h;
359: PetscInt i;
360: const PetscScalar *omega;
363: if (!(bv->nc+j)) return(0);
364: if (!h) { VecGetArray(bv->buffer,&hh); }
365: VecGetArrayRead(bv->omega,&omega);
366: if (inverse) for (i=0;i<bv->nc+j;i++) hh[i] /= PetscRealPart(omega[i]);
367: else for (i=0;i<bv->nc+j;i++) hh[i] *= PetscRealPart(omega[i]);
368: VecRestoreArrayRead(bv->omega,&omega);
369: if (!h) { VecRestoreArray(bv->buffer,&hh); }
370: return(0);
371: }
373: /*
374: BV_SquareRoot_Default - Returns the square root of position j (counted after the constraints)
375: of the coefficients array
376: */
377: PETSC_STATIC_INLINE PetscErrorCode BV_SquareRoot_Default(BV bv,PetscInt j,PetscScalar *h,PetscReal *beta)
378: {
380: PetscScalar *hh=h;
383: if (!h) { VecGetArray(bv->buffer,&hh); }
384: BV_SafeSqrt(bv,hh[bv->nc+j],beta);
385: if (!h) { VecRestoreArray(bv->buffer,&hh); }
386: return(0);
387: }
389: /*
390: BV_StoreCoefficients_Default - Copy the contents of the coefficients array to an array dest
391: provided by the caller (only values from l to j are copied)
392: */
393: PETSC_STATIC_INLINE PetscErrorCode BV_StoreCoefficients_Default(BV bv,PetscInt j,PetscScalar *h,PetscScalar *dest)
394: {
396: PetscScalar *hh=h,*a;
397: PetscInt i;
400: if (!h) {
401: VecGetArray(bv->buffer,&a);
402: hh = a + j*(bv->nc+bv->m);
403: }
404: for (i=bv->l;i<j;i++) dest[i-bv->l] = hh[bv->nc+i];
405: if (!h) { VecRestoreArray(bv->buffer,&a); }
406: return(0);
407: }
409: /*
410: BV_GetEigenvector - retrieves k-th eigenvector from basis vectors V.
411: The argument eigi is the imaginary part of the corresponding eigenvalue.
412: */
413: PETSC_STATIC_INLINE PetscErrorCode BV_GetEigenvector(BV V,PetscInt k,PetscScalar eigi,Vec Vr,Vec Vi)
414: {
418: #if defined(PETSC_USE_COMPLEX)
419: if (Vr) { BVCopyVec(V,k,Vr); }
420: if (Vi) { VecSet(Vi,0.0); }
421: #else
422: if (eigi > 0.0) { /* first value of conjugate pair */
423: if (Vr) { BVCopyVec(V,k,Vr); }
424: if (Vi) { BVCopyVec(V,k+1,Vi); }
425: } else if (eigi < 0.0) { /* second value of conjugate pair */
426: if (Vr) { BVCopyVec(V,k-1,Vr); }
427: if (Vi) {
428: BVCopyVec(V,k,Vi);
429: VecScale(Vi,-1.0);
430: }
431: } else { /* real eigenvalue */
432: if (Vr) { BVCopyVec(V,k,Vr); }
433: if (Vi) { VecSet(Vi,0.0); }
434: }
435: #endif
436: return(0);
437: }
439: /*
440: BV_OrthogonalizeColumn_Safe - this is intended for cases where we know that
441: the resulting vector is going to be numerically zero, so normalization or
442: iterative refinement may cause problems in parallel (collective operation
443: not being called by all processes)
444: */
445: PETSC_STATIC_INLINE PetscErrorCode BV_OrthogonalizeColumn_Safe(BV bv,PetscInt j,PetscScalar *H,PetscReal *norm,PetscBool *lindep)
446: {
447: PetscErrorCode ierr;
448: BVOrthogRefineType orthog_ref;
451: PetscInfo1(bv,"Orthogonalizing column %D without refinement\n",j);
452: orthog_ref = bv->orthog_ref;
453: bv->orthog_ref = BV_ORTHOG_REFINE_NEVER; /* avoid refinement */
454: BVOrthogonalizeColumn(bv,j,H,NULL,NULL);
455: bv->orthog_ref = orthog_ref; /* restore refinement setting */
456: if (norm) *norm = 0.0;
457: if (lindep) *lindep = PETSC_TRUE;
458: return(0);
459: }
461: #if defined(PETSC_HAVE_CUDA)
462: SLEPC_INTERN PetscErrorCode BV_CleanCoefficients_CUDA(BV,PetscInt,PetscScalar*);
463: SLEPC_INTERN PetscErrorCode BV_AddCoefficients_CUDA(BV,PetscInt,PetscScalar*,PetscScalar*);
464: SLEPC_INTERN PetscErrorCode BV_SetValue_CUDA(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar);
465: SLEPC_INTERN PetscErrorCode BV_SquareSum_CUDA(BV,PetscInt,PetscScalar*,PetscReal*);
466: SLEPC_INTERN PetscErrorCode BV_ApplySignature_CUDA(BV,PetscInt,PetscScalar*,PetscBool);
467: SLEPC_INTERN PetscErrorCode BV_SquareRoot_CUDA(BV,PetscInt,PetscScalar*,PetscReal*);
468: SLEPC_INTERN PetscErrorCode BV_StoreCoefficients_CUDA(BV,PetscInt,PetscScalar*,PetscScalar*);
469: #define BV_CleanCoefficients(a,b,c) ((a)->cuda?BV_CleanCoefficients_CUDA:BV_CleanCoefficients_Default)((a),(b),(c))
470: #define BV_AddCoefficients(a,b,c,d) ((a)->cuda?BV_AddCoefficients_CUDA:BV_AddCoefficients_Default)((a),(b),(c),(d))
471: #define BV_SetValue(a,b,c,d,e) ((a)->cuda?BV_SetValue_CUDA:BV_SetValue_Default)((a),(b),(c),(d),(e))
472: #define BV_SquareSum(a,b,c,d) ((a)->cuda?BV_SquareSum_CUDA:BV_SquareSum_Default)((a),(b),(c),(d))
473: #define BV_ApplySignature(a,b,c,d) ((a)->cuda?BV_ApplySignature_CUDA:BV_ApplySignature_Default)((a),(b),(c),(d))
474: #define BV_SquareRoot(a,b,c,d) ((a)->cuda?BV_SquareRoot_CUDA:BV_SquareRoot_Default)((a),(b),(c),(d))
475: #define BV_StoreCoefficients(a,b,c,d) ((a)->cuda?BV_StoreCoefficients_CUDA:BV_StoreCoefficients_Default)((a),(b),(c),(d))
476: #else
477: #define BV_CleanCoefficients(a,b,c) BV_CleanCoefficients_Default((a),(b),(c))
478: #define BV_AddCoefficients(a,b,c,d) BV_AddCoefficients_Default((a),(b),(c),(d))
479: #define BV_SetValue(a,b,c,d,e) BV_SetValue_Default((a),(b),(c),(d),(e))
480: #define BV_SquareSum(a,b,c,d) BV_SquareSum_Default((a),(b),(c),(d))
481: #define BV_ApplySignature(a,b,c,d) BV_ApplySignature_Default((a),(b),(c),(d))
482: #define BV_SquareRoot(a,b,c,d) BV_SquareRoot_Default((a),(b),(c),(d))
483: #define BV_StoreCoefficients(a,b,c,d) BV_StoreCoefficients_Default((a),(b),(c),(d))
484: #endif /* PETSC_HAVE_CUDA */
486: #endif