1080 lines · c
1//===----------------------------------------------------------------------===//2//3// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.4// See https://llvm.org/LICENSE.txt for license information.5// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception6//7//===----------------------------------------------------------------------===//8 9// Copyright (c) Microsoft Corporation.10// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception11 12// This implementation is dedicated to the memory of Mary and Thavatchai.13 14#ifndef _LIBCPP_SRC_INCLUDE_TO_CHARS_FLOATING_POINT_H15#define _LIBCPP_SRC_INCLUDE_TO_CHARS_FLOATING_POINT_H16 17// Avoid formatting to keep the changes with the original code minimal.18// clang-format off19 20#include <__algorithm/find.h>21#include <__algorithm/find_if.h>22#include <__algorithm/lower_bound.h>23#include <__algorithm/min.h>24#include <__assert>25#include <__config>26#include <__functional/operations.h>27#include <__iterator/access.h>28#include <__iterator/size.h>29#include <bit>30#include <cfloat>31#include <climits>32 33#include "include/ryu/ryu.h"34 35_LIBCPP_BEGIN_NAMESPACE_STD36 37namespace __itoa {38inline constexpr char _Charconv_digits[] = {'0', '1', '2', '3', '4', '5', '6', '7', '8', '9', 'a', 'b', 'c', 'd', 'e',39 'f', 'g', 'h', 'i', 'j', 'k', 'l', 'm', 'n', 'o', 'p', 'q', 'r', 's', 't', 'u', 'v', 'w', 'x', 'y', 'z'};40static_assert(std::size(_Charconv_digits) == 36);41} // __itoa42 43// vvvvvvvvvv DERIVED FROM corecrt_internal_fltintrn.h vvvvvvvvvv44 45template <class _FloatingType>46struct _Floating_type_traits;47 48template <>49struct _Floating_type_traits<float> {50 static constexpr int32_t _Mantissa_bits = FLT_MANT_DIG;51 static constexpr int32_t _Exponent_bits = sizeof(float) * CHAR_BIT - FLT_MANT_DIG;52 53 static constexpr int32_t _Maximum_binary_exponent = FLT_MAX_EXP - 1;54 static constexpr int32_t _Minimum_binary_exponent = FLT_MIN_EXP - 1;55 56 static constexpr int32_t _Exponent_bias = 127;57 58 static constexpr int32_t _Sign_shift = _Exponent_bits + _Mantissa_bits - 1;59 static constexpr int32_t _Exponent_shift = _Mantissa_bits - 1;60 61 using _Uint_type = uint32_t;62 63 static constexpr uint32_t _Exponent_mask = (1u << _Exponent_bits) - 1;64 static constexpr uint32_t _Normal_mantissa_mask = (1u << _Mantissa_bits) - 1;65 static constexpr uint32_t _Denormal_mantissa_mask = (1u << (_Mantissa_bits - 1)) - 1;66 static constexpr uint32_t _Special_nan_mantissa_mask = 1u << (_Mantissa_bits - 2);67 static constexpr uint32_t _Shifted_sign_mask = 1u << _Sign_shift;68 static constexpr uint32_t _Shifted_exponent_mask = _Exponent_mask << _Exponent_shift;69};70 71template <>72struct _Floating_type_traits<double> {73 static constexpr int32_t _Mantissa_bits = DBL_MANT_DIG;74 static constexpr int32_t _Exponent_bits = sizeof(double) * CHAR_BIT - DBL_MANT_DIG;75 76 static constexpr int32_t _Maximum_binary_exponent = DBL_MAX_EXP - 1;77 static constexpr int32_t _Minimum_binary_exponent = DBL_MIN_EXP - 1;78 79 static constexpr int32_t _Exponent_bias = 1023;80 81 static constexpr int32_t _Sign_shift = _Exponent_bits + _Mantissa_bits - 1;82 static constexpr int32_t _Exponent_shift = _Mantissa_bits - 1;83 84 using _Uint_type = uint64_t;85 86 static constexpr uint64_t _Exponent_mask = (1ULL << _Exponent_bits) - 1;87 static constexpr uint64_t _Normal_mantissa_mask = (1ULL << _Mantissa_bits) - 1;88 static constexpr uint64_t _Denormal_mantissa_mask = (1ULL << (_Mantissa_bits - 1)) - 1;89 static constexpr uint64_t _Special_nan_mantissa_mask = 1ULL << (_Mantissa_bits - 2);90 static constexpr uint64_t _Shifted_sign_mask = 1ULL << _Sign_shift;91 static constexpr uint64_t _Shifted_exponent_mask = _Exponent_mask << _Exponent_shift;92};93 94// ^^^^^^^^^^ DERIVED FROM corecrt_internal_fltintrn.h ^^^^^^^^^^95 96// FUNCTION to_chars (FLOATING-POINT TO STRING)97template <class _Floating>98[[nodiscard]] _LIBCPP_HIDE_FROM_ABI99to_chars_result _Floating_to_chars_hex_precision(100 char* _First, char* const _Last, const _Floating _Value, int _Precision) noexcept {101 102 // * Determine the effective _Precision.103 // * Later, we'll decrement _Precision when printing each hexit after the decimal point.104 105 // The hexits after the decimal point correspond to the explicitly stored fraction bits.106 // float explicitly stores 23 fraction bits. 23 / 4 == 5.75, which is 6 hexits.107 // double explicitly stores 52 fraction bits. 52 / 4 == 13, which is 13 hexits.108 constexpr int _Full_precision = _IsSame<_Floating, float>::value ? 6 : 13;109 constexpr int _Adjusted_explicit_bits = _Full_precision * 4;110 111 if (_Precision < 0) {112 // C11 7.21.6.1 "The fprintf function"/5: "A negative precision argument is taken as if the precision were113 // omitted." /8: "if the precision is missing and FLT_RADIX is a power of 2, then the precision is sufficient114 // for an exact representation of the value"115 _Precision = _Full_precision;116 }117 118 // * Extract the _Ieee_mantissa and _Ieee_exponent.119 using _Traits = _Floating_type_traits<_Floating>;120 using _Uint_type = typename _Traits::_Uint_type;121 122 const _Uint_type _Uint_value = std::bit_cast<_Uint_type>(_Value);123 const _Uint_type _Ieee_mantissa = _Uint_value & _Traits::_Denormal_mantissa_mask;124 const int32_t _Ieee_exponent = static_cast<int32_t>(_Uint_value >> _Traits::_Exponent_shift);125 126 // * Prepare the _Adjusted_mantissa. This is aligned to hexit boundaries,127 // * with the implicit bit restored (0 for zero values and subnormal values, 1 for normal values).128 // * Also calculate the _Unbiased_exponent. This unifies the processing of zero, subnormal, and normal values.129 _Uint_type _Adjusted_mantissa;130 131 if constexpr (_IsSame<_Floating, float>::value) {132 _Adjusted_mantissa = _Ieee_mantissa << 1; // align to hexit boundary (23 isn't divisible by 4)133 } else {134 _Adjusted_mantissa = _Ieee_mantissa; // already aligned (52 is divisible by 4)135 }136 137 int32_t _Unbiased_exponent;138 139 if (_Ieee_exponent == 0) { // zero or subnormal140 // implicit bit is 0141 142 if (_Ieee_mantissa == 0) { // zero143 // C11 7.21.6.1 "The fprintf function"/8: "If the value is zero, the exponent is zero."144 _Unbiased_exponent = 0;145 } else { // subnormal146 _Unbiased_exponent = 1 - _Traits::_Exponent_bias;147 }148 } else { // normal149 _Adjusted_mantissa |= _Uint_type{1} << _Adjusted_explicit_bits; // implicit bit is 1150 _Unbiased_exponent = _Ieee_exponent - _Traits::_Exponent_bias;151 }152 153 // _Unbiased_exponent is within [-126, 127] for float, [-1022, 1023] for double.154 155 // * Decompose _Unbiased_exponent into _Sign_character and _Absolute_exponent.156 char _Sign_character;157 uint32_t _Absolute_exponent;158 159 if (_Unbiased_exponent < 0) {160 _Sign_character = '-';161 _Absolute_exponent = static_cast<uint32_t>(-_Unbiased_exponent);162 } else {163 _Sign_character = '+';164 _Absolute_exponent = static_cast<uint32_t>(_Unbiased_exponent);165 }166 167 // _Absolute_exponent is within [0, 127] for float, [0, 1023] for double.168 169 // * Perform a single bounds check.170 {171 int32_t _Exponent_length;172 173 if (_Absolute_exponent < 10) {174 _Exponent_length = 1;175 } else if (_Absolute_exponent < 100) {176 _Exponent_length = 2;177 } else if constexpr (_IsSame<_Floating, float>::value) {178 _Exponent_length = 3;179 } else if (_Absolute_exponent < 1000) {180 _Exponent_length = 3;181 } else {182 _Exponent_length = 4;183 }184 185 // _Precision might be enormous; avoid integer overflow by testing it separately.186 ptrdiff_t _Buffer_size = _Last - _First;187 188 if (_Buffer_size < _Precision) {189 return {_Last, errc::value_too_large};190 }191 192 _Buffer_size -= _Precision;193 194 const int32_t _Length_excluding_precision = 1 // leading hexit195 + static_cast<int32_t>(_Precision > 0) // possible decimal point196 // excluding `+ _Precision`, hexits after decimal point197 + 2 // "p+" or "p-"198 + _Exponent_length; // exponent199 200 if (_Buffer_size < _Length_excluding_precision) {201 return {_Last, errc::value_too_large};202 }203 }204 205 // * Perform rounding when we've been asked to omit hexits.206 if (_Precision < _Full_precision) {207 // _Precision is within [0, 5] for float, [0, 12] for double.208 209 // _Dropped_bits is within [4, 24] for float, [4, 52] for double.210 const int _Dropped_bits = (_Full_precision - _Precision) * 4;211 212 // Perform rounding by adding an appropriately-shifted bit.213 214 // This can propagate carries all the way into the leading hexit. Examples:215 // "0.ff9" rounded to a precision of 2 is "1.00".216 // "1.ff9" rounded to a precision of 2 is "2.00".217 218 // Note that the leading hexit participates in the rounding decision. Examples:219 // "0.8" rounded to a precision of 0 is "0".220 // "1.8" rounded to a precision of 0 is "2".221 222 // Reference implementation with suboptimal codegen:223 // bool _Should_round_up(bool _Lsb_bit, bool _Round_bit, bool _Has_tail_bits) {224 // // If there are no insignificant set bits, the value is exactly-representable and should not be rounded.225 // //226 // // If there are insignificant set bits, we need to round according to round_to_nearest.227 // // We need to handle two cases: we round up if either [1] the value is slightly greater228 // // than the midpoint between two exactly-representable values or [2] the value is exactly the midpoint229 // // between two exactly-representable values and the greater of the two is even (this is "round-to-even").230 // return _Round_bit && (_Has_tail_bits || _Lsb_bit);231 //}232 // const bool _Lsb_bit = (_Adjusted_mantissa & (_Uint_type{1} << _Dropped_bits)) != 0;233 // const bool _Round_bit = (_Adjusted_mantissa & (_Uint_type{1} << (_Dropped_bits - 1))) != 0;234 // const bool _Has_tail_bits = (_Adjusted_mantissa & ((_Uint_type{1} << (_Dropped_bits - 1)) - 1)) != 0;235 // const bool _Should_round = _Should_round_up(_Lsb_bit, _Round_bit, _Has_tail_bits);236 // _Adjusted_mantissa += _Uint_type{_Should_round} << _Dropped_bits;237 238 // Example for optimized implementation: Let _Dropped_bits be 8.239 // Bit index: ...[8]76543210240 // _Adjusted_mantissa: ...[L]RTTTTTTT (not depicting known details, like hexit alignment)241 // By focusing on the bit at index _Dropped_bits, we can avoid unnecessary branching and shifting.242 243 // Bit index: ...[8]76543210244 // _Lsb_bit: ...[L]RTTTTTTT245 const _Uint_type _Lsb_bit = _Adjusted_mantissa;246 247 // Bit index: ...9[8]76543210248 // _Round_bit: ...L[R]TTTTTTT0249 const _Uint_type _Round_bit = _Adjusted_mantissa << 1;250 251 // We can detect (without branching) whether any of the trailing bits are set.252 // Due to _Should_round below, this computation will be used if and only if R is 1, so we can assume that here.253 // Bit index: ...9[8]76543210254 // _Round_bit: ...L[1]TTTTTTT0255 // _Has_tail_bits: ....[H]........256 257 // If all of the trailing bits T are 0, then `_Round_bit - 1` will produce 0 for H (due to R being 1).258 // If any of the trailing bits T are 1, then `_Round_bit - 1` will produce 1 for H (due to R being 1).259 const _Uint_type _Has_tail_bits = _Round_bit - 1;260 261 // Finally, we can use _Should_round_up() logic with bitwise-AND and bitwise-OR,262 // selecting just the bit at index _Dropped_bits. This is the appropriately-shifted bit that we want.263 const _Uint_type _Should_round = _Round_bit & (_Has_tail_bits | _Lsb_bit) & (_Uint_type{1} << _Dropped_bits);264 265 // This rounding technique is dedicated to the memory of Peppermint. =^..^=266 _Adjusted_mantissa += _Should_round;267 }268 269 // * Print the leading hexit, then mask it away.270 {271 const uint32_t _Nibble = static_cast<uint32_t>(_Adjusted_mantissa >> _Adjusted_explicit_bits);272 _LIBCPP_ASSERT_INTERNAL(_Nibble < 3, "");273 const char _Leading_hexit = static_cast<char>('0' + _Nibble);274 275 *_First++ = _Leading_hexit;276 277 constexpr _Uint_type _Mask = (_Uint_type{1} << _Adjusted_explicit_bits) - 1;278 _Adjusted_mantissa &= _Mask;279 }280 281 // * Print the decimal point and trailing hexits.282 283 // C11 7.21.6.1 "The fprintf function"/8:284 // "if the precision is zero and the # flag is not specified, no decimal-point character appears."285 if (_Precision > 0) {286 *_First++ = '.';287 288 int32_t _Number_of_bits_remaining = _Adjusted_explicit_bits; // 24 for float, 52 for double289 290 for (;;) {291 _LIBCPP_ASSERT_INTERNAL(_Number_of_bits_remaining >= 4, "");292 _LIBCPP_ASSERT_INTERNAL(_Number_of_bits_remaining % 4 == 0, "");293 _Number_of_bits_remaining -= 4;294 295 const uint32_t _Nibble = static_cast<uint32_t>(_Adjusted_mantissa >> _Number_of_bits_remaining);296 _LIBCPP_ASSERT_INTERNAL(_Nibble < 16, "");297 const char _Hexit = __itoa::_Charconv_digits[_Nibble];298 299 *_First++ = _Hexit;300 301 // _Precision is the number of hexits that still need to be printed.302 --_Precision;303 if (_Precision == 0) {304 break; // We're completely done with this phase.305 }306 // Otherwise, we need to keep printing hexits.307 308 if (_Number_of_bits_remaining == 0) {309 // We've finished printing _Adjusted_mantissa, so all remaining hexits are '0'.310 std::memset(_First, '0', static_cast<size_t>(_Precision));311 _First += _Precision;312 break;313 }314 315 // Mask away the hexit that we just printed, then keep looping.316 // (We skip this when breaking out of the loop above, because _Adjusted_mantissa isn't used later.)317 const _Uint_type _Mask = (_Uint_type{1} << _Number_of_bits_remaining) - 1;318 _Adjusted_mantissa &= _Mask;319 }320 }321 322 // * Print the exponent.323 324 // C11 7.21.6.1 "The fprintf function"/8: "The exponent always contains at least one digit, and only as many more325 // digits as necessary to represent the decimal exponent of 2."326 327 // Performance note: We should take advantage of the known ranges of possible exponents.328 329 *_First++ = 'p';330 *_First++ = _Sign_character;331 332 // We've already printed '-' if necessary, so uint32_t _Absolute_exponent avoids testing that again.333 return std::to_chars(_First, _Last, _Absolute_exponent);334}335 336template <class _Floating>337[[nodiscard]] _LIBCPP_HIDE_FROM_ABI338to_chars_result _Floating_to_chars_hex_shortest(339 char* _First, char* const _Last, const _Floating _Value) noexcept {340 341 // This prints "1.728p+0" instead of "2.e5p-1".342 // This prints "0.000002p-126" instead of "1p-149" for float.343 // This prints "0.0000000000001p-1022" instead of "1p-1074" for double.344 // This prioritizes being consistent with printf's de facto behavior (and hex-precision's behavior)345 // over minimizing the number of characters printed.346 347 using _Traits = _Floating_type_traits<_Floating>;348 using _Uint_type = typename _Traits::_Uint_type;349 350 const _Uint_type _Uint_value = std::bit_cast<_Uint_type>(_Value);351 352 if (_Uint_value == 0) { // zero detected; write "0p+0" and return353 // C11 7.21.6.1 "The fprintf function"/8: "If the value is zero, the exponent is zero."354 // Special-casing zero is necessary because of the exponent.355 const char* const _Str = "0p+0";356 const size_t _Len = 4;357 358 if (_Last - _First < static_cast<ptrdiff_t>(_Len)) {359 return {_Last, errc::value_too_large};360 }361 362 std::memcpy(_First, _Str, _Len);363 364 return {_First + _Len, errc{}};365 }366 367 const _Uint_type _Ieee_mantissa = _Uint_value & _Traits::_Denormal_mantissa_mask;368 const int32_t _Ieee_exponent = static_cast<int32_t>(_Uint_value >> _Traits::_Exponent_shift);369 370 char _Leading_hexit; // implicit bit371 int32_t _Unbiased_exponent;372 373 if (_Ieee_exponent == 0) { // subnormal374 _Leading_hexit = '0';375 _Unbiased_exponent = 1 - _Traits::_Exponent_bias;376 } else { // normal377 _Leading_hexit = '1';378 _Unbiased_exponent = _Ieee_exponent - _Traits::_Exponent_bias;379 }380 381 // Performance note: Consider avoiding per-character bounds checking when there's plenty of space.382 383 if (_First == _Last) {384 return {_Last, errc::value_too_large};385 }386 387 *_First++ = _Leading_hexit;388 389 if (_Ieee_mantissa == 0) {390 // The fraction bits are all 0. Trim them away, including the decimal point.391 } else {392 if (_First == _Last) {393 return {_Last, errc::value_too_large};394 }395 396 *_First++ = '.';397 398 // The hexits after the decimal point correspond to the explicitly stored fraction bits.399 // float explicitly stores 23 fraction bits. 23 / 4 == 5.75, so we'll print at most 6 hexits.400 // double explicitly stores 52 fraction bits. 52 / 4 == 13, so we'll print at most 13 hexits.401 _Uint_type _Adjusted_mantissa;402 int32_t _Number_of_bits_remaining;403 404 if constexpr (_IsSame<_Floating, float>::value) {405 _Adjusted_mantissa = _Ieee_mantissa << 1; // align to hexit boundary (23 isn't divisible by 4)406 _Number_of_bits_remaining = 24; // 23 fraction bits + 1 alignment bit407 } else {408 _Adjusted_mantissa = _Ieee_mantissa; // already aligned (52 is divisible by 4)409 _Number_of_bits_remaining = 52; // 52 fraction bits410 }411 412 // do-while: The condition _Adjusted_mantissa != 0 is initially true - we have nonzero fraction bits and we've413 // printed the decimal point. Each iteration, we print a hexit, mask it away, and keep looping if we still have414 // nonzero fraction bits. If there would be trailing '0' hexits, this trims them. If there wouldn't be trailing415 // '0' hexits, the same condition works (as we print the final hexit and mask it away); we don't need to test416 // _Number_of_bits_remaining.417 do {418 _LIBCPP_ASSERT_INTERNAL(_Number_of_bits_remaining >= 4, "");419 _LIBCPP_ASSERT_INTERNAL(_Number_of_bits_remaining % 4 == 0, "");420 _Number_of_bits_remaining -= 4;421 422 const uint32_t _Nibble = static_cast<uint32_t>(_Adjusted_mantissa >> _Number_of_bits_remaining);423 _LIBCPP_ASSERT_INTERNAL(_Nibble < 16, "");424 const char _Hexit = __itoa::_Charconv_digits[_Nibble];425 426 if (_First == _Last) {427 return {_Last, errc::value_too_large};428 }429 430 *_First++ = _Hexit;431 432 const _Uint_type _Mask = (_Uint_type{1} << _Number_of_bits_remaining) - 1;433 _Adjusted_mantissa &= _Mask;434 435 } while (_Adjusted_mantissa != 0);436 }437 438 // C11 7.21.6.1 "The fprintf function"/8: "The exponent always contains at least one digit, and only as many more439 // digits as necessary to represent the decimal exponent of 2."440 441 // Performance note: We should take advantage of the known ranges of possible exponents.442 443 // float: _Unbiased_exponent is within [-126, 127].444 // double: _Unbiased_exponent is within [-1022, 1023].445 446 if (_Last - _First < 2) {447 return {_Last, errc::value_too_large};448 }449 450 *_First++ = 'p';451 452 if (_Unbiased_exponent < 0) {453 *_First++ = '-';454 _Unbiased_exponent = -_Unbiased_exponent;455 } else {456 *_First++ = '+';457 }458 459 // We've already printed '-' if necessary, so static_cast<uint32_t> avoids testing that again.460 return std::to_chars(_First, _Last, static_cast<uint32_t>(_Unbiased_exponent));461}462 463// For general precision, we can use lookup tables to avoid performing trial formatting.464 465// For a simple example, imagine counting the number of digits D in an integer, and needing to know466// whether D is less than 3, equal to 3/4/5/6, or greater than 6. We could use a lookup table:467// D | Largest integer with D digits468// 2 | 99469// 3 | 999470// 4 | 9'999471// 5 | 99'999472// 6 | 999'999473// 7 | table end474// Looking up an integer in this table with lower_bound() will work:475// * Too-small integers, like 7, 70, and 99, will cause lower_bound() to return the D == 2 row. (If all we care476// about is whether D is less than 3, then it's okay to smash the D == 1 and D == 2 cases together.)477// * Integers in [100, 999] will cause lower_bound() to return the D == 3 row, and so forth.478// * Too-large integers, like 1'000'000 and above, will cause lower_bound() to return the end of the table. If we479// compute D from that index, this will be considered D == 7, which will activate any "greater than 6" logic.480 481// Floating-point is slightly more complicated.482 483// The ordinary lookup tables are for X within [-5, 38] for float, and [-5, 308] for double.484// (-5 absorbs too-negative exponents, outside the P > X >= -4 criterion. 38 and 308 are the maximum exponents.)485// Due to the P > X condition, we can use a subset of the table for X within [-5, P - 1], suitably clamped.486 487// When P is small, rounding can affect X. For example:488// For P == 1, the largest double with X == 0 is: 9.4999999999999982236431605997495353221893310546875489// For P == 2, the largest double with X == 0 is: 9.949999999999999289457264239899814128875732421875490// For P == 3, the largest double with X == 0 is: 9.9949999999999992184029906638897955417633056640625491 492// Exponent adjustment is a concern for P within [1, 7] for float, and [1, 15] for double (determined via493// brute force). While larger values of P still perform rounding, they can't trigger exponent adjustment.494// This is because only values with repeated '9' digits can undergo exponent adjustment during rounding,495// and floating-point granularity limits the number of consecutive '9' digits that can appear.496 497// So, we need special lookup tables for small values of P.498// These tables have varying lengths due to the P > X >= -4 criterion. For example:499// For P == 1, need table entries for X: -5, -4, -3, -2, -1, 0500// For P == 2, need table entries for X: -5, -4, -3, -2, -1, 0, 1501// For P == 3, need table entries for X: -5, -4, -3, -2, -1, 0, 1, 2502// For P == 4, need table entries for X: -5, -4, -3, -2, -1, 0, 1, 2, 3503 504// We can concatenate these tables for compact storage, using triangular numbers to access them.505// The table for P begins at index (P - 1) * (P + 10) / 2 with length P + 5.506 507// For both the ordinary and special lookup tables, after an index I is returned by lower_bound(), X is I - 5.508 509// We need to special-case the floating-point value 0.0, which is considered to have X == 0.510// Otherwise, the lookup tables would consider it to have a highly negative X.511 512// Finally, because we're working with positive floating-point values,513// representation comparisons behave identically to floating-point comparisons.514 515// The following code generated the lookup tables for the scientific exponent X. Don't remove this code.516#if 0517// cl /EHsc /nologo /W4 /MT /O2 /std:c++17 generate_tables.cpp && generate_tables518 519#include <algorithm>520#include <assert.h>521#include <charconv>522#include <cmath>523#include <limits>524#include <map>525#include <stdint.h>526#include <stdio.h>527#include <system_error>528#include <type_traits>529#include <vector>530using namespace std;531 532template <typename UInt, typename Pred>533[[nodiscard]] UInt uint_partition_point(UInt first, const UInt last, Pred pred) {534 // Find the beginning of the false partition in [first, last).535 // [first, last) is partitioned when all of the true values occur before all of the false values.536 537 static_assert(is_unsigned_v<UInt>);538 assert(first <= last);539 540 for (UInt n = last - first; n > 0;) {541 const UInt n2 = n / 2;542 const UInt mid = first + n2;543 544 if (pred(mid)) {545 first = mid + 1;546 n = n - n2 - 1;547 } else {548 n = n2;549 }550 }551 552 return first;553}554 555template <typename Floating>556[[nodiscard]] int scientific_exponent_X(const int P, const Floating flt) {557 char buf[400]; // more than enough558 559 // C11 7.21.6.1 "The fprintf function"/8 performs trial formatting with scientific precision P - 1.560 const auto to_result = to_chars(buf, end(buf), flt, chars_format::scientific, P - 1);561 assert(to_result.ec == errc{});562 563 const char* exp_ptr = find(buf, to_result.ptr, 'e');564 assert(exp_ptr != to_result.ptr);565 566 ++exp_ptr; // advance past 'e'567 568 if (*exp_ptr == '+') {569 ++exp_ptr; // advance past '+' which from_chars() won't parse570 }571 572 int X;573 const auto from_result = from_chars(exp_ptr, to_result.ptr, X);574 assert(from_result.ec == errc{});575 return X;576}577 578template <typename UInt>579void print_table(const vector<UInt>& v, const char* const name) {580 constexpr const char* UIntName = _IsSame<UInt, uint32_t>::value ? "uint32_t" : "uint64_t";581 582 printf("static constexpr %s %s[%zu] = {\n", UIntName, name, v.size());583 584 for (const auto& val : v) {585 if constexpr (_IsSame<UInt, uint32_t>::value) {586 printf("0x%08Xu,\n", val);587 } else {588 printf("0x%016llXu,\n", val);589 }590 }591 592 printf("};\n");593}594 595enum class Mode { Tables, Tests };596 597template <typename Floating>598void generate_tables(const Mode mode) {599 using Limits = numeric_limits<Floating>;600 using UInt = conditional_t<_IsSame<Floating, float>::value, uint32_t, uint64_t>;601 602 map<int, map<int, UInt>> P_X_LargestValWithX;603 604 constexpr int MaxP = Limits::max_exponent10 + 1; // MaxP performs no rounding during trial formatting605 606 for (int P = 1; P <= MaxP; ++P) {607 for (int X = -5; X < P; ++X) {608 constexpr Floating first = static_cast<Floating>(9e-5); // well below 9.5e-5, otherwise arbitrary609 constexpr Floating last = Limits::infinity(); // one bit above Limits::max()610 611 const UInt val_beyond_X = uint_partition_point(reinterpret_cast<const UInt&>(first),612 reinterpret_cast<const UInt&>(last),613 [P, X](const UInt u) { return scientific_exponent_X(P, reinterpret_cast<const Floating&>(u)) <= X; });614 615 P_X_LargestValWithX[P][X] = val_beyond_X - 1;616 }617 }618 619 constexpr const char* FloatingName = _IsSame<Floating, float>::value ? "float" : "double";620 621 constexpr int MaxSpecialP = _IsSame<Floating, float>::value ? 7 : 15; // MaxSpecialP is affected by exponent adjustment622 623 if (mode == Mode::Tables) {624 printf("template <>\n");625 printf("struct _General_precision_tables<%s> {\n", FloatingName);626 627 printf("static constexpr int _Max_special_P = %d;\n", MaxSpecialP);628 629 vector<UInt> special;630 631 for (int P = 1; P <= MaxSpecialP; ++P) {632 for (int X = -5; X < P; ++X) {633 const UInt val = P_X_LargestValWithX[P][X];634 special.push_back(val);635 }636 }637 638 print_table(special, "_Special_X_table");639 640 for (int P = MaxSpecialP + 1; P < MaxP; ++P) {641 for (int X = -5; X < P; ++X) {642 const UInt val = P_X_LargestValWithX[P][X];643 assert(val == P_X_LargestValWithX[MaxP][X]);644 }645 }646 647 printf("static constexpr int _Max_P = %d;\n", MaxP);648 649 vector<UInt> ordinary;650 651 for (int X = -5; X < MaxP; ++X) {652 const UInt val = P_X_LargestValWithX[MaxP][X];653 ordinary.push_back(val);654 }655 656 print_table(ordinary, "_Ordinary_X_table");657 658 printf("};\n");659 } else {660 printf("==========\n");661 printf("Test cases for %s:\n", FloatingName);662 663 constexpr int Hexits = _IsSame<Floating, float>::value ? 6 : 13;664 constexpr const char* Suffix = _IsSame<Floating, float>::value ? "f" : "";665 666 for (int P = 1; P <= MaxP; ++P) {667 for (int X = -5; X < P; ++X) {668 if (P <= MaxSpecialP || P == 25 || P == MaxP || X == P - 1) {669 const UInt val1 = P_X_LargestValWithX[P][X];670 const Floating f1 = reinterpret_cast<const Floating&>(val1);671 const UInt val2 = val1 + 1;672 const Floating f2 = reinterpret_cast<const Floating&>(val2);673 674 printf("{%.*a%s, chars_format::general, %d, \"%.*g\"},\n", Hexits, f1, Suffix, P, P, f1);675 if (isfinite(f2)) {676 printf("{%.*a%s, chars_format::general, %d, \"%.*g\"},\n", Hexits, f2, Suffix, P, P, f2);677 }678 }679 }680 }681 }682}683 684int main() {685 printf("template <class _Floating>\n");686 printf("struct _General_precision_tables;\n");687 generate_tables<float>(Mode::Tables);688 generate_tables<double>(Mode::Tables);689 generate_tables<float>(Mode::Tests);690 generate_tables<double>(Mode::Tests);691}692#endif // 0693 694template <class _Floating>695struct _General_precision_tables;696 697template <>698struct _General_precision_tables<float> {699 static constexpr int _Max_special_P = 7;700 701 static constexpr uint32_t _Special_X_table[63] = {0x38C73ABCu, 0x3A79096Bu, 0x3C1BA5E3u, 0x3DC28F5Cu, 0x3F733333u,702 0x4117FFFFu, 0x38D0AAA7u, 0x3A826AA8u, 0x3C230553u, 0x3DCBC6A7u, 0x3F7EB851u, 0x411F3333u, 0x42C6FFFFu,703 0x38D19C3Fu, 0x3A8301A7u, 0x3C23C211u, 0x3DCCB295u, 0x3F7FDF3Bu, 0x411FEB85u, 0x42C7E666u, 0x4479DFFFu,704 0x38D1B468u, 0x3A8310C1u, 0x3C23D4F1u, 0x3DCCCA2Du, 0x3F7FFCB9u, 0x411FFDF3u, 0x42C7FD70u, 0x4479FCCCu,705 0x461C3DFFu, 0x38D1B6D2u, 0x3A831243u, 0x3C23D6D4u, 0x3DCCCC89u, 0x3F7FFFACu, 0x411FFFCBu, 0x42C7FFBEu,706 0x4479FFAEu, 0x461C3FCCu, 0x47C34FBFu, 0x38D1B710u, 0x3A83126Au, 0x3C23D704u, 0x3DCCCCC6u, 0x3F7FFFF7u,707 0x411FFFFAu, 0x42C7FFF9u, 0x4479FFF7u, 0x461C3FFAu, 0x47C34FF9u, 0x497423F7u, 0x38D1B716u, 0x3A83126Eu,708 0x3C23D709u, 0x3DCCCCCCu, 0x3F7FFFFFu, 0x411FFFFFu, 0x42C7FFFFu, 0x4479FFFFu, 0x461C3FFFu, 0x47C34FFFu,709 0x497423FFu, 0x4B18967Fu};710 711 static constexpr int _Max_P = 39;712 713 static constexpr uint32_t _Ordinary_X_table[44] = {0x38D1B717u, 0x3A83126Eu, 0x3C23D70Au, 0x3DCCCCCCu, 0x3F7FFFFFu,714 0x411FFFFFu, 0x42C7FFFFu, 0x4479FFFFu, 0x461C3FFFu, 0x47C34FFFu, 0x497423FFu, 0x4B18967Fu, 0x4CBEBC1Fu,715 0x4E6E6B27u, 0x501502F8u, 0x51BA43B7u, 0x5368D4A5u, 0x551184E7u, 0x56B5E620u, 0x58635FA9u, 0x5A0E1BC9u,716 0x5BB1A2BCu, 0x5D5E0B6Bu, 0x5F0AC723u, 0x60AD78EBu, 0x6258D726u, 0x64078678u, 0x65A96816u, 0x6753C21Bu,717 0x69045951u, 0x6AA56FA5u, 0x6C4ECB8Fu, 0x6E013F39u, 0x6FA18F07u, 0x7149F2C9u, 0x72FC6F7Cu, 0x749DC5ADu,718 0x76453719u, 0x77F684DFu, 0x799A130Bu, 0x7B4097CEu, 0x7CF0BDC2u, 0x7E967699u, 0x7F7FFFFFu};719};720 721template <>722struct _General_precision_tables<double> {723 static constexpr int _Max_special_P = 15;724 725 static constexpr uint64_t _Special_X_table[195] = {0x3F18E757928E0C9Du, 0x3F4F212D77318FC5u, 0x3F8374BC6A7EF9DBu,726 0x3FB851EB851EB851u, 0x3FEE666666666666u, 0x4022FFFFFFFFFFFFu, 0x3F1A1554FBDAD751u, 0x3F504D551D68C692u,727 0x3F8460AA64C2F837u, 0x3FB978D4FDF3B645u, 0x3FEFD70A3D70A3D7u, 0x4023E66666666666u, 0x4058DFFFFFFFFFFFu,728 0x3F1A3387ECC8EB96u, 0x3F506034F3FD933Eu, 0x3F84784230FCF80Du, 0x3FB99652BD3C3611u, 0x3FEFFBE76C8B4395u,729 0x4023FD70A3D70A3Du, 0x4058FCCCCCCCCCCCu, 0x408F3BFFFFFFFFFFu, 0x3F1A368D04E0BA6Au, 0x3F506218230C7482u,730 0x3F847A9E2BCF91A3u, 0x3FB99945B6C3760Bu, 0x3FEFFF972474538Eu, 0x4023FFBE76C8B439u, 0x4058FFAE147AE147u,731 0x408F3F9999999999u, 0x40C387BFFFFFFFFFu, 0x3F1A36DA54164F19u, 0x3F506248748DF16Fu, 0x3F847ADA91B16DCBu,732 0x3FB99991361DC93Eu, 0x3FEFFFF583A53B8Eu, 0x4023FFF972474538u, 0x4058FFF7CED91687u, 0x408F3FF5C28F5C28u,733 0x40C387F999999999u, 0x40F869F7FFFFFFFFu, 0x3F1A36E20F35445Du, 0x3F50624D49814ABAu, 0x3F847AE09BE19D69u,734 0x3FB99998C2DA04C3u, 0x3FEFFFFEF39085F4u, 0x4023FFFF583A53B8u, 0x4058FFFF2E48E8A7u, 0x408F3FFEF9DB22D0u,735 0x40C387FF5C28F5C2u, 0x40F869FF33333333u, 0x412E847EFFFFFFFFu, 0x3F1A36E2D51EC34Bu, 0x3F50624DC5333A0Eu,736 0x3F847AE136800892u, 0x3FB9999984200AB7u, 0x3FEFFFFFE5280D65u, 0x4023FFFFEF39085Fu, 0x4058FFFFEB074A77u,737 0x408F3FFFE5C91D14u, 0x40C387FFEF9DB22Du, 0x40F869FFEB851EB8u, 0x412E847FE6666666u, 0x416312CFEFFFFFFFu,738 0x3F1A36E2E8E94FFCu, 0x3F50624DD191D1FDu, 0x3F847AE145F6467Du, 0x3FB999999773D81Cu, 0x3FEFFFFFFD50CE23u,739 0x4023FFFFFE5280D6u, 0x4058FFFFFDE7210Bu, 0x408F3FFFFD60E94Eu, 0x40C387FFFE5C91D1u, 0x40F869FFFDF3B645u,740 0x412E847FFD70A3D7u, 0x416312CFFE666666u, 0x4197D783FDFFFFFFu, 0x3F1A36E2EAE3F7A7u, 0x3F50624DD2CE7AC8u,741 0x3F847AE14782197Bu, 0x3FB9999999629FD9u, 0x3FEFFFFFFFBB47D0u, 0x4023FFFFFFD50CE2u, 0x4058FFFFFFCA501Au,742 0x408F3FFFFFBCE421u, 0x40C387FFFFD60E94u, 0x40F869FFFFCB923Au, 0x412E847FFFBE76C8u, 0x416312CFFFD70A3Du,743 0x4197D783FFCCCCCCu, 0x41CDCD64FFBFFFFFu, 0x3F1A36E2EB16A205u, 0x3F50624DD2EE2543u, 0x3F847AE147A9AE94u,744 0x3FB9999999941A39u, 0x3FEFFFFFFFF920C8u, 0x4023FFFFFFFBB47Du, 0x4058FFFFFFFAA19Cu, 0x408F3FFFFFF94A03u,745 0x40C387FFFFFBCE42u, 0x40F869FFFFFAC1D2u, 0x412E847FFFF97247u, 0x416312CFFFFBE76Cu, 0x4197D783FFFAE147u,746 0x41CDCD64FFF99999u, 0x4202A05F1FFBFFFFu, 0x3F1A36E2EB1BB30Fu, 0x3F50624DD2F14FE9u, 0x3F847AE147ADA3E3u,747 0x3FB9999999990CDCu, 0x3FEFFFFFFFFF5014u, 0x4023FFFFFFFF920Cu, 0x4058FFFFFFFF768Fu, 0x408F3FFFFFFF5433u,748 0x40C387FFFFFF94A0u, 0x40F869FFFFFF79C8u, 0x412E847FFFFF583Au, 0x416312CFFFFF9724u, 0x4197D783FFFF7CEDu,749 0x41CDCD64FFFF5C28u, 0x4202A05F1FFF9999u, 0x42374876E7FF7FFFu, 0x3F1A36E2EB1C34C3u, 0x3F50624DD2F1A0FAu,750 0x3F847AE147AE0938u, 0x3FB9999999998B86u, 0x3FEFFFFFFFFFEE68u, 0x4023FFFFFFFFF501u, 0x4058FFFFFFFFF241u,751 0x408F3FFFFFFFEED1u, 0x40C387FFFFFFF543u, 0x40F869FFFFFFF294u, 0x412E847FFFFFEF39u, 0x416312CFFFFFF583u,752 0x4197D783FFFFF2E4u, 0x41CDCD64FFFFEF9Du, 0x4202A05F1FFFF5C2u, 0x42374876E7FFF333u, 0x426D1A94A1FFEFFFu,753 0x3F1A36E2EB1C41BBu, 0x3F50624DD2F1A915u, 0x3F847AE147AE135Au, 0x3FB9999999999831u, 0x3FEFFFFFFFFFFE3Du,754 0x4023FFFFFFFFFEE6u, 0x4058FFFFFFFFFEA0u, 0x408F3FFFFFFFFE48u, 0x40C387FFFFFFFEEDu, 0x40F869FFFFFFFEA8u,755 0x412E847FFFFFFE52u, 0x416312CFFFFFFEF3u, 0x4197D783FFFFFEB0u, 0x41CDCD64FFFFFE5Cu, 0x4202A05F1FFFFEF9u,756 0x42374876E7FFFEB8u, 0x426D1A94A1FFFE66u, 0x42A2309CE53FFEFFu, 0x3F1A36E2EB1C4307u, 0x3F50624DD2F1A9E4u,757 0x3F847AE147AE145Eu, 0x3FB9999999999975u, 0x3FEFFFFFFFFFFFD2u, 0x4023FFFFFFFFFFE3u, 0x4058FFFFFFFFFFDCu,758 0x408F3FFFFFFFFFD4u, 0x40C387FFFFFFFFE4u, 0x40F869FFFFFFFFDDu, 0x412E847FFFFFFFD5u, 0x416312CFFFFFFFE5u,759 0x4197D783FFFFFFDEu, 0x41CDCD64FFFFFFD6u, 0x4202A05F1FFFFFE5u, 0x42374876E7FFFFDFu, 0x426D1A94A1FFFFD7u,760 0x42A2309CE53FFFE6u, 0x42D6BCC41E8FFFDFu, 0x3F1A36E2EB1C4328u, 0x3F50624DD2F1A9F9u, 0x3F847AE147AE1477u,761 0x3FB9999999999995u, 0x3FEFFFFFFFFFFFFBu, 0x4023FFFFFFFFFFFDu, 0x4058FFFFFFFFFFFCu, 0x408F3FFFFFFFFFFBu,762 0x40C387FFFFFFFFFDu, 0x40F869FFFFFFFFFCu, 0x412E847FFFFFFFFBu, 0x416312CFFFFFFFFDu, 0x4197D783FFFFFFFCu,763 0x41CDCD64FFFFFFFBu, 0x4202A05F1FFFFFFDu, 0x42374876E7FFFFFCu, 0x426D1A94A1FFFFFBu, 0x42A2309CE53FFFFDu,764 0x42D6BCC41E8FFFFCu, 0x430C6BF52633FFFBu};765 766 static constexpr int _Max_P = 309;767 768 static constexpr uint64_t _Ordinary_X_table[314] = {0x3F1A36E2EB1C432Cu, 0x3F50624DD2F1A9FBu, 0x3F847AE147AE147Au,769 0x3FB9999999999999u, 0x3FEFFFFFFFFFFFFFu, 0x4023FFFFFFFFFFFFu, 0x4058FFFFFFFFFFFFu, 0x408F3FFFFFFFFFFFu,770 0x40C387FFFFFFFFFFu, 0x40F869FFFFFFFFFFu, 0x412E847FFFFFFFFFu, 0x416312CFFFFFFFFFu, 0x4197D783FFFFFFFFu,771 0x41CDCD64FFFFFFFFu, 0x4202A05F1FFFFFFFu, 0x42374876E7FFFFFFu, 0x426D1A94A1FFFFFFu, 0x42A2309CE53FFFFFu,772 0x42D6BCC41E8FFFFFu, 0x430C6BF52633FFFFu, 0x4341C37937E07FFFu, 0x4376345785D89FFFu, 0x43ABC16D674EC7FFu,773 0x43E158E460913CFFu, 0x4415AF1D78B58C3Fu, 0x444B1AE4D6E2EF4Fu, 0x4480F0CF064DD591u, 0x44B52D02C7E14AF6u,774 0x44EA784379D99DB4u, 0x45208B2A2C280290u, 0x4554ADF4B7320334u, 0x4589D971E4FE8401u, 0x45C027E72F1F1281u,775 0x45F431E0FAE6D721u, 0x46293E5939A08CE9u, 0x465F8DEF8808B024u, 0x4693B8B5B5056E16u, 0x46C8A6E32246C99Cu,776 0x46FED09BEAD87C03u, 0x4733426172C74D82u, 0x476812F9CF7920E2u, 0x479E17B84357691Bu, 0x47D2CED32A16A1B1u,777 0x48078287F49C4A1Du, 0x483D6329F1C35CA4u, 0x48725DFA371A19E6u, 0x48A6F578C4E0A060u, 0x48DCB2D6F618C878u,778 0x4911EFC659CF7D4Bu, 0x49466BB7F0435C9Eu, 0x497C06A5EC5433C6u, 0x49B18427B3B4A05Bu, 0x49E5E531A0A1C872u,779 0x4A1B5E7E08CA3A8Fu, 0x4A511B0EC57E6499u, 0x4A8561D276DDFDC0u, 0x4ABABA4714957D30u, 0x4AF0B46C6CDD6E3Eu,780 0x4B24E1878814C9CDu, 0x4B5A19E96A19FC40u, 0x4B905031E2503DA8u, 0x4BC4643E5AE44D12u, 0x4BF97D4DF19D6057u,781 0x4C2FDCA16E04B86Du, 0x4C63E9E4E4C2F344u, 0x4C98E45E1DF3B015u, 0x4CCF1D75A5709C1Au, 0x4D03726987666190u,782 0x4D384F03E93FF9F4u, 0x4D6E62C4E38FF872u, 0x4DA2FDBB0E39FB47u, 0x4DD7BD29D1C87A19u, 0x4E0DAC74463A989Fu,783 0x4E428BC8ABE49F63u, 0x4E772EBAD6DDC73Cu, 0x4EACFA698C95390Bu, 0x4EE21C81F7DD43A7u, 0x4F16A3A275D49491u,784 0x4F4C4C8B1349B9B5u, 0x4F81AFD6EC0E1411u, 0x4FB61BCCA7119915u, 0x4FEBA2BFD0D5FF5Bu, 0x502145B7E285BF98u,785 0x50559725DB272F7Fu, 0x508AFCEF51F0FB5Eu, 0x50C0DE1593369D1Bu, 0x50F5159AF8044462u, 0x512A5B01B605557Au,786 0x516078E111C3556Cu, 0x5194971956342AC7u, 0x51C9BCDFABC13579u, 0x5200160BCB58C16Cu, 0x52341B8EBE2EF1C7u,787 0x526922726DBAAE39u, 0x529F6B0F092959C7u, 0x52D3A2E965B9D81Cu, 0x53088BA3BF284E23u, 0x533EAE8CAEF261ACu,788 0x53732D17ED577D0Bu, 0x53A7F85DE8AD5C4Eu, 0x53DDF67562D8B362u, 0x5412BA095DC7701Du, 0x5447688BB5394C25u,789 0x547D42AEA2879F2Eu, 0x54B249AD2594C37Cu, 0x54E6DC186EF9F45Cu, 0x551C931E8AB87173u, 0x5551DBF316B346E7u,790 0x558652EFDC6018A1u, 0x55BBE7ABD3781ECAu, 0x55F170CB642B133Eu, 0x5625CCFE3D35D80Eu, 0x565B403DCC834E11u,791 0x569108269FD210CBu, 0x56C54A3047C694FDu, 0x56FA9CBC59B83A3Du, 0x5730A1F5B8132466u, 0x5764CA732617ED7Fu,792 0x5799FD0FEF9DE8DFu, 0x57D03E29F5C2B18Bu, 0x58044DB473335DEEu, 0x583961219000356Au, 0x586FB969F40042C5u,793 0x58A3D3E2388029BBu, 0x58D8C8DAC6A0342Au, 0x590EFB1178484134u, 0x59435CEAEB2D28C0u, 0x59783425A5F872F1u,794 0x59AE412F0F768FADu, 0x59E2E8BD69AA19CCu, 0x5A17A2ECC414A03Fu, 0x5A4D8BA7F519C84Fu, 0x5A827748F9301D31u,795 0x5AB7151B377C247Eu, 0x5AECDA62055B2D9Du, 0x5B22087D4358FC82u, 0x5B568A9C942F3BA3u, 0x5B8C2D43B93B0A8Bu,796 0x5BC19C4A53C4E697u, 0x5BF6035CE8B6203Du, 0x5C2B843422E3A84Cu, 0x5C6132A095CE492Fu, 0x5C957F48BB41DB7Bu,797 0x5CCADF1AEA12525Au, 0x5D00CB70D24B7378u, 0x5D34FE4D06DE5056u, 0x5D6A3DE04895E46Cu, 0x5DA066AC2D5DAEC3u,798 0x5DD4805738B51A74u, 0x5E09A06D06E26112u, 0x5E400444244D7CABu, 0x5E7405552D60DBD6u, 0x5EA906AA78B912CBu,799 0x5EDF485516E7577Eu, 0x5F138D352E5096AFu, 0x5F48708279E4BC5Au, 0x5F7E8CA3185DEB71u, 0x5FB317E5EF3AB327u,800 0x5FE7DDDF6B095FF0u, 0x601DD55745CBB7ECu, 0x6052A5568B9F52F4u, 0x60874EAC2E8727B1u, 0x60BD22573A28F19Du,801 0x60F2357684599702u, 0x6126C2D4256FFCC2u, 0x615C73892ECBFBF3u, 0x6191C835BD3F7D78u, 0x61C63A432C8F5CD6u,802 0x61FBC8D3F7B3340Bu, 0x62315D847AD00087u, 0x6265B4E5998400A9u, 0x629B221EFFE500D3u, 0x62D0F5535FEF2084u,803 0x630532A837EAE8A5u, 0x633A7F5245E5A2CEu, 0x63708F936BAF85C1u, 0x63A4B378469B6731u, 0x63D9E056584240FDu,804 0x64102C35F729689Eu, 0x6444374374F3C2C6u, 0x647945145230B377u, 0x64AF965966BCE055u, 0x64E3BDF7E0360C35u,805 0x6518AD75D8438F43u, 0x654ED8D34E547313u, 0x6583478410F4C7ECu, 0x65B819651531F9E7u, 0x65EE1FBE5A7E7861u,806 0x6622D3D6F88F0B3Cu, 0x665788CCB6B2CE0Cu, 0x668D6AFFE45F818Fu, 0x66C262DFEEBBB0F9u, 0x66F6FB97EA6A9D37u,807 0x672CBA7DE5054485u, 0x6761F48EAF234AD3u, 0x679671B25AEC1D88u, 0x67CC0E1EF1A724EAu, 0x680188D357087712u,808 0x6835EB082CCA94D7u, 0x686B65CA37FD3A0Du, 0x68A11F9E62FE4448u, 0x68D56785FBBDD55Au, 0x690AC1677AAD4AB0u,809 0x6940B8E0ACAC4EAEu, 0x6974E718D7D7625Au, 0x69AA20DF0DCD3AF0u, 0x69E0548B68A044D6u, 0x6A1469AE42C8560Cu,810 0x6A498419D37A6B8Fu, 0x6A7FE52048590672u, 0x6AB3EF342D37A407u, 0x6AE8EB0138858D09u, 0x6B1F25C186A6F04Cu,811 0x6B537798F428562Fu, 0x6B88557F31326BBBu, 0x6BBE6ADEFD7F06AAu, 0x6BF302CB5E6F642Au, 0x6C27C37E360B3D35u,812 0x6C5DB45DC38E0C82u, 0x6C9290BA9A38C7D1u, 0x6CC734E940C6F9C5u, 0x6CFD022390F8B837u, 0x6D3221563A9B7322u,813 0x6D66A9ABC9424FEBu, 0x6D9C5416BB92E3E6u, 0x6DD1B48E353BCE6Fu, 0x6E0621B1C28AC20Bu, 0x6E3BAA1E332D728Eu,814 0x6E714A52DFFC6799u, 0x6EA59CE797FB817Fu, 0x6EDB04217DFA61DFu, 0x6F10E294EEBC7D2Bu, 0x6F451B3A2A6B9C76u,815 0x6F7A6208B5068394u, 0x6FB07D457124123Cu, 0x6FE49C96CD6D16CBu, 0x7019C3BC80C85C7Eu, 0x70501A55D07D39CFu,816 0x708420EB449C8842u, 0x70B9292615C3AA53u, 0x70EF736F9B3494E8u, 0x7123A825C100DD11u, 0x7158922F31411455u,817 0x718EB6BAFD91596Bu, 0x71C33234DE7AD7E2u, 0x71F7FEC216198DDBu, 0x722DFE729B9FF152u, 0x7262BF07A143F6D3u,818 0x72976EC98994F488u, 0x72CD4A7BEBFA31AAu, 0x73024E8D737C5F0Au, 0x7336E230D05B76CDu, 0x736C9ABD04725480u,819 0x73A1E0B622C774D0u, 0x73D658E3AB795204u, 0x740BEF1C9657A685u, 0x74417571DDF6C813u, 0x7475D2CE55747A18u,820 0x74AB4781EAD1989Eu, 0x74E10CB132C2FF63u, 0x75154FDD7F73BF3Bu, 0x754AA3D4DF50AF0Au, 0x7580A6650B926D66u,821 0x75B4CFFE4E7708C0u, 0x75EA03FDE214CAF0u, 0x7620427EAD4CFED6u, 0x7654531E58A03E8Bu, 0x768967E5EEC84E2Eu,822 0x76BFC1DF6A7A61BAu, 0x76F3D92BA28C7D14u, 0x7728CF768B2F9C59u, 0x775F03542DFB8370u, 0x779362149CBD3226u,823 0x77C83A99C3EC7EAFu, 0x77FE494034E79E5Bu, 0x7832EDC82110C2F9u, 0x7867A93A2954F3B7u, 0x789D9388B3AA30A5u,824 0x78D27C35704A5E67u, 0x79071B42CC5CF601u, 0x793CE2137F743381u, 0x79720D4C2FA8A030u, 0x79A6909F3B92C83Du,825 0x79DC34C70A777A4Cu, 0x7A11A0FC668AAC6Fu, 0x7A46093B802D578Bu, 0x7A7B8B8A6038AD6Eu, 0x7AB137367C236C65u,826 0x7AE585041B2C477Eu, 0x7B1AE64521F7595Eu, 0x7B50CFEB353A97DAu, 0x7B8503E602893DD1u, 0x7BBA44DF832B8D45u,827 0x7BF06B0BB1FB384Bu, 0x7C2485CE9E7A065Eu, 0x7C59A742461887F6u, 0x7C9008896BCF54F9u, 0x7CC40AABC6C32A38u,828 0x7CF90D56B873F4C6u, 0x7D2F50AC6690F1F8u, 0x7D63926BC01A973Bu, 0x7D987706B0213D09u, 0x7DCE94C85C298C4Cu,829 0x7E031CFD3999F7AFu, 0x7E37E43C8800759Bu, 0x7E6DDD4BAA009302u, 0x7EA2AA4F4A405BE1u, 0x7ED754E31CD072D9u,830 0x7F0D2A1BE4048F90u, 0x7F423A516E82D9BAu, 0x7F76C8E5CA239028u, 0x7FAC7B1F3CAC7433u, 0x7FE1CCF385EBC89Fu,831 0x7FEFFFFFFFFFFFFFu};832};833 834template <class _Floating>835[[nodiscard]] _LIBCPP_HIDE_FROM_ABI836to_chars_result _Floating_to_chars_general_precision(837 char* _First, char* const _Last, const _Floating _Value, int _Precision) noexcept {838 839 using _Traits = _Floating_type_traits<_Floating>;840 using _Uint_type = typename _Traits::_Uint_type;841 842 const _Uint_type _Uint_value = std::bit_cast<_Uint_type>(_Value);843 844 if (_Uint_value == 0) { // zero detected; write "0" and return; _Precision is irrelevant due to zero-trimming845 if (_First == _Last) {846 return {_Last, errc::value_too_large};847 }848 849 *_First++ = '0';850 851 return {_First, errc{}};852 }853 854 // C11 7.21.6.1 "The fprintf function"/5:855 // "A negative precision argument is taken as if the precision were omitted."856 // /8: "g,G [...] Let P equal the precision if nonzero, 6 if the precision is omitted,857 // or 1 if the precision is zero."858 859 // Performance note: It's possible to rewrite this for branchless codegen,860 // but profiling will be necessary to determine whether that's faster.861 if (_Precision < 0) {862 _Precision = 6;863 } else if (_Precision == 0) {864 _Precision = 1;865 } else if (_Precision < 1'000'000) {866 // _Precision is ok.867 } else {868 // Avoid integer overflow.869 // Due to general notation's zero-trimming behavior, we can simply clamp _Precision.870 // This is further clamped below.871 _Precision = 1'000'000;872 }873 874 // _Precision is now the Standard's P.875 876 // /8: "Then, if a conversion with style E would have an exponent of X:877 // - if P > X >= -4, the conversion is with style f (or F) and precision P - (X + 1).878 // - otherwise, the conversion is with style e (or E) and precision P - 1."879 880 // /8: "Finally, [...] any trailing zeros are removed from the fractional portion of the result881 // and the decimal-point character is removed if there is no fractional portion remaining."882 883 using _Tables = _General_precision_tables<_Floating>;884 885 const _Uint_type* _Table_begin;886 const _Uint_type* _Table_end;887 888 if (_Precision <= _Tables::_Max_special_P) {889 _Table_begin = _Tables::_Special_X_table + (_Precision - 1) * (_Precision + 10) / 2;890 _Table_end = _Table_begin + _Precision + 5;891 } else {892 _Table_begin = _Tables::_Ordinary_X_table;893 _Table_end = _Table_begin + std::min(_Precision, _Tables::_Max_P) + 5;894 }895 896 // Profiling indicates that linear search is faster than binary search for small tables.897 // Performance note: lambda captures may have a small performance cost.898 const _Uint_type* const _Table_lower_bound = [=] {899 if constexpr (!_IsSame<_Floating, float>::value) {900 if (_Precision > 155) { // threshold determined via profiling901 return std::lower_bound(_Table_begin, _Table_end, _Uint_value, less{});902 }903 }904 905 return std::find_if(_Table_begin, _Table_end, [=](const _Uint_type _Elem) { return _Uint_value <= _Elem; });906 }();907 908 const ptrdiff_t _Table_index = _Table_lower_bound - _Table_begin;909 const int _Scientific_exponent_X = static_cast<int>(_Table_index - 5);910 const bool _Use_fixed_notation = _Precision > _Scientific_exponent_X && _Scientific_exponent_X >= -4;911 912 // Performance note: it might (or might not) be faster to modify Ryu Printf to perform zero-trimming.913 // Such modifications would involve a fairly complicated state machine (notably, both '0' and '9' digits would914 // need to be buffered, due to rounding), and that would have performance costs due to increased branching.915 // Here, we're using a simpler approach: writing into a local buffer, manually zero-trimming, and then copying into916 // the output range. The necessary buffer size is reasonably small, the zero-trimming logic is simple and fast,917 // and the final copying is also fast.918 919 constexpr int _Max_output_length =920 _IsSame<_Floating, float>::value ? 117 : 773; // cases: 0x1.fffffep-126f and 0x1.fffffffffffffp-1022921 constexpr int _Max_fixed_precision =922 _IsSame<_Floating, float>::value ? 37 : 66; // cases: 0x1.fffffep-14f and 0x1.fffffffffffffp-14923 constexpr int _Max_scientific_precision =924 _IsSame<_Floating, float>::value ? 111 : 766; // cases: 0x1.fffffep-126f and 0x1.fffffffffffffp-1022925 926 // Note that _Max_output_length is determined by scientific notation and is more than enough for fixed notation.927 // 0x1.fffffep+127f is 39 digits, plus 1 for '.', plus _Max_fixed_precision for '0' digits, equals 77.928 // 0x1.fffffffffffffp+1023 is 309 digits, plus 1 for '.', plus _Max_fixed_precision for '0' digits, equals 376.929 930 char _Buffer[_Max_output_length];931 const char* const _Significand_first = _Buffer; // e.g. "1.234"932 const char* _Significand_last = nullptr;933 const char* _Exponent_first = nullptr; // e.g. "e-05"934 const char* _Exponent_last = nullptr;935 int _Effective_precision; // number of digits printed after the decimal point, before trimming936 937 // Write into the local buffer.938 // Clamping _Effective_precision allows _Buffer to be as small as possible, and increases efficiency.939 if (_Use_fixed_notation) {940 _Effective_precision = std::min(_Precision - (_Scientific_exponent_X + 1), _Max_fixed_precision);941 const to_chars_result _Buf_result =942 _Floating_to_chars_fixed_precision(_Buffer, std::end(_Buffer), _Value, _Effective_precision);943 _LIBCPP_ASSERT_INTERNAL(_Buf_result.ec == errc{}, "");944 _Significand_last = _Buf_result.ptr;945 } else {946 _Effective_precision = std::min(_Precision - 1, _Max_scientific_precision);947 const to_chars_result _Buf_result =948 _Floating_to_chars_scientific_precision(_Buffer, std::end(_Buffer), _Value, _Effective_precision);949 _LIBCPP_ASSERT_INTERNAL(_Buf_result.ec == errc{}, "");950 _Significand_last = std::find(_Buffer, _Buf_result.ptr, 'e');951 _Exponent_first = _Significand_last;952 _Exponent_last = _Buf_result.ptr;953 }954 955 // If we printed a decimal point followed by digits, perform zero-trimming.956 if (_Effective_precision > 0) {957 while (_Significand_last[-1] == '0') { // will stop at '.' or a nonzero digit958 --_Significand_last;959 }960 961 if (_Significand_last[-1] == '.') {962 --_Significand_last;963 }964 }965 966 // Copy the significand to the output range.967 const ptrdiff_t _Significand_distance = _Significand_last - _Significand_first;968 if (_Last - _First < _Significand_distance) {969 return {_Last, errc::value_too_large};970 }971 std::memcpy(_First, _Significand_first, static_cast<size_t>(_Significand_distance));972 _First += _Significand_distance;973 974 // Copy the exponent to the output range.975 if (!_Use_fixed_notation) {976 const ptrdiff_t _Exponent_distance = _Exponent_last - _Exponent_first;977 if (_Last - _First < _Exponent_distance) {978 return {_Last, errc::value_too_large};979 }980 std::memcpy(_First, _Exponent_first, static_cast<size_t>(_Exponent_distance));981 _First += _Exponent_distance;982 }983 984 return {_First, errc{}};985}986 987enum class _Floating_to_chars_overload { _Plain, _Format_only, _Format_precision };988 989template <_Floating_to_chars_overload _Overload, class _Floating>990[[nodiscard]] _LIBCPP_HIDE_FROM_ABI991to_chars_result _Floating_to_chars(992 char* _First, char* const _Last, _Floating _Value, const chars_format _Fmt, const int _Precision) noexcept {993 994 if constexpr (_Overload == _Floating_to_chars_overload::_Plain) {995 _LIBCPP_ASSERT_INTERNAL(_Fmt == chars_format{}, ""); // plain overload must pass chars_format{} internally996 } else {997 _LIBCPP_ASSERT_ARGUMENT_WITHIN_DOMAIN(_Fmt == chars_format::general || _Fmt == chars_format::scientific998 || _Fmt == chars_format::fixed || _Fmt == chars_format::hex,999 "invalid format in to_chars()");1000 }1001 1002 using _Traits = _Floating_type_traits<_Floating>;1003 using _Uint_type = typename _Traits::_Uint_type;1004 1005 _Uint_type _Uint_value = std::bit_cast<_Uint_type>(_Value);1006 1007 const bool _Was_negative = (_Uint_value & _Traits::_Shifted_sign_mask) != 0;1008 1009 if (_Was_negative) { // sign bit detected; write minus sign and clear sign bit1010 if (_First == _Last) {1011 return {_Last, errc::value_too_large};1012 }1013 1014 *_First++ = '-';1015 1016 _Uint_value &= ~_Traits::_Shifted_sign_mask;1017 _Value = std::bit_cast<_Floating>(_Uint_value);1018 }1019 1020 if ((_Uint_value & _Traits::_Shifted_exponent_mask) == _Traits::_Shifted_exponent_mask) {1021 // inf/nan detected; write appropriate string and return1022 const char* _Str;1023 size_t _Len;1024 1025 const _Uint_type _Mantissa = _Uint_value & _Traits::_Denormal_mantissa_mask;1026 1027 if (_Mantissa == 0) {1028 _Str = "inf";1029 _Len = 3;1030 } else if (_Was_negative && _Mantissa == _Traits::_Special_nan_mantissa_mask) {1031 // When a NaN value has the sign bit set, the quiet bit set, and all other mantissa bits cleared,1032 // the UCRT interprets it to mean "indeterminate", and indicates this by printing "-nan(ind)".1033 _Str = "nan(ind)";1034 _Len = 8;1035 } else if ((_Mantissa & _Traits::_Special_nan_mantissa_mask) != 0) {1036 _Str = "nan";1037 _Len = 3;1038 } else {1039 _Str = "nan(snan)";1040 _Len = 9;1041 }1042 1043 if (_Last - _First < static_cast<ptrdiff_t>(_Len)) {1044 return {_Last, errc::value_too_large};1045 }1046 1047 std::memcpy(_First, _Str, _Len);1048 1049 return {_First + _Len, errc{}};1050 }1051 1052 if constexpr (_Overload == _Floating_to_chars_overload::_Plain) {1053 return _Floating_to_chars_ryu(_First, _Last, _Value, chars_format{});1054 } else if constexpr (_Overload == _Floating_to_chars_overload::_Format_only) {1055 if (_Fmt == chars_format::hex) {1056 return _Floating_to_chars_hex_shortest(_First, _Last, _Value);1057 }1058 1059 return _Floating_to_chars_ryu(_First, _Last, _Value, _Fmt);1060 } else if constexpr (_Overload == _Floating_to_chars_overload::_Format_precision) {1061 switch (_Fmt) {1062 case chars_format::scientific:1063 return _Floating_to_chars_scientific_precision(_First, _Last, _Value, _Precision);1064 case chars_format::fixed:1065 return _Floating_to_chars_fixed_precision(_First, _Last, _Value, _Precision);1066 case chars_format::general:1067 return _Floating_to_chars_general_precision(_First, _Last, _Value, _Precision);1068 case chars_format::hex:1069 default: // avoid MSVC warning C4715: not all control paths return a value1070 return _Floating_to_chars_hex_precision(_First, _Last, _Value, _Precision);1071 }1072 }1073}1074 1075// clang-format on1076 1077_LIBCPP_END_NAMESPACE_STD1078 1079#endif // _LIBCPP_SRC_INCLUDE_TO_CHARS_FLOATING_POINT_H1080