Actual source code: xmon.c
petsc-3.3-p7 2013-05-11
2: #include <petsc-private/kspimpl.h> /*I "petscksp.h" I*/
6: /*@C
7: KSPMonitorLGCreate - Creates a line graph context for use with
8: KSP to monitor convergence of preconditioned residual norms.
10: Collective on KSP
12: Input Parameters:
13: + host - the X display to open, or null for the local machine
14: . label - the title to put in the title bar
15: . x, y - the screen coordinates of the upper left coordinate of
16: the window
17: - m, n - the screen width and height in pixels
19: Output Parameter:
20: . draw - the drawing context
22: Options Database Key:
23: . -ksp_monitor_draw - Sets line graph monitor
25: Notes:
26: Use KSPMonitorLGDestroy() to destroy this line graph; do not use PetscDrawLGDestroy().
28: Level: intermediate
30: .keywords: KSP, monitor, line graph, residual, create
32: .seealso: KSPMonitorLGDestroy(), KSPMonitorSet(), KSPMonitorLGTrueResidualCreate()
33: @*/
34: PetscErrorCode KSPMonitorLGCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
35: {
36: PetscDraw win;
40: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
41: PetscDrawSetType(win,PETSC_DRAW_X);
42: PetscDrawLGCreate(win,1,draw);
43: PetscLogObjectParent(*draw,win);
44: return(0);
45: }
49: PetscErrorCode KSPMonitorLG(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
50: {
51: PetscDrawLG lg = (PetscDrawLG)monctx;
53: PetscReal x,y;
54: PetscViewer v;
57: if (!monctx) {
58: MPI_Comm comm;
60: PetscObjectGetComm((PetscObject)ksp,&comm);
61: v = PETSC_VIEWER_DRAW_(comm);
62: PetscViewerDrawGetDrawLG(v,0,&lg);
63: }
65: if (!n) {PetscDrawLGReset(lg);}
66: x = (PetscReal) n;
67: if (rnorm > 0.0) y = log10(rnorm); else y = -15.0;
68: PetscDrawLGAddPoint(lg,&x,&y);
69: if (n < 20 || !(n % 5)) {
70: PetscDrawLGDraw(lg);
71: }
72: return(0);
73: }
74:
77: /*@
78: KSPMonitorLGDestroy - Destroys a line graph context that was created
79: with KSPMonitorLGCreate().
81: Collective on KSP
83: Input Parameter:
84: . draw - the drawing context
86: Level: intermediate
88: .keywords: KSP, monitor, line graph, destroy
90: .seealso: KSPMonitorLGCreate(), KSPMonitorLGTrueResidualDestroy(), KSPMonitorSet()
91: @*/
92: PetscErrorCode KSPMonitorLGDestroy(PetscDrawLG *drawlg)
93: {
94: PetscDraw draw;
98: PetscDrawLGGetDraw(*drawlg,&draw);
99: PetscDrawDestroy(&draw);
100: PetscDrawLGDestroy(drawlg);
101: return(0);
102: }
106: /*@C
107: KSPMonitorLGRangeCreate - Creates a line graph context for use with
108: KSP to monitor convergence of preconditioned residual norms and range of residual element values.
110: Collective on KSP
112: Input Parameters:
113: + host - the X display to open, or null for the local machine
114: . label - the title to put in the title bar
115: . x, y - the screen coordinates of the upper left coordinate of
116: the window
117: - m, n - the screen width and height in pixels
119: Output Parameter:
120: . draw - the drawing context
122: Options Database Key:
123: . -ksp_monitor_range_draw - Sets line graph monitor
125: Notes:
126: Use KSPMonitorLGDestroy() to destroy this line graph; do not use PetscDrawLGDestroy().
128: Level: intermediate
130: .keywords: KSP, monitor, line graph, residual, create
132: .seealso: KSPMonitorLGDestroy(), KSPMonitorSet(), KSPMonitorLGTrueResidualCreate()
133: @*/
134: PetscErrorCode KSPMonitorLGRangeCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
135: {
136: PetscDraw win;
140: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
141: PetscDrawSetType(win,PETSC_DRAW_X);
142: PetscDrawLGCreate(win,2,draw);
143: PetscLogObjectParent(*draw,win);
144: return(0);
145: }
147: extern PetscErrorCode KSPMonitorRange_Private(KSP,PetscInt,PetscReal*);
150: PetscErrorCode KSPMonitorLGRange(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
151: {
152: PetscDrawLG lg;
153: PetscErrorCode ierr;
154: PetscReal x,y,per;
155: PetscViewer v = (PetscViewer)monctx;
156: static PetscReal prev; /* should be in the context */
157: PetscDraw draw;
160: if (!monctx) {
161: MPI_Comm comm;
163: PetscObjectGetComm((PetscObject)ksp,&comm);
164: v = PETSC_VIEWER_DRAW_(comm);
165: }
166: PetscViewerDrawGetDrawLG(v,0,&lg);
167: if (!n) {PetscDrawLGReset(lg);}
168: PetscDrawLGGetDraw(lg,&draw);
169: PetscDrawSetTitle(draw,"Residual norm");
170: x = (PetscReal) n;
171: if (rnorm > 0.0) y = log10(rnorm); else y = -15.0;
172: PetscDrawLGAddPoint(lg,&x,&y);
173: if (n < 20 || !(n % 5)) {
174: PetscDrawLGDraw(lg);
175: }
177: PetscViewerDrawGetDrawLG(v,1,&lg);
178: KSPMonitorRange_Private(ksp,n,&per);
179: if (!n) {PetscDrawLGReset(lg);}
180: PetscDrawLGGetDraw(lg,&draw);
181: PetscDrawSetTitle(draw,"% elemts > .2*max elemt");
182: x = (PetscReal) n;
183: y = 100.0*per;
184: PetscDrawLGAddPoint(lg,&x,&y);
185: if (n < 20 || !(n % 5)) {
186: PetscDrawLGDraw(lg);
187: }
189: PetscViewerDrawGetDrawLG(v,2,&lg);
190: if (!n) {prev = rnorm;PetscDrawLGReset(lg);}
191: PetscDrawLGGetDraw(lg,&draw);
192: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
193: PetscDrawLGGetDraw(lg,&draw);
194: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm");
195: x = (PetscReal) n;
196: y = (prev - rnorm)/prev;
197: PetscDrawLGAddPoint(lg,&x,&y);
198: if (n < 20 || !(n % 5)) {
199: PetscDrawLGDraw(lg);
200: }
202: PetscViewerDrawGetDrawLG(v,3,&lg);
203: if (!n) {PetscDrawLGReset(lg);}
204: PetscDrawLGGetDraw(lg,&draw);
205: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
206: x = (PetscReal) n;
207: y = (prev - rnorm)/(prev*per);
208: if (n > 5) {
209: PetscDrawLGAddPoint(lg,&x,&y);
210: }
211: if (n < 20 || !(n % 5)) {
212: PetscDrawLGDraw(lg);
213: }
214: prev = rnorm;
215: return(0);
216: }
217:
220: /*@
221: KSPMonitorLGRangeDestroy - Destroys a line graph context that was created
222: with KSPMonitorLGRangeCreate().
224: Collective on KSP
226: Input Parameter:
227: . draw - the drawing context
229: Level: intermediate
231: .keywords: KSP, monitor, line graph, destroy
233: .seealso: KSPMonitorLGCreate(), KSPMonitorLGTrueResidualDestroy(), KSPMonitorSet()
234: @*/
235: PetscErrorCode KSPMonitorLGRangeDestroy(PetscDrawLG *drawlg)
236: {
237: PetscDraw draw;
241: PetscDrawLGGetDraw(*drawlg,&draw);
242: PetscDrawDestroy(&draw);
243: PetscDrawLGDestroy(drawlg);
244: return(0);
245: }
249: /*@C
250: KSPMonitorLGTrueResidualNormCreate - Creates a line graph context for use with
251: KSP to monitor convergence of true residual norms (as opposed to
252: preconditioned residual norms).
254: Collective on KSP
256: Input Parameters:
257: + host - the X display to open, or null for the local machine
258: . label - the title to put in the title bar
259: . x, y - the screen coordinates of the upper left coordinate of
260: the window
261: - m, n - the screen width and height in pixels
263: Output Parameter:
264: . draw - the drawing context
266: Options Database Key:
267: . -ksp_monitor_draw_true_residual - Sets true line graph monitor
269: Notes:
270: Use KSPMonitorLGTrueResidualNormDestroy() to destroy this line graph, not
271: PetscDrawLGDestroy().
273: Level: intermediate
275: .keywords: KSP, monitor, line graph, residual, create, true
277: .seealso: KSPMonitorLGDestroy(), KSPMonitorSet(), KSPMonitorDefault()
278: @*/
279: PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
280: {
281: PetscDraw win;
283: PetscMPIInt rank;
286: MPI_Comm_rank(comm,&rank);
287: if (rank) { *draw = 0; return(0);}
289: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
290: PetscDrawSetType(win,PETSC_DRAW_X);
291: PetscDrawLGCreate(win,2,draw);
292: PetscLogObjectParent(*draw,win);
293: return(0);
294: }
298: PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
299: {
300: PetscDrawLG lg = (PetscDrawLG) monctx;
301: PetscReal x[2],y[2],scnorm;
303: PetscMPIInt rank;
304: Vec resid,work;
307: if (!monctx) {
308: MPI_Comm comm;
309: PetscViewer viewer;
311: PetscObjectGetComm((PetscObject)ksp,&comm);
312: viewer = PETSC_VIEWER_DRAW_(comm);
313: PetscViewerDrawGetDrawLG(viewer,0,&lg);
314: }
316: MPI_Comm_rank(((PetscObject)ksp)->comm,&rank);
317: if (!rank) {
318: if (!n) {PetscDrawLGReset(lg);}
319: x[0] = x[1] = (PetscReal) n;
320: if (rnorm > 0.0) y[0] = log10(rnorm); else y[0] = -15.0;
321: }
323: VecDuplicate(ksp->vec_rhs,&work);
324: KSPBuildResidual(ksp,0,work,&resid);
325: VecNorm(resid,NORM_2,&scnorm);
326: VecDestroy(&work);
328: if (!rank) {
329: if (scnorm > 0.0) y[1] = log10(scnorm); else y[1] = -15.0;
330: PetscDrawLGAddPoint(lg,x,y);
331: if (n <= 20 || (n % 3)) {
332: PetscDrawLGDraw(lg);
333: }
334: }
335: return(0);
336: }
337:
340: /*@C
341: KSPMonitorLGTrueResidualNormDestroy - Destroys a line graph context that was created
342: with KSPMonitorLGTrueResidualNormCreate().
344: Collective on KSP
346: Input Parameter:
347: . draw - the drawing context
349: Level: intermediate
351: .keywords: KSP, monitor, line graph, destroy, true
353: .seealso: KSPMonitorLGTrueResidualNormCreate(), KSPMonitorSet()
354: @*/
355: PetscErrorCode KSPMonitorLGTrueResidualNormDestroy(PetscDrawLG *drawlg)
356: {
358: PetscDraw draw;
361: PetscDrawLGGetDraw(*drawlg,&draw);
362: PetscDrawDestroy(&draw);
363: PetscDrawLGDestroy(drawlg);
364: return(0);
365: }