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: . mctx - [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,void *mctx,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,mctx,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] = mctx;
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,void *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: Note:
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: Level: intermediate
308: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorFirst()`, `PEPMonitorAll()`
309: @*/
310: PetscErrorCode PEPMonitorConverged(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
311: {
312: PetscInt i;
313: PetscScalar er,ei;
314: PetscViewer viewer = vf->viewer;
315: SlepcConvMon ctx;
317: PetscFunctionBegin;
320: ctx = (SlepcConvMon)vf->data;
321: if (its==1 && ((PetscObject)pep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)pep)->prefix));
322: if (its==1) ctx->oldnconv = 0;
323: if (ctx->oldnconv!=nconv) {
324: PetscCall(PetscViewerPushFormat(viewer,vf->format));
325: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel));
326: for (i=ctx->oldnconv;i<nconv;i++) {
327: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " PEP converged value (error) #%" PetscInt_FMT,its,i));
328: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
329: er = eigr[i]; ei = eigi[i];
330: PetscCall(PEPMonitorGetTrueEig(pep,&er,&ei));
331: #if defined(PETSC_USE_COMPLEX)
332: PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(er),(double)PetscImaginaryPart(er)));
333: #else
334: PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)er));
335: if (ei!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)ei));
336: #endif
337: PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
338: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
339: }
340: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel));
341: PetscCall(PetscViewerPopFormat(viewer));
342: ctx->oldnconv = nconv;
343: }
344: PetscFunctionReturn(PETSC_SUCCESS);
345: }
347: PetscErrorCode PEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
348: {
349: SlepcConvMon mctx;
351: PetscFunctionBegin;
352: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
353: PetscCall(PetscNew(&mctx));
354: mctx->ctx = ctx;
355: (*vf)->data = (void*)mctx;
356: PetscFunctionReturn(PETSC_SUCCESS);
357: }
359: PetscErrorCode PEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
360: {
361: PetscFunctionBegin;
362: if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
363: PetscCall(PetscFree((*vf)->data));
364: PetscCall(PetscViewerDestroy(&(*vf)->viewer));
365: PetscCall(PetscFree(*vf));
366: PetscFunctionReturn(PETSC_SUCCESS);
367: }
369: /*@C
370: PEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
371: approximation at each iteration of the polynomial eigensolver.
373: Collective
375: Input Parameters:
376: + pep - the polynomial eigensolver context
377: . its - iteration number
378: . nconv - number of converged eigenpairs so far
379: . eigr - real part of the eigenvalues
380: . eigi - imaginary part of the eigenvalues
381: . errest - error estimates
382: . nest - number of error estimates to display
383: - vf - viewer and format for monitoring
385: Options Database Key:
386: . -pep_monitor draw::draw_lg - activates `PEPMonitorFirstDrawLG()`
388: Notes:
389: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
390: function as an argument, to cause the monitor to be used during the `PEP` solve.
392: Call `PEPMonitorFirstDrawLGCreate()` to create the context used with this monitor.
394: Level: intermediate
396: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorFirstDrawLGCreate()`
397: @*/
398: PetscErrorCode PEPMonitorFirstDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
399: {
400: PetscViewer viewer = vf->viewer;
401: PetscDrawLG lg;
402: PetscReal x,y;
404: PetscFunctionBegin;
407: PetscCall(PetscViewerPushFormat(viewer,vf->format));
408: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
409: if (its==1) {
410: PetscCall(PetscDrawLGReset(lg));
411: PetscCall(PetscDrawLGSetDimension(lg,1));
412: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0));
413: }
414: if (nconv<nest) {
415: x = (PetscReal)its;
416: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
417: else y = 0.0;
418: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
419: if (its <= 20 || !(its % 5) || pep->reason) {
420: PetscCall(PetscDrawLGDraw(lg));
421: PetscCall(PetscDrawLGSave(lg));
422: }
423: }
424: PetscCall(PetscViewerPopFormat(viewer));
425: PetscFunctionReturn(PETSC_SUCCESS);
426: }
428: /*@C
429: PEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
431: Collective
433: Input Parameters:
434: + viewer - the viewer
435: . format - the viewer format
436: - ctx - an optional user context
438: Output Parameter:
439: . vf - the viewer and format context
441: Level: intermediate
443: .seealso: [](ch:pep), `PEPMonitorSet()`
444: @*/
445: PetscErrorCode PEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
446: {
447: PetscFunctionBegin;
448: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
449: (*vf)->data = ctx;
450: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
451: PetscFunctionReturn(PETSC_SUCCESS);
452: }
454: /*@C
455: PEPMonitorAllDrawLG - Plots the error estimate of all unconverged
456: approximations at each iteration of the polynomial eigensolver.
458: Collective
460: Input Parameters:
461: + pep - the polynomial eigensolver context
462: . its - iteration number
463: . nconv - number of converged eigenpairs so far
464: . eigr - real part of the eigenvalues
465: . eigi - imaginary part of the eigenvalues
466: . errest - error estimates
467: . nest - number of error estimates to display
468: - vf - viewer and format for monitoring
470: Options Database Key:
471: . -pep_monitor_all draw::draw_lg - activates `PEPMonitorAllDrawLG()`
473: Notes:
474: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
475: function as an argument, to cause the monitor to be used during the `PEP` solve.
477: Call `PEPMonitorAllDrawLGCreate()` to create the context used with this monitor.
479: Level: intermediate
481: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorAllDrawLGCreate()`
482: @*/
483: PetscErrorCode PEPMonitorAllDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
484: {
485: PetscViewer viewer = vf->viewer;
486: PetscDrawLG lg;
487: PetscInt i,n = PetscMin(pep->nev,255);
488: PetscReal *x,*y;
490: PetscFunctionBegin;
493: PetscCall(PetscViewerPushFormat(viewer,vf->format));
494: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
495: if (its==1) {
496: PetscCall(PetscDrawLGReset(lg));
497: PetscCall(PetscDrawLGSetDimension(lg,n));
498: PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(pep->tol)-2,0.0));
499: }
500: PetscCall(PetscMalloc2(n,&x,n,&y));
501: for (i=0;i<n;i++) {
502: x[i] = (PetscReal)its;
503: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
504: else y[i] = 0.0;
505: }
506: PetscCall(PetscDrawLGAddPoint(lg,x,y));
507: if (its <= 20 || !(its % 5) || pep->reason) {
508: PetscCall(PetscDrawLGDraw(lg));
509: PetscCall(PetscDrawLGSave(lg));
510: }
511: PetscCall(PetscFree2(x,y));
512: PetscCall(PetscViewerPopFormat(viewer));
513: PetscFunctionReturn(PETSC_SUCCESS);
514: }
516: /*@C
517: PEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
519: Collective
521: Input Parameters:
522: + viewer - the viewer
523: . format - the viewer format
524: - ctx - an optional user context
526: Output Parameter:
527: . vf - the viewer and format context
529: Level: intermediate
531: .seealso: [](ch:pep), `PEPMonitorSet()`
532: @*/
533: PetscErrorCode PEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
534: {
535: PetscFunctionBegin;
536: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
537: (*vf)->data = ctx;
538: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
539: PetscFunctionReturn(PETSC_SUCCESS);
540: }
542: /*@C
543: PEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
544: at each iteration of the polynomial eigensolver.
546: Collective
548: Input Parameters:
549: + pep - the polynomial eigensolver context
550: . its - iteration number
551: . nconv - number of converged eigenpairs so far
552: . eigr - real part of the eigenvalues
553: . eigi - imaginary part of the eigenvalues
554: . errest - error estimates
555: . nest - number of error estimates to display
556: - vf - viewer and format for monitoring
558: Options Database Key:
559: . -pep_monitor_conv draw::draw_lg - activates `PEPMonitorConvergedDrawLG()`
561: Notes:
562: This is not called directly by users, rather one calls `PEPMonitorSet()`, with this
563: function as an argument, to cause the monitor to be used during the `PEP` solve.
565: Call `PEPMonitorConvergedDrawLGCreate()` to create the context used with this monitor.
567: Level: intermediate
569: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorConvergedDrawLGCreate()`
570: @*/
571: PetscErrorCode PEPMonitorConvergedDrawLG(PEP pep,PetscInt its,PetscInt nconv,PetscScalar eigr[],PetscScalar eigi[],PetscReal errest[],PetscInt nest,PetscViewerAndFormat *vf)
572: {
573: PetscViewer viewer = vf->viewer;
574: PetscDrawLG lg;
575: PetscReal x,y;
577: PetscFunctionBegin;
580: PetscCall(PetscViewerPushFormat(viewer,vf->format));
581: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
582: if (its==1) {
583: PetscCall(PetscDrawLGReset(lg));
584: PetscCall(PetscDrawLGSetDimension(lg,1));
585: PetscCall(PetscDrawLGSetLimits(lg,1,1,0,pep->nev));
586: }
587: x = (PetscReal)its;
588: y = (PetscReal)pep->nconv;
589: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
590: if (its <= 20 || !(its % 5) || pep->reason) {
591: PetscCall(PetscDrawLGDraw(lg));
592: PetscCall(PetscDrawLGSave(lg));
593: }
594: PetscCall(PetscViewerPopFormat(viewer));
595: PetscFunctionReturn(PETSC_SUCCESS);
596: }
598: /*@C
599: PEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
601: Collective
603: Input Parameters:
604: + viewer - the viewer
605: . format - the viewer format
606: - ctx - an optional user context
608: Output Parameter:
609: . vf - the viewer and format context
611: Level: intermediate
613: .seealso: [](ch:pep), `PEPMonitorSet()`
614: @*/
615: PetscErrorCode PEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
616: {
617: SlepcConvMon mctx;
619: PetscFunctionBegin;
620: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
621: PetscCall(PetscNew(&mctx));
622: mctx->ctx = ctx;
623: (*vf)->data = (void*)mctx;
624: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
625: PetscFunctionReturn(PETSC_SUCCESS);
626: }