Actual source code: lmemon.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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> /*I "slepclme.h" I*/
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: {
23: PetscInt i,n = lme->numbermonitors;
26: for (i=0;i<n;i++) {
27: (*lme->monitor[i])(lme,it,errest,lme->monitorcontext[i]);
28: }
29: return(0);
30: }
32: /*@C
33: LMEMonitorSet - Sets an ADDITIONAL function to be called at every
34: iteration to monitor convergence.
36: Logically Collective on LME
38: Input Parameters:
39: + lme - linear matrix equation solver context obtained from LMECreate()
40: . monitor - pointer to function (if this is NULL, it turns off monitoring)
41: . mctx - [optional] context for private data for the
42: monitor routine (use NULL if no context is desired)
43: - monitordestroy - [optional] routine that frees monitor context (may be NULL)
45: Calling Sequence of monitor:
46: $ monitor(LME lme,int its,PetscReal errest,void *mctx)
48: + lme - linear matrix equation solver context obtained from LMECreate()
49: . its - iteration number
50: . errest - error estimate
51: - mctx - optional monitoring context, as set by LMEMonitorSet()
53: Options Database Keys:
54: + -lme_monitor - print the error estimate
55: . -lme_monitor_lg - sets line graph monitor for the error estimate
56: - -lme_monitor_cancel - cancels all monitors that have been hardwired into
57: a code by calls to LMEMonitorSet(), but does not cancel those set via
58: the options database.
60: Notes:
61: Several different monitoring routines may be set by calling
62: LMEMonitorSet() multiple times; all will be called in the
63: order in which they were set.
65: Level: intermediate
67: .seealso: LMEMonitorCancel()
68: @*/
69: PetscErrorCode LMEMonitorSet(LME lme,PetscErrorCode (*monitor)(LME,PetscInt,PetscReal,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
70: {
73: if (lme->numbermonitors >= MAXLMEMONITORS) SETERRQ(PetscObjectComm((PetscObject)lme),PETSC_ERR_ARG_OUTOFRANGE,"Too many LME monitors set");
74: lme->monitor[lme->numbermonitors] = monitor;
75: lme->monitorcontext[lme->numbermonitors] = (void*)mctx;
76: lme->monitordestroy[lme->numbermonitors++] = monitordestroy;
77: return(0);
78: }
80: /*@
81: LMEMonitorCancel - Clears all monitors for an LME object.
83: Logically Collective on LME
85: Input Parameters:
86: . lme - linear matrix equation solver context obtained from LMECreate()
88: Options Database Key:
89: . -lme_monitor_cancel - Cancels all monitors that have been hardwired
90: into a code by calls to LMEMonitorSet(),
91: but does not cancel those set via the options database.
93: Level: intermediate
95: .seealso: LMEMonitorSet()
96: @*/
97: PetscErrorCode LMEMonitorCancel(LME lme)
98: {
100: PetscInt i;
104: for (i=0; i<lme->numbermonitors; i++) {
105: if (lme->monitordestroy[i]) {
106: (*lme->monitordestroy[i])(&lme->monitorcontext[i]);
107: }
108: }
109: lme->numbermonitors = 0;
110: return(0);
111: }
113: /*@C
114: LMEGetMonitorContext - Gets the monitor context, as set by
115: LMEMonitorSet() for the FIRST monitor only.
117: Not Collective
119: Input Parameter:
120: . lme - linear matrix equation solver context obtained from LMECreate()
122: Output Parameter:
123: . ctx - monitor context
125: Level: intermediate
127: .seealso: LMEMonitorSet()
128: @*/
129: PetscErrorCode LMEGetMonitorContext(LME lme,void **ctx)
130: {
133: *ctx = lme->monitorcontext[0];
134: return(0);
135: }
137: /*@C
138: LMEMonitorDefault - Print the error estimate of the current approximation at each
139: iteration of the linear matrix equation solver.
141: Collective on LME
143: Input Parameters:
144: + lme - linear matrix equation solver context
145: . its - iteration number
146: . errest - error estimate
147: - vf - viewer and format for monitoring
149: Level: intermediate
151: .seealso: LMEMonitorSet()
152: @*/
153: PetscErrorCode LMEMonitorDefault(LME lme,PetscInt its,PetscReal errest,PetscViewerAndFormat *vf)
154: {
156: PetscViewer viewer;
161: viewer = vf->viewer;
163: PetscViewerPushFormat(viewer,vf->format);
164: PetscViewerASCIIAddTab(viewer,((PetscObject)lme)->tablevel);
165: if (its == 1 && ((PetscObject)lme)->prefix) {
166: PetscViewerASCIIPrintf(viewer," Error estimates for %s solve.\n",((PetscObject)lme)->prefix);
167: }
168: PetscViewerASCIIPrintf(viewer,"%3D LME Error estimate %14.12e\n",its,(double)errest);
169: PetscViewerASCIISubtractTab(viewer,((PetscObject)lme)->tablevel);
170: PetscViewerPopFormat(viewer);
171: return(0);
172: }
174: /*@C
175: LMEMonitorLGCreate - Creates a line graph context for use with
176: LME to monitor convergence.
178: Collective on MPI_Comm
180: Input Parameters:
181: + comm - communicator context
182: . host - the X display to open, or null for the local machine
183: . label - the title to put in the title bar
184: . x, y - the screen coordinates of the upper left coordinate of
185: the window
186: - m, n - the screen width and height in pixels
188: Output Parameter:
189: . lgctx - the drawing context
191: Options Database Keys:
192: . -lme_monitor_lg - Sets line graph monitor
194: Notes:
195: Use PetscDrawLGDestroy() to destroy this line graph.
197: Level: intermediate
199: .seealso: LMEMonitorSet()
200: @*/
201: PetscErrorCode LMEMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
202: {
203: PetscDraw draw;
204: PetscDrawLG lg;
208: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
209: PetscDrawSetFromOptions(draw);
210: PetscDrawLGCreate(draw,1,&lg);
211: PetscDrawLGSetFromOptions(lg);
212: PetscDrawDestroy(&draw);
213: *lgctx = lg;
214: return(0);
215: }
217: PetscErrorCode LMEMonitorLG(LME lme,PetscInt its,PetscReal errest,void *ctx)
218: {
219: PetscDrawLG lg = (PetscDrawLG)ctx;
220: PetscReal x,y;
225: if (its==1) {
226: PetscDrawLGReset(lg);
227: PetscDrawLGSetDimension(lg,1);
228: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(lme->tol)-2,0.0);
229: }
230: x = (PetscReal)its;
231: if (errest > 0.0) y = PetscLog10Real(errest);
232: else y = 0.0;
233: PetscDrawLGAddPoint(lg,&x,&y);
234: if (its <= 20 || !(its % 5) || lme->reason) {
235: PetscDrawLGDraw(lg);
236: PetscDrawLGSave(lg);
237: }
238: return(0);
239: }