Actual source code: svdmon.c
slepc-main 2024-11-09
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),
60: see PetscCtxDestroyFn for the calling sequence
62: Calling sequence of monitor:
63: $ PetscErrorCode monitor(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *mctx)
64: + svd - singular value solver context obtained from SVDCreate()
65: . its - iteration number
66: . nconv - number of converged singular triplets
67: . sigma - singular values
68: . errest - relative error estimates for each singular triplet
69: . nest - number of error estimates
70: - mctx - optional monitoring context, as set by SVDMonitorSet()
72: Options Database Keys:
73: + -svd_monitor - print only the first error estimate
74: . -svd_monitor_all - print error estimates at each iteration
75: . -svd_monitor_conv - print the singular value approximations only when
76: convergence has been reached
77: . -svd_monitor_conditioning - print the condition number when available
78: . -svd_monitor draw::draw_lg - sets line graph monitor for the first unconverged
79: approximate singular value
80: . -svd_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
81: approximate singular values
82: . -svd_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
83: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
84: a code by calls to SVDMonitorSet(), but does not cancel those set via
85: the options database.
87: Notes:
88: Several different monitoring routines may be set by calling
89: SVDMonitorSet() multiple times; all will be called in the
90: order in which they were set.
92: Level: intermediate
94: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorConditioning(), SVDMonitorCancel()
95: @*/
96: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *mctx),void *mctx,PetscCtxDestroyFn *monitordestroy)
97: {
98: PetscFunctionBegin;
100: PetscCheck(svd->numbermonitors<MAXSVDMONITORS,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
101: svd->monitor[svd->numbermonitors] = monitor;
102: svd->monitorcontext[svd->numbermonitors] = (void*)mctx;
103: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
104: PetscFunctionReturn(PETSC_SUCCESS);
105: }
107: /*@
108: SVDMonitorCancel - Clears all monitors for an SVD object.
110: Logically Collective
112: Input Parameters:
113: . svd - singular value solver context obtained from SVDCreate()
115: Options Database Key:
116: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
117: into a code by calls to SVDMonitorSet(),
118: but does not cancel those set via the options database.
120: Level: intermediate
122: .seealso: SVDMonitorSet()
123: @*/
124: PetscErrorCode SVDMonitorCancel(SVD svd)
125: {
126: PetscInt i;
128: PetscFunctionBegin;
130: for (i=0; i<svd->numbermonitors; i++) {
131: if (svd->monitordestroy[i]) PetscCall((*svd->monitordestroy[i])(&svd->monitorcontext[i]));
132: }
133: svd->numbermonitors = 0;
134: PetscFunctionReturn(PETSC_SUCCESS);
135: }
137: /*@C
138: SVDGetMonitorContext - Gets the monitor context, as set by
139: SVDMonitorSet() for the FIRST monitor only.
141: Not Collective
143: Input Parameter:
144: . svd - singular value solver context obtained from SVDCreate()
146: Output Parameter:
147: . ctx - monitor context
149: Level: intermediate
151: .seealso: SVDMonitorSet()
152: @*/
153: PetscErrorCode SVDGetMonitorContext(SVD svd,void *ctx)
154: {
155: PetscFunctionBegin;
157: *(void**)ctx = svd->monitorcontext[0];
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
161: /*@C
162: SVDMonitorFirst - Print the first unconverged approximate value and
163: error estimate at each iteration of the singular value solver.
165: Collective
167: Input Parameters:
168: + svd - singular value solver context
169: . its - iteration number
170: . nconv - number of converged singular triplets so far
171: . sigma - singular values
172: . errest - error estimates
173: . nest - number of error estimates to display
174: - vf - viewer and format for monitoring
176: Options Database Key:
177: . -svd_monitor - activates SVDMonitorFirst()
179: Level: intermediate
181: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConditioning(), SVDMonitorConverged()
182: @*/
183: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
184: {
185: PetscViewer viewer = vf->viewer;
187: PetscFunctionBegin;
190: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
191: if (nconv<nest) {
192: PetscCall(PetscViewerPushFormat(viewer,vf->format));
193: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
194: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
195: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
196: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]));
197: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
198: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
199: PetscCall(PetscViewerPopFormat(viewer));
200: }
201: PetscFunctionReturn(PETSC_SUCCESS);
202: }
204: /*@C
205: SVDMonitorAll - Print the current approximate values and
206: error estimates at each iteration of the singular value solver.
208: Collective
210: Input Parameters:
211: + svd - singular value solver context
212: . its - iteration number
213: . nconv - number of converged singular triplets so far
214: . sigma - singular values
215: . errest - error estimates
216: . nest - number of error estimates to display
217: - vf - viewer and format for monitoring
219: Options Database Key:
220: . -svd_monitor_all - activates SVDMonitorAll()
222: Level: intermediate
224: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConditioning(), SVDMonitorConverged()
225: @*/
226: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
227: {
228: PetscInt i;
229: PetscViewer viewer = vf->viewer;
231: PetscFunctionBegin;
234: PetscCall(PetscViewerPushFormat(viewer,vf->format));
235: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
236: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
237: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
238: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
239: for (i=0;i<nest;i++) PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]));
240: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
241: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
242: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
243: PetscCall(PetscViewerPopFormat(viewer));
244: PetscFunctionReturn(PETSC_SUCCESS);
245: }
247: /*@C
248: SVDMonitorConverged - Print the approximate values and
249: error estimates as they converge.
251: Collective
253: Input Parameters:
254: + svd - singular value solver context
255: . its - iteration number
256: . nconv - number of converged singular triplets so far
257: . sigma - singular values
258: . errest - error estimates
259: . nest - number of error estimates to display
260: - vf - viewer and format for monitoring
262: Options Database Key:
263: . -svd_monitor_conv - activates SVDMonitorConverged()
265: Level: intermediate
267: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConditioning(), SVDMonitorAll()
268: @*/
269: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
270: {
271: PetscInt i;
272: PetscViewer viewer = vf->viewer;
273: SlepcConvMon ctx;
275: PetscFunctionBegin;
278: ctx = (SlepcConvMon)vf->data;
279: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)svd)->prefix));
280: if (its==1) ctx->oldnconv = 0;
281: if (ctx->oldnconv!=nconv) {
282: PetscCall(PetscViewerPushFormat(viewer,vf->format));
283: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
284: for (i=ctx->oldnconv;i<nconv;i++) {
285: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD converged value (error) #%" PetscInt_FMT,its,i));
286: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
287: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]));
288: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
289: }
290: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
291: PetscCall(PetscViewerPopFormat(viewer));
292: ctx->oldnconv = nconv;
293: }
294: PetscFunctionReturn(PETSC_SUCCESS);
295: }
297: PetscErrorCode SVDMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
298: {
299: SlepcConvMon mctx;
301: PetscFunctionBegin;
302: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
303: PetscCall(PetscNew(&mctx));
304: mctx->ctx = ctx;
305: (*vf)->data = (void*)mctx;
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }
309: PetscErrorCode SVDMonitorConvergedDestroy(PetscViewerAndFormat **vf)
310: {
311: PetscFunctionBegin;
312: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
313: PetscCall(PetscFree((*vf)->data));
314: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
315: PetscCall(PetscDrawLGDestroy(&(*vf)->lg));
316: PetscCall(PetscFree(*vf));
317: PetscFunctionReturn(PETSC_SUCCESS);
318: }
320: /*@C
321: SVDMonitorFirstDrawLG - Plots the error estimate of the first unconverged
322: approximation at each iteration of the singular value solver.
324: Collective
326: Input Parameters:
327: + svd - singular value solver context
328: . its - iteration number
329: . nconv - number of converged singular triplets so far
330: . sigma - singular values
331: . errest - error estimates
332: . nest - number of error estimates to display
333: - vf - viewer and format for monitoring
335: Options Database Key:
336: . -svd_monitor draw::draw_lg - activates SVDMonitorFirstDrawLG()
338: Level: intermediate
340: .seealso: SVDMonitorSet()
341: @*/
342: PetscErrorCode SVDMonitorFirstDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
343: {
344: PetscViewer viewer = vf->viewer;
345: PetscDrawLG lg = vf->lg;
346: PetscReal x,y;
348: PetscFunctionBegin;
352: PetscCall(PetscViewerPushFormat(viewer,vf->format));
353: if (its==1) {
354: PetscCall(PetscDrawLGReset(lg));
355: PetscCall(PetscDrawLGSetDimension(lg,1));
356: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
357: }
358: if (nconv<nest) {
359: x = (PetscReal)its;
360: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
361: else y = 0.0;
362: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
363: if (its <= 20 || !(its % 5) || svd->reason) {
364: PetscCall(PetscDrawLGDraw(lg));
365: PetscCall(PetscDrawLGSave(lg));
366: }
367: }
368: PetscCall(PetscViewerPopFormat(viewer));
369: PetscFunctionReturn(PETSC_SUCCESS);
370: }
372: /*@C
373: SVDMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
375: Collective
377: Input Parameters:
378: + viewer - the viewer
379: . format - the viewer format
380: - ctx - an optional user context
382: Output Parameter:
383: . vf - the viewer and format context
385: Level: intermediate
387: .seealso: SVDMonitorSet()
388: @*/
389: PetscErrorCode SVDMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
390: {
391: PetscFunctionBegin;
392: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
393: (*vf)->data = ctx;
394: PetscCall(SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
395: PetscFunctionReturn(PETSC_SUCCESS);
396: }
398: /*@C
399: SVDMonitorAllDrawLG - Plots the error estimate of all unconverged
400: approximations at each iteration of the singular value solver.
402: Collective
404: Input Parameters:
405: + svd - singular value solver context
406: . its - iteration number
407: . nconv - number of converged singular triplets so far
408: . sigma - singular values
409: . errest - error estimates
410: . nest - number of error estimates to display
411: - vf - viewer and format for monitoring
413: Options Database Key:
414: . -svd_monitor_all draw::draw_lg - activates SVDMonitorAllDrawLG()
416: Level: intermediate
418: .seealso: SVDMonitorSet()
419: @*/
420: PetscErrorCode SVDMonitorAllDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
421: {
422: PetscViewer viewer = vf->viewer;
423: PetscDrawLG lg = vf->lg;
424: PetscInt i,n = PetscMin(svd->nsv,255);
425: PetscReal *x,*y;
427: PetscFunctionBegin;
431: PetscCall(PetscViewerPushFormat(viewer,vf->format));
432: if (its==1) {
433: PetscCall(PetscDrawLGReset(lg));
434: PetscCall(PetscDrawLGSetDimension(lg,n));
435: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
436: }
437: PetscCall(PetscMalloc2(n,&x,n,&y));
438: for (i=0;i<n;i++) {
439: x[i] = (PetscReal)its;
440: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
441: else y[i] = 0.0;
442: }
443: PetscCall(PetscDrawLGAddPoint(lg,x,y));
444: if (its <= 20 || !(its % 5) || svd->reason) {
445: PetscCall(PetscDrawLGDraw(lg));
446: PetscCall(PetscDrawLGSave(lg));
447: }
448: PetscCall(PetscFree2(x,y));
449: PetscCall(PetscViewerPopFormat(viewer));
450: PetscFunctionReturn(PETSC_SUCCESS);
451: }
453: /*@C
454: SVDMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
456: Collective
458: Input Parameters:
459: + viewer - the viewer
460: . format - the viewer format
461: - ctx - an optional user context
463: Output Parameter:
464: . vf - the viewer and format context
466: Level: intermediate
468: .seealso: SVDMonitorSet()
469: @*/
470: PetscErrorCode SVDMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
471: {
472: PetscFunctionBegin;
473: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
474: (*vf)->data = ctx;
475: PetscCall(SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
476: PetscFunctionReturn(PETSC_SUCCESS);
477: }
479: /*@C
480: SVDMonitorConvergedDrawLG - Plots the number of converged eigenvalues
481: at each iteration of the singular value solver.
483: Collective
485: Input Parameters:
486: + svd - singular value solver context
487: . its - iteration number
488: . nconv - number of converged singular triplets so far
489: . sigma - singular values
490: . errest - error estimates
491: . nest - number of error estimates to display
492: - vf - viewer and format for monitoring
494: Options Database Key:
495: . -svd_monitor_conv draw::draw_lg - activates SVDMonitorConvergedDrawLG()
497: Level: intermediate
499: .seealso: SVDMonitorSet()
500: @*/
501: PetscErrorCode SVDMonitorConvergedDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
502: {
503: PetscViewer viewer = vf->viewer;
504: PetscDrawLG lg = vf->lg;
505: PetscReal x,y;
507: PetscFunctionBegin;
511: PetscCall(PetscViewerPushFormat(viewer,vf->format));
512: if (its==1) {
513: PetscCall(PetscDrawLGReset(lg));
514: PetscCall(PetscDrawLGSetDimension(lg,1));
515: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,svd->nsv));
516: }
517: x = (PetscReal)its;
518: y = (PetscReal)svd->nconv;
519: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
520: if (its <= 20 || !(its % 5) || svd->reason) {
521: PetscCall(PetscDrawLGDraw(lg));
522: PetscCall(PetscDrawLGSave(lg));
523: }
524: PetscCall(PetscViewerPopFormat(viewer));
525: PetscFunctionReturn(PETSC_SUCCESS);
526: }
528: /*@C
529: SVDMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
531: Collective
533: Input Parameters:
534: + viewer - the viewer
535: . format - the viewer format
536: - ctx - an optional user context
538: Output Parameter:
539: . vf - the viewer and format context
541: Level: intermediate
543: .seealso: SVDMonitorSet()
544: @*/
545: PetscErrorCode SVDMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
546: {
547: SlepcConvMon mctx;
549: PetscFunctionBegin;
550: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
551: PetscCall(PetscNew(&mctx));
552: mctx->ctx = ctx;
553: (*vf)->data = (void*)mctx;
554: PetscCall(SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
555: PetscFunctionReturn(PETSC_SUCCESS);
556: }
558: /*@C
559: SVDMonitorConditioning - Print the condition number at each iteration of the singular
560: value solver.
562: Collective
564: Input Parameters:
565: + svd - singular value solver context
566: . its - iteration number
567: . nconv - (unused) number of converged singular triplets so far
568: . sigma - (unused) singular values
569: . errest - (unused) error estimates
570: . nest - (unused) number of error estimates to display
571: - vf - viewer and format for monitoring
573: Options Database Key:
574: . -svd_monitor_conditioning - activates SVDMonitorConditioning()
576: Note:
577: Works only for solvers that use a DS of type GSVD. The printed information corresponds
578: to the maximum of the condition number of the two generated bidiagonal matrices.
580: Level: intermediate
582: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorFirst(), SVDMonitorConverged()
583: @*/
584: PetscErrorCode SVDMonitorConditioning(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
585: {
586: PetscViewer viewer = vf->viewer;
587: PetscBool isgsvd;
588: PetscReal cond;
590: PetscFunctionBegin;
593: PetscCall(PetscObjectTypeCompare((PetscObject)svd->ds,DSGSVD,&isgsvd));
594: if (!isgsvd) PetscFunctionReturn(PETSC_SUCCESS);
595: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Condition number of bidiagonal matrices for %s solve.\n",((PetscObject)svd)->prefix));
596: PetscCall(PetscViewerPushFormat(viewer,vf->format));
597: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
598: PetscCall(DSCond(svd->ds,&cond));
599: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD condition number = %g\n",its,(double)cond));
600: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
601: PetscCall(PetscViewerPopFormat(viewer));
602: PetscFunctionReturn(PETSC_SUCCESS);
603: }