Actual source code: bvmat.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 with a dense Mat
 12: */

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

 16: typedef struct {
 17:   Mat       A;
 18:   PetscBool mpi;
 19: } BV_MAT;

 21: PetscErrorCode BVMult_Mat(BV Y,PetscScalar alpha,PetscScalar beta,BV X,Mat Q)
 22: {
 23:   BV_MAT            *y = (BV_MAT*)Y->data,*x = (BV_MAT*)X->data;
 24:   PetscScalar       *py;
 25:   const PetscScalar *px,*q;
 26:   PetscInt          ldq;

 28:   PetscFunctionBegin;
 29:   PetscCall(MatDenseGetArrayRead(x->A,&px));
 30:   PetscCall(MatDenseGetArray(y->A,&py));
 31:   if (Q) {
 32:     PetscCall(MatDenseGetLDA(Q,&ldq));
 33:     PetscCall(MatDenseGetArrayRead(Q,&q));
 34:     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));
 35:     PetscCall(MatDenseRestoreArrayRead(Q,&q));
 36:   } 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));
 37:   PetscCall(MatDenseRestoreArrayRead(x->A,&px));
 38:   PetscCall(MatDenseRestoreArray(y->A,&py));
 39:   PetscFunctionReturn(PETSC_SUCCESS);
 40: }

 42: PetscErrorCode BVMultVec_Mat(BV X,PetscScalar alpha,PetscScalar beta,Vec y,PetscScalar *q)
 43: {
 44:   BV_MAT            *x = (BV_MAT*)X->data;
 45:   PetscScalar       *py,*qq=q;
 46:   const PetscScalar *px;

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

 59: PetscErrorCode BVMultInPlace_Mat(BV V,Mat Q,PetscInt s,PetscInt e)
 60: {
 61:   BV_MAT            *ctx = (BV_MAT*)V->data;
 62:   PetscScalar       *pv;
 63:   const PetscScalar *q;
 64:   PetscInt          ldq;

 66:   PetscFunctionBegin;
 67:   PetscCall(MatDenseGetLDA(Q,&ldq));
 68:   PetscCall(MatDenseGetArray(ctx->A,&pv));
 69:   PetscCall(MatDenseGetArrayRead(Q,&q));
 70:   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));
 71:   PetscCall(MatDenseRestoreArrayRead(Q,&q));
 72:   PetscCall(MatDenseRestoreArray(ctx->A,&pv));
 73:   PetscFunctionReturn(PETSC_SUCCESS);
 74: }

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

 83:   PetscFunctionBegin;
 84:   PetscCall(MatDenseGetLDA(Q,&ldq));
 85:   PetscCall(MatDenseGetArray(ctx->A,&pv));
 86:   PetscCall(MatDenseGetArrayRead(Q,&q));
 87:   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));
 88:   PetscCall(MatDenseRestoreArrayRead(Q,&q));
 89:   PetscCall(MatDenseRestoreArray(ctx->A,&pv));
 90:   PetscFunctionReturn(PETSC_SUCCESS);
 91: }

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

100:   PetscFunctionBegin;
101:   PetscCall(MatDenseGetLDA(M,&ldm));
102:   PetscCall(MatDenseGetArrayRead(x->A,&px));
103:   PetscCall(MatDenseGetArrayRead(y->A,&py));
104:   PetscCall(MatDenseGetArray(M,&m));
105:   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));
106:   PetscCall(MatDenseRestoreArray(M,&m));
107:   PetscCall(MatDenseRestoreArrayRead(x->A,&px));
108:   PetscCall(MatDenseRestoreArrayRead(y->A,&py));
109:   PetscFunctionReturn(PETSC_SUCCESS);
110: }

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

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

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

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

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

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

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

171:   PetscFunctionBegin;
172:   PetscCall(MatDenseGetArrayRead(ctx->A,&array));
173:   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));
174:   else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->n,type,val,ctx->mpi));
175:   PetscCall(MatDenseRestoreArrayRead(ctx->A,&array));
176:   PetscFunctionReturn(PETSC_SUCCESS);
177: }

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

184:   PetscFunctionBegin;
185:   PetscCall(MatDenseGetArrayRead(ctx->A,&array));
186:   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));
187:   else PetscCall(BVNorm_LAPACK_Private(bv,bv->n,1,array+(bv->nc+j)*bv->n,type,val,PETSC_FALSE));
188:   PetscCall(MatDenseRestoreArrayRead(ctx->A,&array));
189:   PetscFunctionReturn(PETSC_SUCCESS);
190: }

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

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

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

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

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

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

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

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

264: PetscErrorCode BVResize_Mat(BV bv,PetscInt m,PetscBool copy)
265: {
266:   BV_MAT            *ctx = (BV_MAT*)bv->data;
267:   PetscScalar       *pnew;
268:   const PetscScalar *pA;
269:   Mat               A;
270:   char              str[50];

272:   PetscFunctionBegin;
273:   PetscCall(MatCreateDense(PetscObjectComm((PetscObject)bv->t),bv->n,PETSC_DECIDE,PETSC_DECIDE,m,NULL,&A));
274:   if (((PetscObject)bv)->name) {
275:     PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
276:     PetscCall(PetscObjectSetName((PetscObject)A,str));
277:   }
278:   if (copy) {
279:     PetscCall(MatDenseGetArrayRead(ctx->A,&pA));
280:     PetscCall(MatDenseGetArrayWrite(A,&pnew));
281:     PetscCall(PetscArraycpy(pnew,pA,PetscMin(m,bv->m)*bv->n));
282:     PetscCall(MatDenseRestoreArrayRead(ctx->A,&pA));
283:     PetscCall(MatDenseRestoreArrayWrite(A,&pnew));
284:   }
285:   PetscCall(MatDestroy(&ctx->A));
286:   ctx->A = A;
287:   PetscFunctionReturn(PETSC_SUCCESS);
288: }

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

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

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

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

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

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

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

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

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

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

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

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

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

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

375: PetscErrorCode BVDestroy_Mat(BV bv)
376: {
377:   BV_MAT         *ctx = (BV_MAT*)bv->data;

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

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

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

400:   PetscCall(PetscObjectTypeCompare((PetscObject)bv->t,VECMPI,&ctx->mpi));
401:   if (!ctx->mpi) {
402:     PetscCall(PetscObjectTypeCompare((PetscObject)bv->t,VECSEQ,&seq));
403:     PetscCheck(seq,PetscObjectComm((PetscObject)bv),PETSC_ERR_SUP,"Cannot create a BVMAT from a non-standard template vector");
404:   }

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

409:   if (PetscUnlikely(bv->issplit)) {
410:     /* split BV: share the memory of the parent BV */
411:     parent = bv->splitparent;
412:     lsplit = parent->lsplit;
413:     PetscCall(MatDenseGetArray(((BV_MAT*)parent->data)->A,&array));
414:     ptr = (bv->issplit==1)? array: array+lsplit*nloc;
415:     PetscCall(MatDenseRestoreArray(((BV_MAT*)parent->data)->A,&array));
416:   } else {
417:     /* regular BV: allocate memory for the BV entries */
418:     ptr = NULL;
419:   }
420:   PetscCall(MatCreateDense(PetscObjectComm((PetscObject)bv->t),nloc,PETSC_DECIDE,PETSC_DECIDE,bv->m,ptr,&ctx->A));
421:   if (((PetscObject)bv)->name) {
422:     PetscCall(PetscSNPrintf(str,sizeof(str),"%s_0",((PetscObject)bv)->name));
423:     PetscCall(PetscObjectSetName((PetscObject)ctx->A,str));
424:   }

426:   if (PetscUnlikely(bv->Acreate)) {
427:     PetscCall(MatCopy(bv->Acreate,ctx->A,SAME_NONZERO_PATTERN));
428:     PetscCall(MatDestroy(&bv->Acreate));
429:   }

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

434:   bv->ops->mult             = BVMult_Mat;
435:   bv->ops->multvec          = BVMultVec_Mat;
436:   bv->ops->multinplace      = BVMultInPlace_Mat;
437:   bv->ops->multinplacetrans = BVMultInPlaceHermitianTranspose_Mat;
438:   bv->ops->dot              = BVDot_Mat;
439:   bv->ops->dotvec           = BVDotVec_Mat;
440:   bv->ops->dotvec_local     = BVDotVec_Local_Mat;
441:   bv->ops->scale            = BVScale_Mat;
442:   bv->ops->norm             = BVNorm_Mat;
443:   bv->ops->norm_local       = BVNorm_Local_Mat;
444:   bv->ops->normalize        = BVNormalize_Mat;
445:   bv->ops->matmult          = BVMatMult_Mat;
446:   bv->ops->copy             = BVCopy_Mat;
447:   bv->ops->copycolumn       = BVCopyColumn_Mat;
448:   bv->ops->resize           = BVResize_Mat;
449:   bv->ops->getcolumn        = BVGetColumn_Mat;
450:   bv->ops->restorecolumn    = BVRestoreColumn_Mat;
451:   bv->ops->getarray         = BVGetArray_Mat;
452:   bv->ops->restorearray     = BVRestoreArray_Mat;
453:   bv->ops->getarrayread     = BVGetArrayRead_Mat;
454:   bv->ops->restorearrayread = BVRestoreArrayRead_Mat;
455:   bv->ops->getmat           = BVGetMat_Default;
456:   bv->ops->restoremat       = BVRestoreMat_Default;
457:   bv->ops->destroy          = BVDestroy_Mat;
458:   if (!ctx->mpi) bv->ops->view = BVView_Mat;
459:   PetscFunctionReturn(PETSC_SUCCESS);
460: }