aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--ChangeLog13
-rw-r--r--NEWS2
-rw-r--r--math/Makefile2
-rw-r--r--sysdeps/i386/fpu/e_log_data.c1
-rw-r--r--sysdeps/ia64/fpu/e_log_data.c1
-rw-r--r--sysdeps/ieee754/dbl-64/e_log.c257
-rw-r--r--sysdeps/ieee754/dbl-64/e_log_data.c347
-rw-r--r--sysdeps/ieee754/dbl-64/math_config.h16
-rw-r--r--sysdeps/ieee754/dbl-64/ulog.h93
-rw-r--r--sysdeps/ieee754/dbl-64/ulog.tbl3326
-rw-r--r--sysdeps/m68k/m680x0/fpu/e_log_data.c1
11 files changed, 492 insertions, 3567 deletions
diff --git a/ChangeLog b/ChangeLog
index 0df8aa58ae..8c990ad7b7 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,16 @@
+2018-09-12 Szabolcs Nagy <szabolcs.nagy@arm.com>
+
+ * NEWS: Mention log improvement.
+ * math/Makefile (type-double-routines): Add e_log_data.
+ * sysdeps/i386/fpu/e_log_data.c: New file.
+ * sysdeps/ia64/fpu/e_log_data.c: New file.
+ * sysdeps/ieee754/dbl-64/e_log.c: Rewrite.
+ * sysdeps/ieee754/dbl-64/e_log_data.c: New file.
+ * sysdeps/ieee754/dbl-64/math_config.h (__log_data): Add.
+ * sysdeps/ieee754/dbl-64/ulog.h: Remove.
+ * sysdeps/ieee754/dbl-64/ulog.tbl: Remove.
+ * sysdeps/m68k/m680x0/fpu/e_log_data.c: New file.
+
2018-09-12 H.J. Lu <hongjiu.lu@intel.com>
Xuepeng Guo <xuepeng.guo@intel.com>
diff --git a/NEWS b/NEWS
index 085325ab87..79bee8ee6b 100644
--- a/NEWS
+++ b/NEWS
@@ -16,7 +16,7 @@ Major new features:
to set the install root if you wish to install into a non-default
configured location.
-* Optimized generic exp, exp2, sinf, cosf, sincosf and tanf.
+* Optimized generic exp, exp2, log, sinf, cosf, sincosf and tanf.
* The reallocarray function is now declared under _DEFAULT_SOURCE, not just
for _GNU_SOURCE, to match BSD environments.
diff --git a/math/Makefile b/math/Makefile
index f1ba1a0c36..8bfbebc4d0 100644
--- a/math/Makefile
+++ b/math/Makefile
@@ -127,7 +127,7 @@ type-ldouble-yes := ldouble
type-double-suffix :=
type-double-routines := branred doasin dosincos mpa mpatan2 \
k_rem_pio2 mpatan mpsqrt mptan sincos32 \
- sincostab math_err e_exp_data
+ sincostab math_err e_exp_data e_log_data
# float support
type-float-suffix := f
diff --git a/sysdeps/i386/fpu/e_log_data.c b/sysdeps/i386/fpu/e_log_data.c
new file mode 100644
index 0000000000..1cc8931700
--- /dev/null
+++ b/sysdeps/i386/fpu/e_log_data.c
@@ -0,0 +1 @@
+/* Not needed. */
diff --git a/sysdeps/ia64/fpu/e_log_data.c b/sysdeps/ia64/fpu/e_log_data.c
new file mode 100644
index 0000000000..1cc8931700
--- /dev/null
+++ b/sysdeps/ia64/fpu/e_log_data.c
@@ -0,0 +1 @@
+/* Not needed. */
diff --git a/sysdeps/ieee754/dbl-64/e_log.c b/sysdeps/ieee754/dbl-64/e_log.c
index 2483dd8551..a56b714fb7 100644
--- a/sysdeps/ieee754/dbl-64/e_log.c
+++ b/sysdeps/ieee754/dbl-64/e_log.c
@@ -1,167 +1,132 @@
-/*
- * IBM Accurate Mathematical Library
- * written by International Business Machines Corp.
- * Copyright (C) 2001-2018 Free Software Foundation, Inc.
- *
- * This program is free software; you can redistribute it and/or modify
- * it under the terms of the GNU Lesser General Public License as published by
- * the Free Software Foundation; either version 2.1 of the License, or
- * (at your option) any later version.
- *
- * This program is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
- * GNU Lesser General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public License
- * along with this program; if not, see <http://www.gnu.org/licenses/>.
- */
-/*********************************************************************/
-/* */
-/* MODULE_NAME:ulog.c */
-/* */
-/* FUNCTION:ulog */
-/* */
-/* FILES NEEDED: dla.h endian.h mpa.h mydefs.h ulog.h */
-/* ulog.tbl */
-/* */
-/* An ultimate log routine. Given an IEEE double machine number x */
-/* it computes the rounded (to nearest) value of log(x). */
-/* Assumption: Machine arithmetic operations are performed in */
-/* round to nearest mode of IEEE 754 standard. */
-/* */
-/*********************************************************************/
+/* Double-precision log(x) function.
+ Copyright (C) 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
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ The GNU C Library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, see
+ <http://www.gnu.org/licenses/>. */
-#include "endian.h"
-#include <dla.h>
-#include "mpa.h"
-#include "MathLib.h"
#include <math.h>
-#include <math_private.h>
+#include <stdint.h>
+#include "math_config.h"
+
+#define T __log_data.tab
+#define T2 __log_data.tab2
+#define B __log_data.poly1
+#define A __log_data.poly
+#define Ln2hi __log_data.ln2hi
+#define Ln2lo __log_data.ln2lo
+#define N (1 << LOG_TABLE_BITS)
+#define OFF 0x3fe6000000000000
+
+/* Top 16 bits of a double. */
+static inline uint32_t
+top16 (double x)
+{
+ return asuint64 (x) >> 48;
+}
#ifndef SECTION
# define SECTION
#endif
-/*********************************************************************/
-/* An ultimate log routine. Given an IEEE double machine number x */
-/* it computes the rounded (to nearest) value of log(x). */
-/*********************************************************************/
double
SECTION
__ieee754_log (double x)
{
- int i, j, n, ux, dx;
- double dbl_n, u, p0, q, r0, w, nln2a, luai, lubi, lvaj, lvbj,
- sij, ssij, ttij, A, B, B0, polI, polII, t8, a, aa, b, bb, c;
-#ifndef DLA_FMS
- double t1, t2, t3, t4, t5;
-#endif
- number num;
-
-#include "ulog.tbl"
-#include "ulog.h"
-
- /* Treating special values of x ( x<=0, x=INF, x=NaN etc.). */
-
- num.d = x;
- ux = num.i[HIGH_HALF];
- dx = num.i[LOW_HALF];
- n = 0;
- if (__glibc_unlikely (ux < 0x00100000))
+ /* double_t for better performance on targets with FLT_EVAL_METHOD==2. */
+ double_t w, z, r, r2, r3, y, invc, logc, kd, hi, lo;
+ uint64_t ix, iz, tmp;
+ uint32_t top;
+ int k, i;
+
+ ix = asuint64 (x);
+ top = top16 (x);
+
+#define LO asuint64 (1.0 - 0x1p-4)
+#define HI asuint64 (1.0 + 0x1.09p-4)
+ if (__glibc_unlikely (ix - LO < HI - LO))
{
- if (__glibc_unlikely (((ux & 0x7fffffff) | dx) == 0))
- return MHALF / 0.0; /* return -INF */
- if (__glibc_unlikely (ux < 0))
- return (x - x) / 0.0; /* return NaN */
- n -= 54;
- x *= two54.d; /* scale x */
- num.d = x;
+ /* Handle close to 1.0 inputs separately. */
+ /* Fix sign of zero with downward rounding when x==1. */
+ if (WANT_ROUNDING && __glibc_unlikely (ix == asuint64 (1.0)))
+ return 0;
+ r = x - 1.0;
+ r2 = r * r;
+ r3 = r * r2;
+ y = r3 * (B[1] + r * B[2] + r2 * B[3]
+ + r3 * (B[4] + r * B[5] + r2 * B[6]
+ + r3 * (B[7] + r * B[8] + r2 * B[9] + r3 * B[10])));
+ /* Worst-case error is around 0.507 ULP. */
+ w = r * 0x1p27;
+ double_t rhi = r + w - w;
+ double_t rlo = r - rhi;
+ w = rhi * rhi * B[0]; /* B[0] == -0.5. */
+ hi = r + w;
+ lo = r - hi + w;
+ lo += B[0] * rlo * (rhi + r);
+ y += lo;
+ y += hi;
+ return y;
}
- if (__glibc_unlikely (ux >= 0x7ff00000))
- return x + x; /* INF or NaN */
-
- /* Regular values of x */
-
- w = x - 1;
- if (__glibc_likely (fabs (w) > U03))
- goto case_03;
-
- /* log (1) is +0 in all rounding modes. */
- if (w == 0.0)
- return 0.0;
-
- /*--- The case abs(x-1) < 0.03 */
-
- t8 = MHALF * w;
- EMULV (t8, w, a, aa, t1, t2, t3, t4, t5);
- EADD (w, a, b, bb);
- /* Evaluate polynomial II */
- polII = b7.d + w * b8.d;
- polII = b6.d + w * polII;
- polII = b5.d + w * polII;
- polII = b4.d + w * polII;
- polII = b3.d + w * polII;
- polII = b2.d + w * polII;
- polII = b1.d + w * polII;
- polII = b0.d + w * polII;
- polII *= w * w * w;
- c = (aa + bb) + polII;
-
- /* Here b contains the high part of the result, and c the low part.
- Maximum error is b * 2.334e-19, so accuracy is >61 bits.
- Therefore max ULP error of b + c is ~0.502. */
- return b + c;
-
- /*--- The case abs(x-1) > 0.03 */
-case_03:
-
- /* Find n,u such that x = u*2**n, 1/sqrt(2) < u < sqrt(2) */
- n += (num.i[HIGH_HALF] >> 20) - 1023;
- num.i[HIGH_HALF] = (num.i[HIGH_HALF] & 0x000fffff) | 0x3ff00000;
- if (num.d > SQRT_2)
+ if (__glibc_unlikely (top - 0x0010 >= 0x7ff0 - 0x0010))
{
- num.d *= HALF;
- n++;
+ /* x < 0x1p-1022 or inf or nan. */
+ if (ix * 2 == 0)
+ return __math_divzero (1);
+ if (ix == asuint64 (INFINITY)) /* log(inf) == inf. */
+ return x;
+ if ((top & 0x8000) || (top & 0x7ff0) == 0x7ff0)
+ return __math_invalid (x);
+ /* x is subnormal, normalize it. */
+ ix = asuint64 (x * 0x1p52);
+ ix -= 52ULL << 52;
}
- u = num.d;
- dbl_n = (double) n;
-
- /* Find i such that ui=1+(i-75)/2**8 is closest to u (i= 0,1,2,...,181) */
- num.d += h1.d;
- i = (num.i[HIGH_HALF] & 0x000fffff) >> 12;
-
- /* Find j such that vj=1+(j-180)/2**16 is closest to v=u/ui (j= 0,...,361) */
- num.d = u * Iu[i].d + h2.d;
- j = (num.i[HIGH_HALF] & 0x000fffff) >> 4;
- /* Compute w=(u-ui*vj)/(ui*vj) */
- p0 = (1 + (i - 75) * DEL_U) * (1 + (j - 180) * DEL_V);
- q = u - p0;
- r0 = Iu[i].d * Iv[j].d;
- w = q * r0;
-
- /* Evaluate polynomial I */
- polI = w + (a2.d + a3.d * w) * w * w;
-
- /* Add up everything */
- nln2a = dbl_n * LN2A;
- luai = Lu[i][0].d;
- lubi = Lu[i][1].d;
- lvaj = Lv[j][0].d;
- lvbj = Lv[j][1].d;
- EADD (luai, lvaj, sij, ssij);
- EADD (nln2a, sij, A, ttij);
- B0 = (((lubi + lvbj) + ssij) + ttij) + dbl_n * LN2B;
- B = polI + B0;
-
- /* Here A contains the high part of the result, and B the low part.
- Maximum abs error is 6.095e-21 and min log (x) is 0.0295 since x > 1.03.
- Therefore max ULP error of A + B is ~0.502. */
- return A + B;
+ /* x = 2^k z; where z is in range [OFF,2*OFF) and exact.
+ The range is split into N subintervals.
+ The ith subinterval contains z and c is near its center. */
+ tmp = ix - OFF;
+ i = (tmp >> (52 - LOG_TABLE_BITS)) % N;
+ k = (int64_t) tmp >> 52; /* arithmetic shift */
+ iz = ix - (tmp & 0xfffULL << 52);
+ invc = T[i].invc;
+ logc = T[i].logc;
+ z = asdouble (iz);
+
+ /* log(x) = log1p(z/c-1) + log(c) + k*Ln2. */
+ /* r ~= z/c - 1, |r| < 1/(2*N). */
+#ifdef __FP_FAST_FMA
+ /* rounding error: 0x1p-55/N. */
+ r = __builtin_fma (z, invc, -1.0);
+#else
+ /* rounding error: 0x1p-55/N + 0x1p-66. */
+ r = (z - T2[i].chi - T2[i].clo) * invc;
+#endif
+ kd = (double_t) k;
+
+ /* hi + lo = r + log(c) + k*Ln2. */
+ w = kd * Ln2hi + logc;
+ hi = w + r;
+ lo = w - hi + r + kd * Ln2lo;
+
+ /* log(x) = lo + (log1p(r) - r) + hi. */
+ r2 = r * r; /* rounding error: 0x1p-54/N^2. */
+ /* Worst case error if |y| > 0x1p-4: 0.519 ULP (0.520 ULP without fma).
+ 0.5 + 2.06/N + abs-poly-error*2^56 ULP (+ 0.001 ULP without fma). */
+ y = lo + r2 * A[0] + r * r2 * (A[1] + r * A[2] + r2 * (A[3] + r * A[4])) + hi;
+ return y;
}
-
#ifndef __ieee754_log
strong_alias (__ieee754_log, __log_finite)
#endif
diff --git a/sysdeps/ieee754/dbl-64/e_log_data.c b/sysdeps/ieee754/dbl-64/e_log_data.c
new file mode 100644
index 0000000000..62f8a09c93
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/e_log_data.c
@@ -0,0 +1,347 @@
+/* Data for log.
+ Copyright (C) 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
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ The GNU C Library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, see
+ <http://www.gnu.org/licenses/>. */
+
+#include "math_config.h"
+
+#define N (1 << LOG_TABLE_BITS)
+
+const struct log_data __log_data = {
+.ln2hi = 0x1.62e42fefa3800p-1,
+.ln2lo = 0x1.ef35793c76730p-45,
+.poly1 = {
+#if LOG_POLY1_ORDER == 12
+// relative error: 0x1.c04d76cp-63
+// in -0x1p-4 0x1.09p-4 (|log(1+x)| > 0x1p-4 outside the interval)
+-0x1p-1,
+0x1.5555555555577p-2,
+-0x1.ffffffffffdcbp-3,
+0x1.999999995dd0cp-3,
+-0x1.55555556745a7p-3,
+0x1.24924a344de3p-3,
+-0x1.fffffa4423d65p-4,
+0x1.c7184282ad6cap-4,
+-0x1.999eb43b068ffp-4,
+0x1.78182f7afd085p-4,
+-0x1.5521375d145cdp-4,
+#endif
+},
+.poly = {
+#if N == 128 && LOG_POLY_ORDER == 6
+// relative error: 0x1.926199e8p-56
+// abs error: 0x1.882ff33p-65
+// in -0x1.fp-9 0x1.fp-9
+-0x1.0000000000001p-1,
+0x1.555555551305bp-2,
+-0x1.fffffffeb459p-3,
+0x1.999b324f10111p-3,
+-0x1.55575e506c89fp-3,
+#endif
+},
+/* Algorithm:
+
+ x = 2^k z
+ log(x) = k ln2 + log(c) + log(z/c)
+ log(z/c) = poly(z/c - 1)
+
+where z is in [1.6p-1; 1.6p0] which is split into N subintervals and z falls
+into the ith one, then table entries are computed as
+
+ tab[i].invc = 1/c
+ tab[i].logc = (double)log(c)
+ tab2[i].chi = (double)c
+ tab2[i].clo = (double)(c - (double)c)
+
+where c is near the center of the subinterval and is chosen by trying +-2^29
+floating point invc candidates around 1/center and selecting one for which
+
+ 1) the rounding error in 0x1.8p9 + logc is 0,
+ 2) the rounding error in z - chi - clo is < 0x1p-66 and
+ 3) the rounding error in (double)log(c) is minimized (< 0x1p-66).
+
+Note: 1) ensures that k*ln2hi + logc can be computed without rounding error,
+2) ensures that z/c - 1 can be computed as (z - chi - clo)*invc with close to
+a single rounding error when there is no fast fma for z*invc - 1, 3) ensures
+that logc + poly(z/c - 1) has small error, however near x == 1 when
+|log(x)| < 0x1p-4, this is not enough so that is special cased. */
+.tab = {
+#if N == 128
+{0x1.734f0c3e0de9fp+0, -0x1.7cc7f79e69000p-2},
+{0x1.713786a2ce91fp+0, -0x1.76feec20d0000p-2},
+{0x1.6f26008fab5a0p+0, -0x1.713e31351e000p-2},
+{0x1.6d1a61f138c7dp+0, -0x1.6b85b38287800p-2},
+{0x1.6b1490bc5b4d1p+0, -0x1.65d5590807800p-2},
+{0x1.69147332f0cbap+0, -0x1.602d076180000p-2},
+{0x1.6719f18224223p+0, -0x1.5a8ca86909000p-2},
+{0x1.6524f99a51ed9p+0, -0x1.54f4356035000p-2},
+{0x1.63356aa8f24c4p+0, -0x1.4f637c36b4000p-2},
+{0x1.614b36b9ddc14p+0, -0x1.49da7fda85000p-2},
+{0x1.5f66452c65c4cp+0, -0x1.445923989a800p-2},
+{0x1.5d867b5912c4fp+0, -0x1.3edf439b0b800p-2},
+{0x1.5babccb5b90dep+0, -0x1.396ce448f7000p-2},
+{0x1.59d61f2d91a78p+0, -0x1.3401e17bda000p-2},
+{0x1.5805612465687p+0, -0x1.2e9e2ef468000p-2},
+{0x1.56397cee76bd3p+0, -0x1.2941b3830e000p-2},
+{0x1.54725e2a77f93p+0, -0x1.23ec58cda8800p-2},
+{0x1.52aff42064583p+0, -0x1.1e9e129279000p-2},
+{0x1.50f22dbb2bddfp+0, -0x1.1956d2b48f800p-2},
+{0x1.4f38f4734ded7p+0, -0x1.141679ab9f800p-2},
+{0x1.4d843cfde2840p+0, -0x1.0edd094ef9800p-2},
+{0x1.4bd3ec078a3c8p+0, -0x1.09aa518db1000p-2},
+{0x1.4a27fc3e0258ap+0, -0x1.047e65263b800p-2},
+{0x1.4880524d48434p+0, -0x1.feb224586f000p-3},
+{0x1.46dce1b192d0bp+0, -0x1.f474a7517b000p-3},
+{0x1.453d9d3391854p+0, -0x1.ea4443d103000p-3},
+{0x1.43a2744b4845ap+0, -0x1.e020d44e9b000p-3},
+{0x1.420b54115f8fbp+0, -0x1.d60a22977f000p-3},
+{0x1.40782da3ef4b1p+0, -0x1.cc00104959000p-3},
+{0x1.3ee8f5d57fe8fp+0, -0x1.c202956891000p-3},
+{0x1.3d5d9a00b4ce9p+0, -0x1.b81178d811000p-3},
+{0x1.3bd60c010c12bp+0, -0x1.ae2c9ccd3d000p-3},
+{0x1.3a5242b75dab8p+0, -0x1.a45402e129000p-3},
+{0x1.38d22cd9fd002p+0, -0x1.9a877681df000p-3},
+{0x1.3755bc5847a1cp+0, -0x1.90c6d69483000p-3},
+{0x1.35dce49ad36e2p+0, -0x1.87120a645c000p-3},
+{0x1.34679984dd440p+0, -0x1.7d68fb4143000p-3},
+{0x1.32f5cceffcb24p+0, -0x1.73cb83c627000p-3},
+{0x1.3187775a10d49p+0, -0x1.6a39a9b376000p-3},
+{0x1.301c8373e3990p+0, -0x1.60b3154b7a000p-3},
+{0x1.2eb4ebb95f841p+0, -0x1.5737d76243000p-3},
+{0x1.2d50a0219a9d1p+0, -0x1.4dc7b8fc23000p-3},
+{0x1.2bef9a8b7fd2ap+0, -0x1.4462c51d20000p-3},
+{0x1.2a91c7a0c1babp+0, -0x1.3b08abc830000p-3},
+{0x1.293726014b530p+0, -0x1.31b996b490000p-3},
+{0x1.27dfa5757a1f5p+0, -0x1.2875490a44000p-3},
+{0x1.268b39b1d3bbfp+0, -0x1.1f3b9f879a000p-3},
+{0x1.2539d838ff5bdp+0, -0x1.160c8252ca000p-3},
+{0x1.23eb7aac9083bp+0, -0x1.0ce7f57f72000p-3},
+{0x1.22a012ba940b6p+0, -0x1.03cdc49fea000p-3},
+{0x1.2157996cc4132p+0, -0x1.f57bdbc4b8000p-4},
+{0x1.201201dd2fc9bp+0, -0x1.e370896404000p-4},
+{0x1.1ecf4494d480bp+0, -0x1.d17983ef94000p-4},
+{0x1.1d8f5528f6569p+0, -0x1.bf9674ed8a000p-4},
+{0x1.1c52311577e7cp+0, -0x1.adc79202f6000p-4},
+{0x1.1b17c74cb26e9p+0, -0x1.9c0c3e7288000p-4},
+{0x1.19e010c2c1ab6p+0, -0x1.8a646b372c000p-4},
+{0x1.18ab07bb670bdp+0, -0x1.78d01b3ac0000p-4},
+{0x1.1778a25efbcb6p+0, -0x1.674f145380000p-4},
+{0x1.1648d354c31dap+0, -0x1.55e0e6d878000p-4},
+{0x1.151b990275fddp+0, -0x1.4485cdea1e000p-4},
+{0x1.13f0ea432d24cp+0, -0x1.333d94d6aa000p-4},
+{0x1.12c8b7210f9dap+0, -0x1.22079f8c56000p-4},
+{0x1.11a3028ecb531p+0, -0x1.10e4698622000p-4},
+{0x1.107fbda8434afp+0, -0x1.ffa6c6ad20000p-5},
+{0x1.0f5ee0f4e6bb3p+0, -0x1.dda8d4a774000p-5},
+{0x1.0e4065d2a9fcep+0, -0x1.bbcece4850000p-5},
+{0x1.0d244632ca521p+0, -0x1.9a1894012c000p-5},
+{0x1.0c0a77ce2981ap+0, -0x1.788583302c000p-5},
+{0x1.0af2f83c636d1p+0, -0x1.5715e67d68000p-5},
+{0x1.09ddb98a01339p+0, -0x1.35c8a49658000p-5},
+{0x1.08cabaf52e7dfp+0, -0x1.149e364154000p-5},
+{0x1.07b9f2f4e28fbp+0, -0x1.e72c082eb8000p-6},
+{0x1.06ab58c358f19p+0, -0x1.a55f152528000p-6},
+{0x1.059eea5ecf92cp+0, -0x1.63d62cf818000p-6},
+{0x1.04949cdd12c90p+0, -0x1.228fb8caa0000p-6},
+{0x1.038c6c6f0ada9p+0, -0x1.c317b20f90000p-7},
+{0x1.02865137932a9p+0, -0x1.419355daa0000p-7},
+{0x1.0182427ea7348p+0, -0x1.81203c2ec0000p-8},
+{0x1.008040614b195p+0, -0x1.0040979240000p-9},
+{0x1.fe01ff726fa1ap-1, 0x1.feff384900000p-9},
+{0x1.fa11cc261ea74p-1, 0x1.7dc41353d0000p-7},
+{0x1.f6310b081992ep-1, 0x1.3cea3c4c28000p-6},
+{0x1.f25f63ceeadcdp-1, 0x1.b9fc114890000p-6},
+{0x1.ee9c8039113e7p-1, 0x1.1b0d8ce110000p-5},
+{0x1.eae8078cbb1abp-1, 0x1.58a5bd001c000p-5},
+{0x1.e741aa29d0c9bp-1, 0x1.95c8340d88000p-5},
+{0x1.e3a91830a99b5p-1, 0x1.d276aef578000p-5},
+{0x1.e01e009609a56p-1, 0x1.07598e598c000p-4},
+{0x1.dca01e577bb98p-1, 0x1.253f5e30d2000p-4},
+{0x1.d92f20b7c9103p-1, 0x1.42edd8b380000p-4},
+{0x1.d5cac66fb5ccep-1, 0x1.606598757c000p-4},
+{0x1.d272caa5ede9dp-1, 0x1.7da76356a0000p-4},
+{0x1.cf26e3e6b2ccdp-1, 0x1.9ab434e1c6000p-4},
+{0x1.cbe6da2a77902p-1, 0x1.b78c7bb0d6000p-4},
+{0x1.c8b266d37086dp-1, 0x1.d431332e72000p-4},
+{0x1.c5894bd5d5804p-1, 0x1.f0a3171de6000p-4},
+{0x1.c26b533bb9f8cp-1, 0x1.067152b914000p-3},
+{0x1.bf583eeece73fp-1, 0x1.147858292b000p-3},
+{0x1.bc4fd75db96c1p-1, 0x1.2266ecdca3000p-3},
+{0x1.b951e0c864a28p-1, 0x1.303d7a6c55000p-3},
+{0x1.b65e2c5ef3e2cp-1, 0x1.3dfc33c331000p-3},
+{0x1.b374867c9888bp-1, 0x1.4ba366b7a8000p-3},
+{0x1.b094b211d304ap-1, 0x1.5933928d1f000p-3},
+{0x1.adbe885f2ef7ep-1, 0x1.66acd2418f000p-3},
+{0x1.aaf1d31603da2p-1, 0x1.740f8ec669000p-3},
+{0x1.a82e63fd358a7p-1, 0x1.815c0f51af000p-3},
+{0x1.a5740ef09738bp-1, 0x1.8e92954f68000p-3},
+{0x1.a2c2a90ab4b27p-1, 0x1.9bb3602f84000p-3},
+{0x1.a01a01393f2d1p-1, 0x1.a8bed1c2c0000p-3},
+{0x1.9d79f24db3c1bp-1, 0x1.b5b515c01d000p-3},
+{0x1.9ae2505c7b190p-1, 0x1.c2967ccbcc000p-3},
+{0x1.9852ef297ce2fp-1, 0x1.cf635d5486000p-3},
+{0x1.95cbaeea44b75p-1, 0x1.dc1bd3446c000p-3},
+{0x1.934c69de74838p-1, 0x1.e8c01b8cfe000p-3},
+{0x1.90d4f2f6752e6p-1, 0x1.f5509c0179000p-3},
+{0x1.8e6528effd79dp-1, 0x1.00e6c121fb800p-2},
+{0x1.8bfce9fcc007cp-1, 0x1.071b80e93d000p-2},
+{0x1.899c0dabec30ep-1, 0x1.0d46b9e867000p-2},
+{0x1.87427aa2317fbp-1, 0x1.13687334bd000p-2},
+{0x1.84f00acb39a08p-1, 0x1.1980d67234800p-2},
+{0x1.82a49e8653e55p-1, 0x1.1f8ffe0cc8000p-2},
+{0x1.8060195f40260p-1, 0x1.2595fd7636800p-2},
+{0x1.7e22563e0a329p-1, 0x1.2b9300914a800p-2},
+{0x1.7beb377dcb5adp-1, 0x1.3187210436000p-2},
+{0x1.79baa679725c2p-1, 0x1.377266dec1800p-2},
+{0x1.77907f2170657p-1, 0x1.3d54ffbaf3000p-2},
+{0x1.756cadbd6130cp-1, 0x1.432eee32fe000p-2},
+#endif
+},
+#ifndef __FP_FAST_FMA
+.tab2 = {
+# if N == 128
+{0x1.61000014fb66bp-1, 0x1.e026c91425b3cp-56},
+{0x1.63000034db495p-1, 0x1.dbfea48005d41p-55},
+{0x1.650000d94d478p-1, 0x1.e7fa786d6a5b7p-55},
+{0x1.67000074e6fadp-1, 0x1.1fcea6b54254cp-57},
+{0x1.68ffffedf0faep-1, -0x1.c7e274c590efdp-56},
+{0x1.6b0000763c5bcp-1, -0x1.ac16848dcda01p-55},
+{0x1.6d0001e5cc1f6p-1, 0x1.33f1c9d499311p-55},
+{0x1.6efffeb05f63ep-1, -0x1.e80041ae22d53p-56},
+{0x1.710000e86978p-1, 0x1.bff6671097952p-56},
+{0x1.72ffffc67e912p-1, 0x1.c00e226bd8724p-55},
+{0x1.74fffdf81116ap-1, -0x1.e02916ef101d2p-57},
+{0x1.770000f679c9p-1, -0x1.7fc71cd549c74p-57},
+{0x1.78ffffa7ec835p-1, 0x1.1bec19ef50483p-55},
+{0x1.7affffe20c2e6p-1, -0x1.07e1729cc6465p-56},
+{0x1.7cfffed3fc9p-1, -0x1.08072087b8b1cp-55},
+{0x1.7efffe9261a76p-1, 0x1.dc0286d9df9aep-55},
+{0x1.81000049ca3e8p-1, 0x1.97fd251e54c33p-55},
+{0x1.8300017932c8fp-1, -0x1.afee9b630f381p-55},
+{0x1.850000633739cp-1, 0x1.9bfbf6b6535bcp-55},
+{0x1.87000204289c6p-1, -0x1.bbf65f3117b75p-55},
+{0x1.88fffebf57904p-1, -0x1.9006ea23dcb57p-55},
+{0x1.8b00022bc04dfp-1, -0x1.d00df38e04b0ap-56},
+{0x1.8cfffe50c1b8ap-1, -0x1.8007146ff9f05p-55},
+{0x1.8effffc918e43p-1, 0x1.3817bd07a7038p-55},
+{0x1.910001efa5fc7p-1, 0x1.93e9176dfb403p-55},
+{0x1.9300013467bb9p-1, 0x1.f804e4b980276p-56},
+{0x1.94fffe6ee076fp-1, -0x1.f7ef0d9ff622ep-55},
+{0x1.96fffde3c12d1p-1, -0x1.082aa962638bap-56},
+{0x1.98ffff4458a0dp-1, -0x1.7801b9164a8efp-55},
+{0x1.9afffdd982e3ep-1, -0x1.740e08a5a9337p-55},
+{0x1.9cfffed49fb66p-1, 0x1.fce08c19bep-60},
+{0x1.9f00020f19c51p-1, -0x1.a3faa27885b0ap-55},
+{0x1.a10001145b006p-1, 0x1.4ff489958da56p-56},
+{0x1.a300007bbf6fap-1, 0x1.cbeab8a2b6d18p-55},
+{0x1.a500010971d79p-1, 0x1.8fecadd78793p-55},
+{0x1.a70001df52e48p-1, -0x1.f41763dd8abdbp-55},
+{0x1.a90001c593352p-1, -0x1.ebf0284c27612p-55},
+{0x1.ab0002a4f3e4bp-1, -0x1.9fd043cff3f5fp-57},
+{0x1.acfffd7ae1ed1p-1, -0x1.23ee7129070b4p-55},
+{0x1.aefffee510478p-1, 0x1.a063ee00edea3p-57},
+{0x1.b0fffdb650d5b