Actual source code: pepmon.c
slepc-3.23.0 2025-03-29
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: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode PEPMonitor(PEP pep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
22: PetscInt i,n = pep->numbermonitors;
24: PetscFunctionBegin;
25: for (i=0;i<n;i++) PetscCall((*pep->monitor[i])(pep,it,nconv,eigr,eigi,errest,nest,pep->monitorcontext[i]));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@C
30: PEPMonitorSet - Sets an ADDITIONAL function to be called at every
31: iteration to monitor the error estimates for each requested eigenpair.
33: Logically Collective
35: Input Parameters:
36: + pep - eigensolver context obtained from PEPCreate()
37: . monitor - pointer to function (if this is NULL, it turns off monitoring)
38: . mctx - [optional] context for private data for the
39: monitor routine (use NULL if no context is desired)
40: - monitordestroy - [optional] routine that frees monitor context (may be NULL),
41: see PetscCtxDestroyFn for the calling sequence
43: Calling sequence of monitor:
44: $ PetscErrorCode monitor(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx)
45: + pep - polynomial eigensolver context obtained from PEPCreate()
46: . its - iteration number
47: . nconv - number of converged eigenpairs
48: . eigr - real part of the eigenvalues
49: . eigi - imaginary part of the eigenvalues
50: . errest - relative error estimates for each eigenpair
51: . nest - number of error estimates
52: - mctx - optional monitoring context, as set by PEPMonitorSet()
54: Options Database Keys:
55: + -pep_monitor - print only the first error estimate
56: . -pep_monitor_all - print error estimates at each iteration
57: . -pep_monitor_conv - print the eigenvalue approximations only when
58: convergence has been reached
59: . -pep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
60: approximate eigenvalue
61: . -pep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
62: approximate eigenvalues
63: . -pep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
64: - -pep_monitor_cancel - cancels all monitors that have been hardwired into
65: a code by calls to PEPMonitorSet(), but does not cancel those set via
66: the options database.
68: Notes:
69: Several different monitoring routines may be set by calling
70: PEPMonitorSet() multiple times; all will be called in the
71: order in which they were set.
73: Level: intermediate
75: .seealso: PEPMonitorFirst(), PEPMonitorAll(), PEPMonitorCancel()
76: @*/
77: 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)
78: {
79: PetscFunctionBegin;
81: PetscCheck(pep->numbermonitors<MAXPEPMONITORS,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Too many PEP monitors set");
82: pep->monitor[pep->numbermonitors] = monitor;
83: pep->monitorcontext[pep->numbermonitors] = (void*)mctx;
84: pep->monitordestroy[pep->numbermonitors++] = monitordestroy;
85: PetscFunctionReturn(PETSC_SUCCESS);
86: }
88: /*@
89: PEPMonitorCancel - Clears all monitors for a PEP object.
91: Logically Collective
93: Input Parameters:
94: . pep - eigensolver context obtained from PEPCreate()
96: Options Database Key:
97: . -pep_monitor_cancel - Cancels all monitors that have been hardwired
98: into a code by calls to PEPMonitorSet(),
99: but does not cancel those set via the options database.
101: Level: intermediate
103: .seealso: PEPMonitorSet()
104: @*/
105: PetscErrorCode PEPMonitorCancel(PEP pep)
106: {
107: PetscInt i;
109: PetscFunctionBegin;
111: for (i=0; i<pep->numbermonitors; i++) {
112: if (pep->monitordestroy[i]) PetscCall((*pep->monitordestroy[i])(&pep->monitorcontext[i]));
113: }
114: pep->numbermonitors = 0;
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }
118: /*@C
119: PEPGetMonitorContext - Gets the monitor context, as set by
120: PEPMonitorSet() for the FIRST monitor only.
122: Not Collective
124: Input Parameter:
125: . pep - eigensolver context obtained from PEPCreate()
127: Output Parameter:
128: . ctx - monitor context
130: Level: intermediate
132: .seealso: PEPMonitorSet(), PEPDefaultMonitor()
133: @*/
134: PetscErrorCode PEPGetMonitorContext(PEP pep,void *ctx)
135: {
136: PetscFunctionBegin;
138: *(void**)ctx = pep->monitorcontext[0];
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: /*
143: Helper function to compute eigenvalue that must be viewed in monitor
144: */
145: static PetscErrorCode PEPMonitorGetTrueEig(PEP pep,PetscScalar *er,PetscScalar *ei)
146: {
147: PetscBool flg,sinv;
149: PetscFunctionBegin;
150: PetscCall(STBackTransform(pep->st,1,er,ei));
151: if (pep->sfactor!=1.0) {
152: PetscCall(STGetTransform(pep->st,&flg));
153: PetscCall(PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv));
154: if (flg && sinv) {
155: *er /= pep->sfactor; *ei /= pep->sfactor;
156: } else {
157: *er *= pep->sfactor; *ei *= pep->sfactor;
158: }
159: }
160: PetscFunctionReturn(PETSC_SUCCESS);
161: }
163: /*@C
164: PEPMonitorFirst - Print the first unconverged approximate value and
165: error estimate at each iteration of the polynomial eigensolver.
167: Collective
169: Input Parameters:
170: + pep - polynomial eigensolver context
171: . its - iteration number
172: . nconv - number of converged eigenpairs so far
173: . eigr - real part of the eigenvalues
174: . eigi - imaginary part of the eigenvalues
175: . errest - error estimates
176: . nest - number of error estimates to display
177: - vf - viewer and format for monitoring
179: Options Database Key:
180: . -pep_monitor - activates PEPMonitorFirst()
182: Level: intermediate
184: .seealso: PEPMonitorSet(), PEPMonitorAll(), PEPMonitorConverged()
185: @*/
186: PetscErrorCode PEPMonitorFirst(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
187: {
188: PetscScalar er,ei;
189: PetscViewer viewer = vf->viewer;
191: PetscFunctionBegin;
194: if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix));
195: if (nconv<nest) {
196: PetscCall(PetscViewerPushFormat(viewer,vf->format));
197: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
198: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
199: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
200: er = eigr[nconv]; ei = eigi[nconv];
201: PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
202: #if defined(PETSC_USE_COMPLEX)
203: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
204: #else
205: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
206: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
207: #endif
208: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
209: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
210: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
211: PetscCall(PetscViewerPopFormat(viewer));
212: }
213: PetscFunctionReturn(PETSC_SUCCESS);
214: }
216: /*@C
217: PEPMonitorAll - Print the current approximate values and
218: error estimates at each iteration of the polynomial eigensolver.
220: Collective
222: Input Parameters:
223: + pep - polynomial eigensolver context
224: . its - iteration number
225: . nconv - number of converged eigenpairs so far
226: . eigr - real part of the eigenvalues
227: . eigi - imaginary part of the eigenvalues
228: . errest - error estimates
229: . nest - number of error estimates to display
230: - vf - viewer and format for monitoring
232: Options Database Key:
233: . -pep_monitor_all - activates PEPMonitorAll()
235: Level: intermediate
237: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorConverged()
238: @*/
239: PetscErrorCode PEPMonitorAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
240: {
241: PetscInt i;
242: PetscScalar er,ei;
243: PetscViewer viewer = vf->viewer;
245: PetscFunctionBegin;
248: PetscCall(PetscViewerPushFormat(viewer,vf->format));
249: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
250: if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix));
251: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
252: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
253: for (i=0;i<nest;i++) {
254: er = eigr[i]; ei = eigi[i];
255: PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
256: #if defined(PETSC_USE_COMPLEX)
257: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
258: #else
259: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
260: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
261: #endif
262: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
263: }
264: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
265: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
266: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
267: PetscCall(PetscViewerPopFormat(viewer));
268: PetscFunctionReturn(PETSC_SUCCESS);
269: }
271: /*@C
272: PEPMonitorConverged - Print the approximate values and
273: error estimates as they converge.
275: Collective
277: Input Parameters:
278: + pep - polynomial eigensolver context
279: . its - iteration number
280: . nconv - number of converged eigenpairs so far
281: . eigr - real part of the eigenvalues
282: . eigi - imaginary part of the eigenvalues
283: . errest - error estimates
284: . nest - number of error estimates to display
285: - vf - viewer and format for monitoring
287: Options Database Key:
288: . -pep_monitor_conv - activates PEPMonitorConverged()
290: Level: intermediate
292: .seealso: PEPMonitorSet(), PEPMonitorFirst(), PEPMonitorAll()
293: @*/
294: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
295: {
296: PetscInt i;
297: PetscScalar er,ei;
298: PetscViewer viewer = vf->viewer;
299: SlepcConvMon ctx;
301: PetscFunctionBegin;
304: ctx = (SlepcConvMon)vf->data;
305: if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)pep)->prefix));
306: if (its==1) ctx->oldnconv = 0;
307: if (ctx->oldnconv!=nconv) {
308: PetscCall(PetscViewerPushFormat(viewer,vf->format));
309: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
310: for (i=ctx->oldnconv;i<nconv;i++) {
311: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP converged value (error) #%" PetscInt_FMT,its,i));
312: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
313: er = eigr[i]; ei = eigi[i];
314: PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
315: #if defined(PETSC_USE_COMPLEX)
316: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
317: #else
318: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
319: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
320: #endif
321: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
322: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
323: }
324: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
325: PetscCall(PetscViewerPopFormat(viewer));
326: ctx->oldnconv = nconv;
327: }
328: PetscFunctionReturn(PETSC_SUCCESS);
329: }
331: PetscErrorCode PEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
332: {
333: SlepcConvMon mctx;
335: PetscFunctionBegin;
336: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
337: PetscCall(PetscNew(&mctx));
338: mctx->ctx = ctx;
339: (*vf)->data = (void*)mctx;
340: PetscFunctionReturn(PETSC_SUCCESS);
341: }
343: PetscErrorCode PEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
344: {
345: PetscFunctionBegin;
346: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
347: PetscCall(PetscFree((*vf)->data));
348: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
349: PetscCall(PetscFree(*vf));
350: PetscFunctionReturn(PETSC_SUCCESS);
351: }
353: /*@C
354: PEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
355: approximation at each iteration of the polynomial eigensolver.
357: Collective
359: Input Parameters:
360: + pep - polynomial eigensolver context
361: . its - iteration number
362: . nconv - number of converged eigenpairs so far
363: . eigr - real part of the eigenvalues
364: . eigi - imaginary part of the eigenvalues
365: . errest - error estimates
366: . nest - number of error estimates to display
367: - vf - viewer and format for monitoring
369: Options Database Key:
370: . -pep_monitor draw::draw_lg - activates PEPMonitorFirstDrawLG()
372: Level: intermediate
374: .seealso: PEPMonitorSet()
375: @*/
376: PetscErrorCode PEPMonitorFirstDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
377: {
378: PetscViewer viewer = vf->viewer;
379: PetscDrawLG lg;
380: PetscReal x,y;
382: PetscFunctionBegin;
385: PetscCall(PetscViewerPushFormat(viewer,vf->format));
386: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
387: if (its==1) {
388: PetscCall(PetscDrawLGReset(lg));
389: PetscCall(PetscDrawLGSetDimension(lg,1));
390: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0));
391: }
392: if (nconv<nest) {
393: x = (PetscReal)its;
394: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
395: else y = 0.0;
396: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
397: if (its <= 20 || !(its % 5) || pep->reason) {
398: PetscCall(PetscDrawLGDraw(lg));
399: PetscCall(PetscDrawLGSave(lg));
400: }
401: }
402: PetscCall(PetscViewerPopFormat(viewer));
403: PetscFunctionReturn(PETSC_SUCCESS);
404: }
406: /*@C
407: PEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
409: Collective
411: Input Parameters:
412: + viewer - the viewer
413: . format - the viewer format
414: - ctx - an optional user context
416: Output Parameter:
417: . vf - the viewer and format context
419: Level: intermediate
421: .seealso: PEPMonitorSet()
422: @*/
423: PetscErrorCode PEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
424: {
425: PetscFunctionBegin;
426: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
427: (*vf)->data = ctx;
428: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
429: PetscFunctionReturn(PETSC_SUCCESS);
430: }
432: /*@C
433: PEPMonitorAllDrawLG - Plots the error estimate of all unconverged
434: approximations at each iteration of the polynomial eigensolver.
436: Collective
438: Input Parameters:
439: + pep - polynomial eigensolver context
440: . its - iteration number
441: . nconv - number of converged eigenpairs so far
442: . eigr - real part of the eigenvalues
443: . eigi - imaginary part of the eigenvalues
444: . errest - error estimates
445: . nest - number of error estimates to display
446: - vf - viewer and format for monitoring
448: Options Database Key:
449: . -pep_monitor_all draw::draw_lg - activates PEPMonitorAllDrawLG()
451: Level: intermediate
453: .seealso: PEPMonitorSet()
454: @*/
455: PetscErrorCode PEPMonitorAllDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
456: {
457: PetscViewer viewer = vf->viewer;
458: PetscDrawLG lg;
459: PetscInt i,n = PetscMin(pep->nev,255);
460: PetscReal *x,*y;
462: PetscFunctionBegin;
465: PetscCall(PetscViewerPushFormat(viewer,vf->format));
466: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
467: if (its==1) {
468: PetscCall(PetscDrawLGReset(lg));
469: PetscCall(PetscDrawLGSetDimension(lg,n));
470: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0));
471: }
472: PetscCall(PetscMalloc2(n,&x,n,&y));
473: for (i=0;i<n;i++) {
474: x[i] = (PetscReal)its;
475: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
476: else y[i] = 0.0;
477: }
478: PetscCall(PetscDrawLGAddPoint(lg,x,y));
479: if (its <= 20 || !(its % 5) || pep->reason) {
480: PetscCall(PetscDrawLGDraw(lg));
481: PetscCall(PetscDrawLGSave(lg));
482: }
483: PetscCall(PetscFree2(x,y));
484: PetscCall(PetscViewerPopFormat(viewer));
485: PetscFunctionReturn(PETSC_SUCCESS);
486: }
488: /*@C
489: PEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
491: Collective
493: Input Parameters:
494: + viewer - the viewer
495: . format - the viewer format
496: - ctx - an optional user context
498: Output Parameter:
499: . vf - the viewer and format context
501: Level: intermediate
503: .seealso: PEPMonitorSet()
504: @*/
505: PetscErrorCode PEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
506: {
507: PetscFunctionBegin;
508: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
509: (*vf)->data = ctx;
510: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
511: PetscFunctionReturn(PETSC_SUCCESS);
512: }
514: /*@C
515: PEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
516: at each iteration of the polynomial eigensolver.
518: Collective
520: Input Parameters:
521: + pep - polynomial eigensolver context
522: . its - iteration number
523: . nconv - number of converged eigenpairs so far
524: . eigr - real part of the eigenvalues
525: . eigi - imaginary part of the eigenvalues
526: . errest - error estimates
527: . nest - number of error estimates to display
528: - vf - viewer and format for monitoring
530: Options Database Key:
531: . -pep_monitor_conv draw::draw_lg - activates PEPMonitorConvergedDrawLG()
533: Level: intermediate
535: .seealso: PEPMonitorSet()
536: @*/
537: PetscErrorCode PEPMonitorConvergedDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
538: {
539: PetscViewer viewer = vf->viewer;
540: PetscDrawLG lg;
541: PetscReal x,y;
543: PetscFunctionBegin;
546: PetscCall(PetscViewerPushFormat(viewer,vf->format));
547: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
548: if (its==1) {
549: PetscCall(PetscDrawLGReset(lg));
550: PetscCall(PetscDrawLGSetDimension(lg,1));
551: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,pep->nev));
552: }
553: x = (PetscReal)its;
554: y = (PetscReal)pep->nconv;
555: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
556: if (its <= 20 || !(its % 5) || pep->reason) {
557: PetscCall(PetscDrawLGDraw(lg));
558: PetscCall(PetscDrawLGSave(lg));
559: }
560: PetscCall(PetscViewerPopFormat(viewer));
561: PetscFunctionReturn(PETSC_SUCCESS);
562: }
564: /*@C
565: PEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
567: Collective
569: Input Parameters:
570: + viewer - the viewer
571: . format - the viewer format
572: - ctx - an optional user context
574: Output Parameter:
575: . vf - the viewer and format context
577: Level: intermediate
579: .seealso: PEPMonitorSet()
580: @*/
581: PetscErrorCode PEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
582: {
583: SlepcConvMon mctx;
585: PetscFunctionBegin;
586: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
587: PetscCall(PetscNew(&mctx));
588: mctx->ctx = ctx;
589: (*vf)->data = (void*)mctx;
590: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
591: PetscFunctionReturn(PETSC_SUCCESS);
592: }