Actual source code: xmon.c
1: #define PETSCKSP_DLL
3: #include private/kspimpl.h
7: /*@C
8: KSPMonitorLGCreate - Creates a line graph context for use with
9: KSP to monitor convergence of preconditioned residual norms.
11: Collective on KSP
13: Input Parameters:
14: + host - the X display to open, or null for the local machine
15: . label - the title to put in the title bar
16: . x, y - the screen coordinates of the upper left coordinate of
17: the window
18: - m, n - the screen width and height in pixels
20: Output Parameter:
21: . draw - the drawing context
23: Options Database Key:
24: . -ksp_monitor_draw - Sets line graph monitor
26: Notes:
27: Use KSPMonitorLGDestroy() to destroy this line graph; do not use PetscDrawLGDestroy().
29: Level: intermediate
31: .keywords: KSP, monitor, line graph, residual, create
33: .seealso: KSPMonitorLGDestroy(), KSPMonitorSet(), KSPMonitorLGTrueResidualCreate()
34: @*/
35: PetscErrorCode KSPMonitorLGCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
36: {
37: PetscDraw win;
41: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
42: PetscDrawSetType(win,PETSC_DRAW_X);
43: PetscDrawLGCreate(win,1,draw);
44: PetscLogObjectParent(*draw,win);
45: return(0);
46: }
50: PetscErrorCode KSPMonitorLG(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
51: {
52: PetscDrawLG lg;
54: PetscReal x,y;
55: PetscViewer v = (PetscViewer)monctx;
58: if (!monctx) {
59: MPI_Comm comm;
61: PetscObjectGetComm((PetscObject)ksp,&comm);
62: v = PETSC_VIEWER_DRAW_(comm);
63: }
64: PetscViewerDrawGetDrawLG(v,0,&lg);
66: if (!n) {PetscDrawLGReset(lg);}
67: x = (PetscReal) n;
68: if (rnorm > 0.0) y = log10(rnorm); else y = -15.0;
69: PetscDrawLGAddPoint(lg,&x,&y);
70: if (n < 20 || !(n % 5)) {
71: PetscDrawLGDraw(lg);
72: }
73: return(0);
74: }
75:
78: /*@
79: KSPMonitorLGDestroy - Destroys a line graph context that was created
80: with KSPMonitorLGCreate().
82: Collective on KSP
84: Input Parameter:
85: . draw - the drawing context
87: Level: intermediate
89: .keywords: KSP, monitor, line graph, destroy
91: .seealso: KSPMonitorLGCreate(), KSPMonitorLGTrueResidualDestroy(), KSPMonitorSet()
92: @*/
93: PetscErrorCode KSPMonitorLGDestroy(PetscDrawLG drawlg)
94: {
95: PetscDraw draw;
99: PetscDrawLGGetDraw(drawlg,&draw);
100: if (draw) { PetscDrawDestroy(draw);}
101: PetscDrawLGDestroy(drawlg);
102: return(0);
103: }
107: /*@C
108: KSPMonitorLGRangeCreate - Creates a line graph context for use with
109: KSP to monitor convergence of preconditioned residual norms and range of residual element values.
111: Collective on KSP
113: Input Parameters:
114: + host - the X display to open, or null for the local machine
115: . label - the title to put in the title bar
116: . x, y - the screen coordinates of the upper left coordinate of
117: the window
118: - m, n - the screen width and height in pixels
120: Output Parameter:
121: . draw - the drawing context
123: Options Database Key:
124: . -ksp_monitor_range_draw - Sets line graph monitor
126: Notes:
127: Use KSPMonitorLGDestroy() to destroy this line graph; do not use PetscDrawLGDestroy().
129: Level: intermediate
131: .keywords: KSP, monitor, line graph, residual, create
133: .seealso: KSPMonitorLGDestroy(), KSPMonitorSet(), KSPMonitorLGTrueResidualCreate()
134: @*/
135: PetscErrorCode KSPMonitorLGRangeCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
136: {
137: PetscDraw win;
141: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
142: PetscDrawSetType(win,PETSC_DRAW_X);
143: PetscDrawLGCreate(win,2,draw);
144: PetscLogObjectParent(*draw,win);
145: return(0);
146: }
151: PetscErrorCode KSPMonitorLGRange(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
152: {
153: PetscDrawLG lg;
154: PetscErrorCode ierr;
155: PetscReal x,y,per;
156: PetscViewer v = (PetscViewer)monctx;
157: static PetscReal prev; /* should be in the context */
158: PetscDraw draw;
161: if (!monctx) {
162: MPI_Comm comm;
164: PetscObjectGetComm((PetscObject)ksp,&comm);
165: v = PETSC_VIEWER_DRAW_(comm);
166: }
167: PetscViewerDrawGetDrawLG(v,0,&lg);
168: if (!n) {PetscDrawLGReset(lg);}
169: PetscDrawLGGetDraw(lg,&draw);
170: PetscDrawSetTitle(draw,"Residual norm");
171: x = (PetscReal) n;
172: if (rnorm > 0.0) y = log10(rnorm); else y = -15.0;
173: PetscDrawLGAddPoint(lg,&x,&y);
174: if (n < 20 || !(n % 5)) {
175: PetscDrawLGDraw(lg);
176: }
178: PetscViewerDrawGetDrawLG(v,1,&lg);
179: KSPMonitorRange_Private(ksp,n,&per);
180: if (!n) {PetscDrawLGReset(lg);}
181: PetscDrawLGGetDraw(lg,&draw);
182: PetscDrawSetTitle(draw,"% elemts > .2*max elemt");
183: x = (PetscReal) n;
184: y = 100.0*per;
185: PetscDrawLGAddPoint(lg,&x,&y);
186: if (n < 20 || !(n % 5)) {
187: PetscDrawLGDraw(lg);
188: }
190: PetscViewerDrawGetDrawLG(v,2,&lg);
191: if (!n) {prev = rnorm;PetscDrawLGReset(lg);}
192: PetscDrawLGGetDraw(lg,&draw);
193: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
194: PetscDrawLGGetDraw(lg,&draw);
195: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm");
196: x = (PetscReal) n;
197: y = (prev - rnorm)/prev;
198: PetscDrawLGAddPoint(lg,&x,&y);
199: if (n < 20 || !(n % 5)) {
200: PetscDrawLGDraw(lg);
201: }
203: PetscViewerDrawGetDrawLG(v,3,&lg);
204: if (!n) {PetscDrawLGReset(lg);}
205: PetscDrawLGGetDraw(lg,&draw);
206: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
207: x = (PetscReal) n;
208: y = (prev - rnorm)/(prev*per);
209: if (n > 5) {
210: PetscDrawLGAddPoint(lg,&x,&y);
211: }
212: if (n < 20 || !(n % 5)) {
213: PetscDrawLGDraw(lg);
214: }
215: prev = rnorm;
216: return(0);
217: }
218:
221: /*@
222: KSPMonitorLGRangeDestroy - Destroys a line graph context that was created
223: with KSPMonitorLGRangeCreate().
225: Collective on KSP
227: Input Parameter:
228: . draw - the drawing context
230: Level: intermediate
232: .keywords: KSP, monitor, line graph, destroy
234: .seealso: KSPMonitorLGCreate(), KSPMonitorLGTrueResidualDestroy(), KSPMonitorSet()
235: @*/
236: PetscErrorCode KSPMonitorLGRangeDestroy(PetscDrawLG drawlg)
237: {
238: PetscDraw draw;
242: PetscDrawLGGetDraw(drawlg,&draw);
243: if (draw) { PetscDrawDestroy(draw);}
244: PetscDrawLGDestroy(drawlg);
245: return(0);
246: }
250: /*@C
251: KSPMonitorLGTrueResidualNormCreate - Creates a line graph context for use with
252: KSP to monitor convergence of true residual norms (as opposed to
253: preconditioned residual norms).
255: Collective on KSP
257: Input Parameters:
258: + host - the X display to open, or null for the local machine
259: . label - the title to put in the title bar
260: . x, y - the screen coordinates of the upper left coordinate of
261: the window
262: - m, n - the screen width and height in pixels
264: Output Parameter:
265: . draw - the drawing context
267: Options Database Key:
268: . -ksp_monitor_draw_true_residual - Sets true line graph monitor
270: Notes:
271: Use KSPMonitorLGTrueResidualNormDestroy() to destroy this line graph, not
272: PetscDrawLGDestroy().
274: Level: intermediate
276: .keywords: KSP, monitor, line graph, residual, create, true
278: .seealso: KSPMonitorLGDestroy(), KSPMonitorSet(), KSPMonitorDefault()
279: @*/
280: PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
281: {
282: PetscDraw win;
284: PetscMPIInt rank;
287: MPI_Comm_rank(comm,&rank);
288: if (rank) { *draw = 0; return(0);}
290: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
291: PetscDrawSetType(win,PETSC_DRAW_X);
292: PetscDrawLGCreate(win,2,draw);
293: PetscLogObjectParent(*draw,win);
294: return(0);
295: }
299: PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
300: {
301: PetscDrawLG lg = (PetscDrawLG) monctx;
302: PetscReal x[2],y[2],scnorm;
304: PetscMPIInt rank;
305: Vec resid,work;
308: if (!monctx) {
309: MPI_Comm comm;
310: PetscViewer viewer;
312: PetscObjectGetComm((PetscObject)ksp,&comm);
313: viewer = PETSC_VIEWER_DRAW_(comm);
314: PetscViewerDrawGetDrawLG(viewer,0,&lg);
315: }
317: MPI_Comm_rank(((PetscObject)ksp)->comm,&rank);
318: if (!rank) {
319: if (!n) {PetscDrawLGReset(lg);}
320: x[0] = x[1] = (PetscReal) n;
321: if (rnorm > 0.0) y[0] = log10(rnorm); else y[0] = -15.0;
322: }
324: VecDuplicate(ksp->vec_rhs,&work);
325: KSPBuildResidual(ksp,0,work,&resid);
326: VecNorm(resid,NORM_2,&scnorm);
327: VecDestroy(work);
329: if (!rank) {
330: if (scnorm > 0.0) y[1] = log10(scnorm); else y[1] = -15.0;
331: PetscDrawLGAddPoint(lg,x,y);
332: if (n <= 20 || (n % 3)) {
333: PetscDrawLGDraw(lg);
334: }
335: }
336: return(0);
337: }
338:
341: /*@C
342: KSPMonitorLGTrueResidualNormDestroy - Destroys a line graph context that was created
343: with KSPMonitorLGTrueResidualNormCreate().
345: Collective on KSP
347: Input Parameter:
348: . draw - the drawing context
350: Level: intermediate
352: .keywords: KSP, monitor, line graph, destroy, true
354: .seealso: KSPMonitorLGTrueResidualNormCreate(), KSPMonitorSet()
355: @*/
356: PetscErrorCode KSPMonitorLGTrueResidualNormDestroy(PetscDrawLG drawlg)
357: {
359: PetscDraw draw;
362: PetscDrawLGGetDraw(drawlg,&draw);
363: PetscDrawDestroy(draw);
364: PetscDrawLGDestroy(drawlg);
365: return(0);
366: }