Actual source code: nepmon.c
slepc-3.16.0 2021-09-30
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: NEP routines related to monitors
12: */
14: #include <slepc/private/nepimpl.h>
15: #include <petscdraw.h>
17: PetscErrorCode NEPMonitorLGCreate(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 NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
41: {
43: PetscInt i,n = nep->numbermonitors;
46: for (i=0;i<n;i++) {
47: (*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]);
48: }
49: return(0);
50: }
52: /*@C
53: NEPMonitorSet - Sets an ADDITIONAL function to be called at every
54: iteration to monitor the error estimates for each requested eigenpair.
56: Logically Collective on nep
58: Input Parameters:
59: + nep - eigensolver context obtained from NEPCreate()
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(NEP nep,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)
68: + nep - nonlinear eigensolver context obtained from NEPCreate()
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 - error estimates for each eigenpair
74: . nest - number of error estimates
75: - mctx - optional monitoring context, as set by NEPMonitorSet()
77: Options Database Keys:
78: + -nep_monitor - print only the first error estimate
79: . -nep_monitor_all - print error estimates at each iteration
80: . -nep_monitor_conv - print the eigenvalue approximations only when
81: convergence has been reached
82: . -nep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
83: approximate eigenvalue
84: . -nep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
85: approximate eigenvalues
86: . -nep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
87: - -nep_monitor_cancel - cancels all monitors that have been hardwired into
88: a code by calls to NEPMonitorSet(), but does not cancel those set via
89: the options database.
91: Notes:
92: Several different monitoring routines may be set by calling
93: NEPMonitorSet() multiple times; all will be called in the
94: order in which they were set.
96: Level: intermediate
98: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
99: @*/
100: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
101: {
104: if (nep->numbermonitors >= MAXNEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
105: nep->monitor[nep->numbermonitors] = monitor;
106: nep->monitorcontext[nep->numbermonitors] = (void*)mctx;
107: nep->monitordestroy[nep->numbermonitors++] = monitordestroy;
108: return(0);
109: }
111: /*@
112: NEPMonitorCancel - Clears all monitors for a NEP object.
114: Logically Collective on nep
116: Input Parameters:
117: . nep - eigensolver context obtained from NEPCreate()
119: Options Database Key:
120: . -nep_monitor_cancel - Cancels all monitors that have been hardwired
121: into a code by calls to NEPMonitorSet(),
122: but does not cancel those set via the options database.
124: Level: intermediate
126: .seealso: NEPMonitorSet()
127: @*/
128: PetscErrorCode NEPMonitorCancel(NEP nep)
129: {
131: PetscInt i;
135: for (i=0; i<nep->numbermonitors; i++) {
136: if (nep->monitordestroy[i]) {
137: (*nep->monitordestroy[i])(&nep->monitorcontext[i]);
138: }
139: }
140: nep->numbermonitors = 0;
141: return(0);
142: }
144: /*@C
145: NEPGetMonitorContext - Gets the monitor context, as set by
146: NEPMonitorSet() for the FIRST monitor only.
148: Not Collective
150: Input Parameter:
151: . nep - eigensolver context obtained from NEPCreate()
153: Output Parameter:
154: . ctx - monitor context
156: Level: intermediate
158: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
159: @*/
160: PetscErrorCode NEPGetMonitorContext(NEP nep,void *ctx)
161: {
164: *(void**)ctx = nep->monitorcontext[0];
165: return(0);
166: }
168: /*@C
169: NEPMonitorFirst - Print the first unconverged approximate value and
170: error estimate at each iteration of the nonlinear eigensolver.
172: Collective on nep
174: Input Parameters:
175: + nep - nonlinear 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: . -nep_monitor - activates NEPMonitorFirst()
187: Level: intermediate
189: .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
190: @*/
191: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
192: {
194: PetscViewer viewer = vf->viewer;
199: if (its==1 && ((PetscObject)nep)->prefix) {
200: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
201: }
202: if (nconv<nest) {
203: PetscViewerPushFormat(viewer,vf->format);
204: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
205: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D first unconverged value (error)",its,nconv);
206: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
207: #if defined(PETSC_USE_COMPLEX)
208: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv]));
209: #else
210: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]);
211: if (eigi[nconv]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]); }
212: #endif
213: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
214: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
215: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
216: PetscViewerPopFormat(viewer);
217: }
218: return(0);
219: }
221: /*@C
222: NEPMonitorAll - Print the current approximate values and
223: error estimates at each iteration of the nonlinear eigensolver.
225: Collective on nep
227: Input Parameters:
228: + nep - nonlinear eigensolver context
229: . its - iteration number
230: . nconv - number of converged eigenpairs so far
231: . eigr - real part of the eigenvalues
232: . eigi - imaginary part of the eigenvalues
233: . errest - error estimates
234: . nest - number of error estimates to display
235: - vf - viewer and format for monitoring
237: Options Database Key:
238: . -nep_monitor_all - activates NEPMonitorAll()
240: Level: intermediate
242: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
243: @*/
244: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
245: {
247: PetscInt i;
248: PetscViewer viewer = vf->viewer;
253: PetscViewerPushFormat(viewer,vf->format);
254: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
255: if (its==1 && ((PetscObject)nep)->prefix) {
256: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
257: }
258: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D Values (Errors)",its,nconv);
259: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
260: for (i=0;i<nest;i++) {
261: #if defined(PETSC_USE_COMPLEX)
262: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
263: #else
264: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
265: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
266: #endif
267: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
268: }
269: PetscViewerASCIIPrintf(viewer,"\n");
270: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
271: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
272: PetscViewerPopFormat(viewer);
273: return(0);
274: }
276: /*@C
277: NEPMonitorConverged - Print the approximate values and
278: error estimates as they converge.
280: Collective on nep
282: Input Parameters:
283: + nep - nonlinear eigensolver context
284: . its - iteration number
285: . nconv - number of converged eigenpairs so far
286: . eigr - real part of the eigenvalues
287: . eigi - imaginary part of the eigenvalues
288: . errest - error estimates
289: . nest - number of error estimates to display
290: - vf - viewer and format for monitoring
292: Options Database Key:
293: . -nep_monitor_conv - activates NEPMonitorConverged()
295: Level: intermediate
297: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
298: @*/
299: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
300: {
302: PetscInt i;
303: PetscViewer viewer = vf->viewer;
304: SlepcConvMon ctx;
309: ctx = (SlepcConvMon)vf->data;
310: if (its==1 && ((PetscObject)nep)->prefix) {
311: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)nep)->prefix);
312: }
313: if (its==1) ctx->oldnconv = 0;
314: if (ctx->oldnconv!=nconv) {
315: PetscViewerPushFormat(viewer,vf->format);
316: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
317: for (i=ctx->oldnconv;i<nconv;i++) {
318: PetscViewerASCIIPrintf(viewer,"%3D NEP converged value (error) #%D",its,i);
319: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
320: #if defined(PETSC_USE_COMPLEX)
321: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
322: #else
323: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
324: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
325: #endif
326: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
327: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
328: }
329: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
330: PetscViewerPopFormat(viewer);
331: ctx->oldnconv = nconv;
332: }
333: return(0);
334: }
336: PetscErrorCode NEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
337: {
339: SlepcConvMon mctx;
342: PetscViewerAndFormatCreate(viewer,format,vf);
343: PetscNew(&mctx);
344: mctx->ctx = ctx;
345: (*vf)->data = (void*)mctx;
346: return(0);
347: }
349: PetscErrorCode NEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
350: {
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: NEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
364: approximation at each iteration of the nonlinear eigensolver.
366: Collective on nep
368: Input Parameters:
369: + nep - nonlinear 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: . -nep_monitor draw::draw_lg - activates NEPMonitorFirstDrawLG()
382: Level: intermediate
384: .seealso: NEPMonitorSet()
385: @*/
386: PetscErrorCode NEPMonitorFirstDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
387: {
389: PetscViewer viewer = vf->viewer;
390: PetscDrawLG lg = vf->lg;
391: PetscReal x,y;
397: PetscViewerPushFormat(viewer,vf->format);
398: if (its==1) {
399: PetscDrawLGReset(lg);
400: PetscDrawLGSetDimension(lg,1);
401: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0);
402: }
403: if (nconv<nest) {
404: x = (PetscReal)its;
405: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
406: else y = 0.0;
407: PetscDrawLGAddPoint(lg,&x,&y);
408: if (its <= 20 || !(its % 5) || nep->reason) {
409: PetscDrawLGDraw(lg);
410: PetscDrawLGSave(lg);
411: }
412: }
413: PetscViewerPopFormat(viewer);
414: return(0);
415: }
417: /*@C
418: NEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
420: Collective on viewer
422: Input Parameters:
423: + viewer - the viewer
424: . format - the viewer format
425: - ctx - an optional user context
427: Output Parameter:
428: . vf - the viewer and format context
430: Level: intermediate
432: .seealso: NEPMonitorSet()
433: @*/
434: PetscErrorCode NEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
435: {
439: PetscViewerAndFormatCreate(viewer,format,vf);
440: (*vf)->data = ctx;
441: NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
442: return(0);
443: }
445: /*@C
446: NEPMonitorAllDrawLG - Plots the error estimate of all unconverged
447: approximations at each iteration of the nonlinear eigensolver.
449: Collective on nep
451: Input Parameters:
452: + nep - nonlinear eigensolver context
453: . its - iteration number
454: . its - iteration number
455: . nconv - number of converged eigenpairs so far
456: . eigr - real part of the eigenvalues
457: . eigi - imaginary part of the eigenvalues
458: . errest - error estimates
459: . nest - number of error estimates to display
460: - vf - viewer and format for monitoring
462: Options Database Key:
463: . -nep_monitor_all draw::draw_lg - activates NEPMonitorAllDrawLG()
465: Level: intermediate
467: .seealso: NEPMonitorSet()
468: @*/
469: PetscErrorCode NEPMonitorAllDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
470: {
472: PetscViewer viewer = vf->viewer;
473: PetscDrawLG lg = vf->lg;
474: PetscInt i,n = PetscMin(nep->nev,255);
475: PetscReal *x,*y;
481: PetscViewerPushFormat(viewer,vf->format);
482: if (its==1) {
483: PetscDrawLGReset(lg);
484: PetscDrawLGSetDimension(lg,n);
485: PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0);
486: }
487: PetscMalloc2(n,&x,n,&y);
488: for (i=0;i<n;i++) {
489: x[i] = (PetscReal)its;
490: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
491: else y[i] = 0.0;
492: }
493: PetscDrawLGAddPoint(lg,x,y);
494: if (its <= 20 || !(its % 5) || nep->reason) {
495: PetscDrawLGDraw(lg);
496: PetscDrawLGSave(lg);
497: }
498: PetscFree2(x,y);
499: PetscViewerPopFormat(viewer);
500: return(0);
501: }
503: /*@C
504: NEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
506: Collective on viewer
508: Input Parameters:
509: + viewer - the viewer
510: . format - the viewer format
511: - ctx - an optional user context
513: Output Parameter:
514: . vf - the viewer and format context
516: Level: intermediate
518: .seealso: NEPMonitorSet()
519: @*/
520: PetscErrorCode NEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
521: {
525: PetscViewerAndFormatCreate(viewer,format,vf);
526: (*vf)->data = ctx;
527: NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
528: return(0);
529: }
531: /*@C
532: NEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
533: at each iteration of the nonlinear eigensolver.
535: Collective on nep
537: Input Parameters:
538: + nep - nonlinear eigensolver context
539: . its - iteration number
540: . its - iteration number
541: . nconv - number of converged eigenpairs so far
542: . eigr - real part of the eigenvalues
543: . eigi - imaginary part of the eigenvalues
544: . errest - error estimates
545: . nest - number of error estimates to display
546: - vf - viewer and format for monitoring
548: Options Database Key:
549: . -nep_monitor_conv draw::draw_lg - activates NEPMonitorConvergedDrawLG()
551: Level: intermediate
553: .seealso: NEPMonitorSet()
554: @*/
555: PetscErrorCode NEPMonitorConvergedDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
556: {
557: PetscErrorCode ierr;
558: PetscViewer viewer = vf->viewer;
559: PetscDrawLG lg = vf->lg;
560: PetscReal x,y;
566: PetscViewerPushFormat(viewer,vf->format);
567: if (its==1) {
568: PetscDrawLGReset(lg);
569: PetscDrawLGSetDimension(lg,1);
570: PetscDrawLGSetLimits(lg,1,1,0,nep->nev);
571: }
572: x = (PetscReal)its;
573: y = (PetscReal)nep->nconv;
574: PetscDrawLGAddPoint(lg,&x,&y);
575: if (its <= 20 || !(its % 5) || nep->reason) {
576: PetscDrawLGDraw(lg);
577: PetscDrawLGSave(lg);
578: }
579: PetscViewerPopFormat(viewer);
580: return(0);
581: }
583: /*@C
584: NEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
586: Collective on viewer
588: Input Parameters:
589: + viewer - the viewer
590: . format - the viewer format
591: - ctx - an optional user context
593: Output Parameter:
594: . vf - the viewer and format context
596: Level: intermediate
598: .seealso: NEPMonitorSet()
599: @*/
600: PetscErrorCode NEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
601: {
603: SlepcConvMon mctx;
606: PetscViewerAndFormatCreate(viewer,format,vf);
607: PetscNew(&mctx);
608: mctx->ctx = ctx;
609: (*vf)->data = (void*)mctx;
610: NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
611: return(0);
612: }