Actual source code: svdview.c
slepc-3.16.1 2021-11-17
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2021, 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: The SVD routines related to various viewers
12: */
14: #include <slepc/private/svdimpl.h>
15: #include <petscdraw.h>
17: /*@C
18: SVDView - Prints the SVD data structure.
20: Collective on svd
22: Input Parameters:
23: + svd - the singular value solver context
24: - viewer - optional visualization context
26: Options Database Key:
27: . -svd_view - Calls SVDView() at end of SVDSolve()
29: Note:
30: The available visualization contexts include
31: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
32: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
33: output where only the first processor opens
34: the file. All other processors send their
35: data to the first processor to print.
37: The user can open an alternative visualization context with
38: PetscViewerASCIIOpen() - output to a specified file.
40: Level: beginner
41: @*/
42: PetscErrorCode SVDView(SVD svd,PetscViewer viewer)
43: {
45: const char *type=NULL;
46: PetscBool isascii,isshell;
50: if (!viewer) {
51: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
52: }
56: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
57: if (isascii) {
58: PetscObjectPrintClassNamePrefixType((PetscObject)svd,viewer);
59: if (svd->ops->view) {
60: PetscViewerASCIIPushTab(viewer);
61: (*svd->ops->view)(svd,viewer);
62: PetscViewerASCIIPopTab(viewer);
63: }
64: if (svd->problem_type) {
65: switch (svd->problem_type) {
66: case SVD_STANDARD: type = "(standard) singular value problem"; break;
67: case SVD_GENERALIZED: type = "generalized singular value problem"; break;
68: }
69: } else type = "not yet set";
70: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
71: PetscViewerASCIIPrintf(viewer," transpose mode: %s\n",svd->impltrans?"implicit":"explicit");
72: if (svd->which == SVD_LARGEST) {
73: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: largest\n");
74: } else {
75: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: smallest\n");
76: }
77: PetscViewerASCIIPrintf(viewer," number of singular values (nsv): %D\n",svd->nsv);
78: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",svd->ncv);
79: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",svd->mpd);
80: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",svd->max_it);
81: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)svd->tol);
82: PetscViewerASCIIPrintf(viewer," convergence test: ");
83: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
84: switch (svd->conv) {
85: case SVD_CONV_ABS:
86: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
87: case SVD_CONV_REL:
88: PetscViewerASCIIPrintf(viewer,"relative to the singular value\n");break;
89: case SVD_CONV_NORM:
90: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
91: PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)svd->nrma);
92: if (svd->isgeneralized) {
93: PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)svd->nrmb);
94: }
95: PetscViewerASCIIPrintf(viewer,"\n");
96: break;
97: case SVD_CONV_MAXIT:
98: PetscViewerASCIIPrintf(viewer,"maximum number of iterations\n");break;
99: case SVD_CONV_USER:
100: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
101: }
102: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
103: if (svd->nini) {
104: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(svd->nini));
105: }
106: if (svd->ninil) {
107: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %D\n",PetscAbs(svd->ninil));
108: }
109: } else {
110: if (svd->ops->view) {
111: (*svd->ops->view)(svd,viewer);
112: }
113: }
114: PetscObjectTypeCompareAny((PetscObject)svd,&isshell,SVDCROSS,SVDCYCLIC,SVDSCALAPACK,SVDELEMENTAL,SVDPRIMME,"");
115: if (!isshell) {
116: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
117: if (!svd->V) { SVDGetBV(svd,&svd->V,NULL); }
118: BVView(svd->V,viewer);
119: if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
120: DSView(svd->ds,viewer);
121: PetscViewerPopFormat(viewer);
122: }
123: return(0);
124: }
126: /*@C
127: SVDViewFromOptions - View from options
129: Collective on SVD
131: Input Parameters:
132: + svd - the singular value solver context
133: . obj - optional object
134: - name - command line option
136: Level: intermediate
138: .seealso: SVDView(), SVDCreate()
139: @*/
140: PetscErrorCode SVDViewFromOptions(SVD svd,PetscObject obj,const char name[])
141: {
146: PetscObjectViewFromOptions((PetscObject)svd,obj,name);
147: return(0);
148: }
150: /*@C
151: SVDConvergedReasonView - Displays the reason an SVD solve converged or diverged.
153: Collective on svd
155: Input Parameters:
156: + svd - the singular value solver context
157: - viewer - the viewer to display the reason
159: Options Database Keys:
160: . -svd_converged_reason - print reason for convergence, and number of iterations
162: Note:
163: To change the format of the output call PetscViewerPushFormat(viewer,format) before
164: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
165: display a reason if it fails. The latter can be set in the command line with
166: -svd_converged_reason ::failed
168: Level: intermediate
170: .seealso: SVDSetTolerances(), SVDGetIterationNumber(), SVDConvergedReasonViewFromOptions()
171: @*/
172: PetscErrorCode SVDConvergedReasonView(SVD svd,PetscViewer viewer)
173: {
174: PetscErrorCode ierr;
175: PetscBool isAscii;
176: PetscViewerFormat format;
179: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)svd));
180: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
181: if (isAscii) {
182: PetscViewerGetFormat(viewer,&format);
183: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
184: if (svd->reason > 0 && format != PETSC_VIEWER_FAILED) {
185: PetscViewerASCIIPrintf(viewer,"%s SVD solve converged (%D singular triplet%s) due to %s; iterations %D\n",((PetscObject)svd)->prefix?((PetscObject)svd)->prefix:"",svd->nconv,(svd->nconv>1)?"s":"",SVDConvergedReasons[svd->reason],svd->its);
186: } else if (svd->reason <= 0) {
187: PetscViewerASCIIPrintf(viewer,"%s SVD solve did not converge due to %s; iterations %D\n",((PetscObject)svd)->prefix?((PetscObject)svd)->prefix:"",SVDConvergedReasons[svd->reason],svd->its);
188: }
189: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
190: }
191: return(0);
192: }
194: /*@
195: SVDConvergedReasonViewFromOptions - Processes command line options to determine if/how
196: the SVD converged reason is to be viewed.
198: Collective on svd
200: Input Parameter:
201: . svd - the singular value solver context
203: Level: developer
205: .seealso: SVDConvergedReasonView()
206: @*/
207: PetscErrorCode SVDConvergedReasonViewFromOptions(SVD svd)
208: {
209: PetscErrorCode ierr;
210: PetscViewer viewer;
211: PetscBool flg;
212: static PetscBool incall = PETSC_FALSE;
213: PetscViewerFormat format;
216: if (incall) return(0);
217: incall = PETSC_TRUE;
218: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_converged_reason",&viewer,&format,&flg);
219: if (flg) {
220: PetscViewerPushFormat(viewer,format);
221: SVDConvergedReasonView(svd,viewer);
222: PetscViewerPopFormat(viewer);
223: PetscViewerDestroy(&viewer);
224: }
225: incall = PETSC_FALSE;
226: return(0);
227: }
229: static PetscErrorCode SVDErrorView_ASCII(SVD svd,SVDErrorType etype,PetscViewer viewer)
230: {
231: PetscReal error,sigma;
232: PetscInt i,j;
236: if (svd->nconv<svd->nsv) {
237: PetscViewerASCIIPrintf(viewer," Problem: less than %D singular values converged\n\n",svd->nsv);
238: return(0);
239: }
240: for (i=0;i<svd->nsv;i++) {
241: SVDComputeError(svd,i,etype,&error);
242: if (error>=5.0*svd->tol) {
243: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",svd->nsv);
244: return(0);
245: }
246: }
247: PetscViewerASCIIPrintf(viewer," All requested singular values computed up to the required tolerance:");
248: for (i=0;i<=(svd->nsv-1)/8;i++) {
249: PetscViewerASCIIPrintf(viewer,"\n ");
250: for (j=0;j<PetscMin(8,svd->nsv-8*i);j++) {
251: SVDGetSingularTriplet(svd,8*i+j,&sigma,NULL,NULL);
252: PetscViewerASCIIPrintf(viewer,"%.5f",(double)sigma);
253: if (8*i+j+1<svd->nsv) { PetscViewerASCIIPrintf(viewer,", "); }
254: }
255: }
256: PetscViewerASCIIPrintf(viewer,"\n\n");
257: return(0);
258: }
260: static PetscErrorCode SVDErrorView_DETAIL(SVD svd,SVDErrorType etype,PetscViewer viewer)
261: {
263: PetscReal error,sigma;
264: PetscInt i;
265: char ex[30],sep[]=" ---------------------- --------------------\n";
268: if (!svd->nconv) return(0);
269: switch (etype) {
270: case SVD_ERROR_ABSOLUTE:
271: PetscSNPrintf(ex,sizeof(ex)," absolute error");
272: break;
273: case SVD_ERROR_RELATIVE:
274: PetscSNPrintf(ex,sizeof(ex)," relative error");
275: break;
276: }
277: PetscViewerASCIIPrintf(viewer,"%s sigma %s\n%s",sep,ex,sep);
278: for (i=0;i<svd->nconv;i++) {
279: SVDGetSingularTriplet(svd,i,&sigma,NULL,NULL);
280: SVDComputeError(svd,i,etype,&error);
281: PetscViewerASCIIPrintf(viewer," % 6f %12g\n",(double)sigma,(double)error);
282: }
283: PetscViewerASCIIPrintf(viewer,"%s",sep);
284: return(0);
285: }
287: static PetscErrorCode SVDErrorView_MATLAB(SVD svd,SVDErrorType etype,PetscViewer viewer)
288: {
290: PetscReal error;
291: PetscInt i;
292: const char *name;
295: PetscObjectGetName((PetscObject)svd,&name);
296: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
297: for (i=0;i<svd->nconv;i++) {
298: SVDComputeError(svd,i,etype,&error);
299: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
300: }
301: PetscViewerASCIIPrintf(viewer,"];\n");
302: return(0);
303: }
305: /*@C
306: SVDErrorView - Displays the errors associated with the computed solution
307: (as well as the singular values).
309: Collective on svd
311: Input Parameters:
312: + svd - the singular value solver context
313: . etype - error type
314: - viewer - optional visualization context
316: Options Database Key:
317: + -svd_error_absolute - print absolute errors of each singular triplet
318: - -svd_error_relative - print relative errors of each singular triplet
320: Notes:
321: By default, this function checks the error of all singular triplets and prints
322: the singular values if all of them are below the requested tolerance.
323: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
324: singular values and corresponding errors is printed.
326: Level: intermediate
328: .seealso: SVDSolve(), SVDValuesView(), SVDVectorsView()
329: @*/
330: PetscErrorCode SVDErrorView(SVD svd,SVDErrorType etype,PetscViewer viewer)
331: {
332: PetscBool isascii;
333: PetscViewerFormat format;
334: PetscErrorCode ierr;
338: if (!viewer) {
339: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
340: }
343: SVDCheckSolved(svd,1);
344: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
345: if (!isascii) return(0);
347: PetscViewerGetFormat(viewer,&format);
348: switch (format) {
349: case PETSC_VIEWER_DEFAULT:
350: case PETSC_VIEWER_ASCII_INFO:
351: SVDErrorView_ASCII(svd,etype,viewer);
352: break;
353: case PETSC_VIEWER_ASCII_INFO_DETAIL:
354: SVDErrorView_DETAIL(svd,etype,viewer);
355: break;
356: case PETSC_VIEWER_ASCII_MATLAB:
357: SVDErrorView_MATLAB(svd,etype,viewer);
358: break;
359: default:
360: PetscInfo1(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
361: }
362: return(0);
363: }
365: /*@
366: SVDErrorViewFromOptions - Processes command line options to determine if/how
367: the errors of the computed solution are to be viewed.
369: Collective on svd
371: Input Parameter:
372: . svd - the singular value solver context
374: Level: developer
375: @*/
376: PetscErrorCode SVDErrorViewFromOptions(SVD svd)
377: {
378: PetscErrorCode ierr;
379: PetscViewer viewer;
380: PetscBool flg;
381: static PetscBool incall = PETSC_FALSE;
382: PetscViewerFormat format;
385: if (incall) return(0);
386: incall = PETSC_TRUE;
387: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_absolute",&viewer,&format,&flg);
388: if (flg) {
389: PetscViewerPushFormat(viewer,format);
390: SVDErrorView(svd,SVD_ERROR_ABSOLUTE,viewer);
391: PetscViewerPopFormat(viewer);
392: PetscViewerDestroy(&viewer);
393: }
394: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_relative",&viewer,&format,&flg);
395: if (flg) {
396: PetscViewerPushFormat(viewer,format);
397: SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer);
398: PetscViewerPopFormat(viewer);
399: PetscViewerDestroy(&viewer);
400: }
401: incall = PETSC_FALSE;
402: return(0);
403: }
405: static PetscErrorCode SVDValuesView_DRAW(SVD svd,PetscViewer viewer)
406: {
408: PetscDraw draw;
409: PetscDrawSP drawsp;
410: PetscReal re,im=0.0;
411: PetscInt i;
414: if (!svd->nconv) return(0);
415: PetscViewerDrawGetDraw(viewer,0,&draw);
416: PetscDrawSetTitle(draw,"Computed singular values");
417: PetscDrawSPCreate(draw,1,&drawsp);
418: for (i=0;i<svd->nconv;i++) {
419: re = svd->sigma[svd->perm[i]];
420: PetscDrawSPAddPoint(drawsp,&re,&im);
421: }
422: PetscDrawSPDraw(drawsp,PETSC_TRUE);
423: PetscDrawSPSave(drawsp);
424: PetscDrawSPDestroy(&drawsp);
425: return(0);
426: }
428: static PetscErrorCode SVDValuesView_BINARY(SVD svd,PetscViewer viewer)
429: {
430: PetscInt i,k;
431: PetscReal *sv;
435: PetscMalloc1(svd->nconv,&sv);
436: for (i=0;i<svd->nconv;i++) {
437: k = svd->perm[i];
438: sv[i] = svd->sigma[k];
439: }
440: PetscViewerBinaryWrite(viewer,sv,svd->nconv,PETSC_REAL);
441: PetscFree(sv);
442: return(0);
443: }
445: #if defined(PETSC_HAVE_HDF5)
446: static PetscErrorCode SVDValuesView_HDF5(SVD svd,PetscViewer viewer)
447: {
449: PetscInt i,k,n,N;
450: PetscMPIInt rank;
451: Vec v;
452: char vname[30];
453: const char *ename;
456: MPI_Comm_rank(PetscObjectComm((PetscObject)svd),&rank);
457: N = svd->nconv;
458: n = rank? 0: N;
459: /* create a vector containing the singular values */
460: VecCreateMPI(PetscObjectComm((PetscObject)svd),n,N,&v);
461: PetscObjectGetName((PetscObject)svd,&ename);
462: PetscSNPrintf(vname,sizeof(vname),"sigma_%s",ename);
463: PetscObjectSetName((PetscObject)v,vname);
464: if (!rank) {
465: for (i=0;i<svd->nconv;i++) {
466: k = svd->perm[i];
467: VecSetValue(v,i,svd->sigma[k],INSERT_VALUES);
468: }
469: }
470: VecAssemblyBegin(v);
471: VecAssemblyEnd(v);
472: VecView(v,viewer);
473: VecDestroy(&v);
474: return(0);
475: }
476: #endif
478: static PetscErrorCode SVDValuesView_ASCII(SVD svd,PetscViewer viewer)
479: {
480: PetscInt i;
484: PetscViewerASCIIPrintf(viewer,"Singular values = \n");
485: for (i=0;i<svd->nconv;i++) {
486: PetscViewerASCIIPrintf(viewer," %.5f\n",(double)svd->sigma[svd->perm[i]]);
487: }
488: PetscViewerASCIIPrintf(viewer,"\n");
489: return(0);
490: }
492: static PetscErrorCode SVDValuesView_MATLAB(SVD svd,PetscViewer viewer)
493: {
495: PetscInt i;
496: const char *name;
499: PetscObjectGetName((PetscObject)svd,&name);
500: PetscViewerASCIIPrintf(viewer,"Sigma_%s = [\n",name);
501: for (i=0;i<svd->nconv;i++) {
502: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)svd->sigma[svd->perm[i]]);
503: }
504: PetscViewerASCIIPrintf(viewer,"];\n");
505: return(0);
506: }
508: /*@C
509: SVDValuesView - Displays the computed singular values in a viewer.
511: Collective on svd
513: Input Parameters:
514: + svd - the singular value solver context
515: - viewer - the viewer
517: Options Database Key:
518: . -svd_view_values - print computed singular values
520: Level: intermediate
522: .seealso: SVDSolve(), SVDVectorsView(), SVDErrorView()
523: @*/
524: PetscErrorCode SVDValuesView(SVD svd,PetscViewer viewer)
525: {
526: PetscBool isascii,isdraw,isbinary;
527: PetscViewerFormat format;
528: PetscErrorCode ierr;
529: #if defined(PETSC_HAVE_HDF5)
530: PetscBool ishdf5;
531: #endif
535: if (!viewer) {
536: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
537: }
540: SVDCheckSolved(svd,1);
541: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
542: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
543: #if defined(PETSC_HAVE_HDF5)
544: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
545: #endif
546: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
547: if (isdraw) {
548: SVDValuesView_DRAW(svd,viewer);
549: } else if (isbinary) {
550: SVDValuesView_BINARY(svd,viewer);
551: #if defined(PETSC_HAVE_HDF5)
552: } else if (ishdf5) {
553: SVDValuesView_HDF5(svd,viewer);
554: #endif
555: } else if (isascii) {
556: PetscViewerGetFormat(viewer,&format);
557: switch (format) {
558: case PETSC_VIEWER_DEFAULT:
559: case PETSC_VIEWER_ASCII_INFO:
560: case PETSC_VIEWER_ASCII_INFO_DETAIL:
561: SVDValuesView_ASCII(svd,viewer);
562: break;
563: case PETSC_VIEWER_ASCII_MATLAB:
564: SVDValuesView_MATLAB(svd,viewer);
565: break;
566: default:
567: PetscInfo1(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
568: }
569: }
570: return(0);
571: }
573: /*@
574: SVDValuesViewFromOptions - Processes command line options to determine if/how
575: the computed singular values are to be viewed.
577: Collective on svd
579: Input Parameter:
580: . svd - the singular value solver context
582: Level: developer
583: @*/
584: PetscErrorCode SVDValuesViewFromOptions(SVD svd)
585: {
586: PetscErrorCode ierr;
587: PetscViewer viewer;
588: PetscBool flg;
589: static PetscBool incall = PETSC_FALSE;
590: PetscViewerFormat format;
593: if (incall) return(0);
594: incall = PETSC_TRUE;
595: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_values",&viewer,&format,&flg);
596: if (flg) {
597: PetscViewerPushFormat(viewer,format);
598: SVDValuesView(svd,viewer);
599: PetscViewerPopFormat(viewer);
600: PetscViewerDestroy(&viewer);
601: }
602: incall = PETSC_FALSE;
603: return(0);
604: }
606: /*@C
607: SVDVectorsView - Outputs computed singular vectors to a viewer.
609: Collective on svd
611: Input Parameters:
612: + svd - the singular value solver context
613: - viewer - the viewer
615: Options Database Keys:
616: . -svd_view_vectors - output singular vectors
618: Note:
619: Right and left singular vectors are interleaved, that is, the vectors are
620: output in the following order: V0, U0, V1, U1, V2, U2, ...
622: Level: intermediate
624: .seealso: SVDSolve(), SVDValuesView(), SVDErrorView()
625: @*/
626: PetscErrorCode SVDVectorsView(SVD svd,PetscViewer viewer)
627: {
629: PetscInt i,k;
630: Vec x;
631: char vname[30];
632: const char *ename;
636: if (!viewer) {
637: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
638: }
641: SVDCheckSolved(svd,1);
642: if (svd->nconv) {
643: PetscObjectGetName((PetscObject)svd,&ename);
644: SVDComputeVectors(svd);
645: for (i=0;i<svd->nconv;i++) {
646: k = svd->perm[i];
647: PetscSNPrintf(vname,sizeof(vname),"V%d_%s",(int)i,ename);
648: BVGetColumn(svd->V,k,&x);
649: PetscObjectSetName((PetscObject)x,vname);
650: VecView(x,viewer);
651: BVRestoreColumn(svd->V,k,&x);
652: PetscSNPrintf(vname,sizeof(vname),"U%d_%s",(int)i,ename);
653: BVGetColumn(svd->U,k,&x);
654: PetscObjectSetName((PetscObject)x,vname);
655: VecView(x,viewer);
656: BVRestoreColumn(svd->U,k,&x);
657: }
658: }
659: return(0);
660: }
662: /*@
663: SVDVectorsViewFromOptions - Processes command line options to determine if/how
664: the computed singular vectors are to be viewed.
666: Collective on svd
668: Input Parameter:
669: . svd - the singular value solver context
671: Level: developer
672: @*/
673: PetscErrorCode SVDVectorsViewFromOptions(SVD svd)
674: {
675: PetscErrorCode ierr;
676: PetscViewer viewer;
677: PetscBool flg = PETSC_FALSE;
678: static PetscBool incall = PETSC_FALSE;
679: PetscViewerFormat format;
682: if (incall) return(0);
683: incall = PETSC_TRUE;
684: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_vectors",&viewer,&format,&flg);
685: if (flg) {
686: PetscViewerPushFormat(viewer,format);
687: SVDVectorsView(svd,viewer);
688: PetscViewerPopFormat(viewer);
689: PetscViewerDestroy(&viewer);
690: }
691: incall = PETSC_FALSE;
692: return(0);
693: }