1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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;
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};
26: /*@C
27: BVFinalizePackage - This function destroys everything in the Slepc interface
28: to the BV package. It is called from SlepcFinalize().
30: Level: developer
32: .seealso: SlepcFinalize()
33: @*/
34: PetscErrorCode BVFinalizePackage(void) 35: {
39: PetscFunctionListDestroy(&BVList);
40: MPI_Op_free(&MPIU_TSQR);
41: MPI_Op_free(&MPIU_LAPY2);
42: BVPackageInitialized = PETSC_FALSE;
43: BVRegisterAllCalled = PETSC_FALSE;
44: return(0);
45: }
47: /*@C
48: BVInitializePackage - This function initializes everything in the BV package.
49: It is called from PetscDLLibraryRegister() when using dynamic libraries, and
50: on the first call to BVCreate() when using static libraries.
52: Level: developer
54: .seealso: SlepcInitialize()
55: @*/
56: PetscErrorCode BVInitializePackage(void) 57: {
58: char logList[256];
59: PetscBool opt,pkg;
60: PetscClassId classids[1];
64: if (BVPackageInitialized) return(0);
65: BVPackageInitialized = PETSC_TRUE;
66: /* Register Classes */
67: PetscClassIdRegister("Basis Vectors",&BV_CLASSID);
68: /* Register Constructors */
69: BVRegisterAll();
70: /* Register Events */
71: PetscLogEventRegister("BVCreate",BV_CLASSID,&BV_Create);
72: PetscLogEventRegister("BVCopy",BV_CLASSID,&BV_Copy);
73: PetscLogEventRegister("BVMult",BV_CLASSID,&BV_Mult);
74: PetscLogEventRegister("BVMultVec",BV_CLASSID,&BV_MultVec);
75: PetscLogEventRegister("BVMultInPlace",BV_CLASSID,&BV_MultInPlace);
76: PetscLogEventRegister("BVDot",BV_CLASSID,&BV_Dot);
77: PetscLogEventRegister("BVDotVec",BV_CLASSID,&BV_DotVec);
78: PetscLogEventRegister("BVOrthogonalize",BV_CLASSID,&BV_Orthogonalize);
79: PetscLogEventRegister("BVOrthogonalizeV",BV_CLASSID,&BV_OrthogonalizeVec);
80: PetscLogEventRegister("BVScale",BV_CLASSID,&BV_Scale);
81: PetscLogEventRegister("BVNorm",BV_CLASSID,&BV_Norm);
82: PetscLogEventRegister("BVNormVec",BV_CLASSID,&BV_NormVec);
83: PetscLogEventRegister("BVNormalize",BV_CLASSID,&BV_Normalize);
84: PetscLogEventRegister("BVSetRandom",BV_CLASSID,&BV_SetRandom);
85: PetscLogEventRegister("BVMatMult",BV_CLASSID,&BV_MatMult);
86: PetscLogEventRegister("BVMatMultVec",BV_CLASSID,&BV_MatMultVec);
87: PetscLogEventRegister("BVMatProject",BV_CLASSID,&BV_MatProject);
88: /* MPI reduction operation used in BVOrthogonalize */
89: MPI_Op_create(SlepcGivensPacked,PETSC_FALSE,&MPIU_TSQR);
90: MPI_Op_create(SlepcPythag,PETSC_TRUE,&MPIU_LAPY2);
91: /* Process Info */
92: classids[0] = BV_CLASSID;
93: PetscInfoProcessClass("bv",1,&classids[0]);
94: /* Process summary exclusions */
95: PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt);
96: if (opt) {
97: PetscStrInList("bv",logList,',',&pkg);
98: if (pkg) { PetscLogEventDeactivateClass(BV_CLASSID); }
99: }
100: /* Register package finalizer */
101: PetscRegisterFinalize(BVFinalizePackage);
102: return(0);
103: }
105: /*@
106: BVDestroy - Destroys BV context that was created with BVCreate().
108: Collective on bv
110: Input Parameter:
111: . bv - the basis vectors context
113: Level: beginner
115: .seealso: BVCreate()
116: @*/
117: PetscErrorCode BVDestroy(BV *bv)118: {
122: if (!*bv) return(0);
124: if ((*bv)->lsplit) SETERRQ(PetscObjectComm((PetscObject)(*bv)),PETSC_ERR_ARG_WRONGSTATE,"Must call BVRestoreSplit before destroying the BV");
125: if (--((PetscObject)(*bv))->refct > 0) { *bv = 0; return(0); }
126: if ((*bv)->ops->destroy) { (*(*bv)->ops->destroy)(*bv); }
127: VecDestroy(&(*bv)->t);
128: MatDestroy(&(*bv)->matrix);
129: VecDestroy(&(*bv)->Bx);
130: VecDestroy(&(*bv)->buffer);
131: BVDestroy(&(*bv)->cached);
132: BVDestroy(&(*bv)->L);
133: BVDestroy(&(*bv)->R);
134: PetscFree((*bv)->work);
135: PetscFree2((*bv)->h,(*bv)->c);
136: VecDestroy(&(*bv)->omega);
137: MatDestroy(&(*bv)->Acreate);
138: MatDestroy(&(*bv)->Aget);
139: MatDestroy(&(*bv)->Abuffer);
140: PetscRandomDestroy(&(*bv)->rand);
141: PetscHeaderDestroy(bv);
142: return(0);
143: }
145: /*@
146: BVCreate - Creates a basis vectors context.
148: Collective
150: Input Parameter:
151: . comm - MPI communicator
153: Output Parameter:
154: . bv - location to put the basis vectors context
156: Level: beginner
158: .seealso: BVSetUp(), BVDestroy(), BV159: @*/
160: PetscErrorCode BVCreate(MPI_Comm comm,BV *newbv)161: {
163: BV bv;
167: *newbv = 0;
168: BVInitializePackage();
169: SlepcHeaderCreate(bv,BV_CLASSID,"BV","Basis Vectors","BV",comm,BVDestroy,BVView);
171: bv->t = NULL;
172: bv->n = -1;
173: bv->N = -1;
174: bv->m = 0;
175: bv->l = 0;
176: bv->k = 0;
177: bv->nc = 0;
178: bv->orthog_type = BV_ORTHOG_CGS;
179: bv->orthog_ref = BV_ORTHOG_REFINE_IFNEEDED;
180: bv->orthog_eta = 0.7071;
181: bv->orthog_block = BV_ORTHOG_BLOCK_GS;
182: bv->matrix = NULL;
183: bv->indef = PETSC_FALSE;
184: bv->vmm = BV_MATMULT_MAT;
186: bv->Bx = NULL;
187: bv->buffer = NULL;
188: bv->Abuffer = 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->rrandom = PETSC_FALSE;
214: bv->Acreate = NULL;
215: bv->Aget = NULL;
216: bv->cuda = PETSC_FALSE;
217: bv->work = NULL;
218: bv->lwork = 0;
219: bv->data = NULL;
221: *newbv = bv;
222: return(0);
223: }
225: /*@
226: BVCreateFromMat - Creates a basis vectors object from a dense Mat object.
228: Collective on A
230: Input Parameter:
231: . A - a dense tall-skinny matrix
233: Output Parameter:
234: . bv - the new basis vectors context
236: Notes:
237: The matrix values are copied to the BV data storage, memory is not shared.
239: The communicator of the BV object will be the same as A, and so will be
240: the dimensions.
242: Level: intermediate
244: .seealso: BVCreate(), BVDestroy(), BVCreateMat()
245: @*/
246: PetscErrorCode BVCreateFromMat(Mat A,BV *bv)247: {
249: PetscBool match;
250: PetscInt n,N,k;
254: PetscObjectTypeCompareAny((PetscObject)A,&match,MATSEQDENSE,MATMPIDENSE,"");
255: if (!match) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_SUP,"Matrix must be of type dense");
257: MatGetSize(A,&N,&k);
258: MatGetLocalSize(A,&n,NULL);
259: BVCreate(PetscObjectComm((PetscObject)A),bv);
260: BVSetSizes(*bv,n,N,k);
262: (*bv)->Acreate = A;
263: PetscObjectReference((PetscObject)A);
264: return(0);
265: }
267: /*@
268: BVInsertVec - Insert a vector into the specified column.
270: Collective on V
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: {
284: PetscInt n,N;
285: Vec v;
292: BVCheckSizes(V,1);
295: VecGetSize(w,&N);
296: VecGetLocalSize(w,&n);
297: 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);
298: 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);
300: BVGetColumn(V,j,&v);
301: VecCopy(w,v);
302: BVRestoreColumn(V,j,&v);
303: PetscObjectStateIncrease((PetscObject)V);
304: return(0);
305: }
307: /*@
308: BVInsertVecs - Insert a set of vectors into the specified columns.
310: Collective on V
312: Input Parameters:
313: + V - basis vectors
314: . s - first column of V to be overwritten
315: . W - set of vectors to be copied
316: - orth - flag indicating if the vectors must be orthogonalized
318: Input/Output Parameter:
319: . m - number of input vectors, on output the number of linearly independent
320: vectors
322: Notes:
323: Copies the contents of vectors W to V(:,s:s+n). If the orthogonalization
324: flag is set, then the vectors are copied one by one and then orthogonalized
325: against the previous ones. If any of them is linearly dependent then it
326: is discarded and the value of m is decreased.
328: Level: intermediate
330: .seealso: BVInsertVec(), BVOrthogonalizeColumn()
331: @*/
332: PetscErrorCode BVInsertVecs(BV V,PetscInt s,PetscInt *m,Vec *W,PetscBool orth)333: {
335: PetscInt n,N,i,ndep;
336: PetscBool lindep;
337: PetscReal norm;
338: Vec v;
345: if (!*m) return(0);
346: if (*m<0) SETERRQ1(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",*m);
351: BVCheckSizes(V,1);
354: VecGetSize(*W,&N);
355: VecGetLocalSize(*W,&n);
356: 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);
357: 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);
358: 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);
360: ndep = 0;
361: for (i=0;i<*m;i++) {
362: BVGetColumn(V,s+i-ndep,&v);
363: VecCopy(W[i],v);
364: BVRestoreColumn(V,s+i-ndep,&v);
365: if (orth) {
366: BVOrthogonalizeColumn(V,s+i-ndep,NULL,&norm,&lindep);
367: if (norm==0.0 || lindep) {
368: PetscInfo1(V,"Removing linearly dependent vector %D\n",i);
369: ndep++;
370: } else {
371: BVScaleColumn(V,s+i-ndep,1.0/norm);
372: }
373: }
374: }
375: *m -= ndep;
376: PetscObjectStateIncrease((PetscObject)V);
377: return(0);
378: }
380: /*@
381: BVInsertConstraints - Insert a set of vectors as constraints.
383: Collective on V
385: Input Parameters:
386: + V - basis vectors
387: - C - set of vectors to be inserted as constraints
389: Input/Output Parameter:
390: . nc - number of input vectors, on output the number of linearly independent
391: vectors
393: Notes:
394: The constraints are relevant only during orthogonalization. Constraint
395: vectors span a subspace that is deflated in every orthogonalization
396: operation, so they are intended for removing those directions from the
397: orthogonal basis computed in regular BV columns.
399: Constraints are not stored in regular BV colums, but in a special part of
400: the storage. They can be accessed with negative indices in BVGetColumn().
402: This operation is DESTRUCTIVE, meaning that all data contained in the
403: columns of V is lost. This is typically invoked just after creating the BV.
404: Once a set of constraints has been set, it is not allowed to call this
405: function again.
407: The vectors are copied one by one and then orthogonalized against the
408: previous ones. If any of them is linearly dependent then it is discarded
409: and the value of nc is decreased. The behaviour is similar to BVInsertVecs().
411: Level: advanced
413: .seealso: BVInsertVecs(), BVOrthogonalizeColumn(), BVGetColumn(), BVGetNumConstraints()
414: @*/
415: PetscErrorCode BVInsertConstraints(BV V,PetscInt *nc,Vec *C)416: {
418: PetscInt msave;
424: if (!*nc) return(0);
425: if (*nc<0) SETERRQ1(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_OUTOFRANGE,"Number of constraints (given %D) cannot be negative",*nc);
429: BVCheckSizes(V,1);
431: if (V->issplit) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Operation not permitted for a BV obtained from BVGetSplit");
432: if (V->nc) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_ARG_WRONGSTATE,"Constraints already present in this BV object");
433: if (V->ci[0]!=-1 || V->ci[1]!=-1) SETERRQ(PetscObjectComm((PetscObject)V),PETSC_ERR_SUP,"Cannot call BVInsertConstraints after BVGetColumn");
435: msave = V->m;
436: BVResize(V,*nc+V->m,PETSC_FALSE);
437: BVInsertVecs(V,0,nc,C,PETSC_TRUE);
438: V->nc = *nc;
439: V->m = msave;
440: V->ci[0] = -V->nc-1;
441: V->ci[1] = -V->nc-1;
442: PetscObjectStateIncrease((PetscObject)V);
443: return(0);
444: }
446: /*@C
447: BVSetOptionsPrefix - Sets the prefix used for searching for all
448: BV options in the database.
450: Logically Collective on bv
452: Input Parameters:
453: + bv - the basis vectors context
454: - prefix - the prefix string to prepend to all BV option requests
456: Notes:
457: A hyphen (-) must NOT be given at the beginning of the prefix name.
458: The first character of all runtime options is AUTOMATICALLY the
459: hyphen.
461: Level: advanced
463: .seealso: BVAppendOptionsPrefix(), BVGetOptionsPrefix()
464: @*/
465: PetscErrorCode BVSetOptionsPrefix(BV bv,const char *prefix)466: {
471: PetscObjectSetOptionsPrefix((PetscObject)bv,prefix);
472: return(0);
473: }
475: /*@C
476: BVAppendOptionsPrefix - Appends to the prefix used for searching for all
477: BV options in the database.
479: Logically Collective on bv
481: Input Parameters:
482: + bv - the basis vectors context
483: - prefix - the prefix string to prepend to all BV option requests
485: Notes:
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: BVSetOptionsPrefix(), BVGetOptionsPrefix()
493: @*/
494: PetscErrorCode BVAppendOptionsPrefix(BV bv,const char *prefix)495: {
500: PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix);
501: return(0);
502: }
504: /*@C
505: BVGetOptionsPrefix - Gets the prefix used for searching for all
506: BV options in the database.
508: Not Collective
510: Input Parameters:
511: . bv - the basis vectors context
513: Output Parameters:
514: . prefix - pointer to the prefix string used, is returned
516: Note:
517: On the Fortran side, the user should pass in a string 'prefix' of
518: sufficient length to hold the prefix.
520: Level: advanced
522: .seealso: BVSetOptionsPrefix(), BVAppendOptionsPrefix()
523: @*/
524: PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])525: {
531: PetscObjectGetOptionsPrefix((PetscObject)bv,prefix);
532: return(0);
533: }
535: static PetscErrorCode BVView_Default(BV bv,PetscViewer viewer)536: {
537: PetscErrorCode ierr;
538: PetscInt j;
539: Vec v;
540: PetscViewerFormat format;
541: PetscBool isascii,ismatlab=PETSC_FALSE;
542: const char *bvname,*name;
545: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
546: if (isascii) {
547: PetscViewerGetFormat(viewer,&format);
548: if (format == PETSC_VIEWER_ASCII_MATLAB) ismatlab = PETSC_TRUE;
549: }
550: if (ismatlab) {
551: PetscObjectGetName((PetscObject)bv,&bvname);
552: PetscViewerASCIIPrintf(viewer,"%s=[];\n",bvname);
553: }
554: for (j=-bv->nc;j<bv->m;j++) {
555: BVGetColumn(bv,j,&v);
556: VecView(v,viewer);
557: if (ismatlab) {
558: PetscObjectGetName((PetscObject)v,&name);
559: PetscViewerASCIIPrintf(viewer,"%s=[%s,%s];clear %s\n",bvname,bvname,name,name);
560: }
561: BVRestoreColumn(bv,j,&v);
562: }
563: return(0);
564: }
566: /*@C
567: BVView - Prints the BV data structure.
569: Collective on bv
571: Input Parameters:
572: + bv - the BV context
573: - viewer - optional visualization context
575: Note:
576: The available visualization contexts include
577: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
578: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
579: output where only the first processor opens
580: the file. All other processors send their
581: data to the first processor to print.
583: The user can open an alternative visualization contexts with
584: PetscViewerASCIIOpen() (output to a specified file).
586: Level: beginner
587: @*/
588: PetscErrorCode BVView(BV bv,PetscViewer viewer)589: {
590: PetscErrorCode ierr;
591: PetscBool isascii;
592: PetscViewerFormat format;
593: const char *orthname[2] = {"classical","modified"};
594: const char *refname[3] = {"if needed","never","always"};
598: if (!viewer) {
599: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer);
600: }
603: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
604: if (isascii) {
605: PetscObjectPrintClassNamePrefixType((PetscObject)bv,viewer);
606: PetscViewerGetFormat(viewer,&format);
607: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
608: PetscViewerASCIIPrintf(viewer," %D columns of global length %D%s\n",bv->m,bv->N,bv->cuda?" (CUDA)":"");
609: if (bv->nc>0) {
610: PetscViewerASCIIPrintf(viewer," number of constraints: %D\n",bv->nc);
611: }
612: PetscViewerASCIIPrintf(viewer," vector orthogonalization method: %s Gram-Schmidt\n",orthname[bv->orthog_type]);
613: switch (bv->orthog_ref) {
614: case BV_ORTHOG_REFINE_IFNEEDED:
615: PetscViewerASCIIPrintf(viewer," orthogonalization refinement: %s (eta: %g)\n",refname[bv->orthog_ref],(double)bv->orthog_eta);
616: break;
617: case BV_ORTHOG_REFINE_NEVER:
618: case BV_ORTHOG_REFINE_ALWAYS:
619: PetscViewerASCIIPrintf(viewer," orthogonalization refinement: %s\n",refname[bv->orthog_ref]);
620: break;
621: }
622: PetscViewerASCIIPrintf(viewer," block orthogonalization method: %s\n",BVOrthogBlockTypes[bv->orthog_block]);
623: if (bv->matrix) {
624: if (bv->indef) {
625: PetscViewerASCIIPrintf(viewer," indefinite inner product\n");
626: } else {
627: PetscViewerASCIIPrintf(viewer," non-standard inner product\n");
628: }
629: PetscViewerASCIIPrintf(viewer," inner product matrix:\n");
630: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
631: PetscViewerASCIIPushTab(viewer);
632: MatView(bv->matrix,viewer);
633: PetscViewerASCIIPopTab(viewer);
634: PetscViewerPopFormat(viewer);
635: }
636: switch (bv->vmm) {
637: case BV_MATMULT_VECS:
638: PetscViewerASCIIPrintf(viewer," doing matmult as matrix-vector products\n");
639: break;
640: case BV_MATMULT_MAT:
641: PetscViewerASCIIPrintf(viewer," doing matmult as a single matrix-matrix product\n");
642: break;
643: case BV_MATMULT_MAT_SAVE:
644: PetscViewerASCIIPrintf(viewer," mat_save is deprecated, use mat\n");
645: break;
646: }
647: if (bv->rrandom) {
648: PetscViewerASCIIPrintf(viewer," generating random vectors independent of the number of processes\n");
649: }
650: if (bv->ops->view) { (*bv->ops->view)(bv,viewer); }
651: } else {
652: if (bv->ops->view) { (*bv->ops->view)(bv,viewer); }
653: else { BVView_Default(bv,viewer); }
654: }
655: } else {
656: (*bv->ops->view)(bv,viewer);
657: }
658: return(0);
659: }
661: /*@C
662: BVViewFromOptions - View from options
664: Collective on BV666: Input Parameters:
667: + bv - the basis vectors context
668: . obj - optional object
669: - name - command line option
671: Level: intermediate
673: .seealso: BVView(), BVCreate()
674: @*/
675: PetscErrorCode BVViewFromOptions(BV bv,PetscObject obj,const char name[])676: {
681: PetscObjectViewFromOptions((PetscObject)bv,obj,name);
682: return(0);
683: }
685: /*@C
686: BVRegister - Adds a new storage format to the BV package.
688: Not collective
690: Input Parameters:
691: + name - name of a new user-defined BV692: - function - routine to create context
694: Notes:
695: BVRegister() may be called multiple times to add several user-defined
696: basis vectors.
698: Level: advanced
700: .seealso: BVRegisterAll()
701: @*/
702: PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))703: {
707: BVInitializePackage();
708: PetscFunctionListAdd(&BVList,name,function);
709: return(0);
710: }
712: PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)713: {
717: if (s>bv->lwork) {
718: PetscFree(bv->work);
719: PetscMalloc1(s,&bv->work);
720: PetscLogObjectMemory((PetscObject)bv,(s-bv->lwork)*sizeof(PetscScalar));
721: bv->lwork = s;
722: }
723: return(0);
724: }
726: #if defined(PETSC_USE_DEBUG)
727: /*
728: SlepcDebugBVView - partially view a BV object, to be used from within a debugger.
730: ini, end: columns to be viewed
731: s: name of Matlab variable
732: filename: optionally write output to a file
733: */
734: PETSC_UNUSED PetscErrorCode SlepcDebugBVView(BV bv,PetscInt ini,PetscInt end,const char *s,const char *filename)735: {
737: PetscInt N,m;
738: PetscScalar *array;
741: BVGetArray(bv,&array);
742: BVGetSizes(bv,NULL,&N,&m);
743: SlepcDebugViewMatrix(N,end-ini+1,array+ini*N,NULL,N,s,filename);
744: BVRestoreArray(bv,&array);
745: return(0);
746: }
747: #endif