1 //===-- High Precision Decimal ----------------------------------*- C++ -*-===// 2 // 3 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. 4 // See httpss//llvm.org/LICENSE.txt for license information. 5 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception 6 // 7 //===----------------------------------------------------------------------===// 8 9 #ifndef LIBC_SRC_SUPPORT_HIGH_PRECISION_DECIMAL_H 10 #define LIBC_SRC_SUPPORT_HIGH_PRECISION_DECIMAL_H 11 12 #include "src/__support/ctype_utils.h" 13 #include "src/__support/str_to_integer.h" 14 #include <stdint.h> 15 16 namespace __llvm_libc { 17 namespace internal { 18 19 struct LShiftTableEntry { 20 uint32_t new_digits; 21 char const *power_of_five; 22 }; 23 24 // This is based on the HPD data structure described as part of the Simple 25 // Decimal Conversion algorithm by Nigel Tao, described at this link: 26 // https://nigeltao.github.io/blog/2020/parse-number-f64-simple.html 27 class HighPrecisionDecimal { 28 29 // This precomputed table speeds up left shifts by having the number of new 30 // digits that will be added by multiplying 5^i by 2^i. If the number is less 31 // than 5^i then it will add one fewer digit. There are only 60 entries since 32 // that's the max shift amount. 33 // This table was generated by the script at 34 // libc/utils/mathtools/GenerateHPDConstants.py 35 static constexpr LShiftTableEntry LEFT_SHIFT_DIGIT_TABLE[] = { 36 {0, ""}, 37 {1, "5"}, 38 {1, "25"}, 39 {1, "125"}, 40 {2, "625"}, 41 {2, "3125"}, 42 {2, "15625"}, 43 {3, "78125"}, 44 {3, "390625"}, 45 {3, "1953125"}, 46 {4, "9765625"}, 47 {4, "48828125"}, 48 {4, "244140625"}, 49 {4, "1220703125"}, 50 {5, "6103515625"}, 51 {5, "30517578125"}, 52 {5, "152587890625"}, 53 {6, "762939453125"}, 54 {6, "3814697265625"}, 55 {6, "19073486328125"}, 56 {7, "95367431640625"}, 57 {7, "476837158203125"}, 58 {7, "2384185791015625"}, 59 {7, "11920928955078125"}, 60 {8, "59604644775390625"}, 61 {8, "298023223876953125"}, 62 {8, "1490116119384765625"}, 63 {9, "7450580596923828125"}, 64 {9, "37252902984619140625"}, 65 {9, "186264514923095703125"}, 66 {10, "931322574615478515625"}, 67 {10, "4656612873077392578125"}, 68 {10, "23283064365386962890625"}, 69 {10, "116415321826934814453125"}, 70 {11, "582076609134674072265625"}, 71 {11, "2910383045673370361328125"}, 72 {11, "14551915228366851806640625"}, 73 {12, "72759576141834259033203125"}, 74 {12, "363797880709171295166015625"}, 75 {12, "1818989403545856475830078125"}, 76 {13, "9094947017729282379150390625"}, 77 {13, "45474735088646411895751953125"}, 78 {13, "227373675443232059478759765625"}, 79 {13, "1136868377216160297393798828125"}, 80 {14, "5684341886080801486968994140625"}, 81 {14, "28421709430404007434844970703125"}, 82 {14, "142108547152020037174224853515625"}, 83 {15, "710542735760100185871124267578125"}, 84 {15, "3552713678800500929355621337890625"}, 85 {15, "17763568394002504646778106689453125"}, 86 {16, "88817841970012523233890533447265625"}, 87 {16, "444089209850062616169452667236328125"}, 88 {16, "2220446049250313080847263336181640625"}, 89 {16, "11102230246251565404236316680908203125"}, 90 {17, "55511151231257827021181583404541015625"}, 91 {17, "277555756156289135105907917022705078125"}, 92 {17, "1387778780781445675529539585113525390625"}, 93 {18, "6938893903907228377647697925567626953125"}, 94 {18, "34694469519536141888238489627838134765625"}, 95 {18, "173472347597680709441192448139190673828125"}, 96 {19, "867361737988403547205962240695953369140625"}, 97 }; 98 99 // The maximum amount we can shift is the number of bits used in the 100 // accumulator, minus the number of bits needed to represent the base (in this 101 // case 4). 102 static constexpr uint32_t MAX_SHIFT_AMOUNT = sizeof(uint64_t) - 4; 103 104 // 800 is an arbitrary number of digits, but should be 105 // large enough for any practical number. 106 static constexpr uint32_t MAX_NUM_DIGITS = 800; 107 108 uint32_t num_digits = 0; 109 int32_t decimal_point = 0; 110 bool truncated = false; 111 uint8_t digits[MAX_NUM_DIGITS]; 112 113 private: should_round_up(int32_t roundToDigit)114 bool should_round_up(int32_t roundToDigit) { 115 if (roundToDigit < 0 || 116 static_cast<uint32_t>(roundToDigit) >= this->num_digits) { 117 return false; 118 } 119 120 // If we're right in the middle and there are no extra digits 121 if (this->digits[roundToDigit] == 5 && 122 static_cast<uint32_t>(roundToDigit + 1) == this->num_digits) { 123 124 // Round up if we've truncated (since that means the result is slightly 125 // higher than what's represented.) 126 if (this->truncated) { 127 return true; 128 } 129 130 // If this exactly halfway, round to even. 131 if (roundToDigit == 0) 132 // When the input is ".5". 133 return false; 134 return this->digits[roundToDigit - 1] % 2 != 0; 135 } 136 // If there are digits after roundToDigit, they must be non-zero since we 137 // trim trailing zeroes after all operations that change digits. 138 return this->digits[roundToDigit] >= 5; 139 } 140 141 // Takes an amount to left shift and returns the number of new digits needed 142 // to store the result based on LEFT_SHIFT_DIGIT_TABLE. get_num_new_digits(uint32_t lShiftAmount)143 uint32_t get_num_new_digits(uint32_t lShiftAmount) { 144 const char *power_of_five = 145 LEFT_SHIFT_DIGIT_TABLE[lShiftAmount].power_of_five; 146 uint32_t new_digits = LEFT_SHIFT_DIGIT_TABLE[lShiftAmount].new_digits; 147 uint32_t digit_index = 0; 148 while (power_of_five[digit_index] != 0) { 149 if (digit_index >= this->num_digits) { 150 return new_digits - 1; 151 } 152 if (this->digits[digit_index] != power_of_five[digit_index] - '0') { 153 return new_digits - 154 ((this->digits[digit_index] < power_of_five[digit_index] - '0') 155 ? 1 156 : 0); 157 } 158 ++digit_index; 159 } 160 return new_digits; 161 } 162 163 // Trim all trailing 0s trim_trailing_zeroes()164 void trim_trailing_zeroes() { 165 while (this->num_digits > 0 && this->digits[this->num_digits - 1] == 0) { 166 --this->num_digits; 167 } 168 if (this->num_digits == 0) { 169 this->decimal_point = 0; 170 } 171 } 172 173 // Perform a digitwise binary non-rounding right shift on this value by 174 // shiftAmount. The shiftAmount can't be more than MAX_SHIFT_AMOUNT to prevent 175 // overflow. right_shift(uint32_t shiftAmount)176 void right_shift(uint32_t shiftAmount) { 177 uint32_t read_index = 0; 178 uint32_t write_index = 0; 179 180 uint64_t accumulator = 0; 181 182 const uint64_t shift_mask = (uint64_t(1) << shiftAmount) - 1; 183 184 // Warm Up phase: we don't have enough digits to start writing, so just 185 // read them into the accumulator. 186 while (accumulator >> shiftAmount == 0) { 187 uint64_t read_digit = 0; 188 // If there are still digits to read, read the next one, else the digit is 189 // assumed to be 0. 190 if (read_index < this->num_digits) { 191 read_digit = this->digits[read_index]; 192 } 193 accumulator = accumulator * 10 + read_digit; 194 ++read_index; 195 } 196 197 // Shift the decimal point by the number of digits it took to fill the 198 // accumulator. 199 this->decimal_point -= read_index - 1; 200 201 // Middle phase: we have enough digits to write, as well as more digits to 202 // read. Keep reading until we run out of digits. 203 while (read_index < this->num_digits) { 204 uint64_t read_digit = this->digits[read_index]; 205 uint64_t write_digit = accumulator >> shiftAmount; 206 accumulator &= shift_mask; 207 this->digits[write_index] = static_cast<uint8_t>(write_digit); 208 accumulator = accumulator * 10 + read_digit; 209 ++read_index; 210 ++write_index; 211 } 212 213 // Cool Down phase: All of the readable digits have been read, so just write 214 // the remainder, while treating any more digits as 0. 215 while (accumulator > 0) { 216 uint64_t write_digit = accumulator >> shiftAmount; 217 accumulator &= shift_mask; 218 if (write_index < MAX_NUM_DIGITS) { 219 this->digits[write_index] = static_cast<uint8_t>(write_digit); 220 ++write_index; 221 } else if (write_digit > 0) { 222 this->truncated = true; 223 } 224 accumulator = accumulator * 10; 225 } 226 this->num_digits = write_index; 227 this->trim_trailing_zeroes(); 228 } 229 230 // Perform a digitwise binary non-rounding left shift on this value by 231 // shiftAmount. The shiftAmount can't be more than MAX_SHIFT_AMOUNT to prevent 232 // overflow. left_shift(uint32_t shiftAmount)233 void left_shift(uint32_t shiftAmount) { 234 uint32_t new_digits = this->get_num_new_digits(shiftAmount); 235 236 int32_t read_index = this->num_digits - 1; 237 uint32_t write_index = this->num_digits + new_digits; 238 239 uint64_t accumulator = 0; 240 241 // No Warm Up phase. Since we're putting digits in at the top and taking 242 // digits from the bottom we don't have to wait for the accumulator to fill. 243 244 // Middle phase: while we have more digits to read, keep reading as well as 245 // writing. 246 while (read_index >= 0) { 247 accumulator += static_cast<uint64_t>(this->digits[read_index]) 248 << shiftAmount; 249 uint64_t next_accumulator = accumulator / 10; 250 uint64_t write_digit = accumulator - (10 * next_accumulator); 251 --write_index; 252 if (write_index < MAX_NUM_DIGITS) { 253 this->digits[write_index] = static_cast<uint8_t>(write_digit); 254 } else if (write_digit != 0) { 255 this->truncated = true; 256 } 257 accumulator = next_accumulator; 258 --read_index; 259 } 260 261 // Cool Down phase: there are no more digits to read, so just write the 262 // remaining digits in the accumulator. 263 while (accumulator > 0) { 264 uint64_t next_accumulator = accumulator / 10; 265 uint64_t write_digit = accumulator - (10 * next_accumulator); 266 --write_index; 267 if (write_index < MAX_NUM_DIGITS) { 268 this->digits[write_index] = static_cast<uint8_t>(write_digit); 269 } else if (write_digit != 0) { 270 this->truncated = true; 271 } 272 accumulator = next_accumulator; 273 } 274 275 this->num_digits += new_digits; 276 if (this->num_digits > MAX_NUM_DIGITS) { 277 this->num_digits = MAX_NUM_DIGITS; 278 } 279 this->decimal_point += new_digits; 280 this->trim_trailing_zeroes(); 281 } 282 283 public: 284 // numString is assumed to be a string of numeric characters. It doesn't 285 // handle leading spaces. HighPrecisionDecimal(const char * __restrict numString)286 HighPrecisionDecimal(const char *__restrict numString) { 287 bool saw_dot = false; 288 // This counts the digits in the number, even if there isn't space to store 289 // them all. 290 uint32_t total_digits = 0; 291 while (isdigit(*numString) || *numString == '.') { 292 if (*numString == '.') { 293 if (saw_dot) { 294 break; 295 } 296 this->decimal_point = total_digits; 297 saw_dot = true; 298 } else { 299 if (*numString == '0' && this->num_digits == 0) { 300 --this->decimal_point; 301 ++numString; 302 continue; 303 } 304 ++total_digits; 305 if (this->num_digits < MAX_NUM_DIGITS) { 306 this->digits[this->num_digits] = 307 static_cast<uint8_t>(*numString - '0'); 308 ++this->num_digits; 309 } else if (*numString != '0') { 310 this->truncated = true; 311 } 312 } 313 ++numString; 314 } 315 316 if (!saw_dot) 317 this->decimal_point = total_digits; 318 319 if ((*numString | 32) == 'e') { 320 ++numString; 321 if (isdigit(*numString) || *numString == '+' || *numString == '-') { 322 int32_t add_to_exp = strtointeger<int32_t>(numString, nullptr, 10); 323 if (add_to_exp > 100000) { 324 add_to_exp = 100000; 325 } else if (add_to_exp < -100000) { 326 add_to_exp = -100000; 327 } 328 this->decimal_point += add_to_exp; 329 } 330 } 331 332 this->trim_trailing_zeroes(); 333 } 334 335 // Binary shift left (shiftAmount > 0) or right (shiftAmount < 0) shift(int shiftAmount)336 void shift(int shiftAmount) { 337 if (shiftAmount == 0) { 338 return; 339 } 340 // Left 341 else if (shiftAmount > 0) { 342 while (static_cast<uint32_t>(shiftAmount) > MAX_SHIFT_AMOUNT) { 343 this->left_shift(MAX_SHIFT_AMOUNT); 344 shiftAmount -= MAX_SHIFT_AMOUNT; 345 } 346 this->left_shift(shiftAmount); 347 } 348 // Right 349 else { 350 while (static_cast<uint32_t>(shiftAmount) < -MAX_SHIFT_AMOUNT) { 351 this->right_shift(MAX_SHIFT_AMOUNT); 352 shiftAmount += MAX_SHIFT_AMOUNT; 353 } 354 this->right_shift(-shiftAmount); 355 } 356 } 357 358 // Round the number represented to the closest value of unsigned int type T. 359 // This is done ignoring overflow. round_to_integer_type()360 template <class T> T round_to_integer_type() { 361 T result = 0; 362 uint32_t cur_digit = 0; 363 364 while (static_cast<int32_t>(cur_digit) < this->decimal_point && 365 cur_digit < this->num_digits) { 366 result = result * 10 + (this->digits[cur_digit]); 367 ++cur_digit; 368 } 369 370 // If there are implicit 0s at the end of the number, include those. 371 while (static_cast<int32_t>(cur_digit) < this->decimal_point) { 372 result *= 10; 373 ++cur_digit; 374 } 375 if (this->should_round_up(this->decimal_point)) { 376 ++result; 377 } 378 return result; 379 } 380 381 // Extra functions for testing. 382 get_digits()383 uint8_t *get_digits() { return this->digits; } get_num_digits()384 uint32_t get_num_digits() { return this->num_digits; } get_decimal_point()385 int32_t get_decimal_point() { return this->decimal_point; } set_truncated(bool trunc)386 void set_truncated(bool trunc) { this->truncated = trunc; } 387 }; 388 389 } // namespace internal 390 } // namespace __llvm_libc 391 392 #endif // LIBC_SRC_SUPPORT_HIGH_PRECISION_DECIMAL_H 393