1: /*
2: NEP routines related to monitors.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/
25: #include <petscdraw.h>
29: /*
30: Runs the user provided monitor routines, if any.
31: */
32: PetscErrorCode NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest) 33: {
35: PetscInt i,n = nep->numbermonitors;
38: for (i=0;i<n;i++) {
39: (*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]);
40: }
41: return(0);
42: }
46: /*@C
47: NEPMonitorSet - Sets an ADDITIONAL function to be called at every
48: iteration to monitor the error estimates for each requested eigenpair.
50: Logically Collective on NEP 52: Input Parameters:
53: + nep - eigensolver context obtained from NEPCreate()
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
58: (may be NULL)
60: Calling Sequence of monitor:
61: $ monitor(NEP nep,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)
63: + nep - nonlinear eigensolver context obtained from NEPCreate()
64: . its - iteration number
65: . nconv - number of converged eigenpairs
66: . eigr - real part of the eigenvalues
67: . eigi - imaginary part of the eigenvalues
68: . errest - error estimates for each eigenpair
69: . nest - number of error estimates
70: - mctx - optional monitoring context, as set by NEPMonitorSet()
72: Options Database Keys:
73: + -nep_monitor - print only the first error estimate
74: . -nep_monitor_all - print error estimates at each iteration
75: . -nep_monitor_conv - print the eigenvalue approximations only when
76: convergence has been reached
77: . -nep_monitor_lg - sets line graph monitor for the first unconverged
78: approximate eigenvalue
79: . -nep_monitor_lg_all - sets line graph monitor for all unconverged
80: approximate eigenvalues
81: - -nep_monitor_cancel - cancels all monitors that have been hardwired into
82: a code by calls to NEPMonitorSet(), but does not cancel those set via
83: the options database.
85: Notes:
86: Several different monitoring routines may be set by calling
87: NEPMonitorSet() multiple times; all will be called in the
88: order in which they were set.
90: Level: intermediate
92: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
93: @*/
94: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**)) 95: {
98: if (nep->numbermonitors >= MAXNEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
99: nep->monitor[nep->numbermonitors] = monitor;
100: nep->monitorcontext[nep->numbermonitors] = (void*)mctx;
101: nep->monitordestroy[nep->numbermonitors++] = monitordestroy;
102: return(0);
103: }
107: /*@
108: NEPMonitorCancel - Clears all monitors for a NEP object.
110: Logically Collective on NEP112: Input Parameters:
113: . nep - eigensolver context obtained from NEPCreate()
115: Options Database Key:
116: . -nep_monitor_cancel - Cancels all monitors that have been hardwired
117: into a code by calls to NEPMonitorSet(),
118: but does not cancel those set via the options database.
120: Level: intermediate
122: .seealso: NEPMonitorSet()
123: @*/
124: PetscErrorCode NEPMonitorCancel(NEP nep)125: {
127: PetscInt i;
131: for (i=0; i<nep->numbermonitors; i++) {
132: if (nep->monitordestroy[i]) {
133: (*nep->monitordestroy[i])(&nep->monitorcontext[i]);
134: }
135: }
136: nep->numbermonitors = 0;
137: return(0);
138: }
142: /*@C
143: NEPGetMonitorContext - Gets the monitor context, as set by
144: NEPMonitorSet() for the FIRST monitor only.
146: Not Collective
148: Input Parameter:
149: . nep - eigensolver context obtained from NEPCreate()
151: Output Parameter:
152: . ctx - monitor context
154: Level: intermediate
156: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
157: @*/
158: PetscErrorCode NEPGetMonitorContext(NEP nep,void **ctx)159: {
162: *ctx = nep->monitorcontext[0];
163: return(0);
164: }
168: /*@C
169: NEPMonitorAll - Print the current approximate values and
170: error estimates at each iteration of the nonlinear eigensolver.
172: Collective on NEP174: 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: Level: intermediate
186: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
187: @*/
188: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)189: {
191: PetscInt i;
192: PetscViewer viewer;
197: viewer = vf->viewer;
199: PetscViewerPushFormat(viewer,vf->format);
200: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
201: if (its==1 && ((PetscObject)nep)->prefix) {
202: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
203: }
204: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D Values (Errors)",its,nconv);
205: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
206: for (i=0;i<nest;i++) {
207: #if defined(PETSC_USE_COMPLEX)
208: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
209: #else
210: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
211: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
212: #endif
213: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
214: }
215: PetscViewerASCIIPrintf(viewer,"\n");
216: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
217: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
218: PetscViewerPopFormat(viewer);
219: return(0);
220: }
224: /*@C
225: NEPMonitorFirst - Print the first unconverged approximate value and
226: error estimate at each iteration of the nonlinear eigensolver.
228: Collective on NEP230: Input Parameters:
231: + nep - nonlinear eigensolver context
232: . its - iteration number
233: . nconv - number of converged eigenpairs so far
234: . eigr - real part of the eigenvalues
235: . eigi - imaginary part of the eigenvalues
236: . errest - error estimates
237: . nest - number of error estimates to display
238: - vf - viewer and format for monitoring
240: Level: intermediate
242: .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
243: @*/
244: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)245: {
247: PetscViewer viewer;
252: viewer = vf->viewer;
254: if (its==1 && ((PetscObject)nep)->prefix) {
255: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
256: }
257: if (nconv<nest) {
258: PetscViewerPushFormat(viewer,vf->format);
259: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
260: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D first unconverged value (error)",its,nconv);
261: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
262: #if defined(PETSC_USE_COMPLEX)
263: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv]));
264: #else
265: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]);
266: if (eigi[nconv]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]); }
267: #endif
268: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
269: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
270: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
271: PetscViewerPopFormat(viewer);
272: }
273: return(0);
274: }
278: /*@C
279: NEPMonitorConverged - Print the approximate values and
280: error estimates as they converge.
282: Collective on NEP284: Input Parameters:
285: + nep - nonlinear eigensolver context
286: . its - iteration number
287: . nconv - number of converged eigenpairs so far
288: . eigr - real part of the eigenvalues
289: . eigi - imaginary part of the eigenvalues
290: . errest - error estimates
291: . nest - number of error estimates to display
292: - ctx - monitor context
294: Level: intermediate
296: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
297: @*/
298: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,SlepcConvMonitor ctx)299: {
301: PetscInt i;
302: PetscViewer viewer;
307: viewer = ctx->viewer;
309: if (its==1 && ((PetscObject)nep)->prefix) {
310: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)nep)->prefix);
311: }
312: if (its==1) ctx->oldnconv = 0;
313: if (ctx->oldnconv!=nconv) {
314: PetscViewerPushFormat(viewer,ctx->format);
315: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
316: for (i=ctx->oldnconv;i<nconv;i++) {
317: PetscViewerASCIIPrintf(viewer,"%3D NEP converged value (error) #%D",its,i);
318: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
319: #if defined(PETSC_USE_COMPLEX)
320: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
321: #else
322: PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
323: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
324: #endif
325: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
326: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
327: }
328: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
329: PetscViewerPopFormat(viewer);
330: ctx->oldnconv = nconv;
331: }
332: return(0);
333: }
337: /*@C
338: NEPMonitorLGCreate - Creates a line graph context for use with
339: NEP to monitor convergence.
341: Collective on MPI_Comm
343: Input Parameters:
344: + comm - communicator context
345: . host - the X display to open, or null for the local machine
346: . label - the title to put in the title bar
347: . x, y - the screen coordinates of the upper left coordinate of
348: the window
349: - m, n - the screen width and height in pixels
351: Output Parameter:
352: . lgctx - the drawing context
354: Options Database Keys:
355: + -nep_monitor_lg - Sets line graph monitor for the first residual
356: - -nep_monitor_lg_all - Sets line graph monitor for all residuals
358: Notes:
359: Use PetscDrawLGDestroy() to destroy this line graph.
361: Level: intermediate
363: .seealso: NEPMonitorSet()
364: @*/
365: PetscErrorCode NEPMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)366: {
367: PetscDraw draw;
368: PetscDrawLG lg;
372: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
373: PetscDrawSetFromOptions(draw);
374: PetscDrawLGCreate(draw,1,&lg);
375: PetscDrawLGSetFromOptions(lg);
376: PetscDrawDestroy(&draw);
377: *lgctx = lg;
378: return(0);
379: }
383: PetscErrorCode NEPMonitorLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)384: {
385: PetscDrawLG lg = (PetscDrawLG)ctx;
386: PetscReal x,y;
391: if (its==1) {
392: PetscDrawLGReset(lg);
393: PetscDrawLGSetDimension(lg,1);
394: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(nep->tol)-2,0.0);
395: }
396: x = (PetscReal)its;
397: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
398: else y = 0.0;
399: PetscDrawLGAddPoint(lg,&x,&y);
400: if (its <= 20 || !(its % 5) || nep->reason) {
401: PetscDrawLGDraw(lg);
402: PetscDrawLGSave(lg);
403: }
404: return(0);
405: }
409: PetscErrorCode NEPMonitorLGAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)410: {
411: PetscDrawLG lg = (PetscDrawLG)ctx;
412: PetscInt i,n = PetscMin(nep->nev,255);
413: PetscReal *x,*y;
418: if (its==1) {
419: PetscDrawLGReset(lg);
420: PetscDrawLGSetDimension(lg,n);
421: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(nep->tol)-2,0.0);
422: }
423: PetscMalloc2(n,&x,n,&y);
424: for (i=0;i<n;i++) {
425: x[i] = (PetscReal)its;
426: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
427: else y[i] = 0.0;
428: }
429: PetscDrawLGAddPoint(lg,x,y);
430: if (its <= 20 || !(its % 5) || nep->reason) {
431: PetscDrawLGDraw(lg);
432: PetscDrawLGSave(lg);
433: }
434: PetscFree2(x,y);
435: return(0);
436: }