1
2 #include <petsc/private/kspimpl.h> /*I "petscksp.h" I*/
3 #include <petscdraw.h>
4
5 /*@C
6 KSPMonitorLGResidualNormCreate - Creates a line graph context for use with
7 KSP to monitor convergence of preconditioned residual norms.
8
9 Collective
10
11 Input Parameters:
12 + comm - communicator context
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
18
19 Output Parameter:
20 . lgctx - the drawing context
21
22 Options Database Key:
23 . -ksp_monitor_lg_residualnorm - Sets line graph monitor
24
25 Notes:
26 Use PetscDrawLGDestroy() to destroy this line graph.
27
28 Level: intermediate
29
30 .seealso: KSPMonitorSet(), KSPMonitorLGTrueResidualCreate()
31 @*/
KSPMonitorLGResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG * lgctx)32 PetscErrorCode KSPMonitorLGResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
33 {
34 PetscDraw draw;
35 PetscErrorCode ierr;
36 PetscDrawAxis axis;
37 PetscDrawLG lg;
38
39 PetscFunctionBegin;
40 ierr = PetscDrawCreate(comm,host,label,x,y,m,n,&draw);CHKERRQ(ierr);
41 ierr = PetscDrawSetFromOptions(draw);CHKERRQ(ierr);
42 ierr = PetscDrawLGCreate(draw,1,&lg);CHKERRQ(ierr);
43 ierr = PetscDrawLGSetFromOptions(lg);CHKERRQ(ierr);
44 ierr = PetscDrawLGGetAxis(lg,&axis);CHKERRQ(ierr);
45 ierr = PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norm");CHKERRQ(ierr);
46 ierr = PetscDrawDestroy(&draw);CHKERRQ(ierr);
47 *lgctx = lg;
48 PetscFunctionReturn(0);
49 }
50
KSPMonitorLGResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void * ctx)51 PetscErrorCode KSPMonitorLGResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *ctx)
52 {
53 PetscDrawLG lg = (PetscDrawLG) ctx;
54 PetscReal x,y;
55 PetscErrorCode ierr;
56
57 PetscFunctionBegin;
58 PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,4);
59 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);}
60 x = (PetscReal) n;
61 if (rnorm > 0.0) y = PetscLog10Real(rnorm);
62 else y = -15.0;
63 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr);
64 if (n <= 20 || !(n % 5) || ksp->reason) {
65 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr);
66 ierr = PetscDrawLGSave(lg);CHKERRQ(ierr);
67 }
68 PetscFunctionReturn(0);
69 }
70
71 extern PetscErrorCode KSPMonitorRange_Private(KSP,PetscInt,PetscReal*);
KSPMonitorLGRange(KSP ksp,PetscInt n,PetscReal rnorm,void * monctx)72 PetscErrorCode KSPMonitorLGRange(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
73 {
74 PetscDrawLG lg;
75 PetscErrorCode ierr;
76 PetscReal x,y,per;
77 PetscViewer v = (PetscViewer)monctx;
78 static PetscReal prev; /* should be in the context */
79 PetscDraw draw;
80
81 PetscFunctionBegin;
82 PetscValidHeaderSpecific(v,PETSC_VIEWER_CLASSID,4);
83
84 ierr = KSPMonitorRange_Private(ksp,n,&per);CHKERRQ(ierr);
85 if (!n) prev = rnorm;
86
87 ierr = PetscViewerDrawGetDrawLG(v,0,&lg);CHKERRQ(ierr);
88 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);}
89 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr);
90 ierr = PetscDrawSetTitle(draw,"Residual norm");CHKERRQ(ierr);
91 x = (PetscReal) n;
92 if (rnorm > 0.0) y = PetscLog10Real(rnorm);
93 else y = -15.0;
94 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr);
95 if (n < 20 || !(n % 5) || ksp->reason) {
96 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr);
97 ierr = PetscDrawLGSave(lg);CHKERRQ(ierr);
98 }
99
100 ierr = PetscViewerDrawGetDrawLG(v,1,&lg);CHKERRQ(ierr);
101 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);}
102 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr);
103 ierr = PetscDrawSetTitle(draw,"% elemts > .2*max elemt");CHKERRQ(ierr);
104 x = (PetscReal) n;
105 y = 100.0*per;
106 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr);
107 if (n < 20 || !(n % 5) || ksp->reason) {
108 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr);
109 ierr = PetscDrawLGSave(lg);CHKERRQ(ierr);
110 }
111
112 ierr = PetscViewerDrawGetDrawLG(v,2,&lg);CHKERRQ(ierr);
113 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);}
114 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr);
115 ierr = PetscDrawSetTitle(draw,"(norm-oldnorm)/oldnorm");CHKERRQ(ierr);
116 x = (PetscReal) n;
117 y = (prev - rnorm)/prev;
118 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr);
119 if (n < 20 || !(n % 5) || ksp->reason) {
120 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr);
121 ierr = PetscDrawLGSave(lg);CHKERRQ(ierr);
122 }
123
124 ierr = PetscViewerDrawGetDrawLG(v,3,&lg);CHKERRQ(ierr);
125 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);}
126 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr);
127 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");CHKERRQ(ierr);
128 x = (PetscReal) n;
129 y = (prev - rnorm)/(prev*per);
130 if (n > 5) {
131 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr);
132 }
133 if (n < 20 || !(n % 5) || ksp->reason) {
134 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr);
135 ierr = PetscDrawLGSave(lg);CHKERRQ(ierr);
136 }
137
138 prev = rnorm;
139 PetscFunctionReturn(0);
140 }
141
142 /*@C
143 KSPMonitorLGTrueResidualNormCreate - Creates a line graph context for use with
144 KSP to monitor convergence of true residual norms (as opposed to
145 preconditioned residual norms).
146
147 Collective
148
149 Input Parameters:
150 + comm - communicator context
151 . host - the X display to open, or null for the local machine
152 . label - the title to put in the title bar
153 . x, y - the screen coordinates of the upper left coordinate of
154 the window
155 - m, n - the screen width and height in pixels
156
157 Output Parameter:
158 . lgctx - the drawing context
159
160 Options Database Key:
161 . -ksp_monitor_lg_true_residualnorm - Sets true line graph monitor
162
163 Notes:
164 Use PetscDrawLGDestroy() to destroy this line graph.
165
166 Level: intermediate
167
168 .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate()
169 @*/
KSPMonitorLGTrueResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG * lgctx)170 PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
171 {
172 PetscDraw draw;
173 PetscErrorCode ierr;
174 PetscDrawAxis axis;
175 PetscDrawLG lg;
176 const char *names[] = {"Preconditioned","True"};
177
178 PetscFunctionBegin;
179 ierr = PetscDrawCreate(comm,host,label,x,y,m,n,&draw);CHKERRQ(ierr);
180 ierr = PetscDrawSetFromOptions(draw);CHKERRQ(ierr);
181 ierr = PetscDrawLGCreate(draw,2,&lg);CHKERRQ(ierr);
182 ierr = PetscDrawLGSetLegend(lg,names);CHKERRQ(ierr);
183 ierr = PetscDrawLGSetFromOptions(lg);CHKERRQ(ierr);
184 ierr = PetscDrawLGGetAxis(lg,&axis);CHKERRQ(ierr);
185 ierr = PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norm");CHKERRQ(ierr);
186 ierr = PetscDrawDestroy(&draw);CHKERRQ(ierr);
187 *lgctx = lg;
188 PetscFunctionReturn(0);
189 }
190
KSPMonitorLGTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void * ctx)191 PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *ctx)
192 {
193 PetscDrawLG lg = (PetscDrawLG) ctx;
194 PetscReal x[2],y[2],scnorm;
195 Vec resid,work;
196 PetscErrorCode ierr;
197
198 PetscFunctionBegin;
199 PetscValidHeaderSpecific(ksp,KSP_CLASSID,1);
200 PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,4);
201 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);}
202 x[0] = x[1] = (PetscReal) n;
203 if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
204 else y[0] = -15.0;
205 ierr = VecDuplicate(ksp->vec_rhs,&work);CHKERRQ(ierr);
206 ierr = KSPBuildResidual(ksp,NULL,work,&resid);CHKERRQ(ierr);
207 ierr = VecNorm(resid,NORM_2,&scnorm);CHKERRQ(ierr);
208 ierr = VecDestroy(&work);CHKERRQ(ierr);
209 if (scnorm > 0.0) y[1] = PetscLog10Real(scnorm);
210 else y[1] = -15.0;
211 ierr = PetscDrawLGAddPoint(lg,x,y);CHKERRQ(ierr);
212 if (n <= 20 || !(n % 5)) {
213 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr);
214 ierr = PetscDrawLGSave(lg);CHKERRQ(ierr);
215 }
216 PetscFunctionReturn(0);
217 }
218