Google Git
Sign in
skia / external / github.com / google / wuffs / 339119a5c18b1dfcd4e3c7c04fa2685650abe5ab / . / internal / cgen / base / floatconv-submodule-code.c
blob: eaa298241f413ee283bd259c800f924c04dd421c [file] [log] [blame]
// After editing this file, run "go generate" in the ../data directory.
// Copyright 2020 The Wuffs Authors.
//
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// https://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
// ---------------- IEEE 754 Floating Point
WUFFS_BASE__MAYBE_STATIC wuffs_base__lossy_value_u16 //
wuffs_base__ieee_754_bit_representation__from_f64_to_u16_truncate(double f) {
uint64_t u = 0;
if (sizeof(uint64_t) == sizeof(double)) {
memcpy(&u, &f, sizeof(uint64_t));
}
uint16_t neg = ((uint16_t)((u >> 63) << 15));
u &= 0x7FFFFFFFFFFFFFFF;
uint64_t exp = u >> 52;
uint64_t man = u & 0x000FFFFFFFFFFFFF;
if (exp == 0x7FF) {
if (man == 0) { // Infinity.
wuffs_base__lossy_value_u16 ret;
ret.value = neg | 0x7C00;
ret.lossy = false;
return ret;
}
// NaN. Shift the 52 mantissa bits to 10 mantissa bits, keeping the most
// significant mantissa bit (quiet vs signaling NaNs). Also set the low 9
// bits of ret.value so that the 10-bit mantissa is non-zero.
wuffs_base__lossy_value_u16 ret;
ret.value = neg | 0x7DFF | ((uint16_t)(man >> 42));
ret.lossy = false;
return ret;
} else if (exp > 0x40E) { // Truncate to the largest finite f16.
wuffs_base__lossy_value_u16 ret;
ret.value = neg | 0x7BFF;
ret.lossy = true;
return ret;
} else if (exp <= 0x3E6) { // Truncate to zero.
wuffs_base__lossy_value_u16 ret;
ret.value = neg;
ret.lossy = (u != 0);
return ret;
} else if (exp <= 0x3F0) { // Normal f64, subnormal f16.
// Convert from a 53-bit mantissa (after realizing the implicit bit) to a
// 10-bit mantissa and then adjust for the exponent.
man |= 0x0010000000000000;
uint32_t shift = ((uint32_t)(1051 - exp)); // 1051 = 0x3F0 + 53 - 10.
uint64_t shifted_man = man >> shift;
wuffs_base__lossy_value_u16 ret;
ret.value = neg | ((uint16_t)shifted_man);
ret.lossy = (shifted_man << shift) != man;
return ret;
}
// Normal f64, normal f16.
// Re-bias from 1023 to 15 and shift above f16's 10 mantissa bits.
exp = (exp - 1008) << 10; // 1008 = 1023 - 15 = 0x3FF - 0xF.
// Convert from a 52-bit mantissa (excluding the implicit bit) to a 10-bit
// mantissa (again excluding the implicit bit). We lose some information if
// any of the bottom 42 bits are non-zero.
wuffs_base__lossy_value_u16 ret;
ret.value = neg | ((uint16_t)exp) | ((uint16_t)(man >> 42));
ret.lossy = (man << 22) != 0;
return ret;
}
WUFFS_BASE__MAYBE_STATIC wuffs_base__lossy_value_u32 //
wuffs_base__ieee_754_bit_representation__from_f64_to_u32_truncate(double f) {
uint64_t u = 0;
if (sizeof(uint64_t) == sizeof(double)) {
memcpy(&u, &f, sizeof(uint64_t));
}
uint32_t neg = ((uint32_t)(u >> 63)) << 31;
u &= 0x7FFFFFFFFFFFFFFF;
uint64_t exp = u >> 52;
uint64_t man = u & 0x000FFFFFFFFFFFFF;
if (exp == 0x7FF) {
if (man == 0) { // Infinity.
wuffs_base__lossy_value_u32 ret;
ret.value = neg | 0x7F800000;
ret.lossy = false;
return ret;
}
// NaN. Shift the 52 mantissa bits to 23 mantissa bits, keeping the most
// significant mantissa bit (quiet vs signaling NaNs). Also set the low 22
// bits of ret.value so that the 23-bit mantissa is non-zero.
wuffs_base__lossy_value_u32 ret;
ret.value = neg | 0x7FBFFFFF | ((uint32_t)(man >> 29));
ret.lossy = false;
return ret;
} else if (exp > 0x47E) { // Truncate to the largest finite f32.
wuffs_base__lossy_value_u32 ret;
ret.value = neg | 0x7F7FFFFF;
ret.lossy = true;
return ret;
} else if (exp <= 0x369) { // Truncate to zero.
wuffs_base__lossy_value_u32 ret;
ret.value = neg;
ret.lossy = (u != 0);
return ret;
} else if (exp <= 0x380) { // Normal f64, subnormal f32.
// Convert from a 53-bit mantissa (after realizing the implicit bit) to a
// 23-bit mantissa and then adjust for the exponent.
man |= 0x0010000000000000;
uint32_t shift = ((uint32_t)(926 - exp)); // 926 = 0x380 + 53 - 23.
uint64_t shifted_man = man >> shift;
wuffs_base__lossy_value_u32 ret;
ret.value = neg | ((uint32_t)shifted_man);
ret.lossy = (shifted_man << shift) != man;
return ret;
}
// Normal f64, normal f32.
// Re-bias from 1023 to 127 and shift above f32's 23 mantissa bits.
exp = (exp - 896) << 23; // 896 = 1023 - 127 = 0x3FF - 0x7F.
// Convert from a 52-bit mantissa (excluding the implicit bit) to a 23-bit
// mantissa (again excluding the implicit bit). We lose some information if
// any of the bottom 29 bits are non-zero.
wuffs_base__lossy_value_u32 ret;
ret.value = neg | ((uint32_t)exp) | ((uint32_t)(man >> 29));
ret.lossy = (man << 35) != 0;
return ret;
}
// --------
#define WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE 2047
#define WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION 800
// WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL is the largest N
// such that ((10 << N) < (1 << 64)).
#define WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL 60
// wuffs_base__private_implementation__high_prec_dec (abbreviated as HPD) is a
// fixed precision floating point decimal number, augmented with ±infinity
// values, but it cannot represent NaN (Not a Number).
//
// "High precision" means that the mantissa holds 800 decimal digits. 800 is
// WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION.
//
// An HPD isn't for general purpose arithmetic, only for conversions to and
// from IEEE 754 double-precision floating point, where the largest and
// smallest positive, finite values are approximately 1.8e+308 and 4.9e-324.
// HPD exponents above +2047 mean infinity, below -2047 mean zero. The ±2047
// bounds are further away from zero than ±(324 + 800), where 800 and 2047 is
// WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION and
// WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE.
//
// digits[.. num_digits] are the number's digits in big-endian order. The
// uint8_t values are in the range [0 ..= 9], not ['0' ..= '9'], where e.g. '7'
// is the ASCII value 0x37.
//
// decimal_point is the index (within digits) of the decimal point. It may be
// negative or be larger than num_digits, in which case the explicit digits are
// padded with implicit zeroes.
//
// For example, if num_digits is 3 and digits is "\x07\x08\x09":
// - A decimal_point of -2 means ".00789"
// - A decimal_point of -1 means ".0789"
// - A decimal_point of +0 means ".789"
// - A decimal_point of +1 means "7.89"
// - A decimal_point of +2 means "78.9"
// - A decimal_point of +3 means "789."
// - A decimal_point of +4 means "7890."
// - A decimal_point of +5 means "78900."
//
// As above, a decimal_point higher than +2047 means that the overall value is
// infinity, lower than -2047 means zero.
//
// negative is a sign bit. An HPD can distinguish positive and negative zero.
//
// truncated is whether there are more than
// WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION digits, and at
// least one of those extra digits are non-zero. The existence of long-tail
// digits can affect rounding.
//
// The "all fields are zero" value is valid, and represents the number +0.
typedef struct {
uint32_t num_digits;
int32_t decimal_point;
bool negative;
bool truncated;
uint8_t digits[WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION];
} wuffs_base__private_implementation__high_prec_dec;
// wuffs_base__private_implementation__high_prec_dec__trim trims trailing
// zeroes from the h->digits[.. h->num_digits] slice. They have no benefit,
// since we explicitly track h->decimal_point.
//
// Preconditions:
// - h is non-NULL.
static inline void //
wuffs_base__private_implementation__high_prec_dec__trim(
wuffs_base__private_implementation__high_prec_dec* h) {
while ((h->num_digits > 0) && (h->digits[h->num_digits - 1] == 0)) {
h->num_digits--;
}
}
// wuffs_base__private_implementation__high_prec_dec__assign sets h to
// represent the number x.
//
// Preconditions:
// - h is non-NULL.
static void //
wuffs_base__private_implementation__high_prec_dec__assign(
wuffs_base__private_implementation__high_prec_dec* h,
uint64_t x,
bool negative) {
uint32_t n = 0;
// Set h->digits.
if (x > 0) {
// Calculate the digits, working right-to-left. After we determine n (how
// many digits there are), copy from buf to h->digits.
//
// UINT64_MAX, 18446744073709551615, is 20 digits long. It can be faster to
// copy a constant number of bytes than a variable number (20 instead of
// n). Make buf large enough (and start writing to it from the middle) so
// that can we always copy 20 bytes: the slice buf[(20-n) .. (40-n)].
uint8_t buf[40] = {0};
uint8_t* ptr = &buf[20];
do {
uint64_t remaining = x / 10;
x -= remaining * 10;
ptr--;
*ptr = (uint8_t)x;
n++;
x = remaining;
} while (x > 0);
memcpy(h->digits, ptr, 20);
}
// Set h's other fields.
h->num_digits = n;
h->decimal_point = (int32_t)n;
h->negative = negative;
h->truncated = false;
wuffs_base__private_implementation__high_prec_dec__trim(h);
}
static wuffs_base__status //
wuffs_base__private_implementation__high_prec_dec__parse(
wuffs_base__private_implementation__high_prec_dec* h,
wuffs_base__slice_u8 s,
uint32_t options) {
if (!h) {
return wuffs_base__make_status(wuffs_base__error__bad_receiver);
}
h->num_digits = 0;
h->decimal_point = 0;
h->negative = false;
h->truncated = false;
uint8_t* p = s.ptr;
uint8_t* q = s.ptr + s.len;
if (options & WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_UNDERSCORES) {
for (;; p++) {
if (p >= q) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
} else if (*p != '_') {
break;
}
}
}
// Parse sign.
do {
if (*p == '+') {
p++;
} else if (*p == '-') {
h->negative = true;
p++;
} else {
break;
}
if (options & WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_UNDERSCORES) {
for (;; p++) {
if (p >= q) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
} else if (*p != '_') {
break;
}
}
}
} while (0);
// Parse digits, up to (and including) a '.', 'E' or 'e'. Examples for each
// limb in this if-else chain:
// - "0.789"
// - "1002.789"
// - ".789"
// - Other (invalid input).
uint32_t nd = 0;
int32_t dp = 0;
bool no_digits_before_separator = false;
if (('0' == *p) &&
!(options &
WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_MULTIPLE_LEADING_ZEROES)) {
p++;
for (;; p++) {
if (p >= q) {
goto after_all;
} else if (*p ==
((options &
WUFFS_BASE__PARSE_NUMBER_FXX__DECIMAL_SEPARATOR_IS_A_COMMA)
? ','
: '.')) {
p++;
goto after_sep;
} else if ((*p == 'E') || (*p == 'e')) {
p++;
goto after_exp;
} else if ((*p != '_') ||
!(options & WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_UNDERSCORES)) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
}
}
} else if (('0' <= *p) && (*p <= '9')) {
if (*p == '0') {
for (; (p < q) && (*p == '0'); p++) {
}
} else {
h->digits[nd++] = (uint8_t)(*p - '0');
dp = (int32_t)nd;
p++;
}
for (;; p++) {
if (p >= q) {
goto after_all;
} else if (('0' <= *p) && (*p <= '9')) {
if (nd < WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION) {
h->digits[nd++] = (uint8_t)(*p - '0');
dp = (int32_t)nd;
} else if ('0' != *p) {
// Long-tail non-zeroes set the truncated bit.
h->truncated = true;
}
} else if (*p ==
((options &
WUFFS_BASE__PARSE_NUMBER_FXX__DECIMAL_SEPARATOR_IS_A_COMMA)
? ','
: '.')) {
p++;
goto after_sep;
} else if ((*p == 'E') || (*p == 'e')) {
p++;
goto after_exp;
} else if ((*p != '_') ||
!(options & WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_UNDERSCORES)) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
}
}
} else if (*p == ((options &
WUFFS_BASE__PARSE_NUMBER_FXX__DECIMAL_SEPARATOR_IS_A_COMMA)
? ','
: '.')) {
p++;
no_digits_before_separator = true;
} else {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
}
after_sep:
for (;; p++) {
if (p >= q) {
goto after_all;
} else if ('0' == *p) {
if (nd == 0) {
// Track leading zeroes implicitly.
dp--;
} else if (nd <
WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION) {
h->digits[nd++] = (uint8_t)(*p - '0');
}
} else if (('0' < *p) && (*p <= '9')) {
if (nd < WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION) {
h->digits[nd++] = (uint8_t)(*p - '0');
} else {
// Long-tail non-zeroes set the truncated bit.
h->truncated = true;
}
} else if ((*p == 'E') || (*p == 'e')) {
p++;
goto after_exp;
} else if ((*p != '_') ||
!(options & WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_UNDERSCORES)) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
}
}
after_exp:
do {
if (options & WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_UNDERSCORES) {
for (;; p++) {
if (p >= q) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
} else if (*p != '_') {
break;
}
}
}
int32_t exp_sign = +1;
if (*p == '+') {
p++;
} else if (*p == '-') {
exp_sign = -1;
p++;
}
int32_t exp = 0;
const int32_t exp_large =
WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE +
WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION;
bool saw_exp_digits = false;
for (; p < q; p++) {
if ((*p == '_') &&
(options & WUFFS_BASE__PARSE_NUMBER_XXX__ALLOW_UNDERSCORES)) {
// No-op.
} else if (('0' <= *p) && (*p <= '9')) {
saw_exp_digits = true;
if (exp < exp_large) {
exp = (10 * exp) + ((int32_t)(*p - '0'));
}
} else {
break;
}
}
if (!saw_exp_digits) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
}
dp += exp_sign * exp;
} while (0);
after_all:
if (p != q) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
}
h->num_digits = nd;
if (nd == 0) {
if (no_digits_before_separator) {
return wuffs_base__make_status(wuffs_base__error__bad_argument);
}
h->decimal_point = 0;
} else if (dp <
-WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE) {
h->decimal_point =
-WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE - 1;
} else if (dp >
+WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE) {
h->decimal_point =
+WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE + 1;
} else {
h->decimal_point = dp;
}
wuffs_base__private_implementation__high_prec_dec__trim(h);
return wuffs_base__make_status(NULL);
}
// --------
// wuffs_base__private_implementation__high_prec_dec__lshift_num_new_digits
// returns the number of additional decimal digits when left-shifting by shift.
//
// See below for preconditions.
static uint32_t //
wuffs_base__private_implementation__high_prec_dec__lshift_num_new_digits(
wuffs_base__private_implementation__high_prec_dec* h,
uint32_t shift) {
// Masking with 0x3F should be unnecessary (assuming the preconditions) but
// it's cheap and ensures that we don't overflow the
// wuffs_base__private_implementation__hpd_left_shift array.
shift &= 63;
uint32_t x_a = wuffs_base__private_implementation__hpd_left_shift[shift];
uint32_t x_b = wuffs_base__private_implementation__hpd_left_shift[shift + 1];
uint32_t num_new_digits = x_a >> 11;
uint32_t pow5_a = 0x7FF & x_a;
uint32_t pow5_b = 0x7FF & x_b;
const uint8_t* pow5 =
&wuffs_base__private_implementation__powers_of_5[pow5_a];
uint32_t i = 0;
uint32_t n = pow5_b - pow5_a;
for (; i < n; i++) {
if (i >= h->num_digits) {
return num_new_digits - 1;
} else if (h->digits[i] == pow5[i]) {
continue;
} else if (h->digits[i] < pow5[i]) {
return num_new_digits - 1;
} else {
return num_new_digits;
}
}
return num_new_digits;
}
// --------
// wuffs_base__private_implementation__high_prec_dec__rounded_integer returns
// the integral (non-fractional) part of h, provided that it is 18 or fewer
// decimal digits. For 19 or more digits, it returns UINT64_MAX. Note that:
// - (1 << 53) is 9007199254740992, which has 16 decimal digits.
// - (1 << 56) is 72057594037927936, which has 17 decimal digits.
// - (1 << 59) is 576460752303423488, which has 18 decimal digits.
// - (1 << 63) is 9223372036854775808, which has 19 decimal digits.
// and that IEEE 754 double precision has 52 mantissa bits.
//
// That integral part is rounded-to-even: rounding 7.5 or 8.5 both give 8.
//
// h's negative bit is ignored: rounding -8.6 returns 9.
//
// See below for preconditions.
static uint64_t //
wuffs_base__private_implementation__high_prec_dec__rounded_integer(
wuffs_base__private_implementation__high_prec_dec* h) {
if ((h->num_digits == 0) || (h->decimal_point < 0)) {
return 0;
} else if (h->decimal_point > 18) {
return UINT64_MAX;
}
uint32_t dp = (uint32_t)(h->decimal_point);
uint64_t n = 0;
uint32_t i = 0;
for (; i < dp; i++) {
n = (10 * n) + ((i < h->num_digits) ? h->digits[i] : 0);
}
bool round_up = false;
if (dp < h->num_digits) {
round_up = h->digits[dp] >= 5;
if ((h->digits[dp] == 5) && (dp + 1 == h->num_digits)) {
// We are exactly halfway. If we're truncated, round up, otherwise round
// to even.
round_up = h->truncated || //
((dp > 0) && (1 & h->digits[dp - 1]));
}
}
if (round_up) {
n++;
}
return n;
}
// wuffs_base__private_implementation__high_prec_dec__small_xshift shifts h's
// number (where 'x' is 'l' or 'r' for left or right) by a small shift value.
//
// Preconditions:
// - h is non-NULL.
// - h->decimal_point is "not extreme".
// - shift is non-zero.
// - shift is "a small shift".
//
// "Not extreme" means within
// ±WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE.
//
// "A small shift" means not more than
// WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL.
//
// wuffs_base__private_implementation__high_prec_dec__rounded_integer and
// wuffs_base__private_implementation__high_prec_dec__lshift_num_new_digits
// have the same preconditions.
//
// wuffs_base__private_implementation__high_prec_dec__lshift keeps the first
// two preconditions but not the last two. Its shift argument is signed and
// does not need to be "small": zero is a no-op, positive means left shift and
// negative means right shift.
static void //
wuffs_base__private_implementation__high_prec_dec__small_lshift(
wuffs_base__private_implementation__high_prec_dec* h,
uint32_t shift) {
if (h->num_digits == 0) {
return;
}
uint32_t num_new_digits =
wuffs_base__private_implementation__high_prec_dec__lshift_num_new_digits(
h, shift);
uint32_t rx = h->num_digits - 1; // Read index.
uint32_t wx = h->num_digits - 1 + num_new_digits; // Write index.
uint64_t n = 0;
// Repeat: pick up a digit, put down a digit, right to left.
while (((int32_t)rx) >= 0) {
n += ((uint64_t)(h->digits[rx])) << shift;
uint64_t quo = n / 10;
uint64_t rem = n - (10 * quo);
if (wx < WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION) {
h->digits[wx] = (uint8_t)rem;
} else if (rem > 0) {
h->truncated = true;
}
n = quo;
wx--;
rx--;
}
// Put down leading digits, right to left.
while (n > 0) {
uint64_t quo = n / 10;
uint64_t rem = n - (10 * quo);
if (wx < WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION) {
h->digits[wx] = (uint8_t)rem;
} else if (rem > 0) {
h->truncated = true;
}
n = quo;
wx--;
}
// Finish.
h->num_digits += num_new_digits;
if (h->num_digits >
WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION) {
h->num_digits = WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION;
}
h->decimal_point += (int32_t)num_new_digits;
wuffs_base__private_implementation__high_prec_dec__trim(h);
}
static void //
wuffs_base__private_implementation__high_prec_dec__small_rshift(
wuffs_base__private_implementation__high_prec_dec* h,
uint32_t shift) {
uint32_t rx = 0; // Read index.
uint32_t wx = 0; // Write index.
uint64_t n = 0;
// Pick up enough leading digits to cover the first shift.
while ((n >> shift) == 0) {
if (rx < h->num_digits) {
// Read a digit.
n = (10 * n) + h->digits[rx++];
} else if (n == 0) {
// h's number used to be zero and remains zero.
return;
} else {
// Read sufficient implicit trailing zeroes.
while ((n >> shift) == 0) {
n = 10 * n;
rx++;
}
break;
}
}
h->decimal_point -= ((int32_t)(rx - 1));
if (h->decimal_point <
-WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE) {
// After the shift, h's number is effectively zero.
h->num_digits = 0;
h->decimal_point = 0;
h->negative = false;
h->truncated = false;
return;
}
// Repeat: pick up a digit, put down a digit, left to right.
uint64_t mask = (((uint64_t)(1)) << shift) - 1;
while (rx < h->num_digits) {
uint8_t new_digit = ((uint8_t)(n >> shift));
n = (10 * (n & mask)) + h->digits[rx++];
h->digits[wx++] = new_digit;
}
// Put down trailing digits, left to right.
while (n > 0) {
uint8_t new_digit = ((uint8_t)(n >> shift));
n = 10 * (n & mask);
if (wx < WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DIGITS_PRECISION) {
h->digits[wx++] = new_digit;
} else if (new_digit > 0) {
h->truncated = true;
}
}
// Finish.
h->num_digits = wx;
wuffs_base__private_implementation__high_prec_dec__trim(h);
}
static void //
wuffs_base__private_implementation__high_prec_dec__lshift(
wuffs_base__private_implementation__high_prec_dec* h,
int32_t shift) {
if (shift > 0) {
while (shift > +WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL) {
wuffs_base__private_implementation__high_prec_dec__small_lshift(
h, WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL);
shift -= WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL;
}
wuffs_base__private_implementation__high_prec_dec__small_lshift(
h, ((uint32_t)(+shift)));
} else if (shift < 0) {
while (shift < -WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL) {
wuffs_base__private_implementation__high_prec_dec__small_rshift(
h, WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL);
shift += WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL;
}
wuffs_base__private_implementation__high_prec_dec__small_rshift(
h, ((uint32_t)(-shift)));
}
}
// --------
// wuffs_base__private_implementation__high_prec_dec__round_etc rounds h's
// number. For those functions that take an n argument, rounding produces at
// most n digits (which is not necessarily at most n decimal places). Negative
// n values are ignored, as well as any n greater than or equal to h's number
// of digits. The etc__round_just_enough function implicitly chooses an n to
// implement WUFFS_BASE__RENDER_NUMBER_FXX__JUST_ENOUGH_PRECISION.
//
// Preconditions:
// - h is non-NULL.
// - h->decimal_point is "not extreme".
//
// "Not extreme" means within
// ±WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE.
static void //
wuffs_base__private_implementation__high_prec_dec__round_down(
wuffs_base__private_implementation__high_prec_dec* h,
int32_t n) {
if ((n < 0) || (h->num_digits <= (uint32_t)n)) {
return;
}
h->num_digits = (uint32_t)(n);
wuffs_base__private_implementation__high_prec_dec__trim(h);
}
static void //
wuffs_base__private_implementation__high_prec_dec__round_up(
wuffs_base__private_implementation__high_prec_dec* h,
int32_t n) {
if ((n < 0) || (h->num_digits <= (uint32_t)n)) {
return;
}
for (n--; n >= 0; n--) {
if (h->digits[n] < 9) {
h->digits[n]++;
h->num_digits = (uint32_t)(n + 1);
return;
}
}
// The number is all 9s. Change to a single 1 and adjust the decimal point.
h->digits[0] = 1;
h->num_digits = 1;
h->decimal_point++;
}
static void //
wuffs_base__private_implementation__high_prec_dec__round_nearest(
wuffs_base__private_implementation__high_prec_dec* h,
int32_t n) {
if ((n < 0) || (h->num_digits <= (uint32_t)n)) {
return;
}
bool up = h->digits[n] >= 5;
if ((h->digits[n] == 5) && ((n + 1) == ((int32_t)(h->num_digits)))) {
up = h->truncated || //
((n > 0) && ((h->digits[n - 1] & 1) != 0));
}
if (up) {
wuffs_base__private_implementation__high_prec_dec__round_up(h, n);
} else {
wuffs_base__private_implementation__high_prec_dec__round_down(h, n);
}
}
static void //
wuffs_base__private_implementation__high_prec_dec__round_just_enough(
wuffs_base__private_implementation__high_prec_dec* h,
int32_t exp2,
uint64_t mantissa) {
// The magic numbers 52 and 53 in this function are because IEEE 754 double
// precision has 52 mantissa bits.
//
// Let f be the floating point number represented by exp2 and mantissa (and
// also the number in h): the number (mantissa * (2 ** (exp2 - 52))).
//
// If f is zero or a small integer, we can return early.
if ((mantissa == 0) ||
((exp2 < 53) && (h->decimal_point >= ((int32_t)(h->num_digits))))) {
return;
}
// The smallest normal f has an exp2 of -1022 and a mantissa of (1 << 52).
// Subnormal numbers have the same exp2 but a smaller mantissa.
static const int32_t min_incl_normal_exp2 = -1022;
static const uint64_t min_incl_normal_mantissa = 0x0010000000000000ul;
// Compute lower and upper bounds such that any number between them (possibly
// inclusive) will round to f. First, the lower bound. Our number f is:
// ((mantissa + 0) * (2 ** ( exp2 - 52)))
//
// The next lowest floating point number is:
// ((mantissa - 1) * (2 ** ( exp2 - 52)))
// unless (mantissa - 1) drops the (1 << 52) bit and exp2 is not the
// min_incl_normal_exp2. Either way, call it:
// ((l_mantissa) * (2 ** (l_exp2 - 52)))
//
// The lower bound is halfway between them (noting that 52 became 53):
// (((2 * l_mantissa) + 1) * (2 ** (l_exp2 - 53)))
int32_t l_exp2 = exp2;
uint64_t l_mantissa = mantissa - 1;
if ((exp2 > min_incl_normal_exp2) && (mantissa <= min_incl_normal_mantissa)) {
l_exp2 = exp2 - 1;
l_mantissa = (2 * mantissa) - 1;
}
wuffs_base__private_implementation__high_prec_dec lower;
wuffs_base__private_implementation__high_prec_dec__assign(
&lower, (2 * l_mantissa) + 1, false);
wuffs_base__private_implementation__high_prec_dec__lshift(&lower,
l_exp2 - 53);
// Next, the upper bound. Our number f is:
// ((mantissa + 0) * (2 ** (exp2 - 52)))
//
// The next highest floating point number is:
// ((mantissa + 1) * (2 ** (exp2 - 52)))
//
// The upper bound is halfway between them (noting that 52 became 53):
// (((2 * mantissa) + 1) * (2 ** (exp2 - 53)))
wuffs_base__private_implementation__high_prec_dec upper;
wuffs_base__private_implementation__high_prec_dec__assign(
&upper, (2 * mantissa) + 1, false);
wuffs_base__private_implementation__high_prec_dec__lshift(&upper, exp2 - 53);
// The lower and upper bounds are possible outputs only if the original
// mantissa is even, so that IEEE round-to-even would round to the original
// mantissa and not its neighbors.
bool inclusive = (mantissa & 1) == 0;
// As we walk the digits, we want to know whether rounding up would fall
// within the upper bound. This is tracked by upper_delta:
// - When -1, the digits of h and upper are the same so far.
// - When +0, we saw a difference of 1 between h and upper on a previous
// digit and subsequently only 9s for h and 0s for upper. Thus, rounding
// up may fall outside of the bound if !inclusive.
// - When +1, the difference is greater than 1 and we know that rounding up
// falls within the bound.
//
// This is a state machine with three states. The numerical value for each
// state (-1, +0 or +1) isn't important, other than their order.
int upper_delta = -1;
// We can now figure out the shortest number of digits required. Walk the
// digits until h has distinguished itself from lower or upper.
//
// The zi and zd variables are indexes and digits, for z in l (lower), h (the
// number) and u (upper).
//
// The lower, h and upper numbers may have their decimal points at different
// places. In this case, upper is the longest, so we iterate ui starting from
// 0 and iterate li and hi starting from either 0 or -1.
int32_t ui = 0;
for (;; ui++) {
// Calculate hd, the middle number's digit.
int32_t hi = ui - upper.decimal_point + h->decimal_point;
if (hi >= ((int32_t)(h->num_digits))) {
break;
}
uint8_t hd = (((uint32_t)hi) < h->num_digits) ? h->digits[hi] : 0;
// Calculate ld, the lower bound's digit.
int32_t li = ui - upper.decimal_point + lower.decimal_point;
uint8_t ld = (((uint32_t)li) < lower.num_digits) ? lower.digits[li] : 0;
// We can round down (truncate) if lower has a different digit than h or if
// lower is inclusive and is exactly the result of rounding down (i.e. we
// have reached the final digit of lower).
bool can_round_down =
(ld != hd) || //
(inclusive && ((li + 1) == ((int32_t)(lower.num_digits))));
// Calculate ud, the upper bound's digit, and update upper_delta.
uint8_t ud = (((uint32_t)ui) < upper.num_digits) ? upper.digits[ui] : 0;
if (upper_delta < 0) {
if ((hd + 1) < ud) {
// For example:
// h = 12345???
// upper = 12347???
upper_delta = +1;
} else if (hd != ud) {
// For example:
// h = 12345???
// upper = 12346???
upper_delta = +0;
}
} else if (upper_delta == 0) {
if ((hd != 9) || (ud != 0)) {
// For example:
// h = 1234598?
// upper = 1234600?
upper_delta = +1;
}
}
// We can round up if upper has a different digit than h and either upper
// is inclusive or upper is bigger than the result of rounding up.
bool can_round_up =
(upper_delta > 0) || //
((upper_delta == 0) && //
(inclusive || ((ui + 1) < ((int32_t)(upper.num_digits)))));
// If we can round either way, round to nearest. If we can round only one
// way, do it. If we can't round, continue the loop.
if (can_round_down) {
if (can_round_up) {
wuffs_base__private_implementation__high_prec_dec__round_nearest(
h, hi + 1);
return;
} else {
wuffs_base__private_implementation__high_prec_dec__round_down(h,
hi + 1);
return;
}
} else {
if (can_round_up) {
wuffs_base__private_implementation__high_prec_dec__round_up(h, hi + 1);
return;
}
}
}
}
// --------
// wuffs_base__private_implementation__parse_number_f64_eisel_lemire produces
// the IEEE 754 double-precision value for an exact mantissa and base-10
// exponent. For example:
// - when parsing "12345.678e+02", man is 12345678 and exp10 is -1.
// - when parsing "-12", man is 12 and exp10 is 0. Processing the leading
// minus sign is the responsibility of the caller, not this function.
//
// On success, it returns a non-negative int64_t such that the low 63 bits hold
// the 11-bit exponent and 52-bit mantissa.
//
// On failure, it returns a negative value.
//
// The algorithm is based on an original idea by Michael Eisel that was refined
// by Daniel Lemire. See
// https://lemire.me/blog/2020/03/10/fast-float-parsing-in-practice/
//
// Preconditions:
// - man is non-zero.
// - exp10 is in the range -326 ..= 310, the same range of the
// wuffs_base__private_implementation__powers_of_10 array.
static int64_t //
wuffs_base__private_implementation__parse_number_f64_eisel_lemire(
uint64_t man,
int32_t exp10) {
// Look up the (possibly truncated) base-2 representation of (10 ** exp10).
// The look-up table was constructed so that it is already normalized: the
// table entry's mantissa's MSB (most significant bit) is on.
const uint32_t* po10 =
&wuffs_base__private_implementation__powers_of_10[5 * (exp10 + 326)];
// Normalize the man argument. The (man != 0) precondition means that a
// non-zero bit exists.
uint32_t clz = wuffs_base__count_leading_zeroes_u64(man);
man <<= clz;
// Calculate the return value's base-2 exponent. We might tweak it by ±1
// later, but its initial value comes from the look-up table and clz.
uint64_t ret_exp2 = ((uint64_t)po10[4]) - ((uint64_t)clz);
// Multiply the two mantissas. Normalization means that both mantissas are at
// least (1<<63), so the 128-bit product must be at least (1<<126). The high
// 64 bits of the product, x_hi, must therefore be at least (1<<62).
//
// As a consequence, x_hi has either 0 or 1 leading zeroes. Shifting x_hi
// right by either 9 or 10 bits (depending on x_hi's MSB) will therefore
// leave the top 10 MSBs (bits 54 ..= 63) off and the 11th MSB (bit 53) on.
wuffs_base__multiply_u64__output x = wuffs_base__multiply_u64(
man, ((uint64_t)po10[2]) | (((uint64_t)po10[3]) << 32));
uint64_t x_hi = x.hi;
uint64_t x_lo = x.lo;
// Before we shift right by at least 9 bits, recall that the look-up table
// entry was possibly truncated. We have so far only calculated a lower bound
// for the product (man * e), where e is (10 ** exp10). The upper bound would
// add a further (man * 1) to the 128-bit product, which overflows the lower
// 64-bit limb if ((x_lo + man) < man).
//
// If overflow occurs, that adds 1 to x_hi. Since we're about to shift right
// by at least 9 bits, that carried 1 can be ignored unless the higher 64-bit
// limb's low 9 bits are all on.
if (((x_hi & 0x1FF) == 0x1FF) && ((x_lo + man) < man)) {
// Refine our calculation of (man * e). Before, our approximation of e used
// a "low resolution" 64-bit mantissa. Now use a "high resolution" 128-bit
// mantissa. We've already calculated x = (man * bits_0_to_63_incl_of_e).
// Now calculate y = (man * bits_64_to_127_incl_of_e).
wuffs_base__multiply_u64__output y = wuffs_base__multiply_u64(
man, ((uint64_t)po10[0]) | (((uint64_t)po10[1]) << 32));
uint64_t y_hi = y.hi;
uint64_t y_lo = y.lo;
// Merge the 128-bit x and 128-bit y, which overlap by 64 bits, to
// calculate the 192-bit product of the 64-bit man by the 128-bit e.
// As we exit this if-block, we only care about the high 128 bits
// (merged_hi and merged_lo) of that 192-bit product.
uint64_t merged_hi = x_hi;
uint64_t merged_lo = x_lo + y_hi;
if (merged_lo < x_lo) {
merged_hi++; // Carry the overflow bit.
}
// The "high resolution" approximation of e is still a lower bound. Once
// again, see if the upper bound is large enough to produce a different
// result. This time, if it does, give up instead of reaching for an even
// more precise approximation to e.
//
// This three-part check is similar to the two-part check that guarded the
// if block that we're now in, but it has an extra term for the middle 64
// bits (checking that adding 1 to merged_lo would overflow).
if (((merged_hi & 0x1FF) == 0x1FF) && ((merged_lo + 1) == 0) &&
(y_lo + man < man)) {
return -1;
}
// Replace the 128-bit x with merged.
x_hi = merged_hi;
x_lo = merged_lo;
}
// As mentioned above, shifting x_hi right by either 9 or 10 bits will leave
// the top 10 MSBs (bits 54 ..= 63) off and the 11th MSB (bit 53) on. If the
// MSB (before shifting) was on, adjust ret_exp2 for the larger shift.
//
// Having bit 53 on (and higher bits off) means that ret_mantissa is a 54-bit
// number.
uint64_t msb = x_hi >> 63;
uint64_t ret_mantissa = x_hi >> (msb + 9);
ret_exp2 -= 1 ^ msb;
// IEEE 754 rounds to-nearest with ties rounded to-even. Rounding to-even can
// be tricky. If we're half-way between two exactly representable numbers
// (x's low 73 bits are zero and the next 2 bits that matter are "01"), give
// up instead of trying to pick the winner.
//
// Technically, we could tighten the condition by changing "73" to "73 or 74,
// depending on msb", but a flat "73" is simpler.
if ((x_lo == 0) && ((x_hi & 0x1FF) == 0) && ((ret_mantissa & 3) == 1)) {
return -1;
}
// If we're not halfway then it's rounding to-nearest. Starting with a 54-bit
// number, carry the lowest bit (bit 0) up if it's on. Regardless of whether
// it was on or off, shifting right by one then produces a 53-bit number. If
// carrying up overflowed, shift again.
ret_mantissa += ret_mantissa & 1;
ret_mantissa >>= 1;
if ((ret_mantissa >> 53) > 0) {
ret_mantissa >>= 1;
ret_exp2++;
}
// Starting with a 53-bit number, IEEE 754 double-precision normal numbers
// have an implicit mantissa bit. Mask that away and keep the low 52 bits.
ret_mantissa &= 0x000FFFFFFFFFFFFF;
// IEEE 754 double-precision floating point has 11 exponent bits. All off (0)
// means subnormal numbers. All on (2047) means infinity or NaN.
if ((ret_exp2 <= 0) || (2047 <= ret_exp2)) {
return -1;
}
// Pack the bits and return.
return ((int64_t)(ret_mantissa | (ret_exp2 << 52)));
}
// --------
static wuffs_base__result_f64 //
wuffs_base__private_implementation__parse_number_f64_special(
wuffs_base__slice_u8 s,
uint32_t options) {
do {
if (options & WUFFS_BASE__PARSE_NUMBER_FXX__REJECT_INF_AND_NAN) {
goto fail;
}
uint8_t* p = s.ptr;
uint8_t* q = s.ptr + s.len;
for (; (p < q) && (*p == '_'); p++) {
}
if (p >= q) {
goto fail;
}
// Parse sign.
bool negative = false;
do {
if (*p == '+') {
p++;
} else if (*p == '-') {
negative = true;
p++;
} else {
break;
}
for (; (p < q) && (*p == '_'); p++) {
}
} while (0);
if (p >= q) {
goto fail;
}
bool nan = false;
switch (p[0]) {
case 'I':
case 'i':
if (((q - p) < 3) || //
((p[1] != 'N') && (p[1] != 'n')) || //
((p[2] != 'F') && (p[2] != 'f'))) {
goto fail;
}
p += 3;
if ((p >= q) || (*p == '_')) {
break;
} else if (((q - p) < 5) || //
((p[0] != 'I') && (p[0] != 'i')) || //
((p[1] != 'N') && (p[1] != 'n')) || //
((p[2] != 'I') && (p[2] != 'i')) || //
((p[3] != 'T') && (p[3] != 't')) || //
((p[4] != 'Y') && (p[4] != 'y'))) {
goto fail;
}
p += 5;
if ((p >= q) || (*p == '_')) {
break;
}
goto fail;
case 'N':
case 'n':
if (((q - p) < 3) || //
((p[1] != 'A') && (p[1] != 'a')) || //
((p[2] != 'N') && (p[2] != 'n'))) {
goto fail;
}
p += 3;
if ((p >= q) || (*p == '_')) {
nan = true;
break;
}
goto fail;
default:
goto fail;
}
// Finish.
for (; (p < q) && (*p == '_'); p++) {
}
if (p != q) {
goto fail;
}
wuffs_base__result_f64 ret;
ret.status.repr = NULL;
ret.value = wuffs_base__ieee_754_bit_representation__from_u64_to_f64(
(nan ? 0x7FFFFFFFFFFFFFFF : 0x7FF0000000000000) |
(negative ? 0x8000000000000000 : 0));
return ret;
} while (0);
fail:
do {
wuffs_base__result_f64 ret;
ret.status.repr = wuffs_base__error__bad_argument;
ret.value = 0;
return ret;
} while (0);
}
WUFFS_BASE__MAYBE_STATIC wuffs_base__result_f64 //
wuffs_base__private_implementation__high_prec_dec__to_f64(
wuffs_base__private_implementation__high_prec_dec* h,
uint32_t options) {
do {
// powers converts decimal powers of 10 to binary powers of 2. For example,
// (10000 >> 13) is 1. It stops before the elements exceed 60, also known
// as WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL.
static const uint32_t num_powers = 19;
static const uint8_t powers[19] = {
0, 3, 6, 9, 13, 16, 19, 23, 26, 29, //
33, 36, 39, 43, 46, 49, 53, 56, 59, //
};
// Handle zero and obvious extremes. The largest and smallest positive
// finite f64 values are approximately 1.8e+308 and 4.9e-324.
if ((h->num_digits == 0) || (h->decimal_point < -326)) {
goto zero;
} else if (h->decimal_point > 310) {
goto infinity;
}
// Try the fast Eisel-Lemire algorithm again. Calculating the (man, exp10)
// pair from the high_prec_dec h is more correct but slower than the
// approach taken in wuffs_base__parse_number_f64. The latter is optimized
// for the common cases (e.g. assuming no underscores or a leading '+'
// sign) rather than the full set of cases allowed by the Wuffs API.
if (h->num_digits <= 19) {
uint64_t man = 0;
uint32_t i;
for (i = 0; i < h->num_digits; i++) {
man = (10 * man) + h->digits[i];
}
int32_t exp10 = h->decimal_point - ((int32_t)(h->num_digits));
if ((man != 0) && (-326 <= exp10) && (exp10 <= 310)) {
int64_t r =
wuffs_base__private_implementation__parse_number_f64_eisel_lemire(
man, exp10);
if (r >= 0) {
wuffs_base__result_f64 ret;
ret.status.repr = NULL;
ret.value = wuffs_base__ieee_754_bit_representation__from_u64_to_f64(
((uint64_t)r) | (((uint64_t)(h->negative)) << 63));
return ret;
}
}
}
// Scale by powers of 2 until we're in the range [½ .. 1], which gives us
// our exponent (in base-2). First we shift right, possibly a little too
// far, ending with a value certainly below 1 and possibly below ½...
const int32_t f64_bias = -1023;
int32_t exp2 = 0;
while (h->decimal_point > 0) {
uint32_t n = (uint32_t)(+h->decimal_point);
uint32_t shift =
(n < num_powers)
? powers[n]
: WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL;
wuffs_base__private_implementation__high_prec_dec__small_rshift(h, shift);
if (h->decimal_point <
-WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE) {
goto zero;
}
exp2 += (int32_t)shift;
}
// ...then we shift left, putting us in [½ .. 1].
while (h->decimal_point <= 0) {
uint32_t shift;
if (h->decimal_point == 0) {
if (h->digits[0] >= 5) {
break;
}
shift = (h->digits[0] < 2) ? 2 : 1;
} else {
uint32_t n = (uint32_t)(-h->decimal_point);
shift = (n < num_powers)
? powers[n]
: WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL;
}
wuffs_base__private_implementation__high_prec_dec__small_lshift(h, shift);
if (h->decimal_point >
+WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__DECIMAL_POINT__RANGE) {
goto infinity;
}
exp2 -= (int32_t)shift;
}
// We're in the range [½ .. 1] but f64 uses [1 .. 2].
exp2--;
// The minimum normal exponent is (f64_bias + 1).
while ((f64_bias + 1) > exp2) {
uint32_t n = (uint32_t)((f64_bias + 1) - exp2);
if (n > WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL) {
n = WUFFS_BASE__PRIVATE_IMPLEMENTATION__HPD__SHIFT__MAX_INCL;
}
wuffs_base__private_implementation__high_prec_dec__small_rshift(h, n);
exp2 += (int32_t)n;
}
// Check for overflow.
if ((exp2 - f64_bias) >= 0x07FF) { // (1 << 11) - 1.
goto infinity;
}
// Extract 53 bits for the mantissa (in base-2).
wuffs_base__private_implementation__high_prec_dec__small_lshift(h, 53);
uint64_t man2 =
wuffs_base__private_implementation__high_prec_dec__rounded_integer(h);
// Rounding might have added one bit. If so, shift and re-check overflow.
if ((man2 >> 53) != 0) {
man2 >>= 1;
exp2++;
if ((exp2 - f64_bias) >= 0x07FF) { // (1 << 11) - 1.
goto infinity;
}
}
// Handle subnormal numbers.
if ((man2 >> 52) == 0) {
exp2 = f64_bias;
}
// Pack the bits and return.
uint64_t exp2_bits =
(uint64_t)((exp2 - f64_bias) & 0x07FF); // (1 << 11) - 1.
uint64_t bits = (man2 & 0x000FFFFFFFFFFFFF) | // (1 << 52) - 1.
(exp2_bits << 52) | //
(h->negative ? 0x8000000000000000 : 0); // (1 << 63).
wuffs_base__result_f64 ret;
ret.status.repr = NULL;
ret.value = wuffs_base__ieee_754_bit_representation__from_u64_to_f64(bits);
return ret;
} while (0);
zero:
do {
uint64_t bits = h->negative ? 0x8000000000000000 : 0;
wuffs_base__result_f64 ret;
ret.status.repr = NULL;
ret.value = wuffs_base__ieee_754_bit_representation__from_u64_to_f64(bits);
return ret;
} while (0);
infinity:
do {
if (options & WUFFS_BASE__PARSE_NUMBER_FXX__REJECT_INF_AND_NAN) {
wuffs_base__result_f64 ret;
ret.status.repr = wuffs_base__error__bad_argument;
ret.value = 0;
return ret;
}
uint64_t bits = h->negative ? 0xFFF0000000000000 : 0x7FF0000000000000;
wuffs_base__result_f64 ret;
ret.status.repr = NULL;
ret.value = wuffs_base__ieee_754_bit_representation__from_u64_to_f64(bits);
return ret;
} while (0);
}
static inline bool //
wuffs_base__private_implementation__is_decimal_digit(uint8_t c) {
return ('0' <= c) && (c <= '9');
}
WUFFS_BASE__MAYBE_STATIC wuffs_base__result_f64 //
wuffs_base__parse_number_f64(wuffs_base__slice_u8 s, uint32_t options) {
// In practice, almost all "dd.ddddE±xxx" numbers can be represented
// losslessly by a uint64_t mantissa "dddddd" and an int32_t base-10
// exponent, adjusting "xxx" for the position (if present) of the decimal
// separator '.' or ','.
//
// This (u64 man, i32 exp10) data structure is superficially similar to the
// "Do It Yourself Floating Point" type from Loitsch (†), but the exponent
// here is base-10, not base-2.
//
// If s's number fits in a (man, exp10), parse that pair with the
// Eisel-Lemire algorithm. If not, or if Eisel-Lemire fails, parsing s with
// the fallback algorithm is slower but comprehensive.
//
// † "Printing Floating-Point Numbers Quickly and Accurately with Integers"
// (https://www.cs.tufts.edu/~nr/cs257/archive/florian-loitsch/printf.pdf).
// Florian Loitsch is also the primary contributor to
// https://github.com/google/double-conversion
do {
// Calculating that (man, exp10) pair needs to stay within s's bounds.
// Provided that s isn't extremely long, work on a NUL-terminated copy of
// s's contents. The NUL byte isn't a valid part of "±dd.ddddE±xxx".
//
// As the pointer p walks the contents, it's faster to repeatedly check "is
// *p a valid digit" than "is p within bounds and *p a valid digit".
if (s.len >= 256) {
goto fallback;
}
uint8_t z[256];
memcpy(&z[0], s.ptr, s.len);
z[s.len] = 0;
const uint8_t* p = &z[0];
// Look for a leading minus sign. Technically, we could also look for an
// optional plus sign, but the "script/process-json-numbers.c with -p"
// benchmark is noticably slower if we do. It's optional and, in practice,
// usually absent. Let the fallback catch it.
bool negative = (*p == '-');
if (negative) {
p++;
}
// After walking "dd.dddd", comparing p later with p now will produce the
// number of "d"s and "."s.
const uint8_t* const start_of_digits_ptr = p;
// Walk the "d"s before a '.', 'E', NUL byte, etc. If it starts with '0',
// it must be a single '0'. If it starts with a non-zero decimal digit, it
// can be a sequence of decimal digits.
//
// Update the man variable during the walk. It's OK if man overflows now.
// We'll detect that later.
uint64_t man;
if (*p == '0') {
man = 0;
p++;
if (wuffs_base__private_implementation__is_decimal_digit(*p)) {
goto fallback;
}
} else if (wuffs_base__private_implementation__is_decimal_digit(*p)) {
man = ((uint8_t)(*p - '0'));
p++;
for (; wuffs_base__private_implementation__is_decimal_digit(*p); p++) {
man = (10 * man) + ((uint8_t)(*p - '0'));
}
} else {
goto fallback;
}
// Walk the "d"s after the optional decimal separator ('.' or ','),
// updating the man and exp10 variables.
int32_t exp10 = 0;
if (*p ==
((options & WUFFS_BASE__PARSE_NUMBER_FXX__DECIMAL_SEPARATOR_IS_A_COMMA)
? ','
: '.')) {
p++;
const uint8_t* first_after_separator_ptr = p;
if (!wuffs_base__private_implementation__is_decimal_digit(*p)) {
goto fallback;
}
man = (10 * man) + ((uint8_t)(*p - '0'));
p++;
for (; wuffs_base__private_implementation__is_decimal_digit(*p); p++) {
man = (10 * man) + ((uint8_t)(*p - '0'));
}
exp10 = ((int32_t)(first_after_separator_ptr - p));
}
// Count the number of digits:
// - for an input of "314159", digit_count is 6.
// - for an input of "3.14159", digit_count is 7.
//
// This is off-by-one if there is a decimal separator. That's OK for now.
// We'll correct for that later. The "script/process-json-numbers.c with
// -p" benchmark is noticably slower if we try to correct for that now.
uint32_t digit_count = (uint32_t)(p - start_of_digits_ptr);
// Update exp10 for the optional exponent, starting with 'E' or 'e'.
if ((*p | 0x20) == 'e') {
p++;
int32_t exp_sign = +1;
if (*p == '-') {
p++;
exp_sign = -1;
} else if (*p == '+') {
p++;
}
if (!wuffs_base__private_implementation__is_decimal_digit(*p)) {
goto fallback;
}
int32_t exp_num = ((uint8_t)(*p - '0'));
p++;
// The rest of the exp_num walking has a peculiar control flow but, once
// again, the "script/process-json-numbers.c with -p" benchmark is
// sensitive to alternative formulations.
if (wuffs_base__private_implementation__is_decimal_digit(*p)) {
exp_num = (10 * exp_num) + ((uint8_t)(*p - '0'));
p++;
}
if (wuffs_base__private_implementation__is_decimal_digit(*p)) {
exp_num = (10 * exp_num) + ((uint8_t)(*p - '0'));
p++;
}
while (wuffs_base__private_implementation__is_decimal_digit(*p)) {
if (exp_num > 0x1000000) {
goto fallback;
}
exp_num = (10 * exp_num) + ((uint8_t)(*p - '0'));
p++;
}
exp10 += exp_sign * exp_num;
}
// The Wuffs API is that the original slice has no trailing data. It also
// allows underscores, which we don't catch here but the fallback should.
if (p != &z[s.len]) {
goto fallback;
}
// Check that the uint64_t typed man variable has not overflowed, based on
// digit_count.
//
// For reference:
// - (1 << 63) is 9223372036854775808, which has 19 decimal digits.
// - (1 << 64) is 18446744073709551616, which has 20 decimal digits.
// - 19 nines, 9999999999999999999, is 0x8AC7230489E7FFFF, which has 64
// bits and 16 hexadecimal digits.
// - 20 nines, 99999999999999999999, is 0x56BC75E2D630FFFFF, which has 67
// bits and 17 hexadecimal digits.
if (digit_count > 19) {
// Even if we have more than 19 pseudo-digits, it's not yet definitely an
// overflow. Recall that digit_count might be off-by-one (too large) if
// there's a decimal separator. It will also over-report the number of
// meaningful digits if the input looks something like "0.000dddExxx".
//
// We adjust by the number of leading '0's and '.'s and re-compare to 19.
// Once again, technically, we could skip ','s too, but that perturbs the
// "script/process-json-numbers.c with -p" benchmark.
const uint8_t* q = start_of_digits_ptr;
for (; (*q == '0') || (*q == '.'); q++) {
}
digit_count -= (uint32_t)(q - start_of_digits_ptr);
if (digit_count > 19) {
goto fallback;
}
}
// The wuffs_base__private_implementation__parse_number_f64_eisel_lemire
// preconditions include that exp10 is in the range -326 ..= 310.
if ((exp10 < -326) || (310 < exp10)) {
goto fallback;
}
// If both man and (10 ** exp10) are exactly representable by a double, we
// don't need to run the Eisel-Lemire algorithm.
if ((-22 <= exp10) && (exp10 <= 22) && ((man >> 53) == 0)) {
double d = (double)man;
if (exp10 >= 0) {
d *= wuffs_base__private_implementation__f64_powers_of_10[+exp10];
} else {
d /= wuffs_base__private_implementation__f64_powers_of_10[-exp10];
}
wuffs_base__result_f64 ret;
ret.status.repr = NULL;
ret.value = negative ? -d : +d;
return ret;
}
// The wuffs_base__private_implementation__parse_number_f64_eisel_lemire
// preconditions include that man is non-zero. Parsing "0" should be caught
// by the "If both man and (10 ** exp10)" above, but "0e99" might not.
if (man == 0) {
goto fallback;
}
// Our man and exp10 are in range. Run the Eisel-Lemire algorithm.
int64_t r =
wuffs_base__private_implementation__parse_number_f64_eisel_lemire(
man, exp10);
if (r < 0) {
goto fallback;
}
wuffs_base__result_f64 ret;
ret.status.repr = NULL;
ret.value = wuffs_base__ieee_754_bit_representation__from_u64_to_f64(
((uint64_t)r) | (((uint64_t)negative) << 63));
return ret;
} while (0);
fallback:
do {
wuffs_base__private_implementation__high_prec_dec h;
wuffs_base__status status =
wuffs_base__private_implementation__high_prec_dec__parse(&h, s,
options);
if (status.repr) {
return wuffs_base__private_implementation__parse_number_f64_special(
s, options);
}
return wuffs_base__private_implementation__high_prec_dec__to_f64(&h,
options);
} while (0);
}
// --------
static inline size_t //
wuffs_base__private_implementation__render_inf(wuffs_base__slice_u8 dst,
bool neg,
uint32_t options) {
if (neg) {
if (dst.len < 4) {
return 0;
}
wuffs_base__store_u32le__no_bounds_check(dst.ptr, 0x666E492D); // '-Inf'le.
return 4;
}
if (options & WUFFS_BASE__RENDER_NUMBER_XXX__LEADING_PLUS_SIGN) {
if (dst.len < 4) {
return 0;
}
wuffs_base__store_u32le__no_bounds_check(dst.ptr, 0x666E492B); // '+Inf'le.
return 4;
}
if (dst.len < 3) {
return 0;
}
wuffs_base__store_u24le__no_bounds_check(dst.ptr, 0x666E49); // 'Inf'le.
return 3;
}
static inline size_t //
wuffs_base__private_implementation__render_nan(wuffs_base__slice_u8 dst) {
if (dst.len < 3) {
return 0;
}
wuffs_base__store_u24le__no_bounds_check(dst.ptr, 0x4E614E); // 'NaN'le.
return 3;
}
static size_t //
wuffs_base__private_implementation__high_prec_dec__render_exponent_absent(
wuffs_base__slice_u8 dst,
wuffs_base__private_implementation__high_prec_dec* h,
uint32_t precision,
uint32_t options) {
size_t n = (h->negative ||
(options & WUFFS_BASE__RENDER_NUMBER_XXX__LEADING_PLUS_SIGN))
? 1
: 0;
if (h->decimal_point <= 0) {
n += 1;
} else {
n += (size_t)(h->decimal_point);
}
if (precision > 0) {
n += precision + 1; // +1 for the '.'.
}
// Don't modify dst if the formatted number won't fit.
if (n > dst.len) {
return 0;
}
// Align-left or align-right.
uint8_t* ptr = (options & WUFFS_BASE__RENDER_NUMBER_XXX__ALIGN_RIGHT)
? &dst.ptr[dst.len - n]
: &dst.ptr[0];
// Leading "±".
if (h->negative) {
*ptr++ = '-';
} else if (options & WUFFS_BASE__RENDER_NUMBER_XXX__LEADING_PLUS_SIGN) {
*ptr++ = '+';
}
// Integral digits.
if (h->decimal_point <= 0) {
*ptr++ = '0';
} else {
uint32_t m =
wuffs_base__u32__min(h->num_digits, (uint32_t)(h->decimal_point));
uint32_t i = 0;
for (; i < m; i++) {
*ptr++ = (uint8_t)('0' | h->digits[i]);
}
for (; i < (uint32_t)(h->decimal_point); i++) {
*ptr++ = '0';
}
}
// Separator and then fractional digits.
if (precision > 0) {
*ptr++ =
(options & WUFFS_BASE__RENDER_NUMBER_FXX__DECIMAL_SEPARATOR_IS_A_COMMA)
? ','
: '.';
uint32_t i = 0;
for (; i < precision; i++) {
uint32_t j = ((uint32_t)(h->decimal_point)) + i;
*ptr++ = (uint8_t)('0' | ((j < h->num_digits) ? h->digits[j] : 0));
}
}
return n;
}
static size_t //
wuffs_base__private_implementation__high_prec_dec__render_exponent_present(
wuffs_base__slice_u8 dst,
wuffs_base__private_implementation__high_prec_dec* h,
uint32_t precision,
uint32_t options) {
int32_t exp = 0;
if (h->num_digits > 0) {
exp = h->decimal_point - 1;
}
bool negative_exp = exp < 0;
if (negative_exp) {
exp = -exp;
}
size_t n = (h->negative ||
(options & WUFFS_BASE__RENDER_NUMBER_XXX__LEADING_PLUS_SIGN))
? 4
: 3; // Mininum 3 bytes: first digit and then "e±".
if (precision > 0) {
n += precision + 1; // +1 for the '.'.
}
n += (exp < 100) ? 2 : 3;
// Don't modify dst if the formatted number won't fit.
if (n > dst.len) {
return 0;
}
// Align-left or align-right.
uint8_t* ptr = (options & WUFFS_BASE__RENDER_NUMBER_XXX__ALIGN_RIGHT)
? &dst.ptr[dst.len - n]
: &dst.ptr[0];
// Leading "±".
if (h->negative) {
*ptr++ = '-';
} else if (options & WUFFS_BASE__RENDER_NUMBER_XXX__LEADING_PLUS_SIGN) {
*ptr++ = '+';
}
// Integral digit.
if (h->num_digits > 0) {
*ptr++ = (uint8_t)('0' | h->digits[0]);
} else {
*ptr++ = '0';
}
// Separator and then fractional digits.
if (precision > 0) {
*ptr++ =
(options & WUFFS_BASE__RENDER_NUMBER_FXX__DECIMAL_SEPARATOR_IS_A_COMMA)
? ','
: '.';
uint32_t i = 1;
uint32_t j = wuffs_base__u32__min(h->num_digits, precision + 1);
for (; i < j; i++) {
*ptr++ = (uint8_t)('0' | h->digits[i]);
}
for (; i <= precision; i++) {
*ptr++ = '0';
}
}
// Exponent: "e±" and then 2 or 3 digits.
*ptr++ = 'e';
*ptr++ = negative_exp ? '-' : '+';
if (exp < 10) {
*ptr++ = '0';
*ptr++ = (uint8_t)('0' | exp);
} else if (exp < 100) {
*ptr++ = (uint8_t)('0' | (exp / 10));
*ptr++ = (uint8_t)('0' | (exp % 10));
} else {
int32_t e = exp / 100;
exp -= e * 100;
*ptr++ = (uint8_t)('0' | e);
*ptr++ = (uint8_t)('0' | (exp / 10));
*ptr++ = (uint8_t)('0' | (exp % 10));
}
return n;
}
WUFFS_BASE__MAYBE_STATIC size_t //
wuffs_base__render_number_f64(wuffs_base__slice_u8 dst,
double x,
uint32_t precision,
uint32_t options) {
// Decompose x (64 bits) into negativity (1 bit), base-2 exponent (11 bits
// with a -1023 bias) and mantissa (52 bits).
uint64_t bits = wuffs_base__ieee_754_bit_representation__from_f64_to_u64(x);
bool neg = (bits >> 63) != 0;
int32_t exp2 = ((int32_t)(bits >> 52)) & 0x7FF;
uint64_t man = bits & 0x000FFFFFFFFFFFFFul;
// Apply the exponent bias and set the implicit top bit of the mantissa,
// unless x is subnormal. Also take care of Inf and NaN.
if (exp2 == 0x7FF) {
if (man != 0) {
return wuffs_base__private_implementation__render_nan(dst);
}
return wuffs_base__private_implementation__render_inf(dst, neg, options);
} else if (exp2 == 0) {
exp2 = -1022;
} else {
exp2 -= 1023;
man |= 0x0010000000000000ul;
}
// Ensure that precision isn't too large.
if (precision > 4095) {
precision = 4095;
}
// Convert from the (neg, exp2, man) tuple to an HPD.
wuffs_base__private_implementation__high_prec_dec h;
wuffs_base__private_implementation__high_prec_dec__assign(&h, man, neg);
if (h.num_digits > 0) {
wuffs_base__private_implementation__high_prec_dec__lshift(
&h, exp2 - 52); // 52 mantissa bits.
}
// Handle the "%e" and "%f" formats.
switch (options & (WUFFS_BASE__RENDER_NUMBER_FXX__EXPONENT_ABSENT |
WUFFS_BASE__RENDER_NUMBER_FXX__EXPONENT_PRESENT)) {
case WUFFS_BASE__RENDER_NUMBER_FXX__EXPONENT_ABSENT: // The "%"f" format.
if (options & WUFFS_BASE__RENDER_NUMBER_FXX__JUST_ENOUGH_PRECISION) {
wuffs_base__private_implementation__high_prec_dec__round_just_enough(
&h, exp2, man);
int32_t p = ((int32_t)(h.num_digits)) - h.decimal_point;
precision = ((uint32_t)(wuffs_base__i32__max(0, p)));
} else {
wuffs_base__private_implementation__high_prec_dec__round_nearest(
&h, ((int32_t)precision) + h.decimal_point);
}
return wuffs_base__private_implementation__high_prec_dec__render_exponent_absent(
dst, &h, precision, options);
case WUFFS_BASE__RENDER_NUMBER_FXX__EXPONENT_PRESENT: // The "%e" format.
if (options & WUFFS_BASE__RENDER_NUMBER_FXX__JUST_ENOUGH_PRECISION) {
wuffs_base__private_implementation__high_prec_dec__round_just_enough(
&h, exp2, man);
precision = (h.num_digits > 0) ? (h.num_digits - 1) : 0;
} else {
wuffs_base__private_implementation__high_prec_dec__round_nearest(
&h, ((int32_t)precision) + 1);
}
return wuffs_base__private_implementation__high_prec_dec__render_exponent_present(
dst, &h, precision, options);
}
// We have the "%g" format and so precision means the number of significant
// digits, not the number of digits after the decimal separator. Perform
// rounding and determine whether to use "%e" or "%f".
int32_t e_threshold = 0;
if (options & WUFFS_BASE__RENDER_NUMBER_FXX__JUST_ENOUGH_PRECISION) {
wuffs_base__private_implementation__high_prec_dec__round_just_enough(
&h, exp2, man);
precision = h.num_digits;
e_threshold = 6;
} else {
if (precision == 0) {
precision = 1;
}
wuffs_base__private_implementation__high_prec_dec__round_nearest(
&h, ((int32_t)precision));
e_threshold = ((int32_t)precision);
int32_t nd = ((int32_t)(h.num_digits));
if ((e_threshold > nd) && (nd >= h.decimal_point)) {
e_threshold = nd;
}
}
// Use the "%e" format if the exponent is large.
int32_t e = h.decimal_point - 1;
if ((e < -4) || (e_threshold <= e)) {
uint32_t p = wuffs_base__u32__min(precision, h.num_digits);
return wuffs_base__private_implementation__high_prec_dec__render_exponent_present(
dst, &h, (p > 0) ? (p - 1) : 0, options);
}
// Use the "%f" format otherwise.
int32_t p = ((int32_t)precision);
if (p > h.decimal_point) {
p = ((int32_t)(h.num_digits));
}
precision = ((uint32_t)(wuffs_base__i32__max(0, p - h.decimal_point)));
return wuffs_base__private_implementation__high_prec_dec__render_exponent_absent(
dst, &h, precision, options);
}