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-2007, Universidad Politecnica de Valencia, Spain
8: This file is part of SLEPc. See the README file for conditions of use
9: and additional information.
10: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
11: */
13: #include src/svd/svdimpl.h
17: /*@C
18: SVDMonitorSet - Sets an ADDITIONAL function to be called at every
19: iteration to monitor the error estimates for each requested singular triplet.
20:
21: Collective on SVD
23: Input Parameters:
24: + svd - singular value solver context obtained from SVDCreate()
25: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring)
26: - mctx - [optional] context for private data for the
27: monitor routine (use PETSC_NULL if no context is desired)
29: Calling Sequence of monitor:
30: $ monitor (SVD svd, int its, int nconv, PetscReal *sigma, PetscReal* errest, int nest, void *mctx)
32: + svd - singular value solver context obtained from SVDCreate()
33: . its - iteration number
34: . nconv - number of converged singular triplets
35: . sigma - singular values
36: . errest - relative error estimates for each singular triplet
37: . nest - number of error estimates
38: - mctx - optional monitoring context, as set by SVDMonitorSet()
40: Options Database Keys:
41: + -svd_monitor - print error estimates at each iteration
42: . -svd_monitor_draw - sets line graph monitor
43: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
44: a code by calls to SVDMonitorSet(), but does not cancel those set via
45: the options database.
47: Notes:
48: Several different monitoring routines may be set by calling
49: SVDMonitorSet() multiple times; all will be called in the
50: order in which they were set.
52: Level: intermediate
54: .seealso: SVDMonitorDefault(), SVDMonitorCancel()
55: @*/
56: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,int,int,PetscReal*,PetscReal*,int,void*),
57: void *mctx,PetscErrorCode (*monitordestroy)(void*))
58: {
61: if (svd->numbermonitors >= MAXSVDMONITORS) {
62: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
63: }
64: svd->monitor[svd->numbermonitors] = monitor;
65: svd->monitorcontext[svd->numbermonitors] = (void*)mctx;
66: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
67: return(0);
68: }
72: /*@
73: SVDMonitorCancel - Clears all monitors for an SVD object.
75: Collective on SVD
77: Input Parameters:
78: . svd - singular value solver context obtained from SVDCreate()
80: Options Database Key:
81: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
82: into a code by calls to SVDMonitorSet(),
83: but does not cancel those set via the options database.
85: Level: intermediate
87: .seealso: SVDMonitorCancel()
88: @*/
89: PetscErrorCode SVDMonitorCancel(SVD svd)
90: {
92: PetscInt i;
96: for (i=0; i<svd->numbermonitors; i++) {
97: if (svd->monitordestroy[i]) {
98: (*svd->monitordestroy[i])(svd->monitorcontext[i]);
99: }
100: }
101: svd->numbermonitors = 0;
102: return(0);
103: }
107: /*@C
108: SVDGetMonitorContext - Gets the monitor context, as set by
109: SVDMonitorSet() for the FIRST monitor only.
111: Not Collective
113: Input Parameter:
114: . svd - singular value solver context obtained from SVDCreate()
116: Output Parameter:
117: . ctx - monitor context
119: Level: intermediate
121: .seealso: SVDMonitorSet(), SVDMonitorDefault()
122: @*/
123: PetscErrorCode SVDGetMonitorContext(SVD svd, void **ctx)
124: {
127: *ctx = (svd->monitorcontext[0]);
128: return(0);
129: }
133: /*@C
134: SVDDefaultMonitor - Print the current approximate values and
135: error estimates at each iteration of the singular value solver.
137: Collective on SVD
139: Input Parameters:
140: + svd - singular value solver context
141: . its - iteration number
142: . nconv - number of converged singular triplets so far
143: . sigma - singular values
144: . errest - error estimates
145: . nest - number of error estimates to display
146: - dummy - unused monitor context
148: Level: intermediate
150: .seealso: SVDMonitorSet()
151: @*/
152: PetscErrorCode SVDMonitorDefault(SVD svd,int its,int nconv,PetscReal *sigma,PetscReal *errest,int nest,void *dummy)
153: {
155: int i;
156: PetscViewer viewer = (PetscViewer) dummy;
159: if (its) {
160: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(svd->comm);
161: PetscViewerASCIIPrintf(viewer,"%3d SVD nconv=%d Values (Errors)",its,nconv);
162: for (i=0;i<nest;i++) {
163: PetscViewerASCIIPrintf(viewer," %g (%10.8e)",sigma[i],errest[i]);
164: }
165: PetscViewerASCIIPrintf(viewer,"\n");
166: }
167: return(0);
168: }
172: PetscErrorCode SVDMonitorLG(SVD svd,int its,int nconv,PetscReal *sigma,PetscReal *errest,int nest,void *monctx)
173: {
174: PetscViewer viewer = (PetscViewer) monctx;
175: PetscDraw draw;
176: PetscDrawLG lg;
178: PetscReal *x,*y;
179: int i,n = svd->nsv;
180: int p;
181: PetscDraw draw1;
182: PetscDrawLG lg1;
186: if (!viewer) { viewer = PETSC_VIEWER_DRAW_(svd->comm); }
188: PetscViewerDrawGetDraw(viewer,0,&draw);
189: PetscViewerDrawGetDrawLG(viewer,0,&lg);
190: PetscViewerDrawGetDraw(viewer,1,&draw1);
191: PetscViewerDrawGetDrawLG(viewer,1,&lg1);
193: if (!its) {
194: PetscDrawSetTitle(draw,"Error estimates");
195: PetscDrawSetDoubleBuffer(draw);
196: PetscDrawLGSetDimension(lg,n);
197: PetscDrawLGReset(lg);
198: PetscDrawLGSetLimits(lg,0,1.0,log10(svd->tol)-2,0.0);
200: PetscDrawSetTitle(draw1,"Approximate singular values");
201: PetscDrawSetDoubleBuffer(draw1);
202: PetscDrawLGSetDimension(lg1,n);
203: PetscDrawLGReset(lg1);
204: PetscDrawLGSetLimits(lg1,0,1.0,1.e20,-1.e20);
205: }
207: PetscMalloc(sizeof(PetscReal)*n,&x);
208: PetscMalloc(sizeof(PetscReal)*n,&y);
209: for (i=0;i<n;i++) {
210: x[i] = (PetscReal) its;
211: if (errest[i] > 0.0) y[i] = log10(errest[i]); else y[i] = 0.0;
212: }
213: PetscDrawLGAddPoint(lg,x,y);
215: PetscDrawLGAddPoint(lg1,x,svd->sigma);
216: PetscDrawGetPause(draw1,&p);
217: PetscDrawSetPause(draw1,0);
218: PetscDrawLGDraw(lg1);
219: PetscDrawSetPause(draw1,p);
220:
221: PetscDrawLGDraw(lg);
222: PetscFree(x);
223: PetscFree(y);
224: return(0);
225: }