Actual source code: bvfunc.c

slepc-3.16.1 2021-11-17
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, 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>

 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_Normalize = 0,BV_SetRandom = 0,BV_MatMult = 0,BV_MatMultVec = 0,BV_MatProject = 0,BV_SVDAndRank = 0;
 18: static PetscBool BVPackageInitialized = PETSC_FALSE;
 19: MPI_Op MPIU_TSQR = 0,MPIU_LAPY2;

 21: const char *BVOrthogTypes[] = {"CGS","MGS","BVOrthogType","BV_ORTHOG_",0};
 22: const char *BVOrthogRefineTypes[] = {"IFNEEDED","NEVER","ALWAYS","BVOrthogRefineType","BV_ORTHOG_REFINE_",0};
 23: const char *BVOrthogBlockTypes[] = {"GS","CHOL","TSQR","TSQRCHOL","SVQB","BVOrthogBlockType","BV_ORTHOG_BLOCK_",0};
 24: const char *BVMatMultTypes[] = {"VECS","MAT","MAT_SAVE","BVMatMultType","BV_MATMULT_",0};
 25: const char *BVSVDMethods[] = {"REFINE","QR","QR_CAA","BVSVDMethod","BV_SVD_METHOD_",0};

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

 31:    Level: developer

 33: .seealso: SlepcFinalize()
 34: @*/
 35: PetscErrorCode BVFinalizePackage(void)
 36: {

 40:   PetscFunctionListDestroy(&BVList);
 41:   MPI_Op_free(&MPIU_TSQR);
 42:   MPI_Op_free(&MPIU_LAPY2);
 43:   BVPackageInitialized = PETSC_FALSE;
 44:   BVRegisterAllCalled  = PETSC_FALSE;
 45:   return(0);
 46: }

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

 53:    Level: developer

 55: .seealso: SlepcInitialize()
 56: @*/
 57: PetscErrorCode BVInitializePackage(void)
 58: {
 59:   char           logList[256];
 60:   PetscBool      opt,pkg;
 61:   PetscClassId   classids[1];

 65:   if (BVPackageInitialized) return(0);
 66:   BVPackageInitialized = PETSC_TRUE;
 67:   /* Register Classes */
 68:   PetscClassIdRegister("Basis Vectors",&BV_CLASSID);
 69:   /* Register Constructors */
 70:   BVRegisterAll();
 71:   /* Register Events */
 72:   PetscLogEventRegister("BVCreate",BV_CLASSID,&BV_Create);
 73:   PetscLogEventRegister("BVCopy",BV_CLASSID,&BV_Copy);
 74:   PetscLogEventRegister("BVMult",BV_CLASSID,&BV_Mult);
 75:   PetscLogEventRegister("BVMultVec",BV_CLASSID,&BV_MultVec);
 76:   PetscLogEventRegister("BVMultInPlace",BV_CLASSID,&BV_MultInPlace);
 77:   PetscLogEventRegister("BVDot",BV_CLASSID,&BV_Dot);
 78:   PetscLogEventRegister("BVDotVec",BV_CLASSID,&BV_DotVec);
 79:   PetscLogEventRegister("BVOrthogonalize",BV_CLASSID,&BV_Orthogonalize);
 80:   PetscLogEventRegister("BVOrthogonalizeV",BV_CLASSID,&BV_OrthogonalizeVec);
 81:   PetscLogEventRegister("BVScale",BV_CLASSID,&BV_Scale);
 82:   PetscLogEventRegister("BVNorm",BV_CLASSID,&BV_Norm);
 83:   PetscLogEventRegister("BVNormVec",BV_CLASSID,&BV_NormVec);
 84:   PetscLogEventRegister("BVNormalize",BV_CLASSID,&BV_Normalize);
 85:   PetscLogEventRegister("BVSetRandom",BV_CLASSID,&BV_SetRandom);
 86:   PetscLogEventRegister("BVMatMult",BV_CLASSID,&BV_MatMult);
 87:   PetscLogEventRegister("BVMatMultVec",BV_CLASSID,&BV_MatMultVec);
 88:   PetscLogEventRegister("BVMatProject",BV_CLASSID,&BV_MatProject);
 89:   PetscLogEventRegister("BVSVDAndRank",BV_CLASSID,&BV_SVDAndRank);
 90:   /* MPI reduction operation used in BVOrthogonalize */
 91:   MPI_Op_create(SlepcGivensPacked,PETSC_FALSE,&MPIU_TSQR);
 92:   MPI_Op_create(SlepcPythag,PETSC_TRUE,&MPIU_LAPY2);
 93:   /* Process Info */
 94:   classids[0] = BV_CLASSID;
 95:   PetscInfoProcessClass("bv",1,&classids[0]);
 96:   /* Process summary exclusions */
 97:   PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt);
 98:   if (opt) {
 99:     PetscStrInList("bv",logList,',',&pkg);
100:     if (pkg) { PetscLogEventDeactivateClass(BV_CLASSID); }
101:   }
102:   /* Register package finalizer */
103:   PetscRegisterFinalize(BVFinalizePackage);
104:   return(0);
105: }

107: /*@C
108:    BVDestroy - Destroys BV context that was created with BVCreate().

110:    Collective on bv

112:    Input Parameter:
113: .  bv - the basis vectors context

115:    Level: beginner

117: .seealso: BVCreate()
118: @*/
119: PetscErrorCode BVDestroy(BV *bv)
120: {

124:   if (!*bv) return(0);
126:   if ((*bv)->lsplit) SETERRQ(PetscObjectComm((PetscObject)(*bv)),PETSC_ERR_ARG_WRONGSTATE,"Must call BVRestoreSplit before destroying the BV");
127:   if (--((PetscObject)(*bv))->refct > 0) { *bv = 0; return(0); }
128:   if ((*bv)->ops->destroy) { (*(*bv)->ops->destroy)(*bv); }
129:   VecDestroy(&(*bv)->t);
130:   MatDestroy(&(*bv)->matrix);
131:   VecDestroy(&(*bv)->Bx);
132:   VecDestroy(&(*bv)->buffer);
133:   BVDestroy(&(*bv)->cached);
134:   BVDestroy(&(*bv)->L);
135:   BVDestroy(&(*bv)->R);
136:   PetscFree((*bv)->work);
137:   PetscFree2((*bv)->h,(*bv)->c);
138:   VecDestroy(&(*bv)->omega);
139:   MatDestroy(&(*bv)->Acreate);
140:   MatDestroy(&(*bv)->Aget);
141:   MatDestroy(&(*bv)->Abuffer);
142:   PetscRandomDestroy(&(*bv)->rand);
143:   PetscHeaderDestroy(bv);
144:   return(0);
145: }

147: /*@
148:    BVCreate - Creates a basis vectors context.

150:    Collective

152:    Input Parameter:
153: .  comm - MPI communicator

155:    Output Parameter:
156: .  bv - location to put the basis vectors context

158:    Level: beginner

160: .seealso: BVSetUp(), BVDestroy(), BV
161: @*/
162: PetscErrorCode BVCreate(MPI_Comm comm,BV *newbv)
163: {
165:   BV             bv;

169:   *newbv = 0;
170:   BVInitializePackage();
171:   SlepcHeaderCreate(bv,BV_CLASSID,"BV","Basis Vectors","BV",comm,BVDestroy,BVView);

173:   bv->t            = NULL;
174:   bv->n            = -1;
175:   bv->N            = -1;
176:   bv->m            = 0;
177:   bv->l            = 0;
178:   bv->k            = 0;
179:   bv->nc           = 0;
180:   bv->orthog_type  = BV_ORTHOG_CGS;
181:   bv->orthog_ref   = BV_ORTHOG_REFINE_IFNEEDED;
182:   bv->orthog_eta   = 0.7071;
183:   bv->orthog_block = BV_ORTHOG_BLOCK_GS;
184:   bv->matrix       = NULL;
185:   bv->indef        = PETSC_FALSE;
186:   bv->vmm          = BV_MATMULT_MAT;
187:   bv->rrandom      = PETSC_FALSE;
188:   bv->deftol       = 10*PETSC_MACHINE_EPSILON;

190:   bv->Bx           = NULL;
191:   bv->buffer       = NULL;
192:   bv->Abuffer      = NULL;
193:   bv->xid          = 0;
194:   bv->xstate       = 0;
195:   bv->cv[0]        = NULL;
196:   bv->cv[1]        = NULL;
197:   bv->ci[0]        = -1;
198:   bv->ci[1]        = -1;
199:   bv->st[0]        = -1;
200:   bv->st[1]        = -1;
201:   bv->id[0]        = 0;
202:   bv->id[1]        = 0;
203:   bv->h            = NULL;
204:   bv->c            = NULL;
205:   bv->omega        = NULL;
206:   bv->defersfo     = PETSC_FALSE;
207:   bv->cached       = NULL;
208:   bv->bvstate      = 0;
209:   bv->L            = NULL;
210:   bv->R            = NULL;
211:   bv->lstate       = 0;
212:   bv->rstate       = 0;
213:   bv->lsplit       = 0;
214:   bv->issplit      = 0;
215:   bv->splitparent  = NULL;
216:   bv->rand         = NULL;
217:   bv->rrandom      = PETSC_FALSE;
218:   bv->Acreate      = NULL;
219:   bv->Aget         = NULL;
220:   bv->cuda         = PETSC_FALSE;
221:   bv->sfocalled    = PETSC_FALSE;
222:   bv->work         = NULL;
223:   bv->lwork        = 0;
224:   bv->data         = NULL;

226:   *newbv = bv;
227:   return(0);
228: }

230: /*@
231:    BVCreateFromMat - Creates a basis vectors object from a dense Mat object.

233:    Collective on A

235:    Input Parameter:
236: .  A - a dense tall-skinny matrix

238:    Output Parameter:
239: .  bv - the new basis vectors context

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

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

247:    Level: intermediate

249: .seealso: BVCreate(), BVDestroy(), BVCreateMat()
250: @*/
251: PetscErrorCode BVCreateFromMat(Mat A,BV *bv)
252: {
254:   PetscInt       n,N,k;


260:   MatGetSize(A,&N,&k);
261:   MatGetLocalSize(A,&n,NULL);
262:   BVCreate(PetscObjectComm((PetscObject)A),bv);
263:   BVSetSizes(*bv,n,N,k);

265:   (*bv)->Acreate = A;
266:   PetscObjectReference((PetscObject)A);
267:   return(0);
268: }

270: /*@
271:    BVInsertVec - Insert a vector into the specified column.

273:    Collective on V

275:    Input Parameters:
276: +  V - basis vectors
277: .  j - the column of V to be overwritten
278: -  w - the vector to be copied

280:    Level: intermediate

282: .seealso: BVInsertVecs()
283: @*/
284: PetscErrorCode BVInsertVec(BV V,PetscInt j,Vec w)
285: {
287:   PetscInt       n,N;
288:   Vec            v;

295:   BVCheckSizes(V,1);

298:   VecGetSize(w,&N);
299:   VecGetLocalSize(w,&n);
300:   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);
301:   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);

303:   BVGetColumn(V,j,&v);
304:   VecCopy(w,v);
305:   BVRestoreColumn(V,j,&v);
306:   PetscObjectStateIncrease((PetscObject)V);
307:   return(0);
308: }

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

313:    Collective on V

315:    Input Parameters:
316: +  V - basis vectors
317: .  s - first column of V to be overwritten
318: .  W - set of vectors to be copied
319: -  orth - flag indicating if the vectors must be orthogonalized

321:    Input/Output Parameter:
322: .  m - number of input vectors, on output the number of linearly independent
323:        vectors

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

331:    Level: intermediate

333: .seealso: BVInsertVec(), BVOrthogonalizeColumn()
334: @*/
335: PetscErrorCode BVInsertVecs(BV V,PetscInt s,PetscInt *m,Vec *W,PetscBool orth)
336: {
338:   PetscInt       n,N,i,ndep;
339:   PetscBool      lindep;
340:   PetscReal      norm;
341:   Vec            v;

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

357:   VecGetSize(*W,&N);
358:   VecGetLocalSize(*W,&n);
359:   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);
360:   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);
361:   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);

363:   ndep = 0;
364:   for (i=0;i<*m;i++) {
365:     BVGetColumn(V,s+i-ndep,&v);
366:     VecCopy(W[i],v);
367:     BVRestoreColumn(V,s+i-ndep,&v);
368:     if (orth) {
369:       BVOrthogonalizeColumn(V,s+i-ndep,NULL,&norm,&lindep);
370:       if (norm==0.0 || lindep) {
371:         PetscInfo1(V,"Removing linearly dependent vector %D\n",i);
372:         ndep++;
373:       } else {
374:         BVScaleColumn(V,s+i-ndep,1.0/norm);
375:       }
376:     }
377:   }
378:   *m -= ndep;
379:   PetscObjectStateIncrease((PetscObject)V);
380:   return(0);
381: }

383: /*@
384:    BVInsertConstraints - Insert a set of vectors as constraints.

386:    Collective on V

388:    Input Parameters:
389: +  V - basis vectors
390: -  C - set of vectors to be inserted as constraints

392:    Input/Output Parameter:
393: .  nc - number of input vectors, on output the number of linearly independent
394:        vectors

396:    Notes:
397:    The constraints are relevant only during orthogonalization. Constraint
398:    vectors span a subspace that is deflated in every orthogonalization
399:    operation, so they are intended for removing those directions from the
400:    orthogonal basis computed in regular BV columns.

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

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

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

414:    Level: advanced

416: .seealso: BVInsertVecs(), BVOrthogonalizeColumn(), BVGetColumn(), BVGetNumConstraints()
417: @*/
418: PetscErrorCode BVInsertConstraints(BV V,PetscInt *nc,Vec *C)
419: {
421:   PetscInt       msave;

427:   if (!*nc) return(0);
428:   if (*nc<0) SETERRQ1(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of constraints (given %D) cannot be negative",*nc);
432:   BVCheckSizes(V,1);
434:   if (V->issplit) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Operation not permitted for a BV obtained from BVGetSplit");
435:   if (V->nc) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Constraints already present in this BV object");
436:   if (V->ci[0]!=-1 || V->ci[1]!=-1) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_SUP,"Cannot call BVInsertConstraints after BVGetColumn");

438:   msave = V->m;
439:   BVResize(V,*nc+V->m,PETSC_FALSE);
440:   BVInsertVecs(V,0,nc,C,PETSC_TRUE);
441:   V->nc = *nc;
442:   V->m  = msave;
443:   V->ci[0] = -V->nc-1;
444:   V->ci[1] = -V->nc-1;
445:   PetscObjectStateIncrease((PetscObject)V);
446:   return(0);
447: }

449: /*@C
450:    BVSetOptionsPrefix - Sets the prefix used for searching for all
451:    BV options in the database.

453:    Logically Collective on bv

455:    Input Parameters:
456: +  bv     - the basis vectors context
457: -  prefix - the prefix string to prepend to all BV option requests

459:    Notes:
460:    A hyphen (-) must NOT be given at the beginning of the prefix name.
461:    The first character of all runtime options is AUTOMATICALLY the
462:    hyphen.

464:    Level: advanced

466: .seealso: BVAppendOptionsPrefix(), BVGetOptionsPrefix()
467: @*/
468: PetscErrorCode BVSetOptionsPrefix(BV bv,const char *prefix)
469: {

474:   PetscObjectSetOptionsPrefix((PetscObject)bv,prefix);
475:   return(0);
476: }

478: /*@C
479:    BVAppendOptionsPrefix - Appends to the prefix used for searching for all
480:    BV options in the database.

482:    Logically Collective on bv

484:    Input Parameters:
485: +  bv     - the basis vectors context
486: -  prefix - the prefix string to prepend to all BV option requests

488:    Notes:
489:    A hyphen (-) must NOT be given at the beginning of the prefix name.
490:    The first character of all runtime options is AUTOMATICALLY the
491:    hyphen.

493:    Level: advanced

495: .seealso: BVSetOptionsPrefix(), BVGetOptionsPrefix()
496: @*/
497: PetscErrorCode BVAppendOptionsPrefix(BV bv,const char *prefix)
498: {

503:   PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix);
504:   return(0);
505: }

507: /*@C
508:    BVGetOptionsPrefix - Gets the prefix used for searching for all
509:    BV options in the database.

511:    Not Collective

513:    Input Parameters:
514: .  bv - the basis vectors context

516:    Output Parameters:
517: .  prefix - pointer to the prefix string used, is returned

519:    Note:
520:    On the Fortran side, the user should pass in a string 'prefix' of
521:    sufficient length to hold the prefix.

523:    Level: advanced

525: .seealso: BVSetOptionsPrefix(), BVAppendOptionsPrefix()
526: @*/
527: PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])
528: {

534:   PetscObjectGetOptionsPrefix((PetscObject)bv,prefix);
535:   return(0);
536: }

538: static PetscErrorCode BVView_Default(BV bv,PetscViewer viewer)
539: {
540:   PetscErrorCode    ierr;
541:   PetscInt          j;
542:   Vec               v;
543:   PetscViewerFormat format;
544:   PetscBool         isascii,ismatlab=PETSC_FALSE;
545:   const char        *bvname,*name;

548:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
549:   if (isascii) {
550:     PetscViewerGetFormat(viewer,&format);
551:     if (format == PETSC_VIEWER_ASCII_MATLAB) ismatlab = PETSC_TRUE;
552:   }
553:   if (ismatlab) {
554:     PetscObjectGetName((PetscObject)bv,&bvname);
555:     PetscViewerASCIIPrintf(viewer,"%s=[];\n",bvname);
556:   }
557:   for (j=-bv->nc;j<bv->m;j++) {
558:     BVGetColumn(bv,j,&v);
559:     VecView(v,viewer);
560:     if (ismatlab) {
561:       PetscObjectGetName((PetscObject)v,&name);
562:       PetscViewerASCIIPrintf(viewer,"%s=[%s,%s];clear %s\n",bvname,bvname,name,name);
563:     }
564:     BVRestoreColumn(bv,j,&v);
565:   }
566:   return(0);
567: }

569: /*@C
570:    BVView - Prints the BV data structure.

572:    Collective on bv

574:    Input Parameters:
575: +  bv     - the BV context
576: -  viewer - optional visualization context

578:    Note:
579:    The available visualization contexts include
580: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
581: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
582:          output where only the first processor opens
583:          the file.  All other processors send their
584:          data to the first processor to print.

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

589:    Level: beginner
590: @*/
591: PetscErrorCode BVView(BV bv,PetscViewer viewer)
592: {
593:   PetscErrorCode    ierr;
594:   PetscBool         isascii;
595:   PetscViewerFormat format;
596:   const char        *orthname[2] = {"classical","modified"};
597:   const char        *refname[3] = {"if needed","never","always"};

601:   if (!viewer) {
602:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer);
603:   }

606:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
607:   if (isascii) {
608:     PetscObjectPrintClassNamePrefixType((PetscObject)bv,viewer);
609:     PetscViewerGetFormat(viewer,&format);
610:     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
611:       PetscViewerASCIIPrintf(viewer,"  %D columns of global length %D%s\n",bv->m,bv->N,bv->cuda?" (CUDA)":"");
612:       if (bv->nc>0) {
613:         PetscViewerASCIIPrintf(viewer,"  number of constraints: %D\n",bv->nc);
614:       }
615:       PetscViewerASCIIPrintf(viewer,"  vector orthogonalization method: %s Gram-Schmidt\n",orthname[bv->orthog_type]);
616:       switch (bv->orthog_ref) {
617:         case BV_ORTHOG_REFINE_IFNEEDED:
618:           PetscViewerASCIIPrintf(viewer,"  orthogonalization refinement: %s (eta: %g)\n",refname[bv->orthog_ref],(double)bv->orthog_eta);
619:           break;
620:         case BV_ORTHOG_REFINE_NEVER:
621:         case BV_ORTHOG_REFINE_ALWAYS:
622:           PetscViewerASCIIPrintf(viewer,"  orthogonalization refinement: %s\n",refname[bv->orthog_ref]);
623:           break;
624:       }
625:       PetscViewerASCIIPrintf(viewer,"  block orthogonalization method: %s\n",BVOrthogBlockTypes[bv->orthog_block]);
626:       if (bv->matrix) {
627:         if (bv->indef) {
628:           PetscViewerASCIIPrintf(viewer,"  indefinite inner product\n");
629:         } else {
630:           PetscViewerASCIIPrintf(viewer,"  non-standard inner product\n");
631:         }
632:         PetscViewerASCIIPrintf(viewer,"  tolerance for definite inner product: %g\n",(double)bv->deftol);
633:         PetscViewerASCIIPrintf(viewer,"  inner product matrix:\n");
634:         PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
635:         PetscViewerASCIIPushTab(viewer);
636:         MatView(bv->matrix,viewer);
637:         PetscViewerASCIIPopTab(viewer);
638:         PetscViewerPopFormat(viewer);
639:       }
640:       switch (bv->vmm) {
641:         case BV_MATMULT_VECS:
642:           PetscViewerASCIIPrintf(viewer,"  doing matmult as matrix-vector products\n");
643:           break;
644:         case BV_MATMULT_MAT:
645:           PetscViewerASCIIPrintf(viewer,"  doing matmult as a single matrix-matrix product\n");
646:           break;
647:         case BV_MATMULT_MAT_SAVE:
648:           PetscViewerASCIIPrintf(viewer,"  mat_save is deprecated, use mat\n");
649:           break;
650:       }
651:       if (bv->rrandom) {
652:         PetscViewerASCIIPrintf(viewer,"  generating random vectors independent of the number of processes\n");
653:       }
654:       if (bv->ops->view) { (*bv->ops->view)(bv,viewer); }
655:     } else {
656:       if (bv->ops->view) { (*bv->ops->view)(bv,viewer); }
657:       else { BVView_Default(bv,viewer); }
658:     }
659:   } else {
660:     (*bv->ops->view)(bv,viewer);
661:   }
662:   return(0);
663: }

665: /*@C
666:    BVViewFromOptions - View from options

668:    Collective on BV

670:    Input Parameters:
671: +  bv   - the basis vectors context
672: .  obj  - optional object
673: -  name - command line option

675:    Level: intermediate

677: .seealso: BVView(), BVCreate()
678: @*/
679: PetscErrorCode BVViewFromOptions(BV bv,PetscObject obj,const char name[])
680: {

685:   PetscObjectViewFromOptions((PetscObject)bv,obj,name);
686:   return(0);
687: }

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

692:    Not collective

694:    Input Parameters:
695: +  name     - name of a new user-defined BV
696: -  function - routine to create context

698:    Notes:
699:    BVRegister() may be called multiple times to add several user-defined
700:    basis vectors.

702:    Level: advanced

704: .seealso: BVRegisterAll()
705: @*/
706: PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))
707: {

711:   BVInitializePackage();
712:   PetscFunctionListAdd(&BVList,name,function);
713:   return(0);
714: }

716: PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)
717: {

721:   if (s>bv->lwork) {
722:     PetscFree(bv->work);
723:     PetscMalloc1(s,&bv->work);
724:     PetscLogObjectMemory((PetscObject)bv,(s-bv->lwork)*sizeof(PetscScalar));
725:     bv->lwork = s;
726:   }
727:   return(0);
728: }

730: #if defined(PETSC_USE_DEBUG)
731: /*
732:    SlepcDebugBVView - partially view a BV object, to be used from within a debugger.

734:      ini, end: columns to be viewed
735:      s: name of Matlab variable
736:      filename: optionally write output to a file
737:  */
738: PETSC_UNUSED PetscErrorCode SlepcDebugBVView(BV bv,PetscInt ini,PetscInt end,const char *s,const char *filename)
739: {
741:   PetscInt       N,m;
742:   PetscScalar    *array;

745:   BVGetArray(bv,&array);
746:   BVGetSizes(bv,NULL,&N,&m);
747:   SlepcDebugViewMatrix(N,end-ini+1,array+ini*N,NULL,N,s,filename);
748:   BVRestoreArray(bv,&array);
749:   return(0);
750: }
751: #endif