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 - the singular value solver context
37: . monitor - pointer to function (if this is `NULL`, it turns off monitoring),
38: see `SVDMonitorFn`
39: . mctx - [optional] context for private data for the monitor routine
40: (use `NULL` if no context is desired)
41: - monitordestroy - [optional] routine that frees monitor context (may be `NULL`),
42: see `PetscCtxDestroyFn` for the calling sequence
44: Options Database Keys:
45: + -svd_monitor - print only the first error estimate
46: . -svd_monitor_all - print error estimates at each iteration
47: . -svd_monitor_conv - print the singular value approximations only when
48: convergence has been reached
49: . -svd_monitor_conditioning - print the condition number when available
50: . -svd_monitor draw::draw_lg - sets line graph monitor for the first unconverged
51: approximate singular value
52: . -svd_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
53: approximate singular values
54: . -svd_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
55: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
56: a code by calls to `SVDMonitorSet()`, but does not cancel
57: those set via the options database.
59: Notes:
60: The options database option `-svd_monitor` and related options are the easiest way
61: to turn on `SVD` iteration monitoring.
63: `SVDMonitorRegister()` provides a way to associate an options database key with `SVD`
64: monitor function.
66: Several different monitoring routines may be set by calling `SVDMonitorSet()` multiple
67: times; all will be called in the order in which they were set.
69: Fortran Note:
70: Only a single monitor function can be set for each `SVD` object.
72: Level: intermediate
74: .seealso: [](ch:svd), `SVDMonitorFirst()`, `SVDMonitorAll()`, `SVDMonitorConverged()`, `SVDMonitorConditioning()`, `SVDMonitorFirstDrawLG()`, `SVDMonitorAllDrawLG()`, `SVDMonitorConvergedDrawLG()`, `SVDMonitorCancel()`
75: @*/
76: PetscErrorCode SVDMonitorSet(SVD svd,SVDMonitorFn *monitor,void *mctx,PetscCtxDestroyFn *monitordestroy)
77: {
78: PetscInt i;
79: PetscBool identical;
81: PetscFunctionBegin;
83: for (i=0;i<svd->numbermonitors;i++) {
84: PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))(PetscVoidFn*)monitor,mctx,monitordestroy,(PetscErrorCode (*)(void))(PetscVoidFn*)svd->monitor[i],svd->monitorcontext[i],svd->monitordestroy[i],&identical));
85: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
86: }
87: PetscCheck(svd->numbermonitors<MAXSVDMONITORS,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
88: svd->monitor[svd->numbermonitors] = monitor;
89: svd->monitorcontext[svd->numbermonitors] = mctx;
90: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
91: PetscFunctionReturn(PETSC_SUCCESS);
92: }
94: /*@
95: SVDMonitorCancel - Clears all monitors for an `SVD` object.
97: Logically Collective
99: Input Parameter:
100: . svd - the singular value solver context
102: Options Database Key:
103: . -svd_monitor_cancel - cancels all monitors that have been hardwired into a code by calls to
104: `SVDMonitorSet()`, but does not cancel those set via the options database.
106: Level: intermediate
108: .seealso: [](ch:svd), `SVDMonitorSet()`
109: @*/
110: PetscErrorCode SVDMonitorCancel(SVD svd)
111: {
112: PetscInt i;
114: PetscFunctionBegin;
116: for (i=0; i<svd->numbermonitors; i++) {
117: if (svd->monitordestroy[i]) PetscCall((*svd->monitordestroy[i])(&svd->monitorcontext[i]));
118: }
119: svd->numbermonitors = 0;
120: PetscFunctionReturn(PETSC_SUCCESS);
121: }
123: /*@C
124: SVDGetMonitorContext - Gets the monitor context, as set by
125: `SVDMonitorSet()` for the FIRST monitor only.
127: Not Collective
129: Input Parameter:
130: . svd - the singular value solver context
132: Output Parameter:
133: . ctx - monitor context
135: Level: intermediate
137: .seealso: [](ch:svd), `SVDMonitorSet()`
138: @*/
139: PetscErrorCode SVDGetMonitorContext(SVD svd,void *ctx)
140: {
141: PetscFunctionBegin;
143: *(void**)ctx = svd->monitorcontext[0];
144: PetscFunctionReturn(PETSC_SUCCESS);
145: }
147: /*@C
148: SVDMonitorFirst - Print the first unconverged approximate value and
149: error estimate at each iteration of the singular value solver.
151: Collective
153: Input Parameters:
154: + svd - the singular value solver context
155: . its - iteration number
156: . nconv - number of converged singular triplets so far
157: . sigma - singular values
158: . errest - error estimates
159: . nest - number of error estimates to display
160: - vf - viewer and format for monitoring
162: Options Database Key:
163: . -svd_monitor - activates `SVDMonitorFirst()`
165: Note:
166: This is not called directly by users, rather one calls `SVDMonitorSet()`, with this
167: function as an argument, to cause the monitor to be used during the `SVD` solve.
169: Level: intermediate
171: .seealso: [](ch:svd), `SVDMonitorSet()`, `SVDMonitorAll()`, `SVDMonitorConditioning()`, `SVDMonitorConverged()`
172: @*/
173: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal sigma[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
174: {
175: PetscViewer viewer = vf->viewer;
177: PetscFunctionBegin;
180: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
181: if (nconv<nest) {
182: PetscCall(PetscViewerPushFormat(viewer,vf->format));
183: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
184: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
185: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
186: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]));
187: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
188: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
189: PetscCall(PetscViewerPopFormat(viewer));
190: }
191: PetscFunctionReturn(PETSC_SUCCESS);
192: }
194: /*@C
195: SVDMonitorAll - Print the current approximate values and
196: error estimates at each iteration of the singular value solver.
198: Collective
200: Input Parameters:
201: + svd - the singular value solver context
202: . its - iteration number
203: . nconv - number of converged singular triplets so far
204: . sigma - singular values
205: . errest - error estimates
206: . nest - number of error estimates to display
207: - vf - viewer and format for monitoring
209: Options Database Key:
210: . -svd_monitor_all - activates `SVDMonitorAll()`
212: Note:
213: This is not called directly by users, rather one calls `SVDMonitorSet()`, with this
214: function as an argument, to cause the monitor to be used during the `SVD` solve.
216: Level: intermediate
218: .seealso: [](ch:svd), `SVDMonitorSet()`, `SVDMonitorFirst()`, `SVDMonitorConditioning()`, `SVDMonitorConverged()`
219: @*/
220: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal sigma[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
221: {
222: PetscInt i;
223: PetscViewer viewer = vf->viewer;
225: PetscFunctionBegin;
228: PetscCall(PetscViewerPushFormat(viewer,vf->format));
229: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
230: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix));
231: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
232: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
233: for (i=0;i<nest;i++) PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]));
234: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
235: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
236: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
237: PetscCall(PetscViewerPopFormat(viewer));
238: PetscFunctionReturn(PETSC_SUCCESS);
239: }
241: /*@C
242: SVDMonitorConverged - Print the approximate values and
243: error estimates as they converge.
245: Collective
247: Input Parameters:
248: + svd - the singular value solver context
249: . its - iteration number
250: . nconv - number of converged singular triplets so far
251: . sigma - singular values
252: . errest - error estimates
253: . nest - number of error estimates to display
254: - vf - viewer and format for monitoring
256: Options Database Key:
257: . -svd_monitor_conv - activates `SVDMonitorConverged()`
259: Note:
260: This is not called directly by users, rather one calls `SVDMonitorSet()`, with this
261: function as an argument, to cause the monitor to be used during the `SVD` solve.
263: Level: intermediate
265: .seealso: [](ch:svd), `SVDMonitorSet()`, `SVDMonitorFirst()`, `SVDMonitorConditioning()`, `SVDMonitorAll()`
266: @*/
267: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal sigma[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
268: {
269: PetscInt i;
270: PetscViewer viewer = vf->viewer;
271: SlepcConvMon ctx;
273: PetscFunctionBegin;
276: ctx = (SlepcConvMon)vf->data;
277: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)svd)->prefix));
278: if (its==1) ctx->oldnconv = 0;
279: if (ctx->oldnconv!=nconv) {
280: PetscCall(PetscViewerPushFormat(viewer,vf->format));
281: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
282: for (i=ctx->oldnconv;i<nconv;i++) {
283: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD converged value (error) #%" PetscInt_FMT,its,i));
284: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
285: PetscCall(PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]));
286: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
287: }
288: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
289: PetscCall(PetscViewerPopFormat(viewer));
290: ctx->oldnconv = nconv;
291: }
292: PetscFunctionReturn(PETSC_SUCCESS);
293: }
295: PetscErrorCode SVDMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
296: {
297: SlepcConvMon mctx;
299: PetscFunctionBegin;
300: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
301: PetscCall(PetscNew(&mctx));
302: mctx->ctx = ctx;
303: (*vf)->data = (void*)mctx;
304: PetscFunctionReturn(PETSC_SUCCESS);
305: }
307: PetscErrorCode SVDMonitorConvergedDestroy(PetscViewerAndFormat **vf)
308: {
309: PetscFunctionBegin;
310: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
311: PetscCall(PetscFree((*vf)->data));
312: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
313: PetscCall(PetscFree(*vf));
314: PetscFunctionReturn(PETSC_SUCCESS);
315: }
317: /*@C
318: SVDMonitorFirstDrawLG - Plots the error estimate of the first unconverged
319: approximation at each iteration of the singular value solver.
321: Collective
323: Input Parameters:
324: + svd - the singular value solver context
325: . its - iteration number
326: . nconv - number of converged singular triplets so far
327: . sigma - singular values
328: . errest - error estimates
329: . nest - number of error estimates to display
330: - vf - viewer and format for monitoring
332: Options Database Key:
333: . -svd_monitor draw::draw_lg - activates `SVDMonitorFirstDrawLG()`
335: Notes:
336: This is not called directly by users, rather one calls `SVDMonitorSet()`, with this
337: function as an argument, to cause the monitor to be used during the `SVD` solve.
339: Call `SVDMonitorFirstDrawLGCreate()` to create the context used with this monitor.
341: Level: intermediate
343: .seealso: [](ch:svd), `SVDMonitorSet()`, `SVDMonitorFirstDrawLGCreate()`
344: @*/
345: PetscErrorCode SVDMonitorFirstDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal sigma[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
346: {
347: PetscViewer viewer = vf->viewer;
348: PetscDrawLG lg;
349: PetscReal x,y;
351: PetscFunctionBegin;
354: PetscCall(PetscViewerPushFormat(viewer,vf->format));
355: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
356: if (its==1) {
357: PetscCall(PetscDrawLGReset(lg));
358: PetscCall(PetscDrawLGSetDimension(lg,1));
359: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
360: }
361: if (nconv<nest) {
362: x = (PetscReal)its;
363: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
364: else y = 0.0;
365: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
366: if (its <= 20 || !(its % 5) || svd->reason) {
367: PetscCall(PetscDrawLGDraw(lg));
368: PetscCall(PetscDrawLGSave(lg));
369: }
370: }
371: PetscCall(PetscViewerPopFormat(viewer));
372: PetscFunctionReturn(PETSC_SUCCESS);
373: }
375: /*@C
376: SVDMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
378: Collective
380: Input Parameters:
381: + viewer - the viewer
382: . format - the viewer format
383: - ctx - an optional user context
385: Output Parameter:
386: . vf - the viewer and format context
388: Level: intermediate
390: .seealso: [](ch:svd), `SVDMonitorSet()`
391: @*/
392: PetscErrorCode SVDMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
393: {
394: PetscFunctionBegin;
395: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
396: (*vf)->data = ctx;
397: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
398: PetscFunctionReturn(PETSC_SUCCESS);
399: }
401: /*@C
402: SVDMonitorAllDrawLG - Plots the error estimate of all unconverged
403: approximations at each iteration of the singular value solver.
405: Collective
407: Input Parameters:
408: + svd - the singular value solver context
409: . its - iteration number
410: . nconv - number of converged singular triplets so far
411: . sigma - singular values
412: . errest - error estimates
413: . nest - number of error estimates to display
414: - vf - viewer and format for monitoring
416: Options Database Key:
417: . -svd_monitor_all draw::draw_lg - activates `SVDMonitorAllDrawLG()`
419: Notes:
420: This is not called directly by users, rather one calls `SVDMonitorSet()`, with this
421: function as an argument, to cause the monitor to be used during the `SVD` solve.
423: Call `SVDMonitorAllDrawLGCreate()` to create the context used with this monitor.
425: Level: intermediate
427: .seealso: [](ch:svd), `SVDMonitorSet()`, `SVDMonitorAllDrawLGCreate()`
428: @*/
429: PetscErrorCode SVDMonitorAllDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal sigma[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
430: {
431: PetscViewer viewer = vf->viewer;
432: PetscDrawLG lg;
433: PetscInt i,n = PetscMin(svd->nsv,255);
434: PetscReal *x,*y;
436: PetscFunctionBegin;
439: PetscCall(PetscViewerPushFormat(viewer,vf->format));
440: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
441: if (its==1) {
442: PetscCall(PetscDrawLGReset(lg));
443: PetscCall(PetscDrawLGSetDimension(lg,n));
444: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0));
445: }
446: PetscCall(PetscMalloc2(n,&x,n,&y));
447: for (i=0;i<n;i++) {
448: x[i] = (PetscReal)its;
449: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
450: else y[i] = 0.0;
451: }
452: PetscCall(PetscDrawLGAddPoint(lg,x,y));
453: if (its <= 20 || !(its % 5) || svd->reason) {
454: PetscCall(PetscDrawLGDraw(lg));
455: PetscCall(PetscDrawLGSave(lg));
456: }
457: PetscCall(PetscFree2(x,y));
458: PetscCall(PetscViewerPopFormat(viewer));
459: PetscFunctionReturn(PETSC_SUCCESS);
460: }
462: /*@C
463: SVDMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
465: Collective
467: Input Parameters:
468: + viewer - the viewer
469: . format - the viewer format
470: - ctx - an optional user context
472: Output Parameter:
473: . vf - the viewer and format context
475: Level: intermediate
477: .seealso: [](ch:svd), `SVDMonitorSet()`
478: @*/
479: PetscErrorCode SVDMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
480: {
481: PetscFunctionBegin;
482: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
483: (*vf)->data = ctx;
484: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
485: PetscFunctionReturn(PETSC_SUCCESS);
486: }
488: /*@C
489: SVDMonitorConvergedDrawLG - Plots the number of converged eigenvalues
490: at each iteration of the singular value solver.
492: Collective
494: Input Parameters:
495: + svd - the singular value solver context
496: . its - iteration number
497: . nconv - number of converged singular triplets so far
498: . sigma - singular values
499: . errest - error estimates
500: . nest - number of error estimates to display
501: - vf - viewer and format for monitoring
503: Options Database Key:
504: . -svd_monitor_conv draw::draw_lg - activates `SVDMonitorConvergedDrawLG()`
506: Notes:
507: This is not called directly by users, rather one calls `SVDMonitorSet()`, with this
508: function as an argument, to cause the monitor to be used during the `SVD` solve.
510: Call `SVDMonitorConvergedDrawLGCreate()` to create the context used with this monitor.
512: Level: intermediate
514: .seealso: [](ch:svd), `SVDMonitorSet()`, `SVDMonitorConvergedDrawLGCreate()`
515: @*/
516: PetscErrorCode SVDMonitorConvergedDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal sigma[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
517: {
518: PetscViewer viewer = vf->viewer;
519: PetscDrawLG lg;
520: PetscReal x,y;
522: PetscFunctionBegin;
525: PetscCall(PetscViewerPushFormat(viewer,vf->format));
526: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
527: if (its==1) {
528: PetscCall(PetscDrawLGReset(lg));
529: PetscCall(PetscDrawLGSetDimension(lg,1));
530: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,svd->nsv));
531: }
532: x = (PetscReal)its;
533: y = (PetscReal)svd->nconv;
534: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
535: if (its <= 20 || !(its % 5) || svd->reason) {
536: PetscCall(PetscDrawLGDraw(lg));
537: PetscCall(PetscDrawLGSave(lg));
538: }
539: PetscCall(PetscViewerPopFormat(viewer));
540: PetscFunctionReturn(PETSC_SUCCESS);
541: }
543: /*@C
544: SVDMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
546: Collective
548: Input Parameters:
549: + viewer - the viewer
550: . format - the viewer format
551: - ctx - an optional user context
553: Output Parameter:
554: . vf - the viewer and format context
556: Level: intermediate
558: .seealso: [](ch:svd), `SVDMonitorSet()`
559: @*/
560: PetscErrorCode SVDMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
561: {
562: SlepcConvMon mctx;
564: PetscFunctionBegin;
565: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
566: PetscCall(PetscNew(&mctx));
567: mctx->ctx = ctx;
568: (*vf)->data = (void*)mctx;
569: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
570: PetscFunctionReturn(PETSC_SUCCESS);
571: }
573: /*@C
574: SVDMonitorConditioning - Print the condition number at each iteration of the singular
575: value solver.
577: Collective
579: Input Parameters:
580: + svd - the singular value solver context
581: . its - iteration number
582: . nconv - (unused) number of converged singular triplets so far
583: . sigma - (unused) singular values
584: . errest - (unused) error estimates
585: . nest - (unused) number of error estimates to display
586: - vf - viewer and format for monitoring
588: Options Database Key:
589: . -svd_monitor_conditioning - activates `SVDMonitorConditioning()`
591: Notes:
592: This is not called directly by users, rather one calls `SVDMonitorSet()`, with this
593: function as an argument, to cause the monitor to be used during the `SVD` solve.
595: Works only for solvers that use a `DS` of type `DSGSVD`. The printed information corresponds
596: to the maximum of the condition number of the two generated bidiagonal matrices.
598: Level: intermediate
600: .seealso: [](ch:svd), `SVDMonitorSet()`, `SVDMonitorAll()`, `SVDMonitorFirst()`, `SVDMonitorConverged()`
601: @*/
602: PetscErrorCode SVDMonitorConditioning(SVD svd,PetscInt its,PetscInt nconv,PetscReal sigma[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
603: {
604: PetscViewer viewer = vf->viewer;
605: PetscBool isgsvd;
606: PetscReal cond;
608: PetscFunctionBegin;
611: PetscCall(PetscObjectTypeCompare((PetscObject)svd->ds,DSGSVD,&isgsvd));
612: if (!isgsvd) PetscFunctionReturn(PETSC_SUCCESS);
613: if (its==1 && ((PetscObject)svd)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Condition number of bidiagonal matrices for %s solve.\n",((PetscObject)svd)->prefix));
614: PetscCall(PetscViewerPushFormat(viewer,vf->format));
615: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel));
616: PetscCall(DSCond(svd->ds,&cond));
617: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD condition number = %g\n",its,(double)cond));
618: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel));
619: PetscCall(PetscViewerPopFormat(viewer));
620: PetscFunctionReturn(PETSC_SUCCESS);
621: }