Actual source code: xmon.c
petsc-3.5.4 2015-05-23
2: #include <petsc-private/kspimpl.h> /*I "petscksp.h" I*/
3: #include <petscdraw.h>
7: /*@C
8: KSPMonitorLGResidualNormCreate - 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_lg_residualnorm - Sets line graph monitor
26: Notes:
27: Use KSPMonitorLGResidualNormDestroy() to destroy this line graph; do not use PetscDrawLGDestroy().
29: Level: intermediate
31: .keywords: KSP, monitor, line graph, residual, create
33: .seealso: KSPMonitorLGResidualNormDestroy(), KSPMonitorSet(), KSPMonitorLGTrueResidualCreate()
34: @*/
35: PetscErrorCode KSPMonitorLGResidualNormCreate(const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
36: {
37: PetscDraw win;
39: PetscDrawAxis axis;
42: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
43: PetscDrawSetFromOptions(win);
44: PetscDrawLGCreate(win,1,draw);
45: PetscDrawLGGetAxis(*draw,&axis);
46: PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norm");
47: PetscLogObjectParent((PetscObject)*draw,(PetscObject)win);
48: return(0);
49: }
53: PetscErrorCode KSPMonitorLGResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
54: {
55: PetscDrawLG lg = (PetscDrawLG)monctx;
57: PetscReal x,y;
60: if (!n) {PetscDrawLGReset(lg);}
61: x = (PetscReal) n;
62: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
63: else y = -15.0;
64: PetscDrawLGAddPoint(lg,&x,&y);
65: if (n < 20 || !(n % 5)) {
66: PetscDrawLGDraw(lg);
67: }
68: return(0);
69: }
73: /*@
74: KSPMonitorLGResidualNormDestroy - Destroys a line graph context that was created
75: with KSPMonitorLGResidualNormCreate().
77: Collective on KSP
79: Input Parameter:
80: . draw - the drawing context
82: Level: intermediate
84: .keywords: KSP, monitor, line graph, destroy
86: .seealso: KSPMonitorLGResidualNormCreate(), KSPMonitorLGTrueResidualDestroy(), KSPMonitorSet()
87: @*/
88: PetscErrorCode KSPMonitorLGResidualNormDestroy(PetscDrawLG *drawlg)
89: {
90: PetscDraw draw;
94: PetscDrawLGGetDraw(*drawlg,&draw);
95: PetscDrawDestroy(&draw);
96: PetscDrawLGDestroy(drawlg);
97: return(0);
98: }
100: extern PetscErrorCode KSPMonitorRange_Private(KSP,PetscInt,PetscReal*);
103: PetscErrorCode KSPMonitorLGRange(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
104: {
105: PetscDrawLG lg;
106: PetscErrorCode ierr;
107: PetscReal x,y,per;
108: PetscViewer v = (PetscViewer)monctx;
109: static PetscReal prev; /* should be in the context */
110: PetscDraw draw;
113: PetscViewerDrawGetDrawLG(v,0,&lg);
114: if (!n) {PetscDrawLGReset(lg);}
115: PetscDrawLGGetDraw(lg,&draw);
116: PetscDrawSetTitle(draw,"Residual norm");
117: x = (PetscReal) n;
118: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
119: else y = -15.0;
120: PetscDrawLGAddPoint(lg,&x,&y);
121: if (n < 20 || !(n % 5)) {
122: PetscDrawLGDraw(lg);
123: }
125: PetscViewerDrawGetDrawLG(v,1,&lg);
126: KSPMonitorRange_Private(ksp,n,&per);
127: if (!n) {PetscDrawLGReset(lg);}
128: PetscDrawLGGetDraw(lg,&draw);
129: PetscDrawSetTitle(draw,"% elemts > .2*max elemt");
130: x = (PetscReal) n;
131: y = 100.0*per;
132: PetscDrawLGAddPoint(lg,&x,&y);
133: if (n < 20 || !(n % 5)) {
134: PetscDrawLGDraw(lg);
135: }
137: PetscViewerDrawGetDrawLG(v,2,&lg);
138: if (!n) {prev = rnorm;PetscDrawLGReset(lg);}
139: PetscDrawLGGetDraw(lg,&draw);
140: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
141: PetscDrawLGGetDraw(lg,&draw);
142: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm");
143: x = (PetscReal) n;
144: y = (prev - rnorm)/prev;
145: PetscDrawLGAddPoint(lg,&x,&y);
146: if (n < 20 || !(n % 5)) {
147: PetscDrawLGDraw(lg);
148: }
150: PetscViewerDrawGetDrawLG(v,3,&lg);
151: if (!n) {PetscDrawLGReset(lg);}
152: PetscDrawLGGetDraw(lg,&draw);
153: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
154: x = (PetscReal) n;
155: y = (prev - rnorm)/(prev*per);
156: if (n > 5) {
157: PetscDrawLGAddPoint(lg,&x,&y);
158: }
159: if (n < 20 || !(n % 5)) {
160: PetscDrawLGDraw(lg);
161: }
162: prev = rnorm;
163: return(0);
164: }
168: /*@C
169: KSPMonitorLGTrueResidualNormCreate - Creates a line graph context for use with
170: KSP to monitor convergence of true residual norms (as opposed to
171: preconditioned residual norms).
173: Collective on KSP
175: Input Parameters:
176: + host - the X display to open, or null for the local machine
177: . label - the title to put in the title bar
178: . x, y - the screen coordinates of the upper left coordinate of
179: the window
180: - m, n - the screen width and height in pixels
182: Output Parameter:
183: . draw - the drawing context
185: Options Database Key:
186: . -ksp_monitor_lg_true_residualnorm - Sets true line graph monitor
188: Notes:
189: Use KSPMonitorLGTrueResidualNormDestroy() to destroy this line graph, not
190: PetscDrawLGDestroy().
192: Level: intermediate
194: .keywords: KSP, monitor, line graph, residual, create, true
196: .seealso: KSPMonitorLGResidualNormDestroy(), KSPMonitorSet(), KSPMonitorDefault()
197: @*/
198: PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *draw)
199: {
200: PetscDraw win;
202: PetscMPIInt rank;
203: PetscDrawAxis axis;
206: MPI_Comm_rank(comm,&rank);
207: if (rank) { *draw = 0; return(0);}
209: PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);
210: PetscDrawSetFromOptions(win);
211: PetscDrawLGCreate(win,2,draw);
212: PetscDrawLGGetAxis(*draw,&axis);
213: PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norms");
214: PetscLogObjectParent((PetscObject)*draw,(PetscObject)win);
215: return(0);
216: }
220: PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
221: {
222: PetscDrawLG lg = (PetscDrawLG) monctx;
223: PetscReal x[2],y[2],scnorm;
225: PetscMPIInt rank;
226: Vec resid,work;
229: MPI_Comm_rank(PetscObjectComm((PetscObject)ksp),&rank);
230: if (!rank) {
231: if (!n) {PetscDrawLGReset(lg);}
232: x[0] = x[1] = (PetscReal) n;
233: if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
234: else y[0] = -15.0;
235: }
237: VecDuplicate(ksp->vec_rhs,&work);
238: KSPBuildResidual(ksp,0,work,&resid);
239: VecNorm(resid,NORM_2,&scnorm);
240: VecDestroy(&work);
242: if (!rank) {
243: if (scnorm > 0.0) y[1] = PetscLog10Real(scnorm);
244: else y[1] = -15.0;
245: PetscDrawLGAddPoint(lg,x,y);
246: if (n <= 20 || (n % 3)) {
247: PetscDrawLGDraw(lg);
248: }
249: }
250: return(0);
251: }
255: /*@C
256: KSPMonitorLGTrueResidualNormDestroy - Destroys a line graph context that was created
257: with KSPMonitorLGTrueResidualNormCreate().
259: Collective on KSP
261: Input Parameter:
262: . draw - the drawing context
264: Level: intermediate
266: .keywords: KSP, monitor, line graph, destroy, true
268: .seealso: KSPMonitorLGTrueResidualNormCreate(), KSPMonitorSet()
269: @*/
270: PetscErrorCode KSPMonitorLGTrueResidualNormDestroy(PetscDrawLG *drawlg)
271: {
273: PetscDraw draw;
276: PetscDrawLGGetDraw(*drawlg,&draw);
277: PetscDrawDestroy(&draw);
278: PetscDrawLGDestroy(drawlg);
279: return(0);
280: }