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     - linear matrix equation solver context obtained from LMECreate()
 37: .  monitor - pointer to function (if this is NULL, it turns off monitoring), see LMEMonitorFn
 38: .  mctx    - [optional] context for private data for the
 39:              monitor routine (use NULL if no context is desired)
 40: -  monitordestroy - [optional] routine that frees monitor context (may be NULL),
 41:              see PetscCtxDestroyFn for the calling sequence

 43:    Options Database Keys:
 44: +    -lme_monitor - print the error estimate
 45: .    -lme_monitor draw::draw_lg - sets line graph monitor for the error estimate
 46: -    -lme_monitor_cancel - cancels all monitors that have been hardwired into
 47:       a code by calls to LMEMonitorSet(), but does not cancel those set via
 48:       the options database.

 50:    Notes:
 51:    The options database option -lme_monitor and related options are the easiest way
 52:    to turn on LME iteration monitoring.

 54:    LMEMonitorRegister() provides a way to associate an options database key with LME
 55:    monitor function.

 57:    Several different monitoring routines may be set by calling LMEMonitorSet() multiple
 58:    times; all will be called in the order in which they were set.

 60:    Level: intermediate

 62: .seealso: `LMEMonitorCancel()`
 63: @*/
 64: PetscErrorCode LMEMonitorSet(LME lme,LMEMonitorFn *monitor,void *mctx,PetscCtxDestroyFn *monitordestroy)
 65: {
 66:   PetscInt  i;
 67:   PetscBool identical;

 69:   PetscFunctionBegin;
 71:   for (i=0;i<lme->numbermonitors;i++) {
 72:     PetscCall(PetscMonitorCompare((PetscErrorCode(*)(void))(PetscVoidFn*)monitor,mctx,monitordestroy,(PetscErrorCode (*)(void))(PetscVoidFn*)lme->monitor[i],lme->monitorcontext[i],lme->monitordestroy[i],&identical));
 73:     if (identical) PetscFunctionReturn(PETSC_SUCCESS);
 74:   }
 75:   PetscCheck(lme->numbermonitors<MAXLMEMONITORS,PetscObjectComm((PetscObject)lme),PETSC_ERR_ARG_OUTOFRANGE,"Too many LME monitors set");
 76:   lme->monitor[lme->numbermonitors]           = monitor;
 77:   lme->monitorcontext[lme->numbermonitors]    = mctx;
 78:   lme->monitordestroy[lme->numbermonitors++]  = monitordestroy;
 79:   PetscFunctionReturn(PETSC_SUCCESS);
 80: }

 82: /*@
 83:    LMEMonitorCancel - Clears all monitors for an LME object.

 85:    Logically Collective

 87:    Input Parameters:
 88: .  lme - linear matrix equation solver context obtained from LMECreate()

 90:    Options Database Key:
 91: .    -lme_monitor_cancel - cancels all monitors that have been hardwired
 92:       into a code by calls to LMEMonitorSet(),
 93:       but does not cancel those set via the options database.

 95:    Level: intermediate

 97: .seealso: `LMEMonitorSet()`
 98: @*/
 99: PetscErrorCode LMEMonitorCancel(LME lme)
100: {
101:   PetscInt       i;

103:   PetscFunctionBegin;
105:   for (i=0; i<lme->numbermonitors; i++) {
106:     if (lme->monitordestroy[i]) PetscCall((*lme->monitordestroy[i])(&lme->monitorcontext[i]));
107:   }
108:   lme->numbermonitors = 0;
109:   PetscFunctionReturn(PETSC_SUCCESS);
110: }

112: /*@C
113:    LMEGetMonitorContext - Gets the monitor context, as set by
114:    LMEMonitorSet() for the FIRST monitor only.

116:    Not Collective

118:    Input Parameter:
119: .  lme - linear matrix equation solver context obtained from LMECreate()

121:    Output Parameter:
122: .  ctx - monitor context

124:    Level: intermediate

126: .seealso: `LMEMonitorSet()`
127: @*/
128: PetscErrorCode LMEGetMonitorContext(LME lme,void *ctx)
129: {
130:   PetscFunctionBegin;
132:   *(void**)ctx = lme->monitorcontext[0];
133:   PetscFunctionReturn(PETSC_SUCCESS);
134: }

136: /*@C
137:    LMEMonitorDefault - Print the error estimate of the current approximation at each
138:    iteration of the linear matrix equation solver.

140:    Collective

142:    Input Parameters:
143: +  lme    - linear matrix equation solver context
144: .  its    - iteration number
145: .  errest - error estimate
146: -  vf     - viewer and format for monitoring

148:    Options Database Key:
149: .  -lme_monitor - activates LMEMonitorDefault()

151:    Level: intermediate

153: .seealso: `LMEMonitorSet()`
154: @*/
155: PetscErrorCode LMEMonitorDefault(LME lme,PetscInt its,PetscReal errest,PetscViewerAndFormat *vf)
156: {
157:   PetscViewer    viewer = vf->viewer;

159:   PetscFunctionBegin;
162:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
163:   PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)lme)->tablevel));
164:   if (its == 1 && ((PetscObject)lme)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer,"  Error estimates for %s solve.\n",((PetscObject)lme)->prefix));
165:   PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " LME Error estimate %14.12e\n",its,(double)errest));
166:   PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)lme)->tablevel));
167:   PetscCall(PetscViewerPopFormat(viewer));
168:   PetscFunctionReturn(PETSC_SUCCESS);
169: }

171: /*@C
172:    LMEMonitorDefaultDrawLG - Plots the error estimate of the current approximation at each
173:    iteration of the linear matrix equation solver.

175:    Collective

177:    Input Parameters:
178: +  lme    - linear matrix equation solver context
179: .  its    - iteration number
180: .  errest - error estimate
181: -  vf     - viewer and format for monitoring

183:    Options Database Key:
184: .  -lme_monitor draw::draw_lg - activates LMEMonitorDefaultDrawLG()

186:    Level: intermediate

188: .seealso: `LMEMonitorSet()`
189: @*/
190: PetscErrorCode LMEMonitorDefaultDrawLG(LME lme,PetscInt its,PetscReal errest,PetscViewerAndFormat *vf)
191: {
192:   PetscViewer    viewer = vf->viewer;
193:   PetscDrawLG    lg;
194:   PetscReal      x,y;

196:   PetscFunctionBegin;
199:   PetscCall(PetscViewerPushFormat(viewer,vf->format));
200:   PetscCall(PetscViewerDrawGetDrawLG(viewer,0,&lg));
201:   if (its==1) {
202:     PetscCall(PetscDrawLGReset(lg));
203:     PetscCall(PetscDrawLGSetDimension(lg,1));
204:     PetscCall(PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(lme->tol)-2,0.0));
205:   }
206:   x = (PetscReal)its;
207:   if (errest > 0.0) y = PetscLog10Real(errest);
208:   else y = 0.0;
209:   PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
210:   if (its <= 20 || !(its % 5) || lme->reason) {
211:     PetscCall(PetscDrawLGDraw(lg));
212:     PetscCall(PetscDrawLGSave(lg));
213:   }
214:   PetscCall(PetscViewerPopFormat(viewer));
215:   PetscFunctionReturn(PETSC_SUCCESS);
216: }

218: /*@C
219:    LMEMonitorDefaultDrawLGCreate - Creates the plotter for the error estimate.

221:    Collective

223:    Input Parameters:
224: +  viewer - the viewer
225: .  format - the viewer format
226: -  ctx    - an optional user context

228:    Output Parameter:
229: .  vf     - the viewer and format context

231:    Level: intermediate

233: .seealso: `LMEMonitorSet()`
234: @*/
235: PetscErrorCode LMEMonitorDefaultDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
236: {
237:   PetscFunctionBegin;
238:   PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
239:   (*vf)->data = ctx;
240:   PetscCall(PetscViewerMonitorLGSetUp(viewer,NULL,"Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300));
241:   PetscFunctionReturn(PETSC_SUCCESS);
242: }