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-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/qepimpl.h
28: /*@C
29: QEPMonitorSet - Sets an ADDITIONAL function to be called at every
30: iteration to monitor the error estimates for each requested eigenpair.
31:
32: Collective on QEP
34: Input Parameters:
35: + qep - eigensolver context obtained from QEPCreate()
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 (QEP qep, int its, int nconv, PetscScalar *eigr, PetscScalar *eigi, PetscReal* errest, int nest, void *mctx)
45: + qep - quadratic eigensolver context obtained from QEPCreate()
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 QEPMonitorSet()
54: Options Database Keys:
55: + -qep_monitor - print only the first error estimate
56: . -qep_monitor_all - print error estimates at each iteration
57: . -qep_monitor_conv - print the eigenvalue approximations only when
58: convergence has been reached
59: . -qep_monitor_draw - sets line graph monitor for the first unconverged
60: approximate eigenvalue
61: . -qep_monitor_draw_all - sets line graph monitor for all unconverged
62: approximate eigenvalue
63: - -qep_monitor_cancel - cancels all monitors that have been hardwired into
64: a code by calls to QEPMonitorSet(), but does not cancel those set via
65: the options database.
67: Notes:
68: Several different monitoring routines may be set by calling
69: QEPMonitorSet() multiple times; all will be called in the
70: order in which they were set.
72: Level: intermediate
74: .seealso: QEPMonitorFirst(), QEPMonitorAll(), QEPMonitorLG(), QEPMonitorLGAll(), QEPMonitorCancel()
75: @*/
76: PetscErrorCode QEPMonitorSet(QEP qep,PetscErrorCode (*monitor)(QEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),
77: void *mctx,PetscErrorCode (*monitordestroy)(void*))
78: {
81: if (qep->numbermonitors >= MAXQEPMONITORS) {
82: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many QEP monitors set");
83: }
84: qep->monitor[qep->numbermonitors] = monitor;
85: qep->monitorcontext[qep->numbermonitors] = (void*)mctx;
86: qep->monitordestroy[qep->numbermonitors++] = monitordestroy;
87: return(0);
88: }
92: /*@
93: QEPMonitorCancel - Clears all monitors for a QEP object.
95: Collective on QEP
97: Input Parameters:
98: . qep - eigensolver context obtained from QEPCreate()
100: Options Database Key:
101: . -qep_monitor_cancel - Cancels all monitors that have been hardwired
102: into a code by calls to QEPMonitorSet(),
103: but does not cancel those set via the options database.
105: Level: intermediate
107: .seealso: QEPMonitorSet()
108: @*/
109: PetscErrorCode QEPMonitorCancel(QEP qep)
110: {
112: PetscInt i;
116: for (i=0; i<qep->numbermonitors; i++) {
117: if (qep->monitordestroy[i]) {
118: (*qep->monitordestroy[i])(qep->monitorcontext[i]);
119: }
120: }
121: qep->numbermonitors = 0;
122: return(0);
123: }
127: /*@C
128: QEPGetMonitorContext - Gets the monitor context, as set by
129: QEPSetMonitor() for the FIRST monitor only.
131: Not Collective
133: Input Parameter:
134: . qep - eigensolver context obtained from QEPCreate()
136: Output Parameter:
137: . ctx - monitor context
139: Level: intermediate
141: .seealso: QEPSetMonitor(), QEPDefaultMonitor()
142: @*/
143: PetscErrorCode QEPGetMonitorContext(QEP qep, void **ctx)
144: {
147: *ctx = (qep->monitorcontext[0]);
148: return(0);
149: }
153: /*@C
154: QEPMonitorAll - Print the current approximate values and
155: error estimates at each iteration of the quadratic eigensolver.
157: Collective on QEP
159: Input Parameters:
160: + qep - quadratic 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: QEPMonitorSet()
172: @*/
173: PetscErrorCode QEPMonitorAll(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
174: {
175: PetscErrorCode ierr;
176: PetscInt i;
177: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
180: if (its) {
181: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)qep)->comm,"stdout",0,&viewer);}
182: PetscViewerASCIIMonitorPrintf(viewer,"%3d QEP nconv=%d Values (Errors)",its,nconv);
183: for (i=0;i<nest;i++) {
184: #if defined(PETSC_USE_COMPLEX)
185: PetscViewerASCIIMonitorPrintf(viewer," %g%+gi",PetscRealPart(eigr[i]),PetscImaginaryPart(eigr[i]));
186: #else
187: PetscViewerASCIIMonitorPrintf(viewer," %g",eigr[i]);
188: if (eigi[i]!=0.0) { PetscViewerASCIIMonitorPrintf(viewer,"%+gi",eigi[i]); }
189: #endif
190: PetscViewerASCIIMonitorPrintf(viewer," (%10.8e)",errest[i]);
191: }
192: PetscViewerASCIIMonitorPrintf(viewer,"\n");
193: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
194: }
195: return(0);
196: }
200: /*@C
201: QEPMonitorFirst - Print the first unconverged approximate value and
202: error estimate at each iteration of the quadratic eigensolver.
204: Collective on QEP
206: Input Parameters:
207: + qep - quadratic eigensolver context
208: . its - iteration number
209: . nconv - number of converged eigenpairs so far
210: . eigr - real part of the eigenvalues
211: . eigi - imaginary part of the eigenvalues
212: . errest - error estimates
213: . nest - number of error estimates to display
214: - dummy - unused monitor context
216: Level: intermediate
218: .seealso: QEPMonitorSet()
219: @*/
220: PetscErrorCode QEPMonitorFirst(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
221: {
222: PetscErrorCode ierr;
223: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
226: if (its && nconv<nest) {
227: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)qep)->comm,"stdout",0,&viewer);}
228: PetscViewerASCIIMonitorPrintf(viewer,"%3d QEP nconv=%d first unconverged value (error)",its,nconv);
229: #if defined(PETSC_USE_COMPLEX)
230: PetscViewerASCIIMonitorPrintf(viewer," %g%+gi",PetscRealPart(eigr[nconv]),PetscImaginaryPart(eigr[nconv]));
231: #else
232: PetscViewerASCIIMonitorPrintf(viewer," %g",eigr[nconv]);
233: if (eigi[nconv]!=0.0) { PetscViewerASCIIMonitorPrintf(viewer,"%+gi",eigi[nconv]); }
234: #endif
235: PetscViewerASCIIMonitorPrintf(viewer," (%10.8e)\n",errest[nconv]);
236: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
237: }
238: return(0);
239: }
243: /*@C
244: QEPMonitorConverged - Print the approximate values and
245: error estimates as they converge.
247: Collective on QEP
249: Input Parameters:
250: + qep - quadratic eigensolver context
251: . its - iteration number
252: . nconv - number of converged eigenpairs so far
253: . eigr - real part of the eigenvalues
254: . eigi - imaginary part of the eigenvalues
255: . errest - error estimates
256: . nest - number of error estimates to display
257: - dummy - unused monitor context
259: Level: intermediate
261: .seealso: QEPMonitorSet()
262: @*/
263: PetscErrorCode QEPMonitorConverged(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *dummy)
264: {
265: PetscErrorCode ierr;
266: PetscInt i;
267: QEPMONITOR_CONV *ctx = (QEPMONITOR_CONV*) dummy;
270: if (!its) {
271: ctx->oldnconv = 0;
272: } else {
273: for (i=ctx->oldnconv;i<nconv;i++) {
274: PetscViewerASCIIMonitorPrintf(ctx->viewer,"%3d QEP converged value (error) #%d",its,i);
275: #if defined(PETSC_USE_COMPLEX)
276: PetscViewerASCIIMonitorPrintf(ctx->viewer," %g%+gi",PetscRealPart(eigr[i]),PetscImaginaryPart(eigr[i]));
277: #else
278: PetscViewerASCIIMonitorPrintf(ctx->viewer," %g",eigr[i]);
279: if (eigi[i]!=0.0) { PetscViewerASCIIMonitorPrintf(ctx->viewer,"%+gi",eigi[i]); }
280: #endif
281: PetscViewerASCIIMonitorPrintf(ctx->viewer," (%10.8e)\n",errest[i]);
282: }
283: ctx->oldnconv = nconv;
284: }
285: return(0);
286: }
288: PetscErrorCode QEPMonitorDestroy_Converged(QEPMONITOR_CONV *ctx)
289: {
290: PetscErrorCode ierr;
292: PetscViewerASCIIMonitorDestroy(ctx->viewer);
293: PetscFree(ctx);
294: return(0);
295: }
299: PetscErrorCode QEPMonitorLG(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
300: {
301: PetscViewer viewer = (PetscViewer) monctx;
302: PetscDraw draw;
303: PetscDrawLG lg;
305: PetscReal x,y;
309: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)qep)->comm); }
311: PetscViewerDrawGetDraw(viewer,0,&draw);
312: PetscViewerDrawGetDrawLG(viewer,0,&lg);
313: if (!its) {
314: PetscDrawSetTitle(draw,"Error estimates");
315: PetscDrawSetDoubleBuffer(draw);
316: PetscDrawLGSetDimension(lg,1);
317: PetscDrawLGReset(lg);
318: PetscDrawLGSetLimits(lg,0,1.0,log10(qep->tol)-2,0.0);
319: }
321: x = (PetscReal) its;
322: if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
323: PetscDrawLGAddPoint(lg,&x,&y);
325: PetscDrawLGDraw(lg);
326: return(0);
327: }
331: PetscErrorCode QEPMonitorLGAll(QEP qep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
332: {
333: PetscViewer viewer = (PetscViewer) monctx;
334: PetscDraw draw;
335: PetscDrawLG lg;
337: PetscReal *x,*y;
338: PetscInt i;
339: int n = PetscMin(qep->nev,255);
343: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)qep)->comm); }
345: PetscViewerDrawGetDraw(viewer,0,&draw);
346: PetscViewerDrawGetDrawLG(viewer,0,&lg);
347: if (!its) {
348: PetscDrawSetTitle(draw,"Error estimates");
349: PetscDrawSetDoubleBuffer(draw);
350: PetscDrawLGSetDimension(lg,n);
351: PetscDrawLGReset(lg);
352: PetscDrawLGSetLimits(lg,0,1.0,log10(qep->tol)-2,0.0);
353: }
355: PetscMalloc(sizeof(PetscReal)*n,&x);
356: PetscMalloc(sizeof(PetscReal)*n,&y);
357: for (i=0;i<n;i++) {
358: x[i] = (PetscReal) its;
359: if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
360: else y[i] = 0.0;
361: }
362: PetscDrawLGAddPoint(lg,x,y);
364: PetscDrawLGDraw(lg);
365: PetscFree(x);
366: PetscFree(y);
367: return(0);
368: }