Actual source code: pepmon.c
slepc-3.16.2 2022-02-01
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: 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;
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 PEPMonitor(PEP pep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
41: {
43: PetscInt i,n = pep->numbermonitors;
46: for (i=0;i<n;i++) {
47: (*pep->monitor[i])(pep,it,nconv,eigr,eigi,errest,nest,pep->monitorcontext[i]);
48: }
49: return(0);
50: }
52: /*@C
53: PEPMonitorSet - Sets an ADDITIONAL function to be called at every
54: iteration to monitor the error estimates for each requested eigenpair.
56: Logically Collective on pep
58: Input Parameters:
59: + pep - eigensolver context obtained from PEPCreate()
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(PEP pep,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)
68: + pep - polynomial eigensolver context obtained from PEPCreate()
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 PEPMonitorSet()
77: Options Database Keys:
78: + -pep_monitor - print only the first error estimate
79: . -pep_monitor_all - print error estimates at each iteration
80: . -pep_monitor_conv - print the eigenvalue approximations only when
81: convergence has been reached
82: . -pep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
83: approximate eigenvalue
84: . -pep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
85: approximate eigenvalues
86: . -pep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
87: - -pep_monitor_cancel - cancels all monitors that have been hardwired into
88: a code by calls to PEPMonitorSet(), but does not cancel those set via
89: the options database.
91: Notes:
92: Several different monitoring routines may be set by calling
93: PEPMonitorSet() multiple times; all will be called in the
94: order in which they were set.
96: Level: intermediate
98: .seealso: PEPMonitorFirst(), PEPMonitorAll(), PEPMonitorCancel()
99: @*/
100: PetscErrorCode PEPMonitorSet(PEP pep,PetscErrorCode (*monitor)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
101: {
104: if (pep->numbermonitors >= MAXPEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Too many PEP monitors set");
105: pep->monitor[pep->numbermonitors] = monitor;
106: pep->monitorcontext[pep->numbermonitors] = (void*)mctx;
107: pep->monitordestroy[pep->numbermonitors++] = monitordestroy;
108: return(0);
109: }
111: /*@
112: PEPMonitorCancel - Clears all monitors for a PEP object.
114: Logically Collective on pep
116: Input Parameters:
117: . pep - eigensolver context obtained from PEPCreate()
119: Options Database Key:
120: . -pep_monitor_cancel - Cancels all monitors that have been hardwired
121: into a code by calls to PEPMonitorSet(),
122: but does not cancel those set via the options database.
124: Level: intermediate
126: .seealso: PEPMonitorSet()
127: @*/
128: PetscErrorCode PEPMonitorCancel(PEP pep)
129: {
131: PetscInt i;
135: for (i=0; i<pep->numbermonitors; i++) {
136: if (pep->monitordestroy[i]) {
137: (*pep->monitordestroy[i])(&pep->monitorcontext[i]);
138: }
139: }
140: pep->numbermonitors = 0;
141: return(0);
142: }
144: /*@C
145: PEPGetMonitorContext - Gets the monitor context, as set by
146: PEPMonitorSet() for the FIRST monitor only.
148: Not Collective
150: Input Parameter:
151: . pep - eigensolver context obtained from PEPCreate()
153: Output Parameter:
154: . ctx - monitor context
156: Level: intermediate
158: .seealso: PEPMonitorSet(), PEPDefaultMonitor()
159: @*/
160: PetscErrorCode PEPGetMonitorContext(PEP pep,void *ctx)
161: {
164: *(void**)ctx = pep->monitorcontext[0];
165: return(0);
166: }
168: /*
169: Helper function to compute eigenvalue that must be viewed in monitor
170: */
171: static PetscErrorCode PEPMonitorGetTrueEig(PEP pep,PetscScalar *er,PetscScalar *ei)
172: {
174: PetscBool flg,sinv;
177: STBackTransform(pep->st,1,er,ei);
178: if (pep->sfactor!=1.0) {
179: STGetTransform(pep->st,&flg);
180: PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
181: if (flg && sinv) {
182: *er /= pep->sfactor; *ei /= pep->sfactor;
183: } else {
184: *er *= pep->sfactor; *ei *= pep->sfactor;
185: }
186: }
187: return(0);
188: }
190: /*@C
191: PEPMonitorFirst - Print the first unconverged approximate value and
192: error estimate at each iteration of the polynomial eigensolver.
194: Collective on pep
196: Input Parameters:
197: + pep - polynomial eigensolver context
198: . its - iteration number
199: . nconv - number of converged eigenpairs so far
200: . eigr - real part of the eigenvalues
201: . eigi - imaginary part of the eigenvalues
202: . errest - error estimates
203: . nest - number of error estimates to display
204: - vf - viewer and format for monitoring
206: Options Database Key:
207: . -pep_monitor - activates PEPMonitorFirst()
209: Level: intermediate
211: .seealso: PEPMonitorSet(), PEPMonitorAll(), PEPMonitorConverged()
212: @*/
213: PetscErrorCode PEPMonitorFirst(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
214: {
216: PetscScalar er,ei;
217: PetscViewer viewer = vf->viewer;
222: if (its==1 && ((PetscObject)pep)->prefix) {
223: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
224: }
225: if (nconv<nest) {
226: PetscViewerPushFormat(viewer,vf->format);
227: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
228: PetscViewerASCIIPrintf(viewer,"%3D PEP nconv=%D first unconverged value (error)",its,nconv);
229: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
230: er = eigr[nconv]; ei = eigi[nconv];
231: PEPMonitorGetTrueEig(pep,&er,&ei);
232: #if defined(PETSC_USE_COMPLEX)
233: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
234: #else
235: PetscViewerASCIIPrintf(viewer," %g",(double)er);
236: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
237: #endif
238: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
239: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
240: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
241: PetscViewerPopFormat(viewer);
242: }
243: return(0);
244: }
246: /*@C
247: PEPMonitorAll - Print the current approximate values and
248: error estimates at each iteration of the polynomial eigensolver.
250: Collective on pep
252: Input Parameters:
253: + pep - polynomial eigensolver context
254: . its - iteration number
255: . nconv - number of converged eigenpairs so far
256: . eigr - real part of the eigenvalues
257: . eigi - imaginary part of the eigenvalues
258: . errest - error estimates
259: . nest - number of error estimates to display
260: - vf - viewer and format for monitoring
262: Options Database Key:
263: . -pep_monitor_all - activates PEPMonitorAll()
265: Level: intermediate
267: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorConverged()
268: @*/
269: PetscErrorCode PEPMonitorAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
270: {
272: PetscInt i;
273: PetscScalar er,ei;
274: PetscViewer viewer = vf->viewer;
279: PetscViewerPushFormat(viewer,vf->format);
280: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
281: if (its==1 && ((PetscObject)pep)->prefix) {
282: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
283: }
284: PetscViewerASCIIPrintf(viewer,"%3D PEP nconv=%D Values (Errors)",its,nconv);
285: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
286: for (i=0;i<nest;i++) {
287: er = eigr[i]; ei = eigi[i];
288: PEPMonitorGetTrueEig(pep,&er,&ei);
289: #if defined(PETSC_USE_COMPLEX)
290: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
291: #else
292: PetscViewerASCIIPrintf(viewer," %g",(double)er);
293: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
294: #endif
295: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
296: }
297: PetscViewerASCIIPrintf(viewer,"\n");
298: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
299: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
300: PetscViewerPopFormat(viewer);
301: return(0);
302: }
304: /*@C
305: PEPMonitorConverged - Print the approximate values and
306: error estimates as they converge.
308: Collective on pep
310: Input Parameters:
311: + pep - polynomial eigensolver context
312: . its - iteration number
313: . nconv - number of converged eigenpairs so far
314: . eigr - real part of the eigenvalues
315: . eigi - imaginary part of the eigenvalues
316: . errest - error estimates
317: . nest - number of error estimates to display
318: - vf - viewer and format for monitoring
320: Options Database Key:
321: . -pep_monitor_conv - activates PEPMonitorConverged()
323: Level: intermediate
325: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorAll()
326: @*/
327: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
328: {
330: PetscInt i;
331: PetscScalar er,ei;
332: PetscViewer viewer = vf->viewer;
333: SlepcConvMon ctx;
338: ctx = (SlepcConvMon)vf->data;
339: if (its==1 && ((PetscObject)pep)->prefix) {
340: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)pep)->prefix);
341: }
342: if (its==1) ctx->oldnconv = 0;
343: if (ctx->oldnconv!=nconv) {
344: PetscViewerPushFormat(viewer,vf->format);
345: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
346: for (i=ctx->oldnconv;i<nconv;i++) {
347: PetscViewerASCIIPrintf(viewer,"%3D PEP converged value (error) #%D",its,i);
348: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
349: er = eigr[i]; ei = eigi[i];
350: PEPMonitorGetTrueEig(pep,&er,&ei);
351: #if defined(PETSC_USE_COMPLEX)
352: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
353: #else
354: PetscViewerASCIIPrintf(viewer," %g",(double)er);
355: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
356: #endif
357: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
358: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
359: }
360: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
361: PetscViewerPopFormat(viewer);
362: ctx->oldnconv = nconv;
363: }
364: return(0);
365: }
367: PetscErrorCode PEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
368: {
370: SlepcConvMon mctx;
373: PetscViewerAndFormatCreate(viewer,format,vf);
374: PetscNew(&mctx);
375: mctx->ctx = ctx;
376: (*vf)->data = (void*)mctx;
377: return(0);
378: }
380: PetscErrorCode PEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
381: {
385: if (!*vf) return(0);
386: PetscFree((*vf)->data);
387: PetscViewerDestroy(&(*vf)->viewer);
388: PetscDrawLGDestroy(&(*vf)->lg);
389: PetscFree(*vf);
390: return(0);
391: }
393: /*@C
394: PEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
395: approximation at each iteration of the polynomial eigensolver.
397: Collective on pep
399: Input Parameters:
400: + pep - polynomial eigensolver context
401: . its - iteration number
402: . its - iteration number
403: . nconv - number of converged eigenpairs so far
404: . eigr - real part of the eigenvalues
405: . eigi - imaginary part of the eigenvalues
406: . errest - error estimates
407: . nest - number of error estimates to display
408: - vf - viewer and format for monitoring
410: Options Database Key:
411: . -pep_monitor draw::draw_lg - activates PEPMonitorFirstDrawLG()
413: Level: intermediate
415: .seealso: PEPMonitorSet()
416: @*/
417: PetscErrorCode PEPMonitorFirstDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
418: {
420: PetscViewer viewer = vf->viewer;
421: PetscDrawLG lg = vf->lg;
422: PetscReal x,y;
428: PetscViewerPushFormat(viewer,vf->format);
429: if (its==1) {
430: PetscDrawLGReset(lg);
431: PetscDrawLGSetDimension(lg,1);
432: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0);
433: }
434: if (nconv<nest) {
435: x = (PetscReal)its;
436: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
437: else y = 0.0;
438: PetscDrawLGAddPoint(lg,&x,&y);
439: if (its <= 20 || !(its % 5) || pep->reason) {
440: PetscDrawLGDraw(lg);
441: PetscDrawLGSave(lg);
442: }
443: }
444: PetscViewerPopFormat(viewer);
445: return(0);
446: }
448: /*@C
449: PEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
451: Collective on viewer
453: Input Parameters:
454: + viewer - the viewer
455: . format - the viewer format
456: - ctx - an optional user context
458: Output Parameter:
459: . vf - the viewer and format context
461: Level: intermediate
463: .seealso: PEPMonitorSet()
464: @*/
465: PetscErrorCode PEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
466: {
470: PetscViewerAndFormatCreate(viewer,format,vf);
471: (*vf)->data = ctx;
472: PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
473: return(0);
474: }
476: /*@C
477: PEPMonitorAllDrawLG - Plots the error estimate of all unconverged
478: approximations at each iteration of the polynomial eigensolver.
480: Collective on pep
482: Input Parameters:
483: + pep - polynomial eigensolver context
484: . its - iteration number
485: . its - iteration number
486: . nconv - number of converged eigenpairs so far
487: . eigr - real part of the eigenvalues
488: . eigi - imaginary part of the eigenvalues
489: . errest - error estimates
490: . nest - number of error estimates to display
491: - vf - viewer and format for monitoring
493: Options Database Key:
494: . -pep_monitor_all draw::draw_lg - activates PEPMonitorAllDrawLG()
496: Level: intermediate
498: .seealso: PEPMonitorSet()
499: @*/
500: PetscErrorCode PEPMonitorAllDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
501: {
503: PetscViewer viewer = vf->viewer;
504: PetscDrawLG lg = vf->lg;
505: PetscInt i,n = PetscMin(pep->nev,255);
506: PetscReal *x,*y;
512: PetscViewerPushFormat(viewer,vf->format);
513: if (its==1) {
514: PetscDrawLGReset(lg);
515: PetscDrawLGSetDimension(lg,n);
516: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0);
517: }
518: PetscMalloc2(n,&x,n,&y);
519: for (i=0;i<n;i++) {
520: x[i] = (PetscReal)its;
521: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
522: else y[i] = 0.0;
523: }
524: PetscDrawLGAddPoint(lg,x,y);
525: if (its <= 20 || !(its % 5) || pep->reason) {
526: PetscDrawLGDraw(lg);
527: PetscDrawLGSave(lg);
528: }
529: PetscFree2(x,y);
530: PetscViewerPopFormat(viewer);
531: return(0);
532: }
534: /*@C
535: PEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
537: Collective on viewer
539: Input Parameters:
540: + viewer - the viewer
541: . format - the viewer format
542: - ctx - an optional user context
544: Output Parameter:
545: . vf - the viewer and format context
547: Level: intermediate
549: .seealso: PEPMonitorSet()
550: @*/
551: PetscErrorCode PEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
552: {
556: PetscViewerAndFormatCreate(viewer,format,vf);
557: (*vf)->data = ctx;
558: PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
559: return(0);
560: }
562: /*@C
563: PEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
564: at each iteration of the polynomial eigensolver.
566: Collective on pep
568: Input Parameters:
569: + pep - polynomial eigensolver context
570: . its - iteration number
571: . its - iteration number
572: . nconv - number of converged eigenpairs so far
573: . eigr - real part of the eigenvalues
574: . eigi - imaginary part of the eigenvalues
575: . errest - error estimates
576: . nest - number of error estimates to display
577: - vf - viewer and format for monitoring
579: Options Database Key:
580: . -pep_monitor_conv draw::draw_lg - activates PEPMonitorConvergedDrawLG()
582: Level: intermediate
584: .seealso: PEPMonitorSet()
585: @*/
586: PetscErrorCode PEPMonitorConvergedDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
587: {
588: PetscErrorCode ierr;
589: PetscViewer viewer = vf->viewer;
590: PetscDrawLG lg = vf->lg;
591: PetscReal x,y;
597: PetscViewerPushFormat(viewer,vf->format);
598: if (its==1) {
599: PetscDrawLGReset(lg);
600: PetscDrawLGSetDimension(lg,1);
601: PetscDrawLGSetLimits(lg,1,1,0,pep->nev);
602: }
603: x = (PetscReal)its;
604: y = (PetscReal)pep->nconv;
605: PetscDrawLGAddPoint(lg,&x,&y);
606: if (its <= 20 || !(its % 5) || pep->reason) {
607: PetscDrawLGDraw(lg);
608: PetscDrawLGSave(lg);
609: }
610: PetscViewerPopFormat(viewer);
611: return(0);
612: }
614: /*@C
615: PEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
617: Collective on viewer
619: Input Parameters:
620: + viewer - the viewer
621: . format - the viewer format
622: - ctx - an optional user context
624: Output Parameter:
625: . vf - the viewer and format context
627: Level: intermediate
629: .seealso: PEPMonitorSet()
630: @*/
631: PetscErrorCode PEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
632: {
634: SlepcConvMon mctx;
637: PetscViewerAndFormatCreate(viewer,format,vf);
638: PetscNew(&mctx);
639: mctx->ctx = ctx;
640: (*vf)->data = (void*)mctx;
641: PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
642: return(0);
643: }