Actual source code: svdmon.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: SVD routines related to monitors
12: */
14: #include <slepc/private/svdimpl.h>
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest)
21: {
22: PetscInt i,n = svd->numbermonitors;
24: PetscFunctionBegin;
25: for (i=0;i<n;i++) PetscCall((*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@C
30: SVDMonitorSet - Sets an ADDITIONAL function to be called at every
31: iteration to monitor the error estimates for each requested singular triplet.
33: Logically Collective
35: Input Parameters:
36: + svd - singular value solver context obtained from SVDCreate()
37: . monitor - pointer to function (if this is NULL, it turns off monitoring), see SVDMonitorFn
38: . mctx - [optional] context for private data for the
39: monitor routine (use NULL if no context is desired)
40: - monitordestroy - [optional] routine that frees monitor context (may be NULL),
41: see PetscCtxDestroyFn for the calling sequence
43: Options Database Keys:
44: + -svd_monitor - print only the first error estimate
45: . -svd_monitor_all - print error estimates at each iteration
46: . -svd_monitor_conv - print the singular value approximations only when
47: convergence has been reached
48: . -svd_monitor_conditioning - print the condition number when available
49: . -svd_monitor draw::draw_lg - sets line graph monitor for the first unconverged
50: approximate singular value
51: . -svd_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
52: approximate singular values
53: . -svd_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
54: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
55: a code by calls to SVDMonitorSet(), but does not cancel those set via
56: the options database.
58: Notes:
59: The options database option -svd_monitor and related options are the easiest way
60: to turn on SVD iteration monitoring.
62: SVDMonitorRegister() provides a way to associate an options database key with SVD
63: monitor function.
65: Several different monitoring routines may be set by calling SVDMonitorSet() multiple
66: times; all will be called in the order in which they were set.
68: Level: intermediate
70: .seealso: `SVDMonitorFirst()`, `SVDMonitorAll()`, `SVDMonitorConditioning()`, `SVDMonitorCancel()`
71: @*/
72: PetscErrorCode SVDMonitorSet(SVD svd,SVDMonitorFn *monitor,void *mctx,PetscCtxDestroyFn *monitordestroy)
73: {
74: PetscInt i;
75: PetscBool identical;
77: PetscFunctionBegin;
79: for (i=0;i<svd->numbermonitors;i++) {
80: PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))(PetscVoidFn*)monitor,mctx,monitordestroy,(PetscErrorCode (*)(void))(PetscVoidFn*)svd->monitor[i],svd->monitorcontext[i],svd->monitordestroy[i],&identical));
81: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
82: }
83: PetscCheck(svd->numbermonitors<MAXSVDMONITORS,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
84: svd->monitor[svd->numbermonitors] = monitor;
85: svd->monitorcontext[svd->numbermonitors] = mctx;
86: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
87: PetscFunctionReturn(PETSC_SUCCESS);
88: }
90: /*@
91: SVDMonitorCancel - Clears all monitors for an SVD object.
93: Logically Collective
95: Input Parameters:
96: . svd - singular value solver context obtained from SVDCreate()
98: Options Database Key:
99: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
100: into a code by calls to SVDMonitorSet(),
101: but does not cancel those set via the options database.
103: Level: intermediate
105: .seealso: `SVDMonitorSet()`
106: @*/
107: PetscErrorCode SVDMonitorCancel(SVD svd)
108: {
109: PetscInt i;
111: PetscFunctionBegin;
113: for (i=0; i<svd->numbermonitors; i++) {
114: if (svd->monitordestroy[i]) PetscCall((*svd->monitordestroy[i])(&svd->monitorcontext[i]));
115: }
116: svd->numbermonitors = 0;
117: PetscFunctionReturn(PETSC_SUCCESS);
118: }
120: /*@C
121: SVDGetMonitorContext - Gets the monitor context, as set by
122: SVDMonitorSet() for the FIRST monitor only.
124: Not Collective
126: Input Parameter:
127: . svd - singular value solver context obtained from SVDCreate()
129: Output Parameter:
130: . ctx - monitor context
132: Level: intermediate
134: .seealso: `SVDMonitorSet()`
135: @*/
136: PetscErrorCode SVDGetMonitorContext(SVD svd,void *ctx)
137: {
138: PetscFunctionBegin;
140: *(void**)ctx = svd->monitorcontext[0];
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
144: /*@C
145: SVDMonitorFirst - Print the first unconverged approximate value and
146: error estimate at each iteration of the singular value solver.
148: Collective
150: Input Parameters:
151: + svd - singular value solver context
152: . its - iteration number
153: . nconv - number of converged singular triplets so far
154: . sigma - singular values
155: . errest - error estimates
156: . nest - number of error estimates to display
157: - vf - viewer and format for monitoring
159: Options Database Key:
160: . -svd_monitor - activates SVDMonitorFirst()
162: Level: intermediate
164: .seealso: `SVDMonitorSet()`, `SVDMonitorAll()`, `SVDMonitorConditioning()`, `SVDMonitorConverged()`
165: @*/
166: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
167: {
168: PetscViewer viewer = vf->viewer;
170: PetscFunctionBegin;
173: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
174: if (nconv<nest) {
175: PetscCall(PetscViewerPushFormat(viewer,vf->format));
176: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
177: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
178: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
179: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]));
180: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
181: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
182: PetscCall(PetscViewerPopFormat(viewer));
183: }
184: PetscFunctionReturn(PETSC_SUCCESS);
185: }
187: /*@C
188: SVDMonitorAll - Print the current approximate values and
189: error estimates at each iteration of the singular value solver.
191: Collective
193: Input Parameters:
194: + svd - singular value solver context
195: . its - iteration number
196: . nconv - number of converged singular triplets so far
197: . sigma - singular values
198: . errest - error estimates
199: . nest - number of error estimates to display
200: - vf - viewer and format for monitoring
202: Options Database Key:
203: . -svd_monitor_all - activates SVDMonitorAll()
205: Level: intermediate
207: .seealso: `SVDMonitorSet()`, `SVDMonitorFirst()`, `SVDMonitorConditioning()`, `SVDMonitorConverged()`
208: @*/
209: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
210: {
211: PetscInt i;
212: PetscViewer viewer = vf->viewer;
214: PetscFunctionBegin;
217: PetscCall(PetscViewerPushFormat(viewer,vf->format));
218: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
219: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
220: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
221: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
222: for (i=0;i<nest;i++) PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]));
223: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
224: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
225: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
226: PetscCall(PetscViewerPopFormat(viewer));
227: PetscFunctionReturn(PETSC_SUCCESS);
228: }
230: /*@C
231: SVDMonitorConverged - Print the approximate values and
232: error estimates as they converge.
234: Collective
236: Input Parameters:
237: + svd - singular value solver context
238: . its - iteration number
239: . nconv - number of converged singular triplets so far
240: . sigma - singular values
241: . errest - error estimates
242: . nest - number of error estimates to display
243: - vf - viewer and format for monitoring
245: Options Database Key:
246: . -svd_monitor_conv - activates SVDMonitorConverged()
248: Level: intermediate
250: .seealso: `SVDMonitorSet()`, `SVDMonitorFirst()`, `SVDMonitorConditioning()`, `SVDMonitorAll()`
251: @*/
252: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
253: {
254: PetscInt i;
255: PetscViewer viewer = vf->viewer;
256: SlepcConvMon ctx;
258: PetscFunctionBegin;
261: ctx = (SlepcConvMon)vf->data;
262: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)svd)->prefix));
263: if (its==1) ctx->oldnconv = 0;
264: if (ctx->oldnconv!=nconv) {
265: PetscCall(PetscViewerPushFormat(viewer,vf->format));
266: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
267: for (i=ctx->oldnconv;i<nconv;i++) {
268: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD converged value (error) #%" PetscInt_FMT,its,i));
269: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
270: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]));
271: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
272: }
273: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
274: PetscCall(PetscViewerPopFormat(viewer));
275: ctx->oldnconv = nconv;
276: }
277: PetscFunctionReturn(PETSC_SUCCESS);
278: }
280: PetscErrorCode SVDMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
281: {
282: SlepcConvMon mctx;
284: PetscFunctionBegin;
285: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
286: PetscCall(PetscNew(&mctx));
287: mctx->ctx = ctx;
288: (*vf)->data = (void*)mctx;
289: PetscFunctionReturn(PETSC_SUCCESS);
290: }
292: PetscErrorCode SVDMonitorConvergedDestroy(PetscViewerAndFormat **vf)
293: {
294: PetscFunctionBegin;
295: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
296: PetscCall(PetscFree((*vf)->data));
297: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
298: PetscCall(PetscFree(*vf));
299: PetscFunctionReturn(PETSC_SUCCESS);
300: }
302: /*@C
303: SVDMonitorFirstDrawLG - Plots the error estimate of the first unconverged
304: approximation at each iteration of the singular value solver.
306: Collective
308: Input Parameters:
309: + svd - singular value solver context
310: . its - iteration number
311: . nconv - number of converged singular triplets so far
312: . sigma - singular values
313: . errest - error estimates
314: . nest - number of error estimates to display
315: - vf - viewer and format for monitoring
317: Options Database Key:
318: . -svd_monitor draw::draw_lg - activates SVDMonitorFirstDrawLG()
320: Level: intermediate
322: .seealso: `SVDMonitorSet()`
323: @*/
324: PetscErrorCode SVDMonitorFirstDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
325: {
326: PetscViewer viewer = vf->viewer;
327: PetscDrawLG lg;
328: PetscReal x,y;
330: PetscFunctionBegin;
333: PetscCall(PetscViewerPushFormat(viewer,vf->format));
334: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
335: if (its==1) {
336: PetscCall(PetscDrawLGReset(lg));
337: PetscCall(PetscDrawLGSetDimension(lg,1));
338: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
339: }
340: if (nconv<nest) {
341: x = (PetscReal)its;
342: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
343: else y = 0.0;
344: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
345: if (its <= 20 || !(its % 5) || svd->reason) {
346: PetscCall(PetscDrawLGDraw(lg));
347: PetscCall(PetscDrawLGSave(lg));
348: }
349: }
350: PetscCall(PetscViewerPopFormat(viewer));
351: PetscFunctionReturn(PETSC_SUCCESS);
352: }
354: /*@C
355: SVDMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
357: Collective
359: Input Parameters:
360: + viewer - the viewer
361: . format - the viewer format
362: - ctx - an optional user context
364: Output Parameter:
365: . vf - the viewer and format context
367: Level: intermediate
369: .seealso: `SVDMonitorSet()`
370: @*/
371: PetscErrorCode SVDMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
372: {
373: PetscFunctionBegin;
374: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
375: (*vf)->data = ctx;
376: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
377: PetscFunctionReturn(PETSC_SUCCESS);
378: }
380: /*@C
381: SVDMonitorAllDrawLG - Plots the error estimate of all unconverged
382: approximations at each iteration of the singular value solver.
384: Collective
386: Input Parameters:
387: + svd - singular value solver context
388: . its - iteration number
389: . nconv - number of converged singular triplets so far
390: . sigma - singular values
391: . errest - error estimates
392: . nest - number of error estimates to display
393: - vf - viewer and format for monitoring
395: Options Database Key:
396: . -svd_monitor_all draw::draw_lg - activates SVDMonitorAllDrawLG()
398: Level: intermediate
400: .seealso: `SVDMonitorSet()`
401: @*/
402: PetscErrorCode SVDMonitorAllDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
403: {
404: PetscViewer viewer = vf->viewer;
405: PetscDrawLG lg;
406: PetscInt i,n = PetscMin(svd->nsv,255);
407: PetscReal *x,*y;
409: PetscFunctionBegin;
412: PetscCall(PetscViewerPushFormat(viewer,vf->format));
413: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
414: if (its==1) {
415: PetscCall(PetscDrawLGReset(lg));
416: PetscCall(PetscDrawLGSetDimension(lg,n));
417: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
418: }
419: PetscCall(PetscMalloc2(n,&x,n,&y));
420: for (i=0;i<n;i++) {
421: x[i] = (PetscReal)its;
422: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
423: else y[i] = 0.0;
424: }
425: PetscCall(PetscDrawLGAddPoint(lg,x,y));
426: if (its <= 20 || !(its % 5) || svd->reason) {
427: PetscCall(PetscDrawLGDraw(lg));
428: PetscCall(PetscDrawLGSave(lg));
429: }
430: PetscCall(PetscFree2(x,y));
431: PetscCall(PetscViewerPopFormat(viewer));
432: PetscFunctionReturn(PETSC_SUCCESS);
433: }
435: /*@C
436: SVDMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
438: Collective
440: Input Parameters:
441: + viewer - the viewer
442: . format - the viewer format
443: - ctx - an optional user context
445: Output Parameter:
446: . vf - the viewer and format context
448: Level: intermediate
450: .seealso: `SVDMonitorSet()`
451: @*/
452: PetscErrorCode SVDMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
453: {
454: PetscFunctionBegin;
455: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
456: (*vf)->data = ctx;
457: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
458: PetscFunctionReturn(PETSC_SUCCESS);
459: }
461: /*@C
462: SVDMonitorConvergedDrawLG - Plots the number of converged eigenvalues
463: at each iteration of the singular value solver.
465: Collective
467: Input Parameters:
468: + svd - singular value solver context
469: . its - iteration number
470: . nconv - number of converged singular triplets so far
471: . sigma - singular values
472: . errest - error estimates
473: . nest - number of error estimates to display
474: - vf - viewer and format for monitoring
476: Options Database Key:
477: . -svd_monitor_conv draw::draw_lg - activates SVDMonitorConvergedDrawLG()
479: Level: intermediate
481: .seealso: `SVDMonitorSet()`
482: @*/
483: PetscErrorCode SVDMonitorConvergedDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
484: {
485: PetscViewer viewer = vf->viewer;
486: PetscDrawLG lg;
487: PetscReal x,y;
489: PetscFunctionBegin;
492: PetscCall(PetscViewerPushFormat(viewer,vf->format));
493: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
494: if (its==1) {
495: PetscCall(PetscDrawLGReset(lg));
496: PetscCall(PetscDrawLGSetDimension(lg,1));
497: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,svd->nsv));
498: }
499: x = (PetscReal)its;
500: y = (PetscReal)svd->nconv;
501: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
502: if (its <= 20 || !(its % 5) || svd->reason) {
503: PetscCall(PetscDrawLGDraw(lg));
504: PetscCall(PetscDrawLGSave(lg));
505: }
506: PetscCall(PetscViewerPopFormat(viewer));
507: PetscFunctionReturn(PETSC_SUCCESS);
508: }
510: /*@C
511: SVDMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
513: Collective
515: Input Parameters:
516: + viewer - the viewer
517: . format - the viewer format
518: - ctx - an optional user context
520: Output Parameter:
521: . vf - the viewer and format context
523: Level: intermediate
525: .seealso: `SVDMonitorSet()`
526: @*/
527: PetscErrorCode SVDMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
528: {
529: SlepcConvMon mctx;
531: PetscFunctionBegin;
532: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
533: PetscCall(PetscNew(&mctx));
534: mctx->ctx = ctx;
535: (*vf)->data = (void*)mctx;
536: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
537: PetscFunctionReturn(PETSC_SUCCESS);
538: }
540: /*@C
541: SVDMonitorConditioning - Print the condition number at each iteration of the singular
542: value solver.
544: Collective
546: Input Parameters:
547: + svd - singular value solver context
548: . its - iteration number
549: . nconv - (unused) number of converged singular triplets so far
550: . sigma - (unused) singular values
551: . errest - (unused) error estimates
552: . nest - (unused) number of error estimates to display
553: - vf - viewer and format for monitoring
555: Options Database Key:
556: . -svd_monitor_conditioning - activates SVDMonitorConditioning()
558: Note:
559: Works only for solvers that use a DS of type GSVD. The printed information corresponds
560: to the maximum of the condition number of the two generated bidiagonal matrices.
562: Level: intermediate
564: .seealso: `SVDMonitorSet()`, `SVDMonitorAll()`, `SVDMonitorFirst()`, `SVDMonitorConverged()`
565: @*/
566: PetscErrorCode SVDMonitorConditioning(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
567: {
568: PetscViewer viewer = vf->viewer;
569: PetscBool isgsvd;
570: PetscReal cond;
572: PetscFunctionBegin;
575: PetscCall(PetscObjectTypeCompare((PetscObject)svd->ds,DSGSVD,&isgsvd));
576: if (!isgsvd) PetscFunctionReturn(PETSC_SUCCESS);
577: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Condition number of bidiagonal matrices for %s solve.\n",((PetscObject)svd)->prefix));
578: PetscCall(PetscViewerPushFormat(viewer,vf->format));
579: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
580: PetscCall(DSCond(svd->ds,&cond));
581: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD condition number = %g\n",its,(double)cond));
582: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
583: PetscCall(PetscViewerPopFormat(viewer));
584: PetscFunctionReturn(PETSC_SUCCESS);
585: }