Actual source code: pepmon.c

slepc-3.15.1 2021-05-28
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:    PEP routines related to monitors
 12: */

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

 17: PetscErrorCode PEPMonitorLGCreate(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 PEPMonitor(PEP pep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
 41: {
 43:   PetscInt       i,n = pep->numbermonitors;

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

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

 56:    Logically Collective on pep

 58:    Input Parameters:
 59: +  pep     - eigensolver context obtained from PEPCreate()
 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(PEP pep,int its,int nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal* errest,int nest,void *mctx)

 68: +  pep    - polynomial eigensolver context obtained from PEPCreate()
 69: .  its    - iteration number
 70: .  nconv  - number of converged eigenpairs
 71: .  eigr   - real part of the eigenvalues
 72: .  eigi   - imaginary part of the eigenvalues
 73: .  errest - relative error estimates for each eigenpair
 74: .  nest   - number of error estimates
 75: -  mctx   - optional monitoring context, as set by PEPMonitorSet()

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

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

 96:    Level: intermediate

 98: .seealso: PEPMonitorFirst(), PEPMonitorAll(), PEPMonitorCancel()
 99: @*/
100: PetscErrorCode PEPMonitorSet(PEP pep,PetscErrorCode (*monitor)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
101: {
104:   if (pep->numbermonitors >= MAXPEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Too many PEP monitors set");
105:   pep->monitor[pep->numbermonitors]           = monitor;
106:   pep->monitorcontext[pep->numbermonitors]    = (void*)mctx;
107:   pep->monitordestroy[pep->numbermonitors++]  = monitordestroy;
108:   return(0);
109: }

111: /*@
112:    PEPMonitorCancel - Clears all monitors for a PEP object.

114:    Logically Collective on pep

116:    Input Parameters:
117: .  pep - eigensolver context obtained from PEPCreate()

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

124:    Level: intermediate

126: .seealso: PEPMonitorSet()
127: @*/
128: PetscErrorCode PEPMonitorCancel(PEP pep)
129: {
131:   PetscInt       i;

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

144: /*@C
145:    PEPGetMonitorContext - Gets the monitor context, as set by
146:    PEPMonitorSet() for the FIRST monitor only.

148:    Not Collective

150:    Input Parameter:
151: .  pep - eigensolver context obtained from PEPCreate()

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

156:    Level: intermediate

158: .seealso: PEPMonitorSet(), PEPDefaultMonitor()
159: @*/
160: PetscErrorCode PEPGetMonitorContext(PEP pep,void **ctx)
161: {
164:   *ctx = pep->monitorcontext[0];
165:   return(0);
166: }

168: /*
169:    Helper function to compute eigenvalue that must be viewed in monitor
170:  */
171: static PetscErrorCode PEPMonitorGetTrueEig(PEP pep,PetscScalar *er,PetscScalar *ei)
172: {
174:   PetscBool      flg,sinv;

177:   STBackTransform(pep->st,1,er,ei);
178:   if (pep->sfactor!=1.0) {
179:     STGetTransform(pep->st,&flg);
180:     PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
181:     if (flg && sinv) {
182:       *er /= pep->sfactor; *ei /= pep->sfactor;
183:     } else {
184:       *er *= pep->sfactor; *ei *= pep->sfactor;
185:     }
186:   }
187:   return(0);
188: }

190: /*@C
191:    PEPMonitorFirst - Print the first unconverged approximate value and
192:    error estimate at each iteration of the polynomial eigensolver.

194:    Collective on pep

196:    Input Parameters:
197: +  pep    - polynomial eigensolver context
198: .  its    - iteration number
199: .  nconv  - number of converged eigenpairs so far
200: .  eigr   - real part of the eigenvalues
201: .  eigi   - imaginary part of the eigenvalues
202: .  errest - error estimates
203: .  nest   - number of error estimates to display
204: -  vf     - viewer and format for monitoring

206:    Options Database Key:
207: .  -pep_monitor - activates PEPMonitorFirst()

209:    Level: intermediate

211: .seealso: PEPMonitorSet(), PEPMonitorAll(), PEPMonitorConverged()
212: @*/
213: PetscErrorCode PEPMonitorFirst(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
214: {
216:   PetscScalar    er,ei;
217:   PetscViewer    viewer;

222:   viewer = vf->viewer;
224:   if (its==1 && ((PetscObject)pep)->prefix) {
225:     PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
226:   }
227:   if (nconv<nest) {
228:     PetscViewerPushFormat(viewer,vf->format);
229:     PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
230:     PetscViewerASCIIPrintf(viewer,"%3D PEP nconv=%D first unconverged value (error)",its,nconv);
231:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
232:     er = eigr[nconv]; ei = eigi[nconv];
233:     PEPMonitorGetTrueEig(pep,&er,&ei);
234: #if defined(PETSC_USE_COMPLEX)
235:     PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
236: #else
237:     PetscViewerASCIIPrintf(viewer," %g",(double)er);
238:     if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
239: #endif
240:     PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
241:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
242:     PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
243:     PetscViewerPopFormat(viewer);
244:   }
245:   return(0);
246: }

248: /*@C
249:    PEPMonitorAll - Print the current approximate values and
250:    error estimates at each iteration of the polynomial eigensolver.

252:    Collective on pep

254:    Input Parameters:
255: +  pep    - polynomial eigensolver context
256: .  its    - iteration number
257: .  nconv  - number of converged eigenpairs so far
258: .  eigr   - real part of the eigenvalues
259: .  eigi   - imaginary part of the eigenvalues
260: .  errest - error estimates
261: .  nest   - number of error estimates to display
262: -  vf     - viewer and format for monitoring

264:    Options Database Key:
265: .  -pep_monitor_all - activates PEPMonitorAll()

267:    Level: intermediate

269: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorConverged()
270: @*/
271: PetscErrorCode PEPMonitorAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
272: {
274:   PetscInt       i;
275:   PetscScalar    er,ei;
276:   PetscViewer    viewer;

281:   viewer = vf->viewer;
283:   PetscViewerPushFormat(viewer,vf->format);
284:   PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
285:   if (its==1 && ((PetscObject)pep)->prefix) {
286:     PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix);
287:   }
288:   PetscViewerASCIIPrintf(viewer,"%3D PEP nconv=%D Values (Errors)",its,nconv);
289:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
290:   for (i=0;i<nest;i++) {
291:     er = eigr[i]; ei = eigi[i];
292:     PEPMonitorGetTrueEig(pep,&er,&ei);
293: #if defined(PETSC_USE_COMPLEX)
294:     PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
295: #else
296:     PetscViewerASCIIPrintf(viewer," %g",(double)er);
297:     if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
298: #endif
299:     PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
300:   }
301:   PetscViewerASCIIPrintf(viewer,"\n");
302:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
303:   PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
304:   PetscViewerPopFormat(viewer);
305:   return(0);
306: }

308: /*@C
309:    PEPMonitorConverged - Print the approximate values and
310:    error estimates as they converge.

312:    Collective on pep

314:    Input Parameters:
315: +  pep    - polynomial eigensolver context
316: .  its    - iteration number
317: .  nconv  - number of converged eigenpairs so far
318: .  eigr   - real part of the eigenvalues
319: .  eigi   - imaginary part of the eigenvalues
320: .  errest - error estimates
321: .  nest   - number of error estimates to display
322: -  vf     - viewer and format for monitoring

324:    Options Database Key:
325: .  -pep_monitor_conv - activates PEPMonitorConverged()

327:    Level: intermediate

329: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorAll()
330: @*/
331: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
332: {
334:   PetscInt       i;
335:   PetscScalar    er,ei;
336:   PetscViewer    viewer;
337:   SlepcConvMon   ctx;

342:   viewer = vf->viewer;
344:   ctx = (SlepcConvMon)vf->data;
345:   if (its==1 && ((PetscObject)pep)->prefix) {
346:     PetscViewerASCIIPrintf(viewer,"  Convergence history for %s solve.\n",((PetscObject)pep)->prefix);
347:   }
348:   if (its==1) ctx->oldnconv = 0;
349:   if (ctx->oldnconv!=nconv) {
350:     PetscViewerPushFormat(viewer,vf->format);
351:     PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
352:     for (i=ctx->oldnconv;i<nconv;i++) {
353:       PetscViewerASCIIPrintf(viewer,"%3D PEP converged value (error) #%D",its,i);
354:       PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
355:       er = eigr[i]; ei = eigi[i];
356:       PEPMonitorGetTrueEig(pep,&er,&ei);
357: #if defined(PETSC_USE_COMPLEX)
358:       PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er));
359: #else
360:       PetscViewerASCIIPrintf(viewer," %g",(double)er);
361:       if (ei!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei); }
362: #endif
363:       PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
364:       PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
365:     }
366:     PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
367:     PetscViewerPopFormat(viewer);
368:     ctx->oldnconv = nconv;
369:   }
370:   return(0);
371: }

373: PetscErrorCode PEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
374: {
376:   SlepcConvMon   mctx;

379:   PetscViewerAndFormatCreate(viewer,format,vf);
380:   PetscNew(&mctx);
381:   mctx->ctx = ctx;
382:   (*vf)->data = (void*)mctx;
383:   return(0);
384: }

386: PetscErrorCode PEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
387: {

391:   if (!*vf) return(0);
392:   PetscFree((*vf)->data);
393:   PetscViewerDestroy(&(*vf)->viewer);
394:   PetscDrawLGDestroy(&(*vf)->lg);
395:   PetscFree(*vf);
396:   return(0);
397: }

399: /*@C
400:    PEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
401:    approximation at each iteration of the polynomial eigensolver.

403:    Collective on pep

405:    Input Parameters:
406: +  pep    - polynomial eigensolver context
407: .  its    - iteration number
408: .  its    - iteration number
409: .  nconv  - number of converged eigenpairs so far
410: .  eigr   - real part of the eigenvalues
411: .  eigi   - imaginary part of the eigenvalues
412: .  errest - error estimates
413: .  nest   - number of error estimates to display
414: -  vf     - viewer and format for monitoring

416:    Options Database Key:
417: .  -pep_monitor draw::draw_lg - activates PEPMonitorFirstDrawLG()

419:    Level: intermediate

421: .seealso: PEPMonitorSet()
422: @*/
423: PetscErrorCode PEPMonitorFirstDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
424: {
426:   PetscViewer    viewer;
427:   PetscDrawLG    lg;
428:   PetscReal      x,y;

433:   viewer = vf->viewer;
435:   lg = vf->lg;
437:   PetscViewerPushFormat(viewer,vf->format);
438:   if (its==1) {
439:     PetscDrawLGReset(lg);
440:     PetscDrawLGSetDimension(lg,1);
441:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0);
442:   }
443:   if (nconv<nest) {
444:     x = (PetscReal)its;
445:     if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
446:     else y = 0.0;
447:     PetscDrawLGAddPoint(lg,&x,&y);
448:     if (its <= 20 || !(its % 5) || pep->reason) {
449:       PetscDrawLGDraw(lg);
450:       PetscDrawLGSave(lg);
451:     }
452:   }
453:   PetscViewerPopFormat(viewer);
454:   return(0);
455: }

457: /*@C
458:    PEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.

460:    Collective on viewer

462:    Input Parameters:
463: +  viewer - the viewer
464: .  format - the viewer format
465: -  ctx    - an optional user context

467:    Output Parameter:
468: .  vf     - the viewer and format context

470:    Level: intermediate

472: .seealso: PEPMonitorSet()
473: @*/
474: PetscErrorCode PEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
475: {

479:   PetscViewerAndFormatCreate(viewer,format,vf);
480:   (*vf)->data = ctx;
481:   PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
482:   return(0);
483: }

485: /*@C
486:    PEPMonitorAllDrawLG - Plots the error estimate of all unconverged
487:    approximations at each iteration of the polynomial eigensolver.

489:    Collective on pep

491:    Input Parameters:
492: +  pep    - polynomial eigensolver context
493: .  its    - iteration number
494: .  its    - iteration number
495: .  nconv  - number of converged eigenpairs so far
496: .  eigr   - real part of the eigenvalues
497: .  eigi   - imaginary part of the eigenvalues
498: .  errest - error estimates
499: .  nest   - number of error estimates to display
500: -  vf     - viewer and format for monitoring

502:    Options Database Key:
503: .  -pep_monitor_all draw::draw_lg - activates PEPMonitorAllDrawLG()

505:    Level: intermediate

507: .seealso: PEPMonitorSet()
508: @*/
509: PetscErrorCode PEPMonitorAllDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
510: {
512:   PetscViewer    viewer;
513:   PetscDrawLG    lg;
514:   PetscInt       i,n = PetscMin(pep->nev,255);
515:   PetscReal      *x,*y;

520:   viewer = vf->viewer;
522:   lg = vf->lg;
524:   PetscViewerPushFormat(viewer,vf->format);
525:   if (its==1) {
526:     PetscDrawLGReset(lg);
527:     PetscDrawLGSetDimension(lg,n);
528:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0);
529:   }
530:   PetscMalloc2(n,&x,n,&y);
531:   for (i=0;i<n;i++) {
532:     x[i] = (PetscReal)its;
533:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
534:     else y[i] = 0.0;
535:   }
536:   PetscDrawLGAddPoint(lg,x,y);
537:   if (its <= 20 || !(its % 5) || pep->reason) {
538:     PetscDrawLGDraw(lg);
539:     PetscDrawLGSave(lg);
540:   }
541:   PetscFree2(x,y);
542:   PetscViewerPopFormat(viewer);
543:   return(0);
544: }

546: /*@C
547:    PEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

549:    Collective on viewer

551:    Input Parameters:
552: +  viewer - the viewer
553: .  format - the viewer format
554: -  ctx    - an optional user context

556:    Output Parameter:
557: .  vf     - the viewer and format context

559:    Level: intermediate

561: .seealso: PEPMonitorSet()
562: @*/
563: PetscErrorCode PEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
564: {

568:   PetscViewerAndFormatCreate(viewer,format,vf);
569:   (*vf)->data = ctx;
570:   PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
571:   return(0);
572: }

574: /*@C
575:    PEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
576:    at each iteration of the polynomial eigensolver.

578:    Collective on pep

580:    Input Parameters:
581: +  pep    - polynomial eigensolver context
582: .  its    - iteration number
583: .  its    - iteration number
584: .  nconv  - number of converged eigenpairs so far
585: .  eigr   - real part of the eigenvalues
586: .  eigi   - imaginary part of the eigenvalues
587: .  errest - error estimates
588: .  nest   - number of error estimates to display
589: -  vf     - viewer and format for monitoring

591:    Options Database Key:
592: .  -pep_monitor_conv draw::draw_lg - activates PEPMonitorConvergedDrawLG()

594:    Level: intermediate

596: .seealso: PEPMonitorSet()
597: @*/
598: PetscErrorCode PEPMonitorConvergedDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
599: {
600:   PetscErrorCode   ierr;
601:   PetscViewer      viewer;
602:   PetscDrawLG      lg;
603:   PetscReal        x,y;

608:   viewer = vf->viewer;
610:   lg = vf->lg;
612:   PetscViewerPushFormat(viewer,vf->format);
613:   if (its==1) {
614:     PetscDrawLGReset(lg);
615:     PetscDrawLGSetDimension(lg,1);
616:     PetscDrawLGSetLimits(lg,1,1,0,pep->nev);
617:   }
618:   x = (PetscReal)its;
619:   y = (PetscReal)pep->nconv;
620:   PetscDrawLGAddPoint(lg,&x,&y);
621:   if (its <= 20 || !(its % 5) || pep->reason) {
622:     PetscDrawLGDraw(lg);
623:     PetscDrawLGSave(lg);
624:   }
625:   PetscViewerPopFormat(viewer);
626:   return(0);
627: }

629: /*@C
630:    PEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

632:    Collective on viewer

634:    Input Parameters:
635: +  viewer - the viewer
636: .  format - the viewer format
637: -  ctx    - an optional user context

639:    Output Parameter:
640: .  vf     - the viewer and format context

642:    Level: intermediate

644: .seealso: PEPMonitorSet()
645: @*/
646: PetscErrorCode PEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
647: {
649:   SlepcConvMon   mctx;

652:   PetscViewerAndFormatCreate(viewer,format,vf);
653:   PetscNew(&mctx);
654:   mctx->ctx = ctx;
655:   (*vf)->data = (void*)mctx;
656:   PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
657:   return(0);
658: }