Actual source code: svdmon.c

slepc-3.18.0 2022-10-01
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:   return 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:   return 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_conditioning - print the condition number when available
 76: .    -svd_monitor draw::draw_lg - sets line graph monitor for the first unconverged
 77:       approximate singular value
 78: .    -svd_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
 79:       approximate singular values
 80: .    -svd_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
 81: -    -svd_monitor_cancel - cancels all monitors that have been hardwired into
 82:       a code by calls to SVDMonitorSet(), but does not cancel those set via
 83:       the options database.

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

 90:    Level: intermediate

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

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

107:    Logically Collective on svd

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

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

117:    Level: intermediate

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

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

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

137:    Not Collective

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

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

145:    Level: intermediate

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

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

160:    Collective on svd

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

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

174:    Level: intermediate

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

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

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

202:    Collective on svd

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

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

216:    Level: intermediate

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

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

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

244:    Collective on svd

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

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

258:    Level: intermediate

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

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

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

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

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

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

314:    Collective on svd

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

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

329:    Level: intermediate

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

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

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

365:    Collective on viewer

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

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

375:    Level: intermediate

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

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

391:    Collective on svd

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

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

406:    Level: intermediate

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

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

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

445:    Collective on viewer

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

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

455:    Level: intermediate

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

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

471:    Collective on svd

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

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

486:    Level: intermediate

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

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

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

519:    Collective on viewer

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

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

529:    Level: intermediate

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

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

545: /*@C
546:    SVDMonitorConditioning - Print the condition number at each iteration of the singular
547:    value solver.

549:    Collective on svd

551:    Input Parameters:
552: +  svd    - singular value solver context
553: .  its    - iteration number
554: .  nconv  - (unused) number of converged singular triplets so far
555: .  sigma  - (unused) singular values
556: .  errest - (unused) error estimates
557: .  nest   - (unused) number of error estimates to display
558: -  vf     - viewer and format for monitoring

560:    Options Database Key:
561: .  -svd_monitor_conditioning - activates SVDMonitorConditioning()

563:    Note:
564:    Works only for solvers that use a DS of type GSVD. The printed information corresponds
565:    to the maximum of the condition number of the two generated bidiagonal matrices.

567:    Level: intermediate

569: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorFirst(), SVDMonitorConverged()
570: @*/
571: PetscErrorCode SVDMonitorConditioning(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
572: {
573:   PetscViewer viewer = vf->viewer;
574:   PetscBool   isgsvd;
575:   PetscReal   cond;

579:   PetscObjectTypeCompare((PetscObject)svd->ds,DSGSVD,&isgsvd);
580:   if (!isgsvd) return 0;
581:   if (its==1 && ((PetscObject)svd)->prefix) PetscViewerASCIIPrintf(viewer,"  Condition number of bidiagonal matrices for %s solve.\n",((PetscObject)svd)->prefix);
582:   PetscViewerPushFormat(viewer,vf->format);
583:   PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
584:   DSCond(svd->ds,&cond);
585:   PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " SVD condition number = %g\n",its,(double)cond);
586:   PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
587:   PetscViewerPopFormat(viewer);
588:   return 0;
589: }