[libc-commits] [libc] [llvm] [libc][math] Impl bfloat16 lgamma function. (PR #199312)

[via libc-commits]

================
@@ -0,0 +1,195 @@
+//===-- Implementation of lgammabf16 ----------------------------*- 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_LGAMMABF16_H
+#define LLVM_LIBC_SRC___SUPPORT_MATH_LGAMMABF16_H
+
+#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/bfloat16.h"
+#include "src/__support/FPUtil/cast.h"
+#include "src/__support/macros/config.h"
+#include "src/__support/macros/optimization.h"
+
+namespace LIBC_NAMESPACE_DECL {
+namespace math {
+
+// Evaluate degree-4 polynomial using Horner's method:
+// p(t) = c0 + t*(c1 + t*(c2 + t*(c3 + t*c4)))
+// coeffs = {c4, c3, c2, c1, c0}
+static inline float poly4(float t, const float c[5]) {
+  return c[4] + t * (c[3] + t * (c[2] + t * (c[1] + t * c[0])));
+}
+
+// Compute natural log of positive float x using range reduction.
+// x = 2^e * m, m in [1, 2)
+// ln(x) = e * ln(2) + ln(m)
+// ln(m) approximated by degree-5 poly centered at 1.5
+static inline float lgamma_logf(float x) {
+  // Coefficients for ln(m) on [1,2), centered at 1.5
+  // Generated by numpy.polyfit, max_err = 1.71e-05 (well within bf16 precision)
+  constexpr float LN_COEFFS[6] = {
+      0x1.ef26300000000p-6f,  // c5
----------------
lntue wrote:

adjust the precision of all the hexadecimal constants in this file.  You overprinted them for float.

https://github.com/llvm/llvm-project/pull/199312