Actual source code: svdview.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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> /*I "slepcsvd.h" I*/
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: STView(), PetscViewerASCIIOpen()
43: @*/
44: PetscErrorCode SVDView(SVD svd,PetscViewer viewer)
45: {
47: PetscBool isascii,isshell;
51: if (!viewer) {
52: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
53: }
57: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
58: if (isascii) {
59: PetscObjectPrintClassNamePrefixType((PetscObject)svd,viewer);
60: if (svd->ops->view) {
61: PetscViewerASCIIPushTab(viewer);
62: (*svd->ops->view)(svd,viewer);
63: PetscViewerASCIIPopTab(viewer);
64: }
65: PetscViewerASCIIPrintf(viewer," transpose mode: %s\n",svd->impltrans?"implicit":"explicit");
66: if (svd->which == SVD_LARGEST) {
67: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: largest\n");
68: } else {
69: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: smallest\n");
70: }
71: PetscViewerASCIIPrintf(viewer," number of singular values (nsv): %D\n",svd->nsv);
72: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",svd->ncv);
73: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",svd->mpd);
74: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\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_USER:
84: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
85: }
86: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
87: if (svd->nini) {
88: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(svd->nini));
89: }
90: if (svd->ninil) {
91: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %D\n",PetscAbs(svd->ninil));
92: }
93: } else {
94: if (svd->ops->view) {
95: (*svd->ops->view)(svd,viewer);
96: }
97: }
98: PetscObjectTypeCompareAny((PetscObject)svd,&isshell,SVDCROSS,SVDCYCLIC,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: SVDReasonView - Displays the reason an SVD solve converged or diverged.
113: Collective on SVD
115: Parameter:
116: + svd - the singular value solver context
117: - viewer - the viewer to display the reason
119: Options Database Keys:
120: . -svd_converged_reason - print reason for convergence, and number of iterations
122: Level: intermediate
124: .seealso: SVDSetTolerances(), SVDGetIterationNumber()
125: @*/
126: PetscErrorCode SVDReasonView(SVD svd,PetscViewer viewer)
127: {
129: PetscBool isAscii;
132: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
133: if (isAscii) {
134: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
135: if (svd->reason > 0) {
136: 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);
137: } else {
138: 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);
139: }
140: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
141: }
142: return(0);
143: }
145: /*@
146: SVDReasonViewFromOptions - Processes command line options to determine if/how
147: the SVD converged reason is to be viewed.
149: Collective on SVD
151: Input Parameters:
152: . svd - the singular value solver context
154: Level: developer
155: @*/
156: PetscErrorCode SVDReasonViewFromOptions(SVD svd)
157: {
158: PetscErrorCode ierr;
159: PetscViewer viewer;
160: PetscBool flg;
161: static PetscBool incall = PETSC_FALSE;
162: PetscViewerFormat format;
165: if (incall) return(0);
166: incall = PETSC_TRUE;
167: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_converged_reason",&viewer,&format,&flg);
168: if (flg) {
169: PetscViewerPushFormat(viewer,format);
170: SVDReasonView(svd,viewer);
171: PetscViewerPopFormat(viewer);
172: PetscViewerDestroy(&viewer);
173: }
174: incall = PETSC_FALSE;
175: return(0);
176: }
178: static PetscErrorCode SVDErrorView_ASCII(SVD svd,SVDErrorType etype,PetscViewer viewer)
179: {
180: PetscReal error,sigma;
181: PetscInt i,j;
185: if (svd->nconv<svd->nsv) {
186: PetscViewerASCIIPrintf(viewer," Problem: less than %D singular values converged\n\n",svd->nsv);
187: return(0);
188: }
189: for (i=0;i<svd->nsv;i++) {
190: SVDComputeError(svd,i,etype,&error);
191: if (error>=5.0*svd->tol) {
192: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",svd->nsv);
193: return(0);
194: }
195: }
196: PetscViewerASCIIPrintf(viewer," All requested singular values computed up to the required tolerance:");
197: for (i=0;i<=(svd->nsv-1)/8;i++) {
198: PetscViewerASCIIPrintf(viewer,"\n ");
199: for (j=0;j<PetscMin(8,svd->nsv-8*i);j++) {
200: SVDGetSingularTriplet(svd,8*i+j,&sigma,NULL,NULL);
201: PetscViewerASCIIPrintf(viewer,"%.5f",(double)sigma);
202: if (8*i+j+1<svd->nsv) { PetscViewerASCIIPrintf(viewer,", "); }
203: }
204: }
205: PetscViewerASCIIPrintf(viewer,"\n\n");
206: return(0);
207: }
209: static PetscErrorCode SVDErrorView_DETAIL(SVD svd,SVDErrorType etype,PetscViewer viewer)
210: {
212: PetscReal error,sigma;
213: PetscInt i;
214: #define EXLEN 30
215: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
218: if (!svd->nconv) return(0);
219: switch (etype) {
220: case SVD_ERROR_ABSOLUTE:
221: PetscSNPrintf(ex,EXLEN," absolute error");
222: break;
223: case SVD_ERROR_RELATIVE:
224: PetscSNPrintf(ex,EXLEN," relative error");
225: break;
226: }
227: PetscViewerASCIIPrintf(viewer,"%s sigma %s\n%s",sep,ex,sep);
228: for (i=0;i<svd->nconv;i++) {
229: SVDGetSingularTriplet(svd,i,&sigma,NULL,NULL);
230: SVDComputeError(svd,i,etype,&error);
231: PetscViewerASCIIPrintf(viewer," % 6f %12g\n",(double)sigma,(double)error);
232: }
233: PetscViewerASCIIPrintf(viewer,sep);
234: return(0);
235: }
237: static PetscErrorCode SVDErrorView_MATLAB(SVD svd,SVDErrorType etype,PetscViewer viewer)
238: {
240: PetscReal error;
241: PetscInt i;
242: const char *name;
245: PetscObjectGetName((PetscObject)svd,&name);
246: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
247: for (i=0;i<svd->nconv;i++) {
248: SVDComputeError(svd,i,etype,&error);
249: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
250: }
251: PetscViewerASCIIPrintf(viewer,"];\n");
252: return(0);
253: }
255: /*@C
256: SVDErrorView - Displays the errors associated with the computed solution
257: (as well as the singular values).
259: Collective on SVD
261: Input Parameters:
262: + svd - the singular value solver context
263: . etype - error type
264: - viewer - optional visualization context
266: Options Database Key:
267: + -svd_error_absolute - print absolute errors of each singular triplet
268: - -svd_error_relative - print relative errors of each singular triplet
270: Notes:
271: By default, this function checks the error of all singular triplets and prints
272: the singular values if all of them are below the requested tolerance.
273: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
274: singular values and corresponding errors is printed.
276: Level: intermediate
278: .seealso: SVDSolve(), SVDValuesView(), SVDVectorsView()
279: @*/
280: PetscErrorCode SVDErrorView(SVD svd,SVDErrorType etype,PetscViewer viewer)
281: {
282: PetscBool isascii;
283: PetscViewerFormat format;
284: PetscErrorCode ierr;
288: if (!viewer) {
289: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
290: }
293: SVDCheckSolved(svd,1);
294: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
295: if (!isascii) return(0);
297: PetscViewerGetFormat(viewer,&format);
298: switch (format) {
299: case PETSC_VIEWER_DEFAULT:
300: case PETSC_VIEWER_ASCII_INFO:
301: SVDErrorView_ASCII(svd,etype,viewer);
302: break;
303: case PETSC_VIEWER_ASCII_INFO_DETAIL:
304: SVDErrorView_DETAIL(svd,etype,viewer);
305: break;
306: case PETSC_VIEWER_ASCII_MATLAB:
307: SVDErrorView_MATLAB(svd,etype,viewer);
308: break;
309: default:
310: PetscInfo1(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
311: }
312: return(0);
313: }
315: /*@
316: SVDErrorViewFromOptions - Processes command line options to determine if/how
317: the errors of the computed solution are to be viewed.
319: Collective on SVD
321: Input Parameters:
322: . svd - the singular value solver context
324: Level: developer
325: @*/
326: PetscErrorCode SVDErrorViewFromOptions(SVD svd)
327: {
328: PetscErrorCode ierr;
329: PetscViewer viewer;
330: PetscBool flg;
331: static PetscBool incall = PETSC_FALSE;
332: PetscViewerFormat format;
335: if (incall) return(0);
336: incall = PETSC_TRUE;
337: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_absolute",&viewer,&format,&flg);
338: if (flg) {
339: PetscViewerPushFormat(viewer,format);
340: SVDErrorView(svd,SVD_ERROR_ABSOLUTE,viewer);
341: PetscViewerPopFormat(viewer);
342: PetscViewerDestroy(&viewer);
343: }
344: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_relative",&viewer,&format,&flg);
345: if (flg) {
346: PetscViewerPushFormat(viewer,format);
347: SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer);
348: PetscViewerPopFormat(viewer);
349: PetscViewerDestroy(&viewer);
350: }
351: incall = PETSC_FALSE;
352: return(0);
353: }
355: static PetscErrorCode SVDValuesView_DRAW(SVD svd,PetscViewer viewer)
356: {
358: PetscDraw draw;
359: PetscDrawSP drawsp;
360: PetscReal re,im=0.0;
361: PetscInt i;
364: if (!svd->nconv) return(0);
365: PetscViewerDrawGetDraw(viewer,0,&draw);
366: PetscDrawSetTitle(draw,"Computed singular values");
367: PetscDrawSPCreate(draw,1,&drawsp);
368: for (i=0;i<svd->nconv;i++) {
369: re = svd->sigma[svd->perm[i]];
370: PetscDrawSPAddPoint(drawsp,&re,&im);
371: }
372: PetscDrawSPDraw(drawsp,PETSC_TRUE);
373: PetscDrawSPSave(drawsp);
374: PetscDrawSPDestroy(&drawsp);
375: return(0);
376: }
378: static PetscErrorCode SVDValuesView_ASCII(SVD svd,PetscViewer viewer)
379: {
380: PetscInt i;
384: PetscViewerASCIIPrintf(viewer,"Singular values = \n");
385: for (i=0;i<svd->nconv;i++) {
386: PetscViewerASCIIPrintf(viewer," %.5f\n",(double)svd->sigma[svd->perm[i]]);
387: }
388: PetscViewerASCIIPrintf(viewer,"\n");
389: return(0);
390: }
392: static PetscErrorCode SVDValuesView_MATLAB(SVD svd,PetscViewer viewer)
393: {
395: PetscInt i;
396: const char *name;
399: PetscObjectGetName((PetscObject)svd,&name);
400: PetscViewerASCIIPrintf(viewer,"Sigma_%s = [\n",name);
401: for (i=0;i<svd->nconv;i++) {
402: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)svd->sigma[svd->perm[i]]);
403: }
404: PetscViewerASCIIPrintf(viewer,"];\n");
405: return(0);
406: }
408: /*@C
409: SVDValuesView - Displays the computed singular values in a viewer.
411: Collective on SVD
413: Input Parameters:
414: + svd - the singular value solver context
415: - viewer - the viewer
417: Options Database Key:
418: . -svd_view_values - print computed singular values
420: Level: intermediate
422: .seealso: SVDSolve(), SVDVectorsView(), SVDErrorView()
423: @*/
424: PetscErrorCode SVDValuesView(SVD svd,PetscViewer viewer)
425: {
426: PetscBool isascii,isdraw;
427: PetscViewerFormat format;
428: PetscErrorCode ierr;
432: if (!viewer) {
433: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
434: }
437: SVDCheckSolved(svd,1);
438: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
439: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
440: if (isdraw) {
441: SVDValuesView_DRAW(svd,viewer);
442: } else if (isascii) {
443: PetscViewerGetFormat(viewer,&format);
444: switch (format) {
445: case PETSC_VIEWER_DEFAULT:
446: case PETSC_VIEWER_ASCII_INFO:
447: case PETSC_VIEWER_ASCII_INFO_DETAIL:
448: SVDValuesView_ASCII(svd,viewer);
449: break;
450: case PETSC_VIEWER_ASCII_MATLAB:
451: SVDValuesView_MATLAB(svd,viewer);
452: break;
453: default:
454: PetscInfo1(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
455: }
456: }
457: return(0);
458: }
460: /*@
461: SVDValuesViewFromOptions - Processes command line options to determine if/how
462: the computed singular values are to be viewed.
464: Collective on SVD
466: Input Parameters:
467: . svd - the singular value solver context
469: Level: developer
470: @*/
471: PetscErrorCode SVDValuesViewFromOptions(SVD svd)
472: {
473: PetscErrorCode ierr;
474: PetscViewer viewer;
475: PetscBool flg;
476: static PetscBool incall = PETSC_FALSE;
477: PetscViewerFormat format;
480: if (incall) return(0);
481: incall = PETSC_TRUE;
482: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_values",&viewer,&format,&flg);
483: if (flg) {
484: PetscViewerPushFormat(viewer,format);
485: SVDValuesView(svd,viewer);
486: PetscViewerPopFormat(viewer);
487: PetscViewerDestroy(&viewer);
488: }
489: incall = PETSC_FALSE;
490: return(0);
491: }
493: /*@C
494: SVDVectorsView - Outputs computed singular vectors to a viewer.
496: Collective on SVD
498: Parameter:
499: + svd - the singular value solver context
500: - viewer - the viewer
502: Options Database Keys:
503: . -svd_view_vectors - output singular vectors
505: Note:
506: Right and left singular vectors are interleaved, that is, the vectors are
507: output in the following order: V0, U0, V1, U1, V2, U2, ...
509: Level: intermediate
511: .seealso: SVDSolve(), SVDValuesView(), SVDErrorView()
512: @*/
513: PetscErrorCode SVDVectorsView(SVD svd,PetscViewer viewer)
514: {
516: PetscInt i,k;
517: Vec x;
518: #define NMLEN 30
519: char vname[NMLEN];
520: const char *ename;
524: if (!viewer) {
525: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer);
526: }
529: SVDCheckSolved(svd,1);
530: if (svd->nconv) {
531: PetscObjectGetName((PetscObject)svd,&ename);
532: SVDComputeVectors(svd);
533: for (i=0;i<svd->nconv;i++) {
534: k = svd->perm[i];
535: PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
536: BVGetColumn(svd->V,k,&x);
537: PetscObjectSetName((PetscObject)x,vname);
538: VecView(x,viewer);
539: BVRestoreColumn(svd->V,k,&x);
540: PetscSNPrintf(vname,NMLEN,"U%d_%s",(int)i,ename);
541: BVGetColumn(svd->U,k,&x);
542: PetscObjectSetName((PetscObject)x,vname);
543: VecView(x,viewer);
544: BVRestoreColumn(svd->U,k,&x);
545: }
546: }
547: return(0);
548: }
550: /*@
551: SVDVectorsViewFromOptions - Processes command line options to determine if/how
552: the computed singular vectors are to be viewed.
554: Collective on SVD
556: Input Parameters:
557: . svd - the singular value solver context
559: Level: developer
560: @*/
561: PetscErrorCode SVDVectorsViewFromOptions(SVD svd)
562: {
563: PetscErrorCode ierr;
564: PetscViewer viewer;
565: PetscBool flg = PETSC_FALSE;
566: static PetscBool incall = PETSC_FALSE;
567: PetscViewerFormat format;
570: if (incall) return(0);
571: incall = PETSC_TRUE;
572: PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_vectors",&viewer,&format,&flg);
573: if (flg) {
574: PetscViewerPushFormat(viewer,format);
575: SVDVectorsView(svd,viewer);
576: PetscViewerPopFormat(viewer);
577: PetscViewerDestroy(&viewer);
578: }
579: incall = PETSC_FALSE;
580: return(0);
581: }