Actual source code: pepmon.c

slepc-main 2024-11-15
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:    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;

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

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

 48: /*@C
 49:    PEPMonitorSet - 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: +  pep     - eigensolver context obtained from PEPCreate()
 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),
 60:              see PetscCtxDestroyFn for the calling sequence

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

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

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

 92:    Level: intermediate

 94: .seealso: PEPMonitorFirst(), PEPMonitorAll(), PEPMonitorCancel()
 95: @*/
 96: PetscErrorCode PEPMonitorSet(PEP pep,PetscErrorCode (*monitor)(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx),void *mctx,PetscCtxDestroyFn *monitordestroy)
 97: {
 98:   PetscFunctionBegin;
100:   PetscCheck(pep->numbermonitors<MAXPEPMONITORS,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Too many PEP monitors set");
101:   pep->monitor[pep->numbermonitors]           = monitor;
102:   pep->monitorcontext[pep->numbermonitors]    = (void*)mctx;
103:   pep->monitordestroy[pep->numbermonitors++]  = monitordestroy;
104:   PetscFunctionReturn(PETSC_SUCCESS);
105: }

107: /*@
108:    PEPMonitorCancel - Clears all monitors for a PEP object.

110:    Logically Collective

112:    Input Parameters:
113: .  pep - eigensolver context obtained from PEPCreate()

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

120:    Level: intermediate

122: .seealso: PEPMonitorSet()
123: @*/
124: PetscErrorCode PEPMonitorCancel(PEP pep)
125: {
126:   PetscInt       i;

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

137: /*@C
138:    PEPGetMonitorContext - Gets the monitor context, as set by
139:    PEPMonitorSet() for the FIRST monitor only.

141:    Not Collective

143:    Input Parameter:
144: .  pep - eigensolver context obtained from PEPCreate()

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

149:    Level: intermediate

151: .seealso: PEPMonitorSet(), PEPDefaultMonitor()
152: @*/
153: PetscErrorCode PEPGetMonitorContext(PEP pep,void *ctx)
154: {
155:   PetscFunctionBegin;
157:   *(void**)ctx = pep->monitorcontext[0];
158:   PetscFunctionReturn(PETSC_SUCCESS);
159: }

161: /*
162:    Helper function to compute eigenvalue that must be viewed in monitor
163:  */
164: static PetscErrorCode PEPMonitorGetTrueEig(PEP pep,PetscScalar *er,PetscScalar *ei)
165: {
166:   PetscBool      flg,sinv;

168:   PetscFunctionBegin;
169:   PetscCall(STBackTransform(pep->st,1,er,ei));
170:   if (pep->sfactor!=1.0) {
171:     PetscCall(STGetTransform(pep->st,&flg));
172:     PetscCall(PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv));
173:     if (flg && sinv) {
174:       *er /= pep->sfactor; *ei /= pep->sfactor;
175:     } else {
176:       *er *= pep->sfactor; *ei *= pep->sfactor;
177:     }
178:   }
179:   PetscFunctionReturn(PETSC_SUCCESS);
180: }

182: /*@C
183:    PEPMonitorFirst - Print the first unconverged approximate value and
184:    error estimate at each iteration of the polynomial eigensolver.

186:    Collective

188:    Input Parameters:
189: +  pep    - polynomial eigensolver context
190: .  its    - iteration number
191: .  nconv  - number of converged eigenpairs so far
192: .  eigr   - real part of the eigenvalues
193: .  eigi   - imaginary part of the eigenvalues
194: .  errest - error estimates
195: .  nest   - number of error estimates to display
196: -  vf     - viewer and format for monitoring

198:    Options Database Key:
199: .  -pep_monitor - activates PEPMonitorFirst()

201:    Level: intermediate

203: .seealso: PEPMonitorSet(), PEPMonitorAll(), PEPMonitorConverged()
204: @*/
205: PetscErrorCode PEPMonitorFirst(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
206: {
207:   PetscScalar    er,ei;
208:   PetscViewer    viewer = vf->viewer;

210:   PetscFunctionBegin;
213:   if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix));
214:   if (nconv<nest) {
215:     PetscCall(PetscViewerPushFormat(viewer,vf->format));
216:     PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
217:     PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
218:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
219:     er = eigr[nconv]; ei = eigi[nconv];
220:     PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
221: #if defined(PETSC_USE_COMPLEX)
222:     PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
223: #else
224:     PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
225:     if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
226: #endif
227:     PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
228:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
229:     PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
230:     PetscCall(PetscViewerPopFormat(viewer));
231:   }
232:   PetscFunctionReturn(PETSC_SUCCESS);
233: }

235: /*@C
236:    PEPMonitorAll - Print the current approximate values and
237:    error estimates at each iteration of the polynomial eigensolver.

239:    Collective

241:    Input Parameters:
242: +  pep    - polynomial eigensolver context
243: .  its    - iteration number
244: .  nconv  - number of converged eigenpairs so far
245: .  eigr   - real part of the eigenvalues
246: .  eigi   - imaginary part of the eigenvalues
247: .  errest - error estimates
248: .  nest   - number of error estimates to display
249: -  vf     - viewer and format for monitoring

251:    Options Database Key:
252: .  -pep_monitor_all - activates PEPMonitorAll()

254:    Level: intermediate

256: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorConverged()
257: @*/
258: PetscErrorCode PEPMonitorAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
259: {
260:   PetscInt       i;
261:   PetscScalar    er,ei;
262:   PetscViewer    viewer = vf->viewer;

264:   PetscFunctionBegin;
267:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
268:   PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
269:   if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer,"  Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix));
270:   PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
271:   PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
272:   for (i=0;i<nest;i++) {
273:     er = eigr[i]; ei = eigi[i];
274:     PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
275: #if defined(PETSC_USE_COMPLEX)
276:     PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
277: #else
278:     PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
279:     if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
280: #endif
281:     PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
282:   }
283:   PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
284:   PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
285:   PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
286:   PetscCall(PetscViewerPopFormat(viewer));
287:   PetscFunctionReturn(PETSC_SUCCESS);
288: }

290: /*@C
291:    PEPMonitorConverged - Print the approximate values and
292:    error estimates as they converge.

294:    Collective

296:    Input Parameters:
297: +  pep    - polynomial eigensolver context
298: .  its    - iteration number
299: .  nconv  - number of converged eigenpairs so far
300: .  eigr   - real part of the eigenvalues
301: .  eigi   - imaginary part of the eigenvalues
302: .  errest - error estimates
303: .  nest   - number of error estimates to display
304: -  vf     - viewer and format for monitoring

306:    Options Database Key:
307: .  -pep_monitor_conv - activates PEPMonitorConverged()

309:    Level: intermediate

311: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorAll()
312: @*/
313: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
314: {
315:   PetscInt       i;
316:   PetscScalar    er,ei;
317:   PetscViewer    viewer = vf->viewer;
318:   SlepcConvMon   ctx;

320:   PetscFunctionBegin;
323:   ctx = (SlepcConvMon)vf->data;
324:   if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer,"  Convergence history for %s solve.\n",((PetscObject)pep)->prefix));
325:   if (its==1) ctx->oldnconv = 0;
326:   if (ctx->oldnconv!=nconv) {
327:     PetscCall(PetscViewerPushFormat(viewer,vf->format));
328:     PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
329:     for (i=ctx->oldnconv;i<nconv;i++) {
330:       PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP converged value (error) #%" PetscInt_FMT,its,i));
331:       PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
332:       er = eigr[i]; ei = eigi[i];
333:       PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
334: #if defined(PETSC_USE_COMPLEX)
335:       PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
336: #else
337:       PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
338:       if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
339: #endif
340:       PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
341:       PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
342:     }
343:     PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
344:     PetscCall(PetscViewerPopFormat(viewer));
345:     ctx->oldnconv = nconv;
346:   }
347:   PetscFunctionReturn(PETSC_SUCCESS);
348: }

350: PetscErrorCode PEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
351: {
352:   SlepcConvMon   mctx;

354:   PetscFunctionBegin;
355:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
356:   PetscCall(PetscNew(&mctx));
357:   mctx->ctx = ctx;
358:   (*vf)->data = (void*)mctx;
359:   PetscFunctionReturn(PETSC_SUCCESS);
360: }

362: PetscErrorCode PEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
363: {
364:   PetscFunctionBegin;
365:   if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
366:   PetscCall(PetscFree((*vf)->data));
367:   PetscCall(PetscViewerDestroy(&(*vf)->viewer));
368:   PetscCall(PetscDrawLGDestroy(&(*vf)->lg));
369:   PetscCall(PetscFree(*vf));
370:   PetscFunctionReturn(PETSC_SUCCESS);
371: }

373: /*@C
374:    PEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
375:    approximation at each iteration of the polynomial eigensolver.

377:    Collective

379:    Input Parameters:
380: +  pep    - polynomial eigensolver context
381: .  its    - iteration number
382: .  nconv  - number of converged eigenpairs so far
383: .  eigr   - real part of the eigenvalues
384: .  eigi   - imaginary part of the eigenvalues
385: .  errest - error estimates
386: .  nest   - number of error estimates to display
387: -  vf     - viewer and format for monitoring

389:    Options Database Key:
390: .  -pep_monitor draw::draw_lg - activates PEPMonitorFirstDrawLG()

392:    Level: intermediate

394: .seealso: PEPMonitorSet()
395: @*/
396: PetscErrorCode PEPMonitorFirstDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
397: {
398:   PetscViewer    viewer = vf->viewer;
399:   PetscDrawLG    lg = vf->lg;
400:   PetscReal      x,y;

402:   PetscFunctionBegin;
406:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
407:   if (its==1) {
408:     PetscCall(PetscDrawLGReset(lg));
409:     PetscCall(PetscDrawLGSetDimension(lg,1));
410:     PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->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:     PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
417:     if (its <= 20 || !(its % 5) || pep->reason) {
418:       PetscCall(PetscDrawLGDraw(lg));
419:       PetscCall(PetscDrawLGSave(lg));
420:     }
421:   }
422:   PetscCall(PetscViewerPopFormat(viewer));
423:   PetscFunctionReturn(PETSC_SUCCESS);
424: }

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

429:    Collective

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: PEPMonitorSet()
442: @*/
443: PetscErrorCode PEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
444: {
445:   PetscFunctionBegin;
446:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
447:   (*vf)->data = ctx;
448:   PetscCall(PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
449:   PetscFunctionReturn(PETSC_SUCCESS);
450: }

452: /*@C
453:    PEPMonitorAllDrawLG - Plots the error estimate of all unconverged
454:    approximations at each iteration of the polynomial eigensolver.

456:    Collective

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

468:    Options Database Key:
469: .  -pep_monitor_all draw::draw_lg - activates PEPMonitorAllDrawLG()

471:    Level: intermediate

473: .seealso: PEPMonitorSet()
474: @*/
475: PetscErrorCode PEPMonitorAllDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
476: {
477:   PetscViewer    viewer = vf->viewer;
478:   PetscDrawLG    lg = vf->lg;
479:   PetscInt       i,n = PetscMin(pep->nev,255);
480:   PetscReal      *x,*y;

482:   PetscFunctionBegin;
486:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
487:   if (its==1) {
488:     PetscCall(PetscDrawLGReset(lg));
489:     PetscCall(PetscDrawLGSetDimension(lg,n));
490:     PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0));
491:   }
492:   PetscCall(PetscMalloc2(n,&x,n,&y));
493:   for (i=0;i<n;i++) {
494:     x[i] = (PetscReal)its;
495:     if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
496:     else y[i] = 0.0;
497:   }
498:   PetscCall(PetscDrawLGAddPoint(lg,x,y));
499:   if (its <= 20 || !(its % 5) || pep->reason) {
500:     PetscCall(PetscDrawLGDraw(lg));
501:     PetscCall(PetscDrawLGSave(lg));
502:   }
503:   PetscCall(PetscFree2(x,y));
504:   PetscCall(PetscViewerPopFormat(viewer));
505:   PetscFunctionReturn(PETSC_SUCCESS);
506: }

508: /*@C
509:    PEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.

511:    Collective

513:    Input Parameters:
514: +  viewer - the viewer
515: .  format - the viewer format
516: -  ctx    - an optional user context

518:    Output Parameter:
519: .  vf     - the viewer and format context

521:    Level: intermediate

523: .seealso: PEPMonitorSet()
524: @*/
525: PetscErrorCode PEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
526: {
527:   PetscFunctionBegin;
528:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
529:   (*vf)->data = ctx;
530:   PetscCall(PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
531:   PetscFunctionReturn(PETSC_SUCCESS);
532: }

534: /*@C
535:    PEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
536:    at each iteration of the polynomial eigensolver.

538:    Collective

540:    Input Parameters:
541: +  pep    - polynomial eigensolver context
542: .  its    - iteration number
543: .  nconv  - number of converged eigenpairs so far
544: .  eigr   - real part of the eigenvalues
545: .  eigi   - imaginary part of the eigenvalues
546: .  errest - error estimates
547: .  nest   - number of error estimates to display
548: -  vf     - viewer and format for monitoring

550:    Options Database Key:
551: .  -pep_monitor_conv draw::draw_lg - activates PEPMonitorConvergedDrawLG()

553:    Level: intermediate

555: .seealso: PEPMonitorSet()
556: @*/
557: PetscErrorCode PEPMonitorConvergedDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
558: {
559:   PetscViewer      viewer = vf->viewer;
560:   PetscDrawLG      lg = vf->lg;
561:   PetscReal        x,y;

563:   PetscFunctionBegin;
567:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
568:   if (its==1) {
569:     PetscCall(PetscDrawLGReset(lg));
570:     PetscCall(PetscDrawLGSetDimension(lg,1));
571:     PetscCall(PetscDrawLGSetLimits(lg,1,1,0,pep->nev));
572:   }
573:   x = (PetscReal)its;
574:   y = (PetscReal)pep->nconv;
575:   PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
576:   if (its <= 20 || !(its % 5) || pep->reason) {
577:     PetscCall(PetscDrawLGDraw(lg));
578:     PetscCall(PetscDrawLGSave(lg));
579:   }
580:   PetscCall(PetscViewerPopFormat(viewer));
581:   PetscFunctionReturn(PETSC_SUCCESS);
582: }

584: /*@C
585:    PEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.

587:    Collective

589:    Input Parameters:
590: +  viewer - the viewer
591: .  format - the viewer format
592: -  ctx    - an optional user context

594:    Output Parameter:
595: .  vf     - the viewer and format context

597:    Level: intermediate

599: .seealso: PEPMonitorSet()
600: @*/
601: PetscErrorCode PEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
602: {
603:   SlepcConvMon   mctx;

605:   PetscFunctionBegin;
606:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
607:   PetscCall(PetscNew(&mctx));
608:   mctx->ctx = ctx;
609:   (*vf)->data = (void*)mctx;
610:   PetscCall(PEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
611:   PetscFunctionReturn(PETSC_SUCCESS);
612: }