Actual source code: bvimpl.h
slepc-3.18.0 2022-10-01
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, 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: /* SUBMANSEC = BV */
19: SLEPC_EXTERN PetscBool BVRegisterAllCalled;
20: SLEPC_EXTERN PetscErrorCode BVRegisterAll(void);
22: 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,BV_SVDAndRank;
24: typedef struct _BVOps *BVOps;
26: struct _BVOps {
27: PetscErrorCode (*mult)(BV,PetscScalar,PetscScalar,BV,Mat);
28: PetscErrorCode (*multvec)(BV,PetscScalar,PetscScalar,Vec,PetscScalar*);
29: PetscErrorCode (*multinplace)(BV,Mat,PetscInt,PetscInt);
30: PetscErrorCode (*multinplacetrans)(BV,Mat,PetscInt,PetscInt);
31: PetscErrorCode (*dot)(BV,BV,Mat);
32: PetscErrorCode (*dotvec)(BV,Vec,PetscScalar*);
33: PetscErrorCode (*dotvec_local)(BV,Vec,PetscScalar*);
34: PetscErrorCode (*dotvec_begin)(BV,Vec,PetscScalar*);
35: PetscErrorCode (*dotvec_end)(BV,Vec,PetscScalar*);
36: PetscErrorCode (*scale)(BV,PetscInt,PetscScalar);
37: PetscErrorCode (*norm)(BV,PetscInt,NormType,PetscReal*);
38: PetscErrorCode (*norm_local)(BV,PetscInt,NormType,PetscReal*);
39: PetscErrorCode (*norm_begin)(BV,PetscInt,NormType,PetscReal*);
40: PetscErrorCode (*norm_end)(BV,PetscInt,NormType,PetscReal*);
41: PetscErrorCode (*normalize)(BV,PetscScalar*);
42: PetscErrorCode (*matmult)(BV,Mat,BV);
43: PetscErrorCode (*copy)(BV,BV);
44: PetscErrorCode (*copycolumn)(BV,PetscInt,PetscInt);
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 (*getarrayread)(BV,const PetscScalar**);
51: PetscErrorCode (*restorearrayread)(BV,const PetscScalar**);
52: PetscErrorCode (*restoresplit)(BV,BV*,BV*);
53: PetscErrorCode (*gramschmidt)(BV,PetscInt,Vec,PetscBool*,PetscScalar*,PetscScalar*,PetscReal*,PetscReal*);
54: PetscErrorCode (*getmat)(BV,Mat*);
55: PetscErrorCode (*restoremat)(BV,Mat*);
56: PetscErrorCode (*duplicate)(BV,BV);
57: PetscErrorCode (*create)(BV);
58: PetscErrorCode (*setfromoptions)(BV,PetscOptionItems*);
59: PetscErrorCode (*view)(BV,PetscViewer);
60: PetscErrorCode (*destroy)(BV);
61: };
63: struct _p_BV {
64: PETSCHEADER(struct _BVOps);
65: /*------------------------- User parameters --------------------------*/
66: Vec t; /* template vector */
67: PetscInt n,N; /* dimensions of vectors (local, global) */
68: PetscInt m; /* number of vectors */
69: PetscInt l; /* number of leading columns */
70: PetscInt k; /* number of active columns */
71: PetscInt nc; /* number of constraints */
72: BVOrthogType orthog_type; /* the method of vector orthogonalization */
73: BVOrthogRefineType orthog_ref; /* refinement method */
74: PetscReal orthog_eta; /* refinement threshold */
75: BVOrthogBlockType orthog_block; /* the method of block orthogonalization */
76: Mat matrix; /* inner product matrix */
77: PetscBool indef; /* matrix is indefinite */
78: BVMatMultType vmm; /* version of matmult operation */
79: PetscBool rrandom; /* reproducible random vectors */
80: PetscReal deftol; /* tolerance for BV_SafeSqrt */
82: /*---------------------- Cached data and workspace -------------------*/
83: Vec buffer; /* buffer vector used in orthogonalization */
84: Mat Abuffer; /* auxiliary seqdense matrix that wraps the buffer */
85: Vec Bx; /* result of matrix times a vector x */
86: PetscObjectId xid; /* object id of vector x */
87: PetscObjectState xstate; /* state of vector x */
88: Vec cv[2]; /* column vectors obtained with BVGetColumn() */
89: PetscInt ci[2]; /* column indices of obtained vectors */
90: PetscObjectState st[2]; /* state of obtained vectors */
91: PetscObjectId id[2]; /* object id of obtained vectors */
92: PetscScalar *h,*c; /* orthogonalization coefficients */
93: Vec omega; /* signature matrix values for indefinite case */
94: PetscBool defersfo; /* deferred call to setfromoptions */
95: BV cached; /* cached BV to store result of matrix times BV */
96: PetscObjectState bvstate; /* state of BV when BVApplyMatrixBV() was called */
97: BV L,R; /* BV objects obtained with BVGetSplit() */
98: PetscObjectState lstate,rstate;/* state of L and R when BVGetSplit() was called */
99: PetscInt lsplit; /* the value of l when BVGetSplit() was called */
100: PetscInt issplit; /* >0 if this BV has been created by splitting (1=left, 2=right) */
101: BV splitparent; /* my parent if I am a split BV */
102: PetscRandom rand; /* random number generator */
103: Mat Acreate; /* matrix given at BVCreateFromMat() */
104: Mat Aget; /* matrix returned for BVGetMat() */
105: PetscBool cuda; /* true if GPU must be used in SVEC */
106: PetscBool sfocalled; /* setfromoptions has been called */
107: PetscScalar *work;
108: PetscInt lwork;
109: void *data;
110: };
112: /*
113: BV_SafeSqrt - Computes the square root of a scalar value alpha, which is
114: assumed to be z'*B*z. The result is
115: if definite inner product: res = sqrt(alpha)
116: if indefinite inner product: res = sgn(alpha)*sqrt(abs(alpha))
117: */
118: static inline PetscErrorCode BV_SafeSqrt(BV bv,PetscScalar alpha,PetscReal *res)
119: {
120: PetscReal absal,realp;
122: absal = PetscAbsScalar(alpha);
123: realp = PetscRealPart(alpha);
124: if (PetscUnlikely(absal<PETSC_MACHINE_EPSILON)) PetscInfo(bv,"Zero norm, either the vector is zero or a semi-inner product is being used\n");
125: #if defined(PETSC_USE_COMPLEX)
127: #endif
128: if (PetscUnlikely(bv->indef)) {
129: *res = (realp<0.0)? -PetscSqrtReal(-realp): PetscSqrtReal(realp);
130: } else {
132: *res = (realp<0.0)? 0.0: PetscSqrtReal(realp);
133: }
134: return 0;
135: }
137: /*
138: BV_IPMatMult - Multiply a vector x by the inner-product matrix, cache the
139: result in Bx.
140: */
141: static inline PetscErrorCode BV_IPMatMult(BV bv,Vec x)
142: {
143: if (((PetscObject)x)->id != bv->xid || ((PetscObject)x)->state != bv->xstate) {
144: if (PetscUnlikely(!bv->Bx)) MatCreateVecs(bv->matrix,&bv->Bx,NULL);
145: MatMult(bv->matrix,x,bv->Bx);
146: PetscObjectGetId((PetscObject)x,&bv->xid);
147: PetscObjectStateGet((PetscObject)x,&bv->xstate);
148: }
149: return 0;
150: }
152: /*
153: BV_IPMatMultBV - Multiply BV by the inner-product matrix, cache the
154: result internally in bv->cached.
155: */
156: static inline PetscErrorCode BV_IPMatMultBV(BV bv)
157: {
158: BVGetCachedBV(bv,&bv->cached);
159: if (((PetscObject)bv)->state != bv->bvstate || bv->l != bv->cached->l || bv->k != bv->cached->k) {
160: BVSetActiveColumns(bv->cached,bv->l,bv->k);
161: if (bv->matrix) BVMatMult(bv,bv->matrix,bv->cached);
162: else BVCopy(bv,bv->cached);
163: bv->bvstate = ((PetscObject)bv)->state;
164: }
165: return 0;
166: }
168: /*
169: BV_AllocateCoeffs - Allocate orthogonalization coefficients if not done already.
170: */
171: static inline PetscErrorCode BV_AllocateCoeffs(BV bv)
172: {
173: if (!bv->h) PetscMalloc2(bv->nc+bv->m,&bv->h,bv->nc+bv->m,&bv->c);
174: return 0;
175: }
177: /*
178: BV_AllocateSignature - Allocate signature coefficients if not done already.
179: */
180: static inline PetscErrorCode BV_AllocateSignature(BV bv)
181: {
182: if (bv->indef && !bv->omega) {
183: if (bv->cuda) {
184: #if defined(PETSC_HAVE_CUDA)
185: VecCreateSeqCUDA(PETSC_COMM_SELF,bv->nc+bv->m,&bv->omega);
186: #else
187: SETERRQ(PetscObjectComm((PetscObject)bv),PETSC_ERR_PLIB,"Something wrong happened");
188: #endif
189: } else VecCreateSeq(PETSC_COMM_SELF,bv->nc+bv->m,&bv->omega);
190: VecSet(bv->omega,1.0);
191: }
192: return 0;
193: }
195: /*
196: BVAvailableVec: First (0) or second (1) vector available for
197: getcolumn operation (or -1 if both vectors already fetched).
198: */
199: #define BVAvailableVec (((bv->ci[0]==-bv->nc-1)? 0: (bv->ci[1]==-bv->nc-1)? 1: -1))
201: /*
202: Macros to test valid BV arguments
203: */
204: #if !defined(PETSC_USE_DEBUG)
206: #define BVCheckSizes(h,arg) do {(void)(h);} while (0)
207: #define BVCheckOp(h,arg,op) do {(void)(h);} while (0)
209: #else
211: #define BVCheckSizes(h,arg) \
212: do { \
214: } while (0)
216: #define BVCheckOp(h,arg,op) \
217: do { \
219: } while (0)
221: #endif
223: SLEPC_INTERN PetscErrorCode BVView_Vecs(BV,PetscViewer);
225: SLEPC_INTERN PetscErrorCode BVAllocateWork_Private(BV,PetscInt);
227: SLEPC_INTERN PetscErrorCode BVMult_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar,const PetscScalar*,const PetscScalar*,PetscScalar,PetscScalar*);
228: SLEPC_INTERN PetscErrorCode BVMultVec_BLAS_Private(BV,PetscInt,PetscInt,PetscScalar,const PetscScalar*,const PetscScalar*,PetscScalar,PetscScalar*);
229: SLEPC_INTERN PetscErrorCode BVMultInPlace_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar*,const PetscScalar*,PetscBool);
230: SLEPC_INTERN PetscErrorCode BVMultInPlace_Vecs_Private(BV,PetscInt,PetscInt,PetscInt,Vec*,const PetscScalar*,PetscBool);
231: SLEPC_INTERN PetscErrorCode BVAXPY_BLAS_Private(BV,PetscInt,PetscInt,PetscScalar,const PetscScalar*,PetscScalar,PetscScalar*);
232: SLEPC_INTERN PetscErrorCode BVDot_BLAS_Private(BV,PetscInt,PetscInt,PetscInt,PetscInt,const PetscScalar*,const PetscScalar*,PetscScalar*,PetscBool);
233: SLEPC_INTERN PetscErrorCode BVDotVec_BLAS_Private(BV,PetscInt,PetscInt,const PetscScalar*,const PetscScalar*,PetscScalar*,PetscBool);
234: SLEPC_INTERN PetscErrorCode BVScale_BLAS_Private(BV,PetscInt,PetscScalar*,PetscScalar);
235: SLEPC_INTERN PetscErrorCode BVNorm_LAPACK_Private(BV,PetscInt,PetscInt,const PetscScalar*,NormType,PetscReal*,PetscBool);
236: SLEPC_INTERN PetscErrorCode BVNormalize_LAPACK_Private(BV,PetscInt,PetscInt,const PetscScalar*,PetscScalar*,PetscBool);
237: SLEPC_INTERN PetscErrorCode BVGetMat_Default(BV,Mat*);
238: SLEPC_INTERN PetscErrorCode BVRestoreMat_Default(BV,Mat*);
239: SLEPC_INTERN PetscErrorCode BVMatCholInv_LAPACK_Private(BV,Mat,Mat);
240: SLEPC_INTERN PetscErrorCode BVMatTriInv_LAPACK_Private(BV,Mat,Mat);
241: SLEPC_INTERN PetscErrorCode BVMatSVQB_LAPACK_Private(BV,Mat,Mat);
242: SLEPC_INTERN PetscErrorCode BVOrthogonalize_LAPACK_TSQR(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscInt);
243: SLEPC_INTERN PetscErrorCode BVOrthogonalize_LAPACK_TSQR_OnlyR(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscInt);
245: /* reduction operations used in BVOrthogonalize and BVNormalize */
246: SLEPC_EXTERN MPI_Op MPIU_TSQR, MPIU_LAPY2;
247: SLEPC_EXTERN void MPIAPI SlepcGivensPacked(void*,void*,PetscMPIInt*,MPI_Datatype*);
248: SLEPC_EXTERN void MPIAPI SlepcPythag(void*,void*,PetscMPIInt*,MPI_Datatype*);
250: /*
251: BV_CleanCoefficients_Default - Sets to zero all entries of column j of the bv buffer
252: */
253: static inline PetscErrorCode BV_CleanCoefficients_Default(BV bv,PetscInt j,PetscScalar *h)
254: {
255: PetscScalar *hh=h,*a;
256: PetscInt i;
258: if (!h) {
259: VecGetArray(bv->buffer,&a);
260: hh = a + j*(bv->nc+bv->m);
261: }
262: for (i=0;i<bv->nc+j;i++) hh[i] = 0.0;
263: if (!h) VecRestoreArray(bv->buffer,&a);
264: return 0;
265: }
267: /*
268: BV_AddCoefficients_Default - Add the contents of the scratch (0-th column) of the bv buffer
269: into column j of the bv buffer
270: */
271: static inline PetscErrorCode BV_AddCoefficients_Default(BV bv,PetscInt j,PetscScalar *h,PetscScalar *c)
272: {
273: PetscScalar *hh=h,*cc=c;
274: PetscInt i;
276: if (!h) {
277: VecGetArray(bv->buffer,&cc);
278: hh = cc + j*(bv->nc+bv->m);
279: }
280: for (i=0;i<bv->nc+j;i++) hh[i] += cc[i];
281: if (!h) VecRestoreArray(bv->buffer,&cc);
282: PetscLogFlops(1.0*(bv->nc+j));
283: return 0;
284: }
286: /*
287: BV_SetValue_Default - Sets value in row j (counted after the constraints) of column k
288: of the coefficients array
289: */
290: static inline PetscErrorCode BV_SetValue_Default(BV bv,PetscInt j,PetscInt k,PetscScalar *h,PetscScalar value)
291: {
292: PetscScalar *hh=h,*a;
294: if (!h) {
295: VecGetArray(bv->buffer,&a);
296: hh = a + k*(bv->nc+bv->m);
297: }
298: hh[bv->nc+j] = value;
299: if (!h) VecRestoreArray(bv->buffer,&a);
300: return 0;
301: }
303: /*
304: BV_SquareSum_Default - Returns the value h'*h, where h represents the contents of the
305: coefficients array (up to position j)
306: */
307: static inline PetscErrorCode BV_SquareSum_Default(BV bv,PetscInt j,PetscScalar *h,PetscReal *sum)
308: {
309: PetscScalar *hh=h;
310: PetscInt i;
312: *sum = 0.0;
313: if (!h) VecGetArray(bv->buffer,&hh);
314: for (i=0;i<bv->nc+j;i++) *sum += PetscRealPart(hh[i]*PetscConj(hh[i]));
315: if (!h) VecRestoreArray(bv->buffer,&hh);
316: PetscLogFlops(2.0*(bv->nc+j));
317: return 0;
318: }
320: /*
321: BV_ApplySignature_Default - Computes the pointwise product h*omega, where h represents
322: the contents of the coefficients array (up to position j) and omega is the signature;
323: if inverse=TRUE then the operation is h/omega
324: */
325: static inline PetscErrorCode BV_ApplySignature_Default(BV bv,PetscInt j,PetscScalar *h,PetscBool inverse)
326: {
327: PetscScalar *hh=h;
328: PetscInt i;
329: const PetscScalar *omega;
331: if (PetscUnlikely(!(bv->nc+j))) return 0;
332: if (!h) VecGetArray(bv->buffer,&hh);
333: VecGetArrayRead(bv->omega,&omega);
334: if (inverse) for (i=0;i<bv->nc+j;i++) hh[i] /= PetscRealPart(omega[i]);
335: else for (i=0;i<bv->nc+j;i++) hh[i] *= PetscRealPart(omega[i]);
336: VecRestoreArrayRead(bv->omega,&omega);
337: if (!h) VecRestoreArray(bv->buffer,&hh);
338: PetscLogFlops(1.0*(bv->nc+j));
339: return 0;
340: }
342: /*
343: BV_SquareRoot_Default - Returns the square root of position j (counted after the constraints)
344: of the coefficients array
345: */
346: static inline PetscErrorCode BV_SquareRoot_Default(BV bv,PetscInt j,PetscScalar *h,PetscReal *beta)
347: {
348: PetscScalar *hh=h;
350: if (!h) VecGetArray(bv->buffer,&hh);
351: BV_SafeSqrt(bv,hh[bv->nc+j],beta);
352: if (!h) VecRestoreArray(bv->buffer,&hh);
353: return 0;
354: }
356: /*
357: BV_StoreCoefficients_Default - Copy the contents of the coefficients array to an array dest
358: provided by the caller (only values from l to j are copied)
359: */
360: static inline PetscErrorCode BV_StoreCoefficients_Default(BV bv,PetscInt j,PetscScalar *h,PetscScalar *dest)
361: {
362: PetscScalar *hh=h,*a;
363: PetscInt i;
365: if (!h) {
366: VecGetArray(bv->buffer,&a);
367: hh = a + j*(bv->nc+bv->m);
368: }
369: for (i=bv->l;i<j;i++) dest[i-bv->l] = hh[bv->nc+i];
370: if (!h) VecRestoreArray(bv->buffer,&a);
371: return 0;
372: }
374: /*
375: BV_GetEigenvector - retrieves k-th eigenvector from basis vectors V.
376: The argument eigi is the imaginary part of the corresponding eigenvalue.
377: */
378: static inline PetscErrorCode BV_GetEigenvector(BV V,PetscInt k,PetscScalar eigi,Vec Vr,Vec Vi)
379: {
380: #if defined(PETSC_USE_COMPLEX)
381: if (Vr) BVCopyVec(V,k,Vr);
382: if (Vi) VecSet(Vi,0.0);
383: #else
384: if (eigi > 0.0) { /* first value of conjugate pair */
385: if (Vr) BVCopyVec(V,k,Vr);
386: if (Vi) BVCopyVec(V,k+1,Vi);
387: } else if (eigi < 0.0) { /* second value of conjugate pair */
388: if (Vr) BVCopyVec(V,k-1,Vr);
389: if (Vi) {
390: BVCopyVec(V,k,Vi);
391: VecScale(Vi,-1.0);
392: }
393: } else { /* real eigenvalue */
394: if (Vr) BVCopyVec(V,k,Vr);
395: if (Vi) VecSet(Vi,0.0);
396: }
397: #endif
398: return 0;
399: }
401: /*
402: BV_OrthogonalizeColumn_Safe - this is intended for cases where we know that
403: the resulting vector is going to be numerically zero, so normalization or
404: iterative refinement may cause problems in parallel (collective operation
405: not being called by all processes)
406: */
407: static inline PetscErrorCode BV_OrthogonalizeColumn_Safe(BV bv,PetscInt j,PetscScalar *H,PetscReal *norm,PetscBool *lindep)
408: {
409: BVOrthogRefineType orthog_ref;
411: PetscInfo(bv,"Orthogonalizing column %" PetscInt_FMT " without refinement\n",j);
412: orthog_ref = bv->orthog_ref;
413: bv->orthog_ref = BV_ORTHOG_REFINE_NEVER; /* avoid refinement */
414: BVOrthogonalizeColumn(bv,j,H,NULL,NULL);
415: bv->orthog_ref = orthog_ref; /* restore refinement setting */
416: if (norm) *norm = 0.0;
417: if (lindep) *lindep = PETSC_TRUE;
418: return 0;
419: }
421: #if defined(PETSC_HAVE_CUDA)
422: SLEPC_INTERN PetscErrorCode BV_CleanCoefficients_CUDA(BV,PetscInt,PetscScalar*);
423: SLEPC_INTERN PetscErrorCode BV_AddCoefficients_CUDA(BV,PetscInt,PetscScalar*,PetscScalar*);
424: SLEPC_INTERN PetscErrorCode BV_SetValue_CUDA(BV,PetscInt,PetscInt,PetscScalar*,PetscScalar);
425: SLEPC_INTERN PetscErrorCode BV_SquareSum_CUDA(BV,PetscInt,PetscScalar*,PetscReal*);
426: SLEPC_INTERN PetscErrorCode BV_ApplySignature_CUDA(BV,PetscInt,PetscScalar*,PetscBool);
427: SLEPC_INTERN PetscErrorCode BV_SquareRoot_CUDA(BV,PetscInt,PetscScalar*,PetscReal*);
428: SLEPC_INTERN PetscErrorCode BV_StoreCoefficients_CUDA(BV,PetscInt,PetscScalar*,PetscScalar*);
429: #define BV_CleanCoefficients(a,b,c) ((a)->cuda?BV_CleanCoefficients_CUDA:BV_CleanCoefficients_Default)((a),(b),(c))
430: #define BV_AddCoefficients(a,b,c,d) ((a)->cuda?BV_AddCoefficients_CUDA:BV_AddCoefficients_Default)((a),(b),(c),(d))
431: #define BV_SetValue(a,b,c,d,e) ((a)->cuda?BV_SetValue_CUDA:BV_SetValue_Default)((a),(b),(c),(d),(e))
432: #define BV_SquareSum(a,b,c,d) ((a)->cuda?BV_SquareSum_CUDA:BV_SquareSum_Default)((a),(b),(c),(d))
433: #define BV_ApplySignature(a,b,c,d) ((a)->cuda?BV_ApplySignature_CUDA:BV_ApplySignature_Default)((a),(b),(c),(d))
434: #define BV_SquareRoot(a,b,c,d) ((a)->cuda?BV_SquareRoot_CUDA:BV_SquareRoot_Default)((a),(b),(c),(d))
435: #define BV_StoreCoefficients(a,b,c,d) ((a)->cuda?BV_StoreCoefficients_CUDA:BV_StoreCoefficients_Default)((a),(b),(c),(d))
436: #else
437: #define BV_CleanCoefficients(a,b,c) BV_CleanCoefficients_Default((a),(b),(c))
438: #define BV_AddCoefficients(a,b,c,d) BV_AddCoefficients_Default((a),(b),(c),(d))
439: #define BV_SetValue(a,b,c,d,e) BV_SetValue_Default((a),(b),(c),(d),(e))
440: #define BV_SquareSum(a,b,c,d) BV_SquareSum_Default((a),(b),(c),(d))
441: #define BV_ApplySignature(a,b,c,d) BV_ApplySignature_Default((a),(b),(c),(d))
442: #define BV_SquareRoot(a,b,c,d) BV_SquareRoot_Default((a),(b),(c),(d))
443: #define BV_StoreCoefficients(a,b,c,d) BV_StoreCoefficients_Default((a),(b),(c),(d))
444: #endif /* PETSC_HAVE_CUDA */
446: #endif