Actual source code: bvfunc.c
slepc-3.19.0 2023-03-31
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: /*@C
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(VecDestroy(&(*bv)->t));
125: PetscCall(MatDestroy(&(*bv)->matrix));
126: PetscCall(VecDestroy(&(*bv)->Bx));
127: PetscCall(VecDestroy(&(*bv)->buffer));
128: PetscCall(BVDestroy(&(*bv)->cached));
129: PetscCall(BVDestroy(&(*bv)->L));
130: PetscCall(BVDestroy(&(*bv)->R));
131: PetscCall(PetscFree((*bv)->work));
132: PetscCall(PetscFree2((*bv)->h,(*bv)->c));
133: PetscCall(VecDestroy(&(*bv)->omega));
134: PetscCall(MatDestroy(&(*bv)->Acreate));
135: PetscCall(MatDestroy(&(*bv)->Aget));
136: PetscCall(MatDestroy(&(*bv)->Abuffer));
137: PetscCall(PetscRandomDestroy(&(*bv)->rand));
138: PetscCall(PetscHeaderDestroy(bv));
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: /*@
143: BVCreate - Creates a basis vectors context.
145: Collective
147: Input Parameter:
148: . comm - MPI communicator
150: Output Parameter:
151: . newbv - location to put the basis vectors context
153: Level: beginner
155: .seealso: BVSetUp(), BVDestroy(), BV
156: @*/
157: PetscErrorCode BVCreate(MPI_Comm comm,BV *newbv)
158: {
159: BV bv;
161: PetscFunctionBegin;
163: *newbv = NULL;
164: PetscCall(BVInitializePackage());
165: PetscCall(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;
181: bv->rrandom = PETSC_FALSE;
182: bv->deftol = 10*PETSC_MACHINE_EPSILON;
184: bv->buffer = NULL;
185: bv->Abuffer = NULL;
186: bv->Bx = NULL;
187: bv->xid = 0;
188: bv->xstate = 0;
189: bv->cv[0] = NULL;
190: bv->cv[1] = NULL;
191: bv->ci[0] = -1;
192: bv->ci[1] = -1;
193: bv->st[0] = -1;
194: bv->st[1] = -1;
195: bv->id[0] = 0;
196: bv->id[1] = 0;
197: bv->h = NULL;
198: bv->c = NULL;
199: bv->omega = NULL;
200: bv->defersfo = PETSC_FALSE;
201: bv->cached = NULL;
202: bv->bvstate = 0;
203: bv->L = NULL;
204: bv->R = NULL;
205: bv->lstate = 0;
206: bv->rstate = 0;
207: bv->lsplit = 0;
208: bv->issplit = 0;
209: bv->splitparent = NULL;
210: bv->rand = NULL;
211: bv->Acreate = NULL;
212: bv->Aget = NULL;
213: bv->cuda = PETSC_FALSE;
214: bv->sfocalled = PETSC_FALSE;
215: bv->work = NULL;
216: bv->lwork = 0;
217: bv->data = NULL;
219: *newbv = bv;
220: PetscFunctionReturn(PETSC_SUCCESS);
221: }
223: /*@
224: BVCreateFromMat - Creates a basis vectors object from a dense Mat object.
226: Collective
228: Input Parameter:
229: . A - a dense tall-skinny matrix
231: Output Parameter:
232: . bv - the new basis vectors context
234: Notes:
235: The matrix values are copied to the BV data storage, memory is not shared.
237: The communicator of the BV object will be the same as A, and so will be
238: the dimensions.
240: Level: intermediate
242: .seealso: BVCreate(), BVDestroy(), BVCreateMat()
243: @*/
244: PetscErrorCode BVCreateFromMat(Mat A,BV *bv)
245: {
246: PetscInt n,N,k;
248: PetscFunctionBegin;
250: PetscCheckTypeNames(A,MATSEQDENSE,MATMPIDENSE);
252: PetscCall(MatGetSize(A,&N,&k));
253: PetscCall(MatGetLocalSize(A,&n,NULL));
254: PetscCall(BVCreate(PetscObjectComm((PetscObject)A),bv));
255: PetscCall(BVSetSizes(*bv,n,N,k));
257: (*bv)->Acreate = A;
258: PetscCall(PetscObjectReference((PetscObject)A));
259: PetscFunctionReturn(PETSC_SUCCESS);
260: }
262: /*@
263: BVInsertVec - Insert a vector into the specified column.
265: Logically Collective
267: Input Parameters:
268: + V - basis vectors
269: . j - the column of V to be overwritten
270: - w - the vector to be copied
272: Level: intermediate
274: .seealso: BVInsertVecs()
275: @*/
276: PetscErrorCode BVInsertVec(BV V,PetscInt j,Vec w)
277: {
278: PetscInt n,N;
279: Vec v;
281: PetscFunctionBegin;
286: BVCheckSizes(V,1);
287: PetscCheckSameComm(V,1,w,3);
289: PetscCall(VecGetSize(w,&N));
290: PetscCall(VecGetLocalSize(w,&n));
291: 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);
292: 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);
294: PetscCall(BVGetColumn(V,j,&v));
295: PetscCall(VecCopy(w,v));
296: PetscCall(BVRestoreColumn(V,j,&v));
297: PetscCall(PetscObjectStateIncrease((PetscObject)V));
298: PetscFunctionReturn(PETSC_SUCCESS);
299: }
301: /*@
302: BVInsertVecs - Insert a set of vectors into the specified columns.
304: Collective
306: Input Parameters:
307: + V - basis vectors
308: . s - first column of V to be overwritten
309: . W - set of vectors to be copied
310: - orth - flag indicating if the vectors must be orthogonalized
312: Input/Output Parameter:
313: . m - number of input vectors, on output the number of linearly independent
314: vectors
316: Notes:
317: Copies the contents of vectors W to V(:,s:s+n). If the orthogonalization
318: flag is set, then the vectors are copied one by one and then orthogonalized
319: against the previous ones. If any of them is linearly dependent then it
320: is discarded and the value of m is decreased.
322: Level: intermediate
324: .seealso: BVInsertVec(), BVOrthogonalizeColumn()
325: @*/
326: PetscErrorCode BVInsertVecs(BV V,PetscInt s,PetscInt *m,Vec *W,PetscBool orth)
327: {
328: PetscInt n,N,i,ndep;
329: PetscBool lindep;
330: PetscReal norm;
331: Vec v;
333: PetscFunctionBegin;
338: if (!*m) PetscFunctionReturn(PETSC_SUCCESS);
339: PetscCheck(*m>0,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %" PetscInt_FMT ") cannot be negative",*m);
344: BVCheckSizes(V,1);
345: PetscCheckSameComm(V,1,*W,4);
347: PetscCall(VecGetSize(*W,&N));
348: PetscCall(VecGetLocalSize(*W,&n));
349: 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);
350: 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);
351: 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);
353: ndep = 0;
354: for (i=0;i<*m;i++) {
355: PetscCall(BVGetColumn(V,s+i-ndep,&v));
356: PetscCall(VecCopy(W[i],v));
357: PetscCall(BVRestoreColumn(V,s+i-ndep,&v));
358: if (orth) {
359: PetscCall(BVOrthogonalizeColumn(V,s+i-ndep,NULL,&norm,&lindep));
360: if (norm==0.0 || lindep) {
361: PetscCall(PetscInfo(V,"Removing linearly dependent vector %" PetscInt_FMT "\n",i));
362: ndep++;
363: } else PetscCall(BVScaleColumn(V,s+i-ndep,1.0/norm));
364: }
365: }
366: *m -= ndep;
367: PetscCall(PetscObjectStateIncrease((PetscObject)V));
368: PetscFunctionReturn(PETSC_SUCCESS);
369: }
371: /*@
372: BVInsertConstraints - Insert a set of vectors as constraints.
374: Collective
376: Input Parameters:
377: + V - basis vectors
378: - C - set of vectors to be inserted as constraints
380: Input/Output Parameter:
381: . nc - number of input vectors, on output the number of linearly independent
382: vectors
384: Notes:
385: The constraints are relevant only during orthogonalization. Constraint
386: vectors span a subspace that is deflated in every orthogonalization
387: operation, so they are intended for removing those directions from the
388: orthogonal basis computed in regular BV columns.
390: Constraints are not stored in regular BV columns, but in a special part of
391: the storage. They can be accessed with negative indices in BVGetColumn().
393: This operation is DESTRUCTIVE, meaning that all data contained in the
394: columns of V is lost. This is typically invoked just after creating the BV.
395: Once a set of constraints has been set, it is not allowed to call this
396: function again.
398: The vectors are copied one by one and then orthogonalized against the
399: previous ones. If any of them is linearly dependent then it is discarded
400: and the value of nc is decreased. The behaviour is similar to BVInsertVecs().
402: Level: advanced
404: .seealso: BVInsertVecs(), BVOrthogonalizeColumn(), BVGetColumn(), BVGetNumConstraints()
405: @*/
406: PetscErrorCode BVInsertConstraints(BV V,PetscInt *nc,Vec *C)
407: {
408: PetscInt msave;
410: PetscFunctionBegin;
414: if (!*nc) PetscFunctionReturn(PETSC_SUCCESS);
415: PetscCheck(*nc>0,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of constraints (given %" PetscInt_FMT ") cannot be negative",*nc);
419: BVCheckSizes(V,1);
420: PetscCheckSameComm(V,1,*C,3);
421: PetscCheck(!V->issplit,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Operation not permitted for a BV obtained from BVGetSplit");
422: PetscCheck(!V->nc,PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Constraints already present in this BV object");
423: PetscCheck(V->ci[0]==-1 && V->ci[1]==-1,PetscObjectComm((PetscObject)V),PETSC_ERR_SUP,"Cannot call BVInsertConstraints after BVGetColumn");
425: msave = V->m;
426: PetscCall(BVResize(V,*nc+V->m,PETSC_FALSE));
427: PetscCall(BVInsertVecs(V,0,nc,C,PETSC_TRUE));
428: V->nc = *nc;
429: V->m = msave;
430: V->ci[0] = -V->nc-1;
431: V->ci[1] = -V->nc-1;
432: PetscCall(PetscObjectStateIncrease((PetscObject)V));
433: PetscFunctionReturn(PETSC_SUCCESS);
434: }
436: /*@C
437: BVSetOptionsPrefix - Sets the prefix used for searching for all
438: BV options in the database.
440: Logically Collective
442: Input Parameters:
443: + bv - the basis vectors context
444: - prefix - the prefix string to prepend to all BV option requests
446: Notes:
447: A hyphen (-) must NOT be given at the beginning of the prefix name.
448: The first character of all runtime options is AUTOMATICALLY the
449: hyphen.
451: Level: advanced
453: .seealso: BVAppendOptionsPrefix(), BVGetOptionsPrefix()
454: @*/
455: PetscErrorCode BVSetOptionsPrefix(BV bv,const char *prefix)
456: {
457: PetscFunctionBegin;
459: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)bv,prefix));
460: PetscFunctionReturn(PETSC_SUCCESS);
461: }
463: /*@C
464: BVAppendOptionsPrefix - Appends to the prefix used for searching for all
465: BV options in the database.
467: Logically Collective
469: Input Parameters:
470: + bv - the basis vectors context
471: - prefix - the prefix string to prepend to all BV option requests
473: Notes:
474: A hyphen (-) must NOT be given at the beginning of the prefix name.
475: The first character of all runtime options is AUTOMATICALLY the
476: hyphen.
478: Level: advanced
480: .seealso: BVSetOptionsPrefix(), BVGetOptionsPrefix()
481: @*/
482: PetscErrorCode BVAppendOptionsPrefix(BV bv,const char *prefix)
483: {
484: PetscFunctionBegin;
486: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix));
487: PetscFunctionReturn(PETSC_SUCCESS);
488: }
490: /*@C
491: BVGetOptionsPrefix - Gets the prefix used for searching for all
492: BV options in the database.
494: Not Collective
496: Input Parameters:
497: . bv - the basis vectors context
499: Output Parameters:
500: . prefix - pointer to the prefix string used, is returned
502: Note:
503: On the Fortran side, the user should pass in a string 'prefix' of
504: sufficient length to hold the prefix.
506: Level: advanced
508: .seealso: BVSetOptionsPrefix(), BVAppendOptionsPrefix()
509: @*/
510: PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])
511: {
512: PetscFunctionBegin;
515: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)bv,prefix));
516: PetscFunctionReturn(PETSC_SUCCESS);
517: }
519: /*@C
520: BVView - Prints the BV data structure.
522: Collective
524: Input Parameters:
525: + bv - the BV context
526: - viewer - optional visualization context
528: Note:
529: The available visualization contexts include
530: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
531: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
532: output where only the first processor opens
533: the file. All other processors send their
534: data to the first processor to print.
536: The user can open an alternative visualization contexts with
537: PetscViewerASCIIOpen() (output to a specified file).
539: Level: beginner
541: .seealso: BVCreate()
542: @*/
543: PetscErrorCode BVView(BV bv,PetscViewer viewer)
544: {
545: PetscBool isascii;
546: PetscViewerFormat format;
547: const char *orthname[2] = {"classical","modified"};
548: const char *refname[3] = {"if needed","never","always"};
550: PetscFunctionBegin;
552: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer));
555: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
556: if (isascii) {
557: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)bv,viewer));
558: PetscCall(PetscViewerGetFormat(viewer,&format));
559: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
560: PetscCall(PetscViewerASCIIPrintf(viewer," %" PetscInt_FMT " columns of global length %" PetscInt_FMT "%s\n",bv->m,bv->N,bv->cuda?" (CUDA)":""));
561: if (bv->nc>0) PetscCall(PetscViewerASCIIPrintf(viewer," number of constraints: %" PetscInt_FMT "\n",bv->nc));
562: PetscCall(PetscViewerASCIIPrintf(viewer," vector orthogonalization method: %s Gram-Schmidt\n",orthname[bv->orthog_type]));
563: switch (bv->orthog_ref) {
564: case BV_ORTHOG_REFINE_IFNEEDED:
565: PetscCall(PetscViewerASCIIPrintf(viewer," orthogonalization refinement: %s (eta: %g)\n",refname[bv->orthog_ref],(double)bv->orthog_eta));
566: break;
567: case BV_ORTHOG_REFINE_NEVER:
568: case BV_ORTHOG_REFINE_ALWAYS:
569: PetscCall(PetscViewerASCIIPrintf(viewer," orthogonalization refinement: %s\n",refname[bv->orthog_ref]));
570: break;
571: }
572: PetscCall(PetscViewerASCIIPrintf(viewer," block orthogonalization method: %s\n",BVOrthogBlockTypes[bv->orthog_block]));
573: if (bv->matrix) {
574: if (bv->indef) PetscCall(PetscViewerASCIIPrintf(viewer," indefinite inner product\n"));
575: else PetscCall(PetscViewerASCIIPrintf(viewer," non-standard inner product\n"));
576: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance for definite inner product: %g\n",(double)bv->deftol));
577: PetscCall(PetscViewerASCIIPrintf(viewer," inner product matrix:\n"));
578: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
579: PetscCall(PetscViewerASCIIPushTab(viewer));
580: PetscCall(MatView(bv->matrix,viewer));
581: PetscCall(PetscViewerASCIIPopTab(viewer));
582: PetscCall(PetscViewerPopFormat(viewer));
583: }
584: switch (bv->vmm) {
585: case BV_MATMULT_VECS:
586: PetscCall(PetscViewerASCIIPrintf(viewer," doing matmult as matrix-vector products\n"));
587: break;
588: case BV_MATMULT_MAT:
589: PetscCall(PetscViewerASCIIPrintf(viewer," doing matmult as a single matrix-matrix product\n"));
590: break;
591: case BV_MATMULT_MAT_SAVE:
592: PetscCall(PetscViewerASCIIPrintf(viewer," mat_save is deprecated, use mat\n"));
593: break;
594: }
595: if (bv->rrandom) PetscCall(PetscViewerASCIIPrintf(viewer," generating random vectors independent of the number of processes\n"));
596: }
597: }
598: PetscTryTypeMethod(bv,view,viewer);
599: PetscFunctionReturn(PETSC_SUCCESS);
600: }
602: /*@C
603: BVViewFromOptions - View from options
605: Collective
607: Input Parameters:
608: + bv - the basis vectors context
609: . obj - optional object
610: - name - command line option
612: Level: intermediate
614: .seealso: BVView(), BVCreate()
615: @*/
616: PetscErrorCode BVViewFromOptions(BV bv,PetscObject obj,const char name[])
617: {
618: PetscFunctionBegin;
620: PetscCall(PetscObjectViewFromOptions((PetscObject)bv,obj,name));
621: PetscFunctionReturn(PETSC_SUCCESS);
622: }
624: /*@C
625: BVRegister - Adds a new storage format to the BV package.
627: Not Collective
629: Input Parameters:
630: + name - name of a new user-defined BV
631: - function - routine to create context
633: Notes:
634: BVRegister() may be called multiple times to add several user-defined
635: basis vectors.
637: Level: advanced
639: .seealso: BVRegisterAll()
640: @*/
641: PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))
642: {
643: PetscFunctionBegin;
644: PetscCall(BVInitializePackage());
645: PetscCall(PetscFunctionListAdd(&BVList,name,function));
646: PetscFunctionReturn(PETSC_SUCCESS);
647: }
649: PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)
650: {
651: PetscFunctionBegin;
652: if (s>bv->lwork) {
653: PetscCall(PetscFree(bv->work));
654: PetscCall(PetscMalloc1(s,&bv->work));
655: bv->lwork = s;
656: }
657: PetscFunctionReturn(PETSC_SUCCESS);
658: }
660: #if defined(PETSC_USE_DEBUG)
661: /*
662: SlepcDebugBVView - partially view a BV object, to be used from within a debugger.
664: ini, end: columns to be viewed
665: s: name of Matlab variable
666: filename: optionally write output to a file
667: */
668: PETSC_UNUSED PetscErrorCode SlepcDebugBVView(BV bv,PetscInt ini,PetscInt end,const char *s,const char *filename)
669: {
670: PetscInt N,m;
671: PetscScalar *array;
673: PetscFunctionBegin;
674: PetscCall(BVGetArray(bv,&array));
675: PetscCall(BVGetSizes(bv,NULL,&N,&m));
676: PetscCall(SlepcDebugViewMatrix(N,end-ini+1,array+ini*N,NULL,N,s,filename));
677: PetscCall(BVRestoreArray(bv,&array));
678: PetscFunctionReturn(PETSC_SUCCESS);
679: }
680: #endif