Actual source code: bvfunc.c

  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: [](sec:bv), `SlepcFinalize()`, `BVInitializePackage()`
 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_slepc()` when using dynamic libraries, and
 49:    on the first call to `BVCreate()` when using shared or static libraries.

 51:    Note:
 52:    This function never needs to be called by SLEPc users.

 54:    Level: developer

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

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

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

110:    Collective

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

115:    Level: beginner

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

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

149:    Collective

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

154:    Output Parameter:
155: .  newbv - location to put the basis vectors context

157:    Note:
158:    A basis vectors object of type `BVTENSOR` should be created using
159:    the helper function `BVCreateTensor()` instead of this one.

161:    Level: beginner

163: .seealso: [](sec:bv), `BV`, `BVDuplicate()`, `BVDestroy()`, `BVCreateTensor()` , `BVCreateFromMat()`
164: @*/
165: PetscErrorCode BVCreate(MPI_Comm comm,BV *newbv)
166: {
167:   BV             bv;

169:   PetscFunctionBegin;
170:   PetscAssertPointer(newbv,2);
171:   PetscCall(BVInitializePackage());
172:   PetscCall(SlepcHeaderCreate(bv,BV_CLASSID,"BV","Basis Vectors","BV",comm,BVDestroy,BVView));

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

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

229:   *newbv = bv;
230:   PetscFunctionReturn(PETSC_SUCCESS);
231: }

233: /*@
234:    BVCreateFromMat - Creates a basis vectors object from a dense `Mat` object.

236:    Collective

238:    Input Parameter:
239: .  A - a dense tall-skinny matrix

241:    Output Parameter:
242: .  bv - the new basis vectors context

244:    Notes:
245:    The matrix values are copied to the `BV` data storage, memory is not shared.

247:    The communicator of the `BV` object will be the same as `A`, and so will be
248:    the dimensions.

250:    Level: intermediate

252: .seealso: [](sec:bv), `BVCreate()`, `BVDestroy()`, `BVCreateMat()`
253: @*/
254: PetscErrorCode BVCreateFromMat(Mat A,BV *bv)
255: {
256:   PetscInt  n,N,k;
257:   VecType   vtype;

259:   PetscFunctionBegin;

262:   PetscCall(MatGetSize(A,&N,&k));
263:   PetscCall(MatGetLocalSize(A,&n,NULL));
264:   PetscCall(MatGetVecType(A,&vtype));
265:   PetscCall(BVCreate(PetscObjectComm((PetscObject)A),bv));
266:   PetscCall(BVSetSizes(*bv,n,N,k));
267:   PetscCall(BVSetVecType(*bv,vtype));

269:   (*bv)->Acreate = A;
270:   PetscCall(PetscObjectReference((PetscObject)A));
271:   PetscFunctionReturn(PETSC_SUCCESS);
272: }

274: /*@
275:    BVInsertVec - Insert a vector into the specified column.

277:    Logically Collective

279:    Input Parameters:
280: +  V - basis vectors
281: .  j - the column of `V` to be overwritten
282: -  w - the vector to be copied

284:    Level: intermediate

286: .seealso: [](sec:bv), `BVInsertVecs()`
287: @*/
288: PetscErrorCode BVInsertVec(BV V,PetscInt j,Vec w)
289: {
290:   PetscInt       n,N;
291:   Vec            v;

293:   PetscFunctionBegin;
298:   BVCheckSizes(V,1);
299:   PetscCheckSameComm(V,1,w,3);

301:   PetscCall(VecGetSize(w,&N));
302:   PetscCall(VecGetLocalSize(w,&n));
303:   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);
304:   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);

306:   PetscCall(BVGetColumn(V,j,&v));
307:   PetscCall(VecCopy(w,v));
308:   PetscCall(BVRestoreColumn(V,j,&v));
309:   PetscCall(PetscObjectStateIncrease((PetscObject)V));
310:   PetscFunctionReturn(PETSC_SUCCESS);
311: }

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

316:    Collective

318:    Input Parameters:
319: +  V - basis vectors
320: .  s - first column of `V` to be overwritten
321: .  W - array of vectors to be copied
322: -  orth - flag indicating if the vectors must be orthogonalized

324:    Input/Output Parameter:
325: .  m - number of input vectors, on output the number of linearly independent
326:        vectors

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

334:    Level: intermediate

336: .seealso: [](sec:bv), `BVInsertVec()`, `BVOrthogonalizeColumn()`
337: @*/
338: PetscErrorCode BVInsertVecs(BV V,PetscInt s,PetscInt *m,Vec W[],PetscBool orth)
339: {
340:   PetscInt       n,N,i,ndep;
341:   PetscBool      lindep;
342:   PetscReal      norm;
343:   Vec            v;

345:   PetscFunctionBegin;
348:   PetscAssertPointer(m,3);
350:   if (!*m) PetscFunctionReturn(PETSC_SUCCESS);
351:   PetscCheck(*m>0,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %" PetscInt_FMT ") cannot be negative",*m);
352:   PetscAssertPointer(W,4);
356:   BVCheckSizes(V,1);
357:   PetscCheckSameComm(V,1,*W,4);

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

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

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

386:    Collective

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 behavior is similar to `BVInsertVecs()`.

414:    Level: advanced

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

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

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

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

452:    Logically Collective

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

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

463:    Level: advanced

465: .seealso: [](sec:bv), `BVAppendOptionsPrefix()`, `BVGetOptionsPrefix()`
466: @*/
467: PetscErrorCode BVSetOptionsPrefix(BV bv,const char prefix[])
468: {
469:   PetscFunctionBegin;
471:   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)bv,prefix));
472:   PetscFunctionReturn(PETSC_SUCCESS);
473: }

475: /*@
476:    BVAppendOptionsPrefix - Appends to the prefix used for searching for all
477:    `BV` options in the database.

479:    Logically Collective

481:    Input Parameters:
482: +  bv     - the basis vectors context
483: -  prefix - the prefix string to prepend to all `BV` option requests

485:    Note:
486:    A hyphen (-) must NOT be given at the beginning of the prefix name.
487:    The first character of all runtime options is AUTOMATICALLY the
488:    hyphen.

490:    Level: advanced

492: .seealso: [](sec:bv), `BVSetOptionsPrefix()`, `BVGetOptionsPrefix()`
493: @*/
494: PetscErrorCode BVAppendOptionsPrefix(BV bv,const char prefix[])
495: {
496:   PetscFunctionBegin;
498:   PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix));
499:   PetscFunctionReturn(PETSC_SUCCESS);
500: }

502: /*@
503:    BVGetOptionsPrefix - Gets the prefix used for searching for all
504:    `BV` options in the database.

506:    Not Collective

508:    Input Parameter:
509: .  bv - the basis vectors context

511:    Output Parameter:
512: .  prefix - pointer to the prefix string used is returned

514:    Level: advanced

516: .seealso: [](sec:bv), `BVSetOptionsPrefix()`, `BVAppendOptionsPrefix()`
517: @*/
518: PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])
519: {
520:   PetscFunctionBegin;
522:   PetscAssertPointer(prefix,2);
523:   PetscCall(PetscObjectGetOptionsPrefix((PetscObject)bv,prefix));
524:   PetscFunctionReturn(PETSC_SUCCESS);
525: }

527: /*@
528:    BVView - Prints the `BV` data structure.

530:    Collective

532:    Input Parameters:
533: +  bv     - the basis vectors context
534: -  viewer - optional visualization context

536:    Notes:
537:    The available visualization contexts include
538: +     `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
539: -     `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
540:          first process opens the file; all other processes send their data to the
541:          first one to print

543:    The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
544:    to output to a specified file.

546:    Use `BVViewFromOptions()` to allow the user to select many different `PetscViewerType`
547:    and formats from the options database.

549:    Level: beginner

551: .seealso: [](sec:bv), `BVCreate()`
552: @*/
553: PetscErrorCode BVView(BV bv,PetscViewer viewer)
554: {
555:   PetscBool         isascii;
556:   PetscViewerFormat format;
557:   const char        *orthname[2] = {"classical","modified"};
558:   const char        *refname[3] = {"if needed","never","always"};

560:   PetscFunctionBegin;
562:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer));

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

612: /*@
613:    BVViewFromOptions - View (print) a `BV` object based on values in the options database.

615:    Collective

617:    Input Parameters:
618: +  bv   - the basis vectors context
619: .  obj  - optional object that provides the options prefix used to query the options database
620: -  name - command line option

622:    Level: intermediate

624: .seealso: [](sec:bv), `BVView()`, `BVCreate()`, `PetscObjectViewFromOptions()`
625: @*/
626: PetscErrorCode BVViewFromOptions(BV bv,PetscObject obj,const char name[])
627: {
628:   PetscFunctionBegin;
630:   PetscCall(PetscObjectViewFromOptions((PetscObject)bv,obj,name));
631:   PetscFunctionReturn(PETSC_SUCCESS);
632: }

634: /*@C
635:    BVRegister - Adds a new storage format to the `BV` package.

637:    Not Collective

639:    Input Parameters:
640: +  name     - name of a new user-defined `BV`
641: -  function - routine to create context

643:    Notes:
644:    `BVRegister()` may be called multiple times to add several user-defined
645:    basis vectors.

647:    Level: advanced

649: .seealso: [](sec:bv), `BVRegisterAll()`
650: @*/
651: PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))
652: {
653:   PetscFunctionBegin;
654:   PetscCall(BVInitializePackage());
655:   PetscCall(PetscFunctionListAdd(&BVList,name,function));
656:   PetscFunctionReturn(PETSC_SUCCESS);
657: }

659: PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)
660: {
661:   PetscFunctionBegin;
662:   if (s>bv->lwork) {
663:     PetscCall(PetscFree(bv->work));
664:     PetscCall(PetscMalloc1(s,&bv->work));
665:     bv->lwork = s;
666:   }
667:   PetscFunctionReturn(PETSC_SUCCESS);
668: }

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

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

683:   PetscFunctionBegin;
684:   PetscCall(BVGetArray(bv,&array));
685:   PetscCall(BVGetSizes(bv,NULL,&N,&m));
686:   PetscCall(SlepcDebugViewMatrix(N,end-ini+1,array+ini*N,NULL,bv->ld,s,filename));
687:   PetscCall(BVRestoreArray(bv,&array));
688:   PetscFunctionReturn(PETSC_SUCCESS);
689: }
690: #endif