Actual source code: svdmon.c

slepc-3.15.2 2021-09-20
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, 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;

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

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

 46:   for (i=0;i<n;i++) {
 47:     (*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]);
 48:   }
 49:   return(0);
 50: }

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

 56:    Logically Collective on svd

 58:    Input Parameters:
 59: +  svd     - singular value solver context obtained from SVDCreate()
 60: .  monitor - pointer to function (if this is NULL, it turns off monitoring)
 61: .  mctx    - [optional] context for private data for the
 62:              monitor routine (use NULL if no context is desired)
 63: -  monitordestroy - [optional] routine that frees monitor context (may be NULL)

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

 68: +  svd    - singular value solver context obtained from SVDCreate()
 69: .  its    - iteration number
 70: .  nconv  - number of converged singular triplets
 71: .  sigma  - singular values
 72: .  errest - relative error estimates for each singular triplet
 73: .  nest   - number of error estimates
 74: -  mctx   - optional monitoring context, as set by SVDMonitorSet()

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

 90:    Notes:
 91:    Several different monitoring routines may be set by calling
 92:    SVDMonitorSet() multiple times; all will be called in the
 93:    order in which they were set.

 95:    Level: intermediate

 97: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorCancel()
 98: @*/
 99: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
100: {
103:   if (svd->numbermonitors >= MAXSVDMONITORS) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
104:   svd->monitor[svd->numbermonitors]           = monitor;
105:   svd->monitorcontext[svd->numbermonitors]    = (void*)mctx;
106:   svd->monitordestroy[svd->numbermonitors++]  = monitordestroy;
107:   return(0);
108: }

110: /*@
111:    SVDMonitorCancel - Clears all monitors for an SVD object.

113:    Logically Collective on svd

115:    Input Parameters:
116: .  svd - singular value solver context obtained from SVDCreate()

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

123:    Level: intermediate

125: .seealso: SVDMonitorSet()
126: @*/
127: PetscErrorCode SVDMonitorCancel(SVD svd)
128: {
130:   PetscInt       i;

134:   for (i=0; i<svd->numbermonitors; i++) {
135:     if (svd->monitordestroy[i]) {
136:       (*svd->monitordestroy[i])(&svd->monitorcontext[i]);
137:     }
138:   }
139:   svd->numbermonitors = 0;
140:   return(0);
141: }

143: /*@C
144:    SVDGetMonitorContext - Gets the monitor context, as set by
145:    SVDMonitorSet() for the FIRST monitor only.

147:    Not Collective

149:    Input Parameter:
150: .  svd - singular value solver context obtained from SVDCreate()

152:    Output Parameter:
153: .  ctx - monitor context

155:    Level: intermediate

157: .seealso: SVDMonitorSet()
158: @*/
159: PetscErrorCode SVDGetMonitorContext(SVD svd,void **ctx)
160: {
163:   *ctx = svd->monitorcontext[0];
164:   return(0);
165: }

167: /*@C
168:    SVDMonitorFirst - Print the first unconverged approximate value and
169:    error estimate at each iteration of the singular value solver.

171:    Collective on svd

173:    Input Parameters:
174: +  svd    - singular value solver context
175: .  its    - iteration number
176: .  nconv  - number of converged singular triplets so far
177: .  sigma  - singular values
178: .  errest - error estimates
179: .  nest   - number of error estimates to display
180: -  vf     - viewer and format for monitoring

182:    Options Database Key:
183: .  -svd_monitor - activates SVDMonitorFirst()

185:    Level: intermediate

187: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConverged()
188: @*/
189: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
190: {
192:   PetscViewer    viewer;

197:   viewer = vf->viewer;
199:   if (its==1 && ((PetscObject)svd)->prefix) {
200:     PetscViewerASCIIPrintf(viewer,"  Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
201:   }
202:   if (nconv<nest) {
203:     PetscViewerPushFormat(viewer,vf->format);
204:     PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
205:     PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D first unconverged value (error)",its,nconv);
206:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
207:     PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]);
208:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
209:     PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
210:     PetscViewerPopFormat(viewer);
211:   }
212:   return(0);
213: }

215: /*@C
216:    SVDMonitorAll - Print the current approximate values and
217:    error estimates at each iteration of the singular value solver.

219:    Collective on svd

221:    Input Parameters:
222: +  svd    - singular value solver context
223: .  its    - iteration number
224: .  nconv  - number of converged singular triplets so far
225: .  sigma  - singular values
226: .  errest - error estimates
227: .  nest   - number of error estimates to display
228: -  vf     - viewer and format for monitoring

230:    Options Database Key:
231: .  -svd_monitor_all - activates SVDMonitorAll()

233:    Level: intermediate

235: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConverged()
236: @*/
237: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
238: {
240:   PetscInt       i;
241:   PetscViewer    viewer;

246:   viewer = vf->viewer;
248:   PetscViewerPushFormat(viewer,vf->format);
249:   PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
250:   if (its==1 && ((PetscObject)svd)->prefix) {
251:     PetscViewerASCIIPrintf(viewer,"  Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
252:   }
253:   PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D Values (Errors)",its,nconv);
254:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
255:   for (i=0;i<nest;i++) {
256:     PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]);
257:   }
258:   PetscViewerASCIIPrintf(viewer,"\n");
259:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
260:   PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
261:   PetscViewerPopFormat(viewer);
262:   return(0);
263: }

265: /*@C
266:    SVDMonitorConverged - Print the approximate values and
267:    error estimates as they converge.

269:    Collective on svd

271:    Input Parameters:
272: +  svd    - singular value solver context
273: .  its    - iteration number
274: .  nconv  - number of converged singular triplets so far
275: .  sigma  - singular values
276: .  errest - error estimates
277: .  nest   - number of error estimates to display
278: -  vf     - viewer and format for monitoring

280:    Options Database Key:
281: .  -svd_monitor_conv - activates SVDMonitorConverged()

283:    Level: intermediate

285: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorAll()
286: @*/
287: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
288: {
290:   PetscInt       i;
291:   PetscViewer    viewer;
292:   SlepcConvMon   ctx;

297:   viewer = vf->viewer;
299:   ctx = (SlepcConvMon)vf->data;
300:   if (its==1 && ((PetscObject)svd)->prefix) {
301:     PetscViewerASCIIPrintf(viewer,"  Convergence history for %s solve.\n",((PetscObject)svd)->prefix);
302:   }
303:   if (its==1) ctx->oldnconv = 0;
304:   if (ctx->oldnconv!=nconv) {
305:     PetscViewerPushFormat(viewer,vf->format);
306:     PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
307:     for (i=ctx->oldnconv;i<nconv;i++) {
308:       PetscViewerASCIIPrintf(viewer,"%3D SVD converged value (error) #%D",its,i);
309:       PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
310:       PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]);
311:       PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
312:     }
313:     PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
314:     PetscViewerPopFormat(viewer);
315:     ctx->oldnconv = nconv;
316:   }
317:   return(0);
318: }

320: PetscErrorCode SVDMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
321: {
323:   SlepcConvMon   mctx;

326:   PetscViewerAndFormatCreate(viewer,format,vf);
327:   PetscNew(&mctx);
328:   mctx->ctx = ctx;
329:   (*vf)->data = (void*)mctx;
330:   return(0);
331: }

333: PetscErrorCode SVDMonitorConvergedDestroy(PetscViewerAndFormat **vf)
334: {

338:   if (!*vf) return(0);
339:   PetscFree((*vf)->data);
340:   PetscViewerDestroy(&(*vf)->viewer);
341:   PetscDrawLGDestroy(&(*vf)->lg);
342:   PetscFree(*vf);
343:   return(0);
344: }

346: /*@C
347:    SVDMonitorFirstDrawLG - Plots the error estimate of the first unconverged
348:    approximation at each iteration of the singular value solver.

350:    Collective on svd

352:    Input Parameters:
353: +  svd    - singular value solver context
354: .  its    - iteration number
355: .  its    - iteration number
356: .  nconv  - number of converged singular triplets so far
357: .  sigma  - singular values
358: .  errest - error estimates
359: .  nest   - number of error estimates to display
360: -  vf     - viewer and format for monitoring

362:    Options Database Key:
363: .  -svd_monitor draw::draw_lg - activates SVDMonitorFirstDrawLG()

365:    Level: intermediate

367: .seealso: SVDMonitorSet()
368: @*/
369: PetscErrorCode SVDMonitorFirstDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
370: {
372:   PetscViewer    viewer;
373:   PetscDrawLG    lg;
374:   PetscReal      x,y;

379:   viewer = vf->viewer;
381:   lg = vf->lg;
383:   PetscViewerPushFormat(viewer,vf->format);
384:   if (its==1) {
385:     PetscDrawLGReset(lg);
386:     PetscDrawLGSetDimension(lg,1);
387:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0);
388:   }
389:   if (nconv<nest) {
390:     x = (PetscReal)its;
391:     if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
392:     else y = 0.0;
393:     PetscDrawLGAddPoint(lg,&x,&y);
394:     if (its <= 20 || !(its % 5) || svd->reason) {
395:       PetscDrawLGDraw(lg);
396:       PetscDrawLGSave(lg);
397:     }
398:   }
399:   PetscViewerPopFormat(viewer);
400:   return(0);
401: }

403: /*@C
404:    SVDMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.

406:    Collective on viewer

408:    Input Parameters:
409: +  viewer - the viewer
410: .  format - the viewer format
411: -  ctx    - an optional user context

413:    Output Parameter:
414: .  vf     - the viewer and format context

416:    Level: intermediate

418: .seealso: SVDMonitorSet()
419: @*/
420: PetscErrorCode SVDMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
421: {

425:   PetscViewerAndFormatCreate(viewer,format,vf);
426:   (*vf)->data = ctx;
427:   SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
428:   return(0);
429: }

431: /*@C
432:    SVDMonitorAllDrawLG - Plots the error estimate of all unconverged
433:    approximations at each iteration of the singular value solver.

435:    Collective on svd

437:    Input Parameters:
438: +  svd    - singular value solver context
439: .  its    - iteration number
440: .  its    - iteration number
441: .  nconv  - number of converged singular triplets so far
442: .  sigma  - singular values
443: .  errest - error estimates
444: .  nest   - number of error estimates to display
445: -  vf     - viewer and format for monitoring

447:    Options Database Key:
448: .  -svd_monitor_all draw::draw_lg - activates SVDMonitorAllDrawLG()

450:    Level: intermediate

452: .seealso: SVDMonitorSet()
453: @*/
454: PetscErrorCode SVDMonitorAllDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
455: {
457:   PetscViewer    viewer;
458:   PetscDrawLG    lg;
459:   PetscInt       i,n = PetscMin(svd->nsv,255);
460:   PetscReal      *x,*y;

465:   viewer = vf->viewer;
467:   lg = vf->lg;
469:   PetscViewerPushFormat(viewer,vf->format);
470:   if (its==1) {
471:     PetscDrawLGReset(lg);
472:     PetscDrawLGSetDimension(lg,n);
473:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(svd->tol)-2,0.0);
474:   }
475:   PetscMalloc2(n,&x,n,&y);
476:   for (i=0;i<n;i++) {
477:     x[i] = (PetscReal)its;
478:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
479:     else y[i] = 0.0;
480:   }
481:   PetscDrawLGAddPoint(lg,x,y);
482:   if (its <= 20 || !(its % 5) || svd->reason) {
483:     PetscDrawLGDraw(lg);
484:     PetscDrawLGSave(lg);
485:   }
486:   PetscFree2(x,y);
487:   PetscViewerPopFormat(viewer);
488:   return(0);
489: }

491: /*@C
492:    SVDMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

494:    Collective on viewer

496:    Input Parameters:
497: +  viewer - the viewer
498: .  format - the viewer format
499: -  ctx    - an optional user context

501:    Output Parameter:
502: .  vf     - the viewer and format context

504:    Level: intermediate

506: .seealso: SVDMonitorSet()
507: @*/
508: PetscErrorCode SVDMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
509: {

513:   PetscViewerAndFormatCreate(viewer,format,vf);
514:   (*vf)->data = ctx;
515:   SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
516:   return(0);
517: }

519: /*@C
520:    SVDMonitorConvergedDrawLG - Plots the number of converged eigenvalues
521:    at each iteration of the singular value solver.

523:    Collective on svd

525:    Input Parameters:
526: +  svd    - singular value solver context
527: .  its    - iteration number
528: .  its    - iteration number
529: .  nconv  - number of converged singular triplets so far
530: .  sigma  - singular values
531: .  errest - error estimates
532: .  nest   - number of error estimates to display
533: -  vf     - viewer and format for monitoring

535:    Options Database Key:
536: .  -svd_monitor_conv draw::draw_lg - activates SVDMonitorConvergedDrawLG()

538:    Level: intermediate

540: .seealso: SVDMonitorSet()
541: @*/
542: PetscErrorCode SVDMonitorConvergedDrawLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
543: {
544:   PetscErrorCode   ierr;
545:   PetscViewer      viewer;
546:   PetscDrawLG      lg;
547:   PetscReal        x,y;

552:   viewer = vf->viewer;
554:   lg = vf->lg;
556:   PetscViewerPushFormat(viewer,vf->format);
557:   if (its==1) {
558:     PetscDrawLGReset(lg);
559:     PetscDrawLGSetDimension(lg,1);
560:     PetscDrawLGSetLimits(lg,1,1,0,svd->nsv);
561:   }
562:   x = (PetscReal)its;
563:   y = (PetscReal)svd->nconv;
564:   PetscDrawLGAddPoint(lg,&x,&y);
565:   if (its <= 20 || !(its % 5) || svd->reason) {
566:     PetscDrawLGDraw(lg);
567:     PetscDrawLGSave(lg);
568:   }
569:   PetscViewerPopFormat(viewer);
570:   return(0);
571: }

573: /*@C
574:    SVDMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

576:    Collective on viewer

578:    Input Parameters:
579: +  viewer - the viewer
580: .  format - the viewer format
581: -  ctx    - an optional user context

583:    Output Parameter:
584: .  vf     - the viewer and format context

586:    Level: intermediate

588: .seealso: SVDMonitorSet()
589: @*/
590: PetscErrorCode SVDMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
591: {
593:   SlepcConvMon   mctx;

596:   PetscViewerAndFormatCreate(viewer,format,vf);
597:   PetscNew(&mctx);
598:   mctx->ctx = ctx;
599:   (*vf)->data = (void*)mctx;
600:   SVDMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
601:   return(0);
602: }