1 2 #include <petsc/private/kspimpl.h> /*I "petscksp.h" I*/ 3 #include <petscdraw.h> 4 5 #undef __FUNCT__ 6 #define __FUNCT__ "KSPMonitorLGResidualNormCreate" 7 /*@C 8 KSPMonitorLGResidualNormCreate - Creates a line graph context for use with 9 KSP to monitor convergence of preconditioned residual norms. 10 11 Collective on KSP 12 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 19 20 Output Parameter: 21 . draw - the drawing context 22 23 Options Database Key: 24 . -ksp_monitor_lg_residualnorm - Sets line graph monitor 25 26 Notes: 27 Use KSPMonitorLGResidualNormDestroy() to destroy this line graph; do not use PetscDrawLGDestroy(). 28 29 Level: intermediate 30 31 .keywords: KSP, monitor, line graph, residual, create 32 33 .seealso: KSPMonitorLGResidualNormDestroy(), KSPMonitorSet(), KSPMonitorLGTrueResidualCreate() 34 @*/ 35 PetscErrorCode KSPMonitorLGResidualNormCreate(const char host[],const char label[],int x,int y,int m,int n,PetscObject **objs) 36 { 37 PetscDraw win; 38 PetscErrorCode ierr; 39 PetscDrawAxis axis; 40 PetscDrawLG draw; 41 42 PetscFunctionBegin; 43 ierr = PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);CHKERRQ(ierr); 44 ierr = PetscDrawSetFromOptions(win);CHKERRQ(ierr); 45 ierr = PetscDrawLGCreate(win,1,&draw);CHKERRQ(ierr); 46 ierr = PetscDrawLGSetFromOptions(draw);CHKERRQ(ierr); 47 ierr = PetscDrawLGGetAxis(draw,&axis);CHKERRQ(ierr); 48 ierr = PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norm");CHKERRQ(ierr); 49 ierr = PetscLogObjectParent((PetscObject)draw,(PetscObject)win);CHKERRQ(ierr); 50 ierr = PetscMalloc1(2,objs);CHKERRQ(ierr); 51 (*objs)[0] = (PetscObject)draw; 52 (*objs)[1] = (PetscObject)win; 53 PetscFunctionReturn(0); 54 } 55 56 #undef __FUNCT__ 57 #define __FUNCT__ "KSPMonitorLGResidualNorm" 58 PetscErrorCode KSPMonitorLGResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,PetscObject *objs) 59 { 60 PetscErrorCode ierr; 61 PetscReal x,y; 62 PetscDrawLG lg = (PetscDrawLG) objs[0]; 63 64 PetscFunctionBegin; 65 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 66 x = (PetscReal) n; 67 if (rnorm > 0.0) y = PetscLog10Real(rnorm); 68 else y = -15.0; 69 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 70 if (n < 20 || !(n % 5)) { 71 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 72 } 73 PetscFunctionReturn(0); 74 } 75 76 #undef __FUNCT__ 77 #define __FUNCT__ "KSPMonitorLGResidualNormDestroy" 78 /*@ 79 KSPMonitorLGResidualNormDestroy - Destroys a line graph context that was created 80 with KSPMonitorLGResidualNormCreate(). 81 82 Collective on KSP 83 84 Input Parameter: 85 . draw - the drawing context 86 87 Level: intermediate 88 89 .keywords: KSP, monitor, line graph, destroy 90 91 .seealso: KSPMonitorLGResidualNormCreate(), KSPMonitorLGTrueResidualDestroy(), KSPMonitorSet() 92 @*/ 93 PetscErrorCode KSPMonitorLGResidualNormDestroy(PetscObject **objs) 94 { 95 PetscErrorCode ierr; 96 PetscDrawLG drawlg = (PetscDrawLG) (*objs)[0]; 97 PetscDraw draw = (PetscDraw) (*objs)[1]; 98 99 PetscFunctionBegin; 100 ierr = PetscDrawDestroy(&draw);CHKERRQ(ierr); 101 ierr = PetscDrawLGDestroy(&drawlg);CHKERRQ(ierr); 102 ierr = PetscFree(*objs);CHKERRQ(ierr); 103 PetscFunctionReturn(0); 104 } 105 106 extern PetscErrorCode KSPMonitorRange_Private(KSP,PetscInt,PetscReal*); 107 #undef __FUNCT__ 108 #define __FUNCT__ "KSPMonitorLGRange" 109 PetscErrorCode KSPMonitorLGRange(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx) 110 { 111 PetscDrawLG lg; 112 PetscErrorCode ierr; 113 PetscReal x,y,per; 114 PetscViewer v = (PetscViewer)monctx; 115 static PetscReal prev; /* should be in the context */ 116 PetscDraw draw; 117 118 PetscFunctionBegin; 119 PetscValidHeaderSpecific(v,PETSC_VIEWER_CLASSID,4); 120 ierr = PetscViewerDrawGetDrawLG(v,0,&lg);CHKERRQ(ierr); 121 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 122 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 123 ierr = PetscDrawSetTitle(draw,"Residual norm");CHKERRQ(ierr); 124 x = (PetscReal) n; 125 if (rnorm > 0.0) y = PetscLog10Real(rnorm); 126 else y = -15.0; 127 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 128 if (n < 20 || !(n % 5)) { 129 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 130 } 131 132 ierr = PetscViewerDrawGetDrawLG(v,1,&lg);CHKERRQ(ierr); 133 ierr = KSPMonitorRange_Private(ksp,n,&per);CHKERRQ(ierr); 134 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 135 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 136 ierr = PetscDrawSetTitle(draw,"% elemts > .2*max elemt");CHKERRQ(ierr); 137 x = (PetscReal) n; 138 y = 100.0*per; 139 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 140 if (n < 20 || !(n % 5)) { 141 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 142 } 143 144 ierr = PetscViewerDrawGetDrawLG(v,2,&lg);CHKERRQ(ierr); 145 if (!n) {prev = rnorm;ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 146 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 147 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");CHKERRQ(ierr); 148 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 149 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm");CHKERRQ(ierr); 150 x = (PetscReal) n; 151 y = (prev - rnorm)/prev; 152 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 153 if (n < 20 || !(n % 5)) { 154 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 155 } 156 157 ierr = PetscViewerDrawGetDrawLG(v,3,&lg);CHKERRQ(ierr); 158 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 159 ierr = PetscDrawLGGetDraw(lg,&draw);CHKERRQ(ierr); 160 ierr = PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");CHKERRQ(ierr); 161 x = (PetscReal) n; 162 y = (prev - rnorm)/(prev*per); 163 if (n > 5) { 164 ierr = PetscDrawLGAddPoint(lg,&x,&y);CHKERRQ(ierr); 165 } 166 if (n < 20 || !(n % 5)) { 167 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 168 } 169 prev = rnorm; 170 PetscFunctionReturn(0); 171 } 172 173 #undef __FUNCT__ 174 #define __FUNCT__ "KSPMonitorLGTrueResidualNormCreate" 175 /*@C 176 KSPMonitorLGTrueResidualNormCreate - Creates a line graph context for use with 177 KSP to monitor convergence of true residual norms (as opposed to 178 preconditioned residual norms). 179 180 Collective on KSP 181 182 Input Parameters: 183 + host - the X display to open, or null for the local machine 184 . label - the title to put in the title bar 185 . x, y - the screen coordinates of the upper left coordinate of 186 the window 187 - m, n - the screen width and height in pixels 188 189 Output Parameter: 190 . draw - the drawing context 191 192 Options Database Key: 193 . -ksp_monitor_lg_true_residualnorm - Sets true line graph monitor 194 195 Notes: 196 Use KSPMonitorLGTrueResidualNormDestroy() to destroy this line graph, not 197 PetscDrawLGDestroy(). 198 199 Level: intermediate 200 201 .keywords: KSP, monitor, line graph, residual, create, true 202 203 .seealso: KSPMonitorLGResidualNormDestroy(), KSPMonitorSet(), KSPMonitorDefault() 204 @*/ 205 PetscErrorCode KSPMonitorLGTrueResidualNormCreate(const char host[],const char label[],int x,int y,int m,int n,PetscObject **objs) 206 { 207 PetscDraw win; 208 PetscErrorCode ierr; 209 PetscDrawAxis axis; 210 PetscDrawLG draw; 211 212 PetscFunctionBegin; 213 ierr = PetscDrawCreate(PETSC_COMM_SELF,host,label,x,y,m,n,&win);CHKERRQ(ierr); 214 ierr = PetscDrawSetFromOptions(win);CHKERRQ(ierr); 215 ierr = PetscDrawLGCreate(win,2,&draw);CHKERRQ(ierr); 216 ierr = PetscDrawLGSetFromOptions(draw);CHKERRQ(ierr); 217 ierr = PetscDrawLGGetAxis(draw,&axis);CHKERRQ(ierr); 218 ierr = PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norms");CHKERRQ(ierr); 219 ierr = PetscLogObjectParent((PetscObject)draw,(PetscObject)win);CHKERRQ(ierr); 220 ierr = PetscMalloc1(2,objs);CHKERRQ(ierr); 221 (*objs)[0] = (PetscObject)draw; 222 (*objs)[1] = (PetscObject)win; 223 PetscFunctionReturn(0); 224 } 225 226 #undef __FUNCT__ 227 #define __FUNCT__ "KSPMonitorLGTrueResidualNorm" 228 PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,PetscObject *objs) 229 { 230 PetscDrawLG lg = (PetscDrawLG) objs[0]; 231 PetscReal x[2],y[2],scnorm; 232 PetscErrorCode ierr; 233 PetscMPIInt rank; 234 Vec resid,work; 235 236 PetscFunctionBegin; 237 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)ksp),&rank);CHKERRQ(ierr); 238 if (!rank) { 239 if (!n) {ierr = PetscDrawLGReset(lg);CHKERRQ(ierr);} 240 x[0] = x[1] = (PetscReal) n; 241 if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm); 242 else y[0] = -15.0; 243 } 244 245 ierr = VecDuplicate(ksp->vec_rhs,&work);CHKERRQ(ierr); 246 ierr = KSPBuildResidual(ksp,0,work,&resid);CHKERRQ(ierr); 247 ierr = VecNorm(resid,NORM_2,&scnorm);CHKERRQ(ierr); 248 ierr = VecDestroy(&work);CHKERRQ(ierr); 249 250 if (!rank) { 251 if (scnorm > 0.0) y[1] = PetscLog10Real(scnorm); 252 else y[1] = -15.0; 253 ierr = PetscDrawLGAddPoint(lg,x,y);CHKERRQ(ierr); 254 if (n <= 20 || (n % 3)) { 255 ierr = PetscDrawLGDraw(lg);CHKERRQ(ierr); 256 } 257 } 258 PetscFunctionReturn(0); 259 } 260 261 #undef __FUNCT__ 262 #define __FUNCT__ "KSPMonitorLGTrueResidualNormDestroy" 263 /*@C 264 KSPMonitorLGTrueResidualNormDestroy - Destroys a line graph context that was created 265 with KSPMonitorLGTrueResidualNormCreate(). 266 267 Collective on KSP 268 269 Input Parameter: 270 . draw - the drawing context 271 272 Level: intermediate 273 274 .keywords: KSP, monitor, line graph, destroy, true 275 276 .seealso: KSPMonitorLGTrueResidualNormCreate(), KSPMonitorSet() 277 @*/ 278 PetscErrorCode KSPMonitorLGTrueResidualNormDestroy(PetscObject **objs) 279 { 280 PetscErrorCode ierr; 281 PetscDrawLG drawlg = (PetscDrawLG) (*objs)[0]; 282 PetscDraw draw = (PetscDraw) (*objs)[1]; 283 284 PetscFunctionBegin; 285 ierr = PetscDrawDestroy(&draw);CHKERRQ(ierr); 286 ierr = PetscDrawLGDestroy(&drawlg);CHKERRQ(ierr); 287 ierr = PetscFree(*objs);CHKERRQ(ierr); 288 PetscFunctionReturn(0); 289 } 290 291 292