Actual source code: bvfunc.c
slepc-3.17.1 2022-04-11
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_",0};
22: const char *BVOrthogRefineTypes[] = {"IFNEEDED","NEVER","ALWAYS","BVOrthogRefineType","BV_ORTHOG_REFINE_",0};
23: const char *BVOrthogBlockTypes[] = {"GS","CHOL","TSQR","TSQRCHOL","SVQB","BVOrthogBlockType","BV_ORTHOG_BLOCK_",0};
24: const char *BVMatMultTypes[] = {"VECS","MAT","MAT_SAVE","BVMatMultType","BV_MATMULT_",0};
25: const char *BVSVDMethods[] = {"REFINE","QR","QR_CAA","BVSVDMethod","BV_SVD_METHOD_",0};
27: /*@C
28: BVFinalizePackage - This function destroys everything in the Slepc interface
29: to the BV package. It is called from SlepcFinalize().
31: Level: developer
33: .seealso: SlepcFinalize()
34: @*/
35: PetscErrorCode BVFinalizePackage(void)
36: {
37: PetscFunctionListDestroy(&BVList);
38: MPI_Op_free(&MPIU_TSQR);
39: MPI_Op_free(&MPIU_LAPY2);
40: BVPackageInitialized = PETSC_FALSE;
41: BVRegisterAllCalled = PETSC_FALSE;
42: PetscFunctionReturn(0);
43: }
45: /*@C
46: BVInitializePackage - This function initializes everything in the BV package.
47: It is called from PetscDLLibraryRegister() when using dynamic libraries, and
48: on the first call to BVCreate() when using static libraries.
50: Level: developer
52: .seealso: SlepcInitialize()
53: @*/
54: PetscErrorCode BVInitializePackage(void)
55: {
56: char logList[256];
57: PetscBool opt,pkg;
58: PetscClassId classids[1];
60: if (BVPackageInitialized) PetscFunctionReturn(0);
61: BVPackageInitialized = PETSC_TRUE;
62: /* Register Classes */
63: PetscClassIdRegister("Basis Vectors",&BV_CLASSID);
64: /* Register Constructors */
65: BVRegisterAll();
66: /* Register Events */
67: PetscLogEventRegister("BVCreate",BV_CLASSID,&BV_Create);
68: PetscLogEventRegister("BVCopy",BV_CLASSID,&BV_Copy);
69: PetscLogEventRegister("BVMult",BV_CLASSID,&BV_Mult);
70: PetscLogEventRegister("BVMultVec",BV_CLASSID,&BV_MultVec);
71: PetscLogEventRegister("BVMultInPlace",BV_CLASSID,&BV_MultInPlace);
72: PetscLogEventRegister("BVDot",BV_CLASSID,&BV_Dot);
73: PetscLogEventRegister("BVDotVec",BV_CLASSID,&BV_DotVec);
74: PetscLogEventRegister("BVOrthogonalize",BV_CLASSID,&BV_Orthogonalize);
75: PetscLogEventRegister("BVOrthogonalizeV",BV_CLASSID,&BV_OrthogonalizeVec);
76: PetscLogEventRegister("BVScale",BV_CLASSID,&BV_Scale);
77: PetscLogEventRegister("BVNorm",BV_CLASSID,&BV_Norm);
78: PetscLogEventRegister("BVNormVec",BV_CLASSID,&BV_NormVec);
79: PetscLogEventRegister("BVNormalize",BV_CLASSID,&BV_Normalize);
80: PetscLogEventRegister("BVSetRandom",BV_CLASSID,&BV_SetRandom);
81: PetscLogEventRegister("BVMatMult",BV_CLASSID,&BV_MatMult);
82: PetscLogEventRegister("BVMatMultVec",BV_CLASSID,&BV_MatMultVec);
83: PetscLogEventRegister("BVMatProject",BV_CLASSID,&BV_MatProject);
84: PetscLogEventRegister("BVSVDAndRank",BV_CLASSID,&BV_SVDAndRank);
85: /* MPI reduction operation used in BVOrthogonalize */
86: MPI_Op_create(SlepcGivensPacked,PETSC_FALSE,&MPIU_TSQR);
87: MPI_Op_create(SlepcPythag,PETSC_TRUE,&MPIU_LAPY2);
88: /* Process Info */
89: classids[0] = BV_CLASSID;
90: PetscInfoProcessClass("bv",1,&classids[0]);
91: /* Process summary exclusions */
92: PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt);
93: if (opt) {
94: PetscStrInList("bv",logList,',',&pkg);
95: if (pkg) PetscLogEventDeactivateClass(BV_CLASSID);
96: }
97: /* Register package finalizer */
98: PetscRegisterFinalize(BVFinalizePackage);
99: PetscFunctionReturn(0);
100: }
102: /*@C
103: BVDestroy - Destroys BV context that was created with BVCreate().
105: Collective on bv
107: Input Parameter:
108: . bv - the basis vectors context
110: Level: beginner
112: .seealso: BVCreate()
113: @*/
114: PetscErrorCode BVDestroy(BV *bv)
115: {
116: if (!*bv) PetscFunctionReturn(0);
119: if (--((PetscObject)(*bv))->refct > 0) { *bv = 0; PetscFunctionReturn(0); }
120: if ((*bv)->ops->destroy) (*(*bv)->ops->destroy)(*bv);
121: VecDestroy(&(*bv)->t);
122: MatDestroy(&(*bv)->matrix);
123: VecDestroy(&(*bv)->Bx);
124: VecDestroy(&(*bv)->buffer);
125: BVDestroy(&(*bv)->cached);
126: BVDestroy(&(*bv)->L);
127: BVDestroy(&(*bv)->R);
128: PetscFree((*bv)->work);
129: PetscFree2((*bv)->h,(*bv)->c);
130: VecDestroy(&(*bv)->omega);
131: MatDestroy(&(*bv)->Acreate);
132: MatDestroy(&(*bv)->Aget);
133: MatDestroy(&(*bv)->Abuffer);
134: PetscRandomDestroy(&(*bv)->rand);
135: PetscHeaderDestroy(bv);
136: PetscFunctionReturn(0);
137: }
139: /*@
140: BVCreate - Creates a basis vectors context.
142: Collective
144: Input Parameter:
145: . comm - MPI communicator
147: Output Parameter:
148: . newbv - location to put the basis vectors context
150: Level: beginner
152: .seealso: BVSetUp(), BVDestroy(), BV
153: @*/
154: PetscErrorCode BVCreate(MPI_Comm comm,BV *newbv)
155: {
156: BV bv;
159: *newbv = 0;
160: BVInitializePackage();
161: SlepcHeaderCreate(bv,BV_CLASSID,"BV","Basis Vectors","BV",comm,BVDestroy,BVView);
163: bv->t = NULL;
164: bv->n = -1;
165: bv->N = -1;
166: bv->m = 0;
167: bv->l = 0;
168: bv->k = 0;
169: bv->nc = 0;
170: bv->orthog_type = BV_ORTHOG_CGS;
171: bv->orthog_ref = BV_ORTHOG_REFINE_IFNEEDED;
172: bv->orthog_eta = 0.7071;
173: bv->orthog_block = BV_ORTHOG_BLOCK_GS;
174: bv->matrix = NULL;
175: bv->indef = PETSC_FALSE;
176: bv->vmm = BV_MATMULT_MAT;
177: bv->rrandom = PETSC_FALSE;
178: bv->deftol = 10*PETSC_MACHINE_EPSILON;
180: bv->Bx = NULL;
181: bv->buffer = NULL;
182: bv->Abuffer = NULL;
183: bv->xid = 0;
184: bv->xstate = 0;
185: bv->cv[0] = NULL;
186: bv->cv[1] = NULL;
187: bv->ci[0] = -1;
188: bv->ci[1] = -1;
189: bv->st[0] = -1;
190: bv->st[1] = -1;
191: bv->id[0] = 0;
192: bv->id[1] = 0;
193: bv->h = NULL;
194: bv->c = NULL;
195: bv->omega = NULL;
196: bv->defersfo = PETSC_FALSE;
197: bv->cached = NULL;
198: bv->bvstate = 0;
199: bv->L = NULL;
200: bv->R = NULL;
201: bv->lstate = 0;
202: bv->rstate = 0;
203: bv->lsplit = 0;
204: bv->issplit = 0;
205: bv->splitparent = NULL;
206: bv->rand = NULL;
207: bv->rrandom = PETSC_FALSE;
208: bv->Acreate = NULL;
209: bv->Aget = NULL;
210: bv->cuda = PETSC_FALSE;
211: bv->sfocalled = PETSC_FALSE;
212: bv->work = NULL;
213: bv->lwork = 0;
214: bv->data = NULL;
216: *newbv = bv;
217: PetscFunctionReturn(0);
218: }
220: /*@
221: BVCreateFromMat - Creates a basis vectors object from a dense Mat object.
223: Collective on A
225: Input Parameter:
226: . A - a dense tall-skinny matrix
228: Output Parameter:
229: . bv - the new basis vectors context
231: Notes:
232: The matrix values are copied to the BV data storage, memory is not shared.
234: The communicator of the BV object will be the same as A, and so will be
235: the dimensions.
237: Level: intermediate
239: .seealso: BVCreate(), BVDestroy(), BVCreateMat()
240: @*/
241: PetscErrorCode BVCreateFromMat(Mat A,BV *bv)
242: {
243: PetscInt n,N,k;
248: MatGetSize(A,&N,&k);
249: MatGetLocalSize(A,&n,NULL);
250: BVCreate(PetscObjectComm((PetscObject)A),bv);
251: BVSetSizes(*bv,n,N,k);
253: (*bv)->Acreate = A;
254: PetscObjectReference((PetscObject)A);
255: PetscFunctionReturn(0);
256: }
258: /*@
259: BVInsertVec - Insert a vector into the specified column.
261: Collective on V
263: Input Parameters:
264: + V - basis vectors
265: . j - the column of V to be overwritten
266: - w - the vector to be copied
268: Level: intermediate
270: .seealso: BVInsertVecs()
271: @*/
272: PetscErrorCode BVInsertVec(BV V,PetscInt j,Vec w)
273: {
274: PetscInt n,N;
275: Vec v;
281: BVCheckSizes(V,1);
284: VecGetSize(w,&N);
285: VecGetLocalSize(w,&n);
289: BVGetColumn(V,j,&v);
290: VecCopy(w,v);
291: BVRestoreColumn(V,j,&v);
292: PetscObjectStateIncrease((PetscObject)V);
293: PetscFunctionReturn(0);
294: }
296: /*@
297: BVInsertVecs - Insert a set of vectors into the specified columns.
299: Collective on V
301: Input Parameters:
302: + V - basis vectors
303: . s - first column of V to be overwritten
304: . W - set of vectors to be copied
305: - orth - flag indicating if the vectors must be orthogonalized
307: Input/Output Parameter:
308: . m - number of input vectors, on output the number of linearly independent
309: vectors
311: Notes:
312: Copies the contents of vectors W to V(:,s:s+n). If the orthogonalization
313: flag is set, then the vectors are copied one by one and then orthogonalized
314: against the previous ones. If any of them is linearly dependent then it
315: is discarded and the value of m is decreased.
317: Level: intermediate
319: .seealso: BVInsertVec(), BVOrthogonalizeColumn()
320: @*/
321: PetscErrorCode BVInsertVecs(BV V,PetscInt s,PetscInt *m,Vec *W,PetscBool orth)
322: {
323: PetscInt n,N,i,ndep;
324: PetscBool lindep;
325: PetscReal norm;
326: Vec v;
332: if (!*m) PetscFunctionReturn(0);
338: BVCheckSizes(V,1);
341: VecGetSize(*W,&N);
342: VecGetLocalSize(*W,&n);
347: ndep = 0;
348: for (i=0;i<*m;i++) {
349: BVGetColumn(V,s+i-ndep,&v);
350: VecCopy(W[i],v);
351: BVRestoreColumn(V,s+i-ndep,&v);
352: if (orth) {
353: BVOrthogonalizeColumn(V,s+i-ndep,NULL,&norm,&lindep);
354: if (norm==0.0 || lindep) {
355: PetscInfo(V,"Removing linearly dependent vector %" PetscInt_FMT "\n",i);
356: ndep++;
357: } else BVScaleColumn(V,s+i-ndep,1.0/norm);
358: }
359: }
360: *m -= ndep;
361: PetscObjectStateIncrease((PetscObject)V);
362: PetscFunctionReturn(0);
363: }
365: /*@
366: BVInsertConstraints - Insert a set of vectors as constraints.
368: Collective on V
370: Input Parameters:
371: + V - basis vectors
372: - C - set of vectors to be inserted as constraints
374: Input/Output Parameter:
375: . nc - number of input vectors, on output the number of linearly independent
376: vectors
378: Notes:
379: The constraints are relevant only during orthogonalization. Constraint
380: vectors span a subspace that is deflated in every orthogonalization
381: operation, so they are intended for removing those directions from the
382: orthogonal basis computed in regular BV columns.
384: Constraints are not stored in regular BV columns, but in a special part of
385: the storage. They can be accessed with negative indices in BVGetColumn().
387: This operation is DESTRUCTIVE, meaning that all data contained in the
388: columns of V is lost. This is typically invoked just after creating the BV.
389: Once a set of constraints has been set, it is not allowed to call this
390: function again.
392: The vectors are copied one by one and then orthogonalized against the
393: previous ones. If any of them is linearly dependent then it is discarded
394: and the value of nc is decreased. The behaviour is similar to BVInsertVecs().
396: Level: advanced
398: .seealso: BVInsertVecs(), BVOrthogonalizeColumn(), BVGetColumn(), BVGetNumConstraints()
399: @*/
400: PetscErrorCode BVInsertConstraints(BV V,PetscInt *nc,Vec *C)
401: {
402: PetscInt msave;
407: if (!*nc) PetscFunctionReturn(0);
412: BVCheckSizes(V,1);
418: msave = V->m;
419: BVResize(V,*nc+V->m,PETSC_FALSE);
420: BVInsertVecs(V,0,nc,C,PETSC_TRUE);
421: V->nc = *nc;
422: V->m = msave;
423: V->ci[0] = -V->nc-1;
424: V->ci[1] = -V->nc-1;
425: PetscObjectStateIncrease((PetscObject)V);
426: PetscFunctionReturn(0);
427: }
429: /*@C
430: BVSetOptionsPrefix - Sets the prefix used for searching for all
431: BV options in the database.
433: Logically Collective on bv
435: Input Parameters:
436: + bv - the basis vectors context
437: - prefix - the prefix string to prepend to all BV option requests
439: Notes:
440: A hyphen (-) must NOT be given at the beginning of the prefix name.
441: The first character of all runtime options is AUTOMATICALLY the
442: hyphen.
444: Level: advanced
446: .seealso: BVAppendOptionsPrefix(), BVGetOptionsPrefix()
447: @*/
448: PetscErrorCode BVSetOptionsPrefix(BV bv,const char *prefix)
449: {
451: PetscObjectSetOptionsPrefix((PetscObject)bv,prefix);
452: PetscFunctionReturn(0);
453: }
455: /*@C
456: BVAppendOptionsPrefix - Appends to the prefix used for searching for all
457: BV options in the database.
459: Logically Collective on bv
461: Input Parameters:
462: + bv - the basis vectors context
463: - prefix - the prefix string to prepend to all BV option requests
465: Notes:
466: A hyphen (-) must NOT be given at the beginning of the prefix name.
467: The first character of all runtime options is AUTOMATICALLY the
468: hyphen.
470: Level: advanced
472: .seealso: BVSetOptionsPrefix(), BVGetOptionsPrefix()
473: @*/
474: PetscErrorCode BVAppendOptionsPrefix(BV bv,const char *prefix)
475: {
477: PetscObjectAppendOptionsPrefix((PetscObject)bv,prefix);
478: PetscFunctionReturn(0);
479: }
481: /*@C
482: BVGetOptionsPrefix - Gets the prefix used for searching for all
483: BV options in the database.
485: Not Collective
487: Input Parameters:
488: . bv - the basis vectors context
490: Output Parameters:
491: . prefix - pointer to the prefix string used, is returned
493: Note:
494: On the Fortran side, the user should pass in a string 'prefix' of
495: sufficient length to hold the prefix.
497: Level: advanced
499: .seealso: BVSetOptionsPrefix(), BVAppendOptionsPrefix()
500: @*/
501: PetscErrorCode BVGetOptionsPrefix(BV bv,const char *prefix[])
502: {
505: PetscObjectGetOptionsPrefix((PetscObject)bv,prefix);
506: PetscFunctionReturn(0);
507: }
509: static PetscErrorCode BVView_Default(BV bv,PetscViewer viewer)
510: {
511: PetscInt j;
512: Vec v;
513: PetscViewerFormat format;
514: PetscBool isascii,ismatlab=PETSC_FALSE;
515: const char *bvname,*name;
517: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
518: if (isascii) {
519: PetscViewerGetFormat(viewer,&format);
520: if (format == PETSC_VIEWER_ASCII_MATLAB) ismatlab = PETSC_TRUE;
521: }
522: if (ismatlab) {
523: PetscObjectGetName((PetscObject)bv,&bvname);
524: PetscViewerASCIIPrintf(viewer,"%s=[];\n",bvname);
525: }
526: for (j=-bv->nc;j<bv->m;j++) {
527: BVGetColumn(bv,j,&v);
528: VecView(v,viewer);
529: if (ismatlab) {
530: PetscObjectGetName((PetscObject)v,&name);
531: PetscViewerASCIIPrintf(viewer,"%s=[%s,%s];clear %s\n",bvname,bvname,name,name);
532: }
533: BVRestoreColumn(bv,j,&v);
534: }
535: PetscFunctionReturn(0);
536: }
538: /*@C
539: BVView - Prints the BV data structure.
541: Collective on bv
543: Input Parameters:
544: + bv - the BV context
545: - viewer - optional visualization context
547: Note:
548: The available visualization contexts include
549: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
550: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
551: output where only the first processor opens
552: the file. All other processors send their
553: data to the first processor to print.
555: The user can open an alternative visualization contexts with
556: PetscViewerASCIIOpen() (output to a specified file).
558: Level: beginner
560: .seealso: BVCreate()
561: @*/
562: PetscErrorCode BVView(BV bv,PetscViewer viewer)
563: {
564: PetscBool isascii;
565: PetscViewerFormat format;
566: const char *orthname[2] = {"classical","modified"};
567: const char *refname[3] = {"if needed","never","always"};
570: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)bv),&viewer);
573: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
574: if (isascii) {
575: PetscObjectPrintClassNamePrefixType((PetscObject)bv,viewer);
576: PetscViewerGetFormat(viewer,&format);
577: if (format == PETSC_VIEWER_ASCII_INFO || format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
578: PetscViewerASCIIPrintf(viewer," %" PetscInt_FMT " columns of global length %" PetscInt_FMT "%s\n",bv->m,bv->N,bv->cuda?" (CUDA)":"");
579: if (bv->nc>0) PetscViewerASCIIPrintf(viewer," number of constraints: %" PetscInt_FMT "\n",bv->nc);
580: PetscViewerASCIIPrintf(viewer," vector orthogonalization method: %s Gram-Schmidt\n",orthname[bv->orthog_type]);
581: switch (bv->orthog_ref) {
582: case BV_ORTHOG_REFINE_IFNEEDED:
583: PetscViewerASCIIPrintf(viewer," orthogonalization refinement: %s (eta: %g)\n",refname[bv->orthog_ref],(double)bv->orthog_eta);
584: break;
585: case BV_ORTHOG_REFINE_NEVER:
586: case BV_ORTHOG_REFINE_ALWAYS:
587: PetscViewerASCIIPrintf(viewer," orthogonalization refinement: %s\n",refname[bv->orthog_ref]);
588: break;
589: }
590: PetscViewerASCIIPrintf(viewer," block orthogonalization method: %s\n",BVOrthogBlockTypes[bv->orthog_block]);
591: if (bv->matrix) {
592: if (bv->indef) PetscViewerASCIIPrintf(viewer," indefinite inner product\n");
593: else PetscViewerASCIIPrintf(viewer," non-standard inner product\n");
594: PetscViewerASCIIPrintf(viewer," tolerance for definite inner product: %g\n",(double)bv->deftol);
595: PetscViewerASCIIPrintf(viewer," inner product matrix:\n");
596: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
597: PetscViewerASCIIPushTab(viewer);
598: MatView(bv->matrix,viewer);
599: PetscViewerASCIIPopTab(viewer);
600: PetscViewerPopFormat(viewer);
601: }
602: switch (bv->vmm) {
603: case BV_MATMULT_VECS:
604: PetscViewerASCIIPrintf(viewer," doing matmult as matrix-vector products\n");
605: break;
606: case BV_MATMULT_MAT:
607: PetscViewerASCIIPrintf(viewer," doing matmult as a single matrix-matrix product\n");
608: break;
609: case BV_MATMULT_MAT_SAVE:
610: PetscViewerASCIIPrintf(viewer," mat_save is deprecated, use mat\n");
611: break;
612: }
613: if (bv->rrandom) PetscViewerASCIIPrintf(viewer," generating random vectors independent of the number of processes\n");
614: if (bv->ops->view) (*bv->ops->view)(bv,viewer);
615: } else {
616: if (bv->ops->view) (*bv->ops->view)(bv,viewer);
617: else BVView_Default(bv,viewer);
618: }
619: } else (*bv->ops->view)(bv,viewer);
620: PetscFunctionReturn(0);
621: }
623: /*@C
624: BVViewFromOptions - View from options
626: Collective on BV
628: Input Parameters:
629: + bv - the basis vectors context
630: . obj - optional object
631: - name - command line option
633: Level: intermediate
635: .seealso: BVView(), BVCreate()
636: @*/
637: PetscErrorCode BVViewFromOptions(BV bv,PetscObject obj,const char name[])
638: {
640: PetscObjectViewFromOptions((PetscObject)bv,obj,name);
641: PetscFunctionReturn(0);
642: }
644: /*@C
645: BVRegister - Adds a new storage format to the BV package.
647: Not collective
649: Input Parameters:
650: + name - name of a new user-defined BV
651: - function - routine to create context
653: Notes:
654: BVRegister() may be called multiple times to add several user-defined
655: basis vectors.
657: Level: advanced
659: .seealso: BVRegisterAll()
660: @*/
661: PetscErrorCode BVRegister(const char *name,PetscErrorCode (*function)(BV))
662: {
663: BVInitializePackage();
664: PetscFunctionListAdd(&BVList,name,function);
665: PetscFunctionReturn(0);
666: }
668: PetscErrorCode BVAllocateWork_Private(BV bv,PetscInt s)
669: {
670: if (s>bv->lwork) {
671: PetscFree(bv->work);
672: PetscMalloc1(s,&bv->work);
673: PetscLogObjectMemory((PetscObject)bv,(s-bv->lwork)*sizeof(PetscScalar));
674: bv->lwork = s;
675: }
676: PetscFunctionReturn(0);
677: }
679: #if defined(PETSC_USE_DEBUG)
680: /*
681: SlepcDebugBVView - partially view a BV object, to be used from within a debugger.
683: ini, end: columns to be viewed
684: s: name of Matlab variable
685: filename: optionally write output to a file
686: */
687: PETSC_UNUSED PetscErrorCode SlepcDebugBVView(BV bv,PetscInt ini,PetscInt end,const char *s,const char *filename)
688: {
689: PetscInt N,m;
690: PetscScalar *array;
692: BVGetArray(bv,&array);
693: BVGetSizes(bv,NULL,&N,&m);
694: SlepcDebugViewMatrix(N,end-ini+1,array+ini*N,NULL,N,s,filename);
695: BVRestoreArray(bv,&array);
696: PetscFunctionReturn(0);
697: }
698: #endif