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-2010, Universidad 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 private/svdimpl.h
28: /*@C
29: SVDMonitorSet - Sets an ADDITIONAL function to be called at every
30: iteration to monitor the error estimates for each requested singular triplet.
31:
32: Collective on SVD
34: Input Parameters:
35: + svd - singular value solver context obtained from SVDCreate()
36: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
37: - mctx - [optional] context for private data for the
38: monitor routine (use PETSC_NULL if no context is desired)
40: Calling Sequence of monitor:
41: $ monitor (SVD svd, PetscInt its, PetscInt nconv, PetscReal *sigma, PetscReal* errest, PetscInt nest, void *mctx)
43: + svd - singular value solver context obtained from SVDCreate()
44: . its - iteration number
45: . nconv - number of converged singular triplets
46: . sigma - singular values
47: . errest - relative error estimates for each singular triplet
48: . nest - number of error estimates
49: - mctx - optional monitoring context, as set by SVDMonitorSet()
51: Options Database Keys:
52: + -svd_monitor - print only the first error estimate
53: . -svd_monitor_all - print error estimates at each iteration
54: . -svd_monitor_conv - print the singular value approximations only when
55: convergence has been reached
56: . -svd_monitor_draw - sets line graph monitor for the first unconverged
57: approximate singular value
58: . -svd_monitor_draw_all - sets line graph monitor for all unconverged
59: approximate singular value
60: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
61: a code by calls to SVDMonitorSet(), but does not cancel those set via
62: the options database.
64: Notes:
65: Several different monitoring routines may be set by calling
66: SVDMonitorSet() multiple times; all will be called in the
67: order in which they were set.
69: Level: intermediate
71: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorLG(), SVDMonitorLGAll(), SVDMonitorCancel()
72: @*/
73: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),
74: void *mctx,PetscErrorCode (*monitordestroy)(void*))
75: {
78: if (svd->numbermonitors >= MAXSVDMONITORS) {
79: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
80: }
81: svd->monitor[svd->numbermonitors] = monitor;
82: svd->monitorcontext[svd->numbermonitors] = (void*)mctx;
83: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
84: return(0);
85: }
89: /*@
90: SVDMonitorCancel - Clears all monitors for an SVD object.
92: Collective on SVD
94: Input Parameters:
95: . svd - singular value solver context obtained from SVDCreate()
97: Options Database Key:
98: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
99: into a code by calls to SVDMonitorSet(),
100: but does not cancel those set via the options database.
102: Level: intermediate
104: .seealso: SVDMonitorSet()
105: @*/
106: PetscErrorCode SVDMonitorCancel(SVD svd)
107: {
109: PetscInt i;
113: for (i=0; i<svd->numbermonitors; i++) {
114: if (svd->monitordestroy[i]) {
115: (*svd->monitordestroy[i])(svd->monitorcontext[i]);
116: }
117: }
118: svd->numbermonitors = 0;
119: return(0);
120: }
124: /*@C
125: SVDGetMonitorContext - Gets the monitor context, as set by
126: SVDMonitorSet() for the FIRST monitor only.
128: Not Collective
130: Input Parameter:
131: . svd - singular value solver context obtained from SVDCreate()
133: Output Parameter:
134: . ctx - monitor context
136: Level: intermediate
138: .seealso: SVDMonitorSet()
139: @*/
140: PetscErrorCode SVDGetMonitorContext(SVD svd, void **ctx)
141: {
144: *ctx = (svd->monitorcontext[0]);
145: return(0);
146: }
150: /*@C
151: SVDMonitorAll - Print the current approximate values and
152: error estimates at each iteration of the singular value solver.
154: Collective on SVD
156: Input Parameters:
157: + svd - singular value solver context
158: . its - iteration number
159: . nconv - number of converged singular triplets so far
160: . sigma - singular values
161: . errest - error estimates
162: . nest - number of error estimates to display
163: - dummy - unused monitor context
165: Level: intermediate
167: .seealso: SVDMonitorSet()
168: @*/
169: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *dummy)
170: {
171: PetscErrorCode ierr;
172: PetscInt i;
173: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
176: if (its) {
177: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)svd)->comm,"stdout",0,&viewer);}
178: PetscViewerASCIIMonitorPrintf(viewer,"%3d SVD nconv=%d Values (Errors)",its,nconv);
179: for (i=0;i<nest;i++) {
180: PetscViewerASCIIMonitorPrintf(viewer," %g (%10.8e)",sigma[i],errest[i]);
181: }
182: PetscViewerASCIIMonitorPrintf(viewer,"\n");
183: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
184: }
185: return(0);
186: }
190: /*@C
191: SVDMonitorFirst - Print the first unconverged approximate values and
192: error estimates at each iteration of the singular value solver.
194: Collective on SVD
196: Input Parameters:
197: + svd - singular value solver context
198: . its - iteration number
199: . nconv - number of converged singular triplets so far
200: . sigma - singular values
201: . errest - error estimates
202: . nest - number of error estimates to display
203: - dummy - unused monitor context
205: Level: intermediate
207: .seealso: SVDMonitorSet()
208: @*/
209: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *dummy)
210: {
211: PetscErrorCode ierr;
212: PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
215: if (its && nconv<nest) {
216: if (!dummy) {PetscViewerASCIIMonitorCreate(((PetscObject)svd)->comm,"stdout",0,&viewer);}
217: PetscViewerASCIIMonitorPrintf(viewer,"%3d SVD nconv=%d first unconverged value (error)",its,nconv);
218: PetscViewerASCIIMonitorPrintf(viewer," %g (%10.8e)\n",sigma[nconv],errest[nconv]);
219: if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
220: }
221: return(0);
222: }
226: /*@C
227: SVDMonitorConverged - Print the approximate values and error estimates as they converge.
229: Collective on SVD
231: Input Parameters:
232: + svd - singular value solver context
233: . its - iteration number
234: . nconv - number of converged singular triplets so far
235: . sigma - singular values
236: . errest - error estimates
237: . nest - number of error estimates to display
238: - dummy - unused monitor context
240: Level: intermediate
242: .seealso: SVDMonitorSet()
243: @*/
244: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *dummy)
245: {
246: PetscErrorCode ierr;
247: PetscInt i;
248: SVDMONITOR_CONV *ctx = (SVDMONITOR_CONV*) dummy;
251: if (!its) {
252: ctx->oldnconv = 0;
253: } else {
254: for (i=ctx->oldnconv;i<nconv;i++) {
255: PetscViewerASCIIMonitorPrintf(ctx->viewer,"%3d SVD converged value (error) #%d",its,i);
256: PetscViewerASCIIMonitorPrintf(ctx->viewer," %g (%10.8e)\n",sigma[i],errest[i]);
257: }
258: ctx->oldnconv = nconv;
259: }
260: return(0);
261: }
263: PetscErrorCode SVDMonitorDestroy_Converged(SVDMONITOR_CONV *ctx)
264: {
265: PetscErrorCode ierr;
267: PetscViewerASCIIMonitorDestroy(ctx->viewer);
268: PetscFree(ctx);
269: return(0);
270: }
274: PetscErrorCode SVDMonitorLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
275: {
276: PetscViewer viewer = (PetscViewer) monctx;
277: PetscDraw draw;
278: PetscDrawLG lg;
280: PetscReal x,y,p;
281: PetscDraw draw1;
282: PetscDrawLG lg1;
286: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); }
288: PetscViewerDrawGetDraw(viewer,0,&draw);
289: PetscViewerDrawGetDrawLG(viewer,0,&lg);
290: PetscViewerDrawGetDraw(viewer,1,&draw1);
291: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
293: if (!its) {
294: PetscDrawSetTitle(draw,"Error estimates");
295: PetscDrawSetDoubleBuffer(draw);
296: PetscDrawLGSetDimension(lg,1);
297: PetscDrawLGReset(lg);
298: PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);
300: PetscDrawSetTitle(draw1,"Approximate singular values");
301: PetscDrawSetDoubleBuffer(draw1);
302: PetscDrawLGSetDimension(lg1,1);
303: PetscDrawLGReset(lg1);
304: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
305: }
307: x = (PetscReal) its;
308: if (errest[nconv] > 0.0) y = log10(errest[nconv]); else y = 0.0;
309: PetscDrawLGAddPoint(lg,&x,&y);
311: PetscDrawLGAddPoint(lg1,&x,svd->sigma);
312: PetscDrawGetPause(draw1,&p);
313: PetscDrawSetPause(draw1,0);
314: PetscDrawLGDraw(lg1);
315: PetscDrawSetPause(draw1,p);
316:
317: PetscDrawLGDraw(lg);
318: return(0);
319: }
323: PetscErrorCode SVDMonitorLGAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *monctx)
324: {
325: PetscViewer viewer = (PetscViewer) monctx;
326: PetscDraw draw;
327: PetscDrawLG lg;
329: PetscReal *x,*y,p;
330: int n = PetscMin(svd->nsv,255);
331: PetscInt i;
332: PetscDraw draw1;
333: PetscDrawLG lg1;
337: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(((PetscObject)svd)->comm); }
339: PetscViewerDrawGetDraw(viewer,0,&draw);
340: PetscViewerDrawGetDrawLG(viewer,0,&lg);
341: PetscViewerDrawGetDraw(viewer,1,&draw1);
342: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
344: if (!its) {
345: PetscDrawSetTitle(draw,"Error estimates");
346: PetscDrawSetDoubleBuffer(draw);
347: PetscDrawLGSetDimension(lg,n);
348: PetscDrawLGReset(lg);
349: PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);
351: PetscDrawSetTitle(draw1,"Approximate singular values");
352: PetscDrawSetDoubleBuffer(draw1);
353: PetscDrawLGSetDimension(lg1,n);
354: PetscDrawLGReset(lg1);
355: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
356: }
358: PetscMalloc(sizeof(PetscReal)*n,&x);
359: PetscMalloc(sizeof(PetscReal)*n,&y);
360: for (i=0;i<n;i++) {
361: x[i] = (PetscReal) its;
362: if (i < nest && errest[i] > 0.0) y[i] = log10(errest[i]);
363: else y[i] = 0.0;
364: }
365: PetscDrawLGAddPoint(lg,x,y);
367: PetscDrawLGAddPoint(lg1,x,svd->sigma);
368: PetscDrawGetPause(draw1,&p);
369: PetscDrawSetPause(draw1,0);
370: PetscDrawLGDraw(lg1);
371: PetscDrawSetPause(draw1,p);
372:
373: PetscDrawLGDraw(lg);
374: PetscFree(x);
375: PetscFree(y);
376: return(0);
377: }