llvm/libc/src/math/generic/pow.cpp

//===-- Double-precision x^y 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/pow.h"
#include "common_constants.h" // Lookup tables EXP_M1 and EXP_M2.
#include "hdr/errno_macros.h"
#include "hdr/fenv_macros.h"
#include "src/__support/CPP/bit.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/multiply_add.h"
#include "src/__support/FPUtil/nearest_integer.h"
#include "src/__support/FPUtil/sqrt.h" // Speedup for pow(x, 1/2) = sqrt(x)
#include "src/__support/common.h"
#include "src/__support/macros/config.h"
#include "src/__support/macros/optimization.h" // LIBC_UNLIKELY

namespace LIBC_NAMESPACE_DECL {

DoubleDouble;

namespace {

// Constants for log2(x) range reduction, generated by Sollya with:
// > for i from 0 to 127 do {
//     r = 2^-8 * ceil( 2^8 * (1 - 2^(-8)) / (1 + i*2^-7) );
//     b = nearestint(log2(r) * 2^41) * 2^-41;
//     c = round(log2(r) - b, D, RN);
//     print("{", -c, ",", -b, "},");
//   };
// This is the same as -log2(RD[i]), with the least significant bits of the
// high part set to be 2^-41, so that the sum of high parts + e_x is exact in
// double precision.
// We also replace the first and the last ones to be 0.
constexpr DoubleDouble LOG2_R_DD[128] = …;

bool is_odd_integer(double x) { … }

bool is_integer(double x) { … }

} // namespace

LLVM_LIBC_FUNCTION(double, pow, (double x, double y)) { … }

} // namespace LIBC_NAMESPACE_DECL