1 2 #include <petsc/private/drawimpl.h> /*I "petscdraw.h" I*/ 3 4 /* 5 val is the label value. sep is the separation to the next (or previous) 6 label; this is useful in determining how many significant figures to 7 keep. 8 */ 9 PetscErrorCode PetscADefLabel(PetscReal val,PetscReal sep,char **p) 10 { 11 static char buf[40]; 12 PetscErrorCode ierr; 13 14 PetscFunctionBegin; 15 /* Find the string */ 16 if (PetscAbsReal(val)/sep < 1.e-4) { 17 buf[0] = '0'; buf[1] = 0; 18 } else { 19 sprintf(buf,"%0.1e",(double)val); 20 ierr = PetscStripZerosPlus(buf);CHKERRQ(ierr); 21 ierr = PetscStripe0(buf);CHKERRQ(ierr); 22 ierr = PetscStripInitialZero(buf);CHKERRQ(ierr); 23 ierr = PetscStripAllZeros(buf);CHKERRQ(ierr); 24 ierr = PetscStripTrailingZeros(buf);CHKERRQ(ierr); 25 } 26 *p = buf; 27 PetscFunctionReturn(0); 28 } 29 30 /* Finds "nice" locations for the ticks */ 31 PetscErrorCode PetscADefTicks(PetscReal low,PetscReal high,int num,int *ntick,PetscReal *tickloc,int maxtick) 32 { 33 PetscErrorCode ierr; 34 int i,power; 35 PetscReal x = 0.0,base=0.0,eps; 36 37 PetscFunctionBegin; 38 ierr = PetscAGetBase(low,high,num,&base,&power);CHKERRQ(ierr); 39 ierr = PetscAGetNice(low,base,-1,&x);CHKERRQ(ierr); 40 41 /* Values are of the form j * base */ 42 /* Find the starting value */ 43 if (x < low) x += base; 44 45 i = 0; eps = base/10; 46 while (i < maxtick && x <= high+eps) { 47 tickloc[i++] = x; x += base; 48 } 49 *ntick = i; 50 tickloc[i-1] = PetscMin(tickloc[i-1],high); 51 52 if (i < 2 && num < 10) { 53 ierr = PetscADefTicks(low,high,num+1,ntick,tickloc,maxtick);CHKERRQ(ierr); 54 } 55 PetscFunctionReturn(0); 56 } 57 58 #define EPS 1.e-6 59 60 PetscErrorCode PetscExp10(PetscReal d,PetscReal *result) 61 { 62 PetscFunctionBegin; 63 *result = PetscPowReal((PetscReal)10.0,d); 64 PetscFunctionReturn(0); 65 } 66 67 PetscErrorCode PetscMod(PetscReal x,PetscReal y,PetscReal *result) 68 { 69 int i; 70 71 PetscFunctionBegin; 72 if (y == 1) { 73 *result = 0.0; 74 PetscFunctionReturn(0); 75 } 76 i = ((int)x) / ((int)y); 77 x = x - i * y; 78 while (x > y) x -= y; 79 *result = x; 80 PetscFunctionReturn(0); 81 } 82 83 PetscErrorCode PetscCopysign(PetscReal a,PetscReal b,PetscReal *result) 84 { 85 PetscFunctionBegin; 86 if (b >= 0) *result = a; 87 else *result = -a; 88 PetscFunctionReturn(0); 89 } 90 91 /* 92 Given a value "in" and a "base", return a nice value. 93 based on "sign", extend up (+1) or down (-1) 94 */ 95 PetscErrorCode PetscAGetNice(PetscReal in,PetscReal base,int sign,PetscReal *result) 96 { 97 PetscReal etmp,s,s2,m; 98 PetscErrorCode ierr; 99 100 PetscFunctionBegin; 101 ierr = PetscCopysign (0.5,(double)sign,&s);CHKERRQ(ierr); 102 etmp = in / base + 0.5 + s; 103 ierr = PetscCopysign (0.5,etmp,&s);CHKERRQ(ierr); 104 ierr = PetscCopysign (EPS * etmp,(double)sign,&s2);CHKERRQ(ierr); 105 etmp = etmp - 0.5 + s - s2; 106 ierr = PetscMod(etmp,1.0,&m);CHKERRQ(ierr); 107 etmp = base * (etmp - m); 108 *result = etmp; 109 PetscFunctionReturn(0); 110 } 111 112 PetscErrorCode PetscAGetBase(PetscReal vmin,PetscReal vmax,int num,PetscReal *Base,int *power) 113 { 114 PetscReal base,ftemp,e10; 115 static PetscReal base_try[5] = {10.0,5.0,2.0,1.0,0.5}; 116 PetscErrorCode ierr; 117 int i; 118 119 PetscFunctionBegin; 120 /* labels of the form n * BASE */ 121 /* get an approximate value for BASE */ 122 base = (vmax - vmin) / (double)(num + 1); 123 124 /* make it of form m x 10^power, m in [1.0, 10) */ 125 if (base <= 0.0) { 126 base = PetscAbsReal(vmin); 127 if (base < 1.0) base = 1.0; 128 } 129 ftemp = PetscLog10Real((1.0 + EPS) * base); 130 if (ftemp < 0.0) ftemp -= 1.0; 131 *power = (int)ftemp; 132 ierr = PetscExp10((double)-*power,&e10);CHKERRQ(ierr); 133 base = base * e10; 134 if (base < 1.0) base = 1.0; 135 /* now reduce it to one of 1, 2, or 5 */ 136 for (i=1; i<5; i++) { 137 if (base >= base_try[i]) { 138 ierr = PetscExp10((double)*power,&e10);CHKERRQ(ierr); 139 base = base_try[i-1] * e10; 140 if (i == 1) *power = *power + 1; 141 break; 142 } 143 } 144 *Base = base; 145 PetscFunctionReturn(0); 146 } 147