Actual source code: xmon.c
2: #include src/ksp/ksp/kspimpl.h
6: /*@C
7: KSPLGMonitorCreate - 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_xmonitor - Sets line graph monitor
25: Notes:
26: Use KSPLGMonitorDestroy() to destroy this line graph; do not use PetscDrawLGDestroy().
28: Level: intermediate
30: .keywords: KSP, monitor, line graph, residual, create
32: .seealso: KSPLGMonitorDestroy(), KSPSetMonitor(), KSPLGTrueMonitorCreate()
33: @*/
34: PetscErrorCode KSPLGMonitorCreate(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 KSPLGMonitor(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
50: {
51: PetscDrawLG lg = (PetscDrawLG) monctx;
53: PetscReal x,y;
56: if (!monctx) {
57: MPI_Comm comm;
58: PetscViewer viewer;
60: PetscObjectGetComm((PetscObject)ksp,&comm);
61: viewer = PETSC_VIEWER_DRAW_(comm);
62: PetscViewerDrawGetDrawLG(viewer,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: /*@C
78: KSPLGMonitorDestroy - Destroys a line graph context that was created
79: with KSPLGMonitorCreate().
81: Collective on KSP
83: Input Parameter:
84: . draw - the drawing context
86: Level: intermediate
88: .keywords: KSP, monitor, line graph, destroy
90: .seealso: KSPLGMonitorCreate(), KSPLGTrueMonitorDestroy(), KSPSetMonitor()
91: @*/
92: PetscErrorCode KSPLGMonitorDestroy(PetscDrawLG drawlg)
93: {
94: PetscDraw draw;
98: PetscDrawLGGetDraw(drawlg,&draw);
99: if (draw) { PetscDrawDestroy(draw);}
100: PetscDrawLGDestroy(drawlg);
101: return(0);
102: }
106: /*@C
107: KSPLGTrueMonitorCreate - Creates a line graph context for use with
108: KSP to monitor convergence of true residual norms (as opposed to
109: preconditioned residual norms).
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_xtruemonitor - Sets true line graph monitor
126: Notes:
127: Use KSPLGTrueMonitorDestroy() to destroy this line graph, not
128: PetscDrawLGDestroy().
130: Level: intermediate
132: .keywords: KSP, monitor, line graph, residual, create, true
134: .seealso: KSPLGMonitorDestroy(), KSPSetMonitor(), KSPDefaultMonitor()
135: @*/
136: PetscErrorCode KSPLGTrueMonitorCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
137: {
138: PetscDraw win;
140: PetscMPIInt rank;
143: MPI_Comm_rank(comm,&rank);
144: if (rank) { *draw = 0; return(0);}
146: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
147: PetscDrawSetType(win,PETSC_DRAW_X);
148: PetscDrawLGCreate(win,2,draw);
149: PetscLogObjectParent(*draw,win);
150: return(0);
151: }
155: PetscErrorCode KSPLGTrueMonitor(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
156: {
157: PetscDrawLG lg = (PetscDrawLG) monctx;
158: PetscReal x[2],y[2],scnorm;
160: PetscMPIInt rank;
161: Vec resid,work;
164: if (!monctx) {
165: MPI_Comm comm;
166: PetscViewer viewer;
168: PetscObjectGetComm((PetscObject)ksp,&comm);
169: viewer = PETSC_VIEWER_DRAW_(comm);
170: PetscViewerDrawGetDrawLG(viewer,0,&lg);
171: }
173: MPI_Comm_rank(ksp->comm,&rank);
174: if (!rank) {
175: if (!n) {PetscDrawLGReset(lg);}
176: x[0] = x[1] = (PetscReal) n;
177: if (rnorm > 0.0) y[0] = log10(rnorm); else y[0] = -15.0;
178: }
180: VecDuplicate(ksp->vec_rhs,&work);
181: KSPBuildResidual(ksp,0,work,&resid);
182: VecNorm(resid,NORM_2,&scnorm);
183: VecDestroy(work);
185: if (!rank) {
186: if (scnorm > 0.0) y[1] = log10(scnorm); else y[1] = -15.0;
187: PetscDrawLGAddPoint(lg,x,y);
188: if (n <= 20 || (n % 3)) {
189: PetscDrawLGDraw(lg);
190: }
191: }
192: return(0);
193: }
194:
197: /*@C
198: KSPLGTrueMonitorDestroy - Destroys a line graph context that was created
199: with KSPLGTrueMonitorCreate().
201: Collective on KSP
203: Input Parameter:
204: . draw - the drawing context
206: Level: intermediate
208: .keywords: KSP, monitor, line graph, destroy, true
210: .seealso: KSPLGTrueMonitorCreate(), KSPSetMonitor()
211: @*/
212: PetscErrorCode KSPLGTrueMonitorDestroy(PetscDrawLG drawlg)
213: {
215: PetscDraw draw;
218: PetscDrawLGGetDraw(drawlg,&draw);
219: PetscDrawDestroy(draw);
220: PetscDrawLGDestroy(drawlg);
221: return(0);
222: }