Actual source code: pepmon.c
slepc-3.18.0 2022-10-01
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: PEP routines related to monitors
12: */
14: #include <slepc/private/pepimpl.h>
15: #include <petscdraw.h>
17: PetscErrorCode PEPMonitorLGCreate(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: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
24: PetscDrawSetFromOptions(draw);
25: PetscDrawLGCreate(draw,l,&lg);
26: if (names) PetscDrawLGSetLegend(lg,names);
27: PetscDrawLGSetFromOptions(lg);
28: PetscDrawLGGetAxis(lg,&axis);
29: PetscDrawAxisSetLabels(axis,"Convergence","Iteration",metric);
30: PetscDrawDestroy(&draw);
31: *lgctx = lg;
32: return 0;
33: }
35: /*
36: Runs the user provided monitor routines, if any.
37: */
38: PetscErrorCode PEPMonitor(PEP pep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
39: {
40: PetscInt i,n = pep->numbermonitors;
42: for (i=0;i<n;i++) (*pep->monitor[i])(pep,it,nconv,eigr,eigi,errest,nest,pep->monitorcontext[i]);
43: return 0;
44: }
46: /*@C
47: PEPMonitorSet - Sets an ADDITIONAL function to be called at every
48: iteration to monitor the error estimates for each requested eigenpair.
50: Logically Collective on pep
52: Input Parameters:
53: + pep - eigensolver context obtained from PEPCreate()
54: . monitor - pointer to function (if this is NULL, it turns off monitoring)
55: . mctx - [optional] context for private data for the
56: monitor routine (use NULL if no context is desired)
57: - monitordestroy - [optional] routine that frees monitor context (may be NULL)
59: Calling Sequence of monitor:
60: $ monitor(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,PetscInt nest,void *mctx)
62: + pep - polynomial eigensolver context obtained from PEPCreate()
63: . its - iteration number
64: . nconv - number of converged eigenpairs
65: . eigr - real part of the eigenvalues
66: . eigi - imaginary part of the eigenvalues
67: . errest - relative error estimates for each eigenpair
68: . nest - number of error estimates
69: - mctx - optional monitoring context, as set by PEPMonitorSet()
71: Options Database Keys:
72: + -pep_monitor - print only the first error estimate
73: . -pep_monitor_all - print error estimates at each iteration
74: . -pep_monitor_conv - print the eigenvalue approximations only when
75: convergence has been reached
76: . -pep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
77: approximate eigenvalue
78: . -pep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
79: approximate eigenvalues
80: . -pep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
81: - -pep_monitor_cancel - cancels all monitors that have been hardwired into
82: a code by calls to PEPMonitorSet(), but does not cancel those set via
83: the options database.
85: Notes:
86: Several different monitoring routines may be set by calling
87: PEPMonitorSet() multiple times; all will be called in the
88: order in which they were set.
90: Level: intermediate
92: .seealso: PEPMonitorFirst(), PEPMonitorAll(), PEPMonitorCancel()
93: @*/
94: PetscErrorCode PEPMonitorSet(PEP pep,PetscErrorCode (*monitor)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
95: {
98: pep->monitor[pep->numbermonitors] = monitor;
99: pep->monitorcontext[pep->numbermonitors] = (void*)mctx;
100: pep->monitordestroy[pep->numbermonitors++] = monitordestroy;
101: return 0;
102: }
104: /*@
105: PEPMonitorCancel - Clears all monitors for a PEP object.
107: Logically Collective on pep
109: Input Parameters:
110: . pep - eigensolver context obtained from PEPCreate()
112: Options Database Key:
113: . -pep_monitor_cancel - Cancels all monitors that have been hardwired
114: into a code by calls to PEPMonitorSet(),
115: but does not cancel those set via the options database.
117: Level: intermediate
119: .seealso: PEPMonitorSet()
120: @*/
121: PetscErrorCode PEPMonitorCancel(PEP pep)
122: {
123: PetscInt i;
126: for (i=0; i<pep->numbermonitors; i++) {
127: if (pep->monitordestroy[i]) (*pep->monitordestroy[i])(&pep->monitorcontext[i]);
128: }
129: pep->numbermonitors = 0;
130: return 0;
131: }
133: /*@C
134: PEPGetMonitorContext - Gets the monitor context, as set by
135: PEPMonitorSet() for the FIRST monitor only.
137: Not Collective
139: Input Parameter:
140: . pep - eigensolver context obtained from PEPCreate()
142: Output Parameter:
143: . ctx - monitor context
145: Level: intermediate
147: .seealso: PEPMonitorSet(), PEPDefaultMonitor()
148: @*/
149: PetscErrorCode PEPGetMonitorContext(PEP pep,void *ctx)
150: {
152: *(void**)ctx = pep->monitorcontext[0];
153: return 0;
154: }
156: /*
157: Helper function to compute eigenvalue that must be viewed in monitor
158: */
159: static PetscErrorCode PEPMonitorGetTrueEig(PEP pep,PetscScalar *er,PetscScalar *ei)
160: {
161: PetscBool flg,sinv;
163: STBackTransform(pep->st,1,er,ei);
164: if (pep->sfactor!=1.0) {
165: STGetTransform(pep->st,&flg);
166: PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
167: if (flg && sinv) {
168: *er /= pep->sfactor; *ei /= pep->sfactor;
169: } else {
170: *er *= pep->sfactor; *ei *= pep->sfactor;
171: }
172: }
173: return 0;
174: }
176: /*@C
177: PEPMonitorFirst - Print the first unconverged approximate value and
178: error estimate at each iteration of the polynomial eigensolver.
180: Collective on pep
182: Input Parameters:
183: + pep - polynomial eigensolver context
184: . its - iteration number
185: . nconv - number of converged eigenpairs so far
186: . eigr - real part of the eigenvalues
187: . eigi - imaginary part of the eigenvalues
188: . errest - error estimates
189: . nest - number of error estimates to display
190: - vf - viewer and format for monitoring
192: Options Database Key:
193: . -pep_monitor - activates PEPMonitorFirst()
195: Level: intermediate
197: .seealso: PEPMonitorSet(), PEPMonitorAll(), PEPMonitorConverged()
198: @*/
199: PetscErrorCode PEPMonitorFirst(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
200: {
201: PetscScalar er,ei;
202: PetscViewer viewer = vf->viewer;
206: if (its==1 && ((PetscObject)pep)->prefix) PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
207: if (nconv<nest) {
208: PetscViewerPushFormat(viewer,vf->format);
209: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
210: PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv);
211: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
212: er = eigr[nconv]; ei = eigi[nconv];
213: PEPMonitorGetTrueEig(pep,&er,&ei);
214: #if defined(PETSC_USE_COMPLEX)
215: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
216: #else
217: PetscViewerASCIIPrintf(viewer," %g",(double)er);
218: if (ei!=0.0) PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei);
219: #endif
220: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
221: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
222: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
223: PetscViewerPopFormat(viewer);
224: }
225: return 0;
226: }
228: /*@C
229: PEPMonitorAll - Print the current approximate values and
230: error estimates at each iteration of the polynomial eigensolver.
232: Collective on pep
234: Input Parameters:
235: + pep - polynomial eigensolver context
236: . its - iteration number
237: . nconv - number of converged eigenpairs so far
238: . eigr - real part of the eigenvalues
239: . eigi - imaginary part of the eigenvalues
240: . errest - error estimates
241: . nest - number of error estimates to display
242: - vf - viewer and format for monitoring
244: Options Database Key:
245: . -pep_monitor_all - activates PEPMonitorAll()
247: Level: intermediate
249: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorConverged()
250: @*/
251: PetscErrorCode PEPMonitorAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
252: {
253: PetscInt i;
254: PetscScalar er,ei;
255: PetscViewer viewer = vf->viewer;
259: PetscViewerPushFormat(viewer,vf->format);
260: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
261: if (its==1 && ((PetscObject)pep)->prefix) PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
262: PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv);
263: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
264: for (i=0;i<nest;i++) {
265: er = eigr[i]; ei = eigi[i];
266: PEPMonitorGetTrueEig(pep,&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)pep)->tablevel);
278: PetscViewerPopFormat(viewer);
279: return 0;
280: }
282: /*@C
283: PEPMonitorConverged - Print the approximate values and
284: error estimates as they converge.
286: Collective on pep
288: Input Parameters:
289: + pep - polynomial 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: . -pep_monitor_conv - activates PEPMonitorConverged()
301: Level: intermediate
303: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorAll()
304: @*/
305: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
306: {
307: PetscInt i;
308: PetscScalar er,ei;
309: PetscViewer viewer = vf->viewer;
310: SlepcConvMon ctx;
314: ctx = (SlepcConvMon)vf->data;
315: if (its==1 && ((PetscObject)pep)->prefix) PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)pep)->prefix);
316: if (its==1) ctx->oldnconv = 0;
317: if (ctx->oldnconv!=nconv) {
318: PetscViewerPushFormat(viewer,vf->format);
319: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
320: for (i=ctx->oldnconv;i<nconv;i++) {
321: PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP converged value (error) #%" PetscInt_FMT,its,i);
322: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
323: er = eigr[i]; ei = eigi[i];
324: PEPMonitorGetTrueEig(pep,&er,&ei);
325: #if defined(PETSC_USE_COMPLEX)
326: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
327: #else
328: PetscViewerASCIIPrintf(viewer," %g",(double)er);
329: if (ei!=0.0) PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei);
330: #endif
331: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
332: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
333: }
334: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
335: PetscViewerPopFormat(viewer);
336: ctx->oldnconv = nconv;
337: }
338: return 0;
339: }
341: PetscErrorCode PEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
342: {
343: SlepcConvMon mctx;
345: PetscViewerAndFormatCreate(viewer,format,vf);
346: PetscNew(&mctx);
347: mctx->ctx = ctx;
348: (*vf)->data = (void*)mctx;
349: return 0;
350: }
352: PetscErrorCode PEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
353: {
354: if (!*vf) return 0;
355: PetscFree((*vf)->data);
356: PetscViewerDestroy(&(*vf)->viewer);
357: PetscDrawLGDestroy(&(*vf)->lg);
358: PetscFree(*vf);
359: return 0;
360: }
362: /*@C
363: PEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
364: approximation at each iteration of the polynomial eigensolver.
366: Collective on pep
368: Input Parameters:
369: + pep - polynomial eigensolver context
370: . its - iteration number
371: . its - iteration number
372: . nconv - number of converged eigenpairs so far
373: . eigr - real part of the eigenvalues
374: . eigi - imaginary part of the eigenvalues
375: . errest - error estimates
376: . nest - number of error estimates to display
377: - vf - viewer and format for monitoring
379: Options Database Key:
380: . -pep_monitor draw::draw_lg - activates PEPMonitorFirstDrawLG()
382: Level: intermediate
384: .seealso: PEPMonitorSet()
385: @*/
386: PetscErrorCode PEPMonitorFirstDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
387: {
388: PetscViewer viewer = vf->viewer;
389: PetscDrawLG lg = vf->lg;
390: PetscReal x,y;
395: PetscViewerPushFormat(viewer,vf->format);
396: if (its==1) {
397: PetscDrawLGReset(lg);
398: PetscDrawLGSetDimension(lg,1);
399: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0);
400: }
401: if (nconv<nest) {
402: x = (PetscReal)its;
403: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
404: else y = 0.0;
405: PetscDrawLGAddPoint(lg,&x,&y);
406: if (its <= 20 || !(its % 5) || pep->reason) {
407: PetscDrawLGDraw(lg);
408: PetscDrawLGSave(lg);
409: }
410: }
411: PetscViewerPopFormat(viewer);
412: return 0;
413: }
415: /*@C
416: PEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
418: Collective on viewer
420: Input Parameters:
421: + viewer - the viewer
422: . format - the viewer format
423: - ctx - an optional user context
425: Output Parameter:
426: . vf - the viewer and format context
428: Level: intermediate
430: .seealso: PEPMonitorSet()
431: @*/
432: PetscErrorCode PEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
433: {
434: PetscViewerAndFormatCreate(viewer,format,vf);
435: (*vf)->data = ctx;
436: PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
437: return 0;
438: }
440: /*@C
441: PEPMonitorAllDrawLG - Plots the error estimate of all unconverged
442: approximations at each iteration of the polynomial eigensolver.
444: Collective on pep
446: Input Parameters:
447: + pep - polynomial eigensolver context
448: . its - iteration number
449: . its - iteration number
450: . nconv - number of converged eigenpairs so far
451: . eigr - real part of the eigenvalues
452: . eigi - imaginary part of the eigenvalues
453: . errest - error estimates
454: . nest - number of error estimates to display
455: - vf - viewer and format for monitoring
457: Options Database Key:
458: . -pep_monitor_all draw::draw_lg - activates PEPMonitorAllDrawLG()
460: Level: intermediate
462: .seealso: PEPMonitorSet()
463: @*/
464: PetscErrorCode PEPMonitorAllDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
465: {
466: PetscViewer viewer = vf->viewer;
467: PetscDrawLG lg = vf->lg;
468: PetscInt i,n = PetscMin(pep->nev,255);
469: PetscReal *x,*y;
474: PetscViewerPushFormat(viewer,vf->format);
475: if (its==1) {
476: PetscDrawLGReset(lg);
477: PetscDrawLGSetDimension(lg,n);
478: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0);
479: }
480: PetscMalloc2(n,&x,n,&y);
481: for (i=0;i<n;i++) {
482: x[i] = (PetscReal)its;
483: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
484: else y[i] = 0.0;
485: }
486: PetscDrawLGAddPoint(lg,x,y);
487: if (its <= 20 || !(its % 5) || pep->reason) {
488: PetscDrawLGDraw(lg);
489: PetscDrawLGSave(lg);
490: }
491: PetscFree2(x,y);
492: PetscViewerPopFormat(viewer);
493: return 0;
494: }
496: /*@C
497: PEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
499: Collective on viewer
501: Input Parameters:
502: + viewer - the viewer
503: . format - the viewer format
504: - ctx - an optional user context
506: Output Parameter:
507: . vf - the viewer and format context
509: Level: intermediate
511: .seealso: PEPMonitorSet()
512: @*/
513: PetscErrorCode PEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
514: {
515: PetscViewerAndFormatCreate(viewer,format,vf);
516: (*vf)->data = ctx;
517: PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
518: return 0;
519: }
521: /*@C
522: PEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
523: at each iteration of the polynomial eigensolver.
525: Collective on pep
527: Input Parameters:
528: + pep - polynomial eigensolver context
529: . its - iteration number
530: . its - iteration number
531: . nconv - number of converged eigenpairs so far
532: . eigr - real part of the eigenvalues
533: . eigi - imaginary part of the eigenvalues
534: . errest - error estimates
535: . nest - number of error estimates to display
536: - vf - viewer and format for monitoring
538: Options Database Key:
539: . -pep_monitor_conv draw::draw_lg - activates PEPMonitorConvergedDrawLG()
541: Level: intermediate
543: .seealso: PEPMonitorSet()
544: @*/
545: PetscErrorCode PEPMonitorConvergedDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
546: {
547: PetscViewer viewer = vf->viewer;
548: PetscDrawLG lg = vf->lg;
549: PetscReal x,y;
554: PetscViewerPushFormat(viewer,vf->format);
555: if (its==1) {
556: PetscDrawLGReset(lg);
557: PetscDrawLGSetDimension(lg,1);
558: PetscDrawLGSetLimits(lg,1,1,0,pep->nev);
559: }
560: x = (PetscReal)its;
561: y = (PetscReal)pep->nconv;
562: PetscDrawLGAddPoint(lg,&x,&y);
563: if (its <= 20 || !(its % 5) || pep->reason) {
564: PetscDrawLGDraw(lg);
565: PetscDrawLGSave(lg);
566: }
567: PetscViewerPopFormat(viewer);
568: return 0;
569: }
571: /*@C
572: PEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
574: Collective on viewer
576: Input Parameters:
577: + viewer - the viewer
578: . format - the viewer format
579: - ctx - an optional user context
581: Output Parameter:
582: . vf - the viewer and format context
584: Level: intermediate
586: .seealso: PEPMonitorSet()
587: @*/
588: PetscErrorCode PEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
589: {
590: SlepcConvMon mctx;
592: PetscViewerAndFormatCreate(viewer,format,vf);
593: PetscNew(&mctx);
594: mctx->ctx = ctx;
595: (*vf)->data = (void*)mctx;
596: PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
597: return 0;
598: }