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