Actual source code: epsmon.c
slepc-3.16.1 2021-11-17
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: EPS routines related to monitors
12: */
14: #include <slepc/private/epsimpl.h>
15: #include <petscdraw.h>
17: PetscErrorCode EPSMonitorLGCreate(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 EPSMonitor(EPS eps,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
41: {
43: PetscInt i,n = eps->numbermonitors;
46: for (i=0;i<n;i++) {
47: (*eps->monitor[i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[i]);
48: }
49: return(0);
50: }
52: /*@C
53: EPSMonitorSet - Sets an ADDITIONAL function to be called at every
54: iteration to monitor the error estimates for each requested eigenpair.
56: Logically Collective on eps
58: Input Parameters:
59: + eps - eigensolver context obtained from EPSCreate()
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(EPS eps,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)
68: + eps - eigensolver context obtained from EPSCreate()
69: . its - iteration number
70: . nconv - number of converged eigenpairs
71: . eigr - real part of the eigenvalues
72: . eigi - imaginary part of the eigenvalues
73: . errest - relative error estimates for each eigenpair
74: . nest - number of error estimates
75: - mctx - optional monitoring context, as set by EPSMonitorSet()
77: Options Database Keys:
78: + -eps_monitor - print only the first error estimate
79: . -eps_monitor_all - print error estimates at each iteration
80: . -eps_monitor_conv - print the eigenvalue approximations only when
81: convergence has been reached
82: . -eps_monitor draw::draw_lg - sets line graph monitor for the first unconverged
83: approximate eigenvalue
84: . -eps_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
85: approximate eigenvalues
86: . -eps_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
87: - -eps_monitor_cancel - cancels all monitors that have been hardwired into
88: a code by calls to EPSMonitorSet(), but does not cancel those set via
89: the options database.
91: Notes:
92: Several different monitoring routines may be set by calling
93: EPSMonitorSet() multiple times; all will be called in the
94: order in which they were set.
96: Level: intermediate
98: .seealso: EPSMonitorFirst(), EPSMonitorAll(), EPSMonitorCancel()
99: @*/
100: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
101: {
104: if (eps->numbermonitors >= MAXEPSMONITORS) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
105: eps->monitor[eps->numbermonitors] = monitor;
106: eps->monitorcontext[eps->numbermonitors] = (void*)mctx;
107: eps->monitordestroy[eps->numbermonitors++] = monitordestroy;
108: return(0);
109: }
111: /*@
112: EPSMonitorCancel - Clears all monitors for an EPS object.
114: Logically Collective on eps
116: Input Parameters:
117: . eps - eigensolver context obtained from EPSCreate()
119: Options Database Key:
120: . -eps_monitor_cancel - Cancels all monitors that have been hardwired
121: into a code by calls to EPSMonitorSet(),
122: but does not cancel those set via the options database.
124: Level: intermediate
126: .seealso: EPSMonitorSet()
127: @*/
128: PetscErrorCode EPSMonitorCancel(EPS eps)
129: {
131: PetscInt i;
135: for (i=0; i<eps->numbermonitors; i++) {
136: if (eps->monitordestroy[i]) {
137: (*eps->monitordestroy[i])(&eps->monitorcontext[i]);
138: }
139: }
140: eps->numbermonitors = 0;
141: return(0);
142: }
144: /*@C
145: EPSGetMonitorContext - Gets the monitor context, as set by
146: EPSMonitorSet() for the FIRST monitor only.
148: Not Collective
150: Input Parameter:
151: . eps - eigensolver context obtained from EPSCreate()
153: Output Parameter:
154: . ctx - monitor context
156: Level: intermediate
158: .seealso: EPSMonitorSet()
159: @*/
160: PetscErrorCode EPSGetMonitorContext(EPS eps,void *ctx)
161: {
164: *(void**)ctx = eps->monitorcontext[0];
165: return(0);
166: }
168: /*@C
169: EPSMonitorFirst - Print the first unconverged approximate value and
170: error estimate at each iteration of the eigensolver.
172: Collective on eps
174: Input Parameters:
175: + eps - eigensolver context
176: . its - iteration number
177: . nconv - number of converged eigenpairs so far
178: . eigr - real part of the eigenvalues
179: . eigi - imaginary part of the eigenvalues
180: . errest - error estimates
181: . nest - number of error estimates to display
182: - vf - viewer and format for monitoring
184: Options Database Key:
185: . -eps_monitor - activates EPSMonitorFirst()
187: Level: intermediate
189: .seealso: EPSMonitorSet(), EPSMonitorAll(), EPSMonitorConverged()
190: @*/
191: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
192: {
194: PetscScalar er,ei;
195: PetscViewer viewer = vf->viewer;
200: if (its==1 && ((PetscObject)eps)->prefix) {
201: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix);
202: }
203: if (nconv<nest) {
204: PetscViewerPushFormat(viewer,vf->format);
205: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
206: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D first unconverged value (error)",its,nconv);
207: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
208: er = eigr[nconv]; ei = eigi[nconv];
209: STBackTransform(eps->st,1,&er,&ei);
210: #if defined(PETSC_USE_COMPLEX)
211: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
212: #else
213: PetscViewerASCIIPrintf(viewer," %g",(double)er);
214: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
215: #endif
216: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
217: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
218: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
219: PetscViewerPopFormat(viewer);
220: }
221: return(0);
222: }
224: /*@C
225: EPSMonitorAll - Print the current approximate values and
226: error estimates at each iteration of the eigensolver.
228: Collective on eps
230: Input Parameters:
231: + eps - eigensolver context
232: . its - iteration number
233: . nconv - number of converged eigenpairs so far
234: . eigr - real part of the eigenvalues
235: . eigi - imaginary part of the eigenvalues
236: . errest - error estimates
237: . nest - number of error estimates to display
238: - vf - viewer and format for monitoring
240: Options Database Key:
241: . -eps_monitor_all - activates EPSMonitorAll()
243: Level: intermediate
245: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorConverged()
246: @*/
247: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
248: {
250: PetscInt i;
251: PetscScalar er,ei;
252: PetscViewer viewer = vf->viewer;
257: PetscViewerPushFormat(viewer,vf->format);
258: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
259: if (its==1 && ((PetscObject)eps)->prefix) {
260: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix);
261: }
262: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D Values (Errors)",its,nconv);
263: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
264: for (i=0;i<nest;i++) {
265: er = eigr[i]; ei = eigi[i];
266: STBackTransform(eps->st,1,&er,&ei);
267: #if defined(PETSC_USE_COMPLEX)
268: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
269: #else
270: PetscViewerASCIIPrintf(viewer," %g",(double)er);
271: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
272: #endif
273: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
274: }
275: PetscViewerASCIIPrintf(viewer,"\n");
276: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
277: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
278: PetscViewerPopFormat(viewer);
279: return(0);
280: }
282: /*@C
283: EPSMonitorConverged - Print the approximate values and
284: error estimates as they converge.
286: Collective on eps
288: Input Parameters:
289: + eps - eigensolver context
290: . its - iteration number
291: . nconv - number of converged eigenpairs so far
292: . eigr - real part of the eigenvalues
293: . eigi - imaginary part of the eigenvalues
294: . errest - error estimates
295: . nest - number of error estimates to display
296: - vf - viewer and format for monitoring
298: Options Database Key:
299: . -eps_monitor_conv - activates EPSMonitorConverged()
301: Level: intermediate
303: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorAll()
304: @*/
305: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
306: {
308: PetscInt i;
309: PetscScalar er,ei;
310: PetscViewer viewer = vf->viewer;
311: SlepcConvMon ctx;
316: ctx = (SlepcConvMon)vf->data;
317: if (its==1 && ((PetscObject)eps)->prefix) {
318: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)eps)->prefix);
319: }
320: if (its==1) ctx->oldnconv = 0;
321: if (ctx->oldnconv!=nconv) {
322: PetscViewerPushFormat(viewer,vf->format);
323: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
324: for (i=ctx->oldnconv;i<nconv;i++) {
325: PetscViewerASCIIPrintf(viewer,"%3D EPS converged value (error) #%D",its,i);
326: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
327: er = eigr[i]; ei = eigi[i];
328: STBackTransform(eps->st,1,&er,&ei);
329: #if defined(PETSC_USE_COMPLEX)
330: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
331: #else
332: PetscViewerASCIIPrintf(viewer," %g",(double)er);
333: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
334: #endif
335: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
336: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
337: }
338: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
339: PetscViewerPopFormat(viewer);
340: ctx->oldnconv = nconv;
341: }
342: return(0);
343: }
345: PetscErrorCode EPSMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
346: {
348: SlepcConvMon mctx;
351: PetscViewerAndFormatCreate(viewer,format,vf);
352: PetscNew(&mctx);
353: mctx->ctx = ctx;
354: (*vf)->data = (void*)mctx;
355: return(0);
356: }
358: PetscErrorCode EPSMonitorConvergedDestroy(PetscViewerAndFormat **vf)
359: {
363: if (!*vf) return(0);
364: PetscFree((*vf)->data);
365: PetscViewerDestroy(&(*vf)->viewer);
366: PetscDrawLGDestroy(&(*vf)->lg);
367: PetscFree(*vf);
368: return(0);
369: }
371: /*@C
372: EPSMonitorFirstDrawLG - Plots the error estimate of the first unconverged
373: approximation at each iteration of the eigensolver.
375: Collective on eps
377: Input Parameters:
378: + eps - eigensolver context
379: . its - iteration number
380: . its - iteration number
381: . nconv - number of converged eigenpairs so far
382: . eigr - real part of the eigenvalues
383: . eigi - imaginary part of the eigenvalues
384: . errest - error estimates
385: . nest - number of error estimates to display
386: - vf - viewer and format for monitoring
388: Options Database Key:
389: . -eps_monitor draw::draw_lg - activates EPSMonitorFirstDrawLG()
391: Level: intermediate
393: .seealso: EPSMonitorSet()
394: @*/
395: PetscErrorCode EPSMonitorFirstDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
396: {
398: PetscViewer viewer = vf->viewer;
399: PetscDrawLG lg = vf->lg;
400: PetscReal x,y;
406: PetscViewerPushFormat(viewer,vf->format);
407: if (its==1) {
408: PetscDrawLGReset(lg);
409: PetscDrawLGSetDimension(lg,1);
410: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0);
411: }
412: if (nconv<nest) {
413: x = (PetscReal)its;
414: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
415: else y = 0.0;
416: PetscDrawLGAddPoint(lg,&x,&y);
417: if (its <= 20 || !(its % 5) || eps->reason) {
418: PetscDrawLGDraw(lg);
419: PetscDrawLGSave(lg);
420: }
421: }
422: PetscViewerPopFormat(viewer);
423: return(0);
424: }
426: /*@C
427: EPSMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
429: Collective on viewer
431: Input Parameters:
432: + viewer - the viewer
433: . format - the viewer format
434: - ctx - an optional user context
436: Output Parameter:
437: . vf - the viewer and format context
439: Level: intermediate
441: .seealso: EPSMonitorSet()
442: @*/
443: PetscErrorCode EPSMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
444: {
448: PetscViewerAndFormatCreate(viewer,format,vf);
449: (*vf)->data = ctx;
450: EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
451: return(0);
452: }
454: /*@C
455: EPSMonitorAllDrawLG - Plots the error estimate of all unconverged
456: approximations at each iteration of the eigensolver.
458: Collective on eps
460: Input Parameters:
461: + eps - eigensolver context
462: . its - iteration number
463: . its - iteration number
464: . nconv - number of converged eigenpairs so far
465: . eigr - real part of the eigenvalues
466: . eigi - imaginary part of the eigenvalues
467: . errest - error estimates
468: . nest - number of error estimates to display
469: - vf - viewer and format for monitoring
471: Options Database Key:
472: . -eps_monitor_all draw::draw_lg - activates EPSMonitorAllDrawLG()
474: Level: intermediate
476: .seealso: EPSMonitorSet()
477: @*/
478: PetscErrorCode EPSMonitorAllDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
479: {
481: PetscViewer viewer = vf->viewer;
482: PetscDrawLG lg = vf->lg;
483: PetscInt i,n = PetscMin(eps->nev,255);
484: PetscReal *x,*y;
490: PetscViewerPushFormat(viewer,vf->format);
491: if (its==1) {
492: PetscDrawLGReset(lg);
493: PetscDrawLGSetDimension(lg,n);
494: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0);
495: }
496: PetscMalloc2(n,&x,n,&y);
497: for (i=0;i<n;i++) {
498: x[i] = (PetscReal)its;
499: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
500: else y[i] = 0.0;
501: }
502: PetscDrawLGAddPoint(lg,x,y);
503: if (its <= 20 || !(its % 5) || eps->reason) {
504: PetscDrawLGDraw(lg);
505: PetscDrawLGSave(lg);
506: }
507: PetscFree2(x,y);
508: PetscViewerPopFormat(viewer);
509: return(0);
510: }
512: /*@C
513: EPSMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
515: Collective on viewer
517: Input Parameters:
518: + viewer - the viewer
519: . format - the viewer format
520: - ctx - an optional user context
522: Output Parameter:
523: . vf - the viewer and format context
525: Level: intermediate
527: .seealso: EPSMonitorSet()
528: @*/
529: PetscErrorCode EPSMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
530: {
534: PetscViewerAndFormatCreate(viewer,format,vf);
535: (*vf)->data = ctx;
536: EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
537: return(0);
538: }
540: /*@C
541: EPSMonitorConvergedDrawLG - Plots the number of converged eigenvalues
542: at each iteration of the eigensolver.
544: Collective on eps
546: Input Parameters:
547: + eps - eigensolver context
548: . its - iteration number
549: . its - iteration number
550: . nconv - number of converged eigenpairs so far
551: . eigr - real part of the eigenvalues
552: . eigi - imaginary part of the eigenvalues
553: . errest - error estimates
554: . nest - number of error estimates to display
555: - vf - viewer and format for monitoring
557: Options Database Key:
558: . -eps_monitor_conv draw::draw_lg - activates EPSMonitorConvergedDrawLG()
560: Level: intermediate
562: .seealso: EPSMonitorSet()
563: @*/
564: PetscErrorCode EPSMonitorConvergedDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
565: {
566: PetscErrorCode ierr;
567: PetscViewer viewer = vf->viewer;
568: PetscDrawLG lg = vf->lg;
569: PetscReal x,y;
575: PetscViewerPushFormat(viewer,vf->format);
576: if (its==1) {
577: PetscDrawLGReset(lg);
578: PetscDrawLGSetDimension(lg,1);
579: PetscDrawLGSetLimits(lg,1,1,0,eps->nev);
580: }
581: x = (PetscReal)its;
582: y = (PetscReal)eps->nconv;
583: PetscDrawLGAddPoint(lg,&x,&y);
584: if (its <= 20 || !(its % 5) || eps->reason) {
585: PetscDrawLGDraw(lg);
586: PetscDrawLGSave(lg);
587: }
588: PetscViewerPopFormat(viewer);
589: return(0);
590: }
592: /*@C
593: EPSMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
595: Collective on viewer
597: Input Parameters:
598: + viewer - the viewer
599: . format - the viewer format
600: - ctx - an optional user context
602: Output Parameter:
603: . vf - the viewer and format context
605: Level: intermediate
607: .seealso: EPSMonitorSet()
608: @*/
609: PetscErrorCode EPSMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
610: {
612: SlepcConvMon mctx;
615: PetscViewerAndFormatCreate(viewer,format,vf);
616: PetscNew(&mctx);
617: mctx->ctx = ctx;
618: (*vf)->data = (void*)mctx;
619: EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
620: return(0);
621: }