Actual source code: xmon.c
1: #define PETSCKSP_DLL
3: #include src/ksp/ksp/kspimpl.h
7: /*@C
8: KSPLGMonitorCreate - 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_xmonitor - Sets line graph monitor
26: Notes:
27: Use KSPLGMonitorDestroy() to destroy this line graph; do not use PetscDrawLGDestroy().
29: Level: intermediate
31: .keywords: KSP, monitor, line graph, residual, create
33: .seealso: KSPLGMonitorDestroy(), KSPSetMonitor(), KSPLGTrueMonitorCreate()
34: @*/
35: PetscErrorCode KSPLGMonitorCreate(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 KSPLGMonitor(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
51: {
52: PetscDrawLG lg = (PetscDrawLG) monctx;
54: PetscReal x,y;
57: if (!monctx) {
58: MPI_Comm comm;
59: PetscViewer viewer;
61: PetscObjectGetComm((PetscObject)ksp,&comm);
62: viewer = PETSC_VIEWER_DRAW_(comm);
63: PetscViewerDrawGetDrawLG(viewer,0,&lg);
64: }
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: KSPLGMonitorDestroy - Destroys a line graph context that was created
80: with KSPLGMonitorCreate().
82: Collective on KSP
84: Input Parameter:
85: . draw - the drawing context
87: Level: intermediate
89: .keywords: KSP, monitor, line graph, destroy
91: .seealso: KSPLGMonitorCreate(), KSPLGTrueMonitorDestroy(), KSPSetMonitor()
92: @*/
93: PetscErrorCode KSPLGMonitorDestroy(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: KSPLGTrueMonitorCreate - Creates a line graph context for use with
109: KSP to monitor convergence of true residual norms (as opposed to
110: preconditioned residual norms).
112: Collective on KSP
114: Input Parameters:
115: + host - the X display to open, or null for the local machine
116: . label - the title to put in the title bar
117: . x, y - the screen coordinates of the upper left coordinate of
118: the window
119: - m, n - the screen width and height in pixels
121: Output Parameter:
122: . draw - the drawing context
124: Options Database Key:
125: . -ksp_xtruemonitor - Sets true line graph monitor
127: Notes:
128: Use KSPLGTrueMonitorDestroy() to destroy this line graph, not
129: PetscDrawLGDestroy().
131: Level: intermediate
133: .keywords: KSP, monitor, line graph, residual, create, true
135: .seealso: KSPLGMonitorDestroy(), KSPSetMonitor(), KSPDefaultMonitor()
136: @*/
137: PetscErrorCode KSPLGTrueMonitorCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
138: {
139: PetscDraw win;
141: PetscMPIInt rank;
144: MPI_Comm_rank(comm,&rank);
145: if (rank) { *draw = 0; return(0);}
147: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
148: PetscDrawSetType(win,PETSC_DRAW_X);
149: PetscDrawLGCreate(win,2,draw);
150: PetscLogObjectParent(*draw,win);
151: return(0);
152: }
156: PetscErrorCode KSPLGTrueMonitor(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
157: {
158: PetscDrawLG lg = (PetscDrawLG) monctx;
159: PetscReal x[2],y[2],scnorm;
161: PetscMPIInt rank;
162: Vec resid,work;
165: if (!monctx) {
166: MPI_Comm comm;
167: PetscViewer viewer;
169: PetscObjectGetComm((PetscObject)ksp,&comm);
170: viewer = PETSC_VIEWER_DRAW_(comm);
171: PetscViewerDrawGetDrawLG(viewer,0,&lg);
172: }
174: MPI_Comm_rank(ksp->comm,&rank);
175: if (!rank) {
176: if (!n) {PetscDrawLGReset(lg);}
177: x[0] = x[1] = (PetscReal) n;
178: if (rnorm > 0.0) y[0] = log10(rnorm); else y[0] = -15.0;
179: }
181: VecDuplicate(ksp->vec_rhs,&work);
182: KSPBuildResidual(ksp,0,work,&resid);
183: VecNorm(resid,NORM_2,&scnorm);
184: VecDestroy(work);
186: if (!rank) {
187: if (scnorm > 0.0) y[1] = log10(scnorm); else y[1] = -15.0;
188: PetscDrawLGAddPoint(lg,x,y);
189: if (n <= 20 || (n % 3)) {
190: PetscDrawLGDraw(lg);
191: }
192: }
193: return(0);
194: }
195:
198: /*@C
199: KSPLGTrueMonitorDestroy - Destroys a line graph context that was created
200: with KSPLGTrueMonitorCreate().
202: Collective on KSP
204: Input Parameter:
205: . draw - the drawing context
207: Level: intermediate
209: .keywords: KSP, monitor, line graph, destroy, true
211: .seealso: KSPLGTrueMonitorCreate(), KSPSetMonitor()
212: @*/
213: PetscErrorCode KSPLGTrueMonitorDestroy(PetscDrawLG drawlg)
214: {
216: PetscDraw draw;
219: PetscDrawLGGetDraw(drawlg,&draw);
220: PetscDrawDestroy(draw);
221: PetscDrawLGDestroy(drawlg);
222: return(0);
223: }