Actual source code: bvfunc.c

slepc-main 2024-11-09
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 (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_",NULL};
 22: const char *BVOrthogRefineTypes[] = {"IFNEEDED","NEVER","ALWAYS","BVOrthogRefineType","BV_ORTHOG_REFINE_",NULL};
 23: const char *BVOrthogBlockTypes[] = {"GS","CHOL","TSQR","TSQRCHOL","SVQB","BVOrthogBlockType","BV_ORTHOG_BLOCK_",NULL};
 24: const char *BVMatMultTypes[] = {"VECS","MAT","MAT_SAVE","BVMatMultType","BV_MATMULT_",NULL};
 25: const char *BVSVDMethods[] = {"REFINE","QR","QR_CAA","BVSVDMethod","BV_SVD_METHOD_",NULL};

 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: {
 37:   PetscFunctionBegin;
 38:   PetscCall(PetscFunctionListDestroy(&BVList));
 39:   PetscCallMPI(MPI_Op_free(&MPIU_TSQR));
 40:   PetscCallMPI(MPI_Op_free(&MPIU_LAPY2));
 41:   BVPackageInitialized = PETSC_FALSE;
 42:   BVRegisterAllCalled  = PETSC_FALSE;
 43:   PetscFunctionReturn(PETSC_SUCCESS);
 44: }

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

 51:    Level: developer

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

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

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

107:    Collective

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

112:    Level: beginner

114: .seealso: BVCreate()
115: @*/
116: PetscErrorCode BVDestroy(BV *bv)
117: {
118:   PetscFunctionBegin;
119:   if (!*bv) PetscFunctionReturn(PETSC_SUCCESS);
121:   PetscCheck(!(*bv)->lsplit,PetscObjectComm((PetscObject)*bv),PETSC_ERR_ARG_WRONGSTATE,"Must call BVRestoreSplit before destroying the BV");
122:   if (--((PetscObject)*bv)->refct > 0) { *bv = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
123:   PetscTryTypeMethod(*bv,destroy);
124:   PetscCall(PetscLayoutDestroy(&(*bv)->map));
125:   PetscCall(PetscFree((*bv)->vtype));
126:   PetscCall(MatDestroy(&(*bv)->matrix));
127:   PetscCall(VecDestroy(&(*bv)->Bx));
128:   PetscCall(VecDestroy(&(*bv)->buffer));
129:   PetscCall(BVDestroy(&(*bv)->cached));
130:   PetscCall(BVDestroy(&(*bv)->L));
131:   PetscCall(BVDestroy(&(*bv)->R));
132:   PetscCall(PetscFree((*bv)->work));
133:   PetscCall(PetscFree2((*bv)->h,(*bv)->c));
134:   PetscCall(VecDestroy(&(*bv)->omega));
135:   PetscCall(MatDestroy(&(*bv)->Acreate));
136:   PetscCall(MatDestroy(&(*bv)->Aget));
137:   PetscCall(MatDestroy(&(*bv)->Abuffer));
138:   PetscCall(PetscRandomDestroy(&(*bv)->rand));
139:   PetscCall(PetscHeaderDestroy(bv));
140:   PetscFunctionReturn(PETSC_SUCCESS);
141: }

143: /*@
144:    BVCreate - Creates a basis vectors context.

146:    Collective

148:    Input Parameter:
149: .  comm - MPI communicator

151:    Output Parameter:
152: .  newbv - location to put the basis vectors context

154:    Level: beginner

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

162:   PetscFunctionBegin;
163:   PetscAssertPointer(newbv,2);
164:   PetscCall(BVInitializePackage());
165:   PetscCall(SlepcHeaderCreate(bv,BV_CLASSID,"BV","Basis Vectors","BV",comm,BVDestroy,BVView));

167:   bv->map          = NULL;
168:   bv->vtype        = NULL;
169:   bv->n            = -1;
170:   bv->N            = -1;
171:   bv->m            = 0;
172:   bv->l            = 0;
173:   bv->k            = 0;
174:   bv->nc           = 0;
175:   bv->ld           = 0;
176:   bv->orthog_type  = BV_ORTHOG_CGS;
177:   bv->orthog_ref   = BV_ORTHOG_REFINE_IFNEEDED;
178:   bv->orthog_eta   = 0.7071;
179:   bv->orthog_block = BV_ORTHOG_BLOCK_GS;
180:   bv->matrix       = NULL;
181:   bv->indef        = PETSC_FALSE;
182:   bv->vmm          = BV_MATMULT_MAT;
183:   bv->rrandom      = PETSC_FALSE;
184:   bv->deftol       = 10*PETSC_MACHINE_EPSILON;

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

222:   *newbv = bv;
223:   PetscFunctionReturn(PETSC_SUCCESS);
224: }

226: /*@
227:    BVCreateFromMat - Creates a basis vectors object from a dense Mat object.

229:    Collective

231:    Input Parameter:
232: .  A - a dense tall-skinny matrix

234:    Output Parameter:
235: .  bv - the new basis vectors context

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

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

243:    Level: intermediate

245: .seealso: BVCreate(), BVDestroy(), BVCreateMat()
246: @*/
247: PetscErrorCode BVCreateFromMat(Mat A,BV *bv)
248: {
249:   PetscInt  n,N,k;
250:   VecType   vtype;

252:   PetscFunctionBegin;

255:   PetscCall(MatGetSize(A,&N,&k));
256:   PetscCall(MatGetLocalSize(A,&n,NULL));
257:   PetscCall(MatGetVecType(A,&vtype));
258:   PetscCall(BVCreate(PetscObjectComm((PetscObject)A),bv));
259:   PetscCall(BVSetSizes(*bv,n,N,k));
260:   PetscCall(BVSetVecType(*bv,vtype));

262:   (*bv)->Acreate = A;
263:   PetscCall(PetscObjectReference((PetscObject)A));
264:   PetscFunctionReturn(PETSC_SUCCESS);
265: }

267: /*@
268:    BVInsertVec - Insert a vector into the specified column.

270:    Logically Collective

272:    Input Parameters:
273: +  V - basis vectors
274: .  j - the column of V to be overwritten
275: -  w - the vector to be copied

277:    Level: intermediate

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

286:   PetscFunctionBegin;
291:   BVCheckSizes(V,1);
292:   PetscCheckSameComm(V,1,w,3);

294:   PetscCall(VecGetSize(w,&N));
295:   PetscCall(VecGetLocalSize(w,&n));
296:   PetscCheck(N==V->N && n==V->n,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_INCOMP,"Vec sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ") do not match BV sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ")",N,n,V->N,V->n);
297:   PetscCheck(j>=-V->nc && j<V->m,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Argument j has wrong value %" PetscInt_FMT ", should be between %" PetscInt_FMT " and %" PetscInt_FMT,j,-V->nc,V->m-1);

299:   PetscCall(BVGetColumn(V,j,&v));
300:   PetscCall(VecCopy(w,v));
301:   PetscCall(BVRestoreColumn(V,j,&v));
302:   PetscCall(PetscObjectStateIncrease((PetscObject)V));
303:   PetscFunctionReturn(PETSC_SUCCESS);
304: }

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

309:    Collective

311:    Input Parameters:
312: +  V - basis vectors
313: .  s - first column of V to be overwritten
314: .  W - set of vectors to be copied
315: -  orth - flag indicating if the vectors must be orthogonalized

317:    Input/Output Parameter:
318: .  m - number of input vectors, on output the number of linearly independent
319:        vectors

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

327:    Level: intermediate

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

338:   PetscFunctionBegin;
341:   PetscAssertPointer(m,3);
343:   if (!*m) PetscFunctionReturn(PETSC_SUCCESS);
344:   PetscCheck(*m>0,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %" PetscInt_FMT ") cannot be negative",*m);
345:   PetscAssertPointer(W,4);
349:   BVCheckSizes(V,1);
350:   PetscCheckSameComm(V,1,*W,4);

352:   PetscCall(VecGetSize(*W,&N));
353:   PetscCall(VecGetLocalSize(*W,&n));
354:   PetscCheck(N==V->N && n==V->n,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_INCOMP,"Vec sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ") do not match BV sizes (global %" PetscInt_FMT ", local %" PetscInt_FMT ")",N,n,V->N,V->n);
355:   PetscCheck(s>=0 && s<V->m,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Argument s has wrong value %" PetscInt_FMT ", should be between 0 and %" PetscInt_FMT,s,V->m-1);
356:   PetscCheck(s+(*m)<=V->m,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Too many vectors provided, there is only room for %" PetscInt_FMT,V->m);

358:   ndep = 0;
359:   for (i=0;i<*m;i++) {
360:     PetscCall(BVGetColumn(V,s+i-ndep,&v));
361:     PetscCall(VecCopy(W[i],v));
362:     PetscCall(BVRestoreColumn(V,s+i-ndep,&v));
363:     if (orth) {
364:       PetscCall(BVOrthogonalizeColumn(V,s+i-ndep,NULL,&norm,&lindep));
365:       if (norm==0.0 || lindep) {
366:         PetscCall(PetscInfo(V,"Removing linearly dependent vector %" PetscInt_FMT "\n",i));
367:         ndep++;
368:       } else PetscCall(BVScaleColumn(V,s+i-ndep,1.0/norm));
369:     }
370:   }
371:   *m -= ndep;
372:   PetscCall(PetscObjectStateIncrease((PetscObject)V));
373:   PetscFunctionReturn(PETSC_SUCCESS);
374: }

376: /*@
377:    BVInsertConstraints - Insert a set of vectors as constraints.

379:    Collective

381:    Input Parameters:
382: +  V - basis vectors
383: -  C - set of vectors to be inserted as constraints

385:    Input/Output Parameter:
386: .  nc - number of input vectors, on output the number of linearly independent
387:        vectors

389:    Notes:
390:    The constraints are relevant only during orthogonalization. Constraint
391:    vectors span a subspace that is deflated in every orthogonalization
392:    operation, so they are intended for removing those directions from the
393:    orthogonal basis computed in regular BV columns.

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

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

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

407:    Level: advanced

409: .seealso: BVInsertVecs(), BVOrthogonalizeColumn(), BVGetColumn(), BVGetNumConstraints()
410: @*/
411: PetscErrorCode BVInsertConstraints(BV V,PetscInt *nc,Vec *C)
412: {
413:   PetscInt       msave;

415:   PetscFunctionBegin;
417:   PetscAssertPointer(nc,2);
419:   if (!*nc) PetscFunctionReturn(PETSC_SUCCESS);
420:   PetscCheck(*nc>0,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of constraints (given %" PetscInt_FMT ") cannot be negative",*nc);
421:   PetscAssertPointer(C,3);
424:   BVCheckSizes(V,1);
425:   PetscCheckSameComm(V,1,*C,3);
426:   PetscCheck(!V->issplit,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Operation not permitted for a BV obtained from BVGetSplit");
427:   PetscCheck(!V->nc,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Constraints already present in this BV object");
428:   PetscCheck(V->ci[0]==-1 && V->ci[1]==-1,PetscObjectComm((PetscObject)V),PETSC_ERR_SUP,"Cannot call BVInsertConstraints after BVGetColumn");

430:   msave = V->m;
431:   PetscCall(BVResize(V,*nc+V->m,PETSC_FALSE));
432:   PetscCall(BVInsertVecs(V,0,nc,C,PETSC_TRUE));
433:   V->nc = *nc;
434:   V->m  = msave;
435:   V->ci[0] = -V->nc-1;
436:   V->ci[1] = -V->nc-1;
437:   PetscCall(PetscObjectStateIncrease((PetscObject)V));
438:   PetscFunctionReturn(PETSC_SUCCESS);
439: }

441: /*@
442:    BVSetOptionsPrefix - Sets the prefix used for searching for all
443:    BV options in the database.

445:    Logically Collective

447:    Input Parameters:
448: +  bv     - the basis vectors context
449: -  prefix - the prefix string to prepend to all BV option requests

451:    Notes:
452:    A hyphen (-) must NOT be given at the beginning of the prefix name.
453:    The first character of all runtime options is AUTOMATICALLY the
454:    hyphen.

456:    Level: advanced

458: .seealso: BVAppendOptionsPrefix(), BVGetOptionsPrefix()
459: @*/
460: PetscErrorCode BVSetOptionsPrefix(BV bv,const char *prefix)
461: {
462:   PetscFunctionBegin;
464:   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)bv,prefix));
465:   PetscFunctionReturn(PETSC_SUCCESS);
466: }

468: /*@
469:    BVAppendOptionsPrefix - Appends to the prefix used for searching for all
470:    BV options in the database.

472:    Logically Collective

474:    Input Parameters:
475: +  bv     - the basis vectors context
476: -  prefix - the prefix string to prepend to all BV option requests

478:    Notes:
479:    A hyphen (-) must NOT be given at the beginning of the prefix name.
480:    The first character of all runtime options is AUTOMATICALLY the
481:    hyphen.

483:    Level: advanced

485: .seealso: BVSetOptionsPrefix(), BVGetOptionsPrefix()
486: @*/
487: PetscErrorCode BVAppendOptionsPrefix(BV bv,const char *prefix)
488: {
489:   PetscFunctionBegin;
491:   PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix));
492:   PetscFunctionReturn(PETSC_SUCCESS);
493: }

495: /*@
496:    BVGetOptionsPrefix - Gets the prefix used for searching for all
497:    BV options in the database.

499:    Not Collective

501:    Input Parameters:
502: .  bv - the basis vectors context

504:    Output Parameters:
505: .  prefix - pointer to the prefix string used, is returned

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

511:    Level: advanced

513: .seealso: BVSetOptionsPrefix(), BVAppendOptionsPrefix()
514: @*/
515: PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])
516: {
517:   PetscFunctionBegin;
519:   PetscAssertPointer(prefix,2);
520:   PetscCall(PetscObjectGetOptionsPrefix((PetscObject)bv,prefix));
521:   PetscFunctionReturn(PETSC_SUCCESS);
522: }

524: /*@
525:    BVView - Prints the BV data structure.

527:    Collective

529:    Input Parameters:
530: +  bv     - the BV context
531: -  viewer - optional visualization context

533:    Note:
534:    The available visualization contexts include
535: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
536: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
537:          output where only the first processor opens
538:          the file.  All other processors send their
539:          data to the first processor to print.

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

544:    Level: beginner

546: .seealso: BVCreate()
547: @*/
548: PetscErrorCode BVView(BV bv,PetscViewer viewer)
549: {
550:   PetscBool         isascii;
551:   PetscViewerFormat format;
552:   const char        *orthname[2] = {"classical","modified"};
553:   const char        *refname[3] = {"if needed","never","always"};

555:   PetscFunctionBegin;
557:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer));

560:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
561:   if (isascii) {
562:     PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)bv,viewer));
563:     PetscCall(PetscViewerGetFormat(viewer,&format));
564:     if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
565:       PetscCall(PetscViewerASCIIPrintf(viewer,"  %" PetscInt_FMT " columns of global length %" PetscInt_FMT "%s\n",bv->m,bv->N,bv->cuda?" (CUDA)":bv->hip?" (HIP)":""));
566:       if (bv->nc>0) PetscCall(PetscViewerASCIIPrintf(viewer,"  number of constraints: %" PetscInt_FMT "\n",bv->nc));
567:       PetscCall(PetscViewerASCIIPrintf(viewer,"  vector orthogonalization method: %s Gram-Schmidt\n",orthname[bv->orthog_type]));
568:       switch (bv->orthog_ref) {
569:         case BV_ORTHOG_REFINE_IFNEEDED:
570:           PetscCall(PetscViewerASCIIPrintf(viewer,"  orthogonalization refinement: %s (eta: %g)\n",refname[bv->orthog_ref],(double)bv->orthog_eta));
571:           break;
572:         case BV_ORTHOG_REFINE_NEVER:
573:         case BV_ORTHOG_REFINE_ALWAYS:
574:           PetscCall(PetscViewerASCIIPrintf(viewer,"  orthogonalization refinement: %s\n",refname[bv->orthog_ref]));
575:           break;
576:       }
577:       PetscCall(PetscViewerASCIIPrintf(viewer,"  block orthogonalization method: %s\n",BVOrthogBlockTypes[bv->orthog_block]));
578:       if (bv->matrix) {
579:         if (bv->indef) PetscCall(PetscViewerASCIIPrintf(viewer,"  indefinite inner product\n"));
580:         else PetscCall(PetscViewerASCIIPrintf(viewer,"  non-standard inner product\n"));
581:         PetscCall(PetscViewerASCIIPrintf(viewer,"  tolerance for definite inner product: %g\n",(double)bv->deftol));
582:         PetscCall(PetscViewerASCIIPrintf(viewer,"  inner product matrix:\n"));
583:         PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
584:         PetscCall(PetscViewerASCIIPushTab(viewer));
585:         PetscCall(MatView(bv->matrix,viewer));
586:         PetscCall(PetscViewerASCIIPopTab(viewer));
587:         PetscCall(PetscViewerPopFormat(viewer));
588:       }
589:       switch (bv->vmm) {
590:         case BV_MATMULT_VECS:
591:           PetscCall(PetscViewerASCIIPrintf(viewer,"  doing matmult as matrix-vector products\n"));
592:           break;
593:         case BV_MATMULT_MAT:
594:           PetscCall(PetscViewerASCIIPrintf(viewer,"  doing matmult as a single matrix-matrix product\n"));
595:           break;
596:         case BV_MATMULT_MAT_SAVE:
597:           PetscCall(PetscViewerASCIIPrintf(viewer,"  mat_save is deprecated, use mat\n"));
598:           break;
599:       }
600:       if (bv->rrandom) PetscCall(PetscViewerASCIIPrintf(viewer,"  generating random vectors independent of the number of processes\n"));
601:     }
602:   }
603:   PetscTryTypeMethod(bv,view,viewer);
604:   PetscFunctionReturn(PETSC_SUCCESS);
605: }

607: /*@
608:    BVViewFromOptions - View from options

610:    Collective

612:    Input Parameters:
613: +  bv   - the basis vectors context
614: .  obj  - optional object
615: -  name - command line option

617:    Level: intermediate

619: .seealso: BVView(), BVCreate()
620: @*/
621: PetscErrorCode BVViewFromOptions(BV bv,PetscObject obj,const char name[])
622: {
623:   PetscFunctionBegin;
625:   PetscCall(PetscObjectViewFromOptions((PetscObject)bv,obj,name));
626:   PetscFunctionReturn(PETSC_SUCCESS);
627: }

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

632:    Not Collective

634:    Input Parameters:
635: +  name     - name of a new user-defined BV
636: -  function - routine to create context

638:    Notes:
639:    BVRegister() may be called multiple times to add several user-defined
640:    basis vectors.

642:    Level: advanced

644: .seealso: BVRegisterAll()
645: @*/
646: PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))
647: {
648:   PetscFunctionBegin;
649:   PetscCall(BVInitializePackage());
650:   PetscCall(PetscFunctionListAdd(&BVList,name,function));
651:   PetscFunctionReturn(PETSC_SUCCESS);
652: }

654: PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)
655: {
656:   PetscFunctionBegin;
657:   if (s>bv->lwork) {
658:     PetscCall(PetscFree(bv->work));
659:     PetscCall(PetscMalloc1(s,&bv->work));
660:     bv->lwork = s;
661:   }
662:   PetscFunctionReturn(PETSC_SUCCESS);
663: }

665: #if defined(PETSC_USE_DEBUG) && !defined(PETSC_CLANG_STATIC_ANALYZER)
666: /*
667:    SlepcDebugBVView - partially view a BV object, to be used from within a debugger.

669:      ini, end: columns to be viewed
670:      s: name of Matlab variable
671:      filename: optionally write output to a file
672:  */
673: PETSC_UNUSED PetscErrorCode SlepcDebugBVView(BV bv,PetscInt ini,PetscInt end,const char *s,const char *filename)
674: {
675:   PetscInt       N,m;
676:   PetscScalar    *array;

678:   PetscFunctionBegin;
679:   PetscCall(BVGetArray(bv,&array));
680:   PetscCall(BVGetSizes(bv,NULL,&N,&m));
681:   PetscCall(SlepcDebugViewMatrix(N,end-ini+1,array+ini*N,NULL,bv->ld,s,filename));
682:   PetscCall(BVRestoreArray(bv,&array));
683:   PetscFunctionReturn(PETSC_SUCCESS);
684: }
685: #endif