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