Actual source code: bvfunc.c

slepc-3.8.2 2017-12-01
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2017, 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 (basis vectors) interface routines, callable by users
 12: */

 14: #include <slepc/private/bvimpl.h>            /*I "slepcbv.h" I*/

 16: PetscClassId     BV_CLASSID = 0;
 17: PetscLogEvent    BV_Create = 0,BV_Copy = 0,BV_Mult = 0,BV_MultVec = 0,BV_MultInPlace = 0,BV_Dot = 0,BV_DotVec = 0,BV_Orthogonalize = 0,BV_OrthogonalizeVec = 0,BV_Scale = 0,BV_Norm = 0,BV_NormVec = 0,BV_SetRandom = 0,BV_MatMult = 0,BV_MatMultVec = 0,BV_MatProject = 0;
 18: static PetscBool BVPackageInitialized = PETSC_FALSE;

 20: const char *BVOrthogTypes[] = {"CGS","MGS","BVOrthogType","BV_ORTHOG_",0};
 21: const char *BVOrthogRefineTypes[] = {"IFNEEDED","NEVER","ALWAYS","BVOrthogRefineType","BV_ORTHOG_REFINE_",0};
 22: const char *BVOrthogBlockTypes[] = {"GS","CHOL","TSQR","BVOrthogBlockType","BV_ORTHOG_BLOCK_",0};
 23: const char *BVMatMultTypes[] = {"VECS","MAT","MAT_SAVE","BVMatMultType","BV_MATMULT_",0};

 25: /*@C
 26:    BVFinalizePackage - This function destroys everything in the Slepc interface
 27:    to the BV package. It is called from SlepcFinalize().

 29:    Level: developer

 31: .seealso: SlepcFinalize()
 32: @*/
 33: PetscErrorCode BVFinalizePackage(void)
 34: {

 38:   PetscFunctionListDestroy(&BVList);
 39:   BVPackageInitialized = PETSC_FALSE;
 40:   BVRegisterAllCalled  = PETSC_FALSE;
 41:   return(0);
 42: }

 44: /*@C
 45:    BVInitializePackage - This function initializes everything in the BV package.
 46:    It is called from PetscDLLibraryRegister() when using dynamic libraries, and
 47:    on the first call to BVCreate() when using static libraries.

 49:    Level: developer

 51: .seealso: SlepcInitialize()
 52: @*/
 53: PetscErrorCode BVInitializePackage(void)
 54: {
 55:   char           logList[256];
 56:   char           *className;
 57:   PetscBool      opt;

 61:   if (BVPackageInitialized) return(0);
 62:   BVPackageInitialized = PETSC_TRUE;
 63:   /* Register Classes */
 64:   PetscClassIdRegister("Basis Vectors",&BV_CLASSID);
 65:   /* Register Constructors */
 66:   BVRegisterAll();
 67:   /* Register Events */
 68:   PetscLogEventRegister("BVCreate",BV_CLASSID,&BV_Create);
 69:   PetscLogEventRegister("BVCopy",BV_CLASSID,&BV_Copy);
 70:   PetscLogEventRegister("BVMult",BV_CLASSID,&BV_Mult);
 71:   PetscLogEventRegister("BVMultVec",BV_CLASSID,&BV_MultVec);
 72:   PetscLogEventRegister("BVMultInPlace",BV_CLASSID,&BV_MultInPlace);
 73:   PetscLogEventRegister("BVDot",BV_CLASSID,&BV_Dot);
 74:   PetscLogEventRegister("BVDotVec",BV_CLASSID,&BV_DotVec);
 75:   PetscLogEventRegister("BVOrthogonalize",BV_CLASSID,&BV_Orthogonalize);
 76:   PetscLogEventRegister("BVOrthogonalizeV",BV_CLASSID,&BV_OrthogonalizeVec);
 77:   PetscLogEventRegister("BVScale",BV_CLASSID,&BV_Scale);
 78:   PetscLogEventRegister("BVNorm",BV_CLASSID,&BV_Norm);
 79:   PetscLogEventRegister("BVNormVec",BV_CLASSID,&BV_NormVec);
 80:   PetscLogEventRegister("BVSetRandom",BV_CLASSID,&BV_SetRandom);
 81:   PetscLogEventRegister("BVMatMult",BV_CLASSID,&BV_MatMult);
 82:   PetscLogEventRegister("BVMatMultVec",BV_CLASSID,&BV_MatMultVec);
 83:   PetscLogEventRegister("BVMatProject",BV_CLASSID,&BV_MatProject);
 84:   /* Process info exclusions */
 85:   PetscOptionsGetString(NULL,NULL,"-info_exclude",logList,256,&opt);
 86:   if (opt) {
 87:     PetscStrstr(logList,"bv",&className);
 88:     if (className) {
 89:       PetscInfoDeactivateClass(BV_CLASSID);
 90:     }
 91:   }
 92:   /* Process summary exclusions */
 93:   PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,256,&opt);
 94:   if (opt) {
 95:     PetscStrstr(logList,"bv",&className);
 96:     if (className) {
 97:       PetscLogEventDeactivateClass(BV_CLASSID);
 98:     }
 99:   }
100:   PetscRegisterFinalize(BVFinalizePackage);
101:   return(0);
102: }

104: /*@
105:    BVDestroy - Destroys BV context that was created with BVCreate().

107:    Collective on BV

109:    Input Parameter:
110: .  bv - the basis vectors context

112:    Level: beginner

114: .seealso: BVCreate()
115: @*/
116: PetscErrorCode BVDestroy(BV *bv)
117: {

121:   if (!*bv) return(0);
123:   if (--((PetscObject)(*bv))->refct > 0) { *bv = 0; return(0); }
124:   if ((*bv)->ops->destroy) { (*(*bv)->ops->destroy)(*bv); }
125:   VecDestroy(&(*bv)->t);
126:   MatDestroy(&(*bv)->matrix);
127:   VecDestroy(&(*bv)->Bx);
128:   VecDestroy(&(*bv)->buffer);
129:   BVDestroy(&(*bv)->cached);
130:   PetscFree((*bv)->work);
131:   PetscFree2((*bv)->h,(*bv)->c);
132:   VecDestroy(&(*bv)->omega);
133:   MatDestroy(&(*bv)->B);
134:   MatDestroy(&(*bv)->C);
135:   MatDestroy(&(*bv)->Acreate);
136:   PetscRandomDestroy(&(*bv)->rand);
137:   PetscHeaderDestroy(bv);
138:   return(0);
139: }

141: /*@
142:    BVCreate - Creates a basis vectors context.

144:    Collective on MPI_Comm

146:    Input Parameter:
147: .  comm - MPI communicator

149:    Output Parameter:
150: .  bv - location to put the basis vectors context

152:    Level: beginner

154: .seealso: BVSetUp(), BVDestroy(), BV
155: @*/
156: PetscErrorCode BVCreate(MPI_Comm comm,BV *newbv)
157: {
159:   BV             bv;

163:   *newbv = 0;
164:   BVInitializePackage();
165:   SlepcHeaderCreate(bv,BV_CLASSID,"BV","Basis Vectors","BV",comm,BVDestroy,BVView);

167:   bv->t            = NULL;
168:   bv->n            = -1;
169:   bv->N            = -1;
170:   bv->m            = 0;
171:   bv->l            = 0;
172:   bv->k            = 0;
173:   bv->nc           = 0;
174:   bv->orthog_type  = BV_ORTHOG_CGS;
175:   bv->orthog_ref   = BV_ORTHOG_REFINE_IFNEEDED;
176:   bv->orthog_eta   = 0.7071;
177:   bv->orthog_block = BV_ORTHOG_BLOCK_GS;
178:   bv->matrix       = NULL;
179:   bv->indef        = PETSC_FALSE;
180:   bv->vmm          = BV_MATMULT_MAT;

182:   bv->Bx           = NULL;
183:   bv->buffer       = NULL;
184:   bv->xid          = 0;
185:   bv->xstate       = 0;
186:   bv->cv[0]        = NULL;
187:   bv->cv[1]        = NULL;
188:   bv->ci[0]        = -1;
189:   bv->ci[1]        = -1;
190:   bv->st[0]        = -1;
191:   bv->st[1]        = -1;
192:   bv->id[0]        = 0;
193:   bv->id[1]        = 0;
194:   bv->h            = NULL;
195:   bv->c            = NULL;
196:   bv->omega        = NULL;
197:   bv->B            = NULL;
198:   bv->C            = NULL;
199:   bv->Aid          = 0;
200:   bv->defersfo     = PETSC_FALSE;
201:   bv->cached       = NULL;
202:   bv->bvstate      = 0;
203:   bv->rand         = NULL;
204:   bv->rrandom      = PETSC_FALSE;
205:   bv->Acreate      = NULL;
206:   bv->cuda         = PETSC_FALSE;
207:   bv->work         = NULL;
208:   bv->lwork        = 0;
209:   bv->data         = NULL;

211:   *newbv = bv;
212:   return(0);
213: }

215: /*@
216:    BVCreateFromMat - Creates a basis vectors object from a dense Mat object.

218:    Collective on Mat

220:    Input Parameter:
221: .  A - a dense tall-skinny matrix

223:    Output Parameter:
224: .  bv - the new basis vectors context

226:    Notes:
227:    The matrix values are copied to the BV data storage, memory is not shared.

229:    The communicator of the BV object will be the same as A, and so will be
230:    the dimensions.

232:    Level: intermediate

234: .seealso: BVCreate(), BVDestroy(), BVCreateMat()
235: @*/
236: PetscErrorCode BVCreateFromMat(Mat A,BV *bv)
237: {
239:   PetscBool      match;
240:   PetscInt       n,N,k;

244:   PetscObjectTypeCompareAny((PetscObject)A,&match,MATSEQDENSE,MATMPIDENSE,"");
245:   if (!match) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_SUP,"Matrix must be of type dense");

247:   MatGetSize(A,&N,&k);
248:   MatGetLocalSize(A,&n,NULL);
249:   BVCreate(PetscObjectComm((PetscObject)A),bv);
250:   BVSetSizes(*bv,n,N,k);

252:   (*bv)->Acreate = A;
253:   PetscObjectReference((PetscObject)A);
254:   return(0);
255: }

257: /*@
258:    BVInsertVec - Insert a vector into the specified column.

260:    Collective on BV

262:    Input Parameters:
263: +  V - basis vectors
264: .  j - the column of V to be overwritten
265: -  w - the vector to be copied

267:    Level: intermediate

269: .seealso: BVInsertVecs()
270: @*/
271: PetscErrorCode BVInsertVec(BV V,PetscInt j,Vec w)
272: {
274:   PetscInt       n,N;
275:   Vec            v;

282:   BVCheckSizes(V,1);

285:   VecGetSize(w,&N);
286:   VecGetLocalSize(w,&n);
287:   if (N!=V->N || n!=V->n) SETERRQ4(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_INCOMP,"Vec sizes (global %D, local %D) do not match BV sizes (global %D, local %D)",N,n,V->N,V->n);
288:   if (j<-V->nc || j>=V->m) SETERRQ3(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Argument j has wrong value %D, should be between %D and %D",j,-V->nc,V->m-1);

290:   BVGetColumn(V,j,&v);
291:   VecCopy(w,v);
292:   BVRestoreColumn(V,j,&v);
293:   PetscObjectStateIncrease((PetscObject)V);
294:   return(0);
295: }

297: /*@
298:    BVInsertVecs - Insert a set of vectors into the specified columns.

300:    Collective on BV

302:    Input Parameters:
303: +  V - basis vectors
304: .  s - first column of V to be overwritten
305: .  W - set of vectors to be copied
306: -  orth - flag indicating if the vectors must be orthogonalized

308:    Input/Output Parameter:
309: .  m - number of input vectors, on output the number of linearly independent
310:        vectors

312:    Notes:
313:    Copies the contents of vectors W to V(:,s:s+n). If the orthogonalization
314:    flag is set, then the vectors are copied one by one and then orthogonalized
315:    against the previous ones. If any of them is linearly dependent then it
316:    is discarded and the value of m is decreased.

318:    Level: intermediate

320: .seealso: BVInsertVec(), BVOrthogonalizeColumn()
321: @*/
322: PetscErrorCode BVInsertVecs(BV V,PetscInt s,PetscInt *m,Vec *W,PetscBool orth)
323: {
325:   PetscInt       n,N,i,ndep;
326:   PetscBool      lindep;
327:   PetscReal      norm;
328:   Vec            v;

335:   if (!*m) return(0);
336:   if (*m<0) SETERRQ1(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",*m);
341:   BVCheckSizes(V,1);

344:   VecGetSize(*W,&N);
345:   VecGetLocalSize(*W,&n);
346:   if (N!=V->N || n!=V->n) SETERRQ4(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_INCOMP,"Vec sizes (global %D, local %D) do not match BV sizes (global %D, local %D)",N,n,V->N,V->n);
347:   if (s<0 || s>=V->m) SETERRQ2(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Argument s has wrong value %D, should be between 0 and %D",s,V->m-1);
348:   if (s+(*m)>V->m) SETERRQ1(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Too many vectors provided, there is only room for %D",V->m);

350:   ndep = 0;
351:   for (i=0;i<*m;i++) {
352:     BVGetColumn(V,s+i-ndep,&v);
353:     VecCopy(W[i],v);
354:     BVRestoreColumn(V,s+i-ndep,&v);
355:     if (orth) {
356:       BVOrthogonalizeColumn(V,s+i-ndep,NULL,&norm,&lindep);
357:       if (norm==0.0 || lindep) {
358:         PetscInfo1(V,"Removing linearly dependent vector %D\n",i);
359:         ndep++;
360:       } else {
361:         BVScaleColumn(V,s+i-ndep,1.0/norm);
362:       }
363:     }
364:   }
365:   *m -= ndep;
366:   PetscObjectStateIncrease((PetscObject)V);
367:   return(0);
368: }

370: /*@
371:    BVInsertConstraints - Insert a set of vectors as constraints.

373:    Collective on BV

375:    Input Parameters:
376: +  V - basis vectors
377: -  C - set of vectors to be inserted as constraints

379:    Input/Output Parameter:
380: .  nc - number of input vectors, on output the number of linearly independent
381:        vectors

383:    Notes:
384:    The constraints are relevant only during orthogonalization. Constraint
385:    vectors span a subspace that is deflated in every orthogonalization
386:    operation, so they are intended for removing those directions from the
387:    orthogonal basis computed in regular BV columns.

389:    Constraints are not stored in regular BV colums, but in a special part of
390:    the storage. They can be accessed with negative indices in BVGetColumn().

392:    This operation is DESTRUCTIVE, meaning that all data contained in the
393:    columns of V is lost. This is typically invoked just after creating the BV.
394:    Once a set of constraints has been set, it is not allowed to call this
395:    function again.

397:    The vectors are copied one by one and then orthogonalized against the
398:    previous ones. If any of them is linearly dependent then it is discarded
399:    and the value of nc is decreased. The behaviour is similar to BVInsertVecs().

401:    Level: advanced

403: .seealso: BVInsertVecs(), BVOrthogonalizeColumn(), BVGetColumn(), BVGetNumConstraints()
404: @*/
405: PetscErrorCode BVInsertConstraints(BV V,PetscInt *nc,Vec *C)
406: {
408:   PetscInt       msave;

414:   if (!*nc) return(0);
415:   if (*nc<0) SETERRQ1(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of constraints (given %D) cannot be negative",*nc);
419:   BVCheckSizes(V,1);
421:   if (V->nc) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Constraints already present in this BV object");
422:   if (V->ci[0]!=-1 || V->ci[1]!=-1) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_SUP,"Cannot call BVInsertConstraints after BVGetColumn");

424:   msave = V->m;
425:   BVResize(V,*nc+V->m,PETSC_FALSE);
426:   BVInsertVecs(V,0,nc,C,PETSC_TRUE);
427:   V->nc = *nc;
428:   V->m  = msave;
429:   V->ci[0] = -V->nc-1;
430:   V->ci[1] = -V->nc-1;
431:   PetscObjectStateIncrease((PetscObject)V);
432:   return(0);
433: }

435: /*@C
436:    BVSetOptionsPrefix - Sets the prefix used for searching for all
437:    BV options in the database.

439:    Logically Collective on BV

441:    Input Parameters:
442: +  bv     - the basis vectors context
443: -  prefix - the prefix string to prepend to all BV option requests

445:    Notes:
446:    A hyphen (-) must NOT be given at the beginning of the prefix name.
447:    The first character of all runtime options is AUTOMATICALLY the
448:    hyphen.

450:    Level: advanced

452: .seealso: BVAppendOptionsPrefix(), BVGetOptionsPrefix()
453: @*/
454: PetscErrorCode BVSetOptionsPrefix(BV bv,const char *prefix)
455: {

460:   PetscObjectSetOptionsPrefix((PetscObject)bv,prefix);
461:   return(0);
462: }

464: /*@C
465:    BVAppendOptionsPrefix - Appends to the prefix used for searching for all
466:    BV options in the database.

468:    Logically Collective on BV

470:    Input Parameters:
471: +  bv     - the basis vectors context
472: -  prefix - the prefix string to prepend to all BV option requests

474:    Notes:
475:    A hyphen (-) must NOT be given at the beginning of the prefix name.
476:    The first character of all runtime options is AUTOMATICALLY the
477:    hyphen.

479:    Level: advanced

481: .seealso: BVSetOptionsPrefix(), BVGetOptionsPrefix()
482: @*/
483: PetscErrorCode BVAppendOptionsPrefix(BV bv,const char *prefix)
484: {

489:   PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix);
490:   return(0);
491: }

493: /*@C
494:    BVGetOptionsPrefix - Gets the prefix used for searching for all
495:    BV options in the database.

497:    Not Collective

499:    Input Parameters:
500: .  bv - the basis vectors context

502:    Output Parameters:
503: .  prefix - pointer to the prefix string used, is returned

505:    Note:
506:    On the Fortran side, the user should pass in a string 'prefix' of
507:    sufficient length to hold the prefix.

509:    Level: advanced

511: .seealso: BVSetOptionsPrefix(), BVAppendOptionsPrefix()
512: @*/
513: PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])
514: {

520:   PetscObjectGetOptionsPrefix((PetscObject)bv,prefix);
521:   return(0);
522: }

524: static PetscErrorCode BVView_Default(BV bv,PetscViewer viewer)
525: {
526:   PetscErrorCode    ierr;
527:   PetscInt          j;
528:   Vec               v;
529:   PetscViewerFormat format;
530:   PetscBool         isascii,ismatlab=PETSC_FALSE;
531:   const char        *bvname,*name;

534:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
535:   if (isascii) {
536:     PetscViewerGetFormat(viewer,&format);
537:     if (format == PETSC_VIEWER_ASCII_MATLAB) ismatlab = PETSC_TRUE;
538:   }
539:   if (ismatlab) {
540:     PetscObjectGetName((PetscObject)bv,&bvname);
541:     PetscViewerASCIIPrintf(viewer,"%s=[];\n",bvname);
542:   }
543:   for (j=-bv->nc;j<bv->m;j++) {
544:     BVGetColumn(bv,j,&v);
545:     VecView(v,viewer);
546:     if (ismatlab) {
547:       PetscObjectGetName((PetscObject)v,&name);
548:       PetscViewerASCIIPrintf(viewer,"%s=[%s,%s];clear %s\n",bvname,bvname,name,name);
549:     }
550:     BVRestoreColumn(bv,j,&v);
551:   }
552:   return(0);
553: }

555: /*@C
556:    BVView - Prints the BV data structure.

558:    Collective on BV

560:    Input Parameters:
561: +  bv     - the BV context
562: -  viewer - optional visualization context

564:    Note:
565:    The available visualization contexts include
566: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
567: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
568:          output where only the first processor opens
569:          the file.  All other processors send their
570:          data to the first processor to print.

572:    The user can open an alternative visualization contexts with
573:    PetscViewerASCIIOpen() (output to a specified file).

575:    Level: beginner

577: .seealso: PetscViewerASCIIOpen()
578: @*/
579: PetscErrorCode BVView(BV bv,PetscViewer viewer)
580: {
581:   PetscErrorCode    ierr;
582:   PetscBool         isascii;
583:   PetscViewerFormat format;
584:   const char        *orthname[2] = {"classical","modified"};
585:   const char        *refname[3] = {"if needed","never","always"};
586:   const char        *borthname[2] = {"Gram-Schmidt","Cholesky"};

590:   if (!viewer) {
591:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer);
592:   }

595:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
596:   if (isascii) {
597:     PetscObjectPrintClassNamePrefixType((PetscObject)bv,viewer);
598:     PetscViewerGetFormat(viewer,&format);
599:     PetscViewerASCIIPushTab(viewer);
600:     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
601:       PetscViewerASCIIPrintf(viewer,"%D columns of global length %D%s\n",bv->m,bv->N,bv->cuda?" (CUDA)":"");
602:       if (bv->nc>0) {
603:         PetscViewerASCIIPrintf(viewer,"number of constraints: %D\n",bv->nc);
604:       }
605:       PetscViewerASCIIPrintf(viewer,"vector orthogonalization method: %s Gram-Schmidt\n",orthname[bv->orthog_type]);
606:       switch (bv->orthog_ref) {
607:         case BV_ORTHOG_REFINE_IFNEEDED:
608:           PetscViewerASCIIPrintf(viewer,"orthogonalization refinement: %s (eta: %g)\n",refname[bv->orthog_ref],(double)bv->orthog_eta);
609:           break;
610:         case BV_ORTHOG_REFINE_NEVER:
611:         case BV_ORTHOG_REFINE_ALWAYS:
612:           PetscViewerASCIIPrintf(viewer,"orthogonalization refinement: %s\n",refname[bv->orthog_ref]);
613:           break;
614:       }
615:       PetscViewerASCIIPrintf(viewer,"block orthogonalization method: %s\n",borthname[bv->orthog_block]);
616:       if (bv->matrix) {
617:         if (bv->indef) {
618:           PetscViewerASCIIPrintf(viewer,"indefinite inner product\n");
619:         } else {
620:           PetscViewerASCIIPrintf(viewer,"non-standard inner product\n");
621:         }
622:         PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
623:         MatView(bv->matrix,viewer);
624:         PetscViewerPopFormat(viewer);
625:       }
626:       switch (bv->vmm) {
627:         case BV_MATMULT_VECS:
628:           PetscViewerASCIIPrintf(viewer,"doing matmult as matrix-vector products\n");
629:           break;
630:         case BV_MATMULT_MAT:
631:           PetscViewerASCIIPrintf(viewer,"doing matmult as a single matrix-matrix product\n");
632:           break;
633:         case BV_MATMULT_MAT_SAVE:
634:           PetscViewerASCIIPrintf(viewer,"doing matmult as a single matrix-matrix product, saving aux matrices\n");
635:           break;
636:       }
637:       if (bv->rrandom) {
638:         PetscViewerASCIIPrintf(viewer,"generating random vectors independent of the number of processes\n");
639:       }
640:     } else {
641:       if (bv->ops->view) { (*bv->ops->view)(bv,viewer); }
642:       else { BVView_Default(bv,viewer); }
643:     }
644:     PetscViewerASCIIPopTab(viewer);
645:   } else {
646:     (*bv->ops->view)(bv,viewer);
647:   }
648:   return(0);
649: }

651: /*@C
652:    BVRegister - Adds a new storage format to the BV package.

654:    Not collective

656:    Input Parameters:
657: +  name     - name of a new user-defined BV
658: -  function - routine to create context

660:    Notes:
661:    BVRegister() may be called multiple times to add several user-defined
662:    basis vectors.

664:    Level: advanced

666: .seealso: BVRegisterAll()
667: @*/
668: PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))
669: {

673:   PetscFunctionListAdd(&BVList,name,function);
674:   return(0);
675: }

677: PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)
678: {

682:   if (s>bv->lwork) {
683:     PetscFree(bv->work);
684:     PetscMalloc1(s,&bv->work);
685:     PetscLogObjectMemory((PetscObject)bv,(s-bv->lwork)*sizeof(PetscScalar));
686:     bv->lwork = s;
687:   }
688:   return(0);
689: }