Actual source code: nepmon.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:    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;

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

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

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

 56:    Logically Collective on nep

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

 68: +  nep    - nonlinear eigensolver context obtained from NEPCreate()
 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 - error estimates for each eigenpair
 74: .  nest   - number of error estimates
 75: -  mctx   - optional monitoring context, as set by NEPMonitorSet()

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

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

 96:    Level: intermediate

 98: .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
 99: @*/
100: PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
101: {
104:   if (nep->numbermonitors >= MAXNEPMONITORS) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
105:   nep->monitor[nep->numbermonitors]           = monitor;
106:   nep->monitorcontext[nep->numbermonitors]    = (void*)mctx;
107:   nep->monitordestroy[nep->numbermonitors++]  = monitordestroy;
108:   return(0);
109: }

111: /*@
112:    NEPMonitorCancel - Clears all monitors for a NEP object.

114:    Logically Collective on nep

116:    Input Parameters:
117: .  nep - eigensolver context obtained from NEPCreate()

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

124:    Level: intermediate

126: .seealso: NEPMonitorSet()
127: @*/
128: PetscErrorCode NEPMonitorCancel(NEP nep)
129: {
131:   PetscInt       i;

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

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

148:    Not Collective

150:    Input Parameter:
151: .  nep - eigensolver context obtained from NEPCreate()

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

156:    Level: intermediate

158: .seealso: NEPMonitorSet(), NEPDefaultMonitor()
159: @*/
160: PetscErrorCode NEPGetMonitorContext(NEP nep,void **ctx)
161: {
164:   *ctx = nep->monitorcontext[0];
165:   return(0);
166: }

168: /*@C
169:    NEPMonitorFirst - Print the first unconverged approximate value and
170:    error estimate at each iteration of the nonlinear eigensolver.

172:    Collective on nep

174:    Input Parameters:
175: +  nep    - nonlinear eigensolver context
176: .  its    - iteration number
177: .  nconv  - number of converged eigenpairs so far
178: .  eigr   - real part of the eigenvalues
179: .  eigi   - imaginary part of the eigenvalues
180: .  errest - error estimates
181: .  nest   - number of error estimates to display
182: -  vf     - viewer and format for monitoring

184:    Options Database Key:
185: .  -nep_monitor - activates NEPMonitorFirst()

187:    Level: intermediate

189: .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
190: @*/
191: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
192: {
194:   PetscViewer    viewer;

199:   viewer = vf->viewer;
201:   if (its==1 && ((PetscObject)nep)->prefix) {
202:     PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
203:   }
204:   if (nconv<nest) {
205:     PetscViewerPushFormat(viewer,vf->format);
206:     PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
207:     PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D first unconverged value (error)",its,nconv);
208:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
209: #if defined(PETSC_USE_COMPLEX)
210:     PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv]));
211: #else
212:     PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]);
213:     if (eigi[nconv]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]); }
214: #endif
215:     PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]);
216:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
217:     PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
218:     PetscViewerPopFormat(viewer);
219:   }
220:   return(0);
221: }

223: /*@C
224:    NEPMonitorAll - Print the current approximate values and
225:    error estimates at each iteration of the nonlinear eigensolver.

227:    Collective on nep

229:    Input Parameters:
230: +  nep    - nonlinear eigensolver context
231: .  its    - iteration number
232: .  nconv  - number of converged eigenpairs so far
233: .  eigr   - real part of the eigenvalues
234: .  eigi   - imaginary part of the eigenvalues
235: .  errest - error estimates
236: .  nest   - number of error estimates to display
237: -  vf     - viewer and format for monitoring

239:    Options Database Key:
240: .  -nep_monitor_all - activates NEPMonitorAll()

242:    Level: intermediate

244: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
245: @*/
246: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
247: {
249:   PetscInt       i;
250:   PetscViewer    viewer;

255:   viewer = vf->viewer;
257:   PetscViewerPushFormat(viewer,vf->format);
258:   PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
259:   if (its==1 && ((PetscObject)nep)->prefix) {
260:     PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix);
261:   }
262:   PetscViewerASCIIPrintf(viewer,"%3D NEP nconv=%D Values (Errors)",its,nconv);
263:   PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
264:   for (i=0;i<nest;i++) {
265: #if defined(PETSC_USE_COMPLEX)
266:     PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
267: #else
268:     PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
269:     if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
270: #endif
271:     PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]);
272:   }
273:   PetscViewerASCIIPrintf(viewer,"\n");
274:   PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
275:   PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
276:   PetscViewerPopFormat(viewer);
277:   return(0);
278: }

280: /*@C
281:    NEPMonitorConverged - Print the approximate values and
282:    error estimates as they converge.

284:    Collective on nep

286:    Input Parameters:
287: +  nep    - nonlinear eigensolver context
288: .  its    - iteration number
289: .  nconv  - number of converged eigenpairs so far
290: .  eigr   - real part of the eigenvalues
291: .  eigi   - imaginary part of the eigenvalues
292: .  errest - error estimates
293: .  nest   - number of error estimates to display
294: -  vf     - viewer and format for monitoring

296:    Options Database Key:
297: .  -nep_monitor_conv - activates NEPMonitorConverged()

299:    Level: intermediate

301: .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
302: @*/
303: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
304: {
306:   PetscInt       i;
307:   PetscViewer    viewer;
308:   SlepcConvMon   ctx;

313:   viewer = vf->viewer;
315:   ctx = (SlepcConvMon)vf->data;
316:   if (its==1 && ((PetscObject)nep)->prefix) {
317:     PetscViewerASCIIPrintf(viewer,"  Convergence history for %s solve.\n",((PetscObject)nep)->prefix);
318:   }
319:   if (its==1) ctx->oldnconv = 0;
320:   if (ctx->oldnconv!=nconv) {
321:     PetscViewerPushFormat(viewer,vf->format);
322:     PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
323:     for (i=ctx->oldnconv;i<nconv;i++) {
324:       PetscViewerASCIIPrintf(viewer,"%3D NEP converged value (error) #%D",its,i);
325:       PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
326: #if defined(PETSC_USE_COMPLEX)
327:       PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i]));
328: #else
329:       PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]);
330:       if (eigi[i]!=0.0) { PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]); }
331: #endif
332:       PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]);
333:       PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
334:     }
335:     PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
336:     PetscViewerPopFormat(viewer);
337:     ctx->oldnconv = nconv;
338:   }
339:   return(0);
340: }

342: PetscErrorCode NEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
343: {
345:   SlepcConvMon   mctx;

348:   PetscViewerAndFormatCreate(viewer,format,vf);
349:   PetscNew(&mctx);
350:   mctx->ctx = ctx;
351:   (*vf)->data = (void*)mctx;
352:   return(0);
353: }

355: PetscErrorCode NEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
356: {

360:   if (!*vf) return(0);
361:   PetscFree((*vf)->data);
362:   PetscViewerDestroy(&(*vf)->viewer);
363:   PetscDrawLGDestroy(&(*vf)->lg);
364:   PetscFree(*vf);
365:   return(0);
366: }

368: /*@C
369:    NEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
370:    approximation at each iteration of the nonlinear eigensolver.

372:    Collective on nep

374:    Input Parameters:
375: +  nep    - nonlinear eigensolver context
376: .  its    - iteration number
377: .  its    - iteration number
378: .  nconv  - number of converged eigenpairs so far
379: .  eigr   - real part of the eigenvalues
380: .  eigi   - imaginary part of the eigenvalues
381: .  errest - error estimates
382: .  nest   - number of error estimates to display
383: -  vf     - viewer and format for monitoring

385:    Options Database Key:
386: .  -nep_monitor draw::draw_lg - activates NEPMonitorFirstDrawLG()

388:    Level: intermediate

390: .seealso: NEPMonitorSet()
391: @*/
392: PetscErrorCode NEPMonitorFirstDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
393: {
395:   PetscViewer    viewer;
396:   PetscDrawLG    lg;
397:   PetscReal      x,y;

402:   viewer = vf->viewer;
404:   lg = vf->lg;
406:   PetscViewerPushFormat(viewer,vf->format);
407:   if (its==1) {
408:     PetscDrawLGReset(lg);
409:     PetscDrawLGSetDimension(lg,1);
410:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0);
411:   }
412:   if (nconv<nest) {
413:     x = (PetscReal)its;
414:     if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
415:     else y = 0.0;
416:     PetscDrawLGAddPoint(lg,&x,&y);
417:     if (its <= 20 || !(its % 5) || nep->reason) {
418:       PetscDrawLGDraw(lg);
419:       PetscDrawLGSave(lg);
420:     }
421:   }
422:   PetscViewerPopFormat(viewer);
423:   return(0);
424: }

426: /*@C
427:    NEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.

429:    Collective on viewer

431:    Input Parameters:
432: +  viewer - the viewer
433: .  format - the viewer format
434: -  ctx    - an optional user context

436:    Output Parameter:
437: .  vf     - the viewer and format context

439:    Level: intermediate

441: .seealso: NEPMonitorSet()
442: @*/
443: PetscErrorCode NEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
444: {

448:   PetscViewerAndFormatCreate(viewer,format,vf);
449:   (*vf)->data = ctx;
450:   NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
451:   return(0);
452: }

454: /*@C
455:    NEPMonitorAllDrawLG - Plots the error estimate of all unconverged
456:    approximations at each iteration of the nonlinear eigensolver.

458:    Collective on nep

460:    Input Parameters:
461: +  nep    - nonlinear eigensolver context
462: .  its    - iteration number
463: .  its    - iteration number
464: .  nconv  - number of converged eigenpairs so far
465: .  eigr   - real part of the eigenvalues
466: .  eigi   - imaginary part of the eigenvalues
467: .  errest - error estimates
468: .  nest   - number of error estimates to display
469: -  vf     - viewer and format for monitoring

471:    Options Database Key:
472: .  -nep_monitor_all draw::draw_lg - activates NEPMonitorAllDrawLG()

474:    Level: intermediate

476: .seealso: NEPMonitorSet()
477: @*/
478: PetscErrorCode NEPMonitorAllDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
479: {
481:   PetscViewer    viewer;
482:   PetscDrawLG    lg;
483:   PetscInt       i,n = PetscMin(nep->nev,255);
484:   PetscReal      *x,*y;

489:   viewer = vf->viewer;
491:   lg = vf->lg;
493:   PetscViewerPushFormat(viewer,vf->format);
494:   if (its==1) {
495:     PetscDrawLGReset(lg);
496:     PetscDrawLGSetDimension(lg,n);
497:     PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0);
498:   }
499:   PetscMalloc2(n,&x,n,&y);
500:   for (i=0;i<n;i++) {
501:     x[i] = (PetscReal)its;
502:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
503:     else y[i] = 0.0;
504:   }
505:   PetscDrawLGAddPoint(lg,x,y);
506:   if (its <= 20 || !(its % 5) || nep->reason) {
507:     PetscDrawLGDraw(lg);
508:     PetscDrawLGSave(lg);
509:   }
510:   PetscFree2(x,y);
511:   PetscViewerPopFormat(viewer);
512:   return(0);
513: }

515: /*@C
516:    NEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

518:    Collective on viewer

520:    Input Parameters:
521: +  viewer - the viewer
522: .  format - the viewer format
523: -  ctx    - an optional user context

525:    Output Parameter:
526: .  vf     - the viewer and format context

528:    Level: intermediate

530: .seealso: NEPMonitorSet()
531: @*/
532: PetscErrorCode NEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
533: {

537:   PetscViewerAndFormatCreate(viewer,format,vf);
538:   (*vf)->data = ctx;
539:   NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
540:   return(0);
541: }

543: /*@C
544:    NEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
545:    at each iteration of the nonlinear eigensolver.

547:    Collective on nep

549:    Input Parameters:
550: +  nep    - nonlinear eigensolver context
551: .  its    - iteration number
552: .  its    - iteration number
553: .  nconv  - number of converged eigenpairs so far
554: .  eigr   - real part of the eigenvalues
555: .  eigi   - imaginary part of the eigenvalues
556: .  errest - error estimates
557: .  nest   - number of error estimates to display
558: -  vf     - viewer and format for monitoring

560:    Options Database Key:
561: .  -nep_monitor_conv draw::draw_lg - activates NEPMonitorConvergedDrawLG()

563:    Level: intermediate

565: .seealso: NEPMonitorSet()
566: @*/
567: PetscErrorCode NEPMonitorConvergedDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
568: {
569:   PetscErrorCode   ierr;
570:   PetscViewer      viewer;
571:   PetscDrawLG      lg;
572:   PetscReal        x,y;

577:   viewer = vf->viewer;
579:   lg = vf->lg;
581:   PetscViewerPushFormat(viewer,vf->format);
582:   if (its==1) {
583:     PetscDrawLGReset(lg);
584:     PetscDrawLGSetDimension(lg,1);
585:     PetscDrawLGSetLimits(lg,1,1,0,nep->nev);
586:   }
587:   x = (PetscReal)its;
588:   y = (PetscReal)nep->nconv;
589:   PetscDrawLGAddPoint(lg,&x,&y);
590:   if (its <= 20 || !(its % 5) || nep->reason) {
591:     PetscDrawLGDraw(lg);
592:     PetscDrawLGSave(lg);
593:   }
594:   PetscViewerPopFormat(viewer);
595:   return(0);
596: }

598: /*@C
599:    NEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

601:    Collective on viewer

603:    Input Parameters:
604: +  viewer - the viewer
605: .  format - the viewer format
606: -  ctx    - an optional user context

608:    Output Parameter:
609: .  vf     - the viewer and format context

611:    Level: intermediate

613: .seealso: NEPMonitorSet()
614: @*/
615: PetscErrorCode NEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
616: {
618:   SlepcConvMon   mctx;

621:   PetscViewerAndFormatCreate(viewer,format,vf);
622:   PetscNew(&mctx);
623:   mctx->ctx = ctx;
624:   (*vf)->data = (void*)mctx;
625:   NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg);
626:   return(0);
627: }