Actual source code: svdmon.c
slepc-3.16.3 2022-04-11
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2021, 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;
25: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
26: PetscDrawSetFromOptions(draw);
27: PetscDrawLGCreate(draw,l,&lg);
28: if (names) { PetscDrawLGSetLegend(lg,names); }
29: PetscDrawLGSetFromOptions(lg);
30: PetscDrawLGGetAxis(lg,&axis);
31: PetscDrawAxisSetLabels(axis,"Convergence","Iteration",metric);
32: PetscDrawDestroy(&draw);
33: *lgctx = lg;
34: return(0);
35: }
37: /*
38: Runs the user provided monitor routines, if any.
39: */
40: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest)
41: {
43: PetscInt i,n = svd->numbermonitors;
46: for (i=0;i<n;i++) {
47: (*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]);
48: }
49: return(0);
50: }
52: /*@C
53: SVDMonitorSet - Sets an ADDITIONAL function to be called at every
54: iteration to monitor the error estimates for each requested singular triplet.
56: Logically Collective on svd
58: Input Parameters:
59: + svd - singular value solver context obtained from SVDCreate()
60: . monitor - pointer to function (if this is NULL, it turns off monitoring)
61: . mctx - [optional] context for private data for the
62: monitor routine (use NULL if no context is desired)
63: - monitordestroy - [optional] routine that frees monitor context (may be NULL)
65: Calling Sequence of monitor:
66: $ monitor(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *mctx)
68: + svd - singular value solver context obtained from SVDCreate()
69: . its - iteration number
70: . nconv - number of converged singular triplets
71: . sigma - singular values
72: . errest - relative error estimates for each singular triplet
73: . nest - number of error estimates
74: - mctx - optional monitoring context, as set by SVDMonitorSet()
76: Options Database Keys:
77: + -svd_monitor - print only the first error estimate
78: . -svd_monitor_all - print error estimates at each iteration
79: . -svd_monitor_conv - print the singular value approximations only when
80: convergence has been reached
81: . -svd_monitor draw::draw_lg - sets line graph monitor for the first unconverged
82: approximate singular value
83: . -svd_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
84: approximate singular values
85: . -svd_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
86: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
87: a code by calls to SVDMonitorSet(), but does not cancel those set via
88: the options database.
90: Notes:
91: Several different monitoring routines may be set by calling
92: SVDMonitorSet() multiple times; all will be called in the
93: order in which they were set.
95: Level: intermediate
97: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorCancel()
98: @*/
99: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
100: {
103: if (svd->numbermonitors >= MAXSVDMONITORS) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
104: svd->monitor[svd->numbermonitors] = monitor;
105: svd->monitorcontext[svd->numbermonitors] = (void*)mctx;
106: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
107: return(0);
108: }
110: /*@
111: SVDMonitorCancel - Clears all monitors for an SVD object.
113: Logically Collective on svd
115: Input Parameters:
116: . svd - singular value solver context obtained from SVDCreate()
118: Options Database Key:
119: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
120: into a code by calls to SVDMonitorSet(),
121: but does not cancel those set via the options database.
123: Level: intermediate
125: .seealso: SVDMonitorSet()
126: @*/
127: PetscErrorCode SVDMonitorCancel(SVD svd)
128: {
130: PetscInt i;
134: for (i=0; i<svd->numbermonitors; i++) {
135: if (svd->monitordestroy[i]) {
136: (*svd->monitordestroy[i])(&svd->monitorcontext[i]);
137: }
138: }
139: svd->numbermonitors = 0;
140: return(0);
141: }
143: /*@C
144: SVDGetMonitorContext - Gets the monitor context, as set by
145: SVDMonitorSet() for the FIRST monitor only.
147: Not Collective
149: Input Parameter:
150: . svd - singular value solver context obtained from SVDCreate()
152: Output Parameter:
153: . ctx - monitor context
155: Level: intermediate
157: .seealso: SVDMonitorSet()
158: @*/
159: PetscErrorCode SVDGetMonitorContext(SVD svd,void *ctx)
160: {
163: *(void**)ctx = svd->monitorcontext[0];
164: return(0);
165: }
167: /*@C
168: SVDMonitorFirst - Print the first unconverged approximate value and
169: error estimate at each iteration of the singular value solver.
171: Collective on svd
173: Input Parameters:
174: + svd - singular value solver context
175: . its - iteration number
176: . nconv - number of converged singular triplets so far
177: . sigma - singular values
178: . errest - error estimates
179: . nest - number of error estimates to display
180: - vf - viewer and format for monitoring
182: Options Database Key:
183: . -svd_monitor - activates SVDMonitorFirst()
185: Level: intermediate
187: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConverged()
188: @*/
189: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
190: {
192: PetscViewer viewer = vf->viewer;
197: if (its==1 && ((PetscObject)svd)->prefix) {
198: PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
199: }
200: if (nconv<nest) {
201: PetscViewerPushFormat(viewer,vf->format);
202: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
203: PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D first unconverged value (error)",its,nconv);
204: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
205: PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]);
206: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
207: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
208: PetscViewerPopFormat(viewer);
209: }
210: return(0);
211: }
213: /*@C
214: SVDMonitorAll - Print the current approximate values and
215: error estimates at each iteration of the singular value solver.
217: Collective on svd
219: Input Parameters:
220: + svd - singular value solver context
221: . its - iteration number
222: . nconv - number of converged singular triplets so far
223: . sigma - singular values
224: . errest - error estimates
225: . nest - number of error estimates to display
226: - vf - viewer and format for monitoring
228: Options Database Key:
229: . -svd_monitor_all - activates SVDMonitorAll()
231: Level: intermediate
233: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConverged()
234: @*/
235: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
236: {
238: PetscInt i;
239: PetscViewer viewer = vf->viewer;
244: PetscViewerPushFormat(viewer,vf->format);
245: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
246: if (its==1 && ((PetscObject)svd)->prefix) {
247: PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
248: }
249: PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D Values (Errors)",its,nconv);
250: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
251: for (i=0;i<nest;i++) {
252: PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]);
253: }
254: PetscViewerASCIIPrintf(viewer,"\n");
255: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
256: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
257: PetscViewerPopFormat(viewer);
258: return(0);
259: }
261: /*@C
262: SVDMonitorConverged - Print the approximate values and
263: error estimates as they converge.
265: Collective on svd
267: Input Parameters:
268: + svd - singular value solver context
269: . its - iteration number
270: . nconv - number of converged singular triplets so far
271: . sigma - singular values
272: . errest - error estimates
273: . nest - number of error estimates to display
274: - vf - viewer and format for monitoring
276: Options Database Key:
277: . -svd_monitor_conv - activates SVDMonitorConverged()
279: Level: intermediate
281: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorAll()
282: @*/
283: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
284: {
286: PetscInt i;
287: PetscViewer viewer = vf->viewer;
288: SlepcConvMon ctx;
293: ctx = (SlepcConvMon)vf->data;
294: if (its==1 && ((PetscObject)svd)->prefix) {
295: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)svd)->prefix);
296: }
297: if (its==1) ctx->oldnconv = 0;
298: if (ctx->oldnconv!=nconv) {
299: PetscViewerPushFormat(viewer,vf->format);
300: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
301: for (i=ctx->oldnconv;i<nconv;i++) {
302: PetscViewerASCIIPrintf(viewer,"%3D SVD converged value (error) #%D",its,i);
303: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
304: PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]);
305: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
306: }
307: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
308: PetscViewerPopFormat(viewer);
309: ctx->oldnconv = nconv;
310: }
311: return(0);
312: }
314: PetscErrorCode SVDMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
315: {
317: SlepcConvMon mctx;
320: PetscViewerAndFormatCreate(viewer,format,vf);
321: PetscNew(&mctx);
322: mctx->ctx = ctx;
323: (*vf)->data = (void*)mctx;
324: return(0);
325: }
327: PetscErrorCode SVDMonitorConvergedDestroy(PetscViewerAndFormat **vf)
328: {
332: if (!*vf) return(0);
333: PetscFree((*vf)->data);
334: PetscViewerDestroy(&(*vf)->viewer);
335: PetscDrawLGDestroy(&(*vf)->lg);
336: PetscFree(*vf);
337: return(0);
338: }
340: /*@C
341: SVDMonitorFirstDrawLG - Plots the error estimate of the first unconverged
342: approximation at each iteration of the singular value solver.
344: Collective on svd
346: Input Parameters:
347: + svd - singular value solver context
348: . its - iteration number
349: . its - iteration number
350: . nconv - number of converged singular triplets so far
351: . sigma - singular values
352: . errest - error estimates
353: . nest - number of error estimates to display
354: - vf - viewer and format for monitoring
356: Options Database Key:
357: . -svd_monitor draw::draw_lg - activates SVDMonitorFirstDrawLG()
359: Level: intermediate
361: .seealso: SVDMonitorSet()
362: @*/
363: PetscErrorCode SVDMonitorFirstDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
364: {
366: PetscViewer viewer = vf->viewer;
367: PetscDrawLG lg = vf->lg;
368: PetscReal x,y;
374: PetscViewerPushFormat(viewer,vf->format);
375: if (its==1) {
376: PetscDrawLGReset(lg);
377: PetscDrawLGSetDimension(lg,1);
378: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0);
379: }
380: if (nconv<nest) {
381: x = (PetscReal)its;
382: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
383: else y = 0.0;
384: PetscDrawLGAddPoint(lg,&x,&y);
385: if (its <= 20 || !(its % 5) || svd->reason) {
386: PetscDrawLGDraw(lg);
387: PetscDrawLGSave(lg);
388: }
389: }
390: PetscViewerPopFormat(viewer);
391: return(0);
392: }
394: /*@C
395: SVDMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
397: Collective on viewer
399: Input Parameters:
400: + viewer - the viewer
401: . format - the viewer format
402: - ctx - an optional user context
404: Output Parameter:
405: . vf - the viewer and format context
407: Level: intermediate
409: .seealso: SVDMonitorSet()
410: @*/
411: PetscErrorCode SVDMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
412: {
416: PetscViewerAndFormatCreate(viewer,format,vf);
417: (*vf)->data = ctx;
418: SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
419: return(0);
420: }
422: /*@C
423: SVDMonitorAllDrawLG - Plots the error estimate of all unconverged
424: approximations at each iteration of the singular value solver.
426: Collective on svd
428: Input Parameters:
429: + svd - singular value solver context
430: . its - iteration number
431: . its - iteration number
432: . nconv - number of converged singular triplets so far
433: . sigma - singular values
434: . errest - error estimates
435: . nest - number of error estimates to display
436: - vf - viewer and format for monitoring
438: Options Database Key:
439: . -svd_monitor_all draw::draw_lg - activates SVDMonitorAllDrawLG()
441: Level: intermediate
443: .seealso: SVDMonitorSet()
444: @*/
445: PetscErrorCode SVDMonitorAllDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
446: {
448: PetscViewer viewer = vf->viewer;
449: PetscDrawLG lg = vf->lg;
450: PetscInt i,n = PetscMin(svd->nsv,255);
451: PetscReal *x,*y;
457: PetscViewerPushFormat(viewer,vf->format);
458: if (its==1) {
459: PetscDrawLGReset(lg);
460: PetscDrawLGSetDimension(lg,n);
461: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0);
462: }
463: PetscMalloc2(n,&x,n,&y);
464: for (i=0;i<n;i++) {
465: x[i] = (PetscReal)its;
466: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
467: else y[i] = 0.0;
468: }
469: PetscDrawLGAddPoint(lg,x,y);
470: if (its <= 20 || !(its % 5) || svd->reason) {
471: PetscDrawLGDraw(lg);
472: PetscDrawLGSave(lg);
473: }
474: PetscFree2(x,y);
475: PetscViewerPopFormat(viewer);
476: return(0);
477: }
479: /*@C
480: SVDMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
482: Collective on viewer
484: Input Parameters:
485: + viewer - the viewer
486: . format - the viewer format
487: - ctx - an optional user context
489: Output Parameter:
490: . vf - the viewer and format context
492: Level: intermediate
494: .seealso: SVDMonitorSet()
495: @*/
496: PetscErrorCode SVDMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
497: {
501: PetscViewerAndFormatCreate(viewer,format,vf);
502: (*vf)->data = ctx;
503: SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
504: return(0);
505: }
507: /*@C
508: SVDMonitorConvergedDrawLG - Plots the number of converged eigenvalues
509: at each iteration of the singular value solver.
511: Collective on svd
513: Input Parameters:
514: + svd - singular value solver context
515: . its - iteration number
516: . its - iteration number
517: . nconv - number of converged singular triplets so far
518: . sigma - singular values
519: . errest - error estimates
520: . nest - number of error estimates to display
521: - vf - viewer and format for monitoring
523: Options Database Key:
524: . -svd_monitor_conv draw::draw_lg - activates SVDMonitorConvergedDrawLG()
526: Level: intermediate
528: .seealso: SVDMonitorSet()
529: @*/
530: PetscErrorCode SVDMonitorConvergedDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
531: {
532: PetscErrorCode ierr;
533: PetscViewer viewer = vf->viewer;
534: PetscDrawLG lg = vf->lg;
535: PetscReal x,y;
541: PetscViewerPushFormat(viewer,vf->format);
542: if (its==1) {
543: PetscDrawLGReset(lg);
544: PetscDrawLGSetDimension(lg,1);
545: PetscDrawLGSetLimits(lg,1,1,0,svd->nsv);
546: }
547: x = (PetscReal)its;
548: y = (PetscReal)svd->nconv;
549: PetscDrawLGAddPoint(lg,&x,&y);
550: if (its <= 20 || !(its % 5) || svd->reason) {
551: PetscDrawLGDraw(lg);
552: PetscDrawLGSave(lg);
553: }
554: PetscViewerPopFormat(viewer);
555: return(0);
556: }
558: /*@C
559: SVDMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
561: Collective on viewer
563: Input Parameters:
564: + viewer - the viewer
565: . format - the viewer format
566: - ctx - an optional user context
568: Output Parameter:
569: . vf - the viewer and format context
571: Level: intermediate
573: .seealso: SVDMonitorSet()
574: @*/
575: PetscErrorCode SVDMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
576: {
578: SlepcConvMon mctx;
581: PetscViewerAndFormatCreate(viewer,format,vf);
582: PetscNew(&mctx);
583: mctx->ctx = ctx;
584: (*vf)->data = (void*)mctx;
585: SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
586: return(0);
587: }