log.cpp | Explore in Territory

//===-- Double-precision log(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/log.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/common.h"
#include "src/__support/integer_literals.h"
#include "src/__support/macros/config.h"
#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY

#include "common_constants.h"
#include "log_range_reduction.h"

namespace LIBC_NAMESPACE_DECL {

// 128-bit precision dyadic floating point numbers.
Float128;

operator""_u128;

namespace {

// A simple upper bound for the error of e_x * log(2) - log(r).
constexpr double HI_ERR = …;

// Extra errors from P is from using x^2 to reduce evaluation latency.
constexpr double P_ERR = …;

// log(2) with 128-bit precision generated by SageMath with:
// def format_hex(value):
//     l = hex(value)[2:]
//     n = 8
//     x = [l[i:i + n] for i in range(0, len(l), n)]
//     return "0x" + "'".join(x) + "_u128"
//  (s, m, e) = RealField(128)(2).log().sign_mantissa_exponent();
//  print(format_hex(m));
constexpr Float128 LOG_2(Sign::POS, /*exponent=*/-128, /*mantissa=*/
                         0xb17217f7'd1cf79ab'c9e3b398'03f2f6af_u128);

alignas(64) constexpr LogRR LOG_TABLE = …;
// > P = fpminimax((log(1 + x) - x)/x^2, 2, [|1, 128...|],
//                 [-0x1.0002143p-29 , 0x1p-29]);
// > P;
// > dirtyinfnorm(log(1 + x)/x - x*P, [-0x1.0002143p-29 , 0x1p-29]);
// 0x1.99a3...p-121
constexpr Float128 BIG_COEFFS[3]{ … };

// Reuse the output of the fast pass range reduction.
// -2^-8 <= m_x < 2^-7
double log_accurate(int e_x, int index, double m_x) { … }

} // namespace

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

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