charconv/detail/compute_float80.hpp

0001 // Copyright 2023 Matt Borland
0002 // Distributed under the Boost Software License, Version 1.0.
0003 // https://www.boost.org/LICENSE_1_0.txt
0004
0005 #ifndef BOOST_CHARCONV_DETAIL_COMPUTE_FLOAT80_HPP
0006 #define BOOST_CHARCONV_DETAIL_COMPUTE_FLOAT80_HPP
0007
0008 #include <boost/charconv/detail/config.hpp>
0009 #include <boost/charconv/detail/emulated128.hpp>
0010 #include <boost/charconv/detail/bit_layouts.hpp>
0011 #include <system_error>
0012 #include <type_traits>
0013 #include <limits>
0014 #include <cstdint>
0015 #include <cmath>
0016 #include <climits>
0017 #include <cfloat>
0018
0019 #ifdef BOOST_CHARCONV_DEBUG_FLOAT128
0020 #include <iostream>
0021 #include <iomanip>
0022 #include <boost/charconv/detail/to_chars_integer_impl.hpp>
0023 #endif
0024
0025 namespace boost { namespace charconv { namespace detail {
0026
0027 #if BOOST_CHARCONV_LDBL_BITS > 64
0028
0029 static constexpr long double powers_of_ten_ld[] = {
0030     1e0L,  1e1L,  1e2L,  1e3L,  1e4L,  1e5L,  1e6L,
0031     1e7L,  1e8L,  1e9L,  1e10L, 1e11L, 1e12L, 1e13L,
0032     1e14L, 1e15L, 1e16L, 1e17L, 1e18L, 1e19L, 1e20L,
0033     1e21L, 1e22L, 1e23L, 1e24L, 1e25L, 1e26L, 1e27L,
0034     1e28L, 1e29L, 1e30L, 1e31L, 1e32L, 1e33L, 1e34L,
0035     1e35L, 1e36L, 1e37L, 1e38L, 1e39L, 1e40L, 1e41L,
0036     1e42L, 1e43L, 1e44L, 1e45L, 1e46L, 1e47L, 1e48L,
0037     1e49L, 1e50L, 1e51L, 1e52L, 1e53L, 1e54L, 1e55L
0038 };
0039
0040 template <typename ResultType, typename Unsigned_Integer, typename ArrayPtr>
0041 inline ResultType fast_path(std::int64_t q, Unsigned_Integer w, bool negative, ArrayPtr table) noexcept
0042 {
0043     // The general idea is as follows.
0044     // if 0 <= s <= 2^64 and if 10^0 <= p <= 10^27
0045     // Both s and p can be represented exactly
0046     // because of this s*p and s/p will produce
0047     // correctly rounded values
0048
0049     auto ld = static_cast<ResultType>(w);
0050
0051     if (q < 0)
0052     {
0053         ld /= table[-q];
0054     }
0055     else
0056     {
0057         ld *= table[q];
0058     }
0059
0060     if (negative)
0061     {
0062         ld = -ld;
0063     }
0064
0065     return ld;
0066 }
0067
0068 template <typename ResultType, typename Unsigned_Integer>
0069 inline ResultType compute_float80(std::int64_t q, Unsigned_Integer w, bool negative, std::errc& success) noexcept
0070 {
0071     // GLIBC uses 2^-16444 but MPFR uses 2^-16445 as the smallest subnormal value for 80 bit
0072     // 39 is the max number of digits in an uint128_t
0073     static constexpr auto smallest_power = -4951 - 39;
0074     static constexpr auto largest_power = 4932;
0075
0076     // We start with a fast path
0077     // It is an extension of what was described in Clinger WD.
0078     // How to read floating point numbers accurately.
0079     // ACM SIGPLAN Notices. 1990
0080     // https://dl.acm.org/doi/pdf/10.1145/93542.93557
0081     static constexpr auto clinger_max_exp = BOOST_CHARCONV_LDBL_BITS == 80 ? 27 : 48;   // NOLINT : Only changes by platform
0082     static constexpr auto clinger_min_exp = BOOST_CHARCONV_LDBL_BITS == 80 ? -34 : -55; // NOLINT
0083
0084     if (clinger_min_exp <= q && q <= clinger_max_exp && w <= static_cast<Unsigned_Integer>(1) << 113)
0085     {
0086         success = std::errc();
0087         return fast_path<ResultType>(q, w, negative, powers_of_ten_ld);
0088     }
0089
0090     if (w == 0)
0091     {
0092         success = std::errc();
0093         return negative ? -0.0L : 0.0L;
0094     }
0095     else if (q > largest_power)
0096     {
0097         success = std::errc::result_out_of_range;
0098         return negative ? -HUGE_VALL : HUGE_VALL;
0099     }
0100     else if (q < smallest_power)
0101     {
0102         success = std::errc::result_out_of_range;
0103         return negative ? -0.0L : 0.0L;
0104     }
0105
0106     success = std::errc::not_supported;
0107     return 0;
0108 }
0109
0110 #endif // BOOST_CHARCONV_LDBL_BITS > 64
0111
0112 }}} // Namespaces
0113
0114 #endif // BOOST_CHARCONV_DETAIL_COMPUTE_FLOAT80_HPP