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;
215: PetscBool requestednsv=PETSC_TRUE; /* user has specified nsv */
216: const char *svstring = svd->isgeneralized?"generalized singular values":"singular values";
218: PetscFunctionBegin;
219: if (svd->stop==SVD_STOP_THRESHOLD && svd->nsv==0) requestednsv = PETSC_FALSE;
220: numsv = requestednsv? svd->nsv: svd->nconv;
221: if (requestednsv && svd->nconv<numsv) {
222: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " %s converged\n\n",numsv,svstring));
223: PetscFunctionReturn(PETSC_SUCCESS);
224: }
225: if (!requestednsv && !numsv) {
226: PetscCall(PetscViewerASCIIPrintf(viewer," No %s have been found\n\n",svstring));
227: PetscFunctionReturn(PETSC_SUCCESS);
228: }
229: for (i=0;i<numsv;i++) {
230: PetscCall(SVDComputeError(svd,i,etype,&error));
231: if (error>=5.0*svd->tol) {
232: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",numsv));
233: PetscFunctionReturn(PETSC_SUCCESS);
234: }
235: }
236: if (!requestednsv) PetscCall(PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " %s, all of them computed up to the required tolerance:",numsv,svstring));
237: else PetscCall(PetscViewerASCIIPrintf(viewer," All requested %s computed up to the required tolerance:",svstring));
238: for (i=0;i<=(numsv-1)/8;i++) {
239: PetscCall(PetscViewerASCIIPrintf(viewer,"\n "));
240: for (j=0;j<PetscMin(8,numsv-8*i);j++) {
241: PetscCall(SVDGetSingularTriplet(svd,8*i+j,&sigma,NULL,NULL));
242: PetscCall(PetscViewerASCIIPrintf(viewer,"%.5f",(double)sigma));
243: if (8*i+j+1<numsv) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
244: }
245: }
246: PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
247: PetscFunctionReturn(PETSC_SUCCESS);
248: }
250: static PetscErrorCode SVDErrorView_DETAIL(SVD svd,SVDErrorType etype,PetscViewer viewer)
251: {
252: PetscReal error,sigma;
253: PetscInt i;
254: char ex[30],sep[]=" ---------------------- --------------------\n";
256: PetscFunctionBegin;
257: if (!svd->nconv) PetscFunctionReturn(PETSC_SUCCESS);
258: switch (etype) {
259: case SVD_ERROR_ABSOLUTE:
260: PetscCall(PetscSNPrintf(ex,sizeof(ex)," absolute error"));
261: break;
262: case SVD_ERROR_RELATIVE:
263: PetscCall(PetscSNPrintf(ex,sizeof(ex)," relative error"));
264: break;
265: case SVD_ERROR_NORM:
266: if (svd->isgeneralized) PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||r||/||[A;B]||"));
267: else PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||r||/||A||"));
268: break;
269: }
270: PetscCall(PetscViewerASCIIPrintf(viewer,"%s sigma %s\n%s",sep,ex,sep));
271: for (i=0;i<svd->nconv;i++) {
272: PetscCall(SVDGetSingularTriplet(svd,i,&sigma,NULL,NULL));
273: PetscCall(SVDComputeError(svd,i,etype,&error));
274: PetscCall(PetscViewerASCIIPrintf(viewer," % 6f %12g\n",(double)sigma,(double)error));
275: }
276: PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
277: PetscFunctionReturn(PETSC_SUCCESS);
278: }
280: static PetscErrorCode SVDErrorView_MATLAB(SVD svd,SVDErrorType etype,PetscViewer viewer)
281: {
282: PetscReal error;
283: PetscInt i;
284: const char *name;
286: PetscFunctionBegin;
287: PetscCall(PetscObjectGetName((PetscObject)svd,&name));
288: PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
289: for (i=0;i<svd->nconv;i++) {
290: PetscCall(SVDComputeError(svd,i,etype,&error));
291: PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
292: }
293: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
294: PetscFunctionReturn(PETSC_SUCCESS);
295: }
297: /*@
298: SVDErrorView - Displays the errors associated with the computed solution
299: (as well as the singular values).
301: Collective
303: Input Parameters:
304: + svd - the singular value solver context
305: . etype - error type
306: - viewer - optional visualization context
308: Options Database Keys:
309: + -svd_error_absolute - print absolute errors of each singular triplet
310: . -svd_error_relative - print relative errors of each singular triplet
311: - -svd_error_norm - print errors relative to the matrix norms of each singular triplet
313: Notes:
314: By default, this function checks the error of all singular triplets and prints
315: the singular values if all of them are below the requested tolerance.
316: If the viewer has format `PETSC_VIEWER_ASCII_INFO_DETAIL` then a table with
317: singular values and corresponding errors is printed.
319: All the command-line options listed above admit an optional argument
320: specifying the viewer type and options. For instance, use
321: `-svd_error_relative :myerr.m:ascii_matlab` to save the errors in a file
322: that can be executed in Matlab.
324: Level: intermediate
326: .seealso: [](ch:svd), `SVDSolve()`, `SVDValuesView()`, `SVDVectorsView()`
327: @*/
328: PetscErrorCode SVDErrorView(SVD svd,SVDErrorType etype,PetscViewer viewer)
329: {
330: PetscBool isascii;
331: PetscViewerFormat format;
333: PetscFunctionBegin;
335: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
337: PetscCheckSameComm(svd,1,viewer,3);
338: SVDCheckSolved(svd,1);
339: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
340: if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);
342: PetscCall(PetscViewerGetFormat(viewer,&format));
343: switch (format) {
344: case PETSC_VIEWER_DEFAULT:
345: case PETSC_VIEWER_ASCII_INFO:
346: PetscCall(SVDErrorView_ASCII(svd,etype,viewer));
347: break;
348: case PETSC_VIEWER_ASCII_INFO_DETAIL:
349: PetscCall(SVDErrorView_DETAIL(svd,etype,viewer));
350: break;
351: case PETSC_VIEWER_ASCII_MATLAB:
352: PetscCall(SVDErrorView_MATLAB(svd,etype,viewer));
353: break;
354: default:
355: PetscCall(PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
356: }
357: PetscFunctionReturn(PETSC_SUCCESS);
358: }
360: /*@
361: SVDErrorViewFromOptions - Processes command line options to determine if/how
362: the errors of the computed solution are to be viewed.
364: Collective
366: Input Parameter:
367: . svd - the singular value solver context
369: Level: developer
371: .seealso: [](ch:svd), `SVDErrorView()`
372: @*/
373: PetscErrorCode SVDErrorViewFromOptions(SVD svd)
374: {
375: PetscViewer viewer;
376: PetscBool flg;
377: static PetscBool incall = PETSC_FALSE;
378: PetscViewerFormat format;
380: PetscFunctionBegin;
381: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
382: incall = PETSC_TRUE;
383: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_absolute",&viewer,&format,&flg));
384: if (flg) {
385: PetscCall(PetscViewerPushFormat(viewer,format));
386: PetscCall(SVDErrorView(svd,SVD_ERROR_ABSOLUTE,viewer));
387: PetscCall(PetscViewerPopFormat(viewer));
388: PetscCall(PetscViewerDestroy(&viewer));
389: }
390: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_relative",&viewer,&format,&flg));
391: if (flg) {
392: PetscCall(PetscViewerPushFormat(viewer,format));
393: PetscCall(SVDErrorView(svd,SVD_ERROR_RELATIVE,viewer));
394: PetscCall(PetscViewerPopFormat(viewer));
395: PetscCall(PetscViewerDestroy(&viewer));
396: }
397: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_error_norm",&viewer,&format,&flg));
398: if (flg) {
399: PetscCall(PetscViewerPushFormat(viewer,format));
400: PetscCall(SVDErrorView(svd,SVD_ERROR_NORM,viewer));
401: PetscCall(PetscViewerPopFormat(viewer));
402: PetscCall(PetscViewerDestroy(&viewer));
403: }
404: incall = PETSC_FALSE;
405: PetscFunctionReturn(PETSC_SUCCESS);
406: }
408: static PetscErrorCode SVDValuesView_DRAW(SVD svd,PetscViewer viewer)
409: {
410: PetscDraw draw;
411: PetscDrawSP drawsp;
412: PetscReal re,im=0.0;
413: PetscInt i;
415: PetscFunctionBegin;
416: if (!svd->nconv) PetscFunctionReturn(PETSC_SUCCESS);
417: PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
418: PetscCall(PetscDrawSetTitle(draw,"Computed singular values"));
419: PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
420: for (i=0;i<svd->nconv;i++) {
421: re = svd->sigma[svd->perm[i]];
422: PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
423: }
424: PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
425: PetscCall(PetscDrawSPSave(drawsp));
426: PetscCall(PetscDrawSPDestroy(&drawsp));
427: PetscFunctionReturn(PETSC_SUCCESS);
428: }
430: static PetscErrorCode SVDValuesView_BINARY(SVD svd,PetscViewer viewer)
431: {
432: PetscInt i,k;
433: PetscReal *sv;
435: PetscFunctionBegin;
436: PetscCall(PetscMalloc1(svd->nconv,&sv));
437: for (i=0;i<svd->nconv;i++) {
438: k = svd->perm[i];
439: sv[i] = svd->sigma[k];
440: }
441: PetscCall(PetscViewerBinaryWrite(viewer,sv,svd->nconv,PETSC_REAL));
442: PetscCall(PetscFree(sv));
443: PetscFunctionReturn(PETSC_SUCCESS);
444: }
446: #if defined(PETSC_HAVE_HDF5)
447: static PetscErrorCode SVDValuesView_HDF5(SVD svd,PetscViewer viewer)
448: {
449: PetscInt i,k,n,N;
450: PetscMPIInt rank;
451: Vec v;
452: char vname[30];
453: const char *ename;
455: PetscFunctionBegin;
456: PetscCallMPI(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: PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)svd),n,N,&v));
461: PetscCall(PetscObjectGetName((PetscObject)svd,&ename));
462: PetscCall(PetscSNPrintf(vname,sizeof(vname),"sigma_%s",ename));
463: PetscCall(PetscObjectSetName((PetscObject)v,vname));
464: if (!rank) {
465: for (i=0;i<svd->nconv;i++) {
466: k = svd->perm[i];
467: PetscCall(VecSetValue(v,i,svd->sigma[k],INSERT_VALUES));
468: }
469: }
470: PetscCall(VecAssemblyBegin(v));
471: PetscCall(VecAssemblyEnd(v));
472: PetscCall(VecView(v,viewer));
473: PetscCall(VecDestroy(&v));
474: PetscFunctionReturn(PETSC_SUCCESS);
475: }
476: #endif
478: static PetscErrorCode SVDValuesView_ASCII(SVD svd,PetscViewer viewer)
479: {
480: PetscInt i;
482: PetscFunctionBegin;
483: PetscCall(PetscViewerASCIIPrintf(viewer,"Singular values = \n"));
484: for (i=0;i<svd->nconv;i++) PetscCall(PetscViewerASCIIPrintf(viewer," %.5f\n",(double)svd->sigma[svd->perm[i]]));
485: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
486: PetscFunctionReturn(PETSC_SUCCESS);
487: }
489: static PetscErrorCode SVDValuesView_MATLAB(SVD svd,PetscViewer viewer)
490: {
491: PetscInt i;
492: const char *name;
494: PetscFunctionBegin;
495: PetscCall(PetscObjectGetName((PetscObject)svd,&name));
496: PetscCall(PetscViewerASCIIPrintf(viewer,"Sigma_%s = [\n",name));
497: for (i=0;i<svd->nconv;i++) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)svd->sigma[svd->perm[i]]));
498: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
499: PetscFunctionReturn(PETSC_SUCCESS);
500: }
502: /*@
503: SVDValuesView - Displays the computed singular values in a viewer.
505: Collective
507: Input Parameters:
508: + svd - the singular value solver context
509: - viewer - the viewer
511: Options Database Key:
512: . -svd_view_values - print computed singular values
514: Note:
515: The command-line option listed above admits an optional argument
516: specifying the viewer type and options. For instance, use
517: `-svd_view_values :evals.m:ascii_matlab` to save the values in a file
518: that can be executed in Matlab.
519: See `PetscObjectViewFromOptions()` for more details.
521: Level: intermediate
523: .seealso: [](ch:svd), `SVDSolve()`, `SVDVectorsView()`, `SVDErrorView()`
524: @*/
525: PetscErrorCode SVDValuesView(SVD svd,PetscViewer viewer)
526: {
527: PetscBool isascii,isdraw,isbinary;
528: PetscViewerFormat format;
529: #if defined(PETSC_HAVE_HDF5)
530: PetscBool ishdf5;
531: #endif
533: PetscFunctionBegin;
535: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
537: PetscCheckSameComm(svd,1,viewer,2);
538: SVDCheckSolved(svd,1);
539: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
540: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
541: #if defined(PETSC_HAVE_HDF5)
542: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
543: #endif
544: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
545: if (isdraw) PetscCall(SVDValuesView_DRAW(svd,viewer));
546: else if (isbinary) PetscCall(SVDValuesView_BINARY(svd,viewer));
547: #if defined(PETSC_HAVE_HDF5)
548: else if (ishdf5) PetscCall(SVDValuesView_HDF5(svd,viewer));
549: #endif
550: else if (isascii) {
551: PetscCall(PetscViewerGetFormat(viewer,&format));
552: switch (format) {
553: case PETSC_VIEWER_DEFAULT:
554: case PETSC_VIEWER_ASCII_INFO:
555: case PETSC_VIEWER_ASCII_INFO_DETAIL:
556: PetscCall(SVDValuesView_ASCII(svd,viewer));
557: break;
558: case PETSC_VIEWER_ASCII_MATLAB:
559: PetscCall(SVDValuesView_MATLAB(svd,viewer));
560: break;
561: default:
562: PetscCall(PetscInfo(svd,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
563: }
564: }
565: PetscFunctionReturn(PETSC_SUCCESS);
566: }
568: /*@
569: SVDValuesViewFromOptions - Processes command line options to determine if/how
570: the computed singular values are to be viewed.
572: Collective
574: Input Parameter:
575: . svd - the singular value solver context
577: Level: developer
579: .seealso: [](ch:svd), `SVDValuesView()`
580: @*/
581: PetscErrorCode SVDValuesViewFromOptions(SVD svd)
582: {
583: PetscViewer viewer;
584: PetscBool flg;
585: static PetscBool incall = PETSC_FALSE;
586: PetscViewerFormat format;
588: PetscFunctionBegin;
589: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
590: incall = PETSC_TRUE;
591: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_values",&viewer,&format,&flg));
592: if (flg) {
593: PetscCall(PetscViewerPushFormat(viewer,format));
594: PetscCall(SVDValuesView(svd,viewer));
595: PetscCall(PetscViewerPopFormat(viewer));
596: PetscCall(PetscViewerDestroy(&viewer));
597: }
598: incall = PETSC_FALSE;
599: PetscFunctionReturn(PETSC_SUCCESS);
600: }
602: /*@
603: SVDVectorsView - Outputs computed singular vectors to a viewer.
605: Collective
607: Input Parameters:
608: + svd - the singular value solver context
609: - viewer - the viewer
611: Options Database Key:
612: . -svd_view_vectors - output singular vectors
614: Notes:
615: Right and left singular vectors are interleaved, that is, the vectors are
616: output in the following order\: `V0, U0, V1, U1, V2, U2, ...`
618: The command-line option listed above admits an optional argument
619: specifying the viewer type and options. For instance, use
620: `-svd_view_vectors binary:svecs.bin` to save the vectors in a binary file.
621: See `PetscObjectViewFromOptions()` for more details.
623: Level: intermediate
625: .seealso: [](ch:svd), `SVDSolve()`, `SVDValuesView()`, `SVDErrorView()`
626: @*/
627: PetscErrorCode SVDVectorsView(SVD svd,PetscViewer viewer)
628: {
629: PetscInt i,k;
630: Vec x;
631: char vname[30];
632: const char *ename;
634: PetscFunctionBegin;
636: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)svd),&viewer));
638: PetscCheckSameComm(svd,1,viewer,2);
639: SVDCheckSolved(svd,1);
640: if (svd->nconv) {
641: PetscCall(PetscObjectGetName((PetscObject)svd,&ename));
642: PetscCall(SVDComputeVectors(svd));
643: for (i=0;i<svd->nconv;i++) {
644: k = svd->perm[i];
645: PetscCall(PetscSNPrintf(vname,sizeof(vname),"V%" PetscInt_FMT "_%s",i,ename));
646: PetscCall(BVGetColumn(svd->V,k,&x));
647: PetscCall(PetscObjectSetName((PetscObject)x,vname));
648: PetscCall(VecView(x,viewer));
649: PetscCall(BVRestoreColumn(svd->V,k,&x));
650: PetscCall(PetscSNPrintf(vname,sizeof(vname),"U%" PetscInt_FMT "_%s",i,ename));
651: PetscCall(BVGetColumn(svd->U,k,&x));
652: PetscCall(PetscObjectSetName((PetscObject)x,vname));
653: PetscCall(VecView(x,viewer));
654: PetscCall(BVRestoreColumn(svd->U,k,&x));
655: }
656: }
657: PetscFunctionReturn(PETSC_SUCCESS);
658: }
660: /*@
661: SVDVectorsViewFromOptions - Processes command line options to determine if/how
662: the computed singular vectors are to be viewed.
664: Collective
666: Input Parameter:
667: . svd - the singular value solver context
669: Level: developer
671: .seealso: [](ch:svd), `SVDVectorsView()`
672: @*/
673: PetscErrorCode SVDVectorsViewFromOptions(SVD svd)
674: {
675: PetscViewer viewer;
676: PetscBool flg = PETSC_FALSE;
677: static PetscBool incall = PETSC_FALSE;
678: PetscViewerFormat format;
680: PetscFunctionBegin;
681: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
682: incall = PETSC_TRUE;
683: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,"-svd_view_vectors",&viewer,&format,&flg));
684: if (flg) {
685: PetscCall(PetscViewerPushFormat(viewer,format));
686: PetscCall(SVDVectorsView(svd,viewer));
687: PetscCall(PetscViewerPopFormat(viewer));
688: PetscCall(PetscViewerDestroy(&viewer));
689: }
690: incall = PETSC_FALSE;
691: PetscFunctionReturn(PETSC_SUCCESS);
692: }