exp10.cpp | Explore in Territory

//===-- Double-precision 10^x function ------------------------------------===//
//
// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
// See https://llvm.org/LICENSE.txt for license information.
// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
//
//===----------------------------------------------------------------------===//

#include "src/math/exp10.h"
#include "common_constants.h" // Lookup tables EXP2_MID1 and EXP_M2.
#include "explogxf.h"         // ziv_test_denorm.
#include "src/__support/CPP/bit.h"
#include "src/__support/CPP/optional.h"
#include "src/__support/FPUtil/FEnvImpl.h"
#include "src/__support/FPUtil/FPBits.h"
#include "src/__support/FPUtil/PolyEval.h"
#include "src/__support/FPUtil/double_double.h"
#include "src/__support/FPUtil/dyadic_float.h"
#include "src/__support/FPUtil/multiply_add.h"
#include "src/__support/FPUtil/nearest_integer.h"
#include "src/__support/FPUtil/rounding_mode.h"
#include "src/__support/FPUtil/triple_double.h"
#include "src/__support/common.h"
#include "src/__support/integer_literals.h"
#include "src/__support/macros/config.h"
#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY

#include <errno.h>

namespace LIBC_NAMESPACE_DECL {

DoubleDouble;
TripleDouble;
Float128;

operator""_u128;

// log2(10)
constexpr double LOG2_10 = …;

// -2^-12 * log10(2)
// > a = -2^-12 * log10(2);
// > b = round(a, 32, RN);
// > c = round(a - b, 32, RN);
// > d = round(a - b - c, D, RN);
// Errors < 1.5 * 2^-144
constexpr double MLOG10_2_EXP2_M12_HI = …;
constexpr double MLOG10_2_EXP2_M12_MID = …;
constexpr double MLOG10_2_EXP2_M12_MID_32 = …;
constexpr double MLOG10_2_EXP2_M12_LO = …;

// Error bounds:
// Errors when using double precision.
constexpr double ERR_D = …;

// Errors when using double-double precision.
constexpr double ERR_DD = …;

namespace {

// Polynomial approximations with double precision.  Generated by Sollya with:
// > P = fpminimax((10^x - 1)/x, 3, [|D...|], [-2^-14, 2^-14]);
// > P;
// Error bounds:
//   | output - (10^dx - 1) / dx | < 2^-52.
LIBC_INLINE double poly_approx_d(double dx) { … }

// Polynomial approximation with double-double precision.  Generated by Solya
// with:
// > P = fpminimax((10^x - 1)/x, 5, [|DD...|], [-2^-14, 2^-14]);
// Error bounds:
//   | output - 10^(dx) | < 2^-101
DoubleDouble poly_approx_dd(const DoubleDouble &dx) { … }

// Polynomial approximation with 128-bit precision:
// Return exp(dx) ~ 1 + a0 * dx + a1 * dx^2 + ... + a6 * dx^7
// For |dx| < 2^-14:
//   | output - 10^dx | < 1.5 * 2^-124.
Float128 poly_approx_f128(const Float128 &dx) { … }

// Compute 10^(x) using 128-bit precision.
// TODO(lntue): investigate triple-double precision implementation for this
// step.
Float128 exp10_f128(double x, double kd, int idx1, int idx2) { … }

// Compute 10^x with double-double precision.
DoubleDouble exp10_double_double(double x, double kd,
                                 const DoubleDouble &exp_mid) { … }

// When output is denormal.
double exp10_denorm(double x) { … }

// Check for exceptional cases when:
//  * log10(1 - 2^-54) < x < log10(1 + 2^-53)
//  * x >= log10(2^1024)
//  * x <= log10(2^-1022)
//  * x is inf or nan
double set_exceptional(double x) { … }

} // namespace

LLVM_LIBC_FUNCTION(double, exp10, (double x)) { … }

} // namespace LIBC_NAMESPACE_DECL
llvm/libc/src/math/generic/exp10.cpp
