summaryrefslogtreecommitdiff
path: root/sysdeps/ieee754/dbl-64/lgamma_product.c
diff options
context:
space:
mode:
Diffstat (limited to 'sysdeps/ieee754/dbl-64/lgamma_product.c')
-rw-r--r--sysdeps/ieee754/dbl-64/lgamma_product.c34
1 files changed, 2 insertions, 32 deletions
diff --git a/sysdeps/ieee754/dbl-64/lgamma_product.c b/sysdeps/ieee754/dbl-64/lgamma_product.c
index d956575bc7..b37ecee73f 100644
--- a/sysdeps/ieee754/dbl-64/lgamma_product.c
+++ b/sysdeps/ieee754/dbl-64/lgamma_product.c
@@ -1,5 +1,5 @@
/* Compute a product of 1 + (T/X), 1 + (T/(X+1)), ....
- Copyright (C) 2015-2016 Free Software Foundation, Inc.
+ Copyright (C) 2015-2018 Free Software Foundation, Inc.
This file is part of the GNU C Library.
The GNU C Library is free software; you can redistribute it and/or
@@ -18,37 +18,7 @@
#include <math.h>
#include <math_private.h>
-#include <float.h>
-
-/* Calculate X * Y exactly and store the result in *HI + *LO. It is
- given that the values are small enough that no overflow occurs and
- large enough (or zero) that no underflow occurs. */
-
-static void
-mul_split (double *hi, double *lo, double x, double y)
-{
-#ifdef __FP_FAST_FMA
- /* Fast built-in fused multiply-add. */
- *hi = x * y;
- *lo = __builtin_fma (x, y, -*hi);
-#elif defined FP_FAST_FMA
- /* Fast library fused multiply-add, compiler before GCC 4.6. */
- *hi = x * y;
- *lo = __fma (x, y, -*hi);
-#else
- /* Apply Dekker's algorithm. */
- *hi = x * y;
-# define C ((1 << (DBL_MANT_DIG + 1) / 2) + 1)
- double x1 = x * C;
- double y1 = y * C;
-# undef C
- x1 = (x - x1) + x1;
- y1 = (y - y1) + y1;
- double x2 = x - x1;
- double y2 = y - y1;
- *lo = (((x1 * y1 - *hi) + x1 * y2) + x2 * y1) + x2 * y2;
-#endif
-}
+#include <mul_split.h>
/* Compute the product of 1 + (T / (X + X_EPS)), 1 + (T / (X + X_EPS +
1)), ..., 1 + (T / (X + X_EPS + N - 1)), minus 1. X is such that