1 2 #include <../src/sys/classes/draw/utils/axisimpl.h> 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "PetscRint" 6 static PetscErrorCode PetscRint(PetscReal x,PetscReal *result) 7 { 8 PetscFunctionBegin; 9 if (x > 0) *result = floor(x + 0.5); 10 else *result = floor(x - 0.5); 11 PetscFunctionReturn(0); 12 } 13 14 #undef __FUNCT__ 15 #define __FUNCT__ "PetscDrawAxisSetLimits" 16 /*@ 17 PetscDrawAxisSetLimits - Sets the limits (in user coords) of the axis 18 19 Not Collective (ignored on all processors except processor 0 of PetscDrawAxis) 20 21 Input Parameters: 22 + axis - the axis 23 . xmin,xmax - limits in x 24 - ymin,ymax - limits in y 25 26 Options Database: 27 . -drawaxis_hold - hold the initial set of axis limits for future plotting 28 29 Level: advanced 30 31 .seealso: PetscDrawAxisSetHoldLimits() 32 33 @*/ 34 PetscErrorCode PetscDrawAxisSetLimits(PetscDrawAxis axis,PetscReal xmin,PetscReal xmax,PetscReal ymin,PetscReal ymax) 35 { 36 PetscErrorCode ierr; 37 38 PetscFunctionBegin; 39 if (!axis) PetscFunctionReturn(0); 40 if (axis->hold) PetscFunctionReturn(0); 41 axis->xlow = xmin; 42 axis->xhigh= xmax; 43 axis->ylow = ymin; 44 axis->yhigh= ymax; 45 ierr = PetscOptionsHasName(((PetscObject)axis)->prefix,"-drawaxis_hold",&axis->hold);CHKERRQ(ierr); 46 PetscFunctionReturn(0); 47 } 48 49 #undef __FUNCT__ 50 #define __FUNCT__ "PetscDrawAxisGetLimits" 51 /*@ 52 PetscDrawAxisGetLimits - Gets the limits (in user coords) of the axis 53 54 Not Collective (ignored on all processors except processor 0 of PetscDrawAxis) 55 56 Input Parameters: 57 + axis - the axis 58 . xmin,xmax - limits in x 59 - ymin,ymax - limits in y 60 61 Level: advanced 62 63 .seealso: PetscDrawAxisSetLimits() 64 65 @*/ 66 PetscErrorCode PetscDrawAxisGetLimits(PetscDrawAxis axis,PetscReal *xmin,PetscReal *xmax,PetscReal *ymin,PetscReal *ymax) 67 { 68 PetscFunctionBegin; 69 if (!axis) PetscFunctionReturn(0); 70 if (axis->hold) PetscFunctionReturn(0); 71 *xmin = axis->xlow; 72 *xmax = axis->xhigh; 73 *ymin = axis->ylow; 74 *ymax = axis->yhigh; 75 PetscFunctionReturn(0); 76 } 77 78 #undef __FUNCT__ 79 #define __FUNCT__ "PetscADefLabel" 80 /* 81 val is the label value. sep is the separation to the next (or previous) 82 label; this is useful in determining how many significant figures to 83 keep. 84 */ 85 PetscErrorCode PetscADefLabel(PetscReal val,PetscReal sep,char **p) 86 { 87 static char buf[40]; 88 char fmat[10]; 89 PetscErrorCode ierr; 90 int w,d; 91 PetscReal rval; 92 93 PetscFunctionBegin; 94 /* Find the string */ 95 if (PetscAbsReal(val)/sep < 1.e-4) { 96 buf[0] = '0'; buf[1] = 0; 97 } else if (PetscAbsReal(val) < 1.0e6 && PetscAbsReal(val) > 1.e-4) { 98 /* Compute the number of digits */ 99 w = 0; 100 d = 0; 101 if (sep > 0.0) { 102 d = (int)ceil(- log10 (sep)); 103 if (d < 0) d = 0; 104 if (PetscAbsReal(val) < 1.0e-6*sep) { 105 /* This is the case where we are near zero and less than a small 106 fraction of the sep. In this case, we use 0 as the value */ 107 val = 0.0; 108 w = d; 109 } 110 else if (val == 0.0) w = d; 111 else w = (int)(ceil(log10(PetscAbsReal(val))) + d); 112 if (w < 1) w ++; 113 if (val < 0) w ++; 114 } 115 116 ierr = PetscRint(val,&rval);CHKERRQ(ierr); 117 if (rval == val) { 118 if (w > 0) sprintf(fmat,"%%%dd",w); 119 else {ierr = PetscStrcpy(fmat,"%d");CHKERRQ(ierr);} 120 sprintf(buf,fmat,(int)val); 121 ierr = PetscStripInitialZero(buf);CHKERRQ(ierr); 122 ierr = PetscStripAllZeros(buf);CHKERRQ(ierr); 123 ierr = PetscStripTrailingZeros(buf);CHKERRQ(ierr); 124 } else { 125 /* The code used here is inappropriate for a val of 0, which 126 tends to print with an excessive numer of digits. In this 127 case, we should look at the next/previous values and 128 use those widths */ 129 if (w > 0) sprintf(fmat,"%%%d.%dlf",w + 1,d); 130 else {ierr = PetscStrcpy(fmat,"%lf");CHKERRQ(ierr);} 131 sprintf(buf,fmat,(double)val); 132 ierr = PetscStripInitialZero(buf);CHKERRQ(ierr); 133 ierr = PetscStripAllZeros(buf);CHKERRQ(ierr); 134 ierr = PetscStripTrailingZeros(buf);CHKERRQ(ierr); 135 } 136 } else { 137 ierr = PetscSNPrintf(buf,40,"%g",(double)val); 138 /* remove the extraneous 0 before the e */ 139 ierr = PetscStripZeros(buf);CHKERRQ(ierr); 140 ierr = PetscStripZerosPlus(buf);CHKERRQ(ierr); 141 ierr = PetscStripInitialZero(buf);CHKERRQ(ierr); 142 ierr = PetscStripAllZeros(buf);CHKERRQ(ierr); 143 ierr = PetscStripTrailingZeros(buf);CHKERRQ(ierr); 144 } 145 *p =buf; 146 PetscFunctionReturn(0); 147 } 148 149 #undef __FUNCT__ 150 #define __FUNCT__ "PetscADefTicks" 151 /* Finds "nice" locations for the ticks */ 152 PetscErrorCode PetscADefTicks(PetscReal low,PetscReal high,int num,int *ntick,PetscReal * tickloc,int maxtick) 153 { 154 PetscErrorCode ierr; 155 int i,power; 156 PetscReal x = 0.0,base=0.0; 157 158 PetscFunctionBegin; 159 /* patch if low == high */ 160 if (low == high) { 161 low -= .01; 162 high += .01; 163 } 164 165 /* if (PetscAbsReal(low-high) < 1.e-8) { 166 low -= .01; 167 high += .01; 168 } */ 169 170 ierr = PetscAGetBase(low,high,num,&base,&power);CHKERRQ(ierr); 171 ierr = PetscAGetNice(low,base,-1,&x);CHKERRQ(ierr); 172 173 /* Values are of the form j * base */ 174 /* Find the starting value */ 175 if (x < low) x += base; 176 177 i = 0; 178 while (i < maxtick && x <= high) { 179 tickloc[i++] = x; 180 x += base; 181 } 182 *ntick = i; 183 184 if (i < 2 && num < 10) { 185 ierr = PetscADefTicks(low,high,num+1,ntick,tickloc,maxtick);CHKERRQ(ierr); 186 } 187 if (num == 10) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"trouble"); 188 PetscFunctionReturn(0); 189 } 190 191 #define EPS 1.e-6 192 193 #undef __FUNCT__ 194 #define __FUNCT__ "PetscExp10" 195 PetscErrorCode PetscExp10(PetscReal d,PetscReal *result) 196 { 197 PetscFunctionBegin; 198 *result = pow((PetscReal)10.0,d); 199 PetscFunctionReturn(0); 200 } 201 202 #undef __FUNCT__ 203 #define __FUNCT__ "PetscMod" 204 PetscErrorCode PetscMod(PetscReal x,PetscReal y,PetscReal *result) 205 { 206 int i; 207 208 PetscFunctionBegin; 209 if (y == 1) { 210 *result = 0.0; 211 PetscFunctionReturn(0); 212 } 213 i = ((int)x) / ((int)y); 214 x = x - i * y; 215 while (x > y) x -= y; 216 *result = x; 217 PetscFunctionReturn(0); 218 } 219 220 #undef __FUNCT__ 221 #define __FUNCT__ "PetscCopysign" 222 PetscErrorCode PetscCopysign(PetscReal a,PetscReal b,PetscReal *result) 223 { 224 PetscFunctionBegin; 225 if (b >= 0) *result = a; 226 else *result = -a; 227 PetscFunctionReturn(0); 228 } 229 230 #undef __FUNCT__ 231 #define __FUNCT__ "PetscAGetNice" 232 /* 233 Given a value "in" and a "base", return a nice value. 234 based on "sign", extend up (+1) or down (-1) 235 */ 236 PetscErrorCode PetscAGetNice(PetscReal in,PetscReal base,int sign,PetscReal *result) 237 { 238 PetscReal etmp,s,s2,m; 239 PetscErrorCode ierr; 240 241 PetscFunctionBegin; 242 ierr = PetscCopysign (0.5,(double)sign,&s);CHKERRQ(ierr); 243 etmp = in / base + 0.5 + s; 244 ierr = PetscCopysign (0.5,etmp,&s);CHKERRQ(ierr); 245 ierr = PetscCopysign (EPS * etmp,(double)sign,&s2);CHKERRQ(ierr); 246 etmp = etmp - 0.5 + s - s2; 247 ierr = PetscMod(etmp,1.0,&m);CHKERRQ(ierr); 248 etmp = base * (etmp - m); 249 *result = etmp; 250 PetscFunctionReturn(0); 251 } 252 253 #undef __FUNCT__ 254 #define __FUNCT__ "PetscAGetBase" 255 PetscErrorCode PetscAGetBase(PetscReal vmin,PetscReal vmax,int num,PetscReal*Base,int*power) 256 { 257 PetscReal base,ftemp,e10; 258 static PetscReal base_try[5] = {10.0,5.0,2.0,1.0,0.5}; 259 PetscErrorCode ierr; 260 int i; 261 262 PetscFunctionBegin; 263 /* labels of the form n * BASE */ 264 /* get an approximate value for BASE */ 265 base = (vmax - vmin) / (double)(num + 1); 266 267 /* make it of form m x 10^power, m in [1.0, 10) */ 268 if (base <= 0.0) { 269 base = PetscAbsReal(vmin); 270 if (base < 1.0) base = 1.0; 271 } 272 ftemp = log10((1.0 + EPS) * base); 273 if (ftemp < 0.0) ftemp -= 1.0; 274 *power = (int)ftemp; 275 ierr = PetscExp10((double)- *power,&e10);CHKERRQ(ierr); 276 base = base * e10; 277 if (base < 1.0) base = 1.0; 278 /* now reduce it to one of 1, 2, or 5 */ 279 for (i=1; i<5; i++) { 280 if (base >= base_try[i]) { 281 ierr = PetscExp10((double)*power,&e10);CHKERRQ(ierr); 282 base = base_try[i-1] * e10; 283 if (i == 1) *power = *power + 1; 284 break; 285 } 286 } 287 *Base = base; 288 PetscFunctionReturn(0); 289 } 290 291 292 293 294 295 296