[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[Commits] r22724 - in /fsf/trunk/libc: ./ ports/ ports/sysdeps/arm/ sysdeps/ieee754/dbl-64/ sysdeps/powerpc/fpu/ sysdeps/powerpc/power...



Author: eglibc
Date: Wed Mar 27 00:01:57 2013
New Revision: 22724

Log:
Import glibc-mainline for 2013-03-27

Added:
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa-arch.h
    fsf/trunk/libc/sysdeps/powerpc/power4/fpu/mpa-arch.h
Modified:
    fsf/trunk/libc/ChangeLog
    fsf/trunk/libc/ports/ChangeLog.arm
    fsf/trunk/libc/ports/sysdeps/arm/preconfigure
    fsf/trunk/libc/ports/sysdeps/arm/preconfigure.in
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/branred.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/dosincos.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_asin.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_atan2.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_exp.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_log.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_pow.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_remainder.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.h
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_atan.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_sin.c
    fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_tan.c
    fsf/trunk/libc/sysdeps/powerpc/fpu/s_llround.c

Modified: fsf/trunk/libc/ChangeLog
==============================================================================
--- fsf/trunk/libc/ChangeLog (original)
+++ fsf/trunk/libc/ChangeLog Wed Mar 27 00:01:57 2013
@@ -1,3 +1,62 @@
+2013-03-26  Siddhesh Poyarekar  <siddhesh@xxxxxxxxxx>
+
+	* sysdeps/ieee754/dbl-64/mpa.c (__acr): Use integral
+	constants.
+	(norm): Likewise.
+	(denorm): Likewise.
+	(__dbl_mp): Likewise.
+	(add_magnitudes): Likewise.
+	(sub_magnitudes): Likewise.
+	(__add): Likewise.
+	(__sub): Likewise.
+	(__mul): Likewise.
+	(__sqr): Likewise.
+	(__inv): Likewise.
+	(__dvd): Likewise.
+
+	* sysdeps/ieee754/dbl-64/branred.c (__branred): Remove
+	commented code.
+	* sysdeps/ieee754/dbl-64/dosincos.c (__dubsin): Likewise.
+	(__dubcos): Likewise.
+	* sysdeps/ieee754/dbl-64/e_asin.c (__ieee754_asin): Likewise.
+	(__ieee754_acos): Likewise.
+	* sysdeps/ieee754/dbl-64/e_atan2.c (__ieee754_atan2): Likewise.
+	* sysdeps/ieee754/dbl-64/e_exp.c (__ieee754_exp): Likewise.
+	(__exp1): Likewise.
+	* sysdeps/ieee754/dbl-64/e_log.c (__ieee754_log): Likewise.
+	* sysdeps/ieee754/dbl-64/e_pow.c (__ieee754_pow): Likewise.
+	(log1): Likewise.
+	(my_log2): Likewise.
+	(checkint): Likewise.
+	* sysdeps/ieee754/dbl-64/e_remainder.c
+	(__ieee754_remainder): Likewise.
+	* sysdeps/ieee754/dbl-64/s_atan.c (atan): Likewise.
+	* sysdeps/ieee754/dbl-64/s_sin.c (__sin): Likewise.
+	(bsloww): Likewise.
+	* sysdeps/ieee754/dbl-64/s_tan.c (tan): Likewise.
+
+	* sysdeps/ieee754/dbl-64/mpa-arch.h: New file.
+	* sysdeps/ieee754/dbl-64/mpa.c (norm): Use MANTISSA_T and
+	MANTISSA_STORE_T to store computations on mantissa.  Use
+	macros for rounding and division.
+	(denorm): Likewise.
+	(__dbl_mp): Likewise.
+	(add_magnitudes): Likewise.
+	(sub_magnitudes): Likewise.
+	(__mul): Likewise.
+	(__sqr): Likewise.
+	* sysdeps/ieee754/dbl-64/mpa.h: Include mpa-arch.h.  Define
+	powers of two in terms of TWOPOW macro.
+	(mp_no): Make type of mantissa as MANTISSA_T.
+	[!RADIXI]: Define RADIXI.
+	[!TWO52]: Define TWO52.
+	* sysdeps/powerpc/power4/fpu/mpa-arch.h: New file.
+
+2013-03-25  Adhemerval Zanella  <azanella@xxxxxxxxxxxxxxxxxx>
+
+	* sysdeps/powerpc/fpu/s_llround.c: Fix libm ABI issue with missing
+	llroundl symbol when building for PPC32.
+
 2013-03-24  Mark H Weaver  <mhw@xxxxxxxxxx>
 
 	* manual/arith.texi (Normalization Functions): Fix prototypes for

Modified: fsf/trunk/libc/ports/ChangeLog.arm
==============================================================================
--- fsf/trunk/libc/ports/ChangeLog.arm (original)
+++ fsf/trunk/libc/ports/ChangeLog.arm Wed Mar 27 00:01:57 2013
@@ -1,3 +1,8 @@
+2013-03-26  Mans Rullgard  <mans@xxxxxxxxx>
+
+	* sysdeps/arm/preconfigure.in: Use "test" instead of [ ].
+	* sysdeps/arm/preconfigure: Regenerated.
+
 2013-03-20  Joseph Myers  <joseph@xxxxxxxxxxxxxxxx>
 
 	* sysdeps/arm/configure.in (default-abi): Set using

Modified: fsf/trunk/libc/ports/sysdeps/arm/preconfigure
==============================================================================
--- fsf/trunk/libc/ports/sysdeps/arm/preconfigure (original)
+++ fsf/trunk/libc/ports/sysdeps/arm/preconfigure Wed Mar 27 00:01:57 2013
@@ -10,7 +10,7 @@
     # avoid this, add -fno-unwind-tables here and remove it in
     # sysdeps/unix/sysv/linux/arm/configure.in after those tests have
     # been run.
-    if  "${CFLAGS+set}" != "set" ; then
+    if test "${CFLAGS+set}" != "set"; then
       CFLAGS="-g -O2"
     fi
     CFLAGS="$CFLAGS -fno-unwind-tables"

Modified: fsf/trunk/libc/ports/sysdeps/arm/preconfigure.in
==============================================================================
--- fsf/trunk/libc/ports/sysdeps/arm/preconfigure.in (original)
+++ fsf/trunk/libc/ports/sysdeps/arm/preconfigure.in Wed Mar 27 00:01:57 2013
@@ -10,7 +10,7 @@
     # avoid this, add -fno-unwind-tables here and remove it in
     # sysdeps/unix/sysv/linux/arm/configure.in after those tests have
     # been run.
-    if [ "${CFLAGS+set}" != "set" ]; then
+    if test "${CFLAGS+set}" != "set"; then
       CFLAGS="-g -O2"
     fi
     CFLAGS="$CFLAGS -fno-unwind-tables"

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/branred.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/branred.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/branred.c Wed Mar 27 00:01:57 2013
@@ -53,13 +53,7 @@
 __branred(double x, double *a, double *aa)
 {
   int i,k;
-#if 0
-  int n;
-#endif
   mynumber  u,gor;
-#if 0
-  mynumber v;
-#endif
   double r[6],s,t,sum,b,bb,sum1,sum2,b1,bb1,b2,bb2,x1,x2,t1,t2;
 
   x*=tm600.x;

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/dosincos.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/dosincos.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/dosincos.c Wed Mar 27 00:01:57 2013
@@ -63,9 +63,6 @@
 #ifndef DLA_FMS
   double p,hx,tx,hy,ty,q;
 #endif
-#if 0
-  double xx,y,yy,z,zz;
-#endif
   mynumber u;
   int4 k;
 
@@ -118,9 +115,6 @@
     sn,ssn,cs,ccs,ds,dss,dc,dcc;
 #ifndef DLA_FMS
   double p,hx,tx,hy,ty,q;
-#endif
-#if 0
-  double xx,y,yy,z,zz;
 #endif
   mynumber u;
   int4 k;

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_asin.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_asin.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_asin.c Wed Mar 27 00:01:57 2013
@@ -62,9 +62,6 @@
   double x1,x2,xx,s1,s2,res1,p,t,res,r,cor,cc,y,c,z,w[2];
   mynumber u,v;
   int4 k,m,n;
-#if 0
-  int4 nn;
-#endif
 
   u.x = x;
   m = u.i[HIGH_HALF];
@@ -344,14 +341,8 @@
 __ieee754_acos(double x)
 {
   double x1,x2,xx,s1,s2,res1,p,t,res,r,cor,cc,y,c,z,w[2],eps;
-#if 0
-  double fc;
-#endif
   mynumber u,v;
   int4 k,m,n;
-#if 0
-  int4 nn;
-#endif
   u.x = x;
   m = u.i[HIGH_HALF];
   k = 0x7fffffff&m;

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_atan2.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_atan2.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_atan2.c Wed Mar 27 00:01:57 2013
@@ -67,22 +67,13 @@
 __ieee754_atan2(double y,double x) {
 
   int i,de,ux,dx,uy,dy;
-#if 0
-  int p;
-#endif
   static const int pr[MM]={6,8,10,20,32};
   double ax,ay,u,du,u9,ua,v,vv,dv,t1,t2,t3,t7,t8,
 	 z,zz,cor,s1,ss1,s2,ss2;
 #ifndef DLA_FMS
   double t4,t5,t6;
 #endif
-#if 0
-  double z1,z2;
-#endif
   number num;
-#if 0
-  mp_no mperr,mpt1,mpx,mpy,mpz,mpz1,mpz2;
-#endif
 
   static const int ep= 59768832,   /*  57*16**5   */
 		   em=-59768832;   /* -57*16**5   */

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_exp.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_exp.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_exp.c Wed Mar 27 00:01:57 2013
@@ -55,9 +55,6 @@
 __ieee754_exp(double x) {
   double bexp, t, eps, del, base, y, al, bet, res, rem, cor;
   mynumber junk1, junk2, binexp  = {{0,0}};
-#if 0
-  int4 k;
-#endif
   int4 i,j,m,n,ex;
   double retval;
 
@@ -174,9 +171,6 @@
 __exp1(double x, double xx, double error) {
   double bexp, t, eps, del, base, y, al, bet, res, rem, cor;
   mynumber junk1, junk2, binexp  = {{0,0}};
-#if 0
-  int4 k;
-#endif
   int4 i,j,m,n,ex;
 
   junk1.x = x;

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_log.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_log.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_log.c Wed Mar 27 00:01:57 2013
@@ -56,9 +56,6 @@
 #define M 4
   static const int pr[M]={8,10,18,32};
   int i,j,n,ux,dx,p;
-#if 0
-  int k;
-#endif
   double dbl_n,u,p0,q,r0,w,nln2a,luai,lubi,lvaj,lvbj,
 	 sij,ssij,ttij,A,B,B0,y,y1,y2,polI,polII,sa,sb,
 	 t1,t2,t7,t8,t,ra,rb,ww,

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_pow.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_pow.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_pow.c Wed Mar 27 00:01:57 2013
@@ -64,9 +64,6 @@
 SECTION
 __ieee754_pow(double x, double y) {
   double z,a,aa,error, t,a1,a2,y1,y2;
-#if 0
-  double gor=1.0;
-#endif
   mynumber u,v;
   int k;
   int4 qx,qy;
@@ -206,13 +203,7 @@
 SECTION
 log1(double x, double *delta, double *error) {
   int i,j,m;
-#if 0
-  int n;
-#endif
   double uu,vv,eps,nx,e,e1,e2,t,t1,t2,res,add=0;
-#if 0
-  double cor;
-#endif
   mynumber u,v;
 #ifdef BIG_ENDI
   mynumber
@@ -300,13 +291,7 @@
 SECTION
 my_log2(double x, double *delta, double *error) {
   int i,j,m;
-#if 0
-  int n;
-#endif
   double uu,vv,eps,nx,e,e1,e2,t,t1,t2,res,add=0;
-#if 0
-  double cor;
-#endif
   double ou1,ou2,lu1,lu2,ov,lv1,lv2,a,a1,a2;
   double y,yy,z,zz,j1,j2,j7,j8;
 #ifndef DLA_FMS
@@ -397,9 +382,6 @@
 checkint(double x) {
   union {int4 i[2]; double x;} u;
   int k,m,n;
-#if 0
-  int l;
-#endif
   u.x = x;
   m = u.i[HIGH_HALF]&0x7fffffff;    /* no sign */
   if (m >= 0x7ff00000) return 0;    /*  x is +/-inf or NaN  */

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_remainder.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_remainder.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/e_remainder.c Wed Mar 27 00:01:57 2013
@@ -42,13 +42,7 @@
 double __ieee754_remainder(double x, double y)
 {
   double z,d,xx;
-#if 0
-  double yy;
-#endif
   int4 kx,ky,n,nn,n1,m1,l;
-#if 0
-  int4 m;
-#endif
   mynumber u,t,w={{0,0}},v={{0,0}},ww={{0,0}},r;
   u.x=x;
   t.x=y;

Added: fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa-arch.h
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa-arch.h (added)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa-arch.h Wed Mar 27 00:01:57 2013
@@ -1,0 +1,47 @@
+/* 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))

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.c Wed Mar 27 00:01:57 2013
@@ -80,14 +80,14 @@
 {
   long i;
 
-  if (X[0] == ZERO)
-    {
-      if (Y[0] == ZERO)
+  if (X[0] == 0)
+    {
+      if (Y[0] == 0)
 	i = 0;
       else
 	i = -1;
     }
-  else if (Y[0] == ZERO)
+  else if (Y[0] == 0)
     i = 1;
   else
     {
@@ -125,7 +125,8 @@
 {
 #define R RADIXI
   long i;
-  double a, c, u, v, z[5];
+  double c;
+  mantissa_t a, u, v, z[5];
   if (p < 5)
     {
       if (p == 1)
@@ -139,44 +140,41 @@
     }
   else
     {
-      for (a = ONE, z[1] = X[1]; z[1] < TWO23;)
-	{
-	  a *= TWO;
-	  z[1] *= TWO;
+      for (a = 1, z[1] = X[1]; z[1] < TWO23;)
+	{
+	  a *= 2;
+	  z[1] *= 2;
 	}
 
       for (i = 2; i < 5; i++)
 	{
-	  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 = (z[3] + TWO71) - TWO71;
-      if (u > z[3])
-	u -= TWO19;
+	  mantissa_store_t d, r;
+	  d = X[i] * (mantissa_store_t) a;
+	  DIV_RADIX (d, r);
+	  z[i] = r;
+	  z[i - 1] += d;
+	}
+
+      u = ALIGN_DOWN_TO (z[3], TWO19);
       v = z[3] - u;
 
       if (v == TWO18)
 	{
-	  if (z[4] == ZERO)
+	  if (z[4] == 0)
 	    {
 	      for (i = 5; i <= p; i++)
 		{
-		  if (X[i] == ZERO)
+		  if (X[i] == 0)
 		    continue;
 		  else
 		    {
-		      z[3] += ONE;
+		      z[3] += 1;
 		      break;
 		    }
 		}
 	    }
 	  else
-	    z[3] += ONE;
+	    z[3] += 1;
 	}
 
       c = (z[1] + R * (z[2] + R * z[3])) / a;
@@ -200,12 +198,13 @@
 {
   long i, k;
   long p2 = p;
-  double c, u, z[5];
+  double c;
+  mantissa_t u, z[5];
 
 #define R RADIXI
   if (EX < -44 || (EX == -44 && X[1] < TWO5))
     {
-      *y = ZERO;
+      *y = 0;
       return;
     }
 
@@ -214,21 +213,21 @@
       if (EX == -42)
 	{
 	  z[1] = X[1] + TWO10;
-	  z[2] = ZERO;
-	  z[3] = ZERO;
+	  z[2] = 0;
+	  z[3] = 0;
 	  k = 3;
 	}
       else if (EX == -43)
 	{
 	  z[1] = TWO10;
 	  z[2] = X[1];
-	  z[3] = ZERO;
+	  z[3] = 0;
 	  k = 2;
 	}
       else
 	{
 	  z[1] = TWO10;
-	  z[2] = ZERO;
+	  z[2] = 0;
 	  z[3] = X[1];
 	  k = 1;
 	}
@@ -239,7 +238,7 @@
 	{
 	  z[1] = X[1] + TWO10;
 	  z[2] = X[2];
-	  z[3] = ZERO;
+	  z[3] = 0;
 	  k = 3;
 	}
       else if (EX == -43)
@@ -252,7 +251,7 @@
       else
 	{
 	  z[1] = TWO10;
-	  z[2] = ZERO;
+	  z[2] = 0;
 	  z[3] = X[1];
 	  k = 1;
 	}
@@ -274,25 +273,23 @@
       else
 	{
 	  z[1] = TWO10;
-	  z[2] = ZERO;
+	  z[2] = 0;
 	  k = 1;
 	}
       z[3] = X[k];
     }
 
-  u = (z[3] + TWO57) - TWO57;
-  if (u > z[3])
-    u -= TWO5;
+  u = ALIGN_DOWN_TO (z[3], TWO5);
 
   if (u == z[3])
     {
       for (i = k + 1; i <= p2; i++)
 	{
-	  if (X[i] == ZERO)
+	  if (X[i] == 0)
 	    continue;
 	  else
 	    {
-	      z[3] += ONE;
+	      z[3] += 1;
 	      break;
 	    }
 	}
@@ -309,9 +306,9 @@
 void
 __mp_dbl (const mp_no *x, double *y, int p)
 {
-  if (X[0] == ZERO)
-    {
-      *y = ZERO;
+  if (X[0] == 0)
+    {
+      *y = 0;
       return;
     }
 
@@ -330,41 +327,36 @@
 {
   long i, n;
   long p2 = p;
-  double u;
 
   /* Sign.  */
-  if (x == ZERO)
-    {
-      Y[0] = ZERO;
-      return;
-    }
-  else if (x > ZERO)
-    Y[0] = ONE;
+  if (x == 0)
+    {
+      Y[0] = 0;
+      return;
+    }
+  else if (x > 0)
+    Y[0] = 1;
   else
     {
-      Y[0] = MONE;
+      Y[0] = -1;
       x = -x;
     }
 
   /* Exponent.  */
-  for (EY = ONE; x >= RADIX; EY += ONE)
+  for (EY = 1; x >= RADIX; EY += 1)
     x *= RADIXI;
-  for (; x < ONE; EY -= ONE)
+  for (; x < 1; EY -= 1)
     x *= RADIX;
 
   /* Digits.  */
   n = MIN (p2, 4);
   for (i = 1; i <= n; i++)
     {
-      u = (x + TWO52) - TWO52;
-      if (u > x)
-	u -= ONE;
-      Y[i] = u;
-      x -= u;
+      INTEGER_OF (x, Y[i]);
       x *= RADIX;
     }
   for (; i <= p2; i++)
-    Y[i] = ZERO;
+    Y[i] = 0;
 }
 
 /* Add magnitudes of *X and *Y assuming that abs (*X) >= abs (*Y) > 0.  The
@@ -377,7 +369,7 @@
 {
   long i, j, k;
   long p2 = p;
-  double zk;
+  mantissa_t zk;
 
   EZ = EX;
 
@@ -391,7 +383,7 @@
       return;
     }
 
-  zk = ZERO;
+  zk = 0;
 
   for (; j > 0; i--, j--)
     {
@@ -399,12 +391,12 @@
       if (zk >= RADIX)
 	{
 	  Z[k--] = zk - RADIX;
-	  zk = ONE;
+	  zk = 1;
 	}
       else
         {
 	  Z[k--] = zk;
-	  zk = ZERO;
+	  zk = 0;
 	}
     }
 
@@ -414,16 +406,16 @@
       if (zk >= RADIX)
 	{
 	  Z[k--] = zk - RADIX;
-	  zk = ONE;
+	  zk = 1;
 	}
       else
         {
 	  Z[k--] = zk;
-	  zk = ZERO;
-	}
-    }
-
-  if (zk == ZERO)
+	  zk = 0;
+	}
+    }
+
+  if (zk == 0)
     {
       for (i = 1; i <= p2; i++)
 	Z[i] = Z[i + 1];
@@ -431,7 +423,7 @@
   else
     {
       Z[1] = zk;
-      EZ += ONE;
+      EZ += 1;
     }
 }
 
@@ -445,7 +437,7 @@
 {
   long i, j, k;
   long p2 = p;
-  double zk;
+  mantissa_t zk;
 
   EZ = EX;
   i = p2;
@@ -461,27 +453,27 @@
 
   /* The relevant least significant digit in Y is non-zero, so we factor it in
      to enhance accuracy.  */
-  if (j < p2 && Y[j + 1] > ZERO)
+  if (j < p2 && Y[j + 1] > 0)
     {
       Z[k + 1] = RADIX - Y[j + 1];
-      zk = MONE;
+      zk = -1;
     }
   else
-    zk = Z[k + 1] = ZERO;
+    zk = Z[k + 1] = 0;
 
   /* Subtract and borrow.  */
   for (; j > 0; i--, j--)
     {
       zk += (X[i] - Y[j]);
-      if (zk < ZERO)
+      if (zk < 0)
 	{
 	  Z[k--] = zk + RADIX;
-	  zk = MONE;
+	  zk = -1;
 	}
       else
         {
 	  Z[k--] = zk;
-	  zk = ZERO;
+	  zk = 0;
 	}
     }
 
@@ -489,25 +481,25 @@
   for (; i > 0; i--)
     {
       zk += X[i];
-      if (zk < ZERO)
+      if (zk < 0)
 	{
 	  Z[k--] = zk + RADIX;
-	  zk = MONE;
+	  zk = -1;
 	}
       else
         {
 	  Z[k--] = zk;
-	  zk = ZERO;
+	  zk = 0;
 	}
     }
 
   /* Normalize.  */
-  for (i = 1; Z[i] == ZERO; i++);
+  for (i = 1; Z[i] == 0; i++);
   EZ = EZ - i + 1;
   for (k = 1; i <= p2 + 1;)
     Z[k++] = Z[i++];
   for (; k <= p2;)
-    Z[k++] = ZERO;
+    Z[k++] = 0;
 }
 
 /* Add *X and *Y and store the result in *Z.  X and Y may overlap, but not X
@@ -519,12 +511,12 @@
 {
   int n;
 
-  if (X[0] == ZERO)
+  if (X[0] == 0)
     {
       __cpy (y, z, p);
       return;
     }
-  else if (Y[0] == ZERO)
+  else if (Y[0] == 0)
     {
       __cpy (x, z, p);
       return;
@@ -556,7 +548,7 @@
 	  Z[0] = Y[0];
 	}
       else
-	Z[0] = ZERO;
+	Z[0] = 0;
     }
 }
 
@@ -569,13 +561,13 @@
 {
   int n;
 
-  if (X[0] == ZERO)
+  if (X[0] == 0)
     {
       __cpy (y, z, p);
       Z[0] = -Z[0];
       return;
     }
-  else if (Y[0] == ZERO)
+  else if (Y[0] == 0)
     {
       __cpy (x, z, p);
       return;
@@ -607,7 +599,7 @@
 	  Z[0] = -Y[0];
 	}
       else
-	Z[0] = ZERO;
+	Z[0] = 0;
     }
 }
 
@@ -621,28 +613,28 @@
 {
   long i, j, k, ip, ip2;
   long p2 = p;
-  double u, zk;
+  mantissa_store_t zk;
   const mp_no *a;
-  double *diag;
+  mantissa_store_t *diag;
 
   /* Is z=0?  */
-  if (__glibc_unlikely (X[0] * Y[0] == ZERO))
-    {
-      Z[0] = ZERO;
+  if (__glibc_unlikely (X[0] * Y[0] == 0))
+    {
+      Z[0] = 0;
       return;
     }
 
   /* We need not iterate through all X's and Y's since it's pointless to
      multiply zeroes.  Here, both are zero...  */
   for (ip2 = p2; ip2 > 0; ip2--)
-    if (X[ip2] != ZERO || Y[ip2] != ZERO)
+    if (X[ip2] != 0 || Y[ip2] != 0)
       break;
 
-  a = X[ip2] != ZERO ? y : x;
+  a = X[ip2] != 0 ? y : x;
 
   /* ... and here, at least one of them is still zero.  */
   for (ip = ip2; ip > 0; ip--)
-    if (a->d[ip] != ZERO)
+    if (a->d[ip] != 0)
       break;
 
   /* The product looks like this for p = 3 (as an example):
@@ -669,22 +661,22 @@
      Another thing that becomes evident is that only the most significant
      ip+ip2 digits of the result are non-zero, where ip and ip2 are the
      'internal precision' of the input numbers, i.e. digits after ip and ip2
-     are all ZERO.  */
+     are all 0.  */
 
   k = (__glibc_unlikely (p2 < 3)) ? p2 + p2 : p2 + 3;
 
   while (k > ip + ip2 + 1)
-    Z[k--] = ZERO;
-
-  zk = ZERO;
+    Z[k--] = 0;
+
+  zk = 0;
 
   /* 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 (double));
-  double d = ZERO;
+  diag = alloca (k * sizeof (mantissa_store_t));
+  mantissa_store_t d = 0;
   for (i = 1; i <= ip; i++)
     {
-      d += X[i] * Y[i];
+      d += X[i] * (mantissa_store_t) Y[i];
       diag[i] = d;
     }
   while (i < k)
@@ -697,18 +689,15 @@
       if (k % 2 == 0)
 	/* We want to add this only once, but since we subtract it in the sum
 	   of products above, we add twice.  */
-	zk += 2 * X[lim] * Y[lim];
+	zk += 2 * X[lim] * (mantissa_store_t) Y[lim];
 
       for (i = k - p2, j = p2; i < j; i++, j--)
-	zk += (X[i] + X[j]) * (Y[i] + Y[j]);
+	zk += (X[i] + X[j]) * (mantissa_store_t) (Y[i] + Y[j]);
 
       zk -= diag[k - 1];
 
-      u = (zk + CUTTER) - CUTTER;
-      if (u > zk)
-	u -= RADIX;
-      Z[k--] = zk - u;
-      zk = u * RADIXI;
+      DIV_RADIX (zk, Z[k]);
+      k--;
     }
 
   /* The real deal.  Mantissa digit Z[k] is the sum of all X[i] * Y[j] where i
@@ -731,18 +720,15 @@
       if (k % 2 == 0)
 	/* We want to add this only once, but since we subtract it in the sum
 	   of products above, we add twice.  */
-        zk += 2 * X[lim] * Y[lim];
+        zk += 2 * X[lim] * (mantissa_store_t) Y[lim];
 
       for (i = 1, j = k - 1; i < j; i++, j--)
-	zk += (X[i] + X[j]) * (Y[i] + Y[j]);
+	zk += (X[i] + X[j]) * (mantissa_store_t) (Y[i] + Y[j]);
 
       zk -= diag[k - 1];
 
-      u = (zk + CUTTER) - CUTTER;
-      if (u > zk)
-	u -= RADIX;
-      Z[k--] = zk - u;
-      zk = u * RADIXI;
+      DIV_RADIX (zk, Z[k]);
+      k--;
     }
   Z[k] = zk;
 
@@ -752,7 +738,7 @@
   int e = EX + EY;
 
   /* Is there a carry beyond the most significant digit?  */
-  if (__glibc_unlikely (Z[1] == ZERO))
+  if (__glibc_unlikely (Z[1] == 0))
     {
       for (i = 1; i <= p2; i++)
 	Z[i] = Z[i + 1];
@@ -774,35 +760,35 @@
 __sqr (const mp_no *x, mp_no *y, int p)
 {
   long i, j, k, ip;
-  double u, yk;
+  mantissa_store_t yk;
 
   /* Is z=0?  */
-  if (__glibc_unlikely (X[0] == ZERO))
-    {
-      Y[0] = ZERO;
+  if (__glibc_unlikely (X[0] == 0))
+    {
+      Y[0] = 0;
       return;
     }
 
   /* We need not iterate through all X's since it's pointless to
      multiply zeroes.  */
   for (ip = p; ip > 0; ip--)
-    if (X[ip] != ZERO)
+    if (X[ip] != 0)
       break;
 
   k = (__glibc_unlikely (p < 3)) ? p + p : p + 3;
 
   while (k > 2 * ip + 1)
-    Y[k--] = ZERO;
-
-  yk = ZERO;
+    Y[k--] = 0;
+
+  yk = 0;
 
   while (k > p)
     {
-      double yk2 = 0.0;
+      mantissa_store_t yk2 = 0;
       long lim = k / 2;
 
       if (k % 2 == 0)
-	yk += X[lim] * X[lim];
+	yk += X[lim] * (mantissa_store_t) X[lim];
 
       /* In __mul, this loop (and the one within the next while loop) run
          between a range to calculate the mantissa as follows:
@@ -814,41 +800,35 @@
 	 result.  For cases where the range size is even, the mid-point needs
 	 to be added separately (above).  */
       for (i = k - p, j = p; i < j; i++, j--)
-	yk2 += X[i] * X[j];
-
-      yk += 2.0 * yk2;
-
-      u = (yk + CUTTER) - CUTTER;
-      if (u > yk)
-	u -= RADIX;
-      Y[k--] = yk - u;
-      yk = u * RADIXI;
+	yk2 += X[i] * (mantissa_store_t) X[j];
+
+      yk += 2 * yk2;
+
+      DIV_RADIX (yk, Y[k]);
+      k--;
     }
 
   while (k > 1)
     {
-      double yk2 = 0.0;
+      mantissa_store_t yk2 = 0;
       long lim = k / 2;
 
       if (k % 2 == 0)
-	yk += X[lim] * X[lim];
+	yk += X[lim] * (mantissa_store_t) X[lim];
 
       /* Likewise for this loop.  */
       for (i = 1, j = k - 1; i < j; i++, j--)
-	yk2 += X[i] * X[j];
-
-      yk += 2.0 * yk2;
-
-      u = (yk + CUTTER) - CUTTER;
-      if (u > yk)
-	u -= RADIX;
-      Y[k--] = yk - u;
-      yk = u * RADIXI;
+	yk2 += X[i] * (mantissa_store_t) X[j];
+
+      yk += 2 * yk2;
+
+      DIV_RADIX (yk, Y[k]);
+      k--;
     }
   Y[k] = yk;
 
   /* Squares are always positive.  */
-  Y[0] = 1.0;
+  Y[0] = 1;
 
   /* Get the exponent sum into an intermediate variable.  This is a subtle
      optimization, where given enough registers, all operations on the exponent
@@ -856,7 +836,7 @@
   int e = EX * 2;
 
   /* Is there a carry beyond the most significant digit?  */
-  if (__glibc_unlikely (Y[1] == ZERO))
+  if (__glibc_unlikely (Y[1] == 0))
     {
       for (i = 1; i <= p; i++)
 	Y[i] = Y[i + 1];
@@ -888,7 +868,7 @@
   __cpy (x, &z, p);
   z.e = 0;
   __mp_dbl (&z, &t, p);
-  t = ONE / t;
+  t = 1 / t;
   __dbl_mp (t, y, p);
   EY -= EX;
 
@@ -914,8 +894,8 @@
 {
   mp_no w;
 
-  if (X[0] == ZERO)
-    Z[0] = ZERO;
+  if (X[0] == 0)
+    Z[0] = 0;
   else
     {
       __inv (y, &w, p);

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.h
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.h (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/mpa.h Wed Mar 27 00:01:57 2013
@@ -35,6 +35,7 @@
 /* Common types and definition                                          */
 /************************************************************************/
 
+#include <mpa-arch.h>
 
 /* The mp_no structure holds the details of a multi-precision floating point
    number.
@@ -61,7 +62,7 @@
 typedef struct
 {
   int e;
-  double d[40];
+  mantissa_t d[40];
 } mp_no;
 
 typedef union
@@ -82,9 +83,13 @@
 
 #define ABS(x)   ((x) <  0  ? -(x) : (x))
 
-#define  RADIX     0x1.0p24		/* 2^24    */
-#define  RADIXI    0x1.0p-24		/* 2^-24   */
-#define  CUTTER    0x1.0p76		/* 2^76    */
+#ifndef RADIXI
+# define  RADIXI    0x1.0p-24		/* 2^-24   */
+#endif
+
+#ifndef TWO52
+# define  TWO52     0x1.0p52		/* 2^52    */
+#endif
 
 #define  ZERO      0.0			/* 0       */
 #define  MZERO     -0.0			/* 0 with the sign bit set */
@@ -92,13 +97,13 @@
 #define  MONE      -1.0			/* -1      */
 #define  TWO       2.0			/*  2      */
 
-#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  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  TWO57     0x1.0p57		/* 2^57    */
 #define  TWO71     0x1.0p71		/* 2^71    */
 #define  TWOM1032  0x1.0p-1032		/* 2^-1032 */

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_atan.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_atan.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_atan.c Wed Mar 27 00:01:57 2013
@@ -62,18 +62,9 @@
 #ifndef DLA_FMS
   double t4,t5,t6;
 #endif
-#if 0
-  double y1,y2;
-#endif
   int i,ux,dx;
-#if 0
-  int p;
-#endif
   static const int pr[M]={6,8,10,32};
   number num;
-#if 0
-  mp_no mpt1,mpx,mpy,mpy1,mpy2,mperr;
-#endif
 
   num.d = x;  ux = num.i[HIGH_HALF];  dx = num.i[LOW_HALF];
 

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_sin.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_sin.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_sin.c Wed Mar 27 00:01:57 2013
@@ -100,14 +100,8 @@
 SECTION
 __sin(double x){
 	double xx,res,t,cor,y,s,c,sn,ssn,cs,ccs,xn,a,da,db,eps,xn1,xn2;
-#if 0
-	double w[2];
-#endif
 	mynumber u,v;
 	int4 k,m,n;
-#if 0
-	int4 nn;
-#endif
 	double retval = 0;
 
 	SET_RESTORE_ROUND_53BIT (FE_TONEAREST);
@@ -902,10 +896,6 @@
 bsloww(double x,double dx, double orig,int n) {
   static const double th2_36 = 206158430208.0;   /*    1.5*2**37   */
   double y,x1,x2,xx,r,t,res,cor,w[2];
-#if 0
-  double a,da,xn;
-  union {int4 i[2]; double x;} v;
-#endif
   x1=(x+th2_36)-th2_36;
   y = aa.x*x1*x1*x1;
   r=x+y;

Modified: fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_tan.c
==============================================================================
--- fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_tan.c (original)
+++ fsf/trunk/libc/sysdeps/ieee754/dbl-64/s_tan.c Wed Mar 27 00:01:57 2013
@@ -64,9 +64,6 @@
   int p;
   number num,v;
   mp_no mpa,mpt1,mpt2;
-#if 0
-  mp_no mpy;
-#endif
 
   double retval;
 

Modified: fsf/trunk/libc/sysdeps/powerpc/fpu/s_llround.c
==============================================================================
--- fsf/trunk/libc/sysdeps/powerpc/fpu/s_llround.c (original)
+++ fsf/trunk/libc/sysdeps/powerpc/fpu/s_llround.c Wed Mar 27 00:01:57 2013
@@ -17,6 +17,7 @@
    <http://www.gnu.org/licenses/>.  */
 
 #include <math.h>
+#include <math_ldbl_opt.h>
 
 /* I think that what this routine is supposed to do is round a value
    to the nearest integer, with values exactly on the boundary rounded
@@ -47,3 +48,6 @@
 strong_alias (__llround, __llroundl)
 weak_alias (__llround, llroundl)
 #endif
+#if LONG_DOUBLE_COMPAT (libm, GLIBC_2_1)
+compat_symbol (libm, __llround, llroundl, GLIBC_2_1);
+#endif

Added: fsf/trunk/libc/sysdeps/powerpc/power4/fpu/mpa-arch.h
==============================================================================
--- fsf/trunk/libc/sysdeps/powerpc/power4/fpu/mpa-arch.h (added)
+++ fsf/trunk/libc/sysdeps/powerpc/power4/fpu/mpa-arch.h Wed Mar 27 00:01:57 2013
@@ -1,0 +1,56 @@
+/* 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;									      \
+  })

_______________________________________________
Commits mailing list
Commits@xxxxxxxxxx
http://eglibc.org/cgi-bin/mailman/listinfo/commits