Actual source code: monitor.c
1: /*
2: EPS routines related to monitors.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2011, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
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 <private/epsimpl.h> /*I "slepceps.h" I*/
28: /*
29: Runs the user provided monitor routines, if any.
30: */
31: PetscErrorCode EPSMonitor(EPS eps,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
32: {
34: PetscInt i,n = eps->numbermonitors;
37: for (i=0;i<n;i++) {
38: (*eps->monitor[i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[i]);
39: }
40: return(0);
41: }
45: /*@C
46: EPSMonitorSet - Sets an ADDITIONAL function to be called at every
47: iteration to monitor the error estimates for each requested eigenpair.
48:
49: Logically Collective on EPS
51: Input Parameters:
52: + eps - eigensolver context obtained from EPSCreate()
53: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
54: . mctx - [optional] context for private data for the
55: monitor routine (use PETSC_NULL if no context is desired)
56: - monitordestroy - [optional] routine that frees monitor context
57: (may be PETSC_NULL)
59: Calling Sequence of monitor:
60: $ monitor (EPS eps, int its, int nconv, PetscScalar *eigr, PetscScalar *eigi, PetscReal* errest, int nest, void *mctx)
62: + eps - eigensolver context obtained from EPSCreate()
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 EPSMonitorSet()
71: Options Database Keys:
72: + -eps_monitor - print only the first error estimate
73: . -eps_monitor_all - print error estimates at each iteration
74: . -eps_monitor_conv - print the eigenvalue approximations only when
75: convergence has been reached
76: . -eps_monitor_draw - sets line graph monitor for the first unconverged
77: approximate eigenvalue
78: . -eps_monitor_draw_all - sets line graph monitor for all unconverged
79: approximate eigenvalue
80: - -eps_monitor_cancel - cancels all monitors that have been hardwired into
81: a code by calls to EPSMonitorSet(), but does not cancel those set via
82: the options database.
84: Notes:
85: Several different monitoring routines may be set by calling
86: EPSMonitorSet() multiple times; all will be called in the
87: order in which they were set.
89: Level: intermediate
91: .seealso: EPSMonitorFirst(), EPSMonitorAll(), EPSMonitorLG(), EPSMonitorLGAll(), EPSMonitorCancel()
92: @*/
93: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
94: {
97: if (eps->numbermonitors >= MAXEPSMONITORS) {
98: SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
99: }
100: eps->monitor[eps->numbermonitors] = monitor;
101: eps->monitorcontext[eps->numbermonitors] = (void*)mctx;
102: eps->monitordestroy[eps->numbermonitors++] = monitordestroy;
103: return(0);
104: }
108: /*@
109: EPSMonitorCancel - Clears all monitors for an EPS object.
111: Logically Collective on EPS
113: Input Parameters:
114: . eps - eigensolver context obtained from EPSCreate()
116: Options Database Key:
117: . -eps_monitor_cancel - Cancels all monitors that have been hardwired
118: into a code by calls to EPSMonitorSet(),
119: but does not cancel those set via the options database.
121: Level: intermediate
123: .seealso: EPSMonitorSet()
124: @*/
125: PetscErrorCode EPSMonitorCancel(EPS eps)
126: {
128: PetscInt i;
132: for (i=0; i<eps->numbermonitors; i++) {
133: if (eps->monitordestroy[i]) {
134: (*eps->monitordestroy[i])(&eps->monitorcontext[i]);
135: }
136: }
137: eps->numbermonitors = 0;
138: return(0);
139: }
143: /*@C
144: EPSGetMonitorContext - Gets the monitor context, as set by
145: EPSMonitorSet() for the FIRST monitor only.
147: Not Collective
149: Input Parameter:
150: . eps - eigensolver context obtained from EPSCreate()
152: Output Parameter:
153: . ctx - monitor context
155: Level: intermediate
157: .seealso: EPSMonitorSet()
158: @*/
159: PetscErrorCode EPSGetMonitorContext(EPS eps,void **ctx)
160: {
163: *ctx = eps->monitorcontext[0];
164: return(0);
165: }
169: /*@C
170: EPSMonitorAll - Print the current approximate values and
171: error estimates at each iteration of the eigensolver.
173: Collective on EPS
175: Input Parameters:
176: + eps - eigensolver context
177: . its - iteration number
178: . nconv - number of converged eigenpairs so far
179: . eigr - real part of the eigenvalues
180: . eigi - imaginary part of the eigenvalues
181: . errest - error estimates
182: . nest - number of error estimates to display
183: - dummy - unused monitor context
185: Level: intermediate
187: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorConverged()
188: @*/
189: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
190: {
192: PetscInt i;
193: PetscScalar er,ei;
194: PetscViewer viewer = dummy? (PetscViewer)dummy: PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);
197: if (its) {
198: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
199: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D Values (Errors)",its,nconv);
200: for (i=0;i<nest;i++) {
201: er = eigr[i]; ei = eigi[i];
202: STBackTransform(eps->OP,1,&er,&ei);
203: #if defined(PETSC_USE_COMPLEX)
204: PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(er),PetscImaginaryPart(er));
205: #else
206: PetscViewerASCIIPrintf(viewer," %G",er);
207: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",ei); }
208: #endif
209: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
210: }
211: PetscViewerASCIIPrintf(viewer,"\n");
212: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
213: }
214: return(0);
215: }
219: /*@C
220: EPSMonitorFirst - Print the first approximate value and
221: error estimate at each iteration of the eigensolver.
223: Collective on EPS
225: Input Parameters:
226: + eps - eigensolver context
227: . its - iteration number
228: . nconv - number of converged eigenpairs so far
229: . eigr - real part of the eigenvalues
230: . eigi - imaginary part of the eigenvalues
231: . errest - error estimates
232: . nest - number of error estimates to display
233: - dummy - unused monitor context
235: Level: intermediate
237: .seealso: EPSMonitorSet(), EPSMonitorAll(), EPSMonitorConverged()
238: @*/
239: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
240: {
242: PetscScalar er,ei;
243: PetscViewer viewer = dummy? (PetscViewer)dummy: PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);
246: if (its && nconv<nest) {
247: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
248: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D first unconverged value (error)",its,nconv);
249: er = eigr[nconv]; ei = eigi[nconv];
250: STBackTransform(eps->OP,1,&er,&ei);
251: #if defined(PETSC_USE_COMPLEX)
252: PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(er),PetscImaginaryPart(er));
253: #else
254: PetscViewerASCIIPrintf(viewer," %G",er);
255: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",ei); }
256: #endif
257: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
258: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
259: }
260: return(0);
261: }
265: /*@C
266: EPSMonitorConverged - Print the approximate values and
267: error estimates as they converge.
269: Collective on EPS
271: Input Parameters:
272: + eps - eigensolver context
273: . its - iteration number
274: . nconv - number of converged eigenpairs so far
275: . eigr - real part of the eigenvalues
276: . eigi - imaginary part of the eigenvalues
277: . errest - error estimates
278: . nest - number of error estimates to display
279: - dummy - unused monitor context
281: Level: intermediate
283: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorAll()
284: @*/
285: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
286: {
287: PetscErrorCode ierr;
288: PetscInt i;
289: PetscScalar er,ei;
290: SlepcConvMonitor ctx = (SlepcConvMonitor)dummy;
293: if (!its) {
294: ctx->oldnconv = 0;
295: } else {
296: for (i=ctx->oldnconv;i<nconv;i++) {
297: PetscViewerASCIIAddTab(ctx->viewer,((PetscObject)eps)->tablevel);
298: PetscViewerASCIIPrintf(ctx->viewer,"%3D EPS converged value (error) #%D",its,i);
299: er = eigr[i]; ei = eigi[i];
300: STBackTransform(eps->OP,1,&er,&ei);
301: #if defined(PETSC_USE_COMPLEX)
302: PetscViewerASCIIPrintf(ctx->viewer," %G%+Gi",PetscRealPart(er),PetscImaginaryPart(er));
303: #else
304: PetscViewerASCIIPrintf(ctx->viewer," %G",er);
305: if (ei!=0.0) { PetscViewerASCIIPrintf(ctx->viewer,"%+Gi",ei); }
306: #endif
307: PetscViewerASCIIPrintf(ctx->viewer," (%10.8e)\n",(double)errest[i]);
308: PetscViewerASCIISubtractTab(ctx->viewer,((PetscObject)eps)->tablevel);
309: }
310: ctx->oldnconv = nconv;
311: }
312: return(0);
313: }
317: PetscErrorCode EPSMonitorLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
318: {
319: PetscViewer viewer = (PetscViewer) monctx;
320: PetscDraw draw,draw1;
321: PetscDrawLG lg,lg1;
323: PetscReal x,y,myeigr,p;
324: PetscScalar er,ei;
327: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); }
328: PetscViewerDrawGetDraw(viewer,0,&draw);
329: PetscViewerDrawGetDrawLG(viewer,0,&lg);
330: if (!its) {
331: PetscDrawSetTitle(draw,"Error estimates");
332: PetscDrawSetDoubleBuffer(draw);
333: PetscDrawLGSetDimension(lg,1);
334: PetscDrawLGReset(lg);
335: PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
336: }
338: /* In the hermitian case, the eigenvalues are real and can be plotted */
339: if (eps->ishermitian) {
340: PetscViewerDrawGetDraw(viewer,1,&draw1);
341: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
342: if (!its) {
343: PetscDrawSetTitle(draw1,"Approximate eigenvalues");
344: PetscDrawSetDoubleBuffer(draw1);
345: PetscDrawLGSetDimension(lg1,1);
346: PetscDrawLGReset(lg1);
347: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
348: }
349: }
351: x = (PetscReal) its;
352: if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
353: PetscDrawLGAddPoint(lg,&x,&y);
354: if (eps->ishermitian) {
355: er = eigr[nconv]; ei = eigi[nconv];
356: STBackTransform(eps->OP,1,&er,&ei);
357: myeigr = PetscRealPart(er);
358: PetscDrawLGAddPoint(lg1,&x,&myeigr);
359: PetscDrawGetPause(draw1,&p);
360: PetscDrawSetPause(draw1,0);
361: PetscDrawLGDraw(lg1);
362: PetscDrawSetPause(draw1,p);
363: }
364: PetscDrawLGDraw(lg);
365: return(0);
366: }
370: PetscErrorCode EPSMonitorLGAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
371: {
372: PetscViewer viewer = (PetscViewer) monctx;
373: PetscDraw draw,draw1;
374: PetscDrawLG lg,lg1;
376: PetscReal *x,*y,*myeigr,p;
377: PetscInt i,n = PetscMin(eps->nev,255);
378: PetscScalar er,ei;
381: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); }
382: PetscViewerDrawGetDraw(viewer,0,&draw);
383: PetscViewerDrawGetDrawLG(viewer,0,&lg);
384: if (!its) {
385: PetscDrawSetTitle(draw,"Error estimates");
386: PetscDrawSetDoubleBuffer(draw);
387: PetscDrawLGSetDimension(lg,n);
388: PetscDrawLGReset(lg);
389: PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
390: }
392: /* In the hermitian case, the eigenvalues are real and can be plotted */
393: if (eps->ishermitian) {
394: PetscViewerDrawGetDraw(viewer,1,&draw1);
395: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
396: if (!its) {
397: PetscDrawSetTitle(draw1,"Approximate eigenvalues");
398: PetscDrawSetDoubleBuffer(draw1);
399: PetscDrawLGSetDimension(lg1,n);
400: PetscDrawLGReset(lg1);
401: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
402: }
403: }
405: PetscMalloc(sizeof(PetscReal)*n,&x);
406: PetscMalloc(sizeof(PetscReal)*n,&y);
407: for (i=0;i<n;i++) {
408: x[i] = (PetscReal) its;
409: if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
410: else y[i] = 0.0;
411: }
412: PetscDrawLGAddPoint(lg,x,y);
413: if (eps->ishermitian) {
414: PetscMalloc(sizeof(PetscReal)*n,&myeigr);
415: for(i=0;i<n;i++) {
416: er = eigr[i]; ei = eigi[i];
417: STBackTransform(eps->OP,1,&er,&ei);
418: if (i < nest) myeigr[i] = PetscRealPart(er);
419: else myeigr[i] = 0.0;
420: }
421: PetscDrawLGAddPoint(lg1,x,myeigr);
422: PetscDrawGetPause(draw1,&p);
423: PetscDrawSetPause(draw1,0);
424: PetscDrawLGDraw(lg1);
425: PetscDrawSetPause(draw1,p);
426: PetscFree(myeigr);
427: }
428: PetscDrawLGDraw(lg);
429: PetscFree(x);
430: PetscFree(y);
431: return(0);
432: }