Actual source code: nepmon.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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> /*I "slepcnep.h" I*/
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
23: PetscInt i,n = nep->numbermonitors;
26: for (i=0;i<n;i++) {
27: (*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]);
28: }
29: return(0);
30: }
32: /*@C
33: NEPMonitorSet - Sets an ADDITIONAL function to be called at every
34: iteration to monitor the error estimates for each requested eigenpair.
36: Logically Collective on NEP
38: Input Parameters:
39: + nep - eigensolver context obtained from NEPCreate()
40: . monitor - pointer to function (if this is NULL, it turns off monitoring)
41: . mctx - [optional] context for private data for the
42: monitor routine (use NULL if no context is desired)
43: - monitordestroy - [optional] routine that frees monitor context
44: (may be NULL)
46: Calling Sequence of monitor:
47: $ monitor(NEP nep,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)
49: + nep - nonlinear eigensolver context obtained from NEPCreate()
50: . its - iteration number
51: . nconv - number of converged eigenpairs
52: . eigr - real part of the eigenvalues
53: . eigi - imaginary part of the eigenvalues
54: . errest - error estimates for each eigenpair
55: . nest - number of error estimates
56: - mctx - optional monitoring context, as set by NEPMonitorSet()
58: Options Database Keys:
59: + -nep_monitor - print only the first error estimate
60: . -nep_monitor_all - print error estimates at each iteration
61: . -nep_monitor_conv - print the eigenvalue approximations only when
62: convergence has been reached
63: . -nep_monitor_lg - sets line graph monitor for the first unconverged
64: approximate eigenvalue
65: . -nep_monitor_lg_all - sets line graph monitor for all unconverged
66: approximate eigenvalues
67: - -nep_monitor_cancel - cancels all monitors that have been hardwired into
68: a code by calls to NEPMonitorSet(), but does not cancel those set via
69: the options database.
71: Notes:
72: Several different monitoring routines may be set by calling
73: NEPMonitorSet() multiple times; all will be called in the
74: order in which they were set.
76: Level: intermediate
78: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
79: @*/
80: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
81: {
84: if (nep->numbermonitors >= MAXNEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
85: nep->monitor[nep->numbermonitors] = monitor;
86: nep->monitorcontext[nep->numbermonitors] = (void*)mctx;
87: nep->monitordestroy[nep->numbermonitors++] = monitordestroy;
88: return(0);
89: }
91: /*@
92: NEPMonitorCancel - Clears all monitors for a NEP object.
94: Logically Collective on NEP
96: Input Parameters:
97: . nep - eigensolver context obtained from NEPCreate()
99: Options Database Key:
100: . -nep_monitor_cancel - Cancels all monitors that have been hardwired
101: into a code by calls to NEPMonitorSet(),
102: but does not cancel those set via the options database.
104: Level: intermediate
106: .seealso: NEPMonitorSet()
107: @*/
108: PetscErrorCode NEPMonitorCancel(NEP nep)
109: {
111: PetscInt i;
115: for (i=0; i<nep->numbermonitors; i++) {
116: if (nep->monitordestroy[i]) {
117: (*nep->monitordestroy[i])(&nep->monitorcontext[i]);
118: }
119: }
120: nep->numbermonitors = 0;
121: return(0);
122: }
124: /*@C
125: NEPGetMonitorContext - Gets the monitor context, as set by
126: NEPMonitorSet() for the FIRST monitor only.
128: Not Collective
130: Input Parameter:
131: . nep - eigensolver context obtained from NEPCreate()
133: Output Parameter:
134: . ctx - monitor context
136: Level: intermediate
138: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
139: @*/
140: PetscErrorCode NEPGetMonitorContext(NEP nep,void **ctx)
141: {
144: *ctx = nep->monitorcontext[0];
145: return(0);
146: }
148: /*@C
149: NEPMonitorAll - Print the current approximate values and
150: error estimates at each iteration of the nonlinear eigensolver.
152: Collective on NEP
154: Input Parameters:
155: + nep - nonlinear eigensolver context
156: . its - iteration number
157: . nconv - number of converged eigenpairs so far
158: . eigr - real part of the eigenvalues
159: . eigi - imaginary part of the eigenvalues
160: . errest - error estimates
161: . nest - number of error estimates to display
162: - vf - viewer and format for monitoring
164: Level: intermediate
166: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
167: @*/
168: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
169: {
171: PetscInt i;
172: PetscViewer viewer;
177: viewer = vf->viewer;
179: PetscViewerPushFormat(viewer,vf->format);
180: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
181: if (its==1 && ((PetscObject)nep)->prefix) {
182: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
183: }
184: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D Values (Errors)",its,nconv);
185: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
186: for (i=0;i<nest;i++) {
187: #if defined(PETSC_USE_COMPLEX)
188: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
189: #else
190: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
191: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
192: #endif
193: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
194: }
195: PetscViewerASCIIPrintf(viewer,"\n");
196: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
197: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
198: PetscViewerPopFormat(viewer);
199: return(0);
200: }
202: /*@C
203: NEPMonitorFirst - Print the first unconverged approximate value and
204: error estimate at each iteration of the nonlinear eigensolver.
206: Collective on NEP
208: Input Parameters:
209: + nep - nonlinear eigensolver context
210: . its - iteration number
211: . nconv - number of converged eigenpairs so far
212: . eigr - real part of the eigenvalues
213: . eigi - imaginary part of the eigenvalues
214: . errest - error estimates
215: . nest - number of error estimates to display
216: - vf - viewer and format for monitoring
218: Level: intermediate
220: .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
221: @*/
222: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
223: {
225: PetscViewer viewer;
230: viewer = vf->viewer;
232: if (its==1 && ((PetscObject)nep)->prefix) {
233: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
234: }
235: if (nconv<nest) {
236: PetscViewerPushFormat(viewer,vf->format);
237: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
238: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D first unconverged value (error)",its,nconv);
239: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
240: #if defined(PETSC_USE_COMPLEX)
241: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv]));
242: #else
243: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]);
244: if (eigi[nconv]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]); }
245: #endif
246: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
247: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
248: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
249: PetscViewerPopFormat(viewer);
250: }
251: return(0);
252: }
254: /*@C
255: NEPMonitorConverged - Print the approximate values and
256: error estimates as they converge.
258: Collective on NEP
260: Input Parameters:
261: + nep - nonlinear eigensolver context
262: . its - iteration number
263: . nconv - number of converged eigenpairs so far
264: . eigr - real part of the eigenvalues
265: . eigi - imaginary part of the eigenvalues
266: . errest - error estimates
267: . nest - number of error estimates to display
268: - ctx - monitor context
270: Level: intermediate
272: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
273: @*/
274: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,SlepcConvMonitor ctx)
275: {
277: PetscInt i;
278: PetscViewer viewer;
283: viewer = ctx->viewer;
285: if (its==1 && ((PetscObject)nep)->prefix) {
286: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)nep)->prefix);
287: }
288: if (its==1) ctx->oldnconv = 0;
289: if (ctx->oldnconv!=nconv) {
290: PetscViewerPushFormat(viewer,ctx->format);
291: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
292: for (i=ctx->oldnconv;i<nconv;i++) {
293: PetscViewerASCIIPrintf(viewer,"%3D NEP converged value (error) #%D",its,i);
294: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
295: #if defined(PETSC_USE_COMPLEX)
296: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
297: #else
298: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
299: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
300: #endif
301: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
302: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
303: }
304: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
305: PetscViewerPopFormat(viewer);
306: ctx->oldnconv = nconv;
307: }
308: return(0);
309: }
311: /*@C
312: NEPMonitorLGCreate - Creates a line graph context for use with
313: NEP to monitor convergence.
315: Collective on MPI_Comm
317: Input Parameters:
318: + comm - communicator context
319: . host - the X display to open, or null for the local machine
320: . label - the title to put in the title bar
321: . x, y - the screen coordinates of the upper left coordinate of
322: the window
323: - m, n - the screen width and height in pixels
325: Output Parameter:
326: . lgctx - the drawing context
328: Options Database Keys:
329: + -nep_monitor_lg - Sets line graph monitor for the first residual
330: - -nep_monitor_lg_all - Sets line graph monitor for all residuals
332: Notes:
333: Use PetscDrawLGDestroy() to destroy this line graph.
335: Level: intermediate
337: .seealso: NEPMonitorSet()
338: @*/
339: PetscErrorCode NEPMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
340: {
341: PetscDraw draw;
342: PetscDrawLG lg;
346: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
347: PetscDrawSetFromOptions(draw);
348: PetscDrawLGCreate(draw,1,&lg);
349: PetscDrawLGSetFromOptions(lg);
350: PetscDrawDestroy(&draw);
351: *lgctx = lg;
352: return(0);
353: }
355: PetscErrorCode NEPMonitorLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
356: {
357: PetscDrawLG lg = (PetscDrawLG)ctx;
358: PetscReal x,y;
363: if (its==1) {
364: PetscDrawLGReset(lg);
365: PetscDrawLGSetDimension(lg,1);
366: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(nep->tol)-2,0.0);
367: }
368: x = (PetscReal)its;
369: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
370: else y = 0.0;
371: PetscDrawLGAddPoint(lg,&x,&y);
372: if (its <= 20 || !(its % 5) || nep->reason) {
373: PetscDrawLGDraw(lg);
374: PetscDrawLGSave(lg);
375: }
376: return(0);
377: }
379: PetscErrorCode NEPMonitorLGAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
380: {
381: PetscDrawLG lg = (PetscDrawLG)ctx;
382: PetscInt i,n = PetscMin(nep->nev,255);
383: PetscReal *x,*y;
388: if (its==1) {
389: PetscDrawLGReset(lg);
390: PetscDrawLGSetDimension(lg,n);
391: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(nep->tol)-2,0.0);
392: }
393: PetscMalloc2(n,&x,n,&y);
394: for (i=0;i<n;i++) {
395: x[i] = (PetscReal)its;
396: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
397: else y[i] = 0.0;
398: }
399: PetscDrawLGAddPoint(lg,x,y);
400: if (its <= 20 || !(its % 5) || nep->reason) {
401: PetscDrawLGDraw(lg);
402: PetscDrawLGSave(lg);
403: }
404: PetscFree2(x,y);
405: return(0);
406: }