Actual source code: svdmon.c
1: /*
2: SVD routines related to monitors.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2012, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc-private/svdimpl.h> /*I "slepcsvd.h" I*/
28: /*
29: Runs the user provided monitor routines, if any.
30: */
31: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest)
32: {
34: PetscInt i,n = svd->numbermonitors;
37: for (i=0;i<n;i++) {
38: (*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]);
39: }
40: return(0);
41: }
45: /*@C
46: SVDMonitorSet - Sets an ADDITIONAL function to be called at every
47: iteration to monitor the error estimates for each requested singular triplet.
48:
49: Collective on SVD
51: Input Parameters:
52: + svd - singular value solver context obtained from SVDCreate()
53: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
54: - mctx - [optional] context for private data for the
55: monitor routine (use PETSC_NULL if no context is desired)
57: Calling Sequence of monitor:
58: $ monitor (SVD svd, PetscInt its, PetscInt nconv, PetscReal *sigma, PetscReal* errest, PetscInt nest, void *mctx)
60: + svd - singular value solver context obtained from SVDCreate()
61: . its - iteration number
62: . nconv - number of converged singular triplets
63: . sigma - singular values
64: . errest - relative error estimates for each singular triplet
65: . nest - number of error estimates
66: - mctx - optional monitoring context, as set by SVDMonitorSet()
68: Options Database Keys:
69: + -svd_monitor - print only the first error estimate
70: . -svd_monitor_all - print error estimates at each iteration
71: . -svd_monitor_conv - print the singular value approximations only when
72: convergence has been reached
73: . -svd_monitor_draw - sets line graph monitor for the first unconverged
74: approximate singular value
75: . -svd_monitor_draw_all - sets line graph monitor for all unconverged
76: approximate singular value
77: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
78: a code by calls to SVDMonitorSet(), but does not cancel those set via
79: the options database.
81: Notes:
82: Several different monitoring routines may be set by calling
83: SVDMonitorSet() multiple times; all will be called in the
84: order in which they were set.
86: Level: intermediate
88: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorCancel()
89: @*/
90: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
91: {
94: if (svd->numbermonitors >= MAXSVDMONITORS) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
95: svd->monitor[svd->numbermonitors] = monitor;
96: svd->monitorcontext[svd->numbermonitors] = (void*)mctx;
97: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
98: return(0);
99: }
103: /*@
104: SVDMonitorCancel - Clears all monitors for an SVD object.
106: Collective on SVD
108: Input Parameters:
109: . svd - singular value solver context obtained from SVDCreate()
111: Options Database Key:
112: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
113: into a code by calls to SVDMonitorSet(),
114: but does not cancel those set via the options database.
116: Level: intermediate
118: .seealso: SVDMonitorSet()
119: @*/
120: PetscErrorCode SVDMonitorCancel(SVD svd)
121: {
123: PetscInt i;
127: for (i=0; i<svd->numbermonitors; i++) {
128: if (svd->monitordestroy[i]) {
129: (*svd->monitordestroy[i])(&svd->monitorcontext[i]);
130: }
131: }
132: svd->numbermonitors = 0;
133: return(0);
134: }
138: /*@C
139: SVDGetMonitorContext - Gets the monitor context, as set by
140: SVDMonitorSet() for the FIRST monitor only.
142: Not Collective
144: Input Parameter:
145: . svd - singular value solver context obtained from SVDCreate()
147: Output Parameter:
148: . ctx - monitor context
150: Level: intermediate
152: .seealso: SVDMonitorSet()
153: @*/
154: PetscErrorCode SVDGetMonitorContext(SVD svd,void **ctx)
155: {
158: *ctx = (svd->monitorcontext[0]);
159: return(0);
160: }
164: /*@C
165: SVDMonitorAll - Print the current approximate values and
166: error estimates at each iteration of the singular value solver.
168: Collective on SVD
170: Input Parameters:
171: + svd - singular value solver context
172: . its - iteration number
173: . nconv - number of converged singular triplets so far
174: . sigma - singular values
175: . errest - error estimates
176: . nest - number of error estimates to display
177: - monctx - monitor context (contains viewer, can be PETSC_NULL)
179: Level: intermediate
181: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConverged()
182: @*/
183: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
184: {
186: PetscInt i;
187: PetscViewer viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)svd)->comm);
190: if (its) {
191: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
192: PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D Values (Errors)",its,nconv);
193: for (i=0;i<nest;i++) {
194: PetscViewerASCIIPrintf(viewer," %G (%10.8e)",sigma[i],(double)errest[i]);
195: }
196: PetscViewerASCIIPrintf(viewer,"\n");
197: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
198: }
199: return(0);
200: }
204: /*@C
205: SVDMonitorFirst - Print the first unconverged approximate values and
206: error estimates at each iteration of the singular value solver.
208: Collective on SVD
210: Input Parameters:
211: + svd - singular value solver context
212: . its - iteration number
213: . nconv - number of converged singular triplets so far
214: . sigma - singular values
215: . errest - error estimates
216: . nest - number of error estimates to display
217: - monctx - monitor context (contains viewer, can be PETSC_NULL)
219: Level: intermediate
221: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConverged()
222: @*/
223: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
224: {
226: PetscViewer viewer = monctx? (PetscViewer)monctx: PETSC_VIEWER_STDOUT_(((PetscObject)svd)->comm);
229: if (its && nconv<nest) {
230: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
231: PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D first unconverged value (error)",its,nconv);
232: PetscViewerASCIIPrintf(viewer," %G (%10.8e)\n",sigma[nconv],(double)errest[nconv]);
233: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
234: }
235: return(0);
236: }
240: /*@C
241: SVDMonitorConverged - Print the approximate values and error estimates as they converge.
243: Collective on SVD
245: Input Parameters:
246: + svd - singular value solver context
247: . its - iteration number
248: . nconv - number of converged singular triplets so far
249: . sigma - singular values
250: . errest - error estimates
251: . nest - number of error estimates to display
252: - monctx - monitor context
254: Note:
255: The monitor context must contain a struct with a PetscViewer and a
256: PetscInt. In Fortran, pass a PETSC_NULL_OBJECT.
258: Level: intermediate
260: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorAll()
261: @*/
262: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
263: {
264: PetscErrorCode ierr;
265: PetscInt i;
266: PetscViewer viewer;
267: SlepcConvMonitor ctx = (SlepcConvMonitor) monctx;
270: if (!monctx) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_WRONG,"Must provide a context for SVDMonitorConverged");
271: if (!its) {
272: ctx->oldnconv = 0;
273: } else {
274: viewer = ctx->viewer? ctx->viewer: PETSC_VIEWER_STDOUT_(((PetscObject)svd)->comm);
275: for (i=ctx->oldnconv;i<nconv;i++) {
276: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
277: PetscViewerASCIIPrintf(viewer,"%3D SVD converged value (error) #%D",its,i);
278: PetscViewerASCIIPrintf(viewer," %G (%10.8e)\n",sigma[i],(double)errest[i]);
279: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
280: }
281: ctx->oldnconv = nconv;
282: }
283: return(0);
284: }
288: PetscErrorCode SVDMonitorLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
289: {
290: PetscViewer viewer = (PetscViewer) monctx;
291: PetscDraw draw,draw1;
292: PetscDrawLG lg,lg1;
294: PetscReal x,y,p;
297: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); }
298: PetscViewerDrawGetDraw(viewer,0,&draw);
299: PetscViewerDrawGetDrawLG(viewer,0,&lg);
300: PetscViewerDrawGetDraw(viewer,1,&draw1);
301: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
303: if (!its) {
304: PetscDrawSetTitle(draw,"Error estimates");
305: PetscDrawSetDoubleBuffer(draw);
306: PetscDrawLGSetDimension(lg,1);
307: PetscDrawLGReset(lg);
308: PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);
310: PetscDrawSetTitle(draw1,"Approximate singular values");
311: PetscDrawSetDoubleBuffer(draw1);
312: PetscDrawLGSetDimension(lg1,1);
313: PetscDrawLGReset(lg1);
314: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
315: }
317: x = (PetscReal) its;
318: if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
319: PetscDrawLGAddPoint(lg,&x,&y);
321: PetscDrawLGAddPoint(lg1,&x,svd->sigma);
322: PetscDrawGetPause(draw1,&p);
323: PetscDrawSetPause(draw1,0);
324: PetscDrawLGDraw(lg1);
325: PetscDrawSetPause(draw1,p);
326:
327: PetscDrawLGDraw(lg);
328: return(0);
329: }
333: PetscErrorCode SVDMonitorLGAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
334: {
335: PetscViewer viewer = (PetscViewer) monctx;
336: PetscDraw draw,draw1;
337: PetscDrawLG lg,lg1;
339: PetscReal *x,*y,p;
340: PetscInt i,n = PetscMin(svd->nsv,255);
343: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); }
344: PetscViewerDrawGetDraw(viewer,0,&draw);
345: PetscViewerDrawGetDrawLG(viewer,0,&lg);
346: PetscViewerDrawGetDraw(viewer,1,&draw1);
347: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
349: if (!its) {
350: PetscDrawSetTitle(draw,"Error estimates");
351: PetscDrawSetDoubleBuffer(draw);
352: PetscDrawLGSetDimension(lg,n);
353: PetscDrawLGReset(lg);
354: PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);
356: PetscDrawSetTitle(draw1,"Approximate singular values");
357: PetscDrawSetDoubleBuffer(draw1);
358: PetscDrawLGSetDimension(lg1,n);
359: PetscDrawLGReset(lg1);
360: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
361: }
363: PetscMalloc(sizeof(PetscReal)*n,&x);
364: PetscMalloc(sizeof(PetscReal)*n,&y);
365: for (i=0;i<n;i++) {
366: x[i] = (PetscReal) its;
367: if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
368: else y[i] = 0.0;
369: }
370: PetscDrawLGAddPoint(lg,x,y);
372: PetscDrawLGAddPoint(lg1,x,svd->sigma);
373: PetscDrawGetPause(draw1,&p);
374: PetscDrawSetPause(draw1,0);
375: PetscDrawLGDraw(lg1);
376: PetscDrawSetPause(draw1,p);
377:
378: PetscDrawLGDraw(lg);
379: PetscFree(x);
380: PetscFree(y);
381: return(0);
382: }