Actual source code: lmemon.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: LME routines related to monitors
12: */
14: #include <slepc/private/lmeimpl.h>
15: #include <petscdraw.h>
17: /*
18: Runs the user provided monitor routines, if any.
19: */
20: PetscErrorCode LMEMonitor(LME lme,PetscInt it,PetscReal errest)
21: {
22: PetscInt i,n = lme->numbermonitors;
24: PetscFunctionBegin;
25: for (i=0;i<n;i++) PetscCall((*lme->monitor[i])(lme,it,errest,lme->monitorcontext[i]));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: /*@C
30: LMEMonitorSet - Sets an ADDITIONAL function to be called at every
31: iteration to monitor convergence.
33: Logically Collective
35: Input Parameters:
36: + lme - the linear matrix equation solver context
37: . monitor - pointer to function (if this is `NULL`, it turns off monitoring),
38: see `LMEMonitorFn`
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: + -lme_monitor - print the error estimate
46: . -lme_monitor draw::draw_lg - sets line graph monitor for the error estimate
47: - -lme_monitor_cancel - cancels all monitors that have been hardwired into
48: a code by calls to `LMEMonitorSet()`, but does not cancel
49: those set via the options database
51: Notes:
52: The options database option `-lme_monitor` and related options are the easiest way
53: to turn on `LME` iteration monitoring.
55: `LMEMonitorRegister()` provides a way to associate an options database key with `LME`
56: monitor function.
58: Several different monitoring routines may be set by calling `LMEMonitorSet()` multiple
59: times; all will be called in the order in which they were set.
61: Fortran Note:
62: Only a single monitor function can be set for each `LME` object.
64: Level: intermediate
66: .seealso: [](ch:lme), `LMEMonitorDefault()`, `LMEMonitorDefaultDrawLG()`, `LMEMonitorCancel()`
67: @*/
68: PetscErrorCode LMEMonitorSet(LME lme,LMEMonitorFn *monitor,void *mctx,PetscCtxDestroyFn *monitordestroy)
69: {
70: PetscInt i;
71: PetscBool identical;
73: PetscFunctionBegin;
75: for (i=0;i<lme->numbermonitors;i++) {
76: PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))(PetscVoidFn*)monitor,mctx,monitordestroy,(PetscErrorCode (*)(void))(PetscVoidFn*)lme->monitor[i],lme->monitorcontext[i],lme->monitordestroy[i],&identical));
77: if (identical) PetscFunctionReturn(PETSC_SUCCESS);
78: }
79: PetscCheck(lme->numbermonitors<MAXLMEMONITORS,PetscObjectComm((PetscObject)lme),PETSC_ERR_ARG_OUTOFRANGE,"Too many LME monitors set");
80: lme->monitor[lme->numbermonitors] = monitor;
81: lme->monitorcontext[lme->numbermonitors] = mctx;
82: lme->monitordestroy[lme->numbermonitors++] = monitordestroy;
83: PetscFunctionReturn(PETSC_SUCCESS);
84: }
86: /*@
87: LMEMonitorCancel - Clears all monitors for an `LME` object.
89: Logically Collective
91: Input Parameter:
92: . lme - the linear matrix equation solver context
94: Options Database Key:
95: . -lme_monitor_cancel - cancels all monitors that have been hardwired into a code by calls to
96: `LMEMonitorSet()`, but does not cancel those set via the options database.
98: Level: intermediate
100: .seealso: [](ch:lme), `LMEMonitorSet()`
101: @*/
102: PetscErrorCode LMEMonitorCancel(LME lme)
103: {
104: PetscInt i;
106: PetscFunctionBegin;
108: for (i=0; i<lme->numbermonitors; i++) {
109: if (lme->monitordestroy[i]) PetscCall((*lme->monitordestroy[i])(&lme->monitorcontext[i]));
110: }
111: lme->numbermonitors = 0;
112: PetscFunctionReturn(PETSC_SUCCESS);
113: }
115: /*@C
116: LMEGetMonitorContext - Gets the monitor context, as set by
117: `LMEMonitorSet()` for the FIRST monitor only.
119: Not Collective
121: Input Parameter:
122: . lme - the linear matrix equation solver context
124: Output Parameter:
125: . ctx - monitor context
127: Level: intermediate
129: .seealso: [](ch:lme), `LMEMonitorSet()`
130: @*/
131: PetscErrorCode LMEGetMonitorContext(LME lme,void *ctx)
132: {
133: PetscFunctionBegin;
135: *(void**)ctx = lme->monitorcontext[0];
136: PetscFunctionReturn(PETSC_SUCCESS);
137: }
139: /*@C
140: LMEMonitorDefault - Print the error estimate of the current approximation at each
141: iteration of the linear matrix equation solver.
143: Collective
145: Input Parameters:
146: + lme - the linear matrix equation solver context
147: . its - iteration number
148: . errest - error estimate
149: - vf - viewer and format for monitoring
151: Options Database Key:
152: . -lme_monitor - activates `LMEMonitorDefault()`
154: Note:
155: This is not called directly by users, rather one calls `LMEMonitorSet()`, with this
156: function as an argument, to cause the monitor to be used during the `LME` solve.
158: Level: intermediate
160: .seealso: [](ch:lme), `LMEMonitorSet()`
161: @*/
162: PetscErrorCode LMEMonitorDefault(LME lme,PetscInt its,PetscReal errest,PetscViewerAndFormat *vf)
163: {
164: PetscViewer viewer = vf->viewer;
166: PetscFunctionBegin;
169: PetscCall(PetscViewerPushFormat(viewer,vf->format));
170: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)lme)->tablevel));
171: if (its == 1 && ((PetscObject)lme)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Error estimates for %s solve.\n",((PetscObject)lme)->prefix));
172: PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " LME Error estimate %14.12e\n",its,(double)errest));
173: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)lme)->tablevel));
174: PetscCall(PetscViewerPopFormat(viewer));
175: PetscFunctionReturn(PETSC_SUCCESS);
176: }
178: /*@C
179: LMEMonitorDefaultDrawLG - Plots the error estimate of the current approximation at each
180: iteration of the linear matrix equation solver.
182: Collective
184: Input Parameters:
185: + lme - the linear matrix equation solver context
186: . its - iteration number
187: . errest - error estimate
188: - vf - viewer and format for monitoring
190: Options Database Key:
191: . -lme_monitor draw::draw_lg - activates `LMEMonitorDefaultDrawLG()`
193: Notes:
194: This is not called directly by users, rather one calls `LMEMonitorSet()`, with this
195: function as an argument, to cause the monitor to be used during the `LME` solve.
197: Call `LMEMonitorDefaultDrawLGCreate()` to create the context used with this monitor.
199: Level: intermediate
201: .seealso: [](ch:lme), `LMEMonitorSet()`, `LMEMonitorDefaultDrawLGCreate()`
202: @*/
203: PetscErrorCode LMEMonitorDefaultDrawLG(LME lme,PetscInt its,PetscReal errest,PetscViewerAndFormat *vf)
204: {
205: PetscViewer viewer = vf->viewer;
206: PetscDrawLG lg;
207: PetscReal x,y;
209: PetscFunctionBegin;
212: PetscCall(PetscViewerPushFormat(viewer,vf->format));
213: PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
214: if (its==1) {
215: PetscCall(PetscDrawLGReset(lg));
216: PetscCall(PetscDrawLGSetDimension(lg,1));
217: PetscCall(PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(lme->tol)-2,0.0));
218: }
219: x = (PetscReal)its;
220: if (errest > 0.0) y = PetscLog10Real(errest);
221: else y = 0.0;
222: PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
223: if (its <= 20 || !(its % 5) || lme->reason) {
224: PetscCall(PetscDrawLGDraw(lg));
225: PetscCall(PetscDrawLGSave(lg));
226: }
227: PetscCall(PetscViewerPopFormat(viewer));
228: PetscFunctionReturn(PETSC_SUCCESS);
229: }
231: /*@C
232: LMEMonitorDefaultDrawLGCreate - Creates the plotter for the error estimate.
234: Collective
236: Input Parameters:
237: + viewer - the viewer
238: . format - the viewer format
239: - ctx - an optional user context
241: Output Parameter:
242: . vf - the viewer and format context
244: Level: intermediate
246: .seealso: [](ch:lme), `LMEMonitorSet()`
247: @*/
248: PetscErrorCode LMEMonitorDefaultDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
249: {
250: PetscFunctionBegin;
251: PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
252: (*vf)->data = ctx;
253: PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
254: PetscFunctionReturn(PETSC_SUCCESS);
255: }