Actual source code: svdmon.c

  1: /*
  2:       SVD 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/svdimpl.h>   /*I "slepcsvd.h" I*/

 28: /*
 29:    Runs the user provided monitor routines, if any.
 30: */
 31: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest)
 32: {
 34:   PetscInt       i,n = svd->numbermonitors;

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

 45: /*@C
 46:    SVDMonitorSet - Sets an ADDITIONAL function to be called at every 
 47:    iteration to monitor the error estimates for each requested singular triplet.
 48:       
 49:    Collective on SVD

 51:    Input Parameters:
 52: +  svd     - singular value solver context obtained from SVDCreate()
 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)

 57:    Calling Sequence of monitor:
 58: $     monitor (SVD svd, PetscInt its, PetscInt nconv, PetscReal *sigma, PetscReal* errest, PetscInt nest, void *mctx)

 60: +  svd    - singular value solver context obtained from SVDCreate()
 61: .  its    - iteration number
 62: .  nconv  - number of converged singular triplets
 63: .  sigma  - singular values
 64: .  errest - relative error estimates for each singular triplet
 65: .  nest   - number of error estimates
 66: -  mctx   - optional monitoring context, as set by SVDMonitorSet()

 68:    Options Database Keys:
 69: +    -svd_monitor          - print only the first error estimate
 70: .    -svd_monitor_all      - print error estimates at each iteration
 71: .    -svd_monitor_conv     - print the singular value approximations only when
 72:       convergence has been reached
 73: .    -svd_monitor_draw     - sets line graph monitor for the first unconverged
 74:       approximate singular value
 75: .    -svd_monitor_draw_all - sets line graph monitor for all unconverged
 76:       approximate singular value
 77: -    -svd_monitor_cancel   - cancels all monitors that have been hardwired into
 78:       a code by calls to SVDMonitorSet(), but does not cancel those set via
 79:       the options database.

 81:    Notes:  
 82:    Several different monitoring routines may be set by calling
 83:    SVDMonitorSet() multiple times; all will be called in the 
 84:    order in which they were set.

 86:    Level: intermediate

 88: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorCancel()
 89: @*/
 90: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
 91: {
 94:   if (svd->numbermonitors >= MAXSVDMONITORS) {
 95:     SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
 96:   }
 97:   svd->monitor[svd->numbermonitors]           = monitor;
 98:   svd->monitorcontext[svd->numbermonitors]    = (void*)mctx;
 99:   svd->monitordestroy[svd->numbermonitors++]  = monitordestroy;
100:   return(0);
101: }

105: /*@
106:    SVDMonitorCancel - Clears all monitors for an SVD object.

108:    Collective on SVD

110:    Input Parameters:
111: .  svd - singular value solver context obtained from SVDCreate()

113:    Options Database Key:
114: .    -svd_monitor_cancel - Cancels all monitors that have been hardwired 
115:       into a code by calls to SVDMonitorSet(),
116:       but does not cancel those set via the options database.

118:    Level: intermediate

120: .seealso: SVDMonitorSet()
121: @*/
122: PetscErrorCode SVDMonitorCancel(SVD svd)
123: {
125:   PetscInt       i;

129:   for (i=0; i<svd->numbermonitors; i++) {
130:     if (svd->monitordestroy[i]) {
131:       (*svd->monitordestroy[i])(&svd->monitorcontext[i]);
132:     }
133:   }
134:   svd->numbermonitors = 0;
135:   return(0);
136: }

140: /*@C
141:    SVDGetMonitorContext - Gets the monitor context, as set by 
142:    SVDMonitorSet() for the FIRST monitor only.

144:    Not Collective

146:    Input Parameter:
147: .  svd - singular value solver context obtained from SVDCreate()

149:    Output Parameter:
150: .  ctx - monitor context

152:    Level: intermediate

154: .seealso: SVDMonitorSet()
155: @*/
156: PetscErrorCode SVDGetMonitorContext(SVD svd,void **ctx)
157: {
160:   *ctx = (svd->monitorcontext[0]);
161:   return(0);
162: }

166: /*@C
167:    SVDMonitorAll - Print the current approximate values and 
168:    error estimates at each iteration of the singular value solver.

170:    Collective on SVD

172:    Input Parameters:
173: +  svd    - singular value solver context
174: .  its    - iteration number
175: .  nconv  - number of converged singular triplets so far
176: .  sigma  - singular values
177: .  errest - error estimates
178: .  nest   - number of error estimates to display
179: -  monctx - monitor context (contains viewer, can be PETSC_NULL)

181:    Level: intermediate

183: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConverged()
184: @*/
185: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
186: {
188:   PetscInt       i;
189:   PetscViewer    viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)svd)->comm);

192:   if (its) {
193:     PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
194:     PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D Values (Errors)",its,nconv);
195:     for (i=0;i<nest;i++) {
196:       PetscViewerASCIIPrintf(viewer," %G (%10.8e)",sigma[i],(double)errest[i]);
197:     }
198:     PetscViewerASCIIPrintf(viewer,"\n");
199:     PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
200:   }
201:   return(0);
202: }

206: /*@C
207:    SVDMonitorFirst - Print the first unconverged approximate values and 
208:    error estimates at each iteration of the singular value solver.

210:    Collective on SVD

212:    Input Parameters:
213: +  svd    - singular value solver context
214: .  its    - iteration number
215: .  nconv  - number of converged singular triplets so far
216: .  sigma  - singular values
217: .  errest - error estimates
218: .  nest   - number of error estimates to display
219: -  monctx - monitor context (contains viewer, can be PETSC_NULL)

221:    Level: intermediate

223: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConverged()
224: @*/
225: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
226: {
228:   PetscViewer    viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)svd)->comm);

231:   if (its && nconv<nest) {
232:     PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
233:     PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D first unconverged value (error)",its,nconv);
234:     PetscViewerASCIIPrintf(viewer," %G (%10.8e)\n",sigma[nconv],(double)errest[nconv]);
235:     PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
236:   }
237:   return(0);
238: }

242: /*@C
243:    SVDMonitorConverged - Print the approximate values and error estimates as they converge.

245:    Collective on SVD

247:    Input Parameters:
248: +  svd    - singular value solver context
249: .  its    - iteration number
250: .  nconv  - number of converged singular triplets so far
251: .  sigma  - singular values
252: .  errest - error estimates
253: .  nest   - number of error estimates to display
254: -  monctx - monitor context 

256:    Note:
257:    The monitor context must contain a struct with a PetscViewer and a
258:    PetscInt. In Fortran, pass a PETSC_NULL_OBJECT.

260:    Level: intermediate

262: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorAll()
263: @*/
264: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
265: {
266:   PetscErrorCode   ierr;
267:   PetscInt         i;
268:   PetscViewer      viewer;
269:   SlepcConvMonitor ctx = (SlepcConvMonitor) monctx;

272:   if (!monctx) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_WRONG,"Must provide a context for SVDMonitorConverged");
273:   if (!its) {
274:     ctx->oldnconv = 0;
275:   } else {
276:     viewer = ctx->viewer? ctx->viewer: PETSC_VIEWER_STDOUT_(((PetscObject)svd)->comm);
277:     for (i=ctx->oldnconv;i<nconv;i++) {
278:       PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
279:       PetscViewerASCIIPrintf(viewer,"%3D SVD converged value (error) #%D",its,i);
280:       PetscViewerASCIIPrintf(viewer," %G (%10.8e)\n",sigma[i],(double)errest[i]);
281:       PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
282:     }
283:     ctx->oldnconv = nconv;
284:   }
285:   return(0);
286: }

290: PetscErrorCode SVDMonitorLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
291: {
292:   PetscViewer    viewer = (PetscViewer) monctx;
293:   PetscDraw      draw,draw1;
294:   PetscDrawLG    lg,lg1;
296:   PetscReal      x,y,p;

299:   if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); }
300:   PetscViewerDrawGetDraw(viewer,0,&draw);
301:   PetscViewerDrawGetDrawLG(viewer,0,&lg);
302:   PetscViewerDrawGetDraw(viewer,1,&draw1);
303:   PetscViewerDrawGetDrawLG(viewer,1,&lg1);

305:   if (!its) {
306:     PetscDrawSetTitle(draw,"Error estimates");
307:     PetscDrawSetDoubleBuffer(draw);
308:     PetscDrawLGSetDimension(lg,1);
309:     PetscDrawLGReset(lg);
310:     PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);

312:     PetscDrawSetTitle(draw1,"Approximate singular values");
313:     PetscDrawSetDoubleBuffer(draw1);
314:     PetscDrawLGSetDimension(lg1,1);
315:     PetscDrawLGReset(lg1);
316:     PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
317:   }

319:   x = (PetscReal) its;
320:   if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
321:   PetscDrawLGAddPoint(lg,&x,&y);

323:   PetscDrawLGAddPoint(lg1,&x,svd->sigma);
324:   PetscDrawGetPause(draw1,&p);
325:   PetscDrawSetPause(draw1,0);
326:   PetscDrawLGDraw(lg1);
327:   PetscDrawSetPause(draw1,p);
328: 
329:   PetscDrawLGDraw(lg);
330:   return(0);
331: }

335: PetscErrorCode SVDMonitorLGAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
336: {
337:   PetscViewer    viewer = (PetscViewer) monctx;
338:   PetscDraw      draw,draw1;
339:   PetscDrawLG    lg,lg1;
341:   PetscReal      *x,*y,p;
342:   PetscInt       i,n = PetscMin(svd->nsv,255);

345:   if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); }
346:   PetscViewerDrawGetDraw(viewer,0,&draw);
347:   PetscViewerDrawGetDrawLG(viewer,0,&lg);
348:   PetscViewerDrawGetDraw(viewer,1,&draw1);
349:   PetscViewerDrawGetDrawLG(viewer,1,&lg1);

351:   if (!its) {
352:     PetscDrawSetTitle(draw,"Error estimates");
353:     PetscDrawSetDoubleBuffer(draw);
354:     PetscDrawLGSetDimension(lg,n);
355:     PetscDrawLGReset(lg);
356:     PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);

358:     PetscDrawSetTitle(draw1,"Approximate singular values");
359:     PetscDrawSetDoubleBuffer(draw1);
360:     PetscDrawLGSetDimension(lg1,n);
361:     PetscDrawLGReset(lg1);
362:     PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
363:   }

365:   PetscMalloc(sizeof(PetscReal)*n,&x);
366:   PetscMalloc(sizeof(PetscReal)*n,&y);
367:   for (i=0;i<n;i++) {
368:     x[i] = (PetscReal) its;
369:     if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
370:     else y[i] = 0.0;
371:   }
372:   PetscDrawLGAddPoint(lg,x,y);

374:   PetscDrawLGAddPoint(lg1,x,svd->sigma);
375:   PetscDrawGetPause(draw1,&p);
376:   PetscDrawSetPause(draw1,0);
377:   PetscDrawLGDraw(lg1);
378:   PetscDrawSetPause(draw1,p);
379: 
380:   PetscDrawLGDraw(lg);
381:   PetscFree(x);
382:   PetscFree(y);
383:   return(0);
384: }