Actual source code: svdmon.c
slepc-3.22.1 2024-10-28
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: SVD routines related to monitors
12: */
14: #include <slepc/private/svdimpl.h>
15: #include <petscdraw.h>
17: PetscErrorCode SVDMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],const char metric[],PetscInt l,const char *names[],int x,int y,int m,int n,PetscDrawLG *lgctx)
18: {
19: PetscDraw draw;
20: PetscDrawAxis axis;
21: PetscDrawLG lg;
23: PetscFunctionBegin;
24: PetscCall(PetscDrawCreate(comm,host,label,x,y,m,n,&draw));
25: PetscCall(PetscDrawSetFromOptions(draw));
26: PetscCall(PetscDrawLGCreate(draw,l,&lg));
27: if (names) PetscCall(PetscDrawLGSetLegend(lg,names));
28: PetscCall(PetscDrawLGSetFromOptions(lg));
29: PetscCall(PetscDrawLGGetAxis(lg,&axis));
30: PetscCall(PetscDrawAxisSetLabels(axis,"Convergence","Iteration",metric));
31: PetscCall(PetscDrawDestroy(&draw));
32: *lgctx = lg;
33: PetscFunctionReturn(PETSC_SUCCESS);
34: }
36: /*
37: Runs the user provided monitor routines, if any.
38: */
39: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest)
40: {
41: PetscInt i,n = svd->numbermonitors;
43: PetscFunctionBegin;
44: for (i=0;i<n;i++) PetscCall((*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]));
45: PetscFunctionReturn(PETSC_SUCCESS);
46: }
48: /*@C
49: SVDMonitorSet - Sets an ADDITIONAL function to be called at every
50: iteration to monitor the error estimates for each requested singular triplet.
52: Logically Collective
54: Input Parameters:
55: + svd - singular value solver context obtained from SVDCreate()
56: . monitor - pointer to function (if this is NULL, it turns off monitoring)
57: . mctx - [optional] context for private data for the
58: monitor routine (use NULL if no context is desired)
59: - monitordestroy - [optional] routine that frees monitor context (may be NULL)
61: Calling sequence of monitor:
62: $ PetscErrorCode monitor(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *mctx)
63: + svd - singular value solver context obtained from SVDCreate()
64: . its - iteration number
65: . nconv - number of converged singular triplets
66: . sigma - singular values
67: . errest - relative error estimates for each singular triplet
68: . nest - number of error estimates
69: - mctx - optional monitoring context, as set by SVDMonitorSet()
71: Options Database Keys:
72: + -svd_monitor - print only the first error estimate
73: . -svd_monitor_all - print error estimates at each iteration
74: . -svd_monitor_conv - print the singular value approximations only when
75: convergence has been reached
76: . -svd_monitor_conditioning - print the condition number when available
77: . -svd_monitor draw::draw_lg - sets line graph monitor for the first unconverged
78: approximate singular value
79: . -svd_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
80: approximate singular values
81: . -svd_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
82: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
83: a code by calls to SVDMonitorSet(), but does not cancel those set via
84: the options database.
86: Notes:
87: Several different monitoring routines may be set by calling
88: SVDMonitorSet() multiple times; all will be called in the
89: order in which they were set.
91: Level: intermediate
93: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorConditioning(), SVDMonitorCancel()
94: @*/
95: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *mctx),void *mctx,PetscErrorCode (*monitordestroy)(void**))
96: {
97: PetscFunctionBegin;
99: PetscCheck(svd->numbermonitors<MAXSVDMONITORS,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
100: svd->monitor[svd->numbermonitors] = monitor;
101: svd->monitorcontext[svd->numbermonitors] = (void*)mctx;
102: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
103: PetscFunctionReturn(PETSC_SUCCESS);
104: }
106: /*@
107: SVDMonitorCancel - Clears all monitors for an SVD object.
109: Logically Collective
111: Input Parameters:
112: . svd - singular value solver context obtained from SVDCreate()
114: Options Database Key:
115: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
116: into a code by calls to SVDMonitorSet(),
117: but does not cancel those set via the options database.
119: Level: intermediate
121: .seealso: SVDMonitorSet()
122: @*/
123: PetscErrorCode SVDMonitorCancel(SVD svd)
124: {
125: PetscInt i;
127: PetscFunctionBegin;
129: for (i=0; i<svd->numbermonitors; i++) {
130: if (svd->monitordestroy[i]) PetscCall((*svd->monitordestroy[i])(&svd->monitorcontext[i]));
131: }
132: svd->numbermonitors = 0;
133: PetscFunctionReturn(PETSC_SUCCESS);
134: }
136: /*@C
137: SVDGetMonitorContext - Gets the monitor context, as set by
138: SVDMonitorSet() for the FIRST monitor only.
140: Not Collective
142: Input Parameter:
143: . svd - singular value solver context obtained from SVDCreate()
145: Output Parameter:
146: . ctx - monitor context
148: Level: intermediate
150: .seealso: SVDMonitorSet()
151: @*/
152: PetscErrorCode SVDGetMonitorContext(SVD svd,void *ctx)
153: {
154: PetscFunctionBegin;
156: *(void**)ctx = svd->monitorcontext[0];
157: PetscFunctionReturn(PETSC_SUCCESS);
158: }
160: /*@C
161: SVDMonitorFirst - Print the first unconverged approximate value and
162: error estimate at each iteration of the singular value solver.
164: Collective
166: Input Parameters:
167: + svd - singular value solver context
168: . its - iteration number
169: . nconv - number of converged singular triplets so far
170: . sigma - singular values
171: . errest - error estimates
172: . nest - number of error estimates to display
173: - vf - viewer and format for monitoring
175: Options Database Key:
176: . -svd_monitor - activates SVDMonitorFirst()
178: Level: intermediate
180: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConditioning(), SVDMonitorConverged()
181: @*/
182: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
183: {
184: PetscViewer viewer = vf->viewer;
186: PetscFunctionBegin;
189: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
190: if (nconv<nest) {
191: PetscCall(PetscViewerPushFormat(viewer,vf->format));
192: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
193: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
194: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
195: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]));
196: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
197: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
198: PetscCall(PetscViewerPopFormat(viewer));
199: }
200: PetscFunctionReturn(PETSC_SUCCESS);
201: }
203: /*@C
204: SVDMonitorAll - Print the current approximate values and
205: error estimates at each iteration of the singular value solver.
207: Collective
209: Input Parameters:
210: + svd - singular value solver context
211: . its - iteration number
212: . nconv - number of converged singular triplets so far
213: . sigma - singular values
214: . errest - error estimates
215: . nest - number of error estimates to display
216: - vf - viewer and format for monitoring
218: Options Database Key:
219: . -svd_monitor_all - activates SVDMonitorAll()
221: Level: intermediate
223: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConditioning(), SVDMonitorConverged()
224: @*/
225: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
226: {
227: PetscInt i;
228: PetscViewer viewer = vf->viewer;
230: PetscFunctionBegin;
233: PetscCall(PetscViewerPushFormat(viewer,vf->format));
234: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
235: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
236: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
237: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
238: for (i=0;i<nest;i++) PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]));
239: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
240: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
241: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
242: PetscCall(PetscViewerPopFormat(viewer));
243: PetscFunctionReturn(PETSC_SUCCESS);
244: }
246: /*@C
247: SVDMonitorConverged - Print the approximate values and
248: error estimates as they converge.
250: Collective
252: Input Parameters:
253: + svd - singular value solver context
254: . its - iteration number
255: . nconv - number of converged singular triplets so far
256: . sigma - singular values
257: . errest - error estimates
258: . nest - number of error estimates to display
259: - vf - viewer and format for monitoring
261: Options Database Key:
262: . -svd_monitor_conv - activates SVDMonitorConverged()
264: Level: intermediate
266: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConditioning(), SVDMonitorAll()
267: @*/
268: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
269: {
270: PetscInt i;
271: PetscViewer viewer = vf->viewer;
272: SlepcConvMon ctx;
274: PetscFunctionBegin;
277: ctx = (SlepcConvMon)vf->data;
278: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)svd)->prefix));
279: if (its==1) ctx->oldnconv = 0;
280: if (ctx->oldnconv!=nconv) {
281: PetscCall(PetscViewerPushFormat(viewer,vf->format));
282: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
283: for (i=ctx->oldnconv;i<nconv;i++) {
284: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD converged value (error) #%" PetscInt_FMT,its,i));
285: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
286: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]));
287: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
288: }
289: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
290: PetscCall(PetscViewerPopFormat(viewer));
291: ctx->oldnconv = nconv;
292: }
293: PetscFunctionReturn(PETSC_SUCCESS);
294: }
296: PetscErrorCode SVDMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
297: {
298: SlepcConvMon mctx;
300: PetscFunctionBegin;
301: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
302: PetscCall(PetscNew(&mctx));
303: mctx->ctx = ctx;
304: (*vf)->data = (void*)mctx;
305: PetscFunctionReturn(PETSC_SUCCESS);
306: }
308: PetscErrorCode SVDMonitorConvergedDestroy(PetscViewerAndFormat **vf)
309: {
310: PetscFunctionBegin;
311: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
312: PetscCall(PetscFree((*vf)->data));
313: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
314: PetscCall(PetscDrawLGDestroy(&(*vf)->lg));
315: PetscCall(PetscFree(*vf));
316: PetscFunctionReturn(PETSC_SUCCESS);
317: }
319: /*@C
320: SVDMonitorFirstDrawLG - Plots the error estimate of the first unconverged
321: approximation at each iteration of the singular value solver.
323: Collective
325: Input Parameters:
326: + svd - singular value solver context
327: . its - iteration number
328: . nconv - number of converged singular triplets so far
329: . sigma - singular values
330: . errest - error estimates
331: . nest - number of error estimates to display
332: - vf - viewer and format for monitoring
334: Options Database Key:
335: . -svd_monitor draw::draw_lg - activates SVDMonitorFirstDrawLG()
337: Level: intermediate
339: .seealso: SVDMonitorSet()
340: @*/
341: PetscErrorCode SVDMonitorFirstDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
342: {
343: PetscViewer viewer = vf->viewer;
344: PetscDrawLG lg = vf->lg;
345: PetscReal x,y;
347: PetscFunctionBegin;
351: PetscCall(PetscViewerPushFormat(viewer,vf->format));
352: if (its==1) {
353: PetscCall(PetscDrawLGReset(lg));
354: PetscCall(PetscDrawLGSetDimension(lg,1));
355: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
356: }
357: if (nconv<nest) {
358: x = (PetscReal)its;
359: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
360: else y = 0.0;
361: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
362: if (its <= 20 || !(its % 5) || svd->reason) {
363: PetscCall(PetscDrawLGDraw(lg));
364: PetscCall(PetscDrawLGSave(lg));
365: }
366: }
367: PetscCall(PetscViewerPopFormat(viewer));
368: PetscFunctionReturn(PETSC_SUCCESS);
369: }
371: /*@C
372: SVDMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
374: Collective
376: Input Parameters:
377: + viewer - the viewer
378: . format - the viewer format
379: - ctx - an optional user context
381: Output Parameter:
382: . vf - the viewer and format context
384: Level: intermediate
386: .seealso: SVDMonitorSet()
387: @*/
388: PetscErrorCode SVDMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
389: {
390: PetscFunctionBegin;
391: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
392: (*vf)->data = ctx;
393: PetscCall(SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
394: PetscFunctionReturn(PETSC_SUCCESS);
395: }
397: /*@C
398: SVDMonitorAllDrawLG - Plots the error estimate of all unconverged
399: approximations at each iteration of the singular value solver.
401: Collective
403: Input Parameters:
404: + svd - singular value solver context
405: . its - iteration number
406: . nconv - number of converged singular triplets so far
407: . sigma - singular values
408: . errest - error estimates
409: . nest - number of error estimates to display
410: - vf - viewer and format for monitoring
412: Options Database Key:
413: . -svd_monitor_all draw::draw_lg - activates SVDMonitorAllDrawLG()
415: Level: intermediate
417: .seealso: SVDMonitorSet()
418: @*/
419: PetscErrorCode SVDMonitorAllDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
420: {
421: PetscViewer viewer = vf->viewer;
422: PetscDrawLG lg = vf->lg;
423: PetscInt i,n = PetscMin(svd->nsv,255);
424: PetscReal *x,*y;
426: PetscFunctionBegin;
430: PetscCall(PetscViewerPushFormat(viewer,vf->format));
431: if (its==1) {
432: PetscCall(PetscDrawLGReset(lg));
433: PetscCall(PetscDrawLGSetDimension(lg,n));
434: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
435: }
436: PetscCall(PetscMalloc2(n,&x,n,&y));
437: for (i=0;i<n;i++) {
438: x[i] = (PetscReal)its;
439: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
440: else y[i] = 0.0;
441: }
442: PetscCall(PetscDrawLGAddPoint(lg,x,y));
443: if (its <= 20 || !(its % 5) || svd->reason) {
444: PetscCall(PetscDrawLGDraw(lg));
445: PetscCall(PetscDrawLGSave(lg));
446: }
447: PetscCall(PetscFree2(x,y));
448: PetscCall(PetscViewerPopFormat(viewer));
449: PetscFunctionReturn(PETSC_SUCCESS);
450: }
452: /*@C
453: SVDMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
455: Collective
457: Input Parameters:
458: + viewer - the viewer
459: . format - the viewer format
460: - ctx - an optional user context
462: Output Parameter:
463: . vf - the viewer and format context
465: Level: intermediate
467: .seealso: SVDMonitorSet()
468: @*/
469: PetscErrorCode SVDMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
470: {
471: PetscFunctionBegin;
472: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
473: (*vf)->data = ctx;
474: PetscCall(SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
475: PetscFunctionReturn(PETSC_SUCCESS);
476: }
478: /*@C
479: SVDMonitorConvergedDrawLG - Plots the number of converged eigenvalues
480: at each iteration of the singular value solver.
482: Collective
484: Input Parameters:
485: + svd - singular value solver context
486: . its - iteration number
487: . nconv - number of converged singular triplets so far
488: . sigma - singular values
489: . errest - error estimates
490: . nest - number of error estimates to display
491: - vf - viewer and format for monitoring
493: Options Database Key:
494: . -svd_monitor_conv draw::draw_lg - activates SVDMonitorConvergedDrawLG()
496: Level: intermediate
498: .seealso: SVDMonitorSet()
499: @*/
500: PetscErrorCode SVDMonitorConvergedDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
501: {
502: PetscViewer viewer = vf->viewer;
503: PetscDrawLG lg = vf->lg;
504: PetscReal x,y;
506: PetscFunctionBegin;
510: PetscCall(PetscViewerPushFormat(viewer,vf->format));
511: if (its==1) {
512: PetscCall(PetscDrawLGReset(lg));
513: PetscCall(PetscDrawLGSetDimension(lg,1));
514: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,svd->nsv));
515: }
516: x = (PetscReal)its;
517: y = (PetscReal)svd->nconv;
518: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
519: if (its <= 20 || !(its % 5) || svd->reason) {
520: PetscCall(PetscDrawLGDraw(lg));
521: PetscCall(PetscDrawLGSave(lg));
522: }
523: PetscCall(PetscViewerPopFormat(viewer));
524: PetscFunctionReturn(PETSC_SUCCESS);
525: }
527: /*@C
528: SVDMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
530: Collective
532: Input Parameters:
533: + viewer - the viewer
534: . format - the viewer format
535: - ctx - an optional user context
537: Output Parameter:
538: . vf - the viewer and format context
540: Level: intermediate
542: .seealso: SVDMonitorSet()
543: @*/
544: PetscErrorCode SVDMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
545: {
546: SlepcConvMon mctx;
548: PetscFunctionBegin;
549: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
550: PetscCall(PetscNew(&mctx));
551: mctx->ctx = ctx;
552: (*vf)->data = (void*)mctx;
553: PetscCall(SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
554: PetscFunctionReturn(PETSC_SUCCESS);
555: }
557: /*@C
558: SVDMonitorConditioning - Print the condition number at each iteration of the singular
559: value solver.
561: Collective
563: Input Parameters:
564: + svd - singular value solver context
565: . its - iteration number
566: . nconv - (unused) number of converged singular triplets so far
567: . sigma - (unused) singular values
568: . errest - (unused) error estimates
569: . nest - (unused) number of error estimates to display
570: - vf - viewer and format for monitoring
572: Options Database Key:
573: . -svd_monitor_conditioning - activates SVDMonitorConditioning()
575: Note:
576: Works only for solvers that use a DS of type GSVD. The printed information corresponds
577: to the maximum of the condition number of the two generated bidiagonal matrices.
579: Level: intermediate
581: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorFirst(), SVDMonitorConverged()
582: @*/
583: PetscErrorCode SVDMonitorConditioning(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
584: {
585: PetscViewer viewer = vf->viewer;
586: PetscBool isgsvd;
587: PetscReal cond;
589: PetscFunctionBegin;
592: PetscCall(PetscObjectTypeCompare((PetscObject)svd->ds,DSGSVD,&isgsvd));
593: if (!isgsvd) PetscFunctionReturn(PETSC_SUCCESS);
594: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Condition number of bidiagonal matrices for %s solve.\n",((PetscObject)svd)->prefix));
595: PetscCall(PetscViewerPushFormat(viewer,vf->format));
596: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
597: PetscCall(DSCond(svd->ds,&cond));
598: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD condition number = %g\n",its,(double)cond));
599: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
600: PetscCall(PetscViewerPopFormat(viewer));
601: PetscFunctionReturn(PETSC_SUCCESS);
602: }