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-2009, 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 error estimates at each iteration
56: . -eps_monitor_draw - sets line graph monitor
57: - -eps_monitor_cancel - cancels all monitors that have been hardwired into
58: a code by calls to EPSMonitorSet(), but does not cancel those set via
59: the options database.
61: Notes:
62: Several different monitoring routines may be set by calling
63: EPSMonitorSet() multiple times; all will be called in the
64: order in which they were set.
66: Level: intermediate
68: .seealso: EPSMonitorDefault(), EPSMonitorCancel()
69: @*/
70: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),
71: void *mctx,PetscErrorCode (*monitordestroy)(void*))
72: {
75: if (eps->numbermonitors >= MAXEPSMONITORS) {
76: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
77: }
78: eps->monitor[eps->numbermonitors] = monitor;
79: eps->monitorcontext[eps->numbermonitors] = (void*)mctx;
80: eps->monitordestroy[eps->numbermonitors++] = monitordestroy;
81: return(0);
82: }
86: /*@
87: EPSMonitorCancel - Clears all monitors for an EPS object.
89: Collective on EPS
91: Input Parameters:
92: . eps - eigensolver context obtained from EPSCreate()
94: Options Database Key:
95: . -eps_monitor_cancel - Cancels all monitors that have been hardwired
96: into a code by calls to EPSMonitorSet(),
97: but does not cancel those set via the options database.
99: Level: intermediate
101: .seealso: EPSMonitorSet()
102: @*/
103: PetscErrorCode EPSMonitorCancel(EPS eps)
104: {
106: PetscInt i;
110: for (i=0; i<eps->numbermonitors; i++) {
111: if (eps->monitordestroy[i]) {
112: (*eps->monitordestroy[i])(eps->monitorcontext[i]);
113: }
114: }
115: eps->numbermonitors = 0;
116: return(0);
117: }
121: /*@C
122: EPSGetMonitorContext - Gets the monitor context, as set by
123: EPSSetMonitor() for the FIRST monitor only.
125: Not Collective
127: Input Parameter:
128: . eps - eigensolver context obtained from EPSCreate()
130: Output Parameter:
131: . ctx - monitor context
133: Level: intermediate
135: .seealso: EPSSetMonitor(), EPSDefaultMonitor()
136: @*/
137: PetscErrorCode EPSGetMonitorContext(EPS eps, void **ctx)
138: {
141: *ctx = (eps->monitorcontext[0]);
142: return(0);
143: }
147: /*@C
148: EPSMonitorDefault - Print the current approximate values and
149: error estimates at each iteration of the eigensolver.
151: Collective on EPS
153: Input Parameters:
154: + eps - eigensolver context
155: . its - iteration number
156: . nconv - number of converged eigenpairs so far
157: . eigr - real part of the eigenvalues
158: . eigi - imaginary part of the eigenvalues
159: . errest - error estimates
160: . nest - number of error estimates to display
161: - dummy - unused monitor context
163: Level: intermediate
165: .seealso: EPSMonitorSet()
166: @*/
167: PetscErrorCode EPSMonitorDefault(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
168: {
169: PetscErrorCode ierr;
170: PetscInt i;
171: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
174: if (its) {
175: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)eps)->comm,"stdout",0,&viewer);}
176: PetscViewerASCIIMonitorPrintf(viewer,"%3d EPS nconv=%d Values (Errors)",its,nconv);
177: for (i=0;i<nest;i++) {
178: #if defined(PETSC_USE_COMPLEX)
179: PetscViewerASCIIMonitorPrintf(viewer," %g%+gi",PetscRealPart(eigr[i]),PetscImaginaryPart(eigr[i]));
180: #else
181: PetscViewerASCIIMonitorPrintf(viewer," %g",eigr[i]);
182: if (eigi[i]!=0.0) { PetscViewerASCIIMonitorPrintf(viewer,"%+gi",eigi[i]); }
183: #endif
184: PetscViewerASCIIMonitorPrintf(viewer," (%10.8e)",errest[i]);
185: }
186: PetscViewerASCIIMonitorPrintf(viewer,"\n");
187: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
188: }
189: return(0);
190: }
194: PetscErrorCode EPSMonitorLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
195: {
196: PetscViewer viewer = (PetscViewer) monctx;
197: PetscDraw draw;
198: PetscDrawLG lg;
200: PetscReal *x,*y;
201: PetscInt i;
202: int n = eps->nev;
203: #if !defined(PETSC_USE_COMPLEX)
204: int p;
205: PetscDraw draw1;
206: PetscDrawLG lg1;
207: #endif
211: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); }
213: PetscViewerDrawGetDraw(viewer,0,&draw);
214: PetscViewerDrawGetDrawLG(viewer,0,&lg);
215: if (!its) {
216: PetscDrawSetTitle(draw,"Error estimates");
217: PetscDrawSetDoubleBuffer(draw);
218: PetscDrawLGSetDimension(lg,n);
219: PetscDrawLGReset(lg);
220: PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
221: }
223: #if !defined(PETSC_USE_COMPLEX)
224: if (eps->ishermitian) {
225: PetscViewerDrawGetDraw(viewer,1,&draw1);
226: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
227: if (!its) {
228: PetscDrawSetTitle(draw1,"Approximate eigenvalues");
229: PetscDrawSetDoubleBuffer(draw1);
230: PetscDrawLGSetDimension(lg1,n);
231: PetscDrawLGReset(lg1);
232: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
233: }
234: }
235: #endif
237: PetscMalloc(sizeof(PetscReal)*n,&x);
238: PetscMalloc(sizeof(PetscReal)*n,&y);
239: for (i=0;i<n;i++) {
240: x[i] = (PetscReal) its;
241: if (errest[i] > 0.0) y[i] = log10(errest[i]); else y[i] = 0.0;
242: }
243: PetscDrawLGAddPoint(lg,x,y);
244: #if !defined(PETSC_USE_COMPLEX)
245: if (eps->ishermitian) {
246: PetscDrawLGAddPoint(lg1,x,eps->eigr);
247: PetscDrawGetPause(draw1,&p);
248: PetscDrawSetPause(draw1,0);
249: PetscDrawLGDraw(lg1);
250: PetscDrawSetPause(draw1,p);
251: }
252: #endif
253: PetscDrawLGDraw(lg);
254: PetscFree(x);
255: PetscFree(y);
256: return(0);
257: }