Actual source code: nepmon.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: NEP routines related to monitors
12: */
14: #include <slepc/private/nepimpl.h>
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
21: {
22: PetscInt i,n = nep->numbermonitors;
24: PetscFunctionBegin;
25: for (i=0;i<n;i++) PetscCall((*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@C
30: NEPMonitorSet - 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: + nep - the nonlinear eigensolver context
37: . monitor - pointer to function (if this is `NULL`, it turns off monitoring),
38: see `NEPMonitorFn`
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: + -nep_monitor - print only the first error estimate
46: . -nep_monitor_all - print error estimates at each iteration
47: . -nep_monitor_conv - print the eigenvalue approximations only when
48: convergence has been reached
49: . -nep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
50: approximate eigenvalue
51: . -nep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
52: approximate eigenvalues
53: . -nep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
54: - -nep_monitor_cancel - cancels all monitors that have been hardwired into
55: a code by calls to `NEPMonitorSet()`, but does not cancel
56: those set via the options database.
58: Notes:
59: The options database option `-nep_monitor` and related options are the easiest way
60: to turn on `NEP` iteration monitoring.
62: `NEPMonitorRegister()` provides a way to associate an options database key with `NEP`
63: monitor function.
65: Several different monitoring routines may be set by calling `NEPMonitorSet()` 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 `NEP` object.
71: Level: intermediate
73: .seealso: [](ch:nep), `NEPMonitorFirst()`, `NEPMonitorAll()`, `NEPMonitorConverged()`, `NEPMonitorFirstDrawLG()`, `NEPMonitorAllDrawLG()`, `NEPMonitorConvergedDrawLG()`, `NEPMonitorCancel()`
74: @*/
75: PetscErrorCode NEPMonitorSet(NEP nep,NEPMonitorFn *monitor,PetscCtx ctx,PetscCtxDestroyFn *monitordestroy)
76: {
77: PetscInt i;
78: PetscBool identical;
80: PetscFunctionBegin;
82: for (i=0;i<nep->numbermonitors;i++) {
83: PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))(PetscVoidFn*)monitor,ctx,monitordestroy,(PetscErrorCode (*)(void))(PetscVoidFn*)nep->monitor[i],nep->monitorcontext[i],nep->monitordestroy[i],&identical));
84: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
85: }
86: PetscCheck(nep->numbermonitors<MAXNEPMONITORS,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
87: nep->monitor[nep->numbermonitors] = monitor;
88: nep->monitorcontext[nep->numbermonitors] = ctx;
89: nep->monitordestroy[nep->numbermonitors++] = monitordestroy;
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: /*@
94: NEPMonitorCancel - Clears all monitors for a `NEP` object.
96: Logically Collective
98: Input Parameter:
99: . nep - the nonlinear eigensolver context
101: Options Database Key:
102: . -nep_monitor_cancel - cancels all monitors that have been hardwired into a code by calls to
103: `NEPMonitorSet()`, but does not cancel those set via the options database.
105: Level: intermediate
107: .seealso: [](ch:nep), `NEPMonitorSet()`
108: @*/
109: PetscErrorCode NEPMonitorCancel(NEP nep)
110: {
111: PetscInt i;
113: PetscFunctionBegin;
115: for (i=0; i<nep->numbermonitors; i++) {
116: if (nep->monitordestroy[i]) PetscCall((*nep->monitordestroy[i])(&nep->monitorcontext[i]));
117: }
118: nep->numbermonitors = 0;
119: PetscFunctionReturn(PETSC_SUCCESS);
120: }
122: /*@C
123: NEPGetMonitorContext - Gets the monitor context, as set by
124: `NEPMonitorSet()` for the FIRST monitor only.
126: Not Collective
128: Input Parameter:
129: . nep - the nonlinear eigensolver context
131: Output Parameter:
132: . ctx - monitor context
134: Level: intermediate
136: .seealso: [](ch:nep), `NEPMonitorSet()`
137: @*/
138: PetscErrorCode NEPGetMonitorContext(NEP nep,PetscCtxRt ctx)
139: {
140: PetscFunctionBegin;
142: *(void**)ctx = nep->monitorcontext[0];
143: PetscFunctionReturn(PETSC_SUCCESS);
144: }
146: /*@C
147: NEPMonitorFirst - Print the first unconverged approximate value and
148: error estimate at each iteration of the nonlinear eigensolver.
150: Collective
152: Input Parameters:
153: + nep - the nonlinear eigensolver context
154: . its - iteration number
155: . nconv - number of converged eigenpairs so far
156: . eigr - real part of the eigenvalues
157: . eigi - imaginary part of the eigenvalues
158: . errest - error estimates
159: . nest - number of error estimates to display
160: - vf - viewer and format for monitoring
162: Options Database Key:
163: . -nep_monitor - activates `NEPMonitorFirst()`
165: Note:
166: This is not called directly by users, rather one calls `NEPMonitorSet()`, with this
167: function as an argument, to cause the monitor to be used during the `NEP` solve.
169: Level: intermediate
171: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorAll()`, `NEPMonitorConverged()`
172: @*/
173: PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
174: {
175: PetscViewer viewer = vf->viewer;
177: PetscFunctionBegin;
180: if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix));
181: if (nconv<nest) {
182: PetscCall(PetscViewerPushFormat(viewer,vf->format));
183: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
184: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
185: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
186: #if defined(PETSC_USE_COMPLEX)
187: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv])));
188: #else
189: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]));
190: if (eigi[nconv]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]));
191: #endif
192: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
193: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
194: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
195: PetscCall(PetscViewerPopFormat(viewer));
196: }
197: PetscFunctionReturn(PETSC_SUCCESS);
198: }
200: /*@C
201: NEPMonitorAll - Print the current approximate values and
202: error estimates at each iteration of the nonlinear eigensolver.
204: Collective
206: Input Parameters:
207: + nep - the nonlinear eigensolver context
208: . its - iteration number
209: . nconv - number of converged eigenpairs so far
210: . eigr - real part of the eigenvalues
211: . eigi - imaginary part of the eigenvalues
212: . errest - error estimates
213: . nest - number of error estimates to display
214: - vf - viewer and format for monitoring
216: Options Database Key:
217: . -nep_monitor_all - activates `NEPMonitorAll()`
219: Note:
220: This is not called directly by users, rather one calls `NEPMonitorSet()`, with this
221: function as an argument, to cause the monitor to be used during the `NEP` solve.
223: Level: intermediate
225: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorFirst()`, `NEPMonitorConverged()`
226: @*/
227: PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
228: {
229: PetscInt i;
230: PetscViewer viewer = vf->viewer;
232: PetscFunctionBegin;
235: PetscCall(PetscViewerPushFormat(viewer,vf->format));
236: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
237: if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix));
238: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
239: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
240: for (i=0;i<nest;i++) {
241: #if defined(PETSC_USE_COMPLEX)
242: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i])));
243: #else
244: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]));
245: if (eigi[i]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]));
246: #endif
247: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
248: }
249: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
250: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
251: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
252: PetscCall(PetscViewerPopFormat(viewer));
253: PetscFunctionReturn(PETSC_SUCCESS);
254: }
256: /*@C
257: NEPMonitorConverged - Print the approximate values and
258: error estimates as they converge.
260: Collective
262: Input Parameters:
263: + nep - the nonlinear eigensolver context
264: . its - iteration number
265: . nconv - number of converged eigenpairs so far
266: . eigr - real part of the eigenvalues
267: . eigi - imaginary part of the eigenvalues
268: . errest - error estimates
269: . nest - number of error estimates to display
270: - vf - viewer and format for monitoring
272: Options Database Key:
273: . -nep_monitor_conv - activates `NEPMonitorConverged()`
275: Notes:
276: This is not called directly by users, rather one calls `NEPMonitorSet()`, with this
277: function as an argument, to cause the monitor to be used during the `NEP` solve.
279: Call `NEPMonitorConvergedCreate()` to create the context used with this monitor.
281: Level: intermediate
283: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorConvergedCreate()`, `NEPMonitorFirst()`, `NEPMonitorAll()`
284: @*/
285: PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
286: {
287: PetscInt i,*oldnconv;
288: PetscViewer viewer = vf->viewer;
290: PetscFunctionBegin;
293: oldnconv = (PetscInt*)vf->data;
294: if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)nep)->prefix));
295: if (its==1) *oldnconv = 0;
296: if (*oldnconv!=nconv) {
297: PetscCall(PetscViewerPushFormat(viewer,vf->format));
298: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
299: for (i=*oldnconv;i<nconv;i++) {
300: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP converged value (error) #%" PetscInt_FMT,its,i));
301: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
302: #if defined(PETSC_USE_COMPLEX)
303: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i])));
304: #else
305: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]));
306: if (eigi[i]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]));
307: #endif
308: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
309: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
310: }
311: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
312: PetscCall(PetscViewerPopFormat(viewer));
313: *oldnconv = nconv;
314: }
315: PetscFunctionReturn(PETSC_SUCCESS);
316: }
318: /*@C
319: NEPMonitorConvergedCreate - Creates the context for the convergence history monitor.
321: Collective
323: Input Parameters:
324: + viewer - the viewer
325: . format - the viewer format
326: - ctx - an optional user context
328: Output Parameter:
329: . vf - the viewer and format context
331: Level: intermediate
333: .seealso: [](ch:nep), `NEPMonitorSet()`
334: @*/
335: PetscErrorCode NEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,PetscCtx ctx,PetscViewerAndFormat **vf)
336: {
337: PetscInt *oldnconv;
339: PetscFunctionBegin;
340: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
341: PetscCall(PetscNew(&oldnconv));
342: (*vf)->data = (void*)oldnconv;
343: (*vf)->data_destroy = PetscCtxDestroyDefault;
344: PetscFunctionReturn(PETSC_SUCCESS);
345: }
347: /*@C
348: NEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
349: approximation at each iteration of the nonlinear eigensolver.
351: Collective
353: Input Parameters:
354: + nep - the nonlinear eigensolver context
355: . its - iteration number
356: . nconv - number of converged eigenpairs so far
357: . eigr - real part of the eigenvalues
358: . eigi - imaginary part of the eigenvalues
359: . errest - error estimates
360: . nest - number of error estimates to display
361: - vf - viewer and format for monitoring
363: Options Database Key:
364: . -nep_monitor draw::draw_lg - activates `NEPMonitorFirstDrawLG()`
366: Notes:
367: This is not called directly by users, rather one calls `NEPMonitorSet()`, with this
368: function as an argument, to cause the monitor to be used during the `NEP` solve.
370: Call `NEPMonitorFirstDrawLGCreate()` to create the context used with this monitor.
372: Level: intermediate
374: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorFirstDrawLGCreate()`
375: @*/
376: PetscErrorCode NEPMonitorFirstDrawLG(NEP nep,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(nep->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) || nep->reason) {
398: PetscCall(PetscDrawLGDraw(lg));
399: PetscCall(PetscDrawLGSave(lg));
400: }
401: }
402: PetscCall(PetscViewerPopFormat(viewer));
403: PetscFunctionReturn(PETSC_SUCCESS);
404: }
406: /*@C
407: NEPMonitorFirstDrawLGCreate - 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: [](ch:nep), `NEPMonitorSet()`
422: @*/
423: PetscErrorCode NEPMonitorFirstDrawLGCreate(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: NEPMonitorAllDrawLG - Plots the error estimate of all unconverged
434: approximations at each iteration of the nonlinear eigensolver.
436: Collective
438: Input Parameters:
439: + nep - the nonlinear 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: . -nep_monitor_all draw::draw_lg - activates `NEPMonitorAllDrawLG()`
451: Notes:
452: This is not called directly by users, rather one calls `NEPMonitorSet()`, with this
453: function as an argument, to cause the monitor to be used during the `NEP` solve.
455: Call `NEPMonitorAllDrawLGCreate()` to create the context used with this monitor.
457: Level: intermediate
459: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorAllDrawLGCreate()`
460: @*/
461: PetscErrorCode NEPMonitorAllDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
462: {
463: PetscViewer viewer = vf->viewer;
464: PetscDrawLG lg;
465: PetscInt i,n = PetscMin(nep->nev,255);
466: PetscReal *x,*y;
468: PetscFunctionBegin;
471: PetscCall(PetscViewerPushFormat(viewer,vf->format));
472: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
473: if (its==1) {
474: PetscCall(PetscDrawLGReset(lg));
475: PetscCall(PetscDrawLGSetDimension(lg,n));
476: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0));
477: }
478: PetscCall(PetscMalloc2(n,&x,n,&y));
479: for (i=0;i<n;i++) {
480: x[i] = (PetscReal)its;
481: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
482: else y[i] = 0.0;
483: }
484: PetscCall(PetscDrawLGAddPoint(lg,x,y));
485: if (its <= 20 || !(its % 5) || nep->reason) {
486: PetscCall(PetscDrawLGDraw(lg));
487: PetscCall(PetscDrawLGSave(lg));
488: }
489: PetscCall(PetscFree2(x,y));
490: PetscCall(PetscViewerPopFormat(viewer));
491: PetscFunctionReturn(PETSC_SUCCESS);
492: }
494: /*@C
495: NEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
497: Collective
499: Input Parameters:
500: + viewer - the viewer
501: . format - the viewer format
502: - ctx - an optional user context
504: Output Parameter:
505: . vf - the viewer and format context
507: Level: intermediate
509: .seealso: [](ch:nep), `NEPMonitorSet()`
510: @*/
511: PetscErrorCode NEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
512: {
513: PetscFunctionBegin;
514: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
515: (*vf)->data = ctx;
516: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
517: PetscFunctionReturn(PETSC_SUCCESS);
518: }
520: /*@C
521: NEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
522: at each iteration of the nonlinear eigensolver.
524: Collective
526: Input Parameters:
527: + nep - the nonlinear eigensolver context
528: . its - iteration number
529: . nconv - number of converged eigenpairs so far
530: . eigr - real part of the eigenvalues
531: . eigi - imaginary part of the eigenvalues
532: . errest - error estimates
533: . nest - number of error estimates to display
534: - vf - viewer and format for monitoring
536: Options Database Key:
537: . -nep_monitor_conv draw::draw_lg - activates `NEPMonitorConvergedDrawLG()`
539: Notes:
540: This is not called directly by users, rather one calls `NEPMonitorSet()`, with this
541: function as an argument, to cause the monitor to be used during the `NEP` solve.
543: Call `NEPMonitorConvergedDrawLGCreate()` to create the context used with this monitor.
545: Level: intermediate
547: .seealso: [](ch:nep), `NEPMonitorSet()`, `NEPMonitorConvergedDrawLGCreate()`
548: @*/
549: PetscErrorCode NEPMonitorConvergedDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
550: {
551: PetscViewer viewer = vf->viewer;
552: PetscDrawLG lg;
553: PetscReal x,y;
555: PetscFunctionBegin;
558: PetscCall(PetscViewerPushFormat(viewer,vf->format));
559: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
560: if (its==1) {
561: PetscCall(PetscDrawLGReset(lg));
562: PetscCall(PetscDrawLGSetDimension(lg,1));
563: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,nep->nev));
564: }
565: x = (PetscReal)its;
566: y = (PetscReal)nep->nconv;
567: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
568: if (its <= 20 || !(its % 5) || nep->reason) {
569: PetscCall(PetscDrawLGDraw(lg));
570: PetscCall(PetscDrawLGSave(lg));
571: }
572: PetscCall(PetscViewerPopFormat(viewer));
573: PetscFunctionReturn(PETSC_SUCCESS);
574: }
576: /*@C
577: NEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
579: Collective
581: Input Parameters:
582: + viewer - the viewer
583: . format - the viewer format
584: - ctx - an optional user context
586: Output Parameter:
587: . vf - the viewer and format context
589: Level: intermediate
591: .seealso: [](ch:nep), `NEPMonitorSet()`
592: @*/
593: PetscErrorCode NEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
594: {
595: PetscInt *oldnconv;
597: PetscFunctionBegin;
598: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
599: PetscCall(PetscNew(&oldnconv));
600: (*vf)->data = (void*)oldnconv;
601: (*vf)->data_destroy = PetscCtxDestroyDefault;
602: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
603: PetscFunctionReturn(PETSC_SUCCESS);
604: }