Actual source code: qepmon.c
1: /*
2: QEP 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/qepimpl.h> /*I "slepcqep.h" I*/
28: /*
29: Runs the user provided monitor routines, if any.
30: */
31: PetscErrorCode QEPMonitor(QEP qep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
32: {
34: PetscInt i,n = qep->numbermonitors;
37: for (i=0;i<n;i++) {
38: (*qep->monitor[i])(qep,it,nconv,eigr,eigi,errest,nest,qep->monitorcontext[i]);
39: }
40: return(0);
41: }
45: /*@C
46: QEPMonitorSet - 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 QEP
51: Input Parameters:
52: + qep - eigensolver context obtained from QEPCreate()
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 (QEP qep, int its, int nconv, PetscScalar *eigr, PetscScalar *eigi, PetscReal* errest, int nest, void *mctx)
62: + qep - quadratic eigensolver context obtained from QEPCreate()
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 QEPMonitorSet()
71: Options Database Keys:
72: + -qep_monitor - print only the first error estimate
73: . -qep_monitor_all - print error estimates at each iteration
74: . -qep_monitor_conv - print the eigenvalue approximations only when
75: convergence has been reached
76: . -qep_monitor_draw - sets line graph monitor for the first unconverged
77: approximate eigenvalue
78: . -qep_monitor_draw_all - sets line graph monitor for all unconverged
79: approximate eigenvalue
80: - -qep_monitor_cancel - cancels all monitors that have been hardwired into
81: a code by calls to QEPMonitorSet(), but does not cancel those set via
82: the options database.
84: Notes:
85: Several different monitoring routines may be set by calling
86: QEPMonitorSet() multiple times; all will be called in the
87: order in which they were set.
89: Level: intermediate
91: .seealso: QEPMonitorFirst(), QEPMonitorAll(), QEPMonitorCancel()
92: @*/
93: PetscErrorCode QEPMonitorSet(QEP qep,PetscErrorCode (*monitor)(QEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
94: {
97: if (qep->numbermonitors >= MAXQEPMONITORS) {
98: SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Too many QEP monitors set");
99: }
100: qep->monitor[qep->numbermonitors] = monitor;
101: qep->monitorcontext[qep->numbermonitors] = (void*)mctx;
102: qep->monitordestroy[qep->numbermonitors++] = monitordestroy;
103: return(0);
104: }
108: /*@
109: QEPMonitorCancel - Clears all monitors for a QEP object.
111: Logically Collective on QEP
113: Input Parameters:
114: . qep - eigensolver context obtained from QEPCreate()
116: Options Database Key:
117: . -qep_monitor_cancel - Cancels all monitors that have been hardwired
118: into a code by calls to QEPMonitorSet(),
119: but does not cancel those set via the options database.
121: Level: intermediate
123: .seealso: QEPMonitorSet()
124: @*/
125: PetscErrorCode QEPMonitorCancel(QEP qep)
126: {
128: PetscInt i;
132: for (i=0; i<qep->numbermonitors; i++) {
133: if (qep->monitordestroy[i]) {
134: (*qep->monitordestroy[i])(&qep->monitorcontext[i]);
135: }
136: }
137: qep->numbermonitors = 0;
138: return(0);
139: }
143: /*@C
144: QEPGetMonitorContext - Gets the monitor context, as set by
145: QEPMonitorSet() for the FIRST monitor only.
147: Not Collective
149: Input Parameter:
150: . qep - eigensolver context obtained from QEPCreate()
152: Output Parameter:
153: . ctx - monitor context
155: Level: intermediate
157: .seealso: QEPMonitorSet(), QEPDefaultMonitor()
158: @*/
159: PetscErrorCode QEPGetMonitorContext(QEP qep,void **ctx)
160: {
163: *ctx = (qep->monitorcontext[0]);
164: return(0);
165: }
169: /*@C
170: QEPMonitorAll - Print the current approximate values and
171: error estimates at each iteration of the quadratic eigensolver.
173: Collective on QEP
175: Input Parameters:
176: + qep - quadratic 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: - monctx - monitor context (contains viewer, can be PETSC_NULL)
185: Level: intermediate
187: .seealso: QEPMonitorSet(), QEPMonitorFirst(), QEPMonitorConverged()
188: @*/
189: PetscErrorCode QEPMonitorAll(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
190: {
192: PetscInt i;
193: PetscViewer viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)qep)->comm);
196: if (its) {
197: PetscViewerASCIIAddTab(viewer,((PetscObject)qep)->tablevel);
198: PetscViewerASCIIPrintf(viewer,"%3D QEP nconv=%D Values (Errors)",its,nconv);
199: for (i=0;i<nest;i++) {
200: #if defined(PETSC_USE_COMPLEX)
201: PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(eigr[i]),PetscImaginaryPart(eigr[i]));
202: #else
203: PetscViewerASCIIPrintf(viewer," %G",eigr[i]);
204: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",eigi[i]); }
205: #endif
206: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
207: }
208: PetscViewerASCIIPrintf(viewer,"\n");
209: PetscViewerASCIISubtractTab(viewer,((PetscObject)qep)->tablevel);
210: }
211: return(0);
212: }
216: /*@C
217: QEPMonitorFirst - Print the first unconverged approximate value and
218: error estimate at each iteration of the quadratic eigensolver.
220: Collective on QEP
222: Input Parameters:
223: + qep - quadratic eigensolver context
224: . its - iteration number
225: . nconv - number of converged eigenpairs so far
226: . eigr - real part of the eigenvalues
227: . eigi - imaginary part of the eigenvalues
228: . errest - error estimates
229: . nest - number of error estimates to display
230: - monctx - monitor context (contains viewer, can be PETSC_NULL)
232: Level: intermediate
234: .seealso: QEPMonitorSet(), QEPMonitorAll(), QEPMonitorConverged()
235: @*/
236: PetscErrorCode QEPMonitorFirst(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
237: {
239: PetscViewer viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)qep)->comm);
242: if (its && nconv<nest) {
243: PetscViewerASCIIAddTab(viewer,((PetscObject)qep)->tablevel);
244: PetscViewerASCIIPrintf(viewer,"%3D QEP nconv=%D first unconverged value (error)",its,nconv);
245: #if defined(PETSC_USE_COMPLEX)
246: PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(eigr[nconv]),PetscImaginaryPart(eigr[nconv]));
247: #else
248: PetscViewerASCIIPrintf(viewer," %G",eigr[nconv]);
249: if (eigi[nconv]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",eigi[nconv]); }
250: #endif
251: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
252: PetscViewerASCIISubtractTab(viewer,((PetscObject)qep)->tablevel);
253: }
254: return(0);
255: }
259: /*@C
260: QEPMonitorConverged - Print the approximate values and
261: error estimates as they converge.
263: Collective on QEP
265: Input Parameters:
266: + qep - quadratic eigensolver context
267: . its - iteration number
268: . nconv - number of converged eigenpairs so far
269: . eigr - real part of the eigenvalues
270: . eigi - imaginary part of the eigenvalues
271: . errest - error estimates
272: . nest - number of error estimates to display
273: - monctx - monitor context
275: Level: intermediate
277: Note:
278: The monitor context must contain a struct with a PetscViewer and a
279: PetscInt. In Fortran, pass a PETSC_NULL_OBJECT.
281: .seealso: QEPMonitorSet(), QEPMonitorFirst(), QEPMonitorAll()
282: @*/
283: PetscErrorCode QEPMonitorConverged(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
284: {
285: PetscErrorCode ierr;
286: PetscInt i;
287: PetscViewer viewer;
288: SlepcConvMonitor ctx = (SlepcConvMonitor)monctx;
291: if (!monctx) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_WRONG,"Must provide a context for QEPMonitorConverged");
292: if (!its) {
293: ctx->oldnconv = 0;
294: } else {
295: viewer = ctx->viewer? ctx->viewer: PETSC_VIEWER_STDOUT_(((PetscObject)qep)->comm);
296: for (i=ctx->oldnconv;i<nconv;i++) {
297: PetscViewerASCIIAddTab(viewer,((PetscObject)qep)->tablevel);
298: PetscViewerASCIIPrintf(viewer,"%3D QEP converged value (error) #%D",its,i);
299: #if defined(PETSC_USE_COMPLEX)
300: PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(eigr[i]),PetscImaginaryPart(eigr[i]));
301: #else
302: PetscViewerASCIIPrintf(viewer," %G",eigr[i]);
303: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",eigi[i]); }
304: #endif
305: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
306: PetscViewerASCIISubtractTab(viewer,((PetscObject)qep)->tablevel);
307: }
308: ctx->oldnconv = nconv;
309: }
310: return(0);
311: }
315: PetscErrorCode QEPMonitorLG(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
316: {
317: PetscViewer viewer = (PetscViewer) monctx;
318: PetscDraw draw;
319: PetscDrawLG lg;
321: PetscReal x,y;
324: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)qep)->comm); }
325: PetscViewerDrawGetDraw(viewer,0,&draw);
326: PetscViewerDrawGetDrawLG(viewer,0,&lg);
327: if (!its) {
328: PetscDrawSetTitle(draw,"Error estimates");
329: PetscDrawSetDoubleBuffer(draw);
330: PetscDrawLGSetDimension(lg,1);
331: PetscDrawLGReset(lg);
332: PetscDrawLGSetLimits(lg,0,1.0,log10(qep->tol)-2,0.0);
333: }
335: x = (PetscReal) its;
336: if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
337: PetscDrawLGAddPoint(lg,&x,&y);
339: PetscDrawLGDraw(lg);
340: return(0);
341: }
345: PetscErrorCode QEPMonitorLGAll(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
346: {
347: PetscViewer viewer = (PetscViewer) monctx;
348: PetscDraw draw;
349: PetscDrawLG lg;
351: PetscReal *x,*y;
352: PetscInt i,n = PetscMin(qep->nev,255);
355: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)qep)->comm); }
356: PetscViewerDrawGetDraw(viewer,0,&draw);
357: PetscViewerDrawGetDrawLG(viewer,0,&lg);
358: if (!its) {
359: PetscDrawSetTitle(draw,"Error estimates");
360: PetscDrawSetDoubleBuffer(draw);
361: PetscDrawLGSetDimension(lg,n);
362: PetscDrawLGReset(lg);
363: PetscDrawLGSetLimits(lg,0,1.0,log10(qep->tol)-2,0.0);
364: }
366: PetscMalloc(sizeof(PetscReal)*n,&x);
367: PetscMalloc(sizeof(PetscReal)*n,&y);
368: for (i=0;i<n;i++) {
369: x[i] = (PetscReal) its;
370: if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
371: else y[i] = 0.0;
372: }
373: PetscDrawLGAddPoint(lg,x,y);
375: PetscDrawLGDraw(lg);
376: PetscFree(x);
377: PetscFree(y);
378: return(0);
379: }