Actual source code: monitor.c
1: /*
2: EPS routines related to monitors.
3: */
4: #include src/eps/epsimpl.h
8: /*@C
9: EPSSetMonitor - Sets an ADDITIONAL function to be called at every
10: iteration to monitor the error estimates for each requested eigenpair.
11:
12: Collective on EPS
14: Input Parameters:
15: + eps - eigensolver context obtained from EPSCreate()
16: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
17: - mctx - [optional] context for private data for the
18: monitor routine (use PETSC_NULL if no context is desired)
20: Calling Sequence of monitor:
21: $ monitor (EPS eps, int its, int nconv, PetscScalar *eigr, PetscScalar *eigi, PetscReal* errest, int nest, void *mctx)
23: + eps - eigensolver context obtained from EPSCreate()
24: . its - iteration number
25: . nconv - number of converged eigenpairs
26: . eigr - real part of the eigenvalues
27: . eigi - imaginary part of the eigenvalues
28: . errest - relative error estimates for each eigenpair
29: . nest - number of error estimates
30: - mctx - optional monitoring context, as set by EPSSetMonitor()
32: Options Database Keys:
33: + -eps_monitor - print error estimates at each iteration
34: - -eps_cancelmonitors - cancels all monitors that have been hardwired into
35: a code by calls to EPSetMonitor(), but does not cancel those set via
36: the options database.
38: Notes:
39: Several different monitoring routines may be set by calling
40: EPSSetMonitor() multiple times; all will be called in the
41: order in which they were set.
43: Level: intermediate
45: .seealso: EPSDefaultMonitor(), EPSClearMonitor()
46: @*/
47: PetscErrorCode EPSSetMonitor(EPS eps, int (*monitor)(EPS,int,int,PetscScalar*,PetscScalar*,PetscReal*,int,void*), void *mctx)
48: {
51: if (eps->numbermonitors >= MAXEPSMONITORS) {
52: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
53: }
54: eps->monitor[eps->numbermonitors] = monitor;
55: eps->monitorcontext[eps->numbermonitors++] = (void*)mctx;
56: return(0);
57: }
61: /*@
62: EPSClearMonitor - Clears all monitors for an EPS object.
64: Collective on EPS
66: Input Parameters:
67: . eps - eigensolver context obtained from EPSCreate()
69: Options Database Key:
70: . -eps_cancelmonitors - Cancels all monitors that have been hardwired
71: into a code by calls to EPSSetMonitor(),
72: but does not cancel those set via the options database.
74: Level: intermediate
76: .seealso: EPSSetMonitor()
77: @*/
78: PetscErrorCode EPSClearMonitor(EPS eps)
79: {
82: eps->numbermonitors = 0;
83: return(0);
84: }
88: /*@C
89: EPSGetMonitorContext - Gets the monitor context, as set by
90: EPSSetMonitor() for the FIRST monitor only.
92: Not Collective
94: Input Parameter:
95: . eps - eigensolver context obtained from EPSCreate()
97: Output Parameter:
98: . ctx - monitor context
100: Level: intermediate
102: .seealso: EPSSetMonitor(), EPSDefaultMonitor()
103: @*/
104: PetscErrorCode EPSGetMonitorContext(EPS eps, void **ctx)
105: {
108: *ctx = (eps->monitorcontext[0]);
109: return(0);
110: }
114: /*@C
115: EPSDefaultMonitor - Print the current approximate values and
116: error estimates at each iteration of the eigensolver.
118: Collective on EPS
120: Input Parameters:
121: + eps - eigensolver context
122: . its - iteration number
123: . nconv - number of converged eigenpairs so far
124: . eigr - real part of the eigenvalues
125: . eigi - imaginary part of the eigenvalues
126: . errest - error estimates
127: . nest - number of error estimates to display
128: - dummy - unused monitor context
130: Level: intermediate
132: .seealso: EPSSetMonitor()
133: @*/
134: PetscErrorCode EPSDefaultMonitor(EPS eps,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,int nest,void *dummy)
135: {
137: int i;
138: PetscViewer viewer = (PetscViewer) dummy;
141: if (its) {
142: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(eps->comm);
143: PetscViewerASCIIPrintf(viewer,"%3d EPS nconv=%d Values (Errors)",its,nconv);
144: for (i=0;i<nest;i++) {
145: #if defined(PETSC_USE_COMPLEX)
146: PetscViewerASCIIPrintf(viewer," %g%+gi",PetscRealPart(eigr[i]),PetscImaginaryPart(eigr[i]));
147: #else
148: PetscViewerASCIIPrintf(viewer," %g",eigr[i]);
149: if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",eigi[i]); }
150: #endif
151: PetscViewerASCIIPrintf(viewer," (%10.8e)",errest[i]);
152: }
153: PetscViewerASCIIPrintf(viewer,"\n");
154: }
155: return(0);
156: }
160: PetscErrorCode EPSLGMonitor(EPS eps,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,int nest,void *monctx)
161: {
162: PetscViewer viewer = (PetscViewer) monctx;
163: PetscDraw draw;
164: PetscDrawLG lg;
166: PetscReal *x,*y;
167: int i,n = eps->nev;
168: #if !defined(PETSC_USE_COMPLEX)
169: int p;
170: PetscDraw draw1;
171: PetscDrawLG lg1;
172: #endif
176: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(eps->comm); }
178: PetscViewerDrawGetDraw(viewer,0,&draw);
179: PetscViewerDrawGetDrawLG(viewer,0,&lg);
180: if (!its) {
181: PetscDrawSetTitle(draw,"Error estimates");
182: PetscDrawSetDoubleBuffer(draw);
183: PetscDrawLGSetDimension(lg,n);
184: PetscDrawLGReset(lg);
185: PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
186: }
188: #if !defined(PETSC_USE_COMPLEX)
189: if (eps->ishermitian) {
190: PetscViewerDrawGetDraw(viewer,1,&draw1);
191: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
192: if (!its) {
193: PetscDrawSetTitle(draw1,"Approximate eigenvalues");
194: PetscDrawSetDoubleBuffer(draw1);
195: PetscDrawLGSetDimension(lg1,n);
196: PetscDrawLGReset(lg1);
197: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
198: }
199: }
200: #endif
202: PetscMalloc(sizeof(PetscReal)*n,&x);
203: PetscMalloc(sizeof(PetscReal)*n,&y);
204: for (i=0;i<n;i++) {
205: x[i] = (PetscReal) its;
206: if (errest[i] > 0.0) y[i] = log10(errest[i]); else y[i] = 0.0;
207: }
208: PetscDrawLGAddPoint(lg,x,y);
209: #if !defined(PETSC_USE_COMPLEX)
210: if (eps->ishermitian) {
211: PetscDrawLGAddPoint(lg1,x,eps->eigr);
212: PetscDrawGetPause(draw1,&p);
213: PetscDrawSetPause(draw1,0);
214: PetscDrawLGDraw(lg1);
215: PetscDrawSetPause(draw1,p);
216: }
217: #endif
218: PetscDrawLGDraw(lg);
219: PetscFree(x);
220: PetscFree(y);
221: return(0);
222: }