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