Actual source code: pepmon.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: 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 - the polynomial eigensolver context
37: . monitor - pointer to function (if this is `NULL`, it turns off monitoring),
38: see `PEPMonitorFn`
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: + -pep_monitor - print only the first error estimate
46: . -pep_monitor_all - print error estimates at each iteration
47: . -pep_monitor_conv - print the eigenvalue approximations only when
48: convergence has been reached
49: . -pep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
50: approximate eigenvalue
51: . -pep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
52: approximate eigenvalues
53: . -pep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
54: - -pep_monitor_cancel - cancels all monitors that have been hardwired into
55: a code by calls to `PEPMonitorSet()`, but does not cancel
56: those set via the options database.
58: Notes:
59: The options database option `-pep_monitor` and related options are the easiest way
60: to turn on `PEP` iteration monitoring.
62: `PEPMonitorRegister()` provides a way to associate an options database key with `PEP`
63: monitor function.
65: Several different monitoring routines may be set by calling `PEPMonitorSet()` 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 `PEP` object.
71: Level: intermediate
73: .seealso: [](ch:pep), `PEPMonitorFirst()`, `PEPMonitorAll()`, `PEPMonitorConverged()`, `PEPMonitorFirstDrawLG()`, `PEPMonitorAllDrawLG()`, `PEPMonitorConvergedDrawLG()`, `PEPMonitorCancel()`
74: @*/
75: PetscErrorCode PEPMonitorSet(PEP pep,PEPMonitorFn *monitor,PetscCtx ctx,PetscCtxDestroyFn *monitordestroy)
76: {
77: PetscInt i;
78: PetscBool identical;
80: PetscFunctionBegin;
82: for (i=0;i<pep->numbermonitors;i++) {
83: PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))(PetscVoidFn*)monitor,ctx,monitordestroy,(PetscErrorCode (*)(void))(PetscVoidFn*)pep->monitor[i],pep->monitorcontext[i],pep->monitordestroy[i],&identical));
84: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
85: }
86: PetscCheck(pep->numbermonitors<MAXPEPMONITORS,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Too many PEP monitors set");
87: pep->monitor[pep->numbermonitors] = monitor;
88: pep->monitorcontext[pep->numbermonitors] = ctx;
89: pep->monitordestroy[pep->numbermonitors++] = monitordestroy;
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: /*@
94: PEPMonitorCancel - Clears all monitors for a `PEP` object.
96: Logically Collective
98: Input Parameter:
99: . pep - the polynomial eigensolver context
101: Options Database Key:
102: . -pep_monitor_cancel - cancels all monitors that have been hardwired into a code by calls to
103: `PEPMonitorSet()`, but does not cancel those set via the options database.
105: Level: intermediate
107: .seealso: [](ch:pep), `PEPMonitorSet()`
108: @*/
109: PetscErrorCode PEPMonitorCancel(PEP pep)
110: {
111: PetscInt i;
113: PetscFunctionBegin;
115: for (i=0; i<pep->numbermonitors; i++) {
116: if (pep->monitordestroy[i]) PetscCall((*pep->monitordestroy[i])(&pep->monitorcontext[i]));
117: }
118: pep->numbermonitors = 0;
119: PetscFunctionReturn(PETSC_SUCCESS);
120: }
122: /*@C
123: PEPGetMonitorContext - Gets the monitor context, as set by
124: `PEPMonitorSet()` for the FIRST monitor only.
126: Not Collective
128: Input Parameter:
129: . pep - the polynomial eigensolver context
131: Output Parameter:
132: . ctx - monitor context
134: Level: intermediate
136: .seealso: [](ch:pep), `PEPMonitorSet()`
137: @*/
138: PetscErrorCode PEPGetMonitorContext(PEP pep,PetscCtxRt ctx)
139: {
140: PetscFunctionBegin;
142: *(void**)ctx = pep->monitorcontext[0];
143: PetscFunctionReturn(PETSC_SUCCESS);
144: }
146: /*
147: Helper function to compute eigenvalue that must be viewed in monitor
148: */
149: static PetscErrorCode PEPMonitorGetTrueEig(PEP pep,PetscScalar *er,PetscScalar *ei)
150: {
151: PetscBool flg,sinv;
153: PetscFunctionBegin;
154: PetscCall(STBackTransform(pep->st,1,er,ei));
155: if (pep->sfactor!=1.0) {
156: PetscCall(STGetTransform(pep->st,&flg));
157: PetscCall(PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv));
158: if (flg && sinv) {
159: *er /= pep->sfactor; *ei /= pep->sfactor;
160: } else {
161: *er *= pep->sfactor; *ei *= pep->sfactor;
162: }
163: }
164: PetscFunctionReturn(PETSC_SUCCESS);
165: }
167: /*@C
168: PEPMonitorFirst - Print the first unconverged approximate value and
169: error estimate at each iteration of the polynomial eigensolver.
171: Collective
173: Input Parameters:
174: + pep - the polynomial eigensolver context
175: . its - iteration number
176: . nconv - number of converged eigenpairs so far
177: . eigr - real part of the eigenvalues
178: . eigi - imaginary part of the eigenvalues
179: . errest - error estimates
180: . nest - number of error estimates to display
181: - vf - viewer and format for monitoring
183: Options Database Key:
184: . -pep_monitor - activates `PEPMonitorFirst()`
186: Note:
187: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
188: function as an argument, to cause the monitor to be used during the `PEP` solve.
190: Level: intermediate
192: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorAll()`, `PEPMonitorConverged()`
193: @*/
194: PetscErrorCode PEPMonitorFirst(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
195: {
196: PetscScalar er,ei;
197: PetscViewer viewer = vf->viewer;
199: PetscFunctionBegin;
202: if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix));
203: if (nconv<nest) {
204: PetscCall(PetscViewerPushFormat(viewer,vf->format));
205: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
206: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
207: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
208: er = eigr[nconv]; ei = eigi[nconv];
209: PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
210: #if defined(PETSC_USE_COMPLEX)
211: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
212: #else
213: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
214: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
215: #endif
216: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
217: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
218: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
219: PetscCall(PetscViewerPopFormat(viewer));
220: }
221: PetscFunctionReturn(PETSC_SUCCESS);
222: }
224: /*@C
225: PEPMonitorAll - Print the current approximate values and
226: error estimates at each iteration of the polynomial eigensolver.
228: Collective
230: Input Parameters:
231: + pep - the polynomial eigensolver context
232: . its - iteration number
233: . nconv - number of converged eigenpairs so far
234: . eigr - real part of the eigenvalues
235: . eigi - imaginary part of the eigenvalues
236: . errest - error estimates
237: . nest - number of error estimates to display
238: - vf - viewer and format for monitoring
240: Options Database Key:
241: . -pep_monitor_all - activates `PEPMonitorAll()`
243: Note:
244: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
245: function as an argument, to cause the monitor to be used during the `PEP` solve.
247: Level: intermediate
249: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorFirst()`, `PEPMonitorConverged()`
250: @*/
251: PetscErrorCode PEPMonitorAll(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
252: {
253: PetscInt i;
254: PetscScalar er,ei;
255: PetscViewer viewer = vf->viewer;
257: PetscFunctionBegin;
260: PetscCall(PetscViewerPushFormat(viewer,vf->format));
261: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
262: if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)pep)->prefix));
263: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
264: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
265: for (i=0;i<nest;i++) {
266: er = eigr[i]; ei = eigi[i];
267: PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
268: #if defined(PETSC_USE_COMPLEX)
269: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
270: #else
271: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
272: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
273: #endif
274: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
275: }
276: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
277: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
278: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
279: PetscCall(PetscViewerPopFormat(viewer));
280: PetscFunctionReturn(PETSC_SUCCESS);
281: }
283: /*@C
284: PEPMonitorConverged - Print the approximate values and
285: error estimates as they converge.
287: Collective
289: Input Parameters:
290: + pep - the polynomial eigensolver context
291: . its - iteration number
292: . nconv - number of converged eigenpairs so far
293: . eigr - real part of the eigenvalues
294: . eigi - imaginary part of the eigenvalues
295: . errest - error estimates
296: . nest - number of error estimates to display
297: - vf - viewer and format for monitoring
299: Options Database Key:
300: . -pep_monitor_conv - activates `PEPMonitorConverged()`
302: Notes:
303: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
304: function as an argument, to cause the monitor to be used during the `PEP` solve.
306: Call `PEPMonitorConvergedCreate()` to create the context used with this monitor.
308: Level: intermediate
310: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorConvergedCreate()`, `PEPMonitorFirst()`, `PEPMonitorAll()`
311: @*/
312: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
313: {
314: PetscInt i,*oldnconv;
315: PetscScalar er,ei;
316: PetscViewer viewer = vf->viewer;
318: PetscFunctionBegin;
321: oldnconv = (PetscInt*)vf->data;
322: if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)pep)->prefix));
323: if (its==1) *oldnconv = 0;
324: if (*oldnconv!=nconv) {
325: PetscCall(PetscViewerPushFormat(viewer,vf->format));
326: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
327: for (i=*oldnconv;i<nconv;i++) {
328: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP converged value (error) #%" PetscInt_FMT,its,i));
329: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
330: er = eigr[i]; ei = eigi[i];
331: PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
332: #if defined(PETSC_USE_COMPLEX)
333: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
334: #else
335: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
336: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
337: #endif
338: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
339: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
340: }
341: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
342: PetscCall(PetscViewerPopFormat(viewer));
343: *oldnconv = nconv;
344: }
345: PetscFunctionReturn(PETSC_SUCCESS);
346: }
348: /*@C
349: PEPMonitorConvergedCreate - Creates the context for the convergence history monitor.
351: Collective
353: Input Parameters:
354: + viewer - the viewer
355: . format - the viewer format
356: - ctx - an optional user context
358: Output Parameter:
359: . vf - the viewer and format context
361: Level: intermediate
363: .seealso: [](ch:pep), `PEPMonitorSet()`
364: @*/
365: PetscErrorCode PEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,PetscCtx ctx,PetscViewerAndFormat **vf)
366: {
367: PetscInt *oldnconv;
369: PetscFunctionBegin;
370: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
371: PetscCall(PetscNew(&oldnconv));
372: (*vf)->data = (void*)oldnconv;
373: (*vf)->data_destroy = PetscCtxDestroyDefault;
374: PetscFunctionReturn(PETSC_SUCCESS);
375: }
377: /*@C
378: PEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
379: approximation at each iteration of the polynomial eigensolver.
381: Collective
383: Input Parameters:
384: + pep - the polynomial eigensolver context
385: . its - iteration number
386: . nconv - number of converged eigenpairs so far
387: . eigr - real part of the eigenvalues
388: . eigi - imaginary part of the eigenvalues
389: . errest - error estimates
390: . nest - number of error estimates to display
391: - vf - viewer and format for monitoring
393: Options Database Key:
394: . -pep_monitor draw::draw_lg - activates `PEPMonitorFirstDrawLG()`
396: Notes:
397: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
398: function as an argument, to cause the monitor to be used during the `PEP` solve.
400: Call `PEPMonitorFirstDrawLGCreate()` to create the context used with this monitor.
402: Level: intermediate
404: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorFirstDrawLGCreate()`
405: @*/
406: PetscErrorCode PEPMonitorFirstDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
407: {
408: PetscViewer viewer = vf->viewer;
409: PetscDrawLG lg;
410: PetscReal x,y;
412: PetscFunctionBegin;
415: PetscCall(PetscViewerPushFormat(viewer,vf->format));
416: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
417: if (its==1) {
418: PetscCall(PetscDrawLGReset(lg));
419: PetscCall(PetscDrawLGSetDimension(lg,1));
420: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0));
421: }
422: if (nconv<nest) {
423: x = (PetscReal)its;
424: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
425: else y = 0.0;
426: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
427: if (its <= 20 || !(its % 5) || pep->reason) {
428: PetscCall(PetscDrawLGDraw(lg));
429: PetscCall(PetscDrawLGSave(lg));
430: }
431: }
432: PetscCall(PetscViewerPopFormat(viewer));
433: PetscFunctionReturn(PETSC_SUCCESS);
434: }
436: /*@C
437: PEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
439: Collective
441: Input Parameters:
442: + viewer - the viewer
443: . format - the viewer format
444: - ctx - an optional user context
446: Output Parameter:
447: . vf - the viewer and format context
449: Level: intermediate
451: .seealso: [](ch:pep), `PEPMonitorSet()`
452: @*/
453: PetscErrorCode PEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
454: {
455: PetscFunctionBegin;
456: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
457: (*vf)->data = ctx;
458: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
459: PetscFunctionReturn(PETSC_SUCCESS);
460: }
462: /*@C
463: PEPMonitorAllDrawLG - Plots the error estimate of all unconverged
464: approximations at each iteration of the polynomial eigensolver.
466: Collective
468: Input Parameters:
469: + pep - the polynomial eigensolver context
470: . its - iteration number
471: . nconv - number of converged eigenpairs so far
472: . eigr - real part of the eigenvalues
473: . eigi - imaginary part of the eigenvalues
474: . errest - error estimates
475: . nest - number of error estimates to display
476: - vf - viewer and format for monitoring
478: Options Database Key:
479: . -pep_monitor_all draw::draw_lg - activates `PEPMonitorAllDrawLG()`
481: Notes:
482: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
483: function as an argument, to cause the monitor to be used during the `PEP` solve.
485: Call `PEPMonitorAllDrawLGCreate()` to create the context used with this monitor.
487: Level: intermediate
489: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorAllDrawLGCreate()`
490: @*/
491: PetscErrorCode PEPMonitorAllDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
492: {
493: PetscViewer viewer = vf->viewer;
494: PetscDrawLG lg;
495: PetscInt i,n = PetscMin(pep->nev,255);
496: PetscReal *x,*y;
498: PetscFunctionBegin;
501: PetscCall(PetscViewerPushFormat(viewer,vf->format));
502: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
503: if (its==1) {
504: PetscCall(PetscDrawLGReset(lg));
505: PetscCall(PetscDrawLGSetDimension(lg,n));
506: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0));
507: }
508: PetscCall(PetscMalloc2(n,&x,n,&y));
509: for (i=0;i<n;i++) {
510: x[i] = (PetscReal)its;
511: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
512: else y[i] = 0.0;
513: }
514: PetscCall(PetscDrawLGAddPoint(lg,x,y));
515: if (its <= 20 || !(its % 5) || pep->reason) {
516: PetscCall(PetscDrawLGDraw(lg));
517: PetscCall(PetscDrawLGSave(lg));
518: }
519: PetscCall(PetscFree2(x,y));
520: PetscCall(PetscViewerPopFormat(viewer));
521: PetscFunctionReturn(PETSC_SUCCESS);
522: }
524: /*@C
525: PEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
527: Collective
529: Input Parameters:
530: + viewer - the viewer
531: . format - the viewer format
532: - ctx - an optional user context
534: Output Parameter:
535: . vf - the viewer and format context
537: Level: intermediate
539: .seealso: [](ch:pep), `PEPMonitorSet()`
540: @*/
541: PetscErrorCode PEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
542: {
543: PetscFunctionBegin;
544: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
545: (*vf)->data = ctx;
546: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
547: PetscFunctionReturn(PETSC_SUCCESS);
548: }
550: /*@C
551: PEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
552: at each iteration of the polynomial eigensolver.
554: Collective
556: Input Parameters:
557: + pep - the polynomial eigensolver context
558: . its - iteration number
559: . nconv - number of converged eigenpairs so far
560: . eigr - real part of the eigenvalues
561: . eigi - imaginary part of the eigenvalues
562: . errest - error estimates
563: . nest - number of error estimates to display
564: - vf - viewer and format for monitoring
566: Options Database Key:
567: . -pep_monitor_conv draw::draw_lg - activates `PEPMonitorConvergedDrawLG()`
569: Notes:
570: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
571: function as an argument, to cause the monitor to be used during the `PEP` solve.
573: Call `PEPMonitorConvergedDrawLGCreate()` to create the context used with this monitor.
575: Level: intermediate
577: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorConvergedDrawLGCreate()`
578: @*/
579: PetscErrorCode PEPMonitorConvergedDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
580: {
581: PetscViewer viewer = vf->viewer;
582: PetscDrawLG lg;
583: PetscReal x,y;
585: PetscFunctionBegin;
588: PetscCall(PetscViewerPushFormat(viewer,vf->format));
589: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
590: if (its==1) {
591: PetscCall(PetscDrawLGReset(lg));
592: PetscCall(PetscDrawLGSetDimension(lg,1));
593: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,pep->nev));
594: }
595: x = (PetscReal)its;
596: y = (PetscReal)pep->nconv;
597: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
598: if (its <= 20 || !(its % 5) || pep->reason) {
599: PetscCall(PetscDrawLGDraw(lg));
600: PetscCall(PetscDrawLGSave(lg));
601: }
602: PetscCall(PetscViewerPopFormat(viewer));
603: PetscFunctionReturn(PETSC_SUCCESS);
604: }
606: /*@C
607: PEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
609: Collective
611: Input Parameters:
612: + viewer - the viewer
613: . format - the viewer format
614: - ctx - an optional user context
616: Output Parameter:
617: . vf - the viewer and format context
619: Level: intermediate
621: .seealso: [](ch:pep), `PEPMonitorSet()`
622: @*/
623: PetscErrorCode PEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
624: {
625: PetscInt *oldnconv;
627: PetscFunctionBegin;
628: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
629: PetscCall(PetscNew(&oldnconv));
630: (*vf)->data = (void*)oldnconv;
631: (*vf)->data_destroy = PetscCtxDestroyDefault;
632: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
633: PetscFunctionReturn(PETSC_SUCCESS);
634: }