1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
|
//===-- Collection of utils for acoshf --------------------------*- C++ -*-===//
//
// 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
//
//===----------------------------------------------------------------------===//
#ifndef LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_UTILS_H
#define LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_UTILS_H
#include "acosh_float_constants.h"
#include "src/__support/FPUtil/FPBits.h"
#include "src/__support/FPUtil/PolyEval.h"
#include "src/__support/FPUtil/multiply_add.h"
namespace LIBC_NAMESPACE_DECL {
namespace acoshf_internal {
// x should be positive, normal finite value
LIBC_INLINE static double log_eval(double x) {
// For x = 2^ex * (1 + mx)
// log(x) = ex * log(2) + log(1 + mx)
using FPB = fputil::FPBits<double>;
FPB bs(x);
double ex = static_cast<double>(bs.get_exponent());
// p1 is the leading 7 bits of mx, i.e.
// p1 * 2^(-7) <= m_x < (p1 + 1) * 2^(-7).
int p1 = static_cast<int>(bs.get_mantissa() >> (FPB::FRACTION_LEN - 7));
// Set bs to (1 + (mx - p1*2^(-7))
bs.set_uintval(bs.uintval() & (FPB::FRACTION_MASK >> 7));
bs.set_biased_exponent(FPB::EXP_BIAS);
// dx = (mx - p1*2^(-7)) / (1 + p1*2^(-7)).
double dx = (bs.get_val() - 1.0) * ONE_OVER_F[p1];
// Minimax polynomial of log(1 + dx) generated by Sollya with:
// > P = fpminimax(log(1 + x)/x, 6, [|D...|], [0, 2^-7]);
const double COEFFS[6] = {-0x1.ffffffffffffcp-2, 0x1.5555555552ddep-2,
-0x1.ffffffefe562dp-3, 0x1.9999817d3a50fp-3,
-0x1.554317b3f67a5p-3, 0x1.1dc5c45e09c18p-3};
double dx2 = dx * dx;
double c1 = fputil::multiply_add(dx, COEFFS[1], COEFFS[0]);
double c2 = fputil::multiply_add(dx, COEFFS[3], COEFFS[2]);
double c3 = fputil::multiply_add(dx, COEFFS[5], COEFFS[4]);
double p = fputil::polyeval(dx2, dx, c1, c2, c3);
double result =
fputil::multiply_add(ex, /*log(2)*/ 0x1.62e42fefa39efp-1, LOG_F[p1] + p);
return result;
}
} // namespace acoshf_internal
} // namespace LIBC_NAMESPACE_DECL
#endif // LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_UTILS_H
|
