Actual source code: pepmon.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: PEP routines related to monitors
12: */
14: #include <slepc/private/pepimpl.h> /*I "slepcpep.h" I*/
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode PEPMonitor(PEP pep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
23: PetscInt i,n = pep->numbermonitors;
26: for (i=0;i<n;i++) {
27: (*pep->monitor[i])(pep,it,nconv,eigr,eigi,errest,nest,pep->monitorcontext[i]);
28: }
29: return(0);
30: }
32: /*@C
33: PEPMonitorSet - Sets an ADDITIONAL function to be called at every
34: iteration to monitor the error estimates for each requested eigenpair.
36: Logically Collective on PEP
38: Input Parameters:
39: + pep - eigensolver context obtained from PEPCreate()
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(PEP pep,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)
48: + pep - polynomial eigensolver context obtained from PEPCreate()
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 PEPMonitorSet()
57: Options Database Keys:
58: + -pep_monitor - print only the first error estimate
59: . -pep_monitor_all - print error estimates at each iteration
60: . -pep_monitor_conv - print the eigenvalue approximations only when
61: convergence has been reached
62: . -pep_monitor_lg - sets line graph monitor for the first unconverged
63: approximate eigenvalue
64: . -pep_monitor_lg_all - sets line graph monitor for all unconverged
65: approximate eigenvalues
66: - -pep_monitor_cancel - cancels all monitors that have been hardwired into
67: a code by calls to PEPMonitorSet(), but does not cancel those set via
68: the options database.
70: Notes:
71: Several different monitoring routines may be set by calling
72: PEPMonitorSet() multiple times; all will be called in the
73: order in which they were set.
75: Level: intermediate
77: .seealso: PEPMonitorFirst(), PEPMonitorAll(), PEPMonitorCancel()
78: @*/
79: PetscErrorCode PEPMonitorSet(PEP pep,PetscErrorCode (*monitor)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
80: {
83: if (pep->numbermonitors >= MAXPEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Too many PEP monitors set");
84: pep->monitor[pep->numbermonitors] = monitor;
85: pep->monitorcontext[pep->numbermonitors] = (void*)mctx;
86: pep->monitordestroy[pep->numbermonitors++] = monitordestroy;
87: return(0);
88: }
90: /*@
91: PEPMonitorCancel - Clears all monitors for a PEP object.
93: Logically Collective on PEP
95: Input Parameters:
96: . pep - eigensolver context obtained from PEPCreate()
98: Options Database Key:
99: . -pep_monitor_cancel - Cancels all monitors that have been hardwired
100: into a code by calls to PEPMonitorSet(),
101: but does not cancel those set via the options database.
103: Level: intermediate
105: .seealso: PEPMonitorSet()
106: @*/
107: PetscErrorCode PEPMonitorCancel(PEP pep)
108: {
110: PetscInt i;
114: for (i=0; i<pep->numbermonitors; i++) {
115: if (pep->monitordestroy[i]) {
116: (*pep->monitordestroy[i])(&pep->monitorcontext[i]);
117: }
118: }
119: pep->numbermonitors = 0;
120: return(0);
121: }
123: /*@C
124: PEPGetMonitorContext - Gets the monitor context, as set by
125: PEPMonitorSet() for the FIRST monitor only.
127: Not Collective
129: Input Parameter:
130: . pep - eigensolver context obtained from PEPCreate()
132: Output Parameter:
133: . ctx - monitor context
135: Level: intermediate
137: .seealso: PEPMonitorSet(), PEPDefaultMonitor()
138: @*/
139: PetscErrorCode PEPGetMonitorContext(PEP pep,void **ctx)
140: {
143: *ctx = pep->monitorcontext[0];
144: return(0);
145: }
147: /*
148: Helper function to compute eigenvalue that must be viewed in monitor
149: */
150: static PetscErrorCode PEPMonitorGetTrueEig(PEP pep,PetscScalar *er,PetscScalar *ei)
151: {
153: PetscBool flg,sinv;
156: STBackTransform(pep->st,1,er,ei);
157: if (pep->sfactor!=1.0) {
158: STGetTransform(pep->st,&flg);
159: PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
160: if (flg && sinv) {
161: *er /= pep->sfactor; *ei /= pep->sfactor;
162: } else {
163: *er *= pep->sfactor; *ei *= pep->sfactor;
164: }
165: }
166: return(0);
167: }
169: /*@C
170: PEPMonitorAll - Print the current approximate values and
171: error estimates at each iteration of the polynomial eigensolver.
173: Collective on PEP
175: Input Parameters:
176: + pep - polynomial 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: - vf - viewer and format for monitoring
185: Level: intermediate
187: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorConverged()
188: @*/
189: PetscErrorCode PEPMonitorAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
190: {
192: PetscInt i;
193: PetscScalar er,ei;
194: PetscViewer viewer;
199: viewer = vf->viewer;
201: PetscViewerPushFormat(viewer,vf->format);
202: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
203: if (its==1 && ((PetscObject)pep)->prefix) {
204: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
205: }
206: PetscViewerASCIIPrintf(viewer,"%3D PEP nconv=%D Values (Errors)",its,nconv);
207: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
208: for (i=0;i<nest;i++) {
209: er = eigr[i]; ei = eigi[i];
210: PEPMonitorGetTrueEig(pep,&er,&ei);
211: #if defined(PETSC_USE_COMPLEX)
212: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
213: #else
214: PetscViewerASCIIPrintf(viewer," %g",(double)er);
215: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
216: #endif
217: PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
218: }
219: PetscViewerASCIIPrintf(viewer,"\n");
220: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
221: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
222: PetscViewerPopFormat(viewer);
223: return(0);
224: }
226: /*@C
227: PEPMonitorFirst - Print the first unconverged approximate value and
228: error estimate at each iteration of the polynomial eigensolver.
230: Collective on PEP
232: Input Parameters:
233: + pep - polynomial eigensolver context
234: . its - iteration number
235: . nconv - number of converged eigenpairs so far
236: . eigr - real part of the eigenvalues
237: . eigi - imaginary part of the eigenvalues
238: . errest - error estimates
239: . nest - number of error estimates to display
240: - vf - viewer and format for monitoring
242: Level: intermediate
244: .seealso: PEPMonitorSet(), PEPMonitorAll(), PEPMonitorConverged()
245: @*/
246: PetscErrorCode PEPMonitorFirst(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
247: {
249: PetscScalar er,ei;
250: PetscViewer viewer;
255: viewer = vf->viewer;
257: if (its==1 && ((PetscObject)pep)->prefix) {
258: PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
259: }
260: if (nconv<nest) {
261: PetscViewerPushFormat(viewer,vf->format);
262: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
263: PetscViewerASCIIPrintf(viewer,"%3D PEP nconv=%D first unconverged value (error)",its,nconv);
264: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
265: er = eigr[nconv]; ei = eigi[nconv];
266: PEPMonitorGetTrueEig(pep,&er,&ei);
267: #if defined(PETSC_USE_COMPLEX)
268: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
269: #else
270: PetscViewerASCIIPrintf(viewer," %g",(double)er);
271: if (eigi[nconv]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
272: #endif
273: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
274: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
275: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
276: PetscViewerPopFormat(viewer);
277: }
278: return(0);
279: }
281: /*@C
282: PEPMonitorConverged - Print the approximate values and
283: error estimates as they converge.
285: Collective on PEP
287: Input Parameters:
288: + pep - polynomial eigensolver context
289: . its - iteration number
290: . nconv - number of converged eigenpairs so far
291: . eigr - real part of the eigenvalues
292: . eigi - imaginary part of the eigenvalues
293: . errest - error estimates
294: . nest - number of error estimates to display
295: - ctx - monitor context
297: Level: intermediate
299: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorAll()
300: @*/
301: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,SlepcConvMonitor ctx)
302: {
304: PetscInt i;
305: PetscScalar er,ei;
306: PetscViewer viewer;
311: viewer = ctx->viewer;
313: if (its==1 && ((PetscObject)pep)->prefix) {
314: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)pep)->prefix);
315: }
316: if (its==1) ctx->oldnconv = 0;
317: if (ctx->oldnconv!=nconv) {
318: PetscViewerPushFormat(viewer,ctx->format);
319: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
320: for (i=ctx->oldnconv;i<nconv;i++) {
321: PetscViewerASCIIPrintf(viewer,"%3D PEP converged value (error) #%D",its,i);
322: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
323: er = eigr[i]; ei = eigi[i];
324: PEPMonitorGetTrueEig(pep,&er,&ei);
325: #if defined(PETSC_USE_COMPLEX)
326: PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
327: #else
328: PetscViewerASCIIPrintf(viewer," %g",(double)er);
329: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
330: #endif
331: PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
332: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
333: }
334: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
335: PetscViewerPopFormat(viewer);
336: ctx->oldnconv = nconv;
337: }
338: return(0);
339: }
341: /*@C
342: PEPMonitorLGCreate - Creates a line graph context for use with
343: PEP to monitor convergence.
345: Collective on MPI_Comm
347: Input Parameters:
348: + comm - communicator context
349: . host - the X display to open, or null for the local machine
350: . label - the title to put in the title bar
351: . x, y - the screen coordinates of the upper left coordinate of
352: the window
353: - m, n - the screen width and height in pixels
355: Output Parameter:
356: . lgctx - the drawing context
358: Options Database Keys:
359: + -pep_monitor_lg - Sets line graph monitor for the first residual
360: - -pep_monitor_lg_all - Sets line graph monitor for all residuals
362: Notes:
363: Use PetscDrawLGDestroy() to destroy this line graph.
365: Level: intermediate
367: .seealso: PEPMonitorSet()
368: @*/
369: PetscErrorCode PEPMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
370: {
371: PetscDraw draw;
372: PetscDrawLG lg;
376: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
377: PetscDrawSetFromOptions(draw);
378: PetscDrawLGCreate(draw,1,&lg);
379: PetscDrawLGSetFromOptions(lg);
380: PetscDrawDestroy(&draw);
381: *lgctx = lg;
382: return(0);
383: }
385: PetscErrorCode PEPMonitorLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
386: {
387: PetscDrawLG lg = (PetscDrawLG)ctx;
388: PetscReal x,y;
393: if (its==1) {
394: PetscDrawLGReset(lg);
395: PetscDrawLGSetDimension(lg,1);
396: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(pep->tol)-2,0.0);
397: }
398: x = (PetscReal)its;
399: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
400: else y = 0.0;
401: PetscDrawLGAddPoint(lg,&x,&y);
402: if (its <= 20 || !(its % 5) || pep->reason) {
403: PetscDrawLGDraw(lg);
404: PetscDrawLGSave(lg);
405: }
406: return(0);
407: }
409: PetscErrorCode PEPMonitorLGAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *ctx)
410: {
411: PetscDrawLG lg = (PetscDrawLG)ctx;
412: PetscInt i,n = PetscMin(pep->nev,255);
413: PetscReal *x,*y;
418: if (its==1) {
419: PetscDrawLGReset(lg);
420: PetscDrawLGSetDimension(lg,n);
421: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(pep->tol)-2,0.0);
422: }
423: PetscMalloc2(n,&x,n,&y);
424: for (i=0;i<n;i++) {
425: x[i] = (PetscReal)its;
426: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
427: else y[i] = 0.0;
428: }
429: PetscDrawLGAddPoint(lg,x,y);
430: if (its <= 20 || !(its % 5) || pep->reason) {
431: PetscDrawLGDraw(lg);
432: PetscDrawLGSave(lg);
433: }
434: PetscFree2(x,y);
435: return(0);
436: }