Actual source code: nepmon.c

slepc-3.22.1 2024-10-28
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:   PetscFunctionBegin;
 24:   PetscCall(PetscDrawCreate(comm,host,label,x,y,m,n,&draw));
 25:   PetscCall(PetscDrawSetFromOptions(draw));
 26:   PetscCall(PetscDrawLGCreate(draw,l,&lg));
 27:   if (names) PetscCall(PetscDrawLGSetLegend(lg,names));
 28:   PetscCall(PetscDrawLGSetFromOptions(lg));
 29:   PetscCall(PetscDrawLGGetAxis(lg,&axis));
 30:   PetscCall(PetscDrawAxisSetLabels(axis,"Convergence","Iteration",metric));
 31:   PetscCall(PetscDrawDestroy(&draw));
 32:   *lgctx = lg;
 33:   PetscFunctionReturn(PETSC_SUCCESS);
 34: }

 36: /*
 37:    Runs the user provided monitor routines, if any.
 38: */
 39: PetscErrorCode NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
 40: {
 41:   PetscInt       i,n = nep->numbermonitors;

 43:   PetscFunctionBegin;
 44:   for (i=0;i<n;i++) PetscCall((*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]));
 45:   PetscFunctionReturn(PETSC_SUCCESS);
 46: }

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

 52:    Logically Collective

 54:    Input Parameters:
 55: +  nep     - eigensolver context obtained from NEPCreate()
 56: .  monitor - pointer to function (if this is NULL, it turns off monitoring)
 57: .  mctx    - [optional] context for private data for the
 58:              monitor routine (use NULL if no context is desired)
 59: -  monitordestroy - [optional] routine that frees monitor context (may be NULL)

 61:    Calling sequence of monitor:
 62: $  PetscErrorCode monitor(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx)
 63: +  nep    - nonlinear eigensolver context obtained from NEPCreate()
 64: .  its    - iteration number
 65: .  nconv  - number of converged eigenpairs
 66: .  eigr   - real part of the eigenvalues
 67: .  eigi   - imaginary part of the eigenvalues
 68: .  errest - error estimates for each eigenpair
 69: .  nest   - number of error estimates
 70: -  mctx   - optional monitoring context, as set by NEPMonitorSet()

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

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

 91:    Level: intermediate

 93: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
 94: @*/
 95: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx),void *mctx,PetscErrorCode (*monitordestroy)(void**))
 96: {
 97:   PetscFunctionBegin;
 99:   PetscCheck(nep->numbermonitors<MAXNEPMONITORS,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
100:   nep->monitor[nep->numbermonitors]           = monitor;
101:   nep->monitorcontext[nep->numbermonitors]    = (void*)mctx;
102:   nep->monitordestroy[nep->numbermonitors++]  = monitordestroy;
103:   PetscFunctionReturn(PETSC_SUCCESS);
104: }

106: /*@
107:    NEPMonitorCancel - Clears all monitors for a NEP object.

109:    Logically Collective

111:    Input Parameters:
112: .  nep - eigensolver context obtained from NEPCreate()

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

119:    Level: intermediate

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

127:   PetscFunctionBegin;
129:   for (i=0; i<nep->numbermonitors; i++) {
130:     if (nep->monitordestroy[i]) PetscCall((*nep->monitordestroy[i])(&nep->monitorcontext[i]));
131:   }
132:   nep->numbermonitors = 0;
133:   PetscFunctionReturn(PETSC_SUCCESS);
134: }

136: /*@C
137:    NEPGetMonitorContext - Gets the monitor context, as set by
138:    NEPMonitorSet() for the FIRST monitor only.

140:    Not Collective

142:    Input Parameter:
143: .  nep - eigensolver context obtained from NEPCreate()

145:    Output Parameter:
146: .  ctx - monitor context

148:    Level: intermediate

150: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
151: @*/
152: PetscErrorCode NEPGetMonitorContext(NEP nep,void *ctx)
153: {
154:   PetscFunctionBegin;
156:   *(void**)ctx = nep->monitorcontext[0];
157:   PetscFunctionReturn(PETSC_SUCCESS);
158: }

160: /*@C
161:    NEPMonitorFirst - Print the first unconverged approximate value and
162:    error estimate at each iteration of the nonlinear eigensolver.

164:    Collective

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

176:    Options Database Key:
177: .  -nep_monitor - activates NEPMonitorFirst()

179:    Level: intermediate

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

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

210: /*@C
211:    NEPMonitorAll - Print the current approximate values and
212:    error estimates at each iteration of the nonlinear eigensolver.

214:    Collective

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

226:    Options Database Key:
227: .  -nep_monitor_all - activates NEPMonitorAll()

229:    Level: intermediate

231: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
232: @*/
233: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
234: {
235:   PetscInt       i;
236:   PetscViewer    viewer = vf->viewer;

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

262: /*@C
263:    NEPMonitorConverged - Print the approximate values and
264:    error estimates as they converge.

266:    Collective

268:    Input Parameters:
269: +  nep    - nonlinear eigensolver context
270: .  its    - iteration number
271: .  nconv  - number of converged eigenpairs so far
272: .  eigr   - real part of the eigenvalues
273: .  eigi   - imaginary part of the eigenvalues
274: .  errest - error estimates
275: .  nest   - number of error estimates to display
276: -  vf     - viewer and format for monitoring

278:    Options Database Key:
279: .  -nep_monitor_conv - activates NEPMonitorConverged()

281:    Level: intermediate

283: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
284: @*/
285: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
286: {
287:   PetscInt       i;
288:   PetscViewer    viewer = vf->viewer;
289:   SlepcConvMon   ctx;

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

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

323:   PetscFunctionBegin;
324:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
325:   PetscCall(PetscNew(&mctx));
326:   mctx->ctx = ctx;
327:   (*vf)->data = (void*)mctx;
328:   PetscFunctionReturn(PETSC_SUCCESS);
329: }

331: PetscErrorCode NEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
332: {
333:   PetscFunctionBegin;
334:   if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
335:   PetscCall(PetscFree((*vf)->data));
336:   PetscCall(PetscViewerDestroy(&(*vf)->viewer));
337:   PetscCall(PetscDrawLGDestroy(&(*vf)->lg));
338:   PetscCall(PetscFree(*vf));
339:   PetscFunctionReturn(PETSC_SUCCESS);
340: }

342: /*@C
343:    NEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
344:    approximation at each iteration of the nonlinear eigensolver.

346:    Collective

348:    Input Parameters:
349: +  nep    - nonlinear eigensolver context
350: .  its    - iteration number
351: .  nconv  - number of converged eigenpairs so far
352: .  eigr   - real part of the eigenvalues
353: .  eigi   - imaginary part of the eigenvalues
354: .  errest - error estimates
355: .  nest   - number of error estimates to display
356: -  vf     - viewer and format for monitoring

358:    Options Database Key:
359: .  -nep_monitor draw::draw_lg - activates NEPMonitorFirstDrawLG()

361:    Level: intermediate

363: .seealso: NEPMonitorSet()
364: @*/
365: PetscErrorCode NEPMonitorFirstDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
366: {
367:   PetscViewer    viewer = vf->viewer;
368:   PetscDrawLG    lg = vf->lg;
369:   PetscReal      x,y;

371:   PetscFunctionBegin;
375:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
376:   if (its==1) {
377:     PetscCall(PetscDrawLGReset(lg));
378:     PetscCall(PetscDrawLGSetDimension(lg,1));
379:     PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0));
380:   }
381:   if (nconv<nest) {
382:     x = (PetscReal)its;
383:     if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
384:     else y = 0.0;
385:     PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
386:     if (its <= 20 || !(its % 5) || nep->reason) {
387:       PetscCall(PetscDrawLGDraw(lg));
388:       PetscCall(PetscDrawLGSave(lg));
389:     }
390:   }
391:   PetscCall(PetscViewerPopFormat(viewer));
392:   PetscFunctionReturn(PETSC_SUCCESS);
393: }

395: /*@C
396:    NEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.

398:    Collective

400:    Input Parameters:
401: +  viewer - the viewer
402: .  format - the viewer format
403: -  ctx    - an optional user context

405:    Output Parameter:
406: .  vf     - the viewer and format context

408:    Level: intermediate

410: .seealso: NEPMonitorSet()
411: @*/
412: PetscErrorCode NEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
413: {
414:   PetscFunctionBegin;
415:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
416:   (*vf)->data = ctx;
417:   PetscCall(NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
418:   PetscFunctionReturn(PETSC_SUCCESS);
419: }

421: /*@C
422:    NEPMonitorAllDrawLG - Plots the error estimate of all unconverged
423:    approximations at each iteration of the nonlinear eigensolver.

425:    Collective

427:    Input Parameters:
428: +  nep    - nonlinear eigensolver context
429: .  its    - iteration number
430: .  nconv  - number of converged eigenpairs so far
431: .  eigr   - real part of the eigenvalues
432: .  eigi   - imaginary part of the eigenvalues
433: .  errest - error estimates
434: .  nest   - number of error estimates to display
435: -  vf     - viewer and format for monitoring

437:    Options Database Key:
438: .  -nep_monitor_all draw::draw_lg - activates NEPMonitorAllDrawLG()

440:    Level: intermediate

442: .seealso: NEPMonitorSet()
443: @*/
444: PetscErrorCode NEPMonitorAllDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
445: {
446:   PetscViewer    viewer = vf->viewer;
447:   PetscDrawLG    lg = vf->lg;
448:   PetscInt       i,n = PetscMin(nep->nev,255);
449:   PetscReal      *x,*y;

451:   PetscFunctionBegin;
455:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
456:   if (its==1) {
457:     PetscCall(PetscDrawLGReset(lg));
458:     PetscCall(PetscDrawLGSetDimension(lg,n));
459:     PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0));
460:   }
461:   PetscCall(PetscMalloc2(n,&x,n,&y));
462:   for (i=0;i<n;i++) {
463:     x[i] = (PetscReal)its;
464:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
465:     else y[i] = 0.0;
466:   }
467:   PetscCall(PetscDrawLGAddPoint(lg,x,y));
468:   if (its <= 20 || !(its % 5) || nep->reason) {
469:     PetscCall(PetscDrawLGDraw(lg));
470:     PetscCall(PetscDrawLGSave(lg));
471:   }
472:   PetscCall(PetscFree2(x,y));
473:   PetscCall(PetscViewerPopFormat(viewer));
474:   PetscFunctionReturn(PETSC_SUCCESS);
475: }

477: /*@C
478:    NEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

480:    Collective

482:    Input Parameters:
483: +  viewer - the viewer
484: .  format - the viewer format
485: -  ctx    - an optional user context

487:    Output Parameter:
488: .  vf     - the viewer and format context

490:    Level: intermediate

492: .seealso: NEPMonitorSet()
493: @*/
494: PetscErrorCode NEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
495: {
496:   PetscFunctionBegin;
497:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
498:   (*vf)->data = ctx;
499:   PetscCall(NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
500:   PetscFunctionReturn(PETSC_SUCCESS);
501: }

503: /*@C
504:    NEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
505:    at each iteration of the nonlinear eigensolver.

507:    Collective

509:    Input Parameters:
510: +  nep    - nonlinear eigensolver context
511: .  its    - iteration number
512: .  nconv  - number of converged eigenpairs so far
513: .  eigr   - real part of the eigenvalues
514: .  eigi   - imaginary part of the eigenvalues
515: .  errest - error estimates
516: .  nest   - number of error estimates to display
517: -  vf     - viewer and format for monitoring

519:    Options Database Key:
520: .  -nep_monitor_conv draw::draw_lg - activates NEPMonitorConvergedDrawLG()

522:    Level: intermediate

524: .seealso: NEPMonitorSet()
525: @*/
526: PetscErrorCode NEPMonitorConvergedDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
527: {
528:   PetscViewer      viewer = vf->viewer;
529:   PetscDrawLG      lg = vf->lg;
530:   PetscReal        x,y;

532:   PetscFunctionBegin;
536:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
537:   if (its==1) {
538:     PetscCall(PetscDrawLGReset(lg));
539:     PetscCall(PetscDrawLGSetDimension(lg,1));
540:     PetscCall(PetscDrawLGSetLimits(lg,1,1,0,nep->nev));
541:   }
542:   x = (PetscReal)its;
543:   y = (PetscReal)nep->nconv;
544:   PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
545:   if (its <= 20 || !(its % 5) || nep->reason) {
546:     PetscCall(PetscDrawLGDraw(lg));
547:     PetscCall(PetscDrawLGSave(lg));
548:   }
549:   PetscCall(PetscViewerPopFormat(viewer));
550:   PetscFunctionReturn(PETSC_SUCCESS);
551: }

553: /*@C
554:    NEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

556:    Collective

558:    Input Parameters:
559: +  viewer - the viewer
560: .  format - the viewer format
561: -  ctx    - an optional user context

563:    Output Parameter:
564: .  vf     - the viewer and format context

566:    Level: intermediate

568: .seealso: NEPMonitorSet()
569: @*/
570: PetscErrorCode NEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
571: {
572:   SlepcConvMon   mctx;

574:   PetscFunctionBegin;
575:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
576:   PetscCall(PetscNew(&mctx));
577:   mctx->ctx = ctx;
578:   (*vf)->data = (void*)mctx;
579:   PetscCall(NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
580:   PetscFunctionReturn(PETSC_SUCCESS);
581: }