Actual source code: svec.c

slepc-3.19.0 2023-03-31
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 as a single Vec
 12: */

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

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

 24:   PetscFunctionBegin;
 25:   PetscCall(VecGetArrayRead(x->v,&px));
 26:   PetscCall(VecGetArray(y->v,&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,ldq,alpha,px+(X->nc+X->l)*X->n,q+Y->l*ldq+X->l,beta,py+(Y->nc+Y->l)*Y->n));
 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->n,beta,py+(Y->nc+Y->l)*Y->n));
 33:   PetscCall(VecRestoreArrayRead(x->v,&px));
 34:   PetscCall(VecRestoreArray(y->v,&py));
 35:   PetscFunctionReturn(PETSC_SUCCESS);
 36: }

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

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

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

 61:   PetscFunctionBegin;
 62:   PetscCall(MatDenseGetLDA(Q,&ldq));
 63:   PetscCall(VecGetArray(ctx->v,&pv));
 64:   PetscCall(MatDenseGetArrayRead(Q,&q));
 65:   PetscCall(BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,ldq,s-V->l,e-V->l,pv+(V->nc+V->l)*V->n,q+V->l*ldq+V->l,PETSC_FALSE));
 66:   PetscCall(MatDenseRestoreArrayRead(Q,&q));
 67:   PetscCall(VecRestoreArray(ctx->v,&pv));
 68:   PetscFunctionReturn(PETSC_SUCCESS);
 69: }

 71: PetscErrorCode BVMultInPlaceHermitianTranspose_Svec(BV V,Mat Q,PetscInt s,PetscInt e)
 72: {
 73:   BV_SVEC           *ctx = (BV_SVEC*)V->data;
 74:   PetscScalar       *pv;
 75:   const PetscScalar *q;
 76:   PetscInt          ldq;

 78:   PetscFunctionBegin;
 79:   PetscCall(MatDenseGetLDA(Q,&ldq));
 80:   PetscCall(VecGetArray(ctx->v,&pv));
 81:   PetscCall(MatDenseGetArrayRead(Q,&q));
 82:   PetscCall(BVMultInPlace_BLAS_Private(V,V->n,V->k-V->l,ldq,s-V->l,e-V->l,pv+(V->nc+V->l)*V->n,q+V->l*ldq+V->l,PETSC_TRUE));
 83:   PetscCall(MatDenseRestoreArrayRead(Q,&q));
 84:   PetscCall(VecRestoreArray(ctx->v,&pv));
 85:   PetscFunctionReturn(PETSC_SUCCESS);
 86: }

 88: PetscErrorCode BVDot_Svec(BV X,BV Y,Mat M)
 89: {
 90:   BV_SVEC           *x = (BV_SVEC*)X->data,*y = (BV_SVEC*)Y->data;
 91:   const PetscScalar *px,*py;
 92:   PetscScalar       *m;
 93:   PetscInt          ldm;

 95:   PetscFunctionBegin;
 96:   PetscCall(MatDenseGetLDA(M,&ldm));
 97:   PetscCall(VecGetArrayRead(x->v,&px));
 98:   PetscCall(VecGetArrayRead(y->v,&py));
 99:   PetscCall(MatDenseGetArray(M,&m));
100:   PetscCall(BVDot_BLAS_Private(X,Y->k-Y->l,X->k-X->l,X->n,ldm,py+(Y->nc+Y->l)*Y->n,px+(X->nc+X->l)*X->n,m+X->l*ldm+Y->l,x->mpi));
101:   PetscCall(MatDenseRestoreArray(M,&m));
102:   PetscCall(VecRestoreArrayRead(x->v,&px));
103:   PetscCall(VecRestoreArrayRead(y->v,&py));
104:   PetscFunctionReturn(PETSC_SUCCESS);
105: }

107: PetscErrorCode BVDotVec_Svec(BV X,Vec y,PetscScalar *q)
108: {
109:   BV_SVEC           *x = (BV_SVEC*)X->data;
110:   const PetscScalar *px,*py;
111:   PetscScalar       *qq=q;
112:   Vec               z = y;

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

129: PetscErrorCode BVDotVec_Local_Svec(BV X,Vec y,PetscScalar *m)
130: {
131:   BV_SVEC        *x = (BV_SVEC*)X->data;
132:   PetscScalar    *px,*py;
133:   Vec            z = y;

135:   PetscFunctionBegin;
136:   if (PetscUnlikely(X->matrix)) {
137:     PetscCall(BV_IPMatMult(X,y));
138:     z = X->Bx;
139:   }
140:   PetscCall(VecGetArray(x->v,&px));
141:   PetscCall(VecGetArray(z,&py));
142:   PetscCall(BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->n,py,m,PETSC_FALSE));
143:   PetscCall(VecRestoreArray(z,&py));
144:   PetscCall(VecRestoreArray(x->v,&px));
145:   PetscFunctionReturn(PETSC_SUCCESS);
146: }

148: PetscErrorCode BVScale_Svec(BV bv,PetscInt j,PetscScalar alpha)
149: {
150:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
151:   PetscScalar    *array;

153:   PetscFunctionBegin;
154:   PetscCall(VecGetArray(ctx->v,&array));
155:   if (PetscUnlikely(j<0)) PetscCall(BVScale_BLAS_Private(bv,(bv->k-bv->l)*bv->n,array+(bv->nc+bv->l)*bv->n,alpha));
156:   else PetscCall(BVScale_BLAS_Private(bv,bv->n,array+(bv->nc+j)*bv->n,alpha));
157:   PetscCall(VecRestoreArray(ctx->v,&array));
158:   PetscFunctionReturn(PETSC_SUCCESS);
159: }

161: PetscErrorCode BVNorm_Svec(BV bv,PetscInt j,NormType type,PetscReal *val)
162: {
163:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
164:   PetscScalar    *array;

166:   PetscFunctionBegin;
167:   PetscCall(VecGetArray(ctx->v,&array));
168:   if (PetscUnlikely(j<0)) PetscCall(BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->n,type,val,ctx->mpi));
169:   else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->n,type,val,ctx->mpi));
170:   PetscCall(VecRestoreArray(ctx->v,&array));
171:   PetscFunctionReturn(PETSC_SUCCESS);
172: }

174: PetscErrorCode BVNorm_Local_Svec(BV bv,PetscInt j,NormType type,PetscReal *val)
175: {
176:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
177:   PetscScalar    *array;

179:   PetscFunctionBegin;
180:   PetscCall(VecGetArray(ctx->v,&array));
181:   if (PetscUnlikely(j<0)) PetscCall(BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->n,type,val,PETSC_FALSE));
182:   else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->n,type,val,PETSC_FALSE));
183:   PetscCall(VecRestoreArray(ctx->v,&array));
184:   PetscFunctionReturn(PETSC_SUCCESS);
185: }

187: PetscErrorCode BVNormalize_Svec(BV bv,PetscScalar *eigi)
188: {
189:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
190:   PetscScalar    *array,*wi=NULL;

192:   PetscFunctionBegin;
193:   PetscCall(VecGetArray(ctx->v,&array));
194:   if (eigi) wi = eigi+bv->l;
195:   PetscCall(BVNormalize_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->n,wi,ctx->mpi));
196:   PetscCall(VecRestoreArray(ctx->v,&array));
197:   PetscFunctionReturn(PETSC_SUCCESS);
198: }

200: PetscErrorCode BVMatMult_Svec(BV V,Mat A,BV W)
201: {
202:   PetscInt       j;
203:   Mat            Vmat,Wmat;
204:   Vec            vv,ww;

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

230: PetscErrorCode BVCopy_Svec(BV V,BV W)
231: {
232:   BV_SVEC        *v = (BV_SVEC*)V->data,*w = (BV_SVEC*)W->data;
233:   PetscScalar    *pv,*pw,*pvc,*pwc;

235:   PetscFunctionBegin;
236:   PetscCall(VecGetArray(v->v,&pv));
237:   PetscCall(VecGetArray(w->v,&pw));
238:   pvc = pv+(V->nc+V->l)*V->n;
239:   pwc = pw+(W->nc+W->l)*W->n;
240:   PetscCall(PetscArraycpy(pwc,pvc,(V->k-V->l)*V->n));
241:   PetscCall(VecRestoreArray(v->v,&pv));
242:   PetscCall(VecRestoreArray(w->v,&pw));
243:   PetscFunctionReturn(PETSC_SUCCESS);
244: }

246: PetscErrorCode BVCopyColumn_Svec(BV V,PetscInt j,PetscInt i)
247: {
248:   BV_SVEC        *v = (BV_SVEC*)V->data;
249:   PetscScalar    *pv;

251:   PetscFunctionBegin;
252:   PetscCall(VecGetArray(v->v,&pv));
253:   PetscCall(PetscArraycpy(pv+(V->nc+i)*V->n,pv+(V->nc+j)*V->n,V->n));
254:   PetscCall(VecRestoreArray(v->v,&pv));
255:   PetscFunctionReturn(PETSC_SUCCESS);
256: }

258: PetscErrorCode BVResize_Svec(BV bv,PetscInt m,PetscBool copy)
259: {
260:   BV_SVEC           *ctx = (BV_SVEC*)bv->data;
261:   PetscScalar       *pnew;
262:   const PetscScalar *pv;
263:   PetscInt          bs;
264:   Vec               vnew;
265:   char              str[50];

267:   PetscFunctionBegin;
268:   PetscCall(VecGetBlockSize(bv->t,&bs));
269:   PetscCall(VecCreate(PetscObjectComm((PetscObject)bv->t),&vnew));
270:   PetscCall(VecSetType(vnew,((PetscObject)bv->t)->type_name));
271:   PetscCall(VecSetSizes(vnew,m*bv->n,PETSC_DECIDE));
272:   PetscCall(VecSetBlockSize(vnew,bs));
273:   if (((PetscObject)bv)->name) {
274:     PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
275:     PetscCall(PetscObjectSetName((PetscObject)vnew,str));
276:   }
277:   if (copy) {
278:     PetscCall(VecGetArrayRead(ctx->v,&pv));
279:     PetscCall(VecGetArray(vnew,&pnew));
280:     PetscCall(PetscArraycpy(pnew,pv,PetscMin(m,bv->m)*bv->n));
281:     PetscCall(VecRestoreArrayRead(ctx->v,&pv));
282:     PetscCall(VecRestoreArray(vnew,&pnew));
283:   }
284:   PetscCall(VecDestroy(&ctx->v));
285:   ctx->v = vnew;
286:   PetscFunctionReturn(PETSC_SUCCESS);
287: }

289: PetscErrorCode BVGetColumn_Svec(BV bv,PetscInt j,Vec *v)
290: {
291:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
292:   PetscScalar    *pv;
293:   PetscInt       l;

295:   PetscFunctionBegin;
296:   l = BVAvailableVec;
297:   PetscCall(VecGetArray(ctx->v,&pv));
298:   PetscCall(VecPlaceArray(bv->cv[l],pv+(bv->nc+j)*bv->n));
299:   PetscFunctionReturn(PETSC_SUCCESS);
300: }

302: PetscErrorCode BVRestoreColumn_Svec(BV bv,PetscInt j,Vec *v)
303: {
304:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
305:   PetscInt       l;

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

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

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

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

327:   PetscFunctionBegin;
328:   PetscCall(VecRestoreArray(ctx->v,a));
329:   PetscFunctionReturn(PETSC_SUCCESS);
330: }

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

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

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

345:   PetscFunctionBegin;
346:   PetscCall(VecRestoreArrayRead(ctx->v,a));
347:   PetscFunctionReturn(PETSC_SUCCESS);
348: }

350: PetscErrorCode BVView_Svec(BV bv,PetscViewer viewer)
351: {
352:   BV_SVEC           *ctx = (BV_SVEC*)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(VecView(ctx->v,viewer));
363:     if (format == PETSC_VIEWER_ASCII_MATLAB) {
364:       PetscCall(PetscObjectGetName((PetscObject)bv,&bvname));
365:       PetscCall(PetscObjectGetName((PetscObject)ctx->v,&name));
366:       PetscCall(PetscViewerASCIIPrintf(viewer,"%s=reshape(%s,%" PetscInt_FMT ",%" PetscInt_FMT ");clear %s\n",bvname,name,bv->N,bv->nc+bv->m,name));
367:       if (bv->nc) PetscCall(PetscViewerASCIIPrintf(viewer,"%s=%s(:,%" PetscInt_FMT ":end);\n",bvname,bvname,bv->nc+1));
368:     }
369:   } else PetscCall(VecView(ctx->v,viewer));
370:   PetscFunctionReturn(PETSC_SUCCESS);
371: }

373: PetscErrorCode BVDestroy_Svec(BV bv)
374: {
375:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;

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

386: SLEPC_EXTERN PetscErrorCode BVCreate_Svec(BV bv)
387: {
388:   BV_SVEC           *ctx;
389:   PetscInt          nloc,N,bs,tglobal=0,tlocal,lsplit;
390:   PetscBool         seq;
391:   PetscScalar       *vv;
392:   const PetscScalar *aa,*array,*ptr;
393:   char              str[50];
394:   BV                parent;
395:   Vec               vpar;
396: #if defined(PETSC_HAVE_CUDA)
397:   PetscScalar       *gpuarray,*gptr;
398: #endif

400:   PetscFunctionBegin;
401:   PetscCall(PetscNew(&ctx));
402:   bv->data = (void*)ctx;

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

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

410:   PetscCall(VecGetLocalSize(bv->t,&nloc));
411:   PetscCall(VecGetSize(bv->t,&N));
412:   PetscCall(VecGetBlockSize(bv->t,&bs));
413:   tlocal = bv->m*nloc;
414:   PetscCall(PetscIntMultError(bv->m,N,&tglobal));

416:   if (PetscUnlikely(bv->issplit)) {
417:     /* split BV: create Vec sharing the memory of the parent BV */
418:     parent = bv->splitparent;
419:     lsplit = parent->lsplit;
420:     vpar = ((BV_SVEC*)parent->data)->v;
421:     if (bv->cuda) {
422: #if defined(PETSC_HAVE_CUDA)
423:       PetscCall(VecCUDAGetArray(vpar,&gpuarray));
424:       gptr = (bv->issplit==1)? gpuarray: gpuarray+lsplit*nloc;
425:       PetscCall(VecCUDARestoreArray(vpar,&gpuarray));
426:       if (ctx->mpi) PetscCall(VecCreateMPICUDAWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,bv->m*N,NULL,&ctx->v));
427:       else PetscCall(VecCreateSeqCUDAWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,NULL,&ctx->v));
428:       PetscCall(VecCUDAPlaceArray(ctx->v,gptr));
429: #endif
430:     } else {
431:       PetscCall(VecGetArrayRead(vpar,&array));
432:       ptr = (bv->issplit==1)? array: array+lsplit*nloc;
433:       PetscCall(VecRestoreArrayRead(vpar,&array));
434:       if (ctx->mpi) PetscCall(VecCreateMPIWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,bv->m*N,NULL,&ctx->v));
435:       else PetscCall(VecCreateSeqWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,NULL,&ctx->v));
436:       PetscCall(VecPlaceArray(ctx->v,ptr));
437:     }
438:   } else {
439:     /* regular BV: create Vec to store the BV entries */
440:     PetscCall(VecCreate(PetscObjectComm((PetscObject)bv->t),&ctx->v));
441:     PetscCall(VecSetType(ctx->v,((PetscObject)bv->t)->type_name));
442:     PetscCall(VecSetSizes(ctx->v,tlocal,tglobal));
443:     PetscCall(VecSetBlockSize(ctx->v,bs));
444:   }
445:   if (((PetscObject)bv)->name) {
446:     PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
447:     PetscCall(PetscObjectSetName((PetscObject)ctx->v,str));
448:   }

450:   if (PetscUnlikely(bv->Acreate)) {
451:     PetscCall(MatDenseGetArrayRead(bv->Acreate,&aa));
452:     PetscCall(VecGetArray(ctx->v,&vv));
453:     PetscCall(PetscArraycpy(vv,aa,tlocal));
454:     PetscCall(VecRestoreArray(ctx->v,&vv));
455:     PetscCall(MatDenseRestoreArrayRead(bv->Acreate,&aa));
456:     PetscCall(MatDestroy(&bv->Acreate));
457:   }

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

462:   if (bv->cuda) {
463: #if defined(PETSC_HAVE_CUDA)
464:     bv->ops->mult             = BVMult_Svec_CUDA;
465:     bv->ops->multvec          = BVMultVec_Svec_CUDA;
466:     bv->ops->multinplace      = BVMultInPlace_Svec_CUDA;
467:     bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Svec_CUDA;
468:     bv->ops->dot              = BVDot_Svec_CUDA;
469:     bv->ops->dotvec           = BVDotVec_Svec_CUDA;
470:     bv->ops->dotvec_local     = BVDotVec_Local_Svec_CUDA;
471:     bv->ops->scale            = BVScale_Svec_CUDA;
472:     bv->ops->matmult          = BVMatMult_Svec_CUDA;
473:     bv->ops->copy             = BVCopy_Svec_CUDA;
474:     bv->ops->copycolumn       = BVCopyColumn_Svec_CUDA;
475:     bv->ops->resize           = BVResize_Svec_CUDA;
476:     bv->ops->getcolumn        = BVGetColumn_Svec_CUDA;
477:     bv->ops->restorecolumn    = BVRestoreColumn_Svec_CUDA;
478:     bv->ops->restoresplit     = BVRestoreSplit_Svec_CUDA;
479:     bv->ops->getmat           = BVGetMat_Svec_CUDA;
480:     bv->ops->restoremat       = BVRestoreMat_Svec_CUDA;
481: #endif
482:   } else {
483:     bv->ops->mult             = BVMult_Svec;
484:     bv->ops->multvec          = BVMultVec_Svec;
485:     bv->ops->multinplace      = BVMultInPlace_Svec;
486:     bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Svec;
487:     bv->ops->dot              = BVDot_Svec;
488:     bv->ops->dotvec           = BVDotVec_Svec;
489:     bv->ops->dotvec_local     = BVDotVec_Local_Svec;
490:     bv->ops->scale            = BVScale_Svec;
491:     bv->ops->matmult          = BVMatMult_Svec;
492:     bv->ops->copy             = BVCopy_Svec;
493:     bv->ops->copycolumn       = BVCopyColumn_Svec;
494:     bv->ops->resize           = BVResize_Svec;
495:     bv->ops->getcolumn        = BVGetColumn_Svec;
496:     bv->ops->restorecolumn    = BVRestoreColumn_Svec;
497:     bv->ops->getmat           = BVGetMat_Default;
498:     bv->ops->restoremat       = BVRestoreMat_Default;
499:   }
500:   bv->ops->norm             = BVNorm_Svec;
501:   bv->ops->norm_local       = BVNorm_Local_Svec;
502:   bv->ops->normalize        = BVNormalize_Svec;
503:   bv->ops->getarray         = BVGetArray_Svec;
504:   bv->ops->restorearray     = BVRestoreArray_Svec;
505:   bv->ops->getarrayread     = BVGetArrayRead_Svec;
506:   bv->ops->restorearrayread = BVRestoreArrayRead_Svec;
507:   bv->ops->destroy          = BVDestroy_Svec;
508:   if (!ctx->mpi) bv->ops->view = BVView_Svec;
509:   PetscFunctionReturn(PETSC_SUCCESS);
510: }