Actual source code: nepmon.c
slepc-3.6.1 2015-09-03
1: /*
2: NEP routines related to monitors.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
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 <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/
25: #include <petscdraw.h>
29: /*
30: Runs the user provided monitor routines, if any.
31: */
32: PetscErrorCode NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eig,PetscReal *errest,PetscInt nest)
33: {
35: PetscInt i,n = nep->numbermonitors;
38: for (i=0;i<n;i++) {
39: (*nep->monitor[i])(nep,it,nconv,eig,errest,nest,nep->monitorcontext[i]);
40: }
41: return(0);
42: }
46: /*@C
47: NEPMonitorSet - Sets an ADDITIONAL function to be called at every
48: iteration to monitor the error estimates for each requested eigenpair.
50: Logically Collective on NEP
52: Input Parameters:
53: + nep - eigensolver context obtained from NEPCreate()
54: . monitor - pointer to function (if this is NULL, it turns off monitoring)
55: . mctx - [optional] context for private data for the
56: monitor routine (use NULL if no context is desired)
57: - monitordestroy - [optional] routine that frees monitor context
58: (may be NULL)
60: Calling Sequence of monitor:
61: $ monitor (NEP nep, int its, int nconv, PetscScalar *eig, PetscReal* errest, int nest, void *mctx)
63: + nep - nonlinear eigensolver context obtained from NEPCreate()
64: . its - iteration number
65: . nconv - number of converged eigenpairs
66: . eig - eigenvalues
67: . errest - error estimates for each eigenpair
68: . nest - number of error estimates
69: - mctx - optional monitoring context, as set by NEPMonitorSet()
71: Options Database Keys:
72: + -nep_monitor - print only the first error estimate
73: . -nep_monitor_all - print error estimates at each iteration
74: . -nep_monitor_conv - print the eigenvalue approximations only when
75: convergence has been reached
76: . -nep_monitor_lg - sets line graph monitor for the first unconverged
77: approximate eigenvalue
78: . -nep_monitor_lg_all - sets line graph monitor for all unconverged
79: approximate eigenvalues
80: - -nep_monitor_cancel - cancels all monitors that have been hardwired into
81: a code by calls to NEPMonitorSet(), but does not cancel those set via
82: the options database.
84: Notes:
85: Several different monitoring routines may be set by calling
86: NEPMonitorSet() multiple times; all will be called in the
87: order in which they were set.
89: Level: intermediate
91: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
92: @*/
93: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
94: {
97: if (nep->numbermonitors >= MAXNEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
98: nep->monitor[nep->numbermonitors] = monitor;
99: nep->monitorcontext[nep->numbermonitors] = (void*)mctx;
100: nep->monitordestroy[nep->numbermonitors++] = monitordestroy;
101: return(0);
102: }
106: /*@
107: NEPMonitorCancel - Clears all monitors for a NEP object.
109: Logically Collective on NEP
111: Input Parameters:
112: . nep - eigensolver context obtained from NEPCreate()
114: Options Database Key:
115: . -nep_monitor_cancel - Cancels all monitors that have been hardwired
116: into a code by calls to NEPMonitorSet(),
117: but does not cancel those set via the options database.
119: Level: intermediate
121: .seealso: NEPMonitorSet()
122: @*/
123: PetscErrorCode NEPMonitorCancel(NEP nep)
124: {
126: PetscInt i;
130: for (i=0; i<nep->numbermonitors; i++) {
131: if (nep->monitordestroy[i]) {
132: (*nep->monitordestroy[i])(&nep->monitorcontext[i]);
133: }
134: }
135: nep->numbermonitors = 0;
136: return(0);
137: }
141: /*@C
142: NEPGetMonitorContext - Gets the monitor context, as set by
143: NEPMonitorSet() for the FIRST monitor only.
145: Not Collective
147: Input Parameter:
148: . nep - eigensolver context obtained from NEPCreate()
150: Output Parameter:
151: . ctx - monitor context
153: Level: intermediate
155: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
156: @*/
157: PetscErrorCode NEPGetMonitorContext(NEP nep,void **ctx)
158: {
161: *ctx = nep->monitorcontext[0];
162: return(0);
163: }
167: /*@C
168: NEPMonitorAll - Print the current approximate values and
169: error estimates at each iteration of the nonlinear eigensolver.
171: Collective on NEP
173: Input Parameters:
174: + nep - nonlinear eigensolver context
175: . its - iteration number
176: . nconv - number of converged eigenpairs so far
177: . eig - eigenvalues
178: . errest - error estimates
179: . nest - number of error estimates to display
180: - monctx - monitor context (contains viewer, can be NULL)
182: Level: intermediate
184: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
185: @*/
186: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eig,PetscReal *errest,PetscInt nest,void *monctx)
187: {
189: PetscInt i;
190: PetscViewer viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
193: if (its) {
194: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
195: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D Values (Errors)",its,nconv);
196: for (i=0;i<nest;i++) {
197: #if defined(PETSC_USE_COMPLEX)
198: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eig[i]),(double)PetscImaginaryPart(eig[i]));
199: #else
200: PetscViewerASCIIPrintf(viewer," %g",(double)eig[i]);
201: #endif
202: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
203: }
204: PetscViewerASCIIPrintf(viewer,"\n");
205: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
206: }
207: return(0);
208: }
212: /*@C
213: NEPMonitorFirst - Print the first unconverged approximate value and
214: error estimate at each iteration of the nonlinear eigensolver.
216: Collective on NEP
218: Input Parameters:
219: + nep - nonlinear eigensolver context
220: . its - iteration number
221: . nconv - number of converged eigenpairs so far
222: . eig - eigenvalues
223: . errest - error estimates
224: . nest - number of error estimates to display
225: - monctx - monitor context (contains viewer, can be NULL)
227: Level: intermediate
229: .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
230: @*/
231: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eig,PetscReal *errest,PetscInt nest,void *monctx)
232: {
234: PetscViewer viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
237: if (its && nconv<nest) {
238: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
239: PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D first unconverged value (error)",its,nconv);
240: #if defined(PETSC_USE_COMPLEX)
241: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eig[nconv]),(double)PetscImaginaryPart(eig[nconv]));
242: #else
243: PetscViewerASCIIPrintf(viewer," %g",(double)eig[nconv]);
244: #endif
245: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
246: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
247: }
248: return(0);
249: }
253: /*@C
254: NEPMonitorConverged - Print the approximate values and
255: error estimates as they converge.
257: Collective on NEP
259: Input Parameters:
260: + nep - nonlinear eigensolver context
261: . its - iteration number
262: . nconv - number of converged eigenpairs so far
263: . eig - eigenvalues
264: . errest - error estimates
265: . nest - number of error estimates to display
266: - monctx - monitor context
268: Level: intermediate
270: Note:
271: The monitor context must contain a struct with a PetscViewer and a
272: PetscInt. In Fortran, pass a PETSC_NULL_OBJECT.
274: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
275: @*/
276: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eig,PetscReal *errest,PetscInt nest,void *monctx)
277: {
278: PetscErrorCode ierr;
279: PetscInt i;
280: PetscViewer viewer;
281: SlepcConvMonitor ctx = (SlepcConvMonitor)monctx;
284: if (!monctx) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"Must provide a context for NEPMonitorConverged");
285: if (!its) {
286: ctx->oldnconv = 0;
287: } else {
288: viewer = ctx->viewer? ctx->viewer: PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
289: for (i=ctx->oldnconv;i<nconv;i++) {
290: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
291: PetscViewerASCIIPrintf(viewer,"%3D NEP converged value (error) #%D",its,i);
292: #if defined(PETSC_USE_COMPLEX)
293: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eig[i]),(double)PetscImaginaryPart(eig[i]));
294: #else
295: PetscViewerASCIIPrintf(viewer," %g",(double)eig[i]);
296: #endif
297: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
298: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
299: }
300: ctx->oldnconv = nconv;
301: }
302: return(0);
303: }
307: PetscErrorCode NEPMonitorLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eig,PetscReal *errest,PetscInt nest,void *monctx)
308: {
309: PetscViewer viewer = (PetscViewer)monctx;
310: PetscDraw draw;
311: PetscDrawLG lg;
313: PetscReal x,y;
316: if (!viewer) viewer = PETSC_VIEWER_DRAW_(PetscObjectComm((PetscObject)nep));
317: PetscViewerDrawGetDraw(viewer,0,&draw);
318: PetscViewerDrawGetDrawLG(viewer,0,&lg);
319: if (!its) {
320: PetscDrawSetTitle(draw,"Error estimates");
321: PetscDrawSetDoubleBuffer(draw);
322: PetscDrawLGSetDimension(lg,1);
323: PetscDrawLGReset(lg);
324: PetscDrawLGSetLimits(lg,0,1.0,PetscLog10Real(nep->rtol)-2,0.0);
325: }
327: x = (PetscReal)its;
328: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]); else y = 0.0;
329: PetscDrawLGAddPoint(lg,&x,&y);
331: PetscDrawLGDraw(lg);
332: return(0);
333: }
337: PetscErrorCode NEPMonitorLGAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eig,PetscReal *errest,PetscInt nest,void *monctx)
338: {
339: PetscViewer viewer = (PetscViewer)monctx;
340: PetscDraw draw;
341: PetscDrawLG lg;
343: PetscReal *x,*y;
344: PetscInt i,n = PetscMin(nep->nev,255);
347: if (!viewer) viewer = PETSC_VIEWER_DRAW_(PetscObjectComm((PetscObject)nep));
348: PetscViewerDrawGetDraw(viewer,0,&draw);
349: PetscViewerDrawGetDrawLG(viewer,0,&lg);
350: if (!its) {
351: PetscDrawSetTitle(draw,"Error estimates");
352: PetscDrawSetDoubleBuffer(draw);
353: PetscDrawLGSetDimension(lg,n);
354: PetscDrawLGReset(lg);
355: PetscDrawLGSetLimits(lg,0,1.0,PetscLog10Real(nep->rtol)-2,0.0);
356: }
358: PetscMalloc2(n,&x,n,&y);
359: for (i=0;i<n;i++) {
360: x[i] = (PetscReal)its;
361: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
362: else y[i] = 0.0;
363: }
364: PetscDrawLGAddPoint(lg,x,y);
366: PetscDrawLGDraw(lg);
367: PetscFree2(x,y);
368: return(0);
369: }