Actual source code: svdmon.c

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    SVD routines related to monitors
 12: */

 14: #include <slepc/private/svdimpl.h>   /*I "slepcsvd.h" I*/
 15: #include <petscdraw.h>

 17: /*
 18:    Runs the user provided monitor routines, if any.
 19: */
 20: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest)
 21: {
 23:   PetscInt       i,n = svd->numbermonitors;

 26:   for (i=0;i<n;i++) {
 27:     (*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]);
 28:   }
 29:   return(0);
 30: }

 32: /*@C
 33:    SVDMonitorSet - Sets an ADDITIONAL function to be called at every
 34:    iteration to monitor the error estimates for each requested singular triplet.

 36:    Collective on SVD

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

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

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

 55:    Options Database Keys:
 56: +    -svd_monitor        - print only the first error estimate
 57: .    -svd_monitor_all    - print error estimates at each iteration
 58: .    -svd_monitor_conv   - print the singular value approximations only when
 59:       convergence has been reached
 60: .    -svd_monitor_lg     - sets line graph monitor for the first unconverged
 61:       approximate singular value
 62: .    -svd_monitor_lg_all - sets line graph monitor for all unconverged
 63:       approximate singular values
 64: -    -svd_monitor_cancel - cancels all monitors that have been hardwired into
 65:       a code by calls to SVDMonitorSet(), but does not cancel those set via
 66:       the options database.

 68:    Notes:
 69:    Several different monitoring routines may be set by calling
 70:    SVDMonitorSet() multiple times; all will be called in the
 71:    order in which they were set.

 73:    Level: intermediate

 75: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorCancel()
 76: @*/
 77: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
 78: {
 81:   if (svd->numbermonitors >= MAXSVDMONITORS) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
 82:   svd->monitor[svd->numbermonitors]           = monitor;
 83:   svd->monitorcontext[svd->numbermonitors]    = (void*)mctx;
 84:   svd->monitordestroy[svd->numbermonitors++]  = monitordestroy;
 85:   return(0);
 86: }

 88: /*@
 89:    SVDMonitorCancel - Clears all monitors for an SVD object.

 91:    Collective on SVD

 93:    Input Parameters:
 94: .  svd - singular value solver context obtained from SVDCreate()

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

101:    Level: intermediate

103: .seealso: SVDMonitorSet()
104: @*/
105: PetscErrorCode SVDMonitorCancel(SVD svd)
106: {
108:   PetscInt       i;

112:   for (i=0; i<svd->numbermonitors; i++) {
113:     if (svd->monitordestroy[i]) {
114:       (*svd->monitordestroy[i])(&svd->monitorcontext[i]);
115:     }
116:   }
117:   svd->numbermonitors = 0;
118:   return(0);
119: }

121: /*@C
122:    SVDGetMonitorContext - Gets the monitor context, as set by
123:    SVDMonitorSet() for the FIRST monitor only.

125:    Not Collective

127:    Input Parameter:
128: .  svd - singular value solver context obtained from SVDCreate()

130:    Output Parameter:
131: .  ctx - monitor context

133:    Level: intermediate

135: .seealso: SVDMonitorSet()
136: @*/
137: PetscErrorCode SVDGetMonitorContext(SVD svd,void **ctx)
138: {
141:   *ctx = svd->monitorcontext[0];
142:   return(0);
143: }

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

149:    Collective on SVD

151:    Input Parameters:
152: +  svd    - singular value solver context
153: .  its    - iteration number
154: .  nconv  - number of converged singular triplets so far
155: .  sigma  - singular values
156: .  errest - error estimates
157: .  nest   - number of error estimates to display
158: -  vf     - viewer and format for monitoring

160:    Level: intermediate

162: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConverged()
163: @*/
164: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
165: {
167:   PetscInt       i;
168:   PetscViewer    viewer;

173:   viewer = vf->viewer;
175:   PetscViewerPushFormat(viewer,vf->format);
176:   PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
177:   if (its==1 && ((PetscObject)svd)->prefix) {
178:     PetscViewerASCIIPrintf(viewer,"  Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
179:   }
180:   PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D Values (Errors)",its,nconv);
181:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
182:   for (i=0;i<nest;i++) {
183:     PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]);
184:   }
185:   PetscViewerASCIIPrintf(viewer,"\n");
186:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
187:   PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
188:   PetscViewerPopFormat(viewer);
189:   return(0);
190: }

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

196:    Collective on SVD

198:    Input Parameters:
199: +  svd    - singular value solver context
200: .  its    - iteration number
201: .  nconv  - number of converged singular triplets so far
202: .  sigma  - singular values
203: .  errest - error estimates
204: .  nest   - number of error estimates to display
205: -  vf     - viewer and format for monitoring

207:    Level: intermediate

209: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConverged()
210: @*/
211: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
212: {
214:   PetscViewer    viewer;

219:   viewer = vf->viewer;
221:   if (its==1 && ((PetscObject)svd)->prefix) {
222:     PetscViewerASCIIPrintf(viewer,"  Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
223:   }
224:   if (nconv<nest) {
225:     PetscViewerPushFormat(viewer,vf->format);
226:     PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
227:     PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D first unconverged value (error)",its,nconv);
228:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
229:     PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]);
230:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
231:     PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
232:     PetscViewerPopFormat(viewer);
233:   }
234:   return(0);
235: }

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

240:    Collective on SVD

242:    Input Parameters:
243: +  svd    - singular value solver context
244: .  its    - iteration number
245: .  nconv  - number of converged singular triplets so far
246: .  sigma  - singular values
247: .  errest - error estimates
248: .  nest   - number of error estimates to display
249: -  ctx    - monitor context

251:    Level: intermediate

253: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorAll()
254: @*/
255: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,SlepcConvMonitor ctx)
256: {
258:   PetscInt       i;
259:   PetscViewer    viewer;

264:   viewer = ctx->viewer;
266:   if (its==1 && ((PetscObject)svd)->prefix) {
267:     PetscViewerASCIIPrintf(viewer,"  Convergence history for %s solve.\n",((PetscObject)svd)->prefix);
268:   }
269:   if (its==1) ctx->oldnconv = 0;
270:   if (ctx->oldnconv!=nconv) {
271:     PetscViewerPushFormat(viewer,ctx->format);
272:     PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
273:     for (i=ctx->oldnconv;i<nconv;i++) {
274:       PetscViewerASCIIPrintf(viewer,"%3D SVD converged value (error) #%D",its,i);
275:       PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
276:       PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]);
277:       PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
278:     }
279:     PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
280:     PetscViewerPopFormat(viewer);
281:     ctx->oldnconv = nconv;
282:   }
283:   return(0);
284: }

286: /*@C
287:    SVDMonitorLGCreate - Creates a line graph context for use with
288:    SVD to monitor convergence.

290:    Collective on MPI_Comm

292:    Input Parameters:
293: +  comm - communicator context
294: .  host - the X display to open, or null for the local machine
295: .  label - the title to put in the title bar
296: .  x, y - the screen coordinates of the upper left coordinate of
297:           the window
298: -  m, n - the screen width and height in pixels

300:    Output Parameter:
301: .  lgctx - the drawing context

303:    Options Database Keys:
304: +  -svd_monitor_lg - Sets line graph monitor for the first residual
305: -  -svd_monitor_lg_all - Sets line graph monitor for all residuals

307:    Notes:
308:    Use PetscDrawLGDestroy() to destroy this line graph.

310:    Level: intermediate

312: .seealso: SVDMonitorSet()
313: @*/
314: PetscErrorCode SVDMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
315: {
316:   PetscDraw      draw;
317:   PetscDrawLG    lg;

321:   PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
322:   PetscDrawSetFromOptions(draw);
323:   PetscDrawLGCreate(draw,1,&lg);
324:   PetscDrawLGSetFromOptions(lg);
325:   PetscDrawDestroy(&draw);
326:   *lgctx = lg;
327:   return(0);
328: }

330: PetscErrorCode SVDMonitorLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *ctx)
331: {
332:   PetscDrawLG    lg = (PetscDrawLG)ctx;
333:   PetscReal      x,y;

338:   if (its==1) {
339:     PetscDrawLGReset(lg);
340:     PetscDrawLGSetDimension(lg,1);
341:     PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(svd->tol)-2,0.0);
342:   }
343:   x = (PetscReal)its;
344:   if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
345:   else y = 0.0;
346:   PetscDrawLGAddPoint(lg,&x,&y);
347:   if (its <= 20 || !(its % 5) || svd->reason) {
348:     PetscDrawLGDraw(lg);
349:     PetscDrawLGSave(lg);
350:   }
351:   return(0);
352: }

354: PetscErrorCode SVDMonitorLGAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *ctx)
355: {
356:   PetscDrawLG    lg = (PetscDrawLG)ctx;
357:   PetscInt       i,n = PetscMin(svd->nsv,255);
358:   PetscReal      *x,*y;

363:   if (its==1) {
364:     PetscDrawLGReset(lg);
365:     PetscDrawLGSetDimension(lg,n);
366:     PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(svd->tol)-2,0.0);
367:   }
368:   PetscMalloc2(n,&x,n,&y);
369:   for (i=0;i<n;i++) {
370:     x[i] = (PetscReal)its;
371:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
372:     else y[i] = 0.0;
373:   }
374:   PetscDrawLGAddPoint(lg,x,y);
375:   if (its <= 20 || !(its % 5) || svd->reason) {
376:     PetscDrawLGDraw(lg);
377:     PetscDrawLGSave(lg);
378:   }
379:   PetscFree2(x,y);
380:   return(0);
381: }