Actual source code: epsmon.c
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: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode EPSMonitor(EPS eps,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
22: PetscInt i,n = eps->numbermonitors;
24: PetscFunctionBegin;
25: for (i=0;i<n;i++) PetscCall((*eps->monitor[i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[i]));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@C
30: EPSMonitorSet - 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: + eps - the linear eigensolver context
37: . monitor - pointer to function (if this is `NULL`, it turns off monitoring),
38: see `EPSMonitorFn`
39: . ctx - [optional] context for private data for the monitor routine
40: (use `NULL` if no context is desired)
41: - monitordestroy - [optional] routine that frees monitor context (may be `NULL`),
42: see `PetscCtxDestroyFn` for the calling sequence
44: Options Database Keys:
45: + -eps_monitor - print only the first error estimate
46: . -eps_monitor_all - print error estimates at each iteration
47: . -eps_monitor_conv - print the eigenvalue approximations only when
48: convergence has been reached
49: . -eps_monitor draw::draw_lg - sets line graph monitor for the first unconverged
50: approximate eigenvalue
51: . -eps_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
52: approximate eigenvalues
53: . -eps_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
54: - -eps_monitor_cancel - cancels all monitors that have been hardwired into
55: a code by calls to `EPSMonitorSet()`, but does not cancel
56: those set via the options database.
58: Notes:
59: The options database option `-eps_monitor` and related options are the easiest way
60: to turn on `EPS` iteration monitoring.
62: `EPSMonitorRegister()` provides a way to associate an options database key with `EPS`
63: monitor function.
65: Several different monitoring routines may be set by calling `EPSMonitorSet()` multiple
66: times; all will be called in the order in which they were set.
68: Fortran Note:
69: Only a single monitor function can be set for each `EPS` object.
71: Level: intermediate
73: .seealso: [](ch:eps), `EPSMonitorFirst()`, `EPSMonitorAll()`, `EPSMonitorConverged()`, `EPSMonitorFirstDrawLG()`, `EPSMonitorAllDrawLG()`, `EPSMonitorConvergedDrawLG()`, `EPSMonitorCancel()`
74: @*/
75: PetscErrorCode EPSMonitorSet(EPS eps,EPSMonitorFn *monitor,PetscCtx ctx,PetscCtxDestroyFn *monitordestroy)
76: {
77: PetscInt i;
78: PetscBool identical;
80: PetscFunctionBegin;
82: for (i=0;i<eps->numbermonitors;i++) {
83: PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))(PetscVoidFn*)monitor,ctx,monitordestroy,(PetscErrorCode (*)(void))(PetscVoidFn*)eps->monitor[i],eps->monitorcontext[i],eps->monitordestroy[i],&identical));
84: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
85: }
86: PetscCheck(eps->numbermonitors<MAXEPSMONITORS,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Too many EPS monitors set");
87: eps->monitor[eps->numbermonitors] = monitor;
88: eps->monitorcontext[eps->numbermonitors] = ctx;
89: eps->monitordestroy[eps->numbermonitors++] = monitordestroy;
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: /*@
94: EPSMonitorCancel - Clears all monitors for an `EPS` object.
96: Logically Collective
98: Input Parameter:
99: . eps - the linear eigensolver context
101: Options Database Key:
102: . -eps_monitor_cancel - cancels all monitors that have been hardwired into a code by calls to
103: `EPSMonitorSet()`, but does not cancel those set via the options database.
105: Level: intermediate
107: .seealso: [](ch:eps), `EPSMonitorSet()`
108: @*/
109: PetscErrorCode EPSMonitorCancel(EPS eps)
110: {
111: PetscInt i;
113: PetscFunctionBegin;
115: for (i=0; i<eps->numbermonitors; i++) {
116: if (eps->monitordestroy[i]) PetscCall((*eps->monitordestroy[i])(&eps->monitorcontext[i]));
117: }
118: eps->numbermonitors = 0;
119: PetscFunctionReturn(PETSC_SUCCESS);
120: }
122: /*@C
123: EPSGetMonitorContext - Gets the monitor context, as set by
124: `EPSMonitorSet()` for the FIRST monitor only.
126: Not Collective
128: Input Parameter:
129: . eps - the linear eigensolver context
131: Output Parameter:
132: . ctx - monitor context
134: Level: intermediate
136: .seealso: [](ch:eps), `EPSMonitorSet()`
137: @*/
138: PetscErrorCode EPSGetMonitorContext(EPS eps,PetscCtxRt ctx)
139: {
140: PetscFunctionBegin;
142: *(void**)ctx = eps->monitorcontext[0];
143: PetscFunctionReturn(PETSC_SUCCESS);
144: }
146: static inline PetscErrorCode EPSMonitorPrintEval(EPS eps,PetscViewer viewer,PetscScalar er,PetscScalar ei)
147: {
148: PetscFunctionBegin;
149: if (eps->problem_type==EPS_HEP || eps->problem_type==EPS_GHEP || eps->problem_type==EPS_BSE) PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)PetscRealPart(er)));
150: else {
151: #if defined(PETSC_USE_COMPLEX)
152: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
153: #else
154: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
155: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
156: #endif
157: }
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
161: /*@C
162: EPSMonitorFirst - Print the first unconverged approximate value and
163: error estimate at each iteration of the eigensolver.
165: Collective
167: Input Parameters:
168: + eps - the linear eigensolver context
169: . its - iteration number
170: . nconv - number of converged eigenpairs so far
171: . eigr - real part of the eigenvalues
172: . eigi - imaginary part of the eigenvalues
173: . errest - error estimates
174: . nest - number of error estimates to display
175: - vf - viewer and format for monitoring
177: Options Database Key:
178: . -eps_monitor - activates `EPSMonitorFirst()`
180: Note:
181: This is not called directly by users, rather one calls `EPSMonitorSet()`, with this
182: function as an argument, to cause the monitor to be used during the `EPS` solve.
184: Level: intermediate
186: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorAll()`, `EPSMonitorConverged()`
187: @*/
188: PetscErrorCode EPSMonitorFirst(EPS eps,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
189: {
190: PetscScalar er,ei;
191: PetscViewer viewer = vf->viewer;
193: PetscFunctionBegin;
196: if (its==1 && ((PetscObject)eps)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix));
197: if (nconv<nest) {
198: PetscCall(PetscViewerPushFormat(viewer,vf->format));
199: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
200: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " EPS nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
201: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
202: er = eigr[nconv]; ei = eigi[nconv];
203: PetscCall(STBackTransform(eps->st,1,&er,&ei));
204: PetscCall(EPSMonitorPrintEval(eps,viewer,er,ei));
205: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
206: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
207: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
208: PetscCall(PetscViewerPopFormat(viewer));
209: }
210: PetscFunctionReturn(PETSC_SUCCESS);
211: }
213: /*@C
214: EPSMonitorAll - Print the current approximate values and
215: error estimates at each iteration of the eigensolver.
217: Collective
219: Input Parameters:
220: + eps - the linear eigensolver context
221: . its - iteration number
222: . nconv - number of converged eigenpairs so far
223: . eigr - real part of the eigenvalues
224: . eigi - imaginary part of the eigenvalues
225: . errest - error estimates
226: . nest - number of error estimates to display
227: - vf - viewer and format for monitoring
229: Options Database Key:
230: . -eps_monitor_all - activates `EPSMonitorAll()`
232: Note:
233: This is not called directly by users, rather one calls `EPSMonitorSet()`, with this
234: function as an argument, to cause the monitor to be used during the `EPS` solve.
236: Level: intermediate
238: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorFirst()`, `EPSMonitorConverged()`
239: @*/
240: PetscErrorCode EPSMonitorAll(EPS eps,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
241: {
242: PetscInt i;
243: PetscScalar er,ei;
244: PetscViewer viewer = vf->viewer;
246: PetscFunctionBegin;
249: PetscCall(PetscViewerPushFormat(viewer,vf->format));
250: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
251: if (its==1 && ((PetscObject)eps)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)eps)->prefix));
252: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " EPS nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
253: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
254: for (i=0;i<nest;i++) {
255: er = eigr[i]; ei = eigi[i];
256: PetscCall(STBackTransform(eps->st,1,&er,&ei));
257: PetscCall(EPSMonitorPrintEval(eps,viewer,er,ei));
258: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
259: }
260: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
261: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
262: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
263: PetscCall(PetscViewerPopFormat(viewer));
264: PetscFunctionReturn(PETSC_SUCCESS);
265: }
267: /*@C
268: EPSMonitorConverged - Print the approximate values and
269: error estimates as they converge.
271: Collective
273: Input Parameters:
274: + eps - the linear eigensolver context
275: . its - iteration number
276: . nconv - number of converged eigenpairs so far
277: . eigr - real part of the eigenvalues
278: . eigi - imaginary part of the eigenvalues
279: . errest - error estimates
280: . nest - number of error estimates to display
281: - vf - viewer and format for monitoring
283: Options Database Key:
284: . -eps_monitor_conv - activates `EPSMonitorConverged()`
286: Notes:
287: This is not called directly by users, rather one calls `EPSMonitorSet()`, with this
288: function as an argument, to cause the monitor to be used during the `EPS` solve.
290: Call `EPSMonitorConvergedCreate()` to create the context used with this monitor.
292: Level: intermediate
294: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorConvergedCreate()`, `EPSMonitorFirst()`, `EPSMonitorAll()`
295: @*/
296: PetscErrorCode EPSMonitorConverged(EPS eps,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
297: {
298: PetscInt i,*oldnconv;
299: PetscScalar er,ei;
300: PetscViewer viewer = vf->viewer;
302: PetscFunctionBegin;
305: oldnconv = (PetscInt*)vf->data;
306: if (its==1 && ((PetscObject)eps)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)eps)->prefix));
307: if (its==1) *oldnconv = 0;
308: if (*oldnconv!=nconv) {
309: PetscCall(PetscViewerPushFormat(viewer,vf->format));
310: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
311: for (i=*oldnconv;i<nconv;i++) {
312: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " EPS converged value (error) #%" PetscInt_FMT,its,i));
313: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
314: er = eigr[i]; ei = eigi[i];
315: PetscCall(STBackTransform(eps->st,1,&er,&ei));
316: PetscCall(EPSMonitorPrintEval(eps,viewer,er,ei));
317: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
318: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
319: }
320: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
321: PetscCall(PetscViewerPopFormat(viewer));
322: *oldnconv = nconv;
323: }
324: PetscFunctionReturn(PETSC_SUCCESS);
325: }
327: /*@C
328: EPSMonitorConvergedCreate - Creates the context for the convergence history monitor.
330: Collective
332: Input Parameters:
333: + viewer - the viewer
334: . format - the viewer format
335: - ctx - an optional user context
337: Output Parameter:
338: . vf - the viewer and format context
340: Level: intermediate
342: .seealso: [](ch:eps), `EPSMonitorSet()`
343: @*/
344: PetscErrorCode EPSMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,PetscCtx ctx,PetscViewerAndFormat **vf)
345: {
346: PetscInt *oldnconv;
348: PetscFunctionBegin;
349: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
350: PetscCall(PetscNew(&oldnconv));
351: (*vf)->data = (void*)oldnconv;
352: (*vf)->data_destroy = PetscCtxDestroyDefault;
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: /*@C
357: EPSMonitorFirstDrawLG - Plots the error estimate of the first unconverged
358: approximation at each iteration of the eigensolver.
360: Collective
362: Input Parameters:
363: + eps - the linear eigensolver context
364: . its - iteration number
365: . nconv - number of converged eigenpairs so far
366: . eigr - real part of the eigenvalues
367: . eigi - imaginary part of the eigenvalues
368: . errest - error estimates
369: . nest - number of error estimates to display
370: - vf - viewer and format for monitoring
372: Options Database Key:
373: . -eps_monitor draw::draw_lg - activates `EPSMonitorFirstDrawLG()`
375: Notes:
376: This is not called directly by users, rather one calls `EPSMonitorSet()`, with this
377: function as an argument, to cause the monitor to be used during the `EPS` solve.
379: Call `EPSMonitorFirstDrawLGCreate()` to create the context used with this monitor.
381: Level: intermediate
383: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorFirstDrawLGCreate()`
384: @*/
385: PetscErrorCode EPSMonitorFirstDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
386: {
387: PetscViewer viewer = vf->viewer;
388: PetscDrawLG lg;
389: PetscReal x,y;
391: PetscFunctionBegin;
394: PetscCall(PetscViewerPushFormat(viewer,vf->format));
395: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
396: if (its==1) {
397: PetscCall(PetscDrawLGReset(lg));
398: PetscCall(PetscDrawLGSetDimension(lg,1));
399: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0));
400: }
401: if (nconv<nest) {
402: x = (PetscReal)its;
403: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
404: else y = 0.0;
405: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
406: if (its <= 20 || !(its % 5) || eps->reason) {
407: PetscCall(PetscDrawLGDraw(lg));
408: PetscCall(PetscDrawLGSave(lg));
409: }
410: }
411: PetscCall(PetscViewerPopFormat(viewer));
412: PetscFunctionReturn(PETSC_SUCCESS);
413: }
415: /*@C
416: EPSMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
418: Collective
420: Input Parameters:
421: + viewer - the viewer
422: . format - the viewer format
423: - ctx - an optional user context
425: Output Parameter:
426: . vf - the viewer and format context
428: Level: intermediate
430: .seealso: [](ch:eps), `EPSMonitorSet()`
431: @*/
432: PetscErrorCode EPSMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
433: {
434: PetscFunctionBegin;
435: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
436: (*vf)->data = ctx;
437: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
438: PetscFunctionReturn(PETSC_SUCCESS);
439: }
441: /*@C
442: EPSMonitorAllDrawLG - Plots the error estimate of all unconverged
443: approximations at each iteration of the eigensolver.
445: Collective
447: Input Parameters:
448: + eps - the linear eigensolver context
449: . its - iteration number
450: . nconv - number of converged eigenpairs so far
451: . eigr - real part of the eigenvalues
452: . eigi - imaginary part of the eigenvalues
453: . errest - error estimates
454: . nest - number of error estimates to display
455: - vf - viewer and format for monitoring
457: Options Database Key:
458: . -eps_monitor_all draw::draw_lg - activates `EPSMonitorAllDrawLG()`
460: Notes:
461: This is not called directly by users, rather one calls `EPSMonitorSet()`, with this
462: function as an argument, to cause the monitor to be used during the `EPS` solve.
464: Call `EPSMonitorAllDrawLGCreate()` to create the context used with this monitor.
466: Level: intermediate
468: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorAllDrawLGCreate()`
469: @*/
470: PetscErrorCode EPSMonitorAllDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
471: {
472: PetscViewer viewer = vf->viewer;
473: PetscDrawLG lg;
474: PetscInt i,n = PetscMin(eps->nev,255);
475: PetscReal *x,*y;
477: PetscFunctionBegin;
480: PetscCall(PetscViewerPushFormat(viewer,vf->format));
481: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
482: if (its==1) {
483: PetscCall(PetscDrawLGReset(lg));
484: PetscCall(PetscDrawLGSetDimension(lg,n));
485: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(eps->tol)-2,0.0));
486: }
487: PetscCall(PetscMalloc2(n,&x,n,&y));
488: for (i=0;i<n;i++) {
489: x[i] = (PetscReal)its;
490: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
491: else y[i] = 0.0;
492: }
493: PetscCall(PetscDrawLGAddPoint(lg,x,y));
494: if (its <= 20 || !(its % 5) || eps->reason) {
495: PetscCall(PetscDrawLGDraw(lg));
496: PetscCall(PetscDrawLGSave(lg));
497: }
498: PetscCall(PetscFree2(x,y));
499: PetscCall(PetscViewerPopFormat(viewer));
500: PetscFunctionReturn(PETSC_SUCCESS);
501: }
503: /*@C
504: EPSMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
506: Collective
508: Input Parameters:
509: + viewer - the viewer
510: . format - the viewer format
511: - ctx - an optional user context
513: Output Parameter:
514: . vf - the viewer and format context
516: Level: intermediate
518: .seealso: [](ch:eps), `EPSMonitorSet()`
519: @*/
520: PetscErrorCode EPSMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
521: {
522: PetscFunctionBegin;
523: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
524: (*vf)->data = ctx;
525: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
526: PetscFunctionReturn(PETSC_SUCCESS);
527: }
529: /*@C
530: EPSMonitorConvergedDrawLG - Plots the number of converged eigenvalues
531: at each iteration of the eigensolver.
533: Collective
535: Input Parameters:
536: + eps - the linear eigensolver context
537: . its - iteration number
538: . nconv - number of converged eigenpairs so far
539: . eigr - real part of the eigenvalues
540: . eigi - imaginary part of the eigenvalues
541: . errest - error estimates
542: . nest - number of error estimates to display
543: - vf - viewer and format for monitoring
545: Options Database Key:
546: . -eps_monitor_conv draw::draw_lg - activates `EPSMonitorConvergedDrawLG()`
548: Notes:
549: This is not called directly by users, rather one calls `EPSMonitorSet()`, with this
550: function as an argument, to cause the monitor to be used during the `EPS` solve.
552: Call `EPSMonitorConvergedDrawLGCreate()` to create the context used with this monitor.
554: Level: intermediate
556: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorConvergedDrawLGCreate()`
557: @*/
558: PetscErrorCode EPSMonitorConvergedDrawLG(EPS eps,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
559: {
560: PetscViewer viewer = vf->viewer;
561: PetscDrawLG lg;
562: PetscReal x,y;
564: PetscFunctionBegin;
567: PetscCall(PetscViewerPushFormat(viewer,vf->format));
568: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
569: if (its==1) {
570: PetscCall(PetscDrawLGReset(lg));
571: PetscCall(PetscDrawLGSetDimension(lg,1));
572: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,eps->nev));
573: }
574: x = (PetscReal)its;
575: y = (PetscReal)eps->nconv;
576: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
577: if (its <= 20 || !(its % 5) || eps->reason) {
578: PetscCall(PetscDrawLGDraw(lg));
579: PetscCall(PetscDrawLGSave(lg));
580: }
581: PetscCall(PetscViewerPopFormat(viewer));
582: PetscFunctionReturn(PETSC_SUCCESS);
583: }
585: /*@C
586: EPSMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
588: Collective
590: Input Parameters:
591: + viewer - the viewer
592: . format - the viewer format
593: - ctx - an optional user context
595: Output Parameter:
596: . vf - the viewer and format context
598: Level: intermediate
600: .seealso: [](ch:eps), `EPSMonitorSet()`
601: @*/
602: PetscErrorCode EPSMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
603: {
604: PetscInt *oldnconv;
606: PetscFunctionBegin;
607: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
608: PetscCall(PetscNew(&oldnconv));
609: (*vf)->data = (void*)oldnconv;
610: (*vf)->data_destroy = PetscCtxDestroyDefault;
611: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
612: PetscFunctionReturn(PETSC_SUCCESS);
613: }