Revert configurable mantissa patch

Reverts d22ca8cdfb

since it is severely broken on 32-bit.
This commit is contained in:
Siddhesh Poyarekar 2013-03-15 23:11:50 +05:30
parent 83a6b66ae9
commit 1e3803454e
5 changed files with 67 additions and 148 deletions

View File

@ -1,3 +1,10 @@
2013-03-15 Siddhesh Poyarekar <siddhesh@redhat.com>
* sysdeps/ieee754/dbl-64/mpa-arch.h: Remove.
* sysdeps/ieee754/dbl-64/mpa.c: Revert last change.
* sysdeps/ieee754/dbl-64/mpa.h: Revert last change.
* sysdeps/powerpc/power4/fpu/mpa-arch.h: Remove.
2013-03-15 Adhemerval Zanella <azanella@linux.vnet.ibm.com>
* sysdeps/unix/sysv/linux/powerpc/bits/libc-vdso.h (VDSO_IFUNC_RET): Add

View File

@ -1,47 +0,0 @@
/* Overridable constants and operations.
Copyright (C) 2013 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/>. */
#include <stdint.h>
typedef long mantissa_t;
typedef int64_t mantissa_store_t;
#define TWOPOW(i) (1L << i)
#define RADIX_EXP 24
#define RADIX TWOPOW (RADIX_EXP) /* 2^24 */
/* Divide D by RADIX and put the remainder in R. D must be a non-negative
integral value. */
#define DIV_RADIX(d, r) \
({ \
r = d & (RADIX - 1); \
d >>= RADIX_EXP; \
})
/* Put the integer component of a double X in R and retain the fraction in
X. This is used in extracting mantissa digits for MP_NO by using the
integer portion of the current value of the number as the current mantissa
digit and then scaling by RADIX to get the next mantissa digit in the same
manner. */
#define INTEGER_OF(x, i) \
({ \
i = (mantissa_t) x; \
x -= i; \
})
/* Align IN down to F. The code assumes that F is a power of two. */
#define ALIGN_DOWN_TO(in, f) ((in) & -(f))

View File

@ -125,8 +125,7 @@ norm (const mp_no *x, double *y, int p)
{
#define R RADIXI
long i;
double c;
mantissa_t a, u, v, z[5];
double a, c, u, v, z[5];
if (p < 5)
{
if (p == 1)
@ -148,14 +147,17 @@ norm (const mp_no *x, double *y, int p)
for (i = 2; i < 5; i++)
{
mantissa_t d, r;
d = X[i] * a;
DIV_RADIX (d, r);
z[i] = r;
z[i - 1] += d;
z[i] = X[i] * a;
u = (z[i] + CUTTER) - CUTTER;
if (u > z[i])
u -= RADIX;
z[i] -= u;
z[i - 1] += u * RADIXI;
}
u = ALIGN_DOWN_TO (z[3], TWO19);
u = (z[3] + TWO71) - TWO71;
if (u > z[3])
u -= TWO19;
v = z[3] - u;
if (v == TWO18)
@ -198,8 +200,7 @@ denorm (const mp_no *x, double *y, int p)
{
long i, k;
long p2 = p;
double c;
mantissa_t u, z[5];
double c, u, z[5];
#define R RADIXI
if (EX < -44 || (EX == -44 && X[1] < TWO5))
@ -279,7 +280,9 @@ denorm (const mp_no *x, double *y, int p)
z[3] = X[k];
}
u = ALIGN_DOWN_TO (z[3], TWO5);
u = (z[3] + TWO57) - TWO57;
if (u > z[3])
u -= TWO5;
if (u == z[3])
{
@ -327,6 +330,7 @@ __dbl_mp (double x, mp_no *y, int p)
{
long i, n;
long p2 = p;
double u;
/* Sign. */
if (x == ZERO)
@ -352,7 +356,11 @@ __dbl_mp (double x, mp_no *y, int p)
n = MIN (p2, 4);
for (i = 1; i <= n; i++)
{
INTEGER_OF (x, Y[i]);
u = (x + TWO52) - TWO52;
if (u > x)
u -= ONE;
Y[i] = u;
x -= u;
x *= RADIX;
}
for (; i <= p2; i++)
@ -369,7 +377,7 @@ add_magnitudes (const mp_no *x, const mp_no *y, mp_no *z, int p)
{
long i, j, k;
long p2 = p;
mantissa_t zk;
double zk;
EZ = EX;
@ -437,7 +445,7 @@ sub_magnitudes (const mp_no *x, const mp_no *y, mp_no *z, int p)
{
long i, j, k;
long p2 = p;
mantissa_t zk;
double zk;
EZ = EX;
i = p2;
@ -613,9 +621,9 @@ __mul (const mp_no *x, const mp_no *y, mp_no *z, int p)
{
long i, j, k, ip, ip2;
long p2 = p;
mantissa_store_t zk;
double u, zk;
const mp_no *a;
mantissa_store_t *diag;
double *diag;
/* Is z=0? */
if (__glibc_unlikely (X[0] * Y[0] == ZERO))
@ -672,8 +680,8 @@ __mul (const mp_no *x, const mp_no *y, mp_no *z, int p)
/* Precompute sums of diagonal elements so that we can directly use them
later. See the next comment to know we why need them. */
diag = alloca (k * sizeof (mantissa_store_t));
mantissa_store_t d = ZERO;
diag = alloca (k * sizeof (double));
double d = ZERO;
for (i = 1; i <= ip; i++)
{
d += X[i] * Y[i];
@ -696,8 +704,11 @@ __mul (const mp_no *x, const mp_no *y, mp_no *z, int p)
zk -= diag[k - 1];
DIV_RADIX (zk, Z[k]);
k--;
u = (zk + CUTTER) - CUTTER;
if (u > zk)
u -= RADIX;
Z[k--] = zk - u;
zk = u * RADIXI;
}
/* The real deal. Mantissa digit Z[k] is the sum of all X[i] * Y[j] where i
@ -727,8 +738,11 @@ __mul (const mp_no *x, const mp_no *y, mp_no *z, int p)
zk -= diag[k - 1];
DIV_RADIX (zk, Z[k]);
k--;
u = (zk + CUTTER) - CUTTER;
if (u > zk)
u -= RADIX;
Z[k--] = zk - u;
zk = u * RADIXI;
}
Z[k] = zk;
@ -760,7 +774,7 @@ SECTION
__sqr (const mp_no *x, mp_no *y, int p)
{
long i, j, k, ip;
mantissa_store_t yk;
double u, yk;
/* Is z=0? */
if (__glibc_unlikely (X[0] == ZERO))
@ -784,7 +798,7 @@ __sqr (const mp_no *x, mp_no *y, int p)
while (k > p)
{
mantissa_store_t yk2 = 0;
double yk2 = 0.0;
long lim = k / 2;
if (k % 2 == 0)
@ -804,13 +818,16 @@ __sqr (const mp_no *x, mp_no *y, int p)
yk += 2.0 * yk2;
DIV_RADIX (yk, Y[k]);
k--;
u = (yk + CUTTER) - CUTTER;
if (u > yk)
u -= RADIX;
Y[k--] = yk - u;
yk = u * RADIXI;
}
while (k > 1)
{
mantissa_store_t yk2 = 0;
double yk2 = 0.0;
long lim = k / 2;
if (k % 2 == 0)
@ -822,8 +839,11 @@ __sqr (const mp_no *x, mp_no *y, int p)
yk += 2.0 * yk2;
DIV_RADIX (yk, Y[k]);
k--;
u = (yk + CUTTER) - CUTTER;
if (u > yk)
u -= RADIX;
Y[k--] = yk - u;
yk = u * RADIXI;
}
Y[k] = yk;

View File

@ -35,7 +35,6 @@
/* Common types and definition */
/************************************************************************/
#include <mpa-arch.h>
/* The mp_no structure holds the details of a multi-precision floating point
number.
@ -62,7 +61,7 @@
typedef struct
{
int e;
mantissa_t d[40];
double d[40];
} mp_no;
typedef union
@ -83,13 +82,9 @@ extern const mp_no mptwo;
#define ABS(x) ((x) < 0 ? -(x) : (x))
#ifndef RADIXI
# define RADIXI 0x1.0p-24 /* 2^-24 */
#endif
#ifndef TWO52
# define TWO52 0x1.0p52 /* 2^52 */
#endif
#define RADIX 0x1.0p24 /* 2^24 */
#define RADIXI 0x1.0p-24 /* 2^-24 */
#define CUTTER 0x1.0p76 /* 2^76 */
#define ZERO 0.0 /* 0 */
#define MZERO -0.0 /* 0 with the sign bit set */
@ -97,13 +92,13 @@ extern const mp_no mptwo;
#define MONE -1.0 /* -1 */
#define TWO 2.0 /* 2 */
#define TWO5 TWOPOW (5) /* 2^5 */
#define TWO8 TWOPOW (8) /* 2^52 */
#define TWO10 TWOPOW (10) /* 2^10 */
#define TWO18 TWOPOW (18) /* 2^18 */
#define TWO19 TWOPOW (19) /* 2^19 */
#define TWO23 TWOPOW (23) /* 2^23 */
#define TWO5 0x1.0p5 /* 2^5 */
#define TWO8 0x1.0p8 /* 2^52 */
#define TWO10 0x1.0p10 /* 2^10 */
#define TWO18 0x1.0p18 /* 2^18 */
#define TWO19 0x1.0p19 /* 2^19 */
#define TWO23 0x1.0p23 /* 2^23 */
#define TWO52 0x1.0p52 /* 2^52 */
#define TWO57 0x1.0p57 /* 2^57 */
#define TWO71 0x1.0p71 /* 2^71 */
#define TWOM1032 0x1.0p-1032 /* 2^-1032 */

View File

@ -1,56 +0,0 @@
/* Overridable constants and operations.
Copyright (C) 2013 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/>. */
typedef double mantissa_t;
typedef double mantissa_store_t;
#define TWOPOW(i) (0x1.0p##i)
#define RADIX TWOPOW (24) /* 2^24 */
#define CUTTER TWOPOW (76) /* 2^76 */
#define RADIXI 0x1.0p-24 /* 2^-24 */
#define TWO52 TWOPOW (52) /* 2^52 */
/* Divide D by RADIX and put the remainder in R. */
#define DIV_RADIX(d,r) \
({ \
double u = ((d) + CUTTER) - CUTTER; \
if (u > (d)) \
u -= RADIX; \
r = (d) - u; \
(d) = u * RADIXI; \
})
/* Put the integer component of a double X in R and retain the fraction in
X. */
#define INTEGER_OF(x, r) \
({ \
double u = ((x) + TWO52) - TWO52; \
if (u > (x)) \
u -= ONE; \
(r) = u; \
(x) -= u; \
})
/* Align IN down to a multiple of F, where F is a power of two. */
#define ALIGN_DOWN_TO(in, f) \
({ \
double factor = f * TWO52; \
double u = (in + factor) - factor; \
if (u > in) \
u -= f; \
u; \
})