Actual source code: epsmon.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: EPS routines related to monitors
12: */
14: #include <slepc/private/epsimpl.h> /*I "slepceps.h" I*/
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode EPSMonitor(EPS eps,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
23: PetscInt i,n = eps->numbermonitors;
26: for (i=0;i<n;i++) {
27: (*eps->monitor[i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[i]);
28: }
29: return(0);
30: }
32: /*@C
33: EPSMonitorSet - Sets an ADDITIONAL function to be called at every
34: iteration to monitor the error estimates for each requested eigenpair.
36: Logically Collective on EPS
38: Input Parameters:
39: + eps - eigensolver context obtained from EPSCreate()
40: . monitor - pointer to function (if this is NULL, it turns off monitoring)
41: . mctx - [optional] context for private data for the
42: monitor routine (use NULL if no context is desired)
43: - monitordestroy - [optional] routine that frees monitor context (may be NULL)
45: Calling Sequence of monitor:
46: $ monitor(EPS eps,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)
48: + eps - eigensolver context obtained from EPSCreate()
49: . its - iteration number
50: . nconv - number of converged eigenpairs
51: . eigr - real part of the eigenvalues
52: . eigi - imaginary part of the eigenvalues
53: . errest - relative error estimates for each eigenpair
54: . nest - number of error estimates
55: - mctx - optional monitoring context, as set by EPSMonitorSet()
57: Options Database Keys:
58: + -eps_monitor - print only the first error estimate
59: . -eps_monitor_all - print error estimates at each iteration
60: . -eps_monitor_conv - print the eigenvalue approximations only when
61: convergence has been reached
62: . -eps_monitor_lg - sets line graph monitor for the first unconverged
63: approximate eigenvalue
64: . -eps_monitor_lg_all - sets line graph monitor for all unconverged
65: approximate eigenvalues
66: - -eps_monitor_cancel - cancels all monitors that have been hardwired into
67: a code by calls to EPSMonitorSet(), but does not cancel those set via
68: the options database.
70: Notes:
71: Several different monitoring routines may be set by calling
72: EPSMonitorSet() multiple times; all will be called in the
73: order in which they were set.
75: Level: intermediate
77: .seealso: EPSMonitorFirst(), EPSMonitorAll(), EPSMonitorCancel()
78: @*/
79: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
80: {
83: if (eps->numbermonitors >= MAXEPSMONITORS) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
84: eps->monitor[eps->numbermonitors] = monitor;
85: eps->monitorcontext[eps->numbermonitors] = (void*)mctx;
86: eps->monitordestroy[eps->numbermonitors++] = monitordestroy;
87: return(0);
88: }
90: /*@
91: EPSMonitorCancel - Clears all monitors for an EPS object.
93: Logically Collective on EPS
95: Input Parameters:
96: . eps - eigensolver context obtained from EPSCreate()
98: Options Database Key:
99: . -eps_monitor_cancel - Cancels all monitors that have been hardwired
100: into a code by calls to EPSMonitorSet(),
101: but does not cancel those set via the options database.
103: Level: intermediate
105: .seealso: EPSMonitorSet()
106: @*/
107: PetscErrorCode EPSMonitorCancel(EPS eps)
108: {
110: PetscInt i;
114: for (i=0; i<eps->numbermonitors; i++) {
115: if (eps->monitordestroy[i]) {
116: (*eps->monitordestroy[i])(&eps->monitorcontext[i]);
117: }
118: }
119: eps->numbermonitors = 0;
120: return(0);
121: }
123: /*@C
124: EPSGetMonitorContext - Gets the monitor context, as set by
125: EPSMonitorSet() for the FIRST monitor only.
127: Not Collective
129: Input Parameter:
130: . eps - eigensolver context obtained from EPSCreate()
132: Output Parameter:
133: . ctx - monitor context
135: Level: intermediate
137: .seealso: EPSMonitorSet()
138: @*/
139: PetscErrorCode EPSGetMonitorContext(EPS eps,void **ctx)
140: {
143: *ctx = eps->monitorcontext[0];
144: return(0);
145: }
147: /*@C
148: EPSMonitorAll - 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: - vf - viewer and format for monitoring
163: Level: intermediate
165: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorConverged()
166: @*/
167: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
168: {
170: PetscInt i;
171: PetscScalar er,ei;
172: PetscViewer viewer;
177: viewer = vf->viewer;
179: PetscViewerPushFormat(viewer,vf->format);
180: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
181: if (its==1 && ((PetscObject)eps)->prefix) {
182: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix);
183: }
184: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D Values (Errors)",its,nconv);
185: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
186: for (i=0;i<nest;i++) {
187: er = eigr[i]; ei = eigi[i];
188: STBackTransform(eps->st,1,&er,&ei);
189: #if defined(PETSC_USE_COMPLEX)
190: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
191: #else
192: PetscViewerASCIIPrintf(viewer," %g",(double)er);
193: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
194: #endif
195: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
196: }
197: PetscViewerASCIIPrintf(viewer,"\n");
198: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
199: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
200: PetscViewerPopFormat(viewer);
201: return(0);
202: }
204: /*@C
205: EPSMonitorFirst - Print the first approximate value and
206: error estimate at each iteration of the eigensolver.
208: Collective on EPS
210: Input Parameters:
211: + eps - eigensolver context
212: . its - iteration number
213: . nconv - number of converged eigenpairs so far
214: . eigr - real part of the eigenvalues
215: . eigi - imaginary part of the eigenvalues
216: . errest - error estimates
217: . nest - number of error estimates to display
218: - vf - viewer and format for monitoring
220: Level: intermediate
222: .seealso: EPSMonitorSet(), EPSMonitorAll(), EPSMonitorConverged()
223: @*/
224: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
225: {
227: PetscScalar er,ei;
228: PetscViewer viewer;
233: viewer = vf->viewer;
235: if (its==1 && ((PetscObject)eps)->prefix) {
236: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix);
237: }
238: if (nconv<nest) {
239: PetscViewerPushFormat(viewer,vf->format);
240: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
241: PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D first unconverged value (error)",its,nconv);
242: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
243: er = eigr[nconv]; ei = eigi[nconv];
244: STBackTransform(eps->st,1,&er,&ei);
245: #if defined(PETSC_USE_COMPLEX)
246: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
247: #else
248: PetscViewerASCIIPrintf(viewer," %g",(double)er);
249: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
250: #endif
251: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
252: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
253: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
254: PetscViewerPopFormat(viewer);
255: }
256: return(0);
257: }
259: /*@C
260: EPSMonitorConverged - Print the approximate values and
261: error estimates as they converge.
263: Collective on EPS
265: Input Parameters:
266: + eps - 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: - ctx - monitor context
275: Level: intermediate
277: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorAll()
278: @*/
279: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,SlepcConvMonitor ctx)
280: {
282: PetscInt i;
283: PetscScalar er,ei;
284: PetscViewer viewer;
289: viewer = ctx->viewer;
291: if (its==1 && ((PetscObject)eps)->prefix) {
292: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)eps)->prefix);
293: }
294: if (its==1) ctx->oldnconv = 0;
295: if (ctx->oldnconv!=nconv) {
296: PetscViewerPushFormat(viewer,ctx->format);
297: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
298: for (i=ctx->oldnconv;i<nconv;i++) {
299: PetscViewerASCIIPrintf(viewer,"%3D EPS converged value (error) #%D",its,i);
300: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
301: er = eigr[i]; ei = eigi[i];
302: STBackTransform(eps->st,1,&er,&ei);
303: #if defined(PETSC_USE_COMPLEX)
304: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
305: #else
306: PetscViewerASCIIPrintf(viewer," %g",(double)er);
307: if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
308: #endif
309: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
310: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
311: }
312: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
313: PetscViewerPopFormat(viewer);
314: ctx->oldnconv = nconv;
315: }
316: return(0);
317: }
319: /*@C
320: EPSMonitorLGCreate - Creates a line graph context for use with
321: EPS to monitor convergence.
323: Collective on MPI_Comm
325: Input Parameters:
326: + comm - communicator context
327: . host - the X display to open, or null for the local machine
328: . label - the title to put in the title bar
329: . x, y - the screen coordinates of the upper left coordinate of
330: the window
331: - m, n - the screen width and height in pixels
333: Output Parameter:
334: . lgctx - the drawing context
336: Options Database Keys:
337: + -eps_monitor_lg - Sets line graph monitor for the first residual
338: - -eps_monitor_lg_all - Sets line graph monitor for all residuals
340: Notes:
341: Use PetscDrawLGDestroy() to destroy this line graph.
343: Level: intermediate
345: .seealso: EPSMonitorSet()
346: @*/
347: PetscErrorCode EPSMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
348: {
349: PetscDraw draw;
350: PetscDrawLG lg;
354: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
355: PetscDrawSetFromOptions(draw);
356: PetscDrawLGCreate(draw,1,&lg);
357: PetscDrawLGSetFromOptions(lg);
358: PetscDrawDestroy(&draw);
359: *lgctx = lg;
360: return(0);
361: }
363: PetscErrorCode EPSMonitorLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
364: {
365: PetscDrawLG lg = (PetscDrawLG)ctx;
366: PetscReal x,y;
371: if (its==1) {
372: PetscDrawLGReset(lg);
373: PetscDrawLGSetDimension(lg,1);
374: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(eps->tol)-2,0.0);
375: }
376: x = (PetscReal)its;
377: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
378: else y = 0.0;
379: PetscDrawLGAddPoint(lg,&x,&y);
380: if (its <= 20 || !(its % 5) || eps->reason) {
381: PetscDrawLGDraw(lg);
382: PetscDrawLGSave(lg);
383: }
384: return(0);
385: }
387: PetscErrorCode EPSMonitorLGAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
388: {
389: PetscDrawLG lg = (PetscDrawLG)ctx;
390: PetscInt i,n = PetscMin(eps->nev,255);
391: PetscReal *x,*y;
396: if (its==1) {
397: PetscDrawLGReset(lg);
398: PetscDrawLGSetDimension(lg,n);
399: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(eps->tol)-2,0.0);
400: }
401: PetscMalloc2(n,&x,n,&y);
402: for (i=0;i<n;i++) {
403: x[i] = (PetscReal)its;
404: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
405: else y[i] = 0.0;
406: }
407: PetscDrawLGAddPoint(lg,x,y);
408: if (its <= 20 || !(its % 5) || eps->reason) {
409: PetscDrawLGDraw(lg);
410: PetscDrawLGSave(lg);
411: }
412: PetscFree2(x,y);
413: return(0);
414: }