Actual source code: svdview.c
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: /*@
18: SVDView - Prints the `SVD` data structure.
20: Collective
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: Notes:
30: The available visualization contexts include
31: + `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
32: - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
33: first process opens the file; all other processes send their data to the
34: first one to print
36: The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
37: to output to a specified file.
39: Level: beginner
41: .seealso: [](ch:svd), `SVDCreate()`, `SVDViewFromOptions()`, `EPSView()`
42: @*/
43: PetscErrorCode SVDView(SVD svd,PetscViewer viewer)
44: {
45: const char *type=NULL;
46: PetscBool isascii,isshell,isexternal;
48: PetscFunctionBegin;
50: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
52: PetscCheckSameComm(svd,1,viewer,2);
54: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
55: if (isascii) {
56: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)svd,viewer));
57: PetscCall(PetscViewerASCIIPushTab(viewer));
58: PetscTryTypeMethod(svd,view,viewer);
59: PetscCall(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: PetscCall(PetscViewerASCIIPrintf(viewer," problem type: %s\n",type));
68: PetscCall(PetscViewerASCIIPrintf(viewer," transpose mode: %s\n",svd->impltrans?"implicit":"explicit"));
69: if (svd->which == SVD_LARGEST) PetscCall(PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: largest\n"));
70: else PetscCall(PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: smallest\n"));
71: if (svd->stop==SVD_STOP_THRESHOLD) PetscCall(PetscViewerASCIIPrintf(viewer," computing singular values %s the threshold: %g%s\n",svd->which==SVD_LARGEST?"above":"below",(double)svd->thres,svd->threlative?" (relative)":""));
72: if (svd->nsv) PetscCall(PetscViewerASCIIPrintf(viewer," number of singular values (nsv): %" PetscInt_FMT "\n",svd->nsv));
73: PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",svd->ncv));
74: PetscCall(PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",svd->mpd));
75: PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",svd->max_it));
76: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)svd->tol));
77: PetscCall(PetscViewerASCIIPrintf(viewer," convergence test: "));
78: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
79: switch (svd->conv) {
80: case SVD_CONV_ABS:
81: PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
82: case SVD_CONV_REL:
83: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the singular value\n"));break;
84: case SVD_CONV_NORM:
85: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n"));
86: PetscCall(PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)svd->nrma));
87: if (svd->isgeneralized) PetscCall(PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)svd->nrmb));
88: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
89: break;
90: case SVD_CONV_MAXIT:
91: PetscCall(PetscViewerASCIIPrintf(viewer,"maximum number of iterations\n"));break;
92: case SVD_CONV_USER:
93: PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
94: }
95: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
96: if (svd->nini) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(svd->nini)));
97: if (svd->ninil) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %" PetscInt_FMT "\n",PetscAbs(svd->ninil)));
98: } else PetscTryTypeMethod(svd,view,viewer);
99: PetscCall(PetscObjectTypeCompareAny((PetscObject)svd,&isshell,SVDCROSS,SVDCYCLIC,""));
100: PetscCall(PetscObjectTypeCompareAny((PetscObject)svd,&isexternal,SVDSCALAPACK,SVDKSVD,SVDELEMENTAL,SVDPRIMME,""));
101: if (!isshell && !isexternal) {
102: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
103: if (!svd->V) PetscCall(SVDGetBV(svd,&svd->V,NULL));
104: PetscCall(BVView(svd->V,viewer));
105: if (!svd->ds) PetscCall(SVDGetDS(svd,&svd->ds));
106: PetscCall(DSView(svd->ds,viewer));
107: PetscCall(PetscViewerPopFormat(viewer));
108: }
109: PetscFunctionReturn(PETSC_SUCCESS);
110: }
112: /*@
113: SVDViewFromOptions - View (print) an `SVD` object based on values in the options database.
115: Collective
117: Input Parameters:
118: + svd - the singular value solver context
119: . obj - optional object that provides the options prefix used to query the options database
120: - name - command line option
122: Level: intermediate
124: .seealso: [](ch:svd), `SVDView()`, `SVDCreate()`, `PetscObjectViewFromOptions()`
125: @*/
126: PetscErrorCode SVDViewFromOptions(SVD svd,PetscObject obj,const char name[])
127: {
128: PetscFunctionBegin;
130: PetscCall(PetscObjectViewFromOptions((PetscObject)svd,obj,name));
131: PetscFunctionReturn(PETSC_SUCCESS);
132: }
134: /*@
135: SVDConvergedReasonView - Displays the reason an `SVD` solve converged or diverged.
137: Collective
139: Input Parameters:
140: + svd - the singular value solver context
141: - viewer - the viewer to display the reason
143: Options Database Key:
144: . -svd_converged_reason - print reason for convergence/divergence, and number of iterations
146: Notes:
147: Use `SVDConvergedReasonViewFromOptions()` to display the reason based on values
148: in the options database.
150: To change the format of the output call `PetscViewerPushFormat()` before this
151: call. Use `PETSC_VIEWER_DEFAULT` for the default, or `PETSC_VIEWER_FAILED` to only
152: display a reason if it fails. The latter can be set in the command line with
153: `-svd_converged_reason ::failed`.
155: Level: intermediate
157: .seealso: [](ch:svd), `SVDSetTolerances()`, `SVDGetIterationNumber()`, `SVDConvergedReasonViewFromOptions()`
158: @*/
159: PetscErrorCode SVDConvergedReasonView(SVD svd,PetscViewer viewer)
160: {
161: PetscBool isAscii;
162: PetscViewerFormat format;
164: PetscFunctionBegin;
165: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)svd));
166: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
167: if (isAscii) {
168: PetscCall(PetscViewerGetFormat(viewer,&format));
169: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
170: if (svd->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(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));
171: else if (svd->reason <= 0) PetscCall(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));
172: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
173: }
174: PetscFunctionReturn(PETSC_SUCCESS);
175: }
177: /*@
178: SVDConvergedReasonViewFromOptions - Processes command line options to determine if/how
179: the `SVD` converged reason is to be viewed.
181: Collective
183: Input Parameter:
184: . svd - the singular value solver context
186: Level: intermediate
188: .seealso: [](ch:svd), `SVDConvergedReasonView()`
189: @*/
190: PetscErrorCode SVDConvergedReasonViewFromOptions(SVD svd)
191: {
192: PetscViewer viewer;
193: PetscBool flg;
194: static PetscBool incall = PETSC_FALSE;
195: PetscViewerFormat format;
197: PetscFunctionBegin;
198: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
199: incall = PETSC_TRUE;
200: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_converged_reason",&viewer,&format,&flg));
201: if (flg) {
202: PetscCall(PetscViewerPushFormat(viewer,format));
203: PetscCall(SVDConvergedReasonView(svd,viewer));
204: PetscCall(PetscViewerPopFormat(viewer));
205: PetscCall(PetscViewerDestroy(&viewer));
206: }
207: incall = PETSC_FALSE;
208: PetscFunctionReturn(PETSC_SUCCESS);
209: }
211: static PetscErrorCode SVDErrorView_ASCII(SVD svd,SVDErrorType etype,PetscViewer viewer)
212: {
213: PetscReal error,sigma;
214: PetscInt i,j,numsv;
216: PetscFunctionBegin;
217: numsv = svd->stop==SVD_STOP_THRESHOLD? svd->nconv: svd->nsv;
218: if (svd->nconv<numsv) {
219: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " singular values converged\n\n",numsv));
220: PetscFunctionReturn(PETSC_SUCCESS);
221: }
222: for (i=0;i<numsv;i++) {
223: PetscCall(SVDComputeError(svd,i,etype,&error));
224: if (error>=5.0*svd->tol) {
225: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",numsv));
226: PetscFunctionReturn(PETSC_SUCCESS);
227: }
228: }
229: PetscCall(PetscViewerASCIIPrintf(viewer," All %s%ssingular values computed up to the required tolerance:",svd->stop==SVD_STOP_THRESHOLD?"":"requested ",svd->isgeneralized?"generalized ":""));
230: for (i=0;i<=(numsv-1)/8;i++) {
231: PetscCall(PetscViewerASCIIPrintf(viewer,"\n "));
232: for (j=0;j<PetscMin(8,numsv-8*i);j++) {
233: PetscCall(SVDGetSingularTriplet(svd,8*i+j,&sigma,NULL,NULL));
234: PetscCall(PetscViewerASCIIPrintf(viewer,"%.5f",(double)sigma));
235: if (8*i+j+1<numsv) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
236: }
237: }
238: PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
239: PetscFunctionReturn(PETSC_SUCCESS);
240: }
242: static PetscErrorCode SVDErrorView_DETAIL(SVD svd,SVDErrorType etype,PetscViewer viewer)
243: {
244: PetscReal error,sigma;
245: PetscInt i;
246: char ex[30],sep[]=" ---------------------- --------------------\n";
248: PetscFunctionBegin;
249: if (!svd->nconv) PetscFunctionReturn(PETSC_SUCCESS);
250: switch (etype) {
251: case SVD_ERROR_ABSOLUTE:
252: PetscCall(PetscSNPrintf(ex,sizeof(ex)," absolute error"));
253: break;
254: case SVD_ERROR_RELATIVE:
255: PetscCall(PetscSNPrintf(ex,sizeof(ex)," relative error"));
256: break;
257: case SVD_ERROR_NORM:
258: if (svd->isgeneralized) PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||r||/||[A;B]||"));
259: else PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||r||/||A||"));
260: break;
261: }
262: PetscCall(PetscViewerASCIIPrintf(viewer,"%s sigma %s\n%s",sep,ex,sep));
263: for (i=0;i<svd->nconv;i++) {
264: PetscCall(SVDGetSingularTriplet(svd,i,&sigma,NULL,NULL));
265: PetscCall(SVDComputeError(svd,i,etype,&error));
266: PetscCall(PetscViewerASCIIPrintf(viewer," % 6f %12g\n",(double)sigma,(double)error));
267: }
268: PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
269: PetscFunctionReturn(PETSC_SUCCESS);
270: }
272: static PetscErrorCode SVDErrorView_MATLAB(SVD svd,SVDErrorType etype,PetscViewer viewer)
273: {
274: PetscReal error;
275: PetscInt i;
276: const char *name;
278: PetscFunctionBegin;
279: PetscCall(PetscObjectGetName((PetscObject)svd,&name));
280: PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
281: for (i=0;i<svd->nconv;i++) {
282: PetscCall(SVDComputeError(svd,i,etype,&error));
283: PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
284: }
285: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
286: PetscFunctionReturn(PETSC_SUCCESS);
287: }
289: /*@
290: SVDErrorView - Displays the errors associated with the computed solution
291: (as well as the singular values).
293: Collective
295: Input Parameters:
296: + svd - the singular value solver context
297: . etype - error type
298: - viewer - optional visualization context
300: Options Database Keys:
301: + -svd_error_absolute - print absolute errors of each singular triplet
302: . -svd_error_relative - print relative errors of each singular triplet
303: - -svd_error_norm - print errors relative to the matrix norms of each singular triplet
305: Notes:
306: By default, this function checks the error of all singular triplets and prints
307: the singular values if all of them are below the requested tolerance.
308: If the viewer has format `PETSC_VIEWER_ASCII_INFO_DETAIL` then a table with
309: singular values and corresponding errors is printed.
311: All the command-line options listed above admit an optional argument
312: specifying the viewer type and options. For instance, use
313: `-svd_error_relative :myerr.m:ascii_matlab` to save the errors in a file
314: that can be executed in Matlab.
316: Level: intermediate
318: .seealso: [](ch:svd), `SVDSolve()`, `SVDValuesView()`, `SVDVectorsView()`
319: @*/
320: PetscErrorCode SVDErrorView(SVD svd,SVDErrorType etype,PetscViewer viewer)
321: {
322: PetscBool isascii;
323: PetscViewerFormat format;
325: PetscFunctionBegin;
327: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
329: PetscCheckSameComm(svd,1,viewer,3);
330: SVDCheckSolved(svd,1);
331: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
332: if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);
334: PetscCall(PetscViewerGetFormat(viewer,&format));
335: switch (format) {
336: case PETSC_VIEWER_DEFAULT:
337: case PETSC_VIEWER_ASCII_INFO:
338: PetscCall(SVDErrorView_ASCII(svd,etype,viewer));
339: break;
340: case PETSC_VIEWER_ASCII_INFO_DETAIL:
341: PetscCall(SVDErrorView_DETAIL(svd,etype,viewer));
342: break;
343: case PETSC_VIEWER_ASCII_MATLAB:
344: PetscCall(SVDErrorView_MATLAB(svd,etype,viewer));
345: break;
346: default:
347: PetscCall(PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
348: }
349: PetscFunctionReturn(PETSC_SUCCESS);
350: }
352: /*@
353: SVDErrorViewFromOptions - Processes command line options to determine if/how
354: the errors of the computed solution are to be viewed.
356: Collective
358: Input Parameter:
359: . svd - the singular value solver context
361: Level: developer
363: .seealso: [](ch:svd), `SVDErrorView()`
364: @*/
365: PetscErrorCode SVDErrorViewFromOptions(SVD svd)
366: {
367: PetscViewer viewer;
368: PetscBool flg;
369: static PetscBool incall = PETSC_FALSE;
370: PetscViewerFormat format;
372: PetscFunctionBegin;
373: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
374: incall = PETSC_TRUE;
375: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_absolute",&viewer,&format,&flg));
376: if (flg) {
377: PetscCall(PetscViewerPushFormat(viewer,format));
378: PetscCall(SVDErrorView(svd,SVD_ERROR_ABSOLUTE,viewer));
379: PetscCall(PetscViewerPopFormat(viewer));
380: PetscCall(PetscViewerDestroy(&viewer));
381: }
382: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_relative",&viewer,&format,&flg));
383: if (flg) {
384: PetscCall(PetscViewerPushFormat(viewer,format));
385: PetscCall(SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer));
386: PetscCall(PetscViewerPopFormat(viewer));
387: PetscCall(PetscViewerDestroy(&viewer));
388: }
389: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_norm",&viewer,&format,&flg));
390: if (flg) {
391: PetscCall(PetscViewerPushFormat(viewer,format));
392: PetscCall(SVDErrorView(svd,SVD_ERROR_NORM,viewer));
393: PetscCall(PetscViewerPopFormat(viewer));
394: PetscCall(PetscViewerDestroy(&viewer));
395: }
396: incall = PETSC_FALSE;
397: PetscFunctionReturn(PETSC_SUCCESS);
398: }
400: static PetscErrorCode SVDValuesView_DRAW(SVD svd,PetscViewer viewer)
401: {
402: PetscDraw draw;
403: PetscDrawSP drawsp;
404: PetscReal re,im=0.0;
405: PetscInt i;
407: PetscFunctionBegin;
408: if (!svd->nconv) PetscFunctionReturn(PETSC_SUCCESS);
409: PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
410: PetscCall(PetscDrawSetTitle(draw,"Computed singular values"));
411: PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
412: for (i=0;i<svd->nconv;i++) {
413: re = svd->sigma[svd->perm[i]];
414: PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
415: }
416: PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
417: PetscCall(PetscDrawSPSave(drawsp));
418: PetscCall(PetscDrawSPDestroy(&drawsp));
419: PetscFunctionReturn(PETSC_SUCCESS);
420: }
422: static PetscErrorCode SVDValuesView_BINARY(SVD svd,PetscViewer viewer)
423: {
424: PetscInt i,k;
425: PetscReal *sv;
427: PetscFunctionBegin;
428: PetscCall(PetscMalloc1(svd->nconv,&sv));
429: for (i=0;i<svd->nconv;i++) {
430: k = svd->perm[i];
431: sv[i] = svd->sigma[k];
432: }
433: PetscCall(PetscViewerBinaryWrite(viewer,sv,svd->nconv,PETSC_REAL));
434: PetscCall(PetscFree(sv));
435: PetscFunctionReturn(PETSC_SUCCESS);
436: }
438: #if defined(PETSC_HAVE_HDF5)
439: static PetscErrorCode SVDValuesView_HDF5(SVD svd,PetscViewer viewer)
440: {
441: PetscInt i,k,n,N;
442: PetscMPIInt rank;
443: Vec v;
444: char vname[30];
445: const char *ename;
447: PetscFunctionBegin;
448: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)svd),&rank));
449: N = svd->nconv;
450: n = rank? 0: N;
451: /* create a vector containing the singular values */
452: PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)svd),n,N,&v));
453: PetscCall(PetscObjectGetName((PetscObject)svd,&ename));
454: PetscCall(PetscSNPrintf(vname,sizeof(vname),"sigma_%s",ename));
455: PetscCall(PetscObjectSetName((PetscObject)v,vname));
456: if (!rank) {
457: for (i=0;i<svd->nconv;i++) {
458: k = svd->perm[i];
459: PetscCall(VecSetValue(v,i,svd->sigma[k],INSERT_VALUES));
460: }
461: }
462: PetscCall(VecAssemblyBegin(v));
463: PetscCall(VecAssemblyEnd(v));
464: PetscCall(VecView(v,viewer));
465: PetscCall(VecDestroy(&v));
466: PetscFunctionReturn(PETSC_SUCCESS);
467: }
468: #endif
470: static PetscErrorCode SVDValuesView_ASCII(SVD svd,PetscViewer viewer)
471: {
472: PetscInt i;
474: PetscFunctionBegin;
475: PetscCall(PetscViewerASCIIPrintf(viewer,"Singular values = \n"));
476: for (i=0;i<svd->nconv;i++) PetscCall(PetscViewerASCIIPrintf(viewer," %.5f\n",(double)svd->sigma[svd->perm[i]]));
477: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
478: PetscFunctionReturn(PETSC_SUCCESS);
479: }
481: static PetscErrorCode SVDValuesView_MATLAB(SVD svd,PetscViewer viewer)
482: {
483: PetscInt i;
484: const char *name;
486: PetscFunctionBegin;
487: PetscCall(PetscObjectGetName((PetscObject)svd,&name));
488: PetscCall(PetscViewerASCIIPrintf(viewer,"Sigma_%s = [\n",name));
489: for (i=0;i<svd->nconv;i++) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)svd->sigma[svd->perm[i]]));
490: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
491: PetscFunctionReturn(PETSC_SUCCESS);
492: }
494: /*@
495: SVDValuesView - Displays the computed singular values in a viewer.
497: Collective
499: Input Parameters:
500: + svd - the singular value solver context
501: - viewer - the viewer
503: Options Database Key:
504: . -svd_view_values - print computed singular values
506: Note:
507: The command-line option listed above admits an optional argument
508: specifying the viewer type and options. For instance, use
509: `-svd_view_values :evals.m:ascii_matlab` to save the values in a file
510: that can be executed in Matlab.
512: Level: intermediate
514: .seealso: [](ch:svd), `SVDSolve()`, `SVDVectorsView()`, `SVDErrorView()`
515: @*/
516: PetscErrorCode SVDValuesView(SVD svd,PetscViewer viewer)
517: {
518: PetscBool isascii,isdraw,isbinary;
519: PetscViewerFormat format;
520: #if defined(PETSC_HAVE_HDF5)
521: PetscBool ishdf5;
522: #endif
524: PetscFunctionBegin;
526: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
528: PetscCheckSameComm(svd,1,viewer,2);
529: SVDCheckSolved(svd,1);
530: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
531: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
532: #if defined(PETSC_HAVE_HDF5)
533: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
534: #endif
535: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
536: if (isdraw) PetscCall(SVDValuesView_DRAW(svd,viewer));
537: else if (isbinary) PetscCall(SVDValuesView_BINARY(svd,viewer));
538: #if defined(PETSC_HAVE_HDF5)
539: else if (ishdf5) PetscCall(SVDValuesView_HDF5(svd,viewer));
540: #endif
541: else if (isascii) {
542: PetscCall(PetscViewerGetFormat(viewer,&format));
543: switch (format) {
544: case PETSC_VIEWER_DEFAULT:
545: case PETSC_VIEWER_ASCII_INFO:
546: case PETSC_VIEWER_ASCII_INFO_DETAIL:
547: PetscCall(SVDValuesView_ASCII(svd,viewer));
548: break;
549: case PETSC_VIEWER_ASCII_MATLAB:
550: PetscCall(SVDValuesView_MATLAB(svd,viewer));
551: break;
552: default:
553: PetscCall(PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
554: }
555: }
556: PetscFunctionReturn(PETSC_SUCCESS);
557: }
559: /*@
560: SVDValuesViewFromOptions - Processes command line options to determine if/how
561: the computed singular values are to be viewed.
563: Collective
565: Input Parameter:
566: . svd - the singular value solver context
568: Level: developer
570: .seealso: [](ch:svd), `SVDValuesView()`
571: @*/
572: PetscErrorCode SVDValuesViewFromOptions(SVD svd)
573: {
574: PetscViewer viewer;
575: PetscBool flg;
576: static PetscBool incall = PETSC_FALSE;
577: PetscViewerFormat format;
579: PetscFunctionBegin;
580: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
581: incall = PETSC_TRUE;
582: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_values",&viewer,&format,&flg));
583: if (flg) {
584: PetscCall(PetscViewerPushFormat(viewer,format));
585: PetscCall(SVDValuesView(svd,viewer));
586: PetscCall(PetscViewerPopFormat(viewer));
587: PetscCall(PetscViewerDestroy(&viewer));
588: }
589: incall = PETSC_FALSE;
590: PetscFunctionReturn(PETSC_SUCCESS);
591: }
593: /*@
594: SVDVectorsView - Outputs computed singular vectors to a viewer.
596: Collective
598: Input Parameters:
599: + svd - the singular value solver context
600: - viewer - the viewer
602: Options Database Key:
603: . -svd_view_vectors - output singular vectors
605: Notes:
606: Right and left singular vectors are interleaved, that is, the vectors are
607: output in the following order\: `V0, U0, V1, U1, V2, U2, ...`
609: The command-line option listed above admits an optional argument
610: specifying the viewer type and options. For instance, use
611: `-svd_view_vectors binary:svecs.bin` to save the vectors in a binary file.
613: Level: intermediate
615: .seealso: [](ch:svd), `SVDSolve()`, `SVDValuesView()`, `SVDErrorView()`
616: @*/
617: PetscErrorCode SVDVectorsView(SVD svd,PetscViewer viewer)
618: {
619: PetscInt i,k;
620: Vec x;
621: char vname[30];
622: const char *ename;
624: PetscFunctionBegin;
626: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
628: PetscCheckSameComm(svd,1,viewer,2);
629: SVDCheckSolved(svd,1);
630: if (svd->nconv) {
631: PetscCall(PetscObjectGetName((PetscObject)svd,&ename));
632: PetscCall(SVDComputeVectors(svd));
633: for (i=0;i<svd->nconv;i++) {
634: k = svd->perm[i];
635: PetscCall(PetscSNPrintf(vname,sizeof(vname),"V%" PetscInt_FMT "_%s",i,ename));
636: PetscCall(BVGetColumn(svd->V,k,&x));
637: PetscCall(PetscObjectSetName((PetscObject)x,vname));
638: PetscCall(VecView(x,viewer));
639: PetscCall(BVRestoreColumn(svd->V,k,&x));
640: PetscCall(PetscSNPrintf(vname,sizeof(vname),"U%" PetscInt_FMT "_%s",i,ename));
641: PetscCall(BVGetColumn(svd->U,k,&x));
642: PetscCall(PetscObjectSetName((PetscObject)x,vname));
643: PetscCall(VecView(x,viewer));
644: PetscCall(BVRestoreColumn(svd->U,k,&x));
645: }
646: }
647: PetscFunctionReturn(PETSC_SUCCESS);
648: }
650: /*@
651: SVDVectorsViewFromOptions - Processes command line options to determine if/how
652: the computed singular vectors are to be viewed.
654: Collective
656: Input Parameter:
657: . svd - the singular value solver context
659: Level: developer
661: .seealso: [](ch:svd), `SVDVectorsView()`
662: @*/
663: PetscErrorCode SVDVectorsViewFromOptions(SVD svd)
664: {
665: PetscViewer viewer;
666: PetscBool flg = PETSC_FALSE;
667: static PetscBool incall = PETSC_FALSE;
668: PetscViewerFormat format;
670: PetscFunctionBegin;
671: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
672: incall = PETSC_TRUE;
673: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_vectors",&viewer,&format,&flg));
674: if (flg) {
675: PetscCall(PetscViewerPushFormat(viewer,format));
676: PetscCall(SVDVectorsView(svd,viewer));
677: PetscCall(PetscViewerPopFormat(viewer));
678: PetscCall(PetscViewerDestroy(&viewer));
679: }
680: incall = PETSC_FALSE;
681: PetscFunctionReturn(PETSC_SUCCESS);
682: }