Actual source code: epsmon.c

slepc-3.17.1 2022-04-11
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:    EPS routines related to monitors
 12: */

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

 17: PetscErrorCode EPSMonitorLGCreate(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:   PetscFunctionReturn(0);
 33: }

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

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

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

 50:    Logically Collective on eps

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

 62: +  eps    - eigensolver context obtained from EPSCreate()
 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 - relative error estimates for each eigenpair
 68: .  nest   - number of error estimates
 69: -  mctx   - optional monitoring context, as set by EPSMonitorSet()

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

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

 90:    Level: intermediate

 92: .seealso: EPSMonitorFirst(), EPSMonitorAll(), EPSMonitorCancel()
 93: @*/
 94: PetscErrorCode EPSMonitorSet(EPS eps,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
 95: {
 98:   eps->monitor[eps->numbermonitors]           = monitor;
 99:   eps->monitorcontext[eps->numbermonitors]    = (void*)mctx;
100:   eps->monitordestroy[eps->numbermonitors++]  = monitordestroy;
101:   PetscFunctionReturn(0);
102: }

104: /*@
105:    EPSMonitorCancel - Clears all monitors for an EPS object.

107:    Logically Collective on eps

109:    Input Parameters:
110: .  eps - eigensolver context obtained from EPSCreate()

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

117:    Level: intermediate

119: .seealso: EPSMonitorSet()
120: @*/
121: PetscErrorCode EPSMonitorCancel(EPS eps)
122: {
123:   PetscInt       i;

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

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

137:    Not Collective

139:    Input Parameter:
140: .  eps - eigensolver context obtained from EPSCreate()

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

145:    Level: intermediate

147: .seealso: EPSMonitorSet()
148: @*/
149: PetscErrorCode EPSGetMonitorContext(EPS eps,void *ctx)
150: {
152:   *(void**)ctx = eps->monitorcontext[0];
153:   PetscFunctionReturn(0);
154: }

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

160:    Collective on eps

162:    Input Parameters:
163: +  eps    - 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: .  -eps_monitor - activates EPSMonitorFirst()

175:    Level: intermediate

177: .seealso: EPSMonitorSet(), EPSMonitorAll(), EPSMonitorConverged()
178: @*/
179: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
180: {
181:   PetscScalar    er,ei;
182:   PetscViewer    viewer = vf->viewer;

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

208: /*@C
209:    EPSMonitorAll - Print the current approximate values and
210:    error estimates at each iteration of the eigensolver.

212:    Collective on eps

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

224:    Options Database Key:
225: .  -eps_monitor_all - activates EPSMonitorAll()

227:    Level: intermediate

229: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorConverged()
230: @*/
231: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
232: {
233:   PetscInt       i;
234:   PetscScalar    er,ei;
235:   PetscViewer    viewer = vf->viewer;

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

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

266:    Collective on eps

268:    Input Parameters:
269: +  eps    - 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: .  -eps_monitor_conv - activates EPSMonitorConverged()

281:    Level: intermediate

283: .seealso: EPSMonitorSet(), EPSMonitorFirst(), EPSMonitorAll()
284: @*/
285: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
286: {
287:   PetscInt       i;
288:   PetscScalar    er,ei;
289:   PetscViewer    viewer = vf->viewer;
290:   SlepcConvMon   ctx;

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

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

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

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

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

346:    Collective on eps

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

359:    Options Database Key:
360: .  -eps_monitor draw::draw_lg - activates EPSMonitorFirstDrawLG()

362:    Level: intermediate

364: .seealso: EPSMonitorSet()
365: @*/
366: PetscErrorCode EPSMonitorFirstDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
367: {
368:   PetscViewer    viewer = vf->viewer;
369:   PetscDrawLG    lg = vf->lg;
370:   PetscReal      x,y;

375:   PetscViewerPushFormat(viewer,vf->format);
376:   if (its==1) {
377:     PetscDrawLGReset(lg);
378:     PetscDrawLGSetDimension(lg,1);
379:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->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:     PetscDrawLGAddPoint(lg,&x,&y);
386:     if (its <= 20 || !(its % 5) || eps->reason) {
387:       PetscDrawLGDraw(lg);
388:       PetscDrawLGSave(lg);
389:     }
390:   }
391:   PetscViewerPopFormat(viewer);
392:   PetscFunctionReturn(0);
393: }

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

398:    Collective on viewer

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: EPSMonitorSet()
411: @*/
412: PetscErrorCode EPSMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
413: {
414:   PetscViewerAndFormatCreate(viewer,format,vf);
415:   (*vf)->data = ctx;
416:   EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
417:   PetscFunctionReturn(0);
418: }

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

424:    Collective on eps

426:    Input Parameters:
427: +  eps    - eigensolver context
428: .  its    - iteration number
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: .  -eps_monitor_all draw::draw_lg - activates EPSMonitorAllDrawLG()

440:    Level: intermediate

442: .seealso: EPSMonitorSet()
443: @*/
444: PetscErrorCode EPSMonitorAllDrawLG(EPS eps,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(eps->nev,255);
449:   PetscReal      *x,*y;

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

476: /*@C
477:    EPSMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

479:    Collective on viewer

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

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

489:    Level: intermediate

491: .seealso: EPSMonitorSet()
492: @*/
493: PetscErrorCode EPSMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
494: {
495:   PetscViewerAndFormatCreate(viewer,format,vf);
496:   (*vf)->data = ctx;
497:   EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
498:   PetscFunctionReturn(0);
499: }

501: /*@C
502:    EPSMonitorConvergedDrawLG - Plots the number of converged eigenvalues
503:    at each iteration of the eigensolver.

505:    Collective on eps

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

518:    Options Database Key:
519: .  -eps_monitor_conv draw::draw_lg - activates EPSMonitorConvergedDrawLG()

521:    Level: intermediate

523: .seealso: EPSMonitorSet()
524: @*/
525: PetscErrorCode EPSMonitorConvergedDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
526: {
527:   PetscViewer      viewer = vf->viewer;
528:   PetscDrawLG      lg = vf->lg;
529:   PetscReal        x,y;

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

551: /*@C
552:    EPSMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

554:    Collective on viewer

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

561:    Output Parameter:
562: .  vf     - the viewer and format context

564:    Level: intermediate

566: .seealso: EPSMonitorSet()
567: @*/
568: PetscErrorCode EPSMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
569: {
570:   SlepcConvMon   mctx;

572:   PetscViewerAndFormatCreate(viewer,format,vf);
573:   PetscNew(&mctx);
574:   mctx->ctx = ctx;
575:   (*vf)->data = (void*)mctx;
576:   EPSMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
577:   PetscFunctionReturn(0);
578: }