1/*
2 * Copyright 2021 ByteDance Inc.
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16
17#include "native.h"
18
19/* decimical shift witout overflow, e.g. 9 << 61 overflow */
20#define MAX_SHIFT 60
21
22/* Decimal represent the integer or float
23 * example 1: 1.1 {"11", 2, 1, 0}
24 * example 2: -0.1 {"1", 1, 0, 1}
25 * example 3: 999 {"999", 3, 3, 0}
26 */
27typedef struct Decimal {
28 char* d;
29 size_t cap;
30 int nd;
31 int dp;
32 int neg;
33 int trunc;
34} Decimal;
35
36/* decimal power of ten to binary power of two.
37 * For example: POW_TAB[1]: 10 ** 1 ~ 2 ** 3
38 */
39static const int POW_TAB[9] = {1, 3, 6, 9, 13, 16, 19, 23, 26};
40
41/* Left shift information for decimal.
42 * For example, {2, "625"}. That means that it will add 2 digits to the new decimal
43 * when the prefix of decimal is from "625" to "999", and 1 digit from "0" to "624".
44 */
45typedef struct lshift_cheat {
46 int delta; // number of added digits when left shift
47 const char cutoff[100]; // minus one digit if under the half(cutoff).
48} lshift_cheat;
49
50/* Look up for the decimal shift information by binary shift bits.
51 * idx is shift bits for binary.
52 * value is the shift information for decimal.
53 * For example, idx is 4, the value is {2, "625"}.
54 * That means the binary shift 4 bits left, will cause add 2 digits to the decimal
55 * if the prefix of decimal is under "625".
56 */
57const static lshift_cheat LSHIFT_TAB[61];
58
59static inline void decimal_init(Decimal *d, char *dbuf, size_t cap) {
60 d->d = dbuf;
61 d->cap = cap;
62 for (int i = 0; i < d->cap; ++i) {
63 d->d[i] = 0;
64 }
65 d->dp = 0;
66 d->nd = 0;
67 d->neg = 0;
68 d->trunc = 0;
69}
70
71static inline void decimal_set(Decimal *d, const char *s, ssize_t len, char *dbuf, ssize_t cap) {
72 int i = 0;
73
74 decimal_init(d, dbuf, cap);
75 if (s[i] == '-') {
76 i++;
77 d->neg = 1;
78 }
79
80 int saw_dot = 0;
81 for (; i < len; i++) {
82 if ('0' <= s[i] && s[i] <= '9') {
83 if (s[i] == '0' && d->nd == 0) { // ignore leading zeros
84 d->dp--;
85 continue;
86 }
87 if (d->nd < d->cap) {
88 d->d[d->nd] = s[i];
89 d->nd++;
90 } else if (s[i] != '0') {
91 /* truncat the remaining digits */
92 d->trunc = 1;
93 }
94 } else if (s[i] == '.') {
95 saw_dot = 1;
96 d->dp = d->nd;
97 } else {
98 break;
99 }
100 }
101
102 /* integer */
103 if (saw_dot == 0) {
104 d->dp = d->nd;
105 }
106
107 /* exponent */
108 if (i < len && (s[i] == 'e' || s[i] == 'E')) {
109 int exp = 0;
110 int esgn = 1;
111
112 i++;
113 if (s[i] == '+') {
114 i++;
115 } else if (s[i] == '-') {
116 i++;
117 esgn = -1;
118 }
119
120 for (; i < len && ('0' <= s[i] && s[i] <= '9') && exp < 10000; i++) {
121 exp = exp * 10 + (s[i] - '0');
122 }
123 d->dp += exp * esgn;
124 }
125
126 return;
127}
128
129/* trim trailing zeros from number */
130static inline void trim(Decimal *d) {
131 while (d->nd > 0 && d->d[d->nd - 1] == '0') {
132 d->nd--;
133 }
134 if (d->nd == 0) {
135 d->dp = 0;
136 }
137}
138
139/* Binary shift right (/ 2) by k bits. k <= maxShift to avoid overflow */
140static inline void right_shift(Decimal *d, uint32_t k) {
141 int r = 0; // read pointer
142 int w = 0; // write pointer
143 uint64_t n = 0;
144
145 /* Pick up enough leading digits to cover first shift */
146 for (; n >> k == 0; r++) {
147 if (r >= d->nd) {
148 if (n == 0) {
149 d->nd = 0; // no digits for this num
150 return;
151 }
152 /* until n has enough bits for right shift */
153 while (n >> k == 0) {
154 n *= 10;
155 r++;
156 }
157 break;
158 }
159 n = n * 10 + d->d[r] - '0'; // read the value from d.d
160 }
161 d->dp -= r - 1; // point shift left
162
163 uint64_t mask = (1ull << k) - 1;
164 uint64_t dig = 0;
165
166 /* Pick up a digit, put down a digit */
167 for (; r < d->nd; r++) {
168 dig = n >> k;
169 n &= mask;
170 d->d[w++] = (char)(dig + '0');
171 n = n * 10 + d->d[r] - '0';
172 }
173
174 /* Put down extra digits */
175 while (n > 0) {
176 dig = n >> k;
177 n &= mask;
178 if (w < d->cap) {
179 d->d[w] = (char)(dig + '0');
180 w++;
181 } else if (dig > 0) {
182 /* truncated */
183 d->trunc = 1;
184 }
185 n *= 10;
186 }
187
188 d->nd = w;
189 trim(d);
190}
191
192/* Compare the leading prefix, if b is lexicographically less, return 0 */
193static inline bool prefix_is_less(const char *b, const char *s, uint64_t bn) {
194 int i = 0;
195 for (; i < bn; i++) {
196 if (s[i] == '\0') {
197 return false;
198 }
199 if (b[i] != s[i]) {
200 return b[i] < s[i];
201 }
202 }
203 return s[i] != '\0';
204}
205
206/* Binary shift left (* 2) by k bits. k <= maxShift to avoid overflow */
207static inline void left_shift(Decimal *d, uint32_t k) {
208 int delta = LSHIFT_TAB[k].delta;
209
210 if (prefix_is_less(d->d, LSHIFT_TAB[k].cutoff, d->nd)){
211 delta--;
212 }
213
214 int r = d->nd; // read index
215 int w = d->nd + delta; // write index
216 uint64_t n = 0;
217 uint64_t quo = 0;
218 uint64_t rem = 0;
219
220 /* Pick up a digit, put down a digit */
221 for (r--; r >= 0; r--) {
222 n += (uint64_t)(d->d[r] - '0') << k;
223 quo = n / 10;
224 rem = n - 10 * quo;
225 w--;
226 if (w < d->cap) {
227 d->d[w] = (char)(rem + '0');
228 } else if (rem != 0) {
229 /* truncated */
230 d->trunc = 1;
231 }
232 n = quo;
233 }
234
235 /* Put down extra digits */
236 while (n > 0) {
237 quo = n / 10;
238 rem = n - 10 * quo;
239 w--;
240 if (w < d->cap) {
241 d->d[w] = (char)(rem + '0');
242 } else if (rem != 0) {
243 /* truncated */
244 d->trunc = 1;
245 }
246 n = quo;
247 }
248
249 d->nd += delta;
250 if (d->nd >= d->cap) {
251 d->nd = d->cap;
252 }
253 d->dp += delta;
254 trim(d);
255}
256
257static inline void decimal_shift(Decimal *d, int k) {
258 if (d->nd == 0 || k == 0) {
259 return;
260 }
261
262 if (k > 0) {
263 while (k > MAX_SHIFT) {
264 left_shift(d, MAX_SHIFT);
265 k -= MAX_SHIFT;
266 }
267 if (k) {
268 left_shift(d, k);
269 }
270 }
271
272 if (k < 0) {
273 while (k < -MAX_SHIFT) {
274 right_shift(d, MAX_SHIFT);
275 k += MAX_SHIFT;
276 }
277 if (k) {
278 right_shift(d, -k);
279 }
280 }
281
282}
283
284static inline int should_roundup(Decimal *d, int nd) {
285 if (nd < 0 || nd >= d->nd) {
286 return 0;
287 }
288
289 /* Exactly halfway - round to even */
290 if (d->d[nd] == '5' && nd+1 == d->nd) {
291 if (d->trunc) {
292 return 1;
293 }
294 return nd > 0 && (d->d[nd-1]-'0')%2 != 0;
295 }
296
297 /* not halfway - round to the nearest */
298 return d->d[nd] >= '5';
299}
300
301/* Extract integer part, rounded appropriately */
302static inline uint64_t rounded_integer(Decimal *d) {
303 if (d->dp > 20) { // overflow
304 return 0xFFFFFFFFFFFFFFFF; //64 bits
305 }
306
307 int i = 0;
308 uint64_t n = 0;
309 for (i = 0; i < d->dp && i < d->nd; i++) {
310 n = n * 10 + (d->d[i] - '0');
311 }
312 for (; i < d->dp; i++) {
313 n *= 10;
314 }
315 if (should_roundup(d, d->dp)) {
316 n++;
317 }
318 return n;
319}
320
321int decimal_to_f64(Decimal *d, double *val) {
322 int exp2 = 0;
323 uint64_t mant = 0;
324 uint64_t bits = 0;
325
326 /* d is zero */
327 if (d->nd == 0) {
328 mant = 0;
329 exp2 = -1023;
330 goto out;
331 }
332
333 /* Overflow, return inf/INF */
334 if (d->dp > 310) {
335 goto overflow;
336 }
337 /* Underflow, return zero */
338 if (d->dp < -330) {
339 mant = 0;
340 exp2 = -1023;
341 goto out;
342 }
343
344 /* Scale by powers of two until in range [0.5, 1.0) */
345 int n = 0;
346 while (d->dp > 0) { // d >= 1
347 if (d->dp >= 9) {
348 n = 27;
349 } else {
350 n = POW_TAB[d->dp];
351 }
352 decimal_shift(d, -n); // shift right
353 exp2 += n;
354 }
355 while ((d->dp < 0) || ((d->dp == 0) && (d->d[0] < '5'))) { // d < 0.5
356 if (-d->dp >= 9) {
357 n = 27;
358 } else {
359 n = POW_TAB[-d->dp];
360 }
361 decimal_shift(d, n); // shift left
362 exp2 -= n;
363 }
364
365 /* Our range is [0.5,1) but floating point range is [1,2) */
366 exp2 --;
367
368 /* Minimum exp2 for doulbe is -1022.
369 * If the exponent is smaller, move it up and
370 * adjust d accordingly.
371 */
372 if (exp2 < -1022) {
373 n = -1022 - exp2;
374 decimal_shift(d, -n); // shift right
375 exp2 += n;
376 }
377
378 /* Exp2 too large */
379 if ((exp2 + 1023) >= 0x7FF) {
380 goto overflow;
381 }
382
383 /* Extract 53 bits. */
384 decimal_shift(d, 53); // shift left
385 mant = rounded_integer(d);
386
387 /* Rounding might have added a bit; shift down. */
388 if (mant == (2ull << 52)) { // mant has 54 bits
389 mant >>= 1;
390 exp2 ++;
391 if ((exp2 + 1023) >= 0x7FF) {
392 goto overflow;
393 }
394 }
395
396 /* Denormalized? */
397 if ((mant & (1ull << 52)) == 0) {
398 exp2 = -1023;
399 }
400 goto out;
401
402overflow:
403 /* ±INF/inf */
404 mant = 0;
405 exp2 = 0x7FF - 1023;
406
407out:
408 /* Assemble bits. */
409 bits = mant & 0x000FFFFFFFFFFFFF;
410 bits |= (uint64_t)((exp2 + 1023) & 0x7FF) << 52;
411 if (d->neg) {
412 bits |= 1ull << 63;
413 }
414 *(uint64_t*)val = bits;
415 return 0;
416}
417
418double atof_native(const char *sp, ssize_t nb, char* dbuf, ssize_t cap) {
419 Decimal d;
420 double val = 0;
421 decimal_set(&d, sp, nb, dbuf, cap);
422 decimal_to_f64(&d, &val);
423 return val;
424}
425
426#undef MAX_SHIFT
427
428const static lshift_cheat LSHIFT_TAB[61] = {
429 // Leading digits of 1/2^i = 5^i.
430 // 5^23 is not an exact 64-bit floating point number,
431 // so have to use bc for the math.
432 // Go up to 60 to be large enough for 32bit and 64bit platforms.
433 /*
434 seq 60 | sed 's/^/5^/' | bc |
435 awk 'BEGIN{ print "\t{ 0, \"\" }," }
436 {
437 log2 = log(2)/log(10)
438 printf("\t{ %d, \"%s\" },\t// * %d\n",
439 int(log2*NR+1), $0, 2**NR)
440 }'
441 */
442 {0, ""},
443 {1, "5"}, // * 2
444 {1, "25"}, // * 4
445 {1, "125"}, // * 8
446 {2, "625"}, // * 16
447 {2, "3125"}, // * 32
448 {2, "15625"}, // * 64
449 {3, "78125"}, // * 128
450 {3, "390625"}, // * 256
451 {3, "1953125"}, // * 512
452 {4, "9765625"}, // * 1024
453 {4, "48828125"}, // * 2048
454 {4, "244140625"}, // * 4096
455 {4, "1220703125"}, // * 8192
456 {5, "6103515625"}, // * 16384
457 {5, "30517578125"}, // * 32768
458 {5, "152587890625"}, // * 65536
459 {6, "762939453125"}, // * 131072
460 {6, "3814697265625"}, // * 262144
461 {6, "19073486328125"}, // * 524288
462 {7, "95367431640625"}, // * 1048576
463 {7, "476837158203125"}, // * 2097152
464 {7, "2384185791015625"}, // * 4194304
465 {7, "11920928955078125"}, // * 8388608
466 {8, "59604644775390625"}, // * 16777216
467 {8, "298023223876953125"}, // * 33554432
468 {8, "1490116119384765625"}, // * 67108864
469 {9, "7450580596923828125"}, // * 134217728
470 {9, "37252902984619140625"}, // * 268435456
471 {9, "186264514923095703125"}, // * 536870912
472 {10, "931322574615478515625"}, // * 1073741824
473 {10, "4656612873077392578125"}, // * 2147483648
474 {10, "23283064365386962890625"}, // * 4294967296
475 {10, "116415321826934814453125"}, // * 8589934592
476 {11, "582076609134674072265625"}, // * 17179869184
477 {11, "2910383045673370361328125"}, // * 34359738368
478 {11, "14551915228366851806640625"}, // * 68719476736
479 {12, "72759576141834259033203125"}, // * 137438953472
480 {12, "363797880709171295166015625"}, // * 274877906944
481 {12, "1818989403545856475830078125"}, // * 549755813888
482 {13, "9094947017729282379150390625"}, // * 1099511627776
483 {13, "45474735088646411895751953125"}, // * 2199023255552
484 {13, "227373675443232059478759765625"}, // * 4398046511104
485 {13, "1136868377216160297393798828125"}, // * 8796093022208
486 {14, "5684341886080801486968994140625"}, // * 17592186044416
487 {14, "28421709430404007434844970703125"}, // * 35184372088832
488 {14, "142108547152020037174224853515625"}, // * 70368744177664
489 {15, "710542735760100185871124267578125"}, // * 140737488355328
490 {15, "3552713678800500929355621337890625"}, // * 281474976710656
491 {15, "17763568394002504646778106689453125"}, // * 562949953421312
492 {16, "88817841970012523233890533447265625"}, // * 1125899906842624
493 {16, "444089209850062616169452667236328125"}, // * 2251799813685248
494 {16, "2220446049250313080847263336181640625"}, // * 4503599627370496
495 {16, "11102230246251565404236316680908203125"}, // * 9007199254740992
496 {17, "55511151231257827021181583404541015625"}, // * 18014398509481984
497 {17, "277555756156289135105907917022705078125"}, // * 36028797018963968
498 {17, "1387778780781445675529539585113525390625"}, // * 72057594037927936
499 {18, "6938893903907228377647697925567626953125"}, // * 144115188075855872
500 {18, "34694469519536141888238489627838134765625"}, // * 288230376151711744
501 {18, "173472347597680709441192448139190673828125"}, // * 576460752303423488
502 {19, "867361737988403547205962240695953369140625"}, // * 1152921504606846976
503};
View as plain text