Actual source code: svec.c

slepc-3.17.1 2022-04-11
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:   VecGetArrayRead(x->v,&px);
 25:   VecGetArray(y->v,&py);
 26:   if (Q) {
 27:     MatDenseGetLDA(Q,&ldq);
 28:     MatDenseGetArrayRead(Q,&q);
 29:     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);
 30:     MatDenseRestoreArrayRead(Q,&q);
 31:   } else 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);
 32:   VecRestoreArrayRead(x->v,&px);
 33:   VecRestoreArray(y->v,&py);
 34:   PetscFunctionReturn(0);
 35: }

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

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

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

 59:   MatDenseGetLDA(Q,&ldq);
 60:   VecGetArray(ctx->v,&pv);
 61:   MatDenseGetArrayRead(Q,&q);
 62:   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);
 63:   MatDenseRestoreArrayRead(Q,&q);
 64:   VecRestoreArray(ctx->v,&pv);
 65:   PetscFunctionReturn(0);
 66: }

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

 75:   MatDenseGetLDA(Q,&ldq);
 76:   VecGetArray(ctx->v,&pv);
 77:   MatDenseGetArrayRead(Q,&q);
 78:   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);
 79:   MatDenseRestoreArrayRead(Q,&q);
 80:   VecRestoreArray(ctx->v,&pv);
 81:   PetscFunctionReturn(0);
 82: }

 84: PetscErrorCode BVDot_Svec(BV X,BV Y,Mat M)
 85: {
 86:   BV_SVEC           *x = (BV_SVEC*)X->data,*y = (BV_SVEC*)Y->data;
 87:   const PetscScalar *px,*py;
 88:   PetscScalar       *m;
 89:   PetscInt          ldm;

 91:   MatDenseGetLDA(M,&ldm);
 92:   VecGetArrayRead(x->v,&px);
 93:   VecGetArrayRead(y->v,&py);
 94:   MatDenseGetArray(M,&m);
 95:   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);
 96:   MatDenseRestoreArray(M,&m);
 97:   VecRestoreArrayRead(x->v,&px);
 98:   VecRestoreArrayRead(y->v,&py);
 99:   PetscFunctionReturn(0);
100: }

102: PetscErrorCode BVDotVec_Svec(BV X,Vec y,PetscScalar *q)
103: {
104:   BV_SVEC           *x = (BV_SVEC*)X->data;
105:   const PetscScalar *px,*py;
106:   PetscScalar       *qq=q;
107:   Vec               z = y;

109:   if (PetscUnlikely(X->matrix)) {
110:     BV_IPMatMult(X,y);
111:     z = X->Bx;
112:   }
113:   VecGetArrayRead(x->v,&px);
114:   VecGetArrayRead(z,&py);
115:   if (!q) VecGetArray(X->buffer,&qq);
116:   BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->n,py,qq,x->mpi);
117:   if (!q) VecRestoreArray(X->buffer,&qq);
118:   VecRestoreArrayRead(z,&py);
119:   VecRestoreArrayRead(x->v,&px);
120:   PetscFunctionReturn(0);
121: }

123: PetscErrorCode BVDotVec_Local_Svec(BV X,Vec y,PetscScalar *m)
124: {
125:   BV_SVEC        *x = (BV_SVEC*)X->data;
126:   PetscScalar    *px,*py;
127:   Vec            z = y;

129:   if (PetscUnlikely(X->matrix)) {
130:     BV_IPMatMult(X,y);
131:     z = X->Bx;
132:   }
133:   VecGetArray(x->v,&px);
134:   VecGetArray(z,&py);
135:   BVDotVec_BLAS_Private(X,X->n,X->k-X->l,px+(X->nc+X->l)*X->n,py,m,PETSC_FALSE);
136:   VecRestoreArray(z,&py);
137:   VecRestoreArray(x->v,&px);
138:   PetscFunctionReturn(0);
139: }

141: PetscErrorCode BVScale_Svec(BV bv,PetscInt j,PetscScalar alpha)
142: {
143:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
144:   PetscScalar    *array;

146:   VecGetArray(ctx->v,&array);
147:   if (PetscUnlikely(j<0)) BVScale_BLAS_Private(bv,(bv->k-bv->l)*bv->n,array+(bv->nc+bv->l)*bv->n,alpha);
148:   else BVScale_BLAS_Private(bv,bv->n,array+(bv->nc+j)*bv->n,alpha);
149:   VecRestoreArray(ctx->v,&array);
150:   PetscFunctionReturn(0);
151: }

153: PetscErrorCode BVNorm_Svec(BV bv,PetscInt j,NormType type,PetscReal *val)
154: {
155:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
156:   PetscScalar    *array;

158:   VecGetArray(ctx->v,&array);
159:   if (PetscUnlikely(j<0)) BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->n,type,val,ctx->mpi);
160:   else BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->n,type,val,ctx->mpi);
161:   VecRestoreArray(ctx->v,&array);
162:   PetscFunctionReturn(0);
163: }

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

170:   VecGetArray(ctx->v,&array);
171:   if (PetscUnlikely(j<0)) BVNorm_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->n,type,val,PETSC_FALSE);
172:   else BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->n,type,val,PETSC_FALSE);
173:   VecRestoreArray(ctx->v,&array);
174:   PetscFunctionReturn(0);
175: }

177: PetscErrorCode BVNormalize_Svec(BV bv,PetscScalar *eigi)
178: {
179:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
180:   PetscScalar    *array,*wi=NULL;

182:   VecGetArray(ctx->v,&array);
183:   if (eigi) wi = eigi+bv->l;
184:   BVNormalize_LAPACK_Private(bv,bv->n,bv->k-bv->l,array+(bv->nc+bv->l)*bv->n,wi,ctx->mpi);
185:   VecRestoreArray(ctx->v,&array);
186:   PetscFunctionReturn(0);
187: }

189: PetscErrorCode BVMatMult_Svec(BV V,Mat A,BV W)
190: {
191:   PetscInt       j;
192:   Mat            Vmat,Wmat;
193:   Vec            vv,ww;

195:   if (V->vmm) {
196:     BVGetMat(V,&Vmat);
197:     BVGetMat(W,&Wmat);
198:     MatProductCreateWithMat(A,Vmat,NULL,Wmat);
199:     MatProductSetType(Wmat,MATPRODUCT_AB);
200:     MatProductSetFromOptions(Wmat);
201:     MatProductSymbolic(Wmat);
202:     MatProductNumeric(Wmat);
203:     MatProductClear(Wmat);
204:     BVRestoreMat(V,&Vmat);
205:     BVRestoreMat(W,&Wmat);
206:   } else {
207:     for (j=0;j<V->k-V->l;j++) {
208:       BVGetColumn(V,V->l+j,&vv);
209:       BVGetColumn(W,W->l+j,&ww);
210:       MatMult(A,vv,ww);
211:       BVRestoreColumn(V,V->l+j,&vv);
212:       BVRestoreColumn(W,W->l+j,&ww);
213:     }
214:   }
215:   PetscFunctionReturn(0);
216: }

218: PetscErrorCode BVCopy_Svec(BV V,BV W)
219: {
220:   BV_SVEC        *v = (BV_SVEC*)V->data,*w = (BV_SVEC*)W->data;
221:   PetscScalar    *pv,*pw,*pvc,*pwc;

223:   VecGetArray(v->v,&pv);
224:   VecGetArray(w->v,&pw);
225:   pvc = pv+(V->nc+V->l)*V->n;
226:   pwc = pw+(W->nc+W->l)*W->n;
227:   PetscArraycpy(pwc,pvc,(V->k-V->l)*V->n);
228:   VecRestoreArray(v->v,&pv);
229:   VecRestoreArray(w->v,&pw);
230:   PetscFunctionReturn(0);
231: }

233: PetscErrorCode BVCopyColumn_Svec(BV V,PetscInt j,PetscInt i)
234: {
235:   BV_SVEC        *v = (BV_SVEC*)V->data;
236:   PetscScalar    *pv;

238:   VecGetArray(v->v,&pv);
239:   PetscArraycpy(pv+(V->nc+i)*V->n,pv+(V->nc+j)*V->n,V->n);
240:   VecRestoreArray(v->v,&pv);
241:   PetscFunctionReturn(0);
242: }

244: PetscErrorCode BVResize_Svec(BV bv,PetscInt m,PetscBool copy)
245: {
246:   BV_SVEC           *ctx = (BV_SVEC*)bv->data;
247:   PetscScalar       *pnew;
248:   const PetscScalar *pv;
249:   PetscInt          bs;
250:   Vec               vnew;
251:   char              str[50];

253:   VecGetBlockSize(bv->t,&bs);
254:   VecCreate(PetscObjectComm((PetscObject)bv->t),&vnew);
255:   VecSetType(vnew,((PetscObject)bv->t)->type_name);
256:   VecSetSizes(vnew,m*bv->n,PETSC_DECIDE);
257:   VecSetBlockSize(vnew,bs);
258:   PetscLogObjectParent((PetscObject)bv,(PetscObject)vnew);
259:   if (((PetscObject)bv)->name) {
260:     PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name);
261:     PetscObjectSetName((PetscObject)vnew,str);
262:   }
263:   if (copy) {
264:     VecGetArrayRead(ctx->v,&pv);
265:     VecGetArray(vnew,&pnew);
266:     PetscArraycpy(pnew,pv,PetscMin(m,bv->m)*bv->n);
267:     VecRestoreArrayRead(ctx->v,&pv);
268:     VecRestoreArray(vnew,&pnew);
269:   }
270:   VecDestroy(&ctx->v);
271:   ctx->v = vnew;
272:   PetscFunctionReturn(0);
273: }

275: PetscErrorCode BVGetColumn_Svec(BV bv,PetscInt j,Vec *v)
276: {
277:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
278:   PetscScalar    *pv;
279:   PetscInt       l;

281:   l = BVAvailableVec;
282:   VecGetArray(ctx->v,&pv);
283:   VecPlaceArray(bv->cv[l],pv+(bv->nc+j)*bv->n);
284:   PetscFunctionReturn(0);
285: }

287: PetscErrorCode BVRestoreColumn_Svec(BV bv,PetscInt j,Vec *v)
288: {
289:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;
290:   PetscInt       l;

292:   l = (j==bv->ci[0])? 0: 1;
293:   VecResetArray(bv->cv[l]);
294:   VecRestoreArray(ctx->v,NULL);
295:   PetscFunctionReturn(0);
296: }

298: PetscErrorCode BVGetArray_Svec(BV bv,PetscScalar **a)
299: {
300:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;

302:   VecGetArray(ctx->v,a);
303:   PetscFunctionReturn(0);
304: }

306: PetscErrorCode BVRestoreArray_Svec(BV bv,PetscScalar **a)
307: {
308:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;

310:   VecRestoreArray(ctx->v,a);
311:   PetscFunctionReturn(0);
312: }

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

318:   VecGetArrayRead(ctx->v,a);
319:   PetscFunctionReturn(0);
320: }

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

326:   VecRestoreArrayRead(ctx->v,a);
327:   PetscFunctionReturn(0);
328: }

330: PetscErrorCode BVView_Svec(BV bv,PetscViewer viewer)
331: {
332:   BV_SVEC           *ctx = (BV_SVEC*)bv->data;
333:   PetscViewerFormat format;
334:   PetscBool         isascii;
335:   const char        *bvname,*name;

337:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
338:   if (isascii) {
339:     PetscViewerGetFormat(viewer,&format);
340:     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) PetscFunctionReturn(0);
341:     VecView(ctx->v,viewer);
342:     if (format == PETSC_VIEWER_ASCII_MATLAB) {
343:       PetscObjectGetName((PetscObject)bv,&bvname);
344:       PetscObjectGetName((PetscObject)ctx->v,&name);
345:       PetscViewerASCIIPrintf(viewer,"%s=reshape(%s,%" PetscInt_FMT ",%" PetscInt_FMT ");clear %s\n",bvname,name,bv->N,bv->nc+bv->m,name);
346:       if (bv->nc) PetscViewerASCIIPrintf(viewer,"%s=%s(:,%" PetscInt_FMT ":end);\n",bvname,bvname,bv->nc+1);
347:     }
348:   } else VecView(ctx->v,viewer);
349:   PetscFunctionReturn(0);
350: }

352: PetscErrorCode BVDestroy_Svec(BV bv)
353: {
354:   BV_SVEC        *ctx = (BV_SVEC*)bv->data;

356:   VecDestroy(&ctx->v);
357:   VecDestroy(&bv->cv[0]);
358:   VecDestroy(&bv->cv[1]);
359:   PetscFree(bv->data);
360:   bv->cuda = PETSC_FALSE;
361:   PetscFunctionReturn(0);
362: }

364: SLEPC_EXTERN PetscErrorCode BVCreate_Svec(BV bv)
365: {
366:   BV_SVEC           *ctx;
367:   PetscInt          nloc,N,bs,tglobal=0,tlocal,lsplit;
368:   PetscBool         seq;
369:   PetscScalar       *vv;
370:   const PetscScalar *aa,*array,*ptr;
371:   char              str[50];
372:   BV                parent;
373:   Vec               vpar;
374: #if defined(PETSC_HAVE_CUDA)
375:   PetscScalar       *gpuarray,*gptr;
376: #endif

378:   PetscNewLog(bv,&ctx);
379:   bv->data = (void*)ctx;

381:   PetscObjectTypeCompareAny((PetscObject)bv->t,&bv->cuda,VECSEQCUDA,VECMPICUDA,"");
382:   PetscObjectTypeCompareAny((PetscObject)bv->t,&ctx->mpi,VECMPI,VECMPICUDA,"");

384:   PetscObjectTypeCompare((PetscObject)bv->t,VECSEQ,&seq);

387:   VecGetLocalSize(bv->t,&nloc);
388:   VecGetSize(bv->t,&N);
389:   VecGetBlockSize(bv->t,&bs);
390:   tlocal = bv->m*nloc;
391:   PetscIntMultError(bv->m,N,&tglobal);

393:   if (PetscUnlikely(bv->issplit)) {
394:     /* split BV: create Vec sharing the memory of the parent BV */
395:     parent = bv->splitparent;
396:     lsplit = parent->lsplit;
397:     vpar = ((BV_SVEC*)parent->data)->v;
398:     if (bv->cuda) {
399: #if defined(PETSC_HAVE_CUDA)
400:       VecCUDAGetArray(vpar,&gpuarray);
401:       gptr = (bv->issplit==1)? gpuarray: gpuarray+lsplit*nloc;
402:       VecCUDARestoreArray(vpar,&gpuarray);
403:       if (ctx->mpi) VecCreateMPICUDAWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,bv->m*N,NULL,&ctx->v);
404:       else VecCreateSeqCUDAWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,NULL,&ctx->v);
405:       VecCUDAPlaceArray(ctx->v,gptr);
406: #endif
407:     } else {
408:       VecGetArrayRead(vpar,&array);
409:       ptr = (bv->issplit==1)? array: array+lsplit*nloc;
410:       VecRestoreArrayRead(vpar,&array);
411:       if (ctx->mpi) VecCreateMPIWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,bv->m*N,NULL,&ctx->v);
412:       else VecCreateSeqWithArray(PetscObjectComm((PetscObject)bv->t),bs,tlocal,NULL,&ctx->v);
413:       VecPlaceArray(ctx->v,ptr);
414:     }
415:   } else {
416:     /* regular BV: create Vec to store the BV entries */
417:     VecCreate(PetscObjectComm((PetscObject)bv->t),&ctx->v);
418:     VecSetType(ctx->v,((PetscObject)bv->t)->type_name);
419:     VecSetSizes(ctx->v,tlocal,tglobal);
420:     VecSetBlockSize(ctx->v,bs);
421:   }
422:   PetscLogObjectParent((PetscObject)bv,(PetscObject)ctx->v);
423:   if (((PetscObject)bv)->name) {
424:     PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name);
425:     PetscObjectSetName((PetscObject)ctx->v,str);
426:   }

428:   if (PetscUnlikely(bv->Acreate)) {
429:     MatDenseGetArrayRead(bv->Acreate,&aa);
430:     VecGetArray(ctx->v,&vv);
431:     PetscArraycpy(vv,aa,tlocal);
432:     VecRestoreArray(ctx->v,&vv);
433:     MatDenseRestoreArrayRead(bv->Acreate,&aa);
434:     MatDestroy(&bv->Acreate);
435:   }

437:   VecDuplicateEmpty(bv->t,&bv->cv[0]);
438:   VecDuplicateEmpty(bv->t,&bv->cv[1]);

440:   if (bv->cuda) {
441: #if defined(PETSC_HAVE_CUDA)
442:     bv->ops->mult             = BVMult_Svec_CUDA;
443:     bv->ops->multvec          = BVMultVec_Svec_CUDA;
444:     bv->ops->multinplace      = BVMultInPlace_Svec_CUDA;
445:     bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Svec_CUDA;
446:     bv->ops->dot              = BVDot_Svec_CUDA;
447:     bv->ops->dotvec           = BVDotVec_Svec_CUDA;
448:     bv->ops->dotvec_local     = BVDotVec_Local_Svec_CUDA;
449:     bv->ops->scale            = BVScale_Svec_CUDA;
450:     bv->ops->matmult          = BVMatMult_Svec_CUDA;
451:     bv->ops->copy             = BVCopy_Svec_CUDA;
452:     bv->ops->copycolumn       = BVCopyColumn_Svec_CUDA;
453:     bv->ops->resize           = BVResize_Svec_CUDA;
454:     bv->ops->getcolumn        = BVGetColumn_Svec_CUDA;
455:     bv->ops->restorecolumn    = BVRestoreColumn_Svec_CUDA;
456:     bv->ops->restoresplit     = BVRestoreSplit_Svec_CUDA;
457:     bv->ops->getmat           = BVGetMat_Svec_CUDA;
458:     bv->ops->restoremat       = BVRestoreMat_Svec_CUDA;
459: #endif
460:   } else {
461:     bv->ops->mult             = BVMult_Svec;
462:     bv->ops->multvec          = BVMultVec_Svec;
463:     bv->ops->multinplace      = BVMultInPlace_Svec;
464:     bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Svec;
465:     bv->ops->dot              = BVDot_Svec;
466:     bv->ops->dotvec           = BVDotVec_Svec;
467:     bv->ops->dotvec_local     = BVDotVec_Local_Svec;
468:     bv->ops->scale            = BVScale_Svec;
469:     bv->ops->matmult          = BVMatMult_Svec;
470:     bv->ops->copy             = BVCopy_Svec;
471:     bv->ops->copycolumn       = BVCopyColumn_Svec;
472:     bv->ops->resize           = BVResize_Svec;
473:     bv->ops->getcolumn        = BVGetColumn_Svec;
474:     bv->ops->restorecolumn    = BVRestoreColumn_Svec;
475:   }
476:   bv->ops->norm             = BVNorm_Svec;
477:   bv->ops->norm_local       = BVNorm_Local_Svec;
478:   bv->ops->normalize        = BVNormalize_Svec;
479:   bv->ops->getarray         = BVGetArray_Svec;
480:   bv->ops->restorearray     = BVRestoreArray_Svec;
481:   bv->ops->getarrayread     = BVGetArrayRead_Svec;
482:   bv->ops->restorearrayread = BVRestoreArrayRead_Svec;
483:   bv->ops->destroy          = BVDestroy_Svec;
484:   if (!ctx->mpi) bv->ops->view = BVView_Svec;
485:   PetscFunctionReturn(0);
486: }