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-2010, Universidad 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
28: /*@C
29: EPSMonitorSet - Sets an ADDITIONAL function to be called at every
30: iteration to monitor the error estimates for each requested eigenpair.
31:
32: Collective on EPS
34: Input Parameters:
35: + eps - eigensolver context obtained from EPSCreate()
36: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
37: . mctx - [optional] context for private data for the
38: monitor routine (use PETSC_NULL if no context is desired)
39: - monitordestroy - [optional] routine that frees monitor context
40: (may be PETSC_NULL)
42: Calling Sequence of monitor:
43: $ monitor (EPS eps, int its, int nconv, PetscScalar *eigr, PetscScalar *eigi, PetscReal* errest, int nest, void *mctx)
45: + eps - eigensolver context obtained from EPSCreate()
46: . its - iteration number
47: . nconv - number of converged eigenpairs
48: . eigr - real part of the eigenvalues
49: . eigi - imaginary part of the eigenvalues
50: . errest - relative error estimates for each eigenpair
51: . nest - number of error estimates
52: - mctx - optional monitoring context, as set by EPSMonitorSet()
54: Options Database Keys:
55: + -eps_monitor - print only the first error estimate
56: . -eps_monitor_all - print error estimates at each iteration
57: . -eps_monitor_conv - print the eigenvalue approximations only when
58: convergence has been reached
59: . -eps_monitor_draw - sets line graph monitor for the first unconverged
60: approximate eigenvalue
61: . -eps_monitor_draw_all - sets line graph monitor for all unconverged
62: approximate eigenvalue
63: - -eps_monitor_cancel - cancels all monitors that have been hardwired into
64: a code by calls to EPSMonitorSet(), but does not cancel those set via
65: the options database.
67: Notes:
68: Several different monitoring routines may be set by calling
69: EPSMonitorSet() multiple times; all will be called in the
70: order in which they were set.
72: Level: intermediate
74: .seealso: EPSMonitorFirst(), EPSMonitorAll(), EPSMonitorLG(), EPSMonitorLGAll(), EPSMonitorCancel()
75: @*/
76: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),
77: void *mctx,PetscErrorCode (*monitordestroy)(void*))
78: {
81: if (eps->numbermonitors >= MAXEPSMONITORS) {
82: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
83: }
84: eps->monitor[eps->numbermonitors] = monitor;
85: eps->monitorcontext[eps->numbermonitors] = (void*)mctx;
86: eps->monitordestroy[eps->numbermonitors++] = monitordestroy;
87: return(0);
88: }
92: /*@
93: EPSMonitorCancel - Clears all monitors for an EPS object.
95: Collective on EPS
97: Input Parameters:
98: . eps - eigensolver context obtained from EPSCreate()
100: Options Database Key:
101: . -eps_monitor_cancel - Cancels all monitors that have been hardwired
102: into a code by calls to EPSMonitorSet(),
103: but does not cancel those set via the options database.
105: Level: intermediate
107: .seealso: EPSMonitorSet()
108: @*/
109: PetscErrorCode EPSMonitorCancel(EPS eps)
110: {
112: PetscInt i;
116: for (i=0; i<eps->numbermonitors; i++) {
117: if (eps->monitordestroy[i]) {
118: (*eps->monitordestroy[i])(eps->monitorcontext[i]);
119: }
120: }
121: eps->numbermonitors = 0;
122: return(0);
123: }
127: /*@C
128: EPSGetMonitorContext - Gets the monitor context, as set by
129: EPSSetMonitor() for the FIRST monitor only.
131: Not Collective
133: Input Parameter:
134: . eps - eigensolver context obtained from EPSCreate()
136: Output Parameter:
137: . ctx - monitor context
139: Level: intermediate
141: .seealso: EPSSetMonitor()
142: @*/
143: PetscErrorCode EPSGetMonitorContext(EPS eps, void **ctx)
144: {
147: *ctx = (eps->monitorcontext[0]);
148: return(0);
149: }
153: /*@C
154: EPSMonitorAll - Print the current approximate values and
155: error estimates at each iteration of the eigensolver.
157: Collective on EPS
159: Input Parameters:
160: + eps - eigensolver context
161: . its - iteration number
162: . nconv - number of converged eigenpairs so far
163: . eigr - real part of the eigenvalues
164: . eigi - imaginary part of the eigenvalues
165: . errest - error estimates
166: . nest - number of error estimates to display
167: - dummy - unused monitor context
169: Level: intermediate
171: .seealso: EPSMonitorSet()
172: @*/
173: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
174: {
175: PetscErrorCode ierr;
176: PetscInt i;
177: PetscScalar er,ei;
178: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
181: if (its) {
182: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)eps)->comm,"stdout",0,&viewer);}
183: PetscViewerASCIIMonitorPrintf(viewer,"%3d EPS nconv=%d Values (Errors)",its,nconv);
184: for (i=0;i<nest;i++) {
185: er = eigr[i]; ei = eigi[i];
186: STBackTransform(eps->OP, 1, &er, &ei);
187: #if defined(PETSC_USE_COMPLEX)
188: PetscViewerASCIIMonitorPrintf(viewer," %g%+gi",PetscRealPart(er),PetscImaginaryPart(er));
189: #else
190: PetscViewerASCIIMonitorPrintf(viewer," %g",er);
191: if (ei!=0.0) { PetscViewerASCIIMonitorPrintf(viewer,"%+gi",ei); }
192: #endif
193: PetscViewerASCIIMonitorPrintf(viewer," (%10.8e)",errest[i]);
194: }
195: PetscViewerASCIIMonitorPrintf(viewer,"\n");
196: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
197: }
198: return(0);
199: }
203: /*@C
204: EPSMonitorFirst - Print the first approximate value and
205: error estimate at each iteration of the eigensolver.
207: Collective on EPS
209: Input Parameters:
210: + eps - eigensolver context
211: . its - iteration number
212: . nconv - number of converged eigenpairs so far
213: . eigr - real part of the eigenvalues
214: . eigi - imaginary part of the eigenvalues
215: . errest - error estimates
216: . nest - number of error estimates to display
217: - dummy - unused monitor context
219: Level: intermediate
221: .seealso: EPSMonitorSet()
222: @*/
223: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
224: {
225: PetscErrorCode ierr;
226: PetscScalar er,ei;
227: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
230: if (its && nconv<nest) {
231: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)eps)->comm,"stdout",0,&viewer);}
232: PetscViewerASCIIMonitorPrintf(viewer,"%3d EPS nconv=%d first unconverged value (error)",its,nconv);
233: er = eigr[nconv]; ei = eigi[nconv];
234: STBackTransform(eps->OP, 1, &er, &ei);
235: #if defined(PETSC_USE_COMPLEX)
236: PetscViewerASCIIMonitorPrintf(viewer," %g%+gi",PetscRealPart(er),PetscImaginaryPart(er));
237: #else
238: PetscViewerASCIIMonitorPrintf(viewer," %g",er);
239: if (ei!=0.0) { PetscViewerASCIIMonitorPrintf(viewer,"%+gi",ei); }
240: #endif
241: PetscViewerASCIIMonitorPrintf(viewer," (%10.8e)\n",errest[nconv]);
242: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
243: }
244: return(0);
245: }
249: /*@C
250: EPSMonitorConverged - Print the approximate values and
251: error estimates as they converge.
253: Collective on EPS
255: Input Parameters:
256: + eps - eigensolver context
257: . its - iteration number
258: . nconv - number of converged eigenpairs so far
259: . eigr - real part of the eigenvalues
260: . eigi - imaginary part of the eigenvalues
261: . errest - error estimates
262: . nest - number of error estimates to display
263: - dummy - unused monitor context
265: Level: intermediate
267: .seealso: EPSMonitorSet()
268: @*/
269: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
270: {
271: PetscErrorCode ierr;
272: PetscInt i;
273: PetscScalar er,ei;
274: EPSMONITOR_CONV *ctx = (EPSMONITOR_CONV*) dummy;
277: if (!its) {
278: ctx->oldnconv = 0;
279: } else {
280: for (i=ctx->oldnconv;i<nconv;i++) {
281: PetscViewerASCIIMonitorPrintf(ctx->viewer,"%3d EPS converged value (error) #%d",its,i);
282: er = eigr[i]; ei = eigi[i];
283: STBackTransform(eps->OP, 1, &er, &ei);
284: #if defined(PETSC_USE_COMPLEX)
285: PetscViewerASCIIMonitorPrintf(ctx->viewer," %g%+gi",PetscRealPart(er),PetscImaginaryPart(er));
286: #else
287: PetscViewerASCIIMonitorPrintf(ctx->viewer," %g",er);
288: if (ei!=0.0) { PetscViewerASCIIMonitorPrintf(ctx->viewer,"%+gi",ei); }
289: #endif
290: PetscViewerASCIIMonitorPrintf(ctx->viewer," (%10.8e)\n",errest[i]);
291: }
292: ctx->oldnconv = nconv;
293: }
294: return(0);
295: }
297: PetscErrorCode EPSMonitorDestroy_Converged(EPSMONITOR_CONV *ctx)
298: {
299: PetscErrorCode ierr;
301: PetscViewerASCIIMonitorDestroy(ctx->viewer);
302: PetscFree(ctx);
303: return(0);
304: }
308: PetscErrorCode EPSMonitorLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
309: {
310: PetscViewer viewer = (PetscViewer) monctx;
311: PetscDraw draw,draw1;
312: PetscDrawLG lg,lg1;
314: PetscReal x,y,myeigr,p;
315: PetscScalar er,ei;
319: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); }
321: PetscViewerDrawGetDraw(viewer,0,&draw);
322: PetscViewerDrawGetDrawLG(viewer,0,&lg);
323: if (!its) {
324: PetscDrawSetTitle(draw,"Error estimates");
325: PetscDrawSetDoubleBuffer(draw);
326: PetscDrawLGSetDimension(lg,1);
327: PetscDrawLGReset(lg);
328: PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
329: }
331: /* In the hermitian case, the eigenvalues are real and can be plotted */
332: if (eps->ishermitian) {
333: PetscViewerDrawGetDraw(viewer,1,&draw1);
334: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
335: if (!its) {
336: PetscDrawSetTitle(draw1,"Approximate eigenvalues");
337: PetscDrawSetDoubleBuffer(draw1);
338: PetscDrawLGSetDimension(lg1,1);
339: PetscDrawLGReset(lg1);
340: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
341: }
342: }
344: x = (PetscReal) its;
345: if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
346: PetscDrawLGAddPoint(lg,&x,&y);
347: if (eps->ishermitian) {
348: er = eigr[nconv]; ei = eigi[nconv];
349: STBackTransform(eps->OP, 1, &er, &ei);
350: myeigr = PetscRealPart(er);
351: PetscDrawLGAddPoint(lg1,&x,&myeigr);
352: PetscDrawGetPause(draw1,&p);
353: PetscDrawSetPause(draw1,0);
354: PetscDrawLGDraw(lg1);
355: PetscDrawSetPause(draw1,p);
356: }
357: PetscDrawLGDraw(lg);
358: return(0);
359: }
364: PetscErrorCode EPSMonitorLGAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
365: {
366: PetscViewer viewer = (PetscViewer) monctx;
367: PetscDraw draw,draw1;
368: PetscDrawLG lg,lg1;
370: PetscReal *x,*y,*myeigr,p;
371: PetscInt i,n = PetscMin(eps->nev,255);
372: PetscScalar er,ei;
376: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); }
378: PetscViewerDrawGetDraw(viewer,0,&draw);
379: PetscViewerDrawGetDrawLG(viewer,0,&lg);
380: if (!its) {
381: PetscDrawSetTitle(draw,"Error estimates");
382: PetscDrawSetDoubleBuffer(draw);
383: PetscDrawLGSetDimension(lg,n);
384: PetscDrawLGReset(lg);
385: PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
386: }
388: /* In the hermitian case, the eigenvalues are real and can be plotted */
389: if (eps->ishermitian) {
390: PetscViewerDrawGetDraw(viewer,1,&draw1);
391: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
392: if (!its) {
393: PetscDrawSetTitle(draw1,"Approximate eigenvalues");
394: PetscDrawSetDoubleBuffer(draw1);
395: PetscDrawLGSetDimension(lg1,n);
396: PetscDrawLGReset(lg1);
397: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
398: }
399: }
401: PetscMalloc(sizeof(PetscReal)*n,&x);
402: PetscMalloc(sizeof(PetscReal)*n,&y);
403: for (i=0;i<n;i++) {
404: x[i] = (PetscReal) its;
405: if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
406: else y[i] = 0.0;
407: }
408: PetscDrawLGAddPoint(lg,x,y);
409: if (eps->ishermitian) {
410: PetscMalloc(sizeof(PetscReal)*n,&myeigr);
411: for(i=0;i<n;i++) {
412: er = eigr[i]; ei = eigi[i];
413: STBackTransform(eps->OP, 1, &er, &ei);
414: if (i < nest) myeigr[i] = PetscRealPart(er);
415: else myeigr[i] = 0.0;
416: }
417: PetscDrawLGAddPoint(lg1,x,myeigr);
418: PetscDrawGetPause(draw1,&p);
419: PetscDrawSetPause(draw1,0);
420: PetscDrawLGDraw(lg1);
421: PetscDrawSetPause(draw1,p);
422: PetscFree(myeigr);
423: }
424: PetscDrawLGDraw(lg);
425: PetscFree(x);
426: PetscFree(y);
427: return(0);
428: }