Actual source code: svdview.c
slepc-3.18.2 2023-01-26
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: 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
42: .seealso: EPSView()
43: @*/
44: PetscErrorCode SVDView(SVD svd,PetscViewer viewer)
45: {
46: const char *type=NULL;
47: PetscBool isascii,isshell;
50: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
54: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
55: if (isascii) {
56: PetscObjectPrintClassNamePrefixType((PetscObject)svd,viewer);
57: PetscViewerASCIIPushTab(viewer);
58: PetscTryTypeMethod(svd,view,viewer);
59: PetscViewerASCIIPopTab(viewer);
60: if (svd->problem_type) {
61: switch (svd->problem_type) {
62: case SVD_STANDARD: type = "(standard) singular value problem"; break;
63: case SVD_GENERALIZED: type = "generalized singular value problem"; break;
64: case SVD_HYPERBOLIC: type = "hyperbolic singular value problem"; break;
65: }
66: } else type = "not yet set";
67: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
68: PetscViewerASCIIPrintf(viewer," transpose mode: %s\n",svd->impltrans?"implicit":"explicit");
69: if (svd->which == SVD_LARGEST) PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: largest\n");
70: else PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: smallest\n");
71: PetscViewerASCIIPrintf(viewer," number of singular values (nsv): %" PetscInt_FMT "\n",svd->nsv);
72: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",svd->ncv);
73: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",svd->mpd);
74: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",svd->max_it);
75: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)svd->tol);
76: PetscViewerASCIIPrintf(viewer," convergence test: ");
77: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
78: switch (svd->conv) {
79: case SVD_CONV_ABS:
80: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
81: case SVD_CONV_REL:
82: PetscViewerASCIIPrintf(viewer,"relative to the singular value\n");break;
83: case SVD_CONV_NORM:
84: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
85: PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)svd->nrma);
86: if (svd->isgeneralized) PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)svd->nrmb);
87: PetscViewerASCIIPrintf(viewer,"\n");
88: break;
89: case SVD_CONV_MAXIT:
90: PetscViewerASCIIPrintf(viewer,"maximum number of iterations\n");break;
91: case SVD_CONV_USER:
92: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
93: }
94: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
95: if (svd->nini) PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(svd->nini));
96: if (svd->ninil) PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %" PetscInt_FMT "\n",PetscAbs(svd->ninil));
97: } else PetscTryTypeMethod(svd,view,viewer);
98: PetscObjectTypeCompareAny((PetscObject)svd,&isshell,SVDCROSS,SVDCYCLIC,SVDSCALAPACK,SVDELEMENTAL,SVDPRIMME,"");
99: if (!isshell) {
100: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
101: if (!svd->V) SVDGetBV(svd,&svd->V,NULL);
102: BVView(svd->V,viewer);
103: if (!svd->ds) SVDGetDS(svd,&svd->ds);
104: DSView(svd->ds,viewer);
105: PetscViewerPopFormat(viewer);
106: }
107: return 0;
108: }
110: /*@C
111: SVDViewFromOptions - View from options
113: Collective on SVD
115: Input Parameters:
116: + svd - the singular value solver context
117: . obj - optional object
118: - name - command line option
120: Level: intermediate
122: .seealso: SVDView(), SVDCreate()
123: @*/
124: PetscErrorCode SVDViewFromOptions(SVD svd,PetscObject obj,const char name[])
125: {
127: PetscObjectViewFromOptions((PetscObject)svd,obj,name);
128: return 0;
129: }
131: /*@C
132: SVDConvergedReasonView - Displays the reason an SVD solve converged or diverged.
134: Collective on svd
136: Input Parameters:
137: + svd - the singular value solver context
138: - viewer - the viewer to display the reason
140: Options Database Keys:
141: . -svd_converged_reason - print reason for convergence, and number of iterations
143: Note:
144: To change the format of the output call PetscViewerPushFormat(viewer,format) before
145: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
146: display a reason if it fails. The latter can be set in the command line with
147: -svd_converged_reason ::failed
149: Level: intermediate
151: .seealso: SVDSetTolerances(), SVDGetIterationNumber(), SVDConvergedReasonViewFromOptions()
152: @*/
153: PetscErrorCode SVDConvergedReasonView(SVD svd,PetscViewer viewer)
154: {
155: PetscBool isAscii;
156: PetscViewerFormat format;
158: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)svd));
159: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
160: if (isAscii) {
161: PetscViewerGetFormat(viewer,&format);
162: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
163: if (svd->reason > 0 && format != PETSC_VIEWER_FAILED) PetscViewerASCIIPrintf(viewer,"%s SVD solve converged (%" PetscInt_FMT " singular triplet%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)svd)->prefix?((PetscObject)svd)->prefix:"",svd->nconv,(svd->nconv>1)?"s":"",SVDConvergedReasons[svd->reason],svd->its);
164: else if (svd->reason <= 0) PetscViewerASCIIPrintf(viewer,"%s SVD solve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)svd)->prefix?((PetscObject)svd)->prefix:"",SVDConvergedReasons[svd->reason],svd->its);
165: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
166: }
167: return 0;
168: }
170: /*@
171: SVDConvergedReasonViewFromOptions - Processes command line options to determine if/how
172: the SVD converged reason is to be viewed.
174: Collective on svd
176: Input Parameter:
177: . svd - the singular value solver context
179: Level: developer
181: .seealso: SVDConvergedReasonView()
182: @*/
183: PetscErrorCode SVDConvergedReasonViewFromOptions(SVD svd)
184: {
185: PetscViewer viewer;
186: PetscBool flg;
187: static PetscBool incall = PETSC_FALSE;
188: PetscViewerFormat format;
190: if (incall) return 0;
191: incall = PETSC_TRUE;
192: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_converged_reason",&viewer,&format,&flg);
193: if (flg) {
194: PetscViewerPushFormat(viewer,format);
195: SVDConvergedReasonView(svd,viewer);
196: PetscViewerPopFormat(viewer);
197: PetscViewerDestroy(&viewer);
198: }
199: incall = PETSC_FALSE;
200: return 0;
201: }
203: static PetscErrorCode SVDErrorView_ASCII(SVD svd,SVDErrorType etype,PetscViewer viewer)
204: {
205: PetscReal error,sigma;
206: PetscInt i,j;
208: if (svd->nconv<svd->nsv) {
209: PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " singular values converged\n\n",svd->nsv);
210: return 0;
211: }
212: for (i=0;i<svd->nsv;i++) {
213: SVDComputeError(svd,i,etype,&error);
214: if (error>=5.0*svd->tol) {
215: PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",svd->nsv);
216: return 0;
217: }
218: }
219: PetscViewerASCIIPrintf(viewer," All requested %ssingular values computed up to the required tolerance:",svd->isgeneralized?"generalized ":"");
220: for (i=0;i<=(svd->nsv-1)/8;i++) {
221: PetscViewerASCIIPrintf(viewer,"\n ");
222: for (j=0;j<PetscMin(8,svd->nsv-8*i);j++) {
223: SVDGetSingularTriplet(svd,8*i+j,&sigma,NULL,NULL);
224: PetscViewerASCIIPrintf(viewer,"%.5f",(double)sigma);
225: if (8*i+j+1<svd->nsv) PetscViewerASCIIPrintf(viewer,", ");
226: }
227: }
228: PetscViewerASCIIPrintf(viewer,"\n\n");
229: return 0;
230: }
232: static PetscErrorCode SVDErrorView_DETAIL(SVD svd,SVDErrorType etype,PetscViewer viewer)
233: {
234: PetscReal error,sigma;
235: PetscInt i;
236: char ex[30],sep[]=" ---------------------- --------------------\n";
238: if (!svd->nconv) return 0;
239: switch (etype) {
240: case SVD_ERROR_ABSOLUTE:
241: PetscSNPrintf(ex,sizeof(ex)," absolute error");
242: break;
243: case SVD_ERROR_RELATIVE:
244: PetscSNPrintf(ex,sizeof(ex)," relative error");
245: break;
246: case SVD_ERROR_NORM:
247: if (svd->isgeneralized) PetscSNPrintf(ex,sizeof(ex)," ||r||/||[A;B]||");
248: else PetscSNPrintf(ex,sizeof(ex)," ||r||/||A||");
249: break;
250: }
251: PetscViewerASCIIPrintf(viewer,"%s sigma %s\n%s",sep,ex,sep);
252: for (i=0;i<svd->nconv;i++) {
253: SVDGetSingularTriplet(svd,i,&sigma,NULL,NULL);
254: SVDComputeError(svd,i,etype,&error);
255: PetscViewerASCIIPrintf(viewer," % 6f %12g\n",(double)sigma,(double)error);
256: }
257: PetscViewerASCIIPrintf(viewer,"%s",sep);
258: return 0;
259: }
261: static PetscErrorCode SVDErrorView_MATLAB(SVD svd,SVDErrorType etype,PetscViewer viewer)
262: {
263: PetscReal error;
264: PetscInt i;
265: const char *name;
267: PetscObjectGetName((PetscObject)svd,&name);
268: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
269: for (i=0;i<svd->nconv;i++) {
270: SVDComputeError(svd,i,etype,&error);
271: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
272: }
273: PetscViewerASCIIPrintf(viewer,"];\n");
274: return 0;
275: }
277: /*@C
278: SVDErrorView - Displays the errors associated with the computed solution
279: (as well as the singular values).
281: Collective on svd
283: Input Parameters:
284: + svd - the singular value solver context
285: . etype - error type
286: - viewer - optional visualization context
288: Options Database Keys:
289: + -svd_error_absolute - print absolute errors of each singular triplet
290: . -svd_error_relative - print relative errors of each singular triplet
291: - -svd_error_norm - print errors relative to the matrix norms of each singular triplet
293: Notes:
294: By default, this function checks the error of all singular triplets and prints
295: the singular values if all of them are below the requested tolerance.
296: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
297: singular values and corresponding errors is printed.
299: Level: intermediate
301: .seealso: SVDSolve(), SVDValuesView(), SVDVectorsView()
302: @*/
303: PetscErrorCode SVDErrorView(SVD svd,SVDErrorType etype,PetscViewer viewer)
304: {
305: PetscBool isascii;
306: PetscViewerFormat format;
309: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
312: SVDCheckSolved(svd,1);
313: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
314: if (!isascii) return 0;
316: PetscViewerGetFormat(viewer,&format);
317: switch (format) {
318: case PETSC_VIEWER_DEFAULT:
319: case PETSC_VIEWER_ASCII_INFO:
320: SVDErrorView_ASCII(svd,etype,viewer);
321: break;
322: case PETSC_VIEWER_ASCII_INFO_DETAIL:
323: SVDErrorView_DETAIL(svd,etype,viewer);
324: break;
325: case PETSC_VIEWER_ASCII_MATLAB:
326: SVDErrorView_MATLAB(svd,etype,viewer);
327: break;
328: default:
329: PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
330: }
331: return 0;
332: }
334: /*@
335: SVDErrorViewFromOptions - Processes command line options to determine if/how
336: the errors of the computed solution are to be viewed.
338: Collective on svd
340: Input Parameter:
341: . svd - the singular value solver context
343: Level: developer
345: .seealso: SVDErrorView()
346: @*/
347: PetscErrorCode SVDErrorViewFromOptions(SVD svd)
348: {
349: PetscViewer viewer;
350: PetscBool flg;
351: static PetscBool incall = PETSC_FALSE;
352: PetscViewerFormat format;
354: if (incall) return 0;
355: incall = PETSC_TRUE;
356: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_absolute",&viewer,&format,&flg);
357: if (flg) {
358: PetscViewerPushFormat(viewer,format);
359: SVDErrorView(svd,SVD_ERROR_ABSOLUTE,viewer);
360: PetscViewerPopFormat(viewer);
361: PetscViewerDestroy(&viewer);
362: }
363: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_relative",&viewer,&format,&flg);
364: if (flg) {
365: PetscViewerPushFormat(viewer,format);
366: SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer);
367: PetscViewerPopFormat(viewer);
368: PetscViewerDestroy(&viewer);
369: }
370: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_norm",&viewer,&format,&flg);
371: if (flg) {
372: PetscViewerPushFormat(viewer,format);
373: SVDErrorView(svd,SVD_ERROR_NORM,viewer);
374: PetscViewerPopFormat(viewer);
375: PetscViewerDestroy(&viewer);
376: }
377: incall = PETSC_FALSE;
378: return 0;
379: }
381: static PetscErrorCode SVDValuesView_DRAW(SVD svd,PetscViewer viewer)
382: {
383: PetscDraw draw;
384: PetscDrawSP drawsp;
385: PetscReal re,im=0.0;
386: PetscInt i;
388: if (!svd->nconv) return 0;
389: PetscViewerDrawGetDraw(viewer,0,&draw);
390: PetscDrawSetTitle(draw,"Computed singular values");
391: PetscDrawSPCreate(draw,1,&drawsp);
392: for (i=0;i<svd->nconv;i++) {
393: re = svd->sigma[svd->perm[i]];
394: PetscDrawSPAddPoint(drawsp,&re,&im);
395: }
396: PetscDrawSPDraw(drawsp,PETSC_TRUE);
397: PetscDrawSPSave(drawsp);
398: PetscDrawSPDestroy(&drawsp);
399: return 0;
400: }
402: static PetscErrorCode SVDValuesView_BINARY(SVD svd,PetscViewer viewer)
403: {
404: PetscInt i,k;
405: PetscReal *sv;
407: PetscMalloc1(svd->nconv,&sv);
408: for (i=0;i<svd->nconv;i++) {
409: k = svd->perm[i];
410: sv[i] = svd->sigma[k];
411: }
412: PetscViewerBinaryWrite(viewer,sv,svd->nconv,PETSC_REAL);
413: PetscFree(sv);
414: return 0;
415: }
417: #if defined(PETSC_HAVE_HDF5)
418: static PetscErrorCode SVDValuesView_HDF5(SVD svd,PetscViewer viewer)
419: {
420: PetscInt i,k,n,N;
421: PetscMPIInt rank;
422: Vec v;
423: char vname[30];
424: const char *ename;
426: MPI_Comm_rank(PetscObjectComm((PetscObject)svd),&rank);
427: N = svd->nconv;
428: n = rank? 0: N;
429: /* create a vector containing the singular values */
430: VecCreateMPI(PetscObjectComm((PetscObject)svd),n,N,&v);
431: PetscObjectGetName((PetscObject)svd,&ename);
432: PetscSNPrintf(vname,sizeof(vname),"sigma_%s",ename);
433: PetscObjectSetName((PetscObject)v,vname);
434: if (!rank) {
435: for (i=0;i<svd->nconv;i++) {
436: k = svd->perm[i];
437: VecSetValue(v,i,svd->sigma[k],INSERT_VALUES);
438: }
439: }
440: VecAssemblyBegin(v);
441: VecAssemblyEnd(v);
442: VecView(v,viewer);
443: VecDestroy(&v);
444: return 0;
445: }
446: #endif
448: static PetscErrorCode SVDValuesView_ASCII(SVD svd,PetscViewer viewer)
449: {
450: PetscInt i;
452: PetscViewerASCIIPrintf(viewer,"Singular values = \n");
453: for (i=0;i<svd->nconv;i++) PetscViewerASCIIPrintf(viewer," %.5f\n",(double)svd->sigma[svd->perm[i]]);
454: PetscViewerASCIIPrintf(viewer,"\n");
455: return 0;
456: }
458: static PetscErrorCode SVDValuesView_MATLAB(SVD svd,PetscViewer viewer)
459: {
460: PetscInt i;
461: const char *name;
463: PetscObjectGetName((PetscObject)svd,&name);
464: PetscViewerASCIIPrintf(viewer,"Sigma_%s = [\n",name);
465: for (i=0;i<svd->nconv;i++) PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)svd->sigma[svd->perm[i]]);
466: PetscViewerASCIIPrintf(viewer,"];\n");
467: return 0;
468: }
470: /*@C
471: SVDValuesView - Displays the computed singular values in a viewer.
473: Collective on svd
475: Input Parameters:
476: + svd - the singular value solver context
477: - viewer - the viewer
479: Options Database Key:
480: . -svd_view_values - print computed singular values
482: Level: intermediate
484: .seealso: SVDSolve(), SVDVectorsView(), SVDErrorView()
485: @*/
486: PetscErrorCode SVDValuesView(SVD svd,PetscViewer viewer)
487: {
488: PetscBool isascii,isdraw,isbinary;
489: PetscViewerFormat format;
490: #if defined(PETSC_HAVE_HDF5)
491: PetscBool ishdf5;
492: #endif
495: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
498: SVDCheckSolved(svd,1);
499: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
500: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
501: #if defined(PETSC_HAVE_HDF5)
502: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
503: #endif
504: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
505: if (isdraw) SVDValuesView_DRAW(svd,viewer);
506: else if (isbinary) SVDValuesView_BINARY(svd,viewer);
507: #if defined(PETSC_HAVE_HDF5)
508: else if (ishdf5) SVDValuesView_HDF5(svd,viewer);
509: #endif
510: else if (isascii) {
511: PetscViewerGetFormat(viewer,&format);
512: switch (format) {
513: case PETSC_VIEWER_DEFAULT:
514: case PETSC_VIEWER_ASCII_INFO:
515: case PETSC_VIEWER_ASCII_INFO_DETAIL:
516: SVDValuesView_ASCII(svd,viewer);
517: break;
518: case PETSC_VIEWER_ASCII_MATLAB:
519: SVDValuesView_MATLAB(svd,viewer);
520: break;
521: default:
522: PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
523: }
524: }
525: return 0;
526: }
528: /*@
529: SVDValuesViewFromOptions - Processes command line options to determine if/how
530: the computed singular values are to be viewed.
532: Collective on svd
534: Input Parameter:
535: . svd - the singular value solver context
537: Level: developer
539: .seealso: SVDValuesView()
540: @*/
541: PetscErrorCode SVDValuesViewFromOptions(SVD svd)
542: {
543: PetscViewer viewer;
544: PetscBool flg;
545: static PetscBool incall = PETSC_FALSE;
546: PetscViewerFormat format;
548: if (incall) return 0;
549: incall = PETSC_TRUE;
550: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_values",&viewer,&format,&flg);
551: if (flg) {
552: PetscViewerPushFormat(viewer,format);
553: SVDValuesView(svd,viewer);
554: PetscViewerPopFormat(viewer);
555: PetscViewerDestroy(&viewer);
556: }
557: incall = PETSC_FALSE;
558: return 0;
559: }
561: /*@C
562: SVDVectorsView - Outputs computed singular vectors to a viewer.
564: Collective on svd
566: Input Parameters:
567: + svd - the singular value solver context
568: - viewer - the viewer
570: Options Database Key:
571: . -svd_view_vectors - output singular vectors
573: Note:
574: Right and left singular vectors are interleaved, that is, the vectors are
575: output in the following order V0, U0, V1, U1, V2, U2, ...
577: Level: intermediate
579: .seealso: SVDSolve(), SVDValuesView(), SVDErrorView()
580: @*/
581: PetscErrorCode SVDVectorsView(SVD svd,PetscViewer viewer)
582: {
583: PetscInt i,k;
584: Vec x;
585: char vname[30];
586: const char *ename;
589: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
592: SVDCheckSolved(svd,1);
593: if (svd->nconv) {
594: PetscObjectGetName((PetscObject)svd,&ename);
595: SVDComputeVectors(svd);
596: for (i=0;i<svd->nconv;i++) {
597: k = svd->perm[i];
598: PetscSNPrintf(vname,sizeof(vname),"V%" PetscInt_FMT "_%s",i,ename);
599: BVGetColumn(svd->V,k,&x);
600: PetscObjectSetName((PetscObject)x,vname);
601: VecView(x,viewer);
602: BVRestoreColumn(svd->V,k,&x);
603: PetscSNPrintf(vname,sizeof(vname),"U%" PetscInt_FMT "_%s",i,ename);
604: BVGetColumn(svd->U,k,&x);
605: PetscObjectSetName((PetscObject)x,vname);
606: VecView(x,viewer);
607: BVRestoreColumn(svd->U,k,&x);
608: }
609: }
610: return 0;
611: }
613: /*@
614: SVDVectorsViewFromOptions - Processes command line options to determine if/how
615: the computed singular vectors are to be viewed.
617: Collective on svd
619: Input Parameter:
620: . svd - the singular value solver context
622: Level: developer
624: .seealso: SVDVectorsView()
625: @*/
626: PetscErrorCode SVDVectorsViewFromOptions(SVD svd)
627: {
628: PetscViewer viewer;
629: PetscBool flg = PETSC_FALSE;
630: static PetscBool incall = PETSC_FALSE;
631: PetscViewerFormat format;
633: if (incall) return 0;
634: incall = PETSC_TRUE;
635: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_vectors",&viewer,&format,&flg);
636: if (flg) {
637: PetscViewerPushFormat(viewer,format);
638: SVDVectorsView(svd,viewer);
639: PetscViewerPopFormat(viewer);
640: PetscViewerDestroy(&viewer);
641: }
642: incall = PETSC_FALSE;
643: return 0;
644: }