Actual source code: epsmon.c
slepc-3.15.0 2021-03-31
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: *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;
200: viewer = vf->viewer;
202: if (its==1 && ((PetscObject)eps)->prefix) {
203: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix);
204: }
205: if (nconv<nest) {
206: PetscViewerPushFormat(viewer,vf->format);
207: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
208: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D first unconverged value (error)",its,nconv);
209: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
210: er = eigr[nconv]; ei = eigi[nconv];
211: STBackTransform(eps->st,1,&er,&ei);
212: #if defined(PETSC_USE_COMPLEX)
213: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
214: #else
215: PetscViewerASCIIPrintf(viewer," %g",(double)er);
216: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
217: #endif
218: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
219: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
220: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
221: PetscViewerPopFormat(viewer);
222: }
223: return(0);
224: }
226: /*@C
227: EPSMonitorAll - Print the current approximate values and
228: error estimates at each iteration of the eigensolver.
230: Collective on eps
232: Input Parameters:
233: + eps - eigensolver context
234: . its - iteration number
235: . nconv - number of converged eigenpairs so far
236: . eigr - real part of the eigenvalues
237: . eigi - imaginary part of the eigenvalues
238: . errest - error estimates
239: . nest - number of error estimates to display
240: - vf - viewer and format for monitoring
242: Options Database Key:
243: . -eps_monitor_all - activates EPSMonitorAll()
245: Level: intermediate
247: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorConverged()
248: @*/
249: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
250: {
252: PetscInt i;
253: PetscScalar er,ei;
254: PetscViewer viewer;
259: viewer = vf->viewer;
261: PetscViewerPushFormat(viewer,vf->format);
262: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
263: if (its==1 && ((PetscObject)eps)->prefix) {
264: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix);
265: }
266: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D Values (Errors)",its,nconv);
267: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
268: for (i=0;i<nest;i++) {
269: er = eigr[i]; ei = eigi[i];
270: STBackTransform(eps->st,1,&er,&ei);
271: #if defined(PETSC_USE_COMPLEX)
272: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
273: #else
274: PetscViewerASCIIPrintf(viewer," %g",(double)er);
275: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
276: #endif
277: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
278: }
279: PetscViewerASCIIPrintf(viewer,"\n");
280: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
281: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
282: PetscViewerPopFormat(viewer);
283: return(0);
284: }
286: /*@C
287: EPSMonitorConverged - Print the approximate values and
288: error estimates as they converge.
290: Collective on eps
292: Input Parameters:
293: + eps - eigensolver context
294: . its - iteration number
295: . nconv - number of converged eigenpairs so far
296: . eigr - real part of the eigenvalues
297: . eigi - imaginary part of the eigenvalues
298: . errest - error estimates
299: . nest - number of error estimates to display
300: - vf - viewer and format for monitoring
302: Options Database Key:
303: . -eps_monitor_conv - activates EPSMonitorConverged()
305: Level: intermediate
307: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorAll()
308: @*/
309: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
310: {
312: PetscInt i;
313: PetscScalar er,ei;
314: PetscViewer viewer;
315: SlepcConvMon ctx;
320: viewer = vf->viewer;
322: ctx = (SlepcConvMon)vf->data;
323: if (its==1 && ((PetscObject)eps)->prefix) {
324: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)eps)->prefix);
325: }
326: if (its==1) ctx->oldnconv = 0;
327: if (ctx->oldnconv!=nconv) {
328: PetscViewerPushFormat(viewer,vf->format);
329: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
330: for (i=ctx->oldnconv;i<nconv;i++) {
331: PetscViewerASCIIPrintf(viewer,"%3D EPS converged value (error) #%D",its,i);
332: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
333: er = eigr[i]; ei = eigi[i];
334: STBackTransform(eps->st,1,&er,&ei);
335: #if defined(PETSC_USE_COMPLEX)
336: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
337: #else
338: PetscViewerASCIIPrintf(viewer," %g",(double)er);
339: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
340: #endif
341: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
342: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
343: }
344: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
345: PetscViewerPopFormat(viewer);
346: ctx->oldnconv = nconv;
347: }
348: return(0);
349: }
351: PetscErrorCode EPSMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
352: {
354: SlepcConvMon mctx;
357: PetscViewerAndFormatCreate(viewer,format,vf);
358: PetscNew(&mctx);
359: mctx->ctx = ctx;
360: (*vf)->data = (void*)mctx;
361: return(0);
362: }
364: PetscErrorCode EPSMonitorConvergedDestroy(PetscViewerAndFormat **vf)
365: {
369: if (!*vf) return(0);
370: PetscFree((*vf)->data);
371: PetscViewerDestroy(&(*vf)->viewer);
372: PetscDrawLGDestroy(&(*vf)->lg);
373: PetscFree(*vf);
374: return(0);
375: }
377: /*@C
378: EPSMonitorFirstDrawLG - Plots the error estimate of the first unconverged
379: approximation at each iteration of the eigensolver.
381: Collective on eps
383: Input Parameters:
384: + eps - eigensolver context
385: . its - iteration number
386: . its - iteration number
387: . nconv - number of converged eigenpairs so far
388: . eigr - real part of the eigenvalues
389: . eigi - imaginary part of the eigenvalues
390: . errest - error estimates
391: . nest - number of error estimates to display
392: - vf - viewer and format for monitoring
394: Options Database Key:
395: . -eps_monitor draw::draw_lg - activates EPSMonitorFirstDrawLG()
397: Level: intermediate
399: .seealso: EPSMonitorSet()
400: @*/
401: PetscErrorCode EPSMonitorFirstDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
402: {
404: PetscViewer viewer;
405: PetscDrawLG lg;
406: PetscReal x,y;
411: viewer = vf->viewer;
413: lg = vf->lg;
415: PetscViewerPushFormat(viewer,vf->format);
416: if (its==1) {
417: PetscDrawLGReset(lg);
418: PetscDrawLGSetDimension(lg,1);
419: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0);
420: }
421: if (nconv<nest) {
422: x = (PetscReal)its;
423: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
424: else y = 0.0;
425: PetscDrawLGAddPoint(lg,&x,&y);
426: if (its <= 20 || !(its % 5) || eps->reason) {
427: PetscDrawLGDraw(lg);
428: PetscDrawLGSave(lg);
429: }
430: }
431: PetscViewerPopFormat(viewer);
432: return(0);
433: }
435: /*@C
436: EPSMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
438: Collective on viewer
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: EPSMonitorSet()
451: @*/
452: PetscErrorCode EPSMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
453: {
457: PetscViewerAndFormatCreate(viewer,format,vf);
458: (*vf)->data = ctx;
459: EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
460: return(0);
461: }
463: /*@C
464: EPSMonitorAllDrawLG - Plots the error estimate of all unconverged
465: approximations at each iteration of the eigensolver.
467: Collective on eps
469: Input Parameters:
470: + eps - eigensolver context
471: . its - iteration number
472: . its - iteration number
473: . nconv - number of converged eigenpairs so far
474: . eigr - real part of the eigenvalues
475: . eigi - imaginary part of the eigenvalues
476: . errest - error estimates
477: . nest - number of error estimates to display
478: - vf - viewer and format for monitoring
480: Options Database Key:
481: . -eps_monitor_all draw::draw_lg - activates EPSMonitorAllDrawLG()
483: Level: intermediate
485: .seealso: EPSMonitorSet()
486: @*/
487: PetscErrorCode EPSMonitorAllDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
488: {
490: PetscViewer viewer;
491: PetscDrawLG lg;
492: PetscInt i,n = PetscMin(eps->nev,255);
493: PetscReal *x,*y;
498: viewer = vf->viewer;
500: lg = vf->lg;
502: PetscViewerPushFormat(viewer,vf->format);
503: if (its==1) {
504: PetscDrawLGReset(lg);
505: PetscDrawLGSetDimension(lg,n);
506: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0);
507: }
508: PetscMalloc2(n,&x,n,&y);
509: for (i=0;i<n;i++) {
510: x[i] = (PetscReal)its;
511: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
512: else y[i] = 0.0;
513: }
514: PetscDrawLGAddPoint(lg,x,y);
515: if (its <= 20 || !(its % 5) || eps->reason) {
516: PetscDrawLGDraw(lg);
517: PetscDrawLGSave(lg);
518: }
519: PetscFree2(x,y);
520: PetscViewerPopFormat(viewer);
521: return(0);
522: }
524: /*@C
525: EPSMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
527: Collective on viewer
529: Input Parameters:
530: + viewer - the viewer
531: . format - the viewer format
532: - ctx - an optional user context
534: Output Parameter:
535: . vf - the viewer and format context
537: Level: intermediate
539: .seealso: EPSMonitorSet()
540: @*/
541: PetscErrorCode EPSMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
542: {
546: PetscViewerAndFormatCreate(viewer,format,vf);
547: (*vf)->data = ctx;
548: EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
549: return(0);
550: }
552: /*@C
553: EPSMonitorConvergedDrawLG - Plots the number of converged eigenvalues
554: at each iteration of the eigensolver.
556: Collective on eps
558: Input Parameters:
559: + eps - eigensolver context
560: . its - iteration number
561: . its - iteration number
562: . nconv - number of converged eigenpairs so far
563: . eigr - real part of the eigenvalues
564: . eigi - imaginary part of the eigenvalues
565: . errest - error estimates
566: . nest - number of error estimates to display
567: - vf - viewer and format for monitoring
569: Options Database Key:
570: . -eps_monitor_conv draw::draw_lg - activates EPSMonitorConvergedDrawLG()
572: Level: intermediate
574: .seealso: EPSMonitorSet()
575: @*/
576: PetscErrorCode EPSMonitorConvergedDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
577: {
578: PetscErrorCode ierr;
579: PetscViewer viewer;
580: PetscDrawLG lg;
581: PetscReal x,y;
586: viewer = vf->viewer;
588: lg = vf->lg;
590: PetscViewerPushFormat(viewer,vf->format);
591: if (its==1) {
592: PetscDrawLGReset(lg);
593: PetscDrawLGSetDimension(lg,1);
594: PetscDrawLGSetLimits(lg,1,1,0,eps->nev);
595: }
596: x = (PetscReal)its;
597: y = (PetscReal)eps->nconv;
598: PetscDrawLGAddPoint(lg,&x,&y);
599: if (its <= 20 || !(its % 5) || eps->reason) {
600: PetscDrawLGDraw(lg);
601: PetscDrawLGSave(lg);
602: }
603: PetscViewerPopFormat(viewer);
604: return(0);
605: }
607: /*@C
608: EPSMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
610: Collective on viewer
612: Input Parameters:
613: + viewer - the viewer
614: . format - the viewer format
615: - ctx - an optional user context
617: Output Parameter:
618: . vf - the viewer and format context
620: Level: intermediate
622: .seealso: EPSMonitorSet()
623: @*/
624: PetscErrorCode EPSMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
625: {
627: SlepcConvMon mctx;
630: PetscViewerAndFormatCreate(viewer,format,vf);
631: PetscNew(&mctx);
632: mctx->ctx = ctx;
633: (*vf)->data = (void*)mctx;
634: EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
635: return(0);
636: }