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