Actual source code: svdmon.c

slepc-3.17.2 2022-08-09
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, 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>
 15: #include <petscdraw.h>

 17: PetscErrorCode SVDMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],const char metric[],PetscInt l,const char *names[],int x,int y,int m,int n,PetscDrawLG *lgctx)
 18: {
 19:   PetscDraw      draw;
 20:   PetscDrawAxis  axis;
 21:   PetscDrawLG    lg;

 23:   PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
 24:   PetscDrawSetFromOptions(draw);
 25:   PetscDrawLGCreate(draw,l,&lg);
 26:   if (names) PetscDrawLGSetLegend(lg,names);
 27:   PetscDrawLGSetFromOptions(lg);
 28:   PetscDrawLGGetAxis(lg,&axis);
 29:   PetscDrawAxisSetLabels(axis,"Convergence","Iteration",metric);
 30:   PetscDrawDestroy(&draw);
 31:   *lgctx = lg;
 32:   PetscFunctionReturn(0);
 33: }

 35: /*
 36:    Runs the user provided monitor routines, if any.
 37: */
 38: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest)
 39: {
 40:   PetscInt       i,n = svd->numbermonitors;

 42:   for (i=0;i<n;i++) (*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]);
 43:   PetscFunctionReturn(0);
 44: }

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

 50:    Logically Collective on svd

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

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

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

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

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

 89:    Level: intermediate

 91: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorCancel()
 92: @*/
 93: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
 94: {
 97:   svd->monitor[svd->numbermonitors]           = monitor;
 98:   svd->monitorcontext[svd->numbermonitors]    = (void*)mctx;
 99:   svd->monitordestroy[svd->numbermonitors++]  = monitordestroy;
100:   PetscFunctionReturn(0);
101: }

103: /*@
104:    SVDMonitorCancel - Clears all monitors for an SVD object.

106:    Logically Collective on svd

108:    Input Parameters:
109: .  svd - singular value solver context obtained from SVDCreate()

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

116:    Level: intermediate

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

125:   for (i=0; i<svd->numbermonitors; i++) {
126:     if (svd->monitordestroy[i]) (*svd->monitordestroy[i])(&svd->monitorcontext[i]);
127:   }
128:   svd->numbermonitors = 0;
129:   PetscFunctionReturn(0);
130: }

132: /*@C
133:    SVDGetMonitorContext - Gets the monitor context, as set by
134:    SVDMonitorSet() for the FIRST monitor only.

136:    Not Collective

138:    Input Parameter:
139: .  svd - singular value solver context obtained from SVDCreate()

141:    Output Parameter:
142: .  ctx - monitor context

144:    Level: intermediate

146: .seealso: SVDMonitorSet()
147: @*/
148: PetscErrorCode SVDGetMonitorContext(SVD svd,void *ctx)
149: {
151:   *(void**)ctx = svd->monitorcontext[0];
152:   PetscFunctionReturn(0);
153: }

155: /*@C
156:    SVDMonitorFirst - Print the first unconverged approximate value and
157:    error estimate at each iteration of the singular value solver.

159:    Collective on svd

161:    Input Parameters:
162: +  svd    - singular value solver context
163: .  its    - iteration number
164: .  nconv  - number of converged singular triplets so far
165: .  sigma  - singular values
166: .  errest - error estimates
167: .  nest   - number of error estimates to display
168: -  vf     - viewer and format for monitoring

170:    Options Database Key:
171: .  -svd_monitor - activates SVDMonitorFirst()

173:    Level: intermediate

175: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConverged()
176: @*/
177: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
178: {
179:   PetscViewer    viewer = vf->viewer;

183:   if (its==1 && ((PetscObject)svd)->prefix) PetscViewerASCIIPrintf(viewer,"  Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
184:   if (nconv<nest) {
185:     PetscViewerPushFormat(viewer,vf->format);
186:     PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
187:     PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv);
188:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
189:     PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]);
190:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
191:     PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
192:     PetscViewerPopFormat(viewer);
193:   }
194:   PetscFunctionReturn(0);
195: }

197: /*@C
198:    SVDMonitorAll - Print the current approximate values and
199:    error estimates at each iteration of the singular value solver.

201:    Collective on svd

203:    Input Parameters:
204: +  svd    - singular value solver context
205: .  its    - iteration number
206: .  nconv  - number of converged singular triplets so far
207: .  sigma  - singular values
208: .  errest - error estimates
209: .  nest   - number of error estimates to display
210: -  vf     - viewer and format for monitoring

212:    Options Database Key:
213: .  -svd_monitor_all - activates SVDMonitorAll()

215:    Level: intermediate

217: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConverged()
218: @*/
219: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
220: {
221:   PetscInt       i;
222:   PetscViewer    viewer = vf->viewer;

226:   PetscViewerPushFormat(viewer,vf->format);
227:   PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
228:   if (its==1 && ((PetscObject)svd)->prefix) PetscViewerASCIIPrintf(viewer,"  Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
229:   PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD nconv=%" PetscInt_FMT " Values (Errors)",its,nconv);
230:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
231:   for (i=0;i<nest;i++) PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]);
232:   PetscViewerASCIIPrintf(viewer,"\n");
233:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
234:   PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
235:   PetscViewerPopFormat(viewer);
236:   PetscFunctionReturn(0);
237: }

239: /*@C
240:    SVDMonitorConverged - Print the approximate values and
241:    error estimates as they converge.

243:    Collective on svd

245:    Input Parameters:
246: +  svd    - singular value solver context
247: .  its    - iteration number
248: .  nconv  - number of converged singular triplets so far
249: .  sigma  - singular values
250: .  errest - error estimates
251: .  nest   - number of error estimates to display
252: -  vf     - viewer and format for monitoring

254:    Options Database Key:
255: .  -svd_monitor_conv - activates SVDMonitorConverged()

257:    Level: intermediate

259: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorAll()
260: @*/
261: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
262: {
263:   PetscInt       i;
264:   PetscViewer    viewer = vf->viewer;
265:   SlepcConvMon   ctx;

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

288: PetscErrorCode SVDMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
289: {
290:   SlepcConvMon   mctx;

292:   PetscViewerAndFormatCreate(viewer,format,vf);
293:   PetscNew(&mctx);
294:   mctx->ctx = ctx;
295:   (*vf)->data = (void*)mctx;
296:   PetscFunctionReturn(0);
297: }

299: PetscErrorCode SVDMonitorConvergedDestroy(PetscViewerAndFormat **vf)
300: {
301:   if (!*vf) PetscFunctionReturn(0);
302:   PetscFree((*vf)->data);
303:   PetscViewerDestroy(&(*vf)->viewer);
304:   PetscDrawLGDestroy(&(*vf)->lg);
305:   PetscFree(*vf);
306:   PetscFunctionReturn(0);
307: }

309: /*@C
310:    SVDMonitorFirstDrawLG - Plots the error estimate of the first unconverged
311:    approximation at each iteration of the singular value solver.

313:    Collective on svd

315:    Input Parameters:
316: +  svd    - singular value solver context
317: .  its    - iteration number
318: .  its    - iteration number
319: .  nconv  - number of converged singular triplets so far
320: .  sigma  - singular values
321: .  errest - error estimates
322: .  nest   - number of error estimates to display
323: -  vf     - viewer and format for monitoring

325:    Options Database Key:
326: .  -svd_monitor draw::draw_lg - activates SVDMonitorFirstDrawLG()

328:    Level: intermediate

330: .seealso: SVDMonitorSet()
331: @*/
332: PetscErrorCode SVDMonitorFirstDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
333: {
334:   PetscViewer    viewer = vf->viewer;
335:   PetscDrawLG    lg = vf->lg;
336:   PetscReal      x,y;

341:   PetscViewerPushFormat(viewer,vf->format);
342:   if (its==1) {
343:     PetscDrawLGReset(lg);
344:     PetscDrawLGSetDimension(lg,1);
345:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0);
346:   }
347:   if (nconv<nest) {
348:     x = (PetscReal)its;
349:     if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
350:     else y = 0.0;
351:     PetscDrawLGAddPoint(lg,&x,&y);
352:     if (its <= 20 || !(its % 5) || svd->reason) {
353:       PetscDrawLGDraw(lg);
354:       PetscDrawLGSave(lg);
355:     }
356:   }
357:   PetscViewerPopFormat(viewer);
358:   PetscFunctionReturn(0);
359: }

361: /*@C
362:    SVDMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.

364:    Collective on viewer

366:    Input Parameters:
367: +  viewer - the viewer
368: .  format - the viewer format
369: -  ctx    - an optional user context

371:    Output Parameter:
372: .  vf     - the viewer and format context

374:    Level: intermediate

376: .seealso: SVDMonitorSet()
377: @*/
378: PetscErrorCode SVDMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
379: {
380:   PetscViewerAndFormatCreate(viewer,format,vf);
381:   (*vf)->data = ctx;
382:   SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
383:   PetscFunctionReturn(0);
384: }

386: /*@C
387:    SVDMonitorAllDrawLG - Plots the error estimate of all unconverged
388:    approximations at each iteration of the singular value solver.

390:    Collective on svd

392:    Input Parameters:
393: +  svd    - singular value solver context
394: .  its    - iteration number
395: .  its    - iteration number
396: .  nconv  - number of converged singular triplets so far
397: .  sigma  - singular values
398: .  errest - error estimates
399: .  nest   - number of error estimates to display
400: -  vf     - viewer and format for monitoring

402:    Options Database Key:
403: .  -svd_monitor_all draw::draw_lg - activates SVDMonitorAllDrawLG()

405:    Level: intermediate

407: .seealso: SVDMonitorSet()
408: @*/
409: PetscErrorCode SVDMonitorAllDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
410: {
411:   PetscViewer    viewer = vf->viewer;
412:   PetscDrawLG    lg = vf->lg;
413:   PetscInt       i,n = PetscMin(svd->nsv,255);
414:   PetscReal      *x,*y;

419:   PetscViewerPushFormat(viewer,vf->format);
420:   if (its==1) {
421:     PetscDrawLGReset(lg);
422:     PetscDrawLGSetDimension(lg,n);
423:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0);
424:   }
425:   PetscMalloc2(n,&x,n,&y);
426:   for (i=0;i<n;i++) {
427:     x[i] = (PetscReal)its;
428:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
429:     else y[i] = 0.0;
430:   }
431:   PetscDrawLGAddPoint(lg,x,y);
432:   if (its <= 20 || !(its % 5) || svd->reason) {
433:     PetscDrawLGDraw(lg);
434:     PetscDrawLGSave(lg);
435:   }
436:   PetscFree2(x,y);
437:   PetscViewerPopFormat(viewer);
438:   PetscFunctionReturn(0);
439: }

441: /*@C
442:    SVDMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

444:    Collective on viewer

446:    Input Parameters:
447: +  viewer - the viewer
448: .  format - the viewer format
449: -  ctx    - an optional user context

451:    Output Parameter:
452: .  vf     - the viewer and format context

454:    Level: intermediate

456: .seealso: SVDMonitorSet()
457: @*/
458: PetscErrorCode SVDMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
459: {
460:   PetscViewerAndFormatCreate(viewer,format,vf);
461:   (*vf)->data = ctx;
462:   SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
463:   PetscFunctionReturn(0);
464: }

466: /*@C
467:    SVDMonitorConvergedDrawLG - Plots the number of converged eigenvalues
468:    at each iteration of the singular value solver.

470:    Collective on svd

472:    Input Parameters:
473: +  svd    - singular value solver context
474: .  its    - iteration number
475: .  its    - iteration number
476: .  nconv  - number of converged singular triplets so far
477: .  sigma  - singular values
478: .  errest - error estimates
479: .  nest   - number of error estimates to display
480: -  vf     - viewer and format for monitoring

482:    Options Database Key:
483: .  -svd_monitor_conv draw::draw_lg - activates SVDMonitorConvergedDrawLG()

485:    Level: intermediate

487: .seealso: SVDMonitorSet()
488: @*/
489: PetscErrorCode SVDMonitorConvergedDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
490: {
491:   PetscViewer      viewer = vf->viewer;
492:   PetscDrawLG      lg = vf->lg;
493:   PetscReal        x,y;

498:   PetscViewerPushFormat(viewer,vf->format);
499:   if (its==1) {
500:     PetscDrawLGReset(lg);
501:     PetscDrawLGSetDimension(lg,1);
502:     PetscDrawLGSetLimits(lg,1,1,0,svd->nsv);
503:   }
504:   x = (PetscReal)its;
505:   y = (PetscReal)svd->nconv;
506:   PetscDrawLGAddPoint(lg,&x,&y);
507:   if (its <= 20 || !(its % 5) || svd->reason) {
508:     PetscDrawLGDraw(lg);
509:     PetscDrawLGSave(lg);
510:   }
511:   PetscViewerPopFormat(viewer);
512:   PetscFunctionReturn(0);
513: }

515: /*@C
516:    SVDMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

518:    Collective on viewer

520:    Input Parameters:
521: +  viewer - the viewer
522: .  format - the viewer format
523: -  ctx    - an optional user context

525:    Output Parameter:
526: .  vf     - the viewer and format context

528:    Level: intermediate

530: .seealso: SVDMonitorSet()
531: @*/
532: PetscErrorCode SVDMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
533: {
534:   SlepcConvMon   mctx;

536:   PetscViewerAndFormatCreate(viewer,format,vf);
537:   PetscNew(&mctx);
538:   mctx->ctx = ctx;
539:   (*vf)->data = (void*)mctx;
540:   SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
541:   PetscFunctionReturn(0);
542: }