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