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-2007, Universidad Politecnica de Valencia, Spain

  8:       This file is part of SLEPc. See the README file for conditions of use
  9:       and additional information.
 10:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 11: */

 13:  #include src/eps/epsimpl.h

 17: /*@C
 18:    EPSMonitorSet - Sets an ADDITIONAL function to be called at every 
 19:    iteration to monitor the error estimates for each requested eigenpair.
 20:       
 21:    Collective on EPS

 23:    Input Parameters:
 24: +  eps     - eigensolver context obtained from EPSCreate()
 25: .  monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
 26: .  mctx    - [optional] context for private data for the
 27:              monitor routine (use PETSC_NULL if no context is desired)
 28: -  monitordestroy - [optional] routine that frees monitor context
 29:           (may be PETSC_NULL)

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

 34: +  eps    - eigensolver context obtained from EPSCreate()
 35: .  its    - iteration number
 36: .  nconv  - number of converged eigenpairs
 37: .  eigr   - real part of the eigenvalues
 38: .  eigi   - imaginary part of the eigenvalues
 39: .  errest - relative error estimates for each eigenpair
 40: .  nest   - number of error estimates
 41: -  mctx   - optional monitoring context, as set by EPSMonitorSet()

 43:    Options Database Keys:
 44: +    -eps_monitor        - print error estimates at each iteration
 45: .    -eps_monitor_draw   - sets line graph monitor
 46: -    -eps_monitor_cancel - cancels all monitors that have been hardwired into
 47:       a code by calls to EPSMonitorSet(), but does not cancel those set via
 48:       the options database.

 50:    Notes:  
 51:    Several different monitoring routines may be set by calling
 52:    EPSMonitorSet() multiple times; all will be called in the 
 53:    order in which they were set.

 55:    Level: intermediate

 57: .seealso: EPSMonitorDefault(), EPSMonitorCancel()
 58: @*/
 59: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,int,int,PetscScalar*,PetscScalar*,PetscReal*,int,void*),
 60:                              void *mctx,PetscErrorCode (*monitordestroy)(void*))
 61: {
 64:   if (eps->numbermonitors >= MAXEPSMONITORS) {
 65:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
 66:   }
 67:   eps->monitor[eps->numbermonitors]           = monitor;
 68:   eps->monitorcontext[eps->numbermonitors]    = (void*)mctx;
 69:   eps->monitordestroy[eps->numbermonitors++]  = monitordestroy;
 70:   return(0);
 71: }

 75: /*@
 76:    EPSMonitorCancel - Clears all monitors for an EPS object.

 78:    Collective on EPS

 80:    Input Parameters:
 81: .  eps - eigensolver context obtained from EPSCreate()

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

 88:    Level: intermediate

 90: .seealso: EPSMonitorSet()
 91: @*/
 92: PetscErrorCode EPSMonitorCancel(EPS eps)
 93: {
 95:   PetscInt       i;

 99:   for (i=0; i<eps->numbermonitors; i++) {
100:     if (eps->monitordestroy[i]) {
101:       (*eps->monitordestroy[i])(eps->monitorcontext[i]);
102:     }
103:   }
104:   eps->numbermonitors = 0;
105:   return(0);
106: }

110: /*@C
111:    EPSGetMonitorContext - Gets the monitor context, as set by 
112:    EPSSetMonitor() for the FIRST monitor only.

114:    Not Collective

116:    Input Parameter:
117: .  eps - eigensolver context obtained from EPSCreate()

119:    Output Parameter:
120: .  ctx - monitor context

122:    Level: intermediate

124: .seealso: EPSSetMonitor(), EPSDefaultMonitor()
125: @*/
126: PetscErrorCode EPSGetMonitorContext(EPS eps, void **ctx)
127: {
130:   *ctx =      (eps->monitorcontext[0]);
131:   return(0);
132: }

136: /*@C
137:    EPSMonitorDefault - Print the current approximate values and 
138:    error estimates at each iteration of the eigensolver.

140:    Collective on EPS

142:    Input Parameters:
143: +  eps    - eigensolver context
144: .  its    - iteration number
145: .  nconv  - number of converged eigenpairs so far
146: .  eigr   - real part of the eigenvalues
147: .  eigi   - imaginary part of the eigenvalues
148: .  errest - error estimates
149: .  nest   - number of error estimates to display
150: -  dummy  - unused monitor context 

152:    Level: intermediate

154: .seealso: EPSMonitorSet()
155: @*/
156: PetscErrorCode EPSMonitorDefault(EPS eps,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,int nest,void *dummy)
157: {
159:   int            i;
160:   PetscViewer    viewer = (PetscViewer) dummy;

163:   if (its) {
164:     if (!viewer) viewer = PETSC_VIEWER_STDOUT_(eps->comm);
165:     PetscViewerASCIIPrintf(viewer,"%3d EPS nconv=%d Values (Errors)",its,nconv);
166:     for (i=0;i<nest;i++) {
167: #if defined(PETSC_USE_COMPLEX)
168:       PetscViewerASCIIPrintf(viewer," %g%+gi",PetscRealPart(eigr[i]),PetscImaginaryPart(eigr[i]));
169: #else
170:       PetscViewerASCIIPrintf(viewer," %g",eigr[i]);
171:       if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",eigi[i]); }
172: #endif
173:       PetscViewerASCIIPrintf(viewer," (%10.8e)",errest[i]);
174:     }
175:     PetscViewerASCIIPrintf(viewer,"\n");
176:   }
177:   return(0);
178: }

182: PetscErrorCode EPSMonitorLG(EPS eps,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,int nest,void *monctx)
183: {
184:   PetscViewer    viewer = (PetscViewer) monctx;
185:   PetscDraw      draw;
186:   PetscDrawLG    lg;
188:   PetscReal      *x,*y;
189:   int            i,n = eps->nev;
190: #if !defined(PETSC_USE_COMPLEX)
191:   int            p;
192:   PetscDraw      draw1;
193:   PetscDrawLG    lg1;
194: #endif


198:   if (!viewer) { viewer = PETSC_VIEWER_DRAW_(eps->comm); }

200:   PetscViewerDrawGetDraw(viewer,0,&draw);
201:   PetscViewerDrawGetDrawLG(viewer,0,&lg);
202:   if (!its) {
203:     PetscDrawSetTitle(draw,"Error estimates");
204:     PetscDrawSetDoubleBuffer(draw);
205:     PetscDrawLGSetDimension(lg,n);
206:     PetscDrawLGReset(lg);
207:     PetscDrawLGSetLimits(lg,0,1.0,log10(eps->tol)-2,0.0);
208:   }

210: #if !defined(PETSC_USE_COMPLEX)
211:   if (eps->ishermitian) {
212:     PetscViewerDrawGetDraw(viewer,1,&draw1);
213:     PetscViewerDrawGetDrawLG(viewer,1,&lg1);
214:     if (!its) {
215:       PetscDrawSetTitle(draw1,"Approximate eigenvalues");
216:       PetscDrawSetDoubleBuffer(draw1);
217:       PetscDrawLGSetDimension(lg1,n);
218:       PetscDrawLGReset(lg1);
219:       PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
220:     }
221:   }
222: #endif

224:   PetscMalloc(sizeof(PetscReal)*n,&x);
225:   PetscMalloc(sizeof(PetscReal)*n,&y);
226:   for (i=0;i<n;i++) {
227:     x[i] = (PetscReal) its;
228:     if (errest[i] > 0.0) y[i] = log10(errest[i]); else y[i] = 0.0;
229:   }
230:   PetscDrawLGAddPoint(lg,x,y);
231: #if !defined(PETSC_USE_COMPLEX)
232:   if (eps->ishermitian) {
233:     PetscDrawLGAddPoint(lg1,x,eps->eigr);
234:     PetscDrawGetPause(draw1,&p);
235:     PetscDrawSetPause(draw1,0);
236:     PetscDrawLGDraw(lg1);
237:     PetscDrawSetPause(draw1,p);
238:   }
239: #endif  
240:   PetscDrawLGDraw(lg);
241:   PetscFree(x);
242:   PetscFree(y);
243:   return(0);
244: }