Actual source code: bvmat.c

slepc-3.20.1 2023-11-27
Report Typos and Errors
  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: */
 10: /*
 11:    BV implemented with a dense Mat
 12: */

 14: #include <slepc/private/bvimpl.h>
 15: #include "bvmat.h"

 17: static PetscErrorCode BVMult_Mat(BV Y,PetscScalar alpha,PetscScalar beta,BV X,Mat Q)
 18: {
 19:   BV_MAT            *y = (BV_MAT*)Y->data,*x = (BV_MAT*)X->data;
 20:   PetscScalar       *py;
 21:   const PetscScalar *px,*q;
 22:   PetscInt          ldq;

 24:   PetscFunctionBegin;
 25:   PetscCall(MatDenseGetArrayRead(x->A,&px));
 26:   PetscCall(MatDenseGetArray(y->A,&py));
 27:   if (Q) {
 28:     PetscCall(MatDenseGetLDA(Q,&ldq));
 29:     PetscCall(MatDenseGetArrayRead(Q,&q));
 30:     PetscCall(BVMult_BLAS_Private(Y,Y->n,Y->k-Y->l,X->k-X->l,alpha,px+(X->nc+X->l)*X->ld,X->ld,q+Y->l*ldq+X->l,ldq,beta,py+(Y->nc+Y->l)*Y->ld,Y->ld));
 31:     PetscCall(MatDenseRestoreArrayRead(Q,&q));
 32:   } else PetscCall(BVAXPY_BLAS_Private(Y,Y->n,Y->k-Y->l,alpha,px+(X->nc+X->l)*X->ld,X->ld,beta,py+(Y->nc+Y->l)*Y->ld,Y->ld));
 33:   PetscCall(MatDenseRestoreArrayRead(x->A,&px));
 34:   PetscCall(MatDenseRestoreArray(y->A,&py));
 35:   PetscFunctionReturn(PETSC_SUCCESS);
 36: }

 38: static PetscErrorCode BVMultVec_Mat(BV X,PetscScalar alpha,PetscScalar beta,Vec y,PetscScalar *q)
 39: {
 40:   BV_MAT            *x = (BV_MAT*)X->data;
 41:   PetscScalar       *py,*qq=q;
 42:   const PetscScalar *px;

 44:   PetscFunctionBegin;
 45:   PetscCall(MatDenseGetArrayRead(x->A,&px));
 46:   PetscCall(VecGetArray(y,&py));
 47:   if (!q) PetscCall(VecGetArray(X->buffer,&qq));
 48:   PetscCall(BVMultVec_BLAS_Private(X,X->n,X->k-X->l,alpha,px+(X->nc+X->l)*X->ld,X->ld,qq,beta,py));
 49:   if (!q) PetscCall(VecRestoreArray(X->buffer,&qq));
 50:   PetscCall(MatDenseRestoreArrayRead(x->A,&px));
 51:   PetscCall(VecRestoreArray(y,&py));
 52:   PetscFunctionReturn(PETSC_SUCCESS);
 53: }

 55: static PetscErrorCode BVMultInPlace_Mat(BV V,Mat Q,PetscInt s,PetscInt e)
 56: {
 57:   BV_MAT            *ctx = (BV_MAT*)V->data;
 58:   PetscScalar       *pv;
 59:   const PetscScalar *q;
 60:   PetscInt          ldq;

 62:   PetscFunctionBegin;
 63:   if (s>=e || !V->n) PetscFunctionReturn(PETSC_SUCCESS);
 64:   PetscCall(MatDenseGetLDA(Q,&ldq));
 65:   PetscCall(MatDenseGetArray(ctx->A,&pv));
 66:   PetscCall(MatDenseGetArrayRead(Q,&q));
 67:   PetscCall(BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,s-V->l,e-V->l,pv+(V->nc+V->l)*V->ld,V->ld,q+V->l*ldq+V->l,ldq,PETSC_FALSE));
 68:   PetscCall(MatDenseRestoreArrayRead(Q,&q));
 69:   PetscCall(MatDenseRestoreArray(ctx->A,&pv));
 70:   PetscFunctionReturn(PETSC_SUCCESS);
 71: }

 73: static PetscErrorCode BVMultInPlaceHermitianTranspose_Mat(BV V,Mat Q,PetscInt s,PetscInt e)
 74: {
 75:   BV_MAT            *ctx = (BV_MAT*)V->data;
 76:   PetscScalar       *pv;
 77:   const PetscScalar *q;
 78:   PetscInt          ldq;

 80:   PetscFunctionBegin;
 81:   if (s>=e || !V->n) PetscFunctionReturn(PETSC_SUCCESS);
 82:   PetscCall(MatDenseGetLDA(Q,&ldq));
 83:   PetscCall(MatDenseGetArray(ctx->A,&pv));
 84:   PetscCall(MatDenseGetArrayRead(Q,&q));
 85:   PetscCall(BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,s-V->l,e-V->l,pv+(V->nc+V->l)*V->ld,V->ld,q+V->l*ldq+V->l,ldq,PETSC_TRUE));
 86:   PetscCall(MatDenseRestoreArrayRead(Q,&q));
 87:   PetscCall(MatDenseRestoreArray(ctx->A,&pv));
 88:   PetscFunctionReturn(PETSC_SUCCESS);
 89: }

 91: static PetscErrorCode BVDot_Mat(BV X,BV Y,Mat M)
 92: {
 93:   BV_MAT            *x = (BV_MAT*)X->data,*y = (BV_MAT*)Y->data;
 94:   PetscScalar       *m;
 95:   const PetscScalar *px,*py;
 96:   PetscInt          ldm;

 98:   PetscFunctionBegin;
 99:   PetscCall(MatDenseGetLDA(M,&ldm));
100:   PetscCall(MatDenseGetArrayRead(x->A,&px));
101:   PetscCall(MatDenseGetArrayRead(y->A,&py));
102:   PetscCall(MatDenseGetArray(M,&m));
103:   PetscCall(BVDot_BLAS_Private(X,Y->k-Y->l,X->k-X->l,X->n,py+(Y->nc+Y->l)*Y->ld,Y->ld,px+(X->nc+X->l)*X->ld,X->ld,m+X->l*ldm+Y->l,ldm,x->mpi));
104:   PetscCall(MatDenseRestoreArray(M,&m));
105:   PetscCall(MatDenseRestoreArrayRead(x->A,&px));
106:   PetscCall(MatDenseRestoreArrayRead(y->A,&py));
107:   PetscFunctionReturn(PETSC_SUCCESS);
108: }

110: static PetscErrorCode BVDotVec_Mat(BV X,Vec y,PetscScalar *q)
111: {
112:   BV_MAT            *x = (BV_MAT*)X->data;
113:   PetscScalar       *qq=q;
114:   const PetscScalar *px,*py;
115:   Vec               z = y;

117:   PetscFunctionBegin;
118:   if (PetscUnlikely(X->matrix)) {
119:     PetscCall(BV_IPMatMult(X,y));
120:     z = X->Bx;
121:   }
122:   PetscCall(MatDenseGetArrayRead(x->A,&px));
123:   PetscCall(VecGetArrayRead(z,&py));
124:   if (!q) PetscCall(VecGetArray(X->buffer,&qq));
125:   PetscCall(BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->ld,X->ld,py,qq,x->mpi));
126:   if (!q) PetscCall(VecRestoreArray(X->buffer,&qq));
127:   PetscCall(VecRestoreArrayRead(z,&py));
128:   PetscCall(MatDenseRestoreArrayRead(x->A,&px));
129:   PetscFunctionReturn(PETSC_SUCCESS);
130: }

132: static PetscErrorCode BVDotVec_Local_Mat(BV X,Vec y,PetscScalar *m)
133: {
134:   BV_MAT            *x = (BV_MAT*)X->data;
135:   const PetscScalar *px,*py;
136:   Vec               z = y;

138:   PetscFunctionBegin;
139:   if (PetscUnlikely(X->matrix)) {
140:     PetscCall(BV_IPMatMult(X,y));
141:     z = X->Bx;
142:   }
143:   PetscCall(MatDenseGetArrayRead(x->A,&px));
144:   PetscCall(VecGetArrayRead(z,&py));
145:   PetscCall(BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->ld,X->ld,py,m,PETSC_FALSE));
146:   PetscCall(VecRestoreArrayRead(z,&py));
147:   PetscCall(MatDenseRestoreArrayRead(x->A,&px));
148:   PetscFunctionReturn(PETSC_SUCCESS);
149: }

151: static PetscErrorCode BVScale_Mat(BV bv,PetscInt j,PetscScalar alpha)
152: {
153:   BV_MAT         *ctx = (BV_MAT*)bv->data;
154:   PetscScalar    *array;

156:   PetscFunctionBegin;
157:   if (!bv->n) PetscFunctionReturn(PETSC_SUCCESS);
158:   PetscCall(MatDenseGetArray(ctx->A,&array));
159:   if (PetscUnlikely(j<0)) PetscCall(BVScale_BLAS_Private(bv,(bv->k-bv->l)*bv->ld,array+(bv->nc+bv->l)*bv->ld,alpha));
160:   else PetscCall(BVScale_BLAS_Private(bv,bv->n,array+(bv->nc+j)*bv->ld,alpha));
161:   PetscCall(MatDenseRestoreArray(ctx->A,&array));
162:   PetscFunctionReturn(PETSC_SUCCESS);
163: }

165: static PetscErrorCode BVNorm_Mat(BV bv,PetscInt j,NormType type,PetscReal *val)
166: {
167:   BV_MAT            *ctx = (BV_MAT*)bv->data;
168:   const PetscScalar *array;

170:   PetscFunctionBegin;
171:   PetscCall(MatDenseGetArrayRead(ctx->A,&array));
172:   if (PetscUnlikely(j<0)) PetscCall(BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->ld,bv->ld,type,val,ctx->mpi));
173:   else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->ld,bv->ld,type,val,ctx->mpi));
174:   PetscCall(MatDenseRestoreArrayRead(ctx->A,&array));
175:   PetscFunctionReturn(PETSC_SUCCESS);
176: }

178: static PetscErrorCode BVNorm_Local_Mat(BV bv,PetscInt j,NormType type,PetscReal *val)
179: {
180:   BV_MAT            *ctx = (BV_MAT*)bv->data;
181:   const PetscScalar *array;

183:   PetscFunctionBegin;
184:   PetscCall(MatDenseGetArrayRead(ctx->A,&array));
185:   if (PetscUnlikely(j<0)) PetscCall(BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->ld,bv->ld,type,val,PETSC_FALSE));
186:   else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->ld,bv->ld,type,val,PETSC_FALSE));
187:   PetscCall(MatDenseRestoreArrayRead(ctx->A,&array));
188:   PetscFunctionReturn(PETSC_SUCCESS);
189: }

191: static PetscErrorCode BVNormalize_Mat(BV bv,PetscScalar *eigi)
192: {
193:   BV_MAT         *ctx = (BV_MAT*)bv->data;
194:   PetscScalar    *array,*wi=NULL;

196:   PetscFunctionBegin;
197:   PetscCall(MatDenseGetArray(ctx->A,&array));
198:   if (eigi) wi = eigi+bv->l;
199:   PetscCall(BVNormalize_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->ld,bv->ld,wi,ctx->mpi));
200:   PetscCall(MatDenseRestoreArray(ctx->A,&array));
201:   PetscFunctionReturn(PETSC_SUCCESS);
202: }

204: static PetscErrorCode BVMatMult_Mat(BV V,Mat A,BV W)
205: {
206:   PetscInt       j;
207:   Mat            Vmat,Wmat;
208:   Vec            vv,ww;

210:   PetscFunctionBegin;
211:   if (V->vmm) {
212:     PetscCall(BVGetMat(V,&Vmat));
213:     PetscCall(BVGetMat(W,&Wmat));
214:     PetscCall(MatProductCreateWithMat(A,Vmat,NULL,Wmat));
215:     PetscCall(MatProductSetType(Wmat,MATPRODUCT_AB));
216:     PetscCall(MatProductSetFromOptions(Wmat));
217:     PetscCall(MatProductSymbolic(Wmat));
218:     PetscCall(MatProductNumeric(Wmat));
219:     PetscCall(MatProductClear(Wmat));
220:     PetscCall(BVRestoreMat(V,&Vmat));
221:     PetscCall(BVRestoreMat(W,&Wmat));
222:   } else {
223:     for (j=0;j<V->k-V->l;j++) {
224:       PetscCall(BVGetColumn(V,V->l+j,&vv));
225:       PetscCall(BVGetColumn(W,W->l+j,&ww));
226:       PetscCall(MatMult(A,vv,ww));
227:       PetscCall(BVRestoreColumn(V,V->l+j,&vv));
228:       PetscCall(BVRestoreColumn(W,W->l+j,&ww));
229:     }
230:   }
231:   PetscFunctionReturn(PETSC_SUCCESS);
232: }

234: static PetscErrorCode BVCopy_Mat(BV V,BV W)
235: {
236:   BV_MAT            *v = (BV_MAT*)V->data,*w = (BV_MAT*)W->data;
237:   const PetscScalar *pv;
238:   PetscScalar       *pw;
239:   PetscInt          j;

241:   PetscFunctionBegin;
242:   PetscCall(MatDenseGetArrayRead(v->A,&pv));
243:   PetscCall(MatDenseGetArray(w->A,&pw));
244:   for (j=0;j<V->k-V->l;j++) PetscCall(PetscArraycpy(pw+(W->nc+W->l+j)*W->ld,pv+(V->nc+V->l+j)*V->ld,V->n));
245:   PetscCall(MatDenseRestoreArrayRead(v->A,&pv));
246:   PetscCall(MatDenseRestoreArray(w->A,&pw));
247:   PetscFunctionReturn(PETSC_SUCCESS);
248: }

250: static PetscErrorCode BVCopyColumn_Mat(BV V,PetscInt j,PetscInt i)
251: {
252:   BV_MAT         *v = (BV_MAT*)V->data;
253:   PetscScalar    *pv;

255:   PetscFunctionBegin;
256:   PetscCall(MatDenseGetArray(v->A,&pv));
257:   PetscCall(PetscArraycpy(pv+(V->nc+i)*V->ld,pv+(V->nc+j)*V->ld,V->n));
258:   PetscCall(MatDenseRestoreArray(v->A,&pv));
259:   PetscFunctionReturn(PETSC_SUCCESS);
260: }

262: static PetscErrorCode BVResize_Mat(BV bv,PetscInt m,PetscBool copy)
263: {
264:   BV_MAT            *ctx = (BV_MAT*)bv->data;
265:   Mat               A,Msrc,Mdst;
266:   VecType           vtype;
267:   char              str[50];

269:   PetscFunctionBegin;
270:   PetscCall(VecGetType(bv->t,&vtype));
271:   PetscCall(MatCreateDenseFromVecType(PetscObjectComm((PetscObject)bv->t),vtype,bv->n,PETSC_DECIDE,bv->N,m,bv->ld,NULL,&A));
272:   if (((PetscObject)bv)->name) {
273:     PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
274:     PetscCall(PetscObjectSetName((PetscObject)A,str));
275:   }
276:   if (copy) {
277:     PetscCall(MatDenseGetSubMatrix(ctx->A,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PetscMin(m,bv->m),&Msrc));
278:     PetscCall(MatDenseGetSubMatrix(A,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PetscMin(m,bv->m),&Mdst));
279:     PetscCall(MatCopy(Msrc,Mdst,SAME_NONZERO_PATTERN));
280:     PetscCall(MatDenseRestoreSubMatrix(ctx->A,&Msrc));
281:     PetscCall(MatDenseRestoreSubMatrix(A,&Mdst));
282:   }
283:   PetscCall(MatDestroy(&ctx->A));
284:   ctx->A = A;
285:   PetscFunctionReturn(PETSC_SUCCESS);
286: }

288: static PetscErrorCode BVGetColumn_Mat(BV bv,PetscInt j,Vec *v)
289: {
290:   BV_MAT         *ctx = (BV_MAT*)bv->data;
291:   PetscScalar    *pA;
292:   PetscInt       l;

294:   PetscFunctionBegin;
295:   l = BVAvailableVec;
296:   PetscCall(MatDenseGetArray(ctx->A,&pA));
297:   PetscCall(VecPlaceArray(bv->cv[l],pA+(bv->nc+j)*bv->ld));
298:   PetscFunctionReturn(PETSC_SUCCESS);
299: }

301: static PetscErrorCode BVRestoreColumn_Mat(BV bv,PetscInt j,Vec *v)
302: {
303:   BV_MAT         *ctx = (BV_MAT*)bv->data;
304:   PetscScalar    *pA;
305:   PetscInt       l;

307:   PetscFunctionBegin;
308:   l = (j==bv->ci[0])? 0: 1;
309:   PetscCall(VecResetArray(bv->cv[l]));
310:   PetscCall(MatDenseRestoreArray(ctx->A,&pA));
311:   PetscFunctionReturn(PETSC_SUCCESS);
312: }

314: static PetscErrorCode BVGetArray_Mat(BV bv,PetscScalar **a)
315: {
316:   BV_MAT         *ctx = (BV_MAT*)bv->data;

318:   PetscFunctionBegin;
319:   PetscCall(MatDenseGetArray(ctx->A,a));
320:   PetscFunctionReturn(PETSC_SUCCESS);
321: }

323: static PetscErrorCode BVRestoreArray_Mat(BV bv,PetscScalar **a)
324: {
325:   BV_MAT         *ctx = (BV_MAT*)bv->data;

327:   PetscFunctionBegin;
328:   if (a) PetscCall(MatDenseRestoreArray(ctx->A,a));
329:   PetscFunctionReturn(PETSC_SUCCESS);
330: }

332: static PetscErrorCode BVGetArrayRead_Mat(BV bv,const PetscScalar **a)
333: {
334:   BV_MAT         *ctx = (BV_MAT*)bv->data;

336:   PetscFunctionBegin;
337:   PetscCall(MatDenseGetArrayRead(ctx->A,a));
338:   PetscFunctionReturn(PETSC_SUCCESS);
339: }

341: static PetscErrorCode BVRestoreArrayRead_Mat(BV bv,const PetscScalar **a)
342: {
343:   BV_MAT         *ctx = (BV_MAT*)bv->data;

345:   PetscFunctionBegin;
346:   if (a) PetscCall(MatDenseRestoreArrayRead(ctx->A,a));
347:   PetscFunctionReturn(PETSC_SUCCESS);
348: }

350: static PetscErrorCode BVView_Mat(BV bv,PetscViewer viewer)
351: {
352:   BV_MAT            *ctx = (BV_MAT*)bv->data;
353:   PetscViewerFormat format;
354:   PetscBool         isascii;
355:   const char        *bvname,*name;

357:   PetscFunctionBegin;
358:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
359:   if (isascii) {
360:     PetscCall(PetscViewerGetFormat(viewer,&format));
361:     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) PetscFunctionReturn(PETSC_SUCCESS);
362:     PetscCall(MatView(ctx->A,viewer));
363:     if (format == PETSC_VIEWER_ASCII_MATLAB) {
364:       PetscCall(PetscObjectGetName((PetscObject)bv,&bvname));
365:       PetscCall(PetscObjectGetName((PetscObject)ctx->A,&name));
366:       PetscCall(PetscViewerASCIIPrintf(viewer,"%s=%s;clear %s\n",bvname,name,name));
367:       if (bv->nc) PetscCall(PetscViewerASCIIPrintf(viewer,"%s=%s(:,%" PetscInt_FMT ":end);\n",bvname,bvname,bv->nc+1));
368:     }
369:   } else PetscCall(MatView(ctx->A,viewer));
370:   PetscFunctionReturn(PETSC_SUCCESS);
371: }

373: static PetscErrorCode BVDestroy_Mat(BV bv)
374: {
375:   BV_MAT         *ctx = (BV_MAT*)bv->data;

377:   PetscFunctionBegin;
378:   PetscCall(MatDestroy(&ctx->A));
379:   PetscCall(VecDestroy(&bv->cv[0]));
380:   PetscCall(VecDestroy(&bv->cv[1]));
381:   PetscCall(PetscFree(bv->data));
382:   PetscFunctionReturn(PETSC_SUCCESS);
383: }

385: SLEPC_EXTERN PetscErrorCode BVCreate_Mat(BV bv)
386: {
387:   BV_MAT         *ctx;
388:   PetscInt       nloc,bs,lsplit;
389:   PetscBool      seq;
390:   char           str[50];
391:   PetscScalar    *array,*ptr=NULL;
392:   BV             parent;
393:   Mat            Apar;
394:   VecType        vtype;

396:   PetscFunctionBegin;
397:   PetscCall(PetscNew(&ctx));
398:   bv->data = (void*)ctx;

400:   PetscCall(PetscObjectTypeCompareAny((PetscObject)bv->t,&bv->cuda,VECSEQCUDA,VECMPICUDA,""));
401:   PetscCall(PetscObjectTypeCompareAny((PetscObject)bv->t,&ctx->mpi,VECMPI,VECMPICUDA,""));

403:   PetscCall(PetscObjectTypeCompare((PetscObject)bv->t,VECSEQ,&seq));
404:   PetscCheck(seq || ctx->mpi || bv->cuda,PetscObjectComm((PetscObject)bv),PETSC_ERR_SUP,"BVMAT does not support the type of the provided template vector");

406:   PetscCall(VecGetLocalSize(bv->t,&nloc));
407:   PetscCall(VecGetBlockSize(bv->t,&bs));
408:   PetscCall(BV_SetDefaultLD(bv,nloc));

410:   if (PetscUnlikely(bv->issplit)) {
411:     /* split BV: share the memory of the parent BV */
412:     parent = bv->splitparent;
413:     lsplit = parent->lsplit;
414:     Apar = ((BV_MAT*)parent->data)->A;
415:     if (bv->cuda) {
416: #if defined(PETSC_HAVE_CUDA)
417:       PetscCall(MatDenseCUDAGetArray(Apar,&array));
418:       ptr = (bv->issplit==1)? array: array+lsplit*bv->ld;
419:       PetscCall(MatDenseCUDARestoreArray(Apar,&array));
420: #endif
421:     } else {
422:       PetscCall(MatDenseGetArray(Apar,&array));
423:       ptr = (bv->issplit==1)? array: array+lsplit*bv->ld;
424:       PetscCall(MatDenseRestoreArray(Apar,&array));
425:     }
426:   }

428:   PetscCall(VecGetType(bv->t,&vtype));
429:   PetscCall(MatCreateDenseFromVecType(PetscObjectComm((PetscObject)bv->t),vtype,nloc,PETSC_DECIDE,bv->N,bv->m,bv->ld,ptr,&ctx->A));
430:   if (((PetscObject)bv)->name) {
431:     PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
432:     PetscCall(PetscObjectSetName((PetscObject)ctx->A,str));
433:   }

435:   if (PetscUnlikely(bv->Acreate)) {
436:     PetscCall(MatConvert(bv->Acreate,MATDENSE,MAT_REUSE_MATRIX,&ctx->A));
437:     PetscCall(MatDestroy(&bv->Acreate));
438:   }

440:   PetscCall(VecDuplicateEmpty(bv->t,&bv->cv[0]));
441:   PetscCall(VecDuplicateEmpty(bv->t,&bv->cv[1]));

443:   if (bv->cuda) {
444: #if defined(PETSC_HAVE_CUDA)
445:     bv->ops->mult             = BVMult_Mat_CUDA;
446:     bv->ops->multvec          = BVMultVec_Mat_CUDA;
447:     bv->ops->multinplace      = BVMultInPlace_Mat_CUDA;
448:     bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Mat_CUDA;
449:     bv->ops->dot              = BVDot_Mat_CUDA;
450:     bv->ops->dotvec           = BVDotVec_Mat_CUDA;
451:     bv->ops->dotvec_local     = BVDotVec_Local_Mat_CUDA;
452:     bv->ops->scale            = BVScale_Mat_CUDA;
453:     bv->ops->matmult          = BVMatMult_Mat_CUDA;
454:     bv->ops->copy             = BVCopy_Mat_CUDA;
455:     bv->ops->copycolumn       = BVCopyColumn_Mat_CUDA;
456:     bv->ops->getcolumn        = BVGetColumn_Mat_CUDA;
457:     bv->ops->restorecolumn    = BVRestoreColumn_Mat_CUDA;
458:     bv->ops->restoresplit     = BVRestoreSplit_Mat_CUDA;
459:     bv->ops->getmat           = BVGetMat_Mat_CUDA;
460:     bv->ops->restoremat       = BVRestoreMat_Mat_CUDA;
461: #endif
462:   } else {
463:     bv->ops->mult             = BVMult_Mat;
464:     bv->ops->multvec          = BVMultVec_Mat;
465:     bv->ops->multinplace      = BVMultInPlace_Mat;
466:     bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Mat;
467:     bv->ops->dot              = BVDot_Mat;
468:     bv->ops->dotvec           = BVDotVec_Mat;
469:     bv->ops->dotvec_local     = BVDotVec_Local_Mat;
470:     bv->ops->scale            = BVScale_Mat;
471:     bv->ops->matmult          = BVMatMult_Mat;
472:     bv->ops->copy             = BVCopy_Mat;
473:     bv->ops->copycolumn       = BVCopyColumn_Mat;
474:     bv->ops->getcolumn        = BVGetColumn_Mat;
475:     bv->ops->restorecolumn    = BVRestoreColumn_Mat;
476:     bv->ops->getmat           = BVGetMat_Default;
477:     bv->ops->restoremat       = BVRestoreMat_Default;
478:   }
479:   bv->ops->norm             = BVNorm_Mat;
480:   bv->ops->norm_local       = BVNorm_Local_Mat;
481:   bv->ops->normalize        = BVNormalize_Mat;
482:   bv->ops->resize           = BVResize_Mat;
483:   bv->ops->getarray         = BVGetArray_Mat;
484:   bv->ops->restorearray     = BVRestoreArray_Mat;
485:   bv->ops->getarrayread     = BVGetArrayRead_Mat;
486:   bv->ops->restorearrayread = BVRestoreArrayRead_Mat;
487:   bv->ops->destroy          = BVDestroy_Mat;
488:   if (!ctx->mpi) bv->ops->view = BVView_Mat;
489:   PetscFunctionReturn(PETSC_SUCCESS);
490: }