Actual source code: nepmon.c

slepc-3.18.3 2023-03-24
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:    NEP routines related to monitors
 12: */

 14: #include <slepc/private/nepimpl.h>
 15: #include <petscdraw.h>

 17: PetscErrorCode NEPMonitorLGCreate(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 NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
 39: {
 40:   PetscInt       i,n = nep->numbermonitors;

 42:   for (i=0;i<n;i++) (*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]);
 43:   return 0;
 44: }

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

 50:    Logically Collective on nep

 52:    Input Parameters:
 53: +  nep     - eigensolver context obtained from NEPCreate()
 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(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,PetscInt nest,void *mctx)

 62: +  nep    - nonlinear eigensolver context obtained from NEPCreate()
 63: .  its    - iteration number
 64: .  nconv  - number of converged eigenpairs
 65: .  eigr   - real part of the eigenvalues
 66: .  eigi   - imaginary part of the eigenvalues
 67: .  errest - error estimates for each eigenpair
 68: .  nest   - number of error estimates
 69: -  mctx   - optional monitoring context, as set by NEPMonitorSet()

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

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

 90:    Level: intermediate

 92: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
 93: @*/
 94: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
 95: {
 98:   nep->monitor[nep->numbermonitors]           = monitor;
 99:   nep->monitorcontext[nep->numbermonitors]    = (void*)mctx;
100:   nep->monitordestroy[nep->numbermonitors++]  = monitordestroy;
101:   return 0;
102: }

104: /*@
105:    NEPMonitorCancel - Clears all monitors for a NEP object.

107:    Logically Collective on nep

109:    Input Parameters:
110: .  nep - eigensolver context obtained from NEPCreate()

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

117:    Level: intermediate

119: .seealso: NEPMonitorSet()
120: @*/
121: PetscErrorCode NEPMonitorCancel(NEP nep)
122: {
123:   PetscInt       i;

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

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

137:    Not Collective

139:    Input Parameter:
140: .  nep - eigensolver context obtained from NEPCreate()

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

145:    Level: intermediate

147: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
148: @*/
149: PetscErrorCode NEPGetMonitorContext(NEP nep,void *ctx)
150: {
152:   *(void**)ctx = nep->monitorcontext[0];
153:   return 0;
154: }

156: /*@C
157:    NEPMonitorFirst - Print the first unconverged approximate value and
158:    error estimate at each iteration of the nonlinear eigensolver.

160:    Collective on nep

162:    Input Parameters:
163: +  nep    - nonlinear eigensolver context
164: .  its    - iteration number
165: .  nconv  - number of converged eigenpairs so far
166: .  eigr   - real part of the eigenvalues
167: .  eigi   - imaginary part of the eigenvalues
168: .  errest - error estimates
169: .  nest   - number of error estimates to display
170: -  vf     - viewer and format for monitoring

172:    Options Database Key:
173: .  -nep_monitor - activates NEPMonitorFirst()

175:    Level: intermediate

177: .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
178: @*/
179: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
180: {
181:   PetscViewer    viewer = vf->viewer;

185:   if (its==1 && ((PetscObject)nep)->prefix) PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
186:   if (nconv<nest) {
187:     PetscViewerPushFormat(viewer,vf->format);
188:     PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
189:     PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv);
190:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
191: #if defined(PETSC_USE_COMPLEX)
192:     PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv]));
193: #else
194:     PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]);
195:     if (eigi[nconv]!=0.0) PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]);
196: #endif
197:     PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
198:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
199:     PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
200:     PetscViewerPopFormat(viewer);
201:   }
202:   return 0;
203: }

205: /*@C
206:    NEPMonitorAll - Print the current approximate values and
207:    error estimates at each iteration of the nonlinear eigensolver.

209:    Collective on nep

211:    Input Parameters:
212: +  nep    - nonlinear eigensolver context
213: .  its    - iteration number
214: .  nconv  - number of converged eigenpairs so far
215: .  eigr   - real part of the eigenvalues
216: .  eigi   - imaginary part of the eigenvalues
217: .  errest - error estimates
218: .  nest   - number of error estimates to display
219: -  vf     - viewer and format for monitoring

221:    Options Database Key:
222: .  -nep_monitor_all - activates NEPMonitorAll()

224:    Level: intermediate

226: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
227: @*/
228: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
229: {
230:   PetscInt       i;
231:   PetscViewer    viewer = vf->viewer;

235:   PetscViewerPushFormat(viewer,vf->format);
236:   PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
237:   if (its==1 && ((PetscObject)nep)->prefix) PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
238:   PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv);
239:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
240:   for (i=0;i<nest;i++) {
241: #if defined(PETSC_USE_COMPLEX)
242:     PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
243: #else
244:     PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
245:     if (eigi[i]!=0.0) PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]);
246: #endif
247:     PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
248:   }
249:   PetscViewerASCIIPrintf(viewer,"\n");
250:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
251:   PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
252:   PetscViewerPopFormat(viewer);
253:   return 0;
254: }

256: /*@C
257:    NEPMonitorConverged - Print the approximate values and
258:    error estimates as they converge.

260:    Collective on nep

262:    Input Parameters:
263: +  nep    - nonlinear eigensolver context
264: .  its    - iteration number
265: .  nconv  - number of converged eigenpairs so far
266: .  eigr   - real part of the eigenvalues
267: .  eigi   - imaginary part of the eigenvalues
268: .  errest - error estimates
269: .  nest   - number of error estimates to display
270: -  vf     - viewer and format for monitoring

272:    Options Database Key:
273: .  -nep_monitor_conv - activates NEPMonitorConverged()

275:    Level: intermediate

277: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
278: @*/
279: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
280: {
281:   PetscInt       i;
282:   PetscViewer    viewer = vf->viewer;
283:   SlepcConvMon   ctx;

287:   ctx = (SlepcConvMon)vf->data;
288:   if (its==1 && ((PetscObject)nep)->prefix) PetscViewerASCIIPrintf(viewer,"  Convergence history for %s solve.\n",((PetscObject)nep)->prefix);
289:   if (its==1) ctx->oldnconv = 0;
290:   if (ctx->oldnconv!=nconv) {
291:     PetscViewerPushFormat(viewer,vf->format);
292:     PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
293:     for (i=ctx->oldnconv;i<nconv;i++) {
294:       PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP converged value (error) #%" PetscInt_FMT,its,i);
295:       PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
296: #if defined(PETSC_USE_COMPLEX)
297:       PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
298: #else
299:       PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
300:       if (eigi[i]!=0.0) PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]);
301: #endif
302:       PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
303:       PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
304:     }
305:     PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
306:     PetscViewerPopFormat(viewer);
307:     ctx->oldnconv = nconv;
308:   }
309:   return 0;
310: }

312: PetscErrorCode NEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
313: {
314:   SlepcConvMon   mctx;

316:   PetscViewerAndFormatCreate(viewer,format,vf);
317:   PetscNew(&mctx);
318:   mctx->ctx = ctx;
319:   (*vf)->data = (void*)mctx;
320:   return 0;
321: }

323: PetscErrorCode NEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
324: {
325:   if (!*vf) return 0;
326:   PetscFree((*vf)->data);
327:   PetscViewerDestroy(&(*vf)->viewer);
328:   PetscDrawLGDestroy(&(*vf)->lg);
329:   PetscFree(*vf);
330:   return 0;
331: }

333: /*@C
334:    NEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
335:    approximation at each iteration of the nonlinear eigensolver.

337:    Collective on nep

339:    Input Parameters:
340: +  nep    - nonlinear eigensolver context
341: .  its    - iteration number
342: .  its    - iteration number
343: .  nconv  - number of converged eigenpairs so far
344: .  eigr   - real part of the eigenvalues
345: .  eigi   - imaginary part of the eigenvalues
346: .  errest - error estimates
347: .  nest   - number of error estimates to display
348: -  vf     - viewer and format for monitoring

350:    Options Database Key:
351: .  -nep_monitor draw::draw_lg - activates NEPMonitorFirstDrawLG()

353:    Level: intermediate

355: .seealso: NEPMonitorSet()
356: @*/
357: PetscErrorCode NEPMonitorFirstDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
358: {
359:   PetscViewer    viewer = vf->viewer;
360:   PetscDrawLG    lg = vf->lg;
361:   PetscReal      x,y;

366:   PetscViewerPushFormat(viewer,vf->format);
367:   if (its==1) {
368:     PetscDrawLGReset(lg);
369:     PetscDrawLGSetDimension(lg,1);
370:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0);
371:   }
372:   if (nconv<nest) {
373:     x = (PetscReal)its;
374:     if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
375:     else y = 0.0;
376:     PetscDrawLGAddPoint(lg,&x,&y);
377:     if (its <= 20 || !(its % 5) || nep->reason) {
378:       PetscDrawLGDraw(lg);
379:       PetscDrawLGSave(lg);
380:     }
381:   }
382:   PetscViewerPopFormat(viewer);
383:   return 0;
384: }

386: /*@C
387:    NEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.

389:    Collective on viewer

391:    Input Parameters:
392: +  viewer - the viewer
393: .  format - the viewer format
394: -  ctx    - an optional user context

396:    Output Parameter:
397: .  vf     - the viewer and format context

399:    Level: intermediate

401: .seealso: NEPMonitorSet()
402: @*/
403: PetscErrorCode NEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
404: {
405:   PetscViewerAndFormatCreate(viewer,format,vf);
406:   (*vf)->data = ctx;
407:   NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
408:   return 0;
409: }

411: /*@C
412:    NEPMonitorAllDrawLG - Plots the error estimate of all unconverged
413:    approximations at each iteration of the nonlinear eigensolver.

415:    Collective on nep

417:    Input Parameters:
418: +  nep    - nonlinear eigensolver context
419: .  its    - iteration number
420: .  its    - iteration number
421: .  nconv  - number of converged eigenpairs so far
422: .  eigr   - real part of the eigenvalues
423: .  eigi   - imaginary part of the eigenvalues
424: .  errest - error estimates
425: .  nest   - number of error estimates to display
426: -  vf     - viewer and format for monitoring

428:    Options Database Key:
429: .  -nep_monitor_all draw::draw_lg - activates NEPMonitorAllDrawLG()

431:    Level: intermediate

433: .seealso: NEPMonitorSet()
434: @*/
435: PetscErrorCode NEPMonitorAllDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
436: {
437:   PetscViewer    viewer = vf->viewer;
438:   PetscDrawLG    lg = vf->lg;
439:   PetscInt       i,n = PetscMin(nep->nev,255);
440:   PetscReal      *x,*y;

445:   PetscViewerPushFormat(viewer,vf->format);
446:   if (its==1) {
447:     PetscDrawLGReset(lg);
448:     PetscDrawLGSetDimension(lg,n);
449:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0);
450:   }
451:   PetscMalloc2(n,&x,n,&y);
452:   for (i=0;i<n;i++) {
453:     x[i] = (PetscReal)its;
454:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
455:     else y[i] = 0.0;
456:   }
457:   PetscDrawLGAddPoint(lg,x,y);
458:   if (its <= 20 || !(its % 5) || nep->reason) {
459:     PetscDrawLGDraw(lg);
460:     PetscDrawLGSave(lg);
461:   }
462:   PetscFree2(x,y);
463:   PetscViewerPopFormat(viewer);
464:   return 0;
465: }

467: /*@C
468:    NEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

470:    Collective on viewer

472:    Input Parameters:
473: +  viewer - the viewer
474: .  format - the viewer format
475: -  ctx    - an optional user context

477:    Output Parameter:
478: .  vf     - the viewer and format context

480:    Level: intermediate

482: .seealso: NEPMonitorSet()
483: @*/
484: PetscErrorCode NEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
485: {
486:   PetscViewerAndFormatCreate(viewer,format,vf);
487:   (*vf)->data = ctx;
488:   NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
489:   return 0;
490: }

492: /*@C
493:    NEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
494:    at each iteration of the nonlinear eigensolver.

496:    Collective on nep

498:    Input Parameters:
499: +  nep    - nonlinear eigensolver context
500: .  its    - iteration number
501: .  its    - iteration number
502: .  nconv  - number of converged eigenpairs so far
503: .  eigr   - real part of the eigenvalues
504: .  eigi   - imaginary part of the eigenvalues
505: .  errest - error estimates
506: .  nest   - number of error estimates to display
507: -  vf     - viewer and format for monitoring

509:    Options Database Key:
510: .  -nep_monitor_conv draw::draw_lg - activates NEPMonitorConvergedDrawLG()

512:    Level: intermediate

514: .seealso: NEPMonitorSet()
515: @*/
516: PetscErrorCode NEPMonitorConvergedDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
517: {
518:   PetscViewer      viewer = vf->viewer;
519:   PetscDrawLG      lg = vf->lg;
520:   PetscReal        x,y;

525:   PetscViewerPushFormat(viewer,vf->format);
526:   if (its==1) {
527:     PetscDrawLGReset(lg);
528:     PetscDrawLGSetDimension(lg,1);
529:     PetscDrawLGSetLimits(lg,1,1,0,nep->nev);
530:   }
531:   x = (PetscReal)its;
532:   y = (PetscReal)nep->nconv;
533:   PetscDrawLGAddPoint(lg,&x,&y);
534:   if (its <= 20 || !(its % 5) || nep->reason) {
535:     PetscDrawLGDraw(lg);
536:     PetscDrawLGSave(lg);
537:   }
538:   PetscViewerPopFormat(viewer);
539:   return 0;
540: }

542: /*@C
543:    NEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

545:    Collective on viewer

547:    Input Parameters:
548: +  viewer - the viewer
549: .  format - the viewer format
550: -  ctx    - an optional user context

552:    Output Parameter:
553: .  vf     - the viewer and format context

555:    Level: intermediate

557: .seealso: NEPMonitorSet()
558: @*/
559: PetscErrorCode NEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
560: {
561:   SlepcConvMon   mctx;

563:   PetscViewerAndFormatCreate(viewer,format,vf);
564:   PetscNew(&mctx);
565:   mctx->ctx = ctx;
566:   (*vf)->data = (void*)mctx;
567:   NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
568:   return 0;
569: }