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-2009, Universidad 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

 28: /*@C
 29:    SVDMonitorSet - Sets an ADDITIONAL function to be called at every 
 30:    iteration to monitor the error estimates for each requested singular triplet.
 31:       
 32:    Collective on SVD

 34:    Input Parameters:
 35: +  svd     - singular value solver context obtained from SVDCreate()
 36: .  monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
 37: -  mctx    - [optional] context for private data for the
 38:              monitor routine (use PETSC_NULL if no context is desired)

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

 43: +  svd    - singular value solver context obtained from SVDCreate()
 44: .  its    - iteration number
 45: .  nconv  - number of converged singular triplets
 46: .  sigma  - singular values
 47: .  errest - relative error estimates for each singular triplet
 48: .  nest   - number of error estimates
 49: -  mctx   - optional monitoring context, as set by SVDMonitorSet()

 51:    Options Database Keys:
 52: +    -svd_monitor        - print error estimates at each iteration
 53: .    -svd_monitor_draw   - sets line graph monitor
 54: -    -svd_monitor_cancel - cancels all monitors that have been hardwired into
 55:       a code by calls to SVDMonitorSet(), but does not cancel those set via
 56:       the options database.

 58:    Notes:  
 59:    Several different monitoring routines may be set by calling
 60:    SVDMonitorSet() multiple times; all will be called in the 
 61:    order in which they were set.

 63:    Level: intermediate

 65: .seealso: SVDMonitorDefault(), SVDMonitorCancel()
 66: @*/
 67: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),
 68:                              void *mctx,PetscErrorCode (*monitordestroy)(void*))
 69: {
 72:   if (svd->numbermonitors >= MAXSVDMONITORS) {
 73:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
 74:   }
 75:   svd->monitor[svd->numbermonitors]           = monitor;
 76:   svd->monitorcontext[svd->numbermonitors]    = (void*)mctx;
 77:   svd->monitordestroy[svd->numbermonitors++]  = monitordestroy;
 78:   return(0);
 79: }

 83: /*@
 84:    SVDMonitorCancel - Clears all monitors for an SVD object.

 86:    Collective on SVD

 88:    Input Parameters:
 89: .  svd - singular value solver context obtained from SVDCreate()

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

 96:    Level: intermediate

 98: .seealso: SVDMonitorCancel()
 99: @*/
100: PetscErrorCode SVDMonitorCancel(SVD svd)
101: {
103:   PetscInt       i;

107:   for (i=0; i<svd->numbermonitors; i++) {
108:     if (svd->monitordestroy[i]) {
109:       (*svd->monitordestroy[i])(svd->monitorcontext[i]);
110:     }
111:   }
112:   svd->numbermonitors = 0;
113:   return(0);
114: }

118: /*@C
119:    SVDGetMonitorContext - Gets the monitor context, as set by 
120:    SVDMonitorSet() for the FIRST monitor only.

122:    Not Collective

124:    Input Parameter:
125: .  svd - singular value solver context obtained from SVDCreate()

127:    Output Parameter:
128: .  ctx - monitor context

130:    Level: intermediate

132: .seealso: SVDMonitorSet(), SVDMonitorDefault()
133: @*/
134: PetscErrorCode SVDGetMonitorContext(SVD svd, void **ctx)
135: {
138:   *ctx =      (svd->monitorcontext[0]);
139:   return(0);
140: }

144: /*@C
145:    SVDDefaultMonitor - Print the current approximate values and 
146:    error estimates at each iteration of the singular value solver.

148:    Collective on SVD

150:    Input Parameters:
151: +  svd    - singular value solver context
152: .  its    - iteration number
153: .  nconv  - number of converged singular triplets so far
154: .  sigma  - singular values
155: .  errest - error estimates
156: .  nest   - number of error estimates to display
157: -  dummy  - unused monitor context 

159:    Level: intermediate

161: .seealso: SVDMonitorSet()
162: @*/
163: PetscErrorCode SVDMonitorDefault(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *dummy)
164: {
165:   PetscErrorCode          ierr;
166:   PetscInt                i;
167:   PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;

170:   if (its) {
171:     if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)svd)->comm,"stdout",0,&viewer);}
172:     PetscViewerASCIIMonitorPrintf(viewer,"%3d SVD nconv=%d Values (Errors)",its,nconv);
173:     for (i=0;i<nest;i++) {
174:       PetscViewerASCIIMonitorPrintf(viewer," %g (%10.8e)",sigma[i],errest[i]);
175:     }
176:     PetscViewerASCIIMonitorPrintf(viewer,"\n");
177:     if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
178:   }
179:   return(0);
180: }

184: PetscErrorCode SVDMonitorLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
185: {
186:   PetscViewer    viewer = (PetscViewer) monctx;
187:   PetscDraw      draw;
188:   PetscDrawLG    lg;
190:   PetscReal      *x,*y;
191:   int            p,n = svd->nsv;
192:   PetscInt       i;
193:   PetscDraw      draw1;
194:   PetscDrawLG    lg1;


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

200:   PetscViewerDrawGetDraw(viewer,0,&draw);
201:   PetscViewerDrawGetDrawLG(viewer,0,&lg);
202:   PetscViewerDrawGetDraw(viewer,1,&draw1);
203:   PetscViewerDrawGetDrawLG(viewer,1,&lg1);

205:   if (!its) {
206:     PetscDrawSetTitle(draw,"Error estimates");
207:     PetscDrawSetDoubleBuffer(draw);
208:     PetscDrawLGSetDimension(lg,n);
209:     PetscDrawLGReset(lg);
210:     PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);

212:     PetscDrawSetTitle(draw1,"Approximate singular values");
213:     PetscDrawSetDoubleBuffer(draw1);
214:     PetscDrawLGSetDimension(lg1,n);
215:     PetscDrawLGReset(lg1);
216:     PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
217:   }

219:   PetscMalloc(sizeof(PetscReal)*n,&x);
220:   PetscMalloc(sizeof(PetscReal)*n,&y);
221:   for (i=0;i<n;i++) {
222:     x[i] = (PetscReal) its;
223:     if (errest[i] > 0.0) y[i] = log10(errest[i]); else y[i] = 0.0;
224:   }
225:   PetscDrawLGAddPoint(lg,x,y);

227:   PetscDrawLGAddPoint(lg1,x,svd->sigma);
228:   PetscDrawGetPause(draw1,&p);
229:   PetscDrawSetPause(draw1,0);
230:   PetscDrawLGDraw(lg1);
231:   PetscDrawSetPause(draw1,p);
232: 
233:   PetscDrawLGDraw(lg);
234:   PetscFree(x);
235:   PetscFree(y);
236:   return(0);
237: }