Actual source code: monitor.c

  1: /*
  2:       EPS routines related to monitors.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2011, Universitat 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/epsimpl.h>   /*I "slepceps.h" I*/

 28: /*
 29:    Runs the user provided monitor routines, if any.
 30: */
 31: PetscErrorCode EPSMonitor(EPS eps,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
 32: {
 34:   PetscInt       i,n = eps->numbermonitors;

 37:   for (i=0;i<n;i++) {
 38:     (*eps->monitor[i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[i]);
 39:   }
 40:   return(0);
 41: }

 45: /*@C
 46:    EPSMonitorSet - Sets an ADDITIONAL function to be called at every 
 47:    iteration to monitor the error estimates for each requested eigenpair.
 48:       
 49:    Logically Collective on EPS

 51:    Input Parameters:
 52: +  eps     - eigensolver context obtained from EPSCreate()
 53: .  monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
 54: .  mctx    - [optional] context for private data for the
 55:              monitor routine (use PETSC_NULL if no context is desired)
 56: -  monitordestroy - [optional] routine that frees monitor context
 57:           (may be PETSC_NULL)

 59:    Calling Sequence of monitor:
 60: $     monitor (EPS eps, int its, int nconv, PetscScalar *eigr, PetscScalar *eigi, PetscReal* errest, int nest, void *mctx)

 62: +  eps    - eigensolver context obtained from EPSCreate()
 63: .  its    - iteration number
 64: .  nconv  - number of converged eigenpairs
 65: .  eigr   - real part of the eigenvalues
 66: .  eigi   - imaginary part of the eigenvalues
 67: .  errest - relative error estimates for each eigenpair
 68: .  nest   - number of error estimates
 69: -  mctx   - optional monitoring context, as set by EPSMonitorSet()

 71:    Options Database Keys:
 72: +    -eps_monitor          - print only the first error estimate
 73: .    -eps_monitor_all      - print error estimates at each iteration
 74: .    -eps_monitor_conv     - print the eigenvalue approximations only when
 75:       convergence has been reached
 76: .    -eps_monitor_draw     - sets line graph monitor for the first unconverged
 77:       approximate eigenvalue
 78: .    -eps_monitor_draw_all - sets line graph monitor for all unconverged
 79:       approximate eigenvalue
 80: -    -eps_monitor_cancel   - cancels all monitors that have been hardwired into
 81:       a code by calls to EPSMonitorSet(), but does not cancel those set via
 82:       the options database.

 84:    Notes:  
 85:    Several different monitoring routines may be set by calling
 86:    EPSMonitorSet() multiple times; all will be called in the 
 87:    order in which they were set.

 89:    Level: intermediate

 91: .seealso: EPSMonitorFirst(), EPSMonitorAll(), EPSMonitorCancel()
 92: @*/
 93: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
 94: {
 97:   if (eps->numbermonitors >= MAXEPSMONITORS) {
 98:     SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
 99:   }
100:   eps->monitor[eps->numbermonitors]           = monitor;
101:   eps->monitorcontext[eps->numbermonitors]    = (void*)mctx;
102:   eps->monitordestroy[eps->numbermonitors++]  = monitordestroy;
103:   return(0);
104: }

108: /*@
109:    EPSMonitorCancel - Clears all monitors for an EPS object.

111:    Logically Collective on EPS

113:    Input Parameters:
114: .  eps - eigensolver context obtained from EPSCreate()

116:    Options Database Key:
117: .    -eps_monitor_cancel - Cancels all monitors that have been hardwired 
118:       into a code by calls to EPSMonitorSet(),
119:       but does not cancel those set via the options database.

121:    Level: intermediate

123: .seealso: EPSMonitorSet()
124: @*/
125: PetscErrorCode EPSMonitorCancel(EPS eps)
126: {
128:   PetscInt       i;

132:   for (i=0; i<eps->numbermonitors; i++) {
133:     if (eps->monitordestroy[i]) {
134:       (*eps->monitordestroy[i])(&eps->monitorcontext[i]);
135:     }
136:   }
137:   eps->numbermonitors = 0;
138:   return(0);
139: }

143: /*@C
144:    EPSGetMonitorContext - Gets the monitor context, as set by 
145:    EPSMonitorSet() for the FIRST monitor only.

147:    Not Collective

149:    Input Parameter:
150: .  eps - eigensolver context obtained from EPSCreate()

152:    Output Parameter:
153: .  ctx - monitor context

155:    Level: intermediate

157: .seealso: EPSMonitorSet()
158: @*/
159: PetscErrorCode EPSGetMonitorContext(EPS eps,void **ctx)
160: {
163:   *ctx = eps->monitorcontext[0];
164:   return(0);
165: }

169: /*@C
170:    EPSMonitorAll - Print the current approximate values and 
171:    error estimates at each iteration of the eigensolver.

173:    Collective on EPS

175:    Input Parameters:
176: +  eps    - 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: -  monctx - monitor context (contains viewer, can be PETSC_NULL)

185:    Level: intermediate

187: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorConverged()
188: @*/
189: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
190: {
192:   PetscInt       i;
193:   PetscScalar    er,ei;
194:   PetscViewer    viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);

197:   if (its) {
198:     PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
199:     PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D Values (Errors)",its,nconv);
200:     for (i=0;i<nest;i++) {
201:       er = eigr[i]; ei = eigi[i];
202:       STBackTransform(eps->OP,1,&er,&ei);
203: #if defined(PETSC_USE_COMPLEX)
204:       PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(er),PetscImaginaryPart(er));
205: #else
206:       PetscViewerASCIIPrintf(viewer," %G",er);
207:       if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",ei); }
208: #endif
209:       PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
210:     }
211:     PetscViewerASCIIPrintf(viewer,"\n");
212:     PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
213:   }
214:   return(0);
215: }

219: /*@C
220:    EPSMonitorFirst - Print the first approximate value and 
221:    error estimate at each iteration of the eigensolver.

223:    Collective on EPS

225:    Input Parameters:
226: +  eps    - eigensolver context
227: .  its    - iteration number
228: .  nconv  - number of converged eigenpairs so far
229: .  eigr   - real part of the eigenvalues
230: .  eigi   - imaginary part of the eigenvalues
231: .  errest - error estimates
232: .  nest   - number of error estimates to display
233: -  monctx - monitor context (contains viewer, can be PETSC_NULL)

235:    Level: intermediate

237: .seealso: EPSMonitorSet(), EPSMonitorAll(), EPSMonitorConverged()
238: @*/
239: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
240: {
242:   PetscScalar    er,ei;
243:   PetscViewer    viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);

246:   if (its && nconv<nest) {
247:     PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
248:     PetscViewerASCIIPrintf(viewer,"%3D EPS nconv=%D first unconverged value (error)",its,nconv);
249:     er = eigr[nconv]; ei = eigi[nconv];
250:     STBackTransform(eps->OP,1,&er,&ei);
251: #if defined(PETSC_USE_COMPLEX)
252:     PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(er),PetscImaginaryPart(er));
253: #else
254:     PetscViewerASCIIPrintf(viewer," %G",er);
255:     if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",ei); }
256: #endif
257:     PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
258:     PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
259:   }
260:   return(0);
261: }

265: /*@C
266:    EPSMonitorConverged - Print the approximate values and 
267:    error estimates as they converge.

269:    Collective on EPS

271:    Input Parameters:
272: +  eps    - eigensolver context
273: .  its    - iteration number
274: .  nconv  - number of converged eigenpairs so far
275: .  eigr   - real part of the eigenvalues
276: .  eigi   - imaginary part of the eigenvalues
277: .  errest - error estimates
278: .  nest   - number of error estimates to display
279: -  monctx - monitor context

281:    Note:
282:    The monitor context must contain a struct with a PetscViewer and a
283:    PetscInt. In Fortran, pass a PETSC_NULL_OBJECT.

285:    Level: intermediate

287: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorAll()
288: @*/
289: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
290: {
291:   PetscErrorCode   ierr;
292:   PetscInt         i;
293:   PetscScalar      er,ei;
294:   PetscViewer      viewer;
295:   SlepcConvMonitor ctx = (SlepcConvMonitor)monctx;

298:   if (!monctx) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Must provide a context for EPSMonitorConverged");
299:   if (!its) {
300:     ctx->oldnconv = 0;
301:   } else {
302:     viewer = ctx->viewer? ctx->viewer: PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);
303:     for (i=ctx->oldnconv;i<nconv;i++) {
304:       PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
305:       PetscViewerASCIIPrintf(viewer,"%3D EPS converged value (error) #%D",its,i);
306:       er = eigr[i]; ei = eigi[i];
307:       STBackTransform(eps->OP,1,&er,&ei);
308: #if defined(PETSC_USE_COMPLEX)
309:       PetscViewerASCIIPrintf(viewer," %G%+Gi",PetscRealPart(er),PetscImaginaryPart(er));
310: #else
311:       PetscViewerASCIIPrintf(viewer," %G",er);
312:       if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+Gi",ei); }
313: #endif
314:       PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
315:       PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
316:     }
317:     ctx->oldnconv = nconv;
318:   }
319:   return(0);
320: }

324: PetscErrorCode EPSMonitorLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
325: {
326:   PetscViewer    viewer = (PetscViewer) monctx;
327:   PetscDraw      draw,draw1;
328:   PetscDrawLG    lg,lg1;
330:   PetscReal      x,y,myeigr,p;
331:   PetscScalar    er,ei;

334:   if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); }
335:   PetscViewerDrawGetDraw(viewer,0,&draw);
336:   PetscViewerDrawGetDrawLG(viewer,0,&lg);
337:   if (!its) {
338:     PetscDrawSetTitle(draw,"Error estimates");
339:     PetscDrawSetDoubleBuffer(draw);
340:     PetscDrawLGSetDimension(lg,1);
341:     PetscDrawLGReset(lg);
342:     PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
343:   }

345:   /* In the hermitian case, the eigenvalues are real and can be plotted */
346:   if (eps->ishermitian) {
347:     PetscViewerDrawGetDraw(viewer,1,&draw1);
348:     PetscViewerDrawGetDrawLG(viewer,1,&lg1);
349:     if (!its) {
350:       PetscDrawSetTitle(draw1,"Approximate eigenvalues");
351:       PetscDrawSetDoubleBuffer(draw1);
352:       PetscDrawLGSetDimension(lg1,1);
353:       PetscDrawLGReset(lg1);
354:       PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
355:     }
356:   }

358:   x = (PetscReal) its;
359:   if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
360:   PetscDrawLGAddPoint(lg,&x,&y);
361:   if (eps->ishermitian) {
362:     er = eigr[nconv]; ei = eigi[nconv];
363:     STBackTransform(eps->OP,1,&er,&ei);
364:     myeigr = PetscRealPart(er);
365:     PetscDrawLGAddPoint(lg1,&x,&myeigr);
366:     PetscDrawGetPause(draw1,&p);
367:     PetscDrawSetPause(draw1,0);
368:     PetscDrawLGDraw(lg1);
369:     PetscDrawSetPause(draw1,p);
370:   }
371:   PetscDrawLGDraw(lg);
372:   return(0);
373: }

377: PetscErrorCode EPSMonitorLGAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *monctx)
378: {
379:   PetscViewer    viewer = (PetscViewer) monctx;
380:   PetscDraw      draw,draw1;
381:   PetscDrawLG    lg,lg1;
383:   PetscReal      *x,*y,*myeigr,p;
384:   PetscInt       i,n = PetscMin(eps->nev,255);
385:   PetscScalar    er,ei;

388:   if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)eps)->comm); }
389:   PetscViewerDrawGetDraw(viewer,0,&draw);
390:   PetscViewerDrawGetDrawLG(viewer,0,&lg);
391:   if (!its) {
392:     PetscDrawSetTitle(draw,"Error estimates");
393:     PetscDrawSetDoubleBuffer(draw);
394:     PetscDrawLGSetDimension(lg,n);
395:     PetscDrawLGReset(lg);
396:     PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
397:   }

399:   /* In the hermitian case, the eigenvalues are real and can be plotted */
400:   if (eps->ishermitian) {
401:     PetscViewerDrawGetDraw(viewer,1,&draw1);
402:     PetscViewerDrawGetDrawLG(viewer,1,&lg1);
403:     if (!its) {
404:       PetscDrawSetTitle(draw1,"Approximate eigenvalues");
405:       PetscDrawSetDoubleBuffer(draw1);
406:       PetscDrawLGSetDimension(lg1,n);
407:       PetscDrawLGReset(lg1);
408:       PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
409:     }
410:   }

412:   PetscMalloc(sizeof(PetscReal)*n,&x);
413:   PetscMalloc(sizeof(PetscReal)*n,&y);
414:   for (i=0;i<n;i++) {
415:     x[i] = (PetscReal) its;
416:     if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
417:     else y[i] = 0.0;
418:   }
419:   PetscDrawLGAddPoint(lg,x,y);
420:   if (eps->ishermitian) {
421:     PetscMalloc(sizeof(PetscReal)*n,&myeigr);
422:     for(i=0;i<n;i++) {
423:       er = eigr[i]; ei = eigi[i];
424:       STBackTransform(eps->OP,1,&er,&ei);
425:       if (i < nest) myeigr[i] = PetscRealPart(er);
426:       else myeigr[i] = 0.0;
427:     }
428:     PetscDrawLGAddPoint(lg1,x,myeigr);
429:     PetscDrawGetPause(draw1,&p);
430:     PetscDrawSetPause(draw1,0);
431:     PetscDrawLGDraw(lg1);
432:     PetscDrawSetPause(draw1,p);
433:     PetscFree(myeigr);
434:   }
435:   PetscDrawLGDraw(lg);
436:   PetscFree(x);
437:   PetscFree(y);
438:   return(0);
439: }