1: /*
2: SVD routines related to monitors.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
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*/
25: #include <petscdraw.h>
29: /*
30: Runs the user provided monitor routines, if any.
31: */
32: PetscErrorCode SVDMonitor(SVD svd,PetscInt it,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest) 33: {
35: PetscInt i,n = svd->numbermonitors;
38: for (i=0;i<n;i++) {
39: (*svd->monitor[i])(svd,it,nconv,sigma,errest,nest,svd->monitorcontext[i]);
40: }
41: return(0);
42: }
46: /*@C
47: SVDMonitorSet - Sets an ADDITIONAL function to be called at every
48: iteration to monitor the error estimates for each requested singular triplet.
50: Collective on SVD 52: Input Parameters:
53: + svd - singular value solver context obtained from SVDCreate()
54: . monitor - pointer to function (if this is NULL, it turns off monitoring)
55: - mctx - [optional] context for private data for the
56: monitor routine (use NULL if no context is desired)
58: Calling Sequence of monitor:
59: $ monitor(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *mctx)
61: + svd - singular value solver context obtained from SVDCreate()
62: . its - iteration number
63: . nconv - number of converged singular triplets
64: . sigma - singular values
65: . errest - relative error estimates for each singular triplet
66: . nest - number of error estimates
67: - mctx - optional monitoring context, as set by SVDMonitorSet()
69: Options Database Keys:
70: + -svd_monitor - print only the first error estimate
71: . -svd_monitor_all - print error estimates at each iteration
72: . -svd_monitor_conv - print the singular value approximations only when
73: convergence has been reached
74: . -svd_monitor_lg - sets line graph monitor for the first unconverged
75: approximate singular value
76: . -svd_monitor_lg_all - sets line graph monitor for all unconverged
77: approximate singular values
78: - -svd_monitor_cancel - cancels all monitors that have been hardwired into
79: a code by calls to SVDMonitorSet(), but does not cancel those set via
80: the options database.
82: Notes:
83: Several different monitoring routines may be set by calling
84: SVDMonitorSet() multiple times; all will be called in the
85: order in which they were set.
87: Level: intermediate
89: .seealso: SVDMonitorFirst(), SVDMonitorAll(), SVDMonitorCancel()
90: @*/
91: PetscErrorCode SVDMonitorSet(SVD svd,PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**)) 92: {
95: if (svd->numbermonitors >= MAXSVDMONITORS) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Too many SVD monitors set");
96: svd->monitor[svd->numbermonitors] = monitor;
97: svd->monitorcontext[svd->numbermonitors] = (void*)mctx;
98: svd->monitordestroy[svd->numbermonitors++] = monitordestroy;
99: return(0);
100: }
104: /*@
105: SVDMonitorCancel - Clears all monitors for an SVD object.
107: Collective on SVD109: Input Parameters:
110: . svd - singular value solver context obtained from SVDCreate()
112: Options Database Key:
113: . -svd_monitor_cancel - Cancels all monitors that have been hardwired
114: into a code by calls to SVDMonitorSet(),
115: but does not cancel those set via the options database.
117: Level: intermediate
119: .seealso: SVDMonitorSet()
120: @*/
121: PetscErrorCode SVDMonitorCancel(SVD svd)122: {
124: PetscInt i;
128: for (i=0; i<svd->numbermonitors; i++) {
129: if (svd->monitordestroy[i]) {
130: (*svd->monitordestroy[i])(&svd->monitorcontext[i]);
131: }
132: }
133: svd->numbermonitors = 0;
134: return(0);
135: }
139: /*@C
140: SVDGetMonitorContext - Gets the monitor context, as set by
141: SVDMonitorSet() for the FIRST monitor only.
143: Not Collective
145: Input Parameter:
146: . svd - singular value solver context obtained from SVDCreate()
148: Output Parameter:
149: . ctx - monitor context
151: Level: intermediate
153: .seealso: SVDMonitorSet()
154: @*/
155: PetscErrorCode SVDGetMonitorContext(SVD svd,void **ctx)156: {
159: *ctx = svd->monitorcontext[0];
160: return(0);
161: }
165: /*@C
166: SVDMonitorAll - Print the current approximate values and
167: error estimates at each iteration of the singular value solver.
169: Collective on SVD171: Input Parameters:
172: + svd - singular value solver context
173: . its - iteration number
174: . nconv - number of converged singular triplets so far
175: . sigma - singular values
176: . errest - error estimates
177: . nest - number of error estimates to display
178: - vf - viewer and format for monitoring
180: Level: intermediate
182: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorConverged()
183: @*/
184: PetscErrorCode SVDMonitorAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)185: {
187: PetscInt i;
188: PetscViewer viewer;
193: viewer = vf->viewer;
195: PetscViewerPushFormat(viewer,vf->format);
196: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
197: if (its==1 && ((PetscObject)svd)->prefix) {
198: PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
199: }
200: PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D Values (Errors)",its,nconv);
201: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
202: for (i=0;i<nest;i++) {
203: PetscViewerASCIIPrintf(viewer," %g (%10.8e)",(double)sigma[i],(double)errest[i]);
204: }
205: PetscViewerASCIIPrintf(viewer,"\n");
206: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
207: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
208: PetscViewerPopFormat(viewer);
209: return(0);
210: }
214: /*@C
215: SVDMonitorFirst - Print the first unconverged approximate values and
216: error estimates at each iteration of the singular value solver.
218: Collective on SVD220: Input Parameters:
221: + svd - singular value solver context
222: . its - iteration number
223: . nconv - number of converged singular triplets so far
224: . sigma - singular values
225: . errest - error estimates
226: . nest - number of error estimates to display
227: - vf - viewer and format for monitoring
229: Level: intermediate
231: .seealso: SVDMonitorSet(), SVDMonitorAll(), SVDMonitorConverged()
232: @*/
233: PetscErrorCode SVDMonitorFirst(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)234: {
236: PetscViewer viewer;
241: viewer = vf->viewer;
243: if (its==1 && ((PetscObject)svd)->prefix) {
244: PetscViewerASCIIPrintf(viewer," Singular value approximations and residual norms for %s solve.\n",((PetscObject)svd)->prefix);
245: }
246: if (nconv<nest) {
247: PetscViewerPushFormat(viewer,vf->format);
248: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
249: PetscViewerASCIIPrintf(viewer,"%3D SVD nconv=%D first unconverged value (error)",its,nconv);
250: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
251: PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[nconv],(double)errest[nconv]);
252: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
253: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
254: PetscViewerPopFormat(viewer);
255: }
256: return(0);
257: }
261: /*@C
262: SVDMonitorConverged - Print the approximate values and error estimates as they converge.
264: Collective on SVD266: Input Parameters:
267: + svd - singular value solver context
268: . its - iteration number
269: . nconv - number of converged singular triplets so far
270: . sigma - singular values
271: . errest - error estimates
272: . nest - number of error estimates to display
273: - ctx - monitor context
275: Level: intermediate
277: .seealso: SVDMonitorSet(), SVDMonitorFirst(), SVDMonitorAll()
278: @*/
279: PetscErrorCode SVDMonitorConverged(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,SlepcConvMonitor ctx)280: {
282: PetscInt i;
283: PetscViewer viewer;
288: viewer = ctx->viewer;
290: if (its==1 && ((PetscObject)svd)->prefix) {
291: PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)svd)->prefix);
292: }
293: if (its==1) ctx->oldnconv = 0;
294: if (ctx->oldnconv!=nconv) {
295: PetscViewerPushFormat(viewer,ctx->format);
296: PetscViewerASCIIAddTab(viewer,((PetscObject)svd)->tablevel);
297: for (i=ctx->oldnconv;i<nconv;i++) {
298: PetscViewerASCIIPrintf(viewer,"%3D SVD converged value (error) #%D",its,i);
299: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
300: PetscViewerASCIIPrintf(viewer," %g (%10.8e)\n",(double)sigma[i],(double)errest[i]);
301: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
302: }
303: PetscViewerASCIISubtractTab(viewer,((PetscObject)svd)->tablevel);
304: PetscViewerPopFormat(viewer);
305: ctx->oldnconv = nconv;
306: }
307: return(0);
308: }
312: /*@C
313: SVDMonitorLGCreate - Creates a line graph context for use with
314: SVD to monitor convergence.
316: Collective on MPI_Comm
318: Input Parameters:
319: + comm - communicator context
320: . host - the X display to open, or null for the local machine
321: . label - the title to put in the title bar
322: . x, y - the screen coordinates of the upper left coordinate of
323: the window
324: - m, n - the screen width and height in pixels
326: Output Parameter:
327: . lgctx - the drawing context
329: Options Database Keys:
330: + -svd_monitor_lg - Sets line graph monitor for the first residual
331: - -svd_monitor_lg_all - Sets line graph monitor for all residuals
333: Notes:
334: Use PetscDrawLGDestroy() to destroy this line graph.
336: Level: intermediate
338: .seealso: SVDMonitorSet()
339: @*/
340: PetscErrorCode SVDMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)341: {
342: PetscDraw draw;
343: PetscDrawLG lg;
347: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
348: PetscDrawSetFromOptions(draw);
349: PetscDrawLGCreate(draw,1,&lg);
350: PetscDrawLGSetFromOptions(lg);
351: PetscDrawDestroy(&draw);
352: *lgctx = lg;
353: return(0);
354: }
358: PetscErrorCode SVDMonitorLG(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *ctx)359: {
360: PetscDrawLG lg = (PetscDrawLG)ctx;
361: PetscReal x,y;
366: if (its==1) {
367: PetscDrawLGReset(lg);
368: PetscDrawLGSetDimension(lg,1);
369: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(svd->tol)-2,0.0);
370: }
371: x = (PetscReal)its;
372: if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
373: else y = 0.0;
374: PetscDrawLGAddPoint(lg,&x,&y);
375: if (its <= 20 || !(its % 5) || svd->reason) {
376: PetscDrawLGDraw(lg);
377: PetscDrawLGSave(lg);
378: }
379: return(0);
380: }
384: PetscErrorCode SVDMonitorLGAll(SVD svd,PetscInt its,PetscInt nconv,PetscReal *sigma,PetscReal *errest,PetscInt nest,void *ctx)385: {
386: PetscDrawLG lg = (PetscDrawLG)ctx;
387: PetscInt i,n = PetscMin(svd->nsv,255);
388: PetscReal *x,*y;
393: if (its==1) {
394: PetscDrawLGReset(lg);
395: PetscDrawLGSetDimension(lg,n);
396: PetscDrawLGSetLimits(lg,1,1.0,PetscLog10Real(svd->tol)-2,0.0);
397: }
398: PetscMalloc2(n,&x,n,&y);
399: for (i=0;i<n;i++) {
400: x[i] = (PetscReal)its;
401: if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
402: else y[i] = 0.0;
403: }
404: PetscDrawLGAddPoint(lg,x,y);
405: if (its <= 20 || !(its % 5) || svd->reason) {
406: PetscDrawLGDraw(lg);
407: PetscDrawLGSave(lg);
408: }
409: PetscFree2(x,y);
410: return(0);
411: }