This is the mail archive of the libc-alpha@sourceware.org mailing list for the glibc project.

Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

Re: Fix sin, cos, tan in non-default rounding modes (bug 3976)

On 03/02/2012 06:39 PM, Joseph S. Myers wrote: 
Similar to the patch for exp, this patch fixes problems with sin, cos
and tan producing wild results in non-default rounding modes by making
them save and restore the rounding mode and do computations in
round-to-nearest.

Tested on x86_64 and x86 and ULPs updated based on those tests.

2012-03-02 Joseph Myers<joseph@codesourcery.com>

	[BZ #3976]
	* sysdeps/ieee754/dbl-64/s_sin.c: Include<fenv.h>
	(__sin): Save and restore rounding mode and use round-to-nearest
	for all computations.
	(__cos): Save and restore rounding mode and use round-to-nearest
	for all computations.
	* sysdeps/ieee754/dbl-64/s_tan.c: Include "math_private.h" and
	<fenv.h>.
	(tan): Save and restore rounding mode and use round-to-nearest for
	all computations.
	* math/libm-test.inc (cos_test_tonearest): New function.
	(cos_test_towardzero): Likewise.
	(cos_test_downward): Likewise.
	(cos_test_upward): Likewise.
	(sin_test_tonearest): Likewise.
	(sin_test_towardzero): Likewise.
	(sin_test_downward): Likewise.
	(sin_test_upward): Likewise.
	(tan_test_tonearest): Likewise.
	(tan_test_towardzero): Likewise.
	(tan_test_downward): Likewise.
	(tan_test_upward): Likewise.
	(main): Call the new functions.
	* sysdeps/i386/fpu/libm-test-ulps: Update.
	* sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
> [...]

thanks, this looks fine, I just have a few suggestions for minor nits below.



diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index 5b79854..32ba66d 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -53,6 +53,7 @@
  #include "usncs.h"
  #include "MathLib.h"
  #include "math_private.h"
+#include<fenv.h>

  #ifndef SECTION
  # define SECTION
@@ -107,12 +108,16 @@ __sin(double x){
  #if 0
  	int4 nn;
  #endif
+	fenv_t env;
+	double retval = 0;
+
+	libc_feholdexcept_setround (&env, FE_TONEAREST);

  	u.x = x;
  	m = u.i[HIGH_HALF];
  	k = 0x7fffffff&m;              /* no sign           */
  	if (k<  0x3e500000)            /* if x->0 =>sin(x)=x */
-	 return x;
+	  { retval = x; goto ret; }
   /*---------------------------- 2^-26<  |x|<  0.25 ----------------------*/
  	else  if (k<  0x3fd00000){
  	  xx = x*x;
@@ -120,7 +125,8 @@ __sin(double x){
  	  t = ((((s5.x*xx + s4.x)*xx + s3.x)*xx + s2.x)*xx + s1.x)*(xx*x);
  	  res = x+t;
  	  cor = (x-res)+t;
-	  return (res == res + 1.07*cor)? res : slow(x);
+	  retval = (res == res + 1.07*cor)? res : slow(x);
+	  goto ret;
  	}    /*  else  if (k<  0x3fd00000)    */
  /*---------------------------- 0.25<|x|<  0.855469---------------------- */
  	else if (k<  0x3feb6000)  {
@@ -137,7 +143,8 @@ __sin(double x){
  	  cor=(ssn+s*ccs-sn*c)+cs*s;
  	  res=sn+cor;
  	  cor=(sn-res)+cor;
-	  return (res==res+1.096*cor)? res : slow1(x);
+	  retval = (res==res+1.096*cor)? res : slow1(x);
+	  goto ret;
  	}    /*   else  if (k<  0x3feb6000)    */

  /*----------------------- 0.855469<|x|<2.426265  ----------------------*/
@@ -163,7 +170,8 @@ __sin(double x){
  	  cor=(ccs-s*ssn-cs*c)-sn*s;
  	  res=cs+cor;
  	  cor=(cs-res)+cor;
-	  return (res==res+1.020*cor)? ((m>0)?res:-res) : slow2(x);
+	  retval = (res==res+1.020*cor)? ((m>0)?res:-res) : slow2(x);
+	  goto ret;
  	} /*   else  if (k<  0x400368fd)    */

  /*-------------------------- 2.426265<|x|<  105414350 ----------------------*/
@@ -189,7 +197,8 @@ __sin(double x){
  	      res = a+t;
  	      cor = (a-res)+t;
  	      cor = (cor>0)? 1.02*cor+eps : 1.02*cor -eps;
-	      return (res == res + cor)? res : sloww(a,da,x);
+	      retval = (res == res + cor)? res : sloww(a,da,x);
+	      goto ret;
  	    }
  	    else  {
  	      if (a>0)
@@ -210,7 +219,8 @@ __sin(double x){
  	      res=sn+cor;
  	      cor=(sn-res)+cor;
  	      cor = (cor>0)? 1.035*cor+eps : 1.035*cor-eps;
-	      return (res==res+cor)? ((m)?res:-res) : sloww1(a,da,x);
+	      retval = (res==res+cor)? ((m)?res:-res) : sloww1(a,da,x);
+	      goto ret;
  	    }
  	    break;

@@ -232,7 +242,8 @@ __sin(double x){
  	    res=cs+cor;
  	    cor=(cs-res)+cor;
  	    cor = (cor>0)? 1.025*cor+eps : 1.025*cor-eps;
-	    return (res==res+cor)? ((n&2)?-res:res) : sloww2(a,da,x,n);
+	    retval = (res==res+cor)? ((n&2)?-res:res) : sloww2(a,da,x,n);
+	    goto ret;

break;

@@ -268,7 +279,8 @@ __sin(double x){
  	      res = a+t;
  	      cor = (a-res)+t;
  	      cor = (cor>0)? 1.02*cor+eps : 1.02*cor -eps;
-	      return (res == res + cor)? res : bsloww(a,da,x,n);
+	      retval = (res == res + cor)? res : bsloww(a,da,x,n);
+	      goto ret;
  	    }
  	    else  {
  	      if (a>0) {m=1;t=a;db=da;}
@@ -287,7 +299,8 @@ __sin(double x){
  	      res=sn+cor;
  	      cor=(sn-res)+cor;
  	      cor = (cor>0)? 1.035*cor+eps : 1.035*cor-eps;
-	      return (res==res+cor)? ((m)?res:-res) : bsloww1(a,da,x,n);
+	      retval = (res==res+cor)? ((m)?res:-res) : bsloww1(a,da,x,n);
+	      goto ret;
  		   }
  	    break;

@@ -309,7 +322,8 @@ __sin(double x){
  	    res=cs+cor;
  	    cor=(cs-res)+cor;
  	    cor = (cor>0)? 1.025*cor+eps : 1.025*cor-eps;
-	    return (res==res+cor)? ((n&2)?-res:res) : bsloww2(a,da,x,n);
+	    retval = (res==res+cor)? ((n&2)?-res:res) : bsloww2(a,da,x,n);
+	    goto ret;

break;

@@ -323,17 +337,20 @@ __sin(double x){
  	  n = __branred(x,&a,&da);
  	  switch (n) {
  	  case 0:
-	    if (a*a<  0.01588) return bsloww(a,da,x,n);
-	    else return bsloww1(a,da,x,n);
+	    if (a*a<  0.01588) retval = bsloww(a,da,x,n);
+	    else retval = bsloww1(a,da,x,n);
+	    goto ret;
  	    break;
  	  case 2:
-	    if (a*a<  0.01588) return bsloww(-a,-da,x,n);
-	    else return bsloww1(-a,-da,x,n);
+	    if (a*a<  0.01588) retval = bsloww(-a,-da,x,n);
+	    else retval = bsloww1(-a,-da,x,n);
+	    goto ret;
  	    break;

  	  case 1:
  	  case 3:
-	    return  bsloww2(a,da,x,n);
+	    retval = bsloww2(a,da,x,n);
+	    goto ret;
  	    break;
  	  }

@@ -343,9 +360,13 @@ __sin(double x){
  	else {
  	  if (k == 0x7ff00000&&  u.i[LOW_HALF] == 0)
  	    __set_errno (EDOM);
-	  return x / x;
+	  retval = x / x;
+	  goto ret;

I suggest to remove the goto ret.


  	}
-	return 0;         /* unreachable */
+
+ ret:
+	libc_feupdateenv (&env);
+	return retval;
  }


@@ -362,11 +383,16 @@ __cos(double x)
    mynumber u,v;
    int4 k,m,n;

+  fenv_t env;
+  double retval = 0;
+
+  libc_feholdexcept_setround (&env, FE_TONEAREST);
+
    u.x = x;
    m = u.i[HIGH_HALF];
    k = 0x7fffffff&m;

-  if (k<  0x3e400000 ) return 1.0; /* |x|<2^-27 =>  cos(x)=1 */
+  if (k<  0x3e400000 ) { retval = 1.0; goto ret; } /* |x|<2^-27 =>  cos(x)=1 */

    else if (k<  0x3feb6000 ) {/* 2^-27<  |x|<  0.855469 */
      y=ABS(x);
@@ -383,7 +409,8 @@ __cos(double x)
      cor=(ccs-s*ssn-cs*c)-sn*s;
      res=cs+cor;
      cor=(cs-res)+cor;
-    return (res==res+1.020*cor)? res : cslow2(x);
+    retval = (res==res+1.020*cor)? res : cslow2(x);
+    goto ret;

} /* else if (k< 0x3feb6000) */

@@ -397,7 +424,8 @@ __cos(double x)
        res = a+t;
        cor = (a-res)+t;
        cor = (cor>0)? 1.02*cor+1.0e-31 : 1.02*cor -1.0e-31;
-      return (res == res + cor)? res : csloww(a,da,x);
+      retval = (res == res + cor)? res : csloww(a,da,x);
+      goto ret;
      }
      else  {
        if (a>0) {m=1;t=a;db=da;}
@@ -416,7 +444,8 @@ __cos(double x)
        res=sn+cor;
        cor=(sn-res)+cor;
        cor = (cor>0)? 1.035*cor+1.0e-31 : 1.035*cor-1.0e-31;
-      return (res==res+cor)? ((m)?res:-res) : csloww1(a,da,x);
+      retval = (res==res+cor)? ((m)?res:-res) : csloww1(a,da,x);
+      goto ret;
  }

  }    /*   else  if (k<  0x400368fd)    */
@@ -443,7 +472,8 @@ __cos(double x)
  	res = a+t;
  	cor = (a-res)+t;
  	cor = (cor>0)? 1.02*cor+eps : 1.02*cor -eps;
-	return (res == res + cor)? res : csloww(a,da,x);
+	retval = (res == res + cor)? res : csloww(a,da,x);
+	goto ret;
        }
        else  {
  	if (a>0) {m=1;t=a;db=da;}
@@ -462,7 +492,8 @@ __cos(double x)
  	res=sn+cor;
  	cor=(sn-res)+cor;
  	cor = (cor>0)? 1.035*cor+eps : 1.035*cor-eps;
-	return (res==res+cor)? ((m)?res:-res) : csloww1(a,da,x);
+	retval = (res==res+cor)? ((m)?res:-res) : csloww1(a,da,x);
+	goto ret;
        }
        break;

@@ -483,7 +514,8 @@ __cos(double x)
        res=cs+cor;
        cor=(cs-res)+cor;
        cor = (cor>0)? 1.025*cor+eps : 1.025*cor-eps;
-      return (res==res+cor)? ((n)?-res:res) : csloww2(a,da,x,n);
+      retval = (res==res+cor)? ((n)?-res:res) : csloww2(a,da,x,n);
+      goto ret;

break;

@@ -518,7 +550,8 @@ __cos(double x)
  	res = a+t;
  	cor = (a-res)+t;
  	cor = (cor>0)? 1.02*cor+eps : 1.02*cor -eps;
-	return (res == res + cor)? res : bsloww(a,da,x,n);
+	retval = (res == res + cor)? res : bsloww(a,da,x,n);
+	goto ret;
        }
        else  {
  	if (a>0) {m=1;t=a;db=da;}
@@ -537,7 +570,8 @@ __cos(double x)
  	res=sn+cor;
  	cor=(sn-res)+cor;
  	cor = (cor>0)? 1.035*cor+eps : 1.035*cor-eps;
-	return (res==res+cor)? ((m)?res:-res) : bsloww1(a,da,x,n);
+	retval = (res==res+cor)? ((m)?res:-res) : bsloww1(a,da,x,n);
+	goto ret;
        }
        break;

@@ -558,7 +592,8 @@ __cos(double x)
        res=cs+cor;
        cor=(cs-res)+cor;
        cor = (cor>0)? 1.025*cor+eps : 1.025*cor-eps;
-      return (res==res+cor)? ((n)?-res:res) : bsloww2(a,da,x,n);
+      retval = (res==res+cor)? ((n)?-res:res) : bsloww2(a,da,x,n);
+      goto ret;
        break;

      }
@@ -570,17 +605,20 @@ __cos(double x)
      n = __branred(x,&a,&da);
      switch (n) {
      case 1:
-      if (a*a<  0.01588) return bsloww(-a,-da,x,n);
-      else return bsloww1(-a,-da,x,n);
+      if (a*a<  0.01588) retval = bsloww(-a,-da,x,n);
+      else retval = bsloww1(-a,-da,x,n);
+      goto ret;
        break;
  		case 3:
-		  if (a*a<  0.01588) return bsloww(a,da,x,n);
-		  else return bsloww1(a,da,x,n);
+		  if (a*a<  0.01588) retval = bsloww(a,da,x,n);
+		  else retval = bsloww1(a,da,x,n);
+		  goto ret;
  		  break;

      case 0:
      case 2:
-      return  bsloww2(a,da,x,n);
+      retval = bsloww2(a,da,x,n);
+      goto ret;
        break;
      }

@@ -592,10 +630,13 @@ __cos(double x)
    else {
      if (k == 0x7ff00000&&  u.i[LOW_HALF] == 0)
        __set_errno (EDOM);
-    return x / x; /* |x|>  2^1024 */
+    retval = x / x; /* |x|>  2^1024 */
+    goto ret;


I suggest to omit the goto ret.

    }
-  return 0;

+ ret:
+  libc_feupdateenv (&env);
+  return retval;
  }

  /************************************************************************/
diff --git a/sysdeps/ieee754/dbl-64/s_tan.c b/sysdeps/ieee754/dbl-64/s_tan.c
index 962a4eb..2c26756 100644
--- a/sysdeps/ieee754/dbl-64/s_tan.c
+++ b/sysdeps/ieee754/dbl-64/s_tan.c
@@ -39,6 +39,8 @@
  #include "mpa.h"
  #include "MathLib.h"
  #include "math.h"
+#include "math_private.h"
+#include<fenv.h>

  #ifndef SECTION
  # define SECTION
@@ -66,21 +68,27 @@ tan(double x) {
    mp_no mpy;
  #endif

+  fenv_t env;
+  double retval;
+
    int __branred(double, double *, double *);
    int __mpranred(double, mp_no *, int);

+  libc_feholdexcept_setround (&env, FE_TONEAREST);
+
    /* x=+-INF, x=NaN */
    num.d = x;  ux = num.i[HIGH_HALF];
    if ((ux&0x7ff00000)==0x7ff00000) {
      if ((ux&0x7fffffff)==0x7ff00000)
        __set_errno (EDOM);
-    return x-x;
+    retval = x-x;
+    goto ret;
    }

w=(x<ZERO) ? -x : x;

    /* (I) The case abs(x)<= 1.259e-8 */
-  if (w<=g1.d)  return x;
+  if (w<=g1.d) { retval = x; goto ret; }

    /* (II) The case 1.259e-8<  abs(x)<= 0.0608 */
    if (w<=g2.d) {
@@ -88,7 +96,7 @@ tan(double x) {
      /* First stage */
      x2 = x*x;
      t2 = x*x2*(d3.d+x2*(d5.d+x2*(d7.d+x2*(d9.d+x2*d11.d))));
-    if ((y=x+(t2-u1.d*t2)) == x+(t2+u1.d*t2))  return y;
+    if ((y=x+(t2-u1.d*t2)) == x+(t2+u1.d*t2)) { retval = y; goto ret; }

      /* Second stage */
      c1 = x2*(a15.d+x2*(a17.d+x2*(a19.d+x2*(a21.d+x2*(a23.d+x2*(a25.d+
@@ -108,8 +116,9 @@ tan(double x) {
      MUL2(x2,xx2,c2,cc2,c1,cc1,t1,t2,t3,t4,t5,t6,t7,t8)
      MUL2(x ,zero.d,c1,cc1,c2,cc2,t1,t2,t3,t4,t5,t6,t7,t8)
      ADD2(x    ,zero.d,c2,cc2,c1,cc1,t1,t2)
-    if ((y=c1+(cc1-u2.d*c1)) == c1+(cc1+u2.d*c1))  return y;
-    return tanMp(x);
+    if ((y=c1+(cc1-u2.d*c1)) == c1+(cc1+u2.d*c1)) { retval = y; goto ret; }
+    retval = tanMp(x);
+    goto ret;
    }

    /* (III) The case 0.0608<  abs(x)<= 0.787 */
@@ -120,10 +129,10 @@ tan(double x) {
      z = w-xfg[i][0].d;  z2 = z*z;   s = (x<ZERO) ? MONE : ONE;
      pz = z+z*z2*(e0.d+z2*e1.d);
      fi = xfg[i][1].d;   gi = xfg[i][2].d;   t2 = pz*(gi+fi)/(gi-pz);
-    if ((y=fi+(t2-fi*u3.d))==fi+(t2+fi*u3.d))  return (s*y);
+    if ((y=fi+(t2-fi*u3.d))==fi+(t2+fi*u3.d)) { retval = (s*y); goto ret; }
      t3 = (t2<ZERO) ? -t2 : t2;
      t4 = fi*ua3.d+t3*ub3.d;
-    if ((y=fi+(t2-t4))==fi+(t2+t4))  return (s*y);
+    if ((y=fi+(t2-t4))==fi+(t2+t4)) { retval = (s*y); goto ret; }

      /* Second stage */
      ffi = xfg[i][3].d;
@@ -141,8 +150,9 @@ tan(double x) {
      SUB2(one.d,zero.d,c3,cc3,c1,cc1,t1,t2)
      DIV2(c2,cc2,c1,cc1,c3,cc3,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)

-    if ((y=c3+(cc3-u4.d*c3))==c3+(cc3+u4.d*c3))  return (s*y);
-    return tanMp(x);
+      if ((y=c3+(cc3-u4.d*c3))==c3+(cc3+u4.d*c3)) { retval = (s*y); goto ret; }
+    retval = tanMp(x);
+    goto ret;
    }

    /* (---) The case 0.787<  abs(x)<= 25 */
@@ -160,7 +170,7 @@ tan(double x) {
      else         {ya= a;  yya= da;  sy= ONE;}

      /* (IV),(V) The case 0.787<  abs(x)<= 25,    abs(y)<= 1e-7 */
-    if (ya<=gy1.d)  return tanMp(x);
+    if (ya<=gy1.d) { retval = tanMp(x); goto ret; }

      /* (VI) The case 0.787<  abs(x)<= 25,    1e-7<  abs(y)<= 0.0608 */
      if (ya<=gy2.d) {
@@ -170,10 +180,10 @@ tan(double x) {
  	/* First stage -cot */
  	EADD(a,t2,b,db)
  	DIV2(one.d,zero.d,b,db,c,dc,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-	if ((y=c+(dc-u6.d*c))==c+(dc+u6.d*c))  return (-y); }
+	if ((y=c+(dc-u6.d*c))==c+(dc+u6.d*c)) { retval = (-y); goto ret; } }
        else {
  	/* First stage tan */
-	if ((y=a+(t2-u5.d*a))==a+(t2+u5.d*a))  return y; }
+	if ((y=a+(t2-u5.d*a))==a+(t2+u5.d*a)) { retval = y; goto ret; } }
        /* Second stage */
        /* Range reduction by algorithm ii */
        t = (x*hpinv.d + toint.d);
@@ -211,11 +221,12 @@ tan(double x) {
        if (n) {
  	/* Second stage -cot */
  	DIV2(one.d,zero.d,c1,cc1,c2,cc2,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-	if ((y=c2+(cc2-u8.d*c2)) == c2+(cc2+u8.d*c2))  return (-y); }
+	if ((y=c2+(cc2-u8.d*c2)) == c2+(cc2+u8.d*c2)) { retval = (-y); goto ret; } }
        else {
  	/* Second stage tan */
-	if ((y=c1+(cc1-u7.d*c1)) == c1+(cc1+u7.d*c1))  return y; }
-      return tanMp(x);
+	if ((y=c1+(cc1-u7.d*c1)) == c1+(cc1+u7.d*c1)) { retval = y; goto ret; } }
+      retval = tanMp(x);
+      goto ret;
      }

      /* (VII) The case 0.787<  abs(x)<= 25,    0.0608<  abs(y)<= 0.787 */
@@ -229,17 +240,17 @@ tan(double x) {
      if (n) {
        /* -cot */
        t2 = pz*(fi+gi)/(fi+pz);
-      if ((y=gi-(t2-gi*u10.d))==gi-(t2+gi*u10.d))  return (-sy*y);
+      if ((y=gi-(t2-gi*u10.d))==gi-(t2+gi*u10.d)) { retval = (-sy*y); goto ret; }
        t3 = (t2<ZERO) ? -t2 : t2;
        t4 = gi*ua10.d+t3*ub10.d;
-      if ((y=gi-(t2-t4))==gi-(t2+t4))  return (-sy*y); }
+      if ((y=gi-(t2-t4))==gi-(t2+t4)) { retval = (-sy*y); goto ret; } }
      else   {
        /* tan */
        t2 = pz*(gi+fi)/(gi-pz);
-      if ((y=fi+(t2-fi*u9.d))==fi+(t2+fi*u9.d))  return (sy*y);
+      if ((y=fi+(t2-fi*u9.d))==fi+(t2+fi*u9.d)) { retval = (sy*y); goto ret; }
        t3 = (t2<ZERO) ? -t2 : t2;
        t4 = fi*ua9.d+t3*ub9.d;
-      if ((y=fi+(t2-t4))==fi+(t2+t4))  return (sy*y); }
+      if ((y=fi+(t2-t4))==fi+(t2+t4)) { retval = (sy*y); goto ret; } }

      /* Second stage */
      ffi = xfg[i][3].d;
@@ -260,13 +271,14 @@ tan(double x) {
      if (n) {
        /* -cot */
        DIV2(c1,cc1,c2,cc2,c3,cc3,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-      if ((y=c3+(cc3-u12.d*c3))==c3+(cc3+u12.d*c3))  return (-sy*y); }
+      if ((y=c3+(cc3-u12.d*c3))==c3+(cc3+u12.d*c3)) { retval = (-sy*y); goto ret; } }
      else {
        /* tan */
        DIV2(c2,cc2,c1,cc1,c3,cc3,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-      if ((y=c3+(cc3-u11.d*c3))==c3+(cc3+u11.d*c3))  return (sy*y); }
+      if ((y=c3+(cc3-u11.d*c3))==c3+(cc3+u11.d*c3)) { retval = (sy*y); goto ret; } }

-    return tanMp(x);
+    retval = tanMp(x);
+    goto ret;
    }

    /* (---) The case 25<  abs(x)<= 1e8 */
@@ -288,7 +300,7 @@ tan(double x) {
      else         {ya= a;  yya= da;  sy= ONE;}

      /* (+++) The case 25<  abs(x)<= 1e8,    abs(y)<= 1e-7 */
-    if (ya<=gy1.d)  return tanMp(x);
+    if (ya<=gy1.d) { retval = tanMp(x); goto ret; }

      /* (VIII) The case 25<  abs(x)<= 1e8,    1e-7<  abs(y)<= 0.0608 */
      if (ya<=gy2.d) {
@@ -298,10 +310,10 @@ tan(double x) {
  	/* First stage -cot */
  	EADD(a,t2,b,db)
  	DIV2(one.d,zero.d,b,db,c,dc,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-	if ((y=c+(dc-u14.d*c))==c+(dc+u14.d*c))  return (-y); }
+	if ((y=c+(dc-u14.d*c))==c+(dc+u14.d*c)) { retval = (-y); goto ret; } }
        else {
  	/* First stage tan */
-	if ((y=a+(t2-u13.d*a))==a+(t2+u13.d*a))  return y; }
+	if ((y=a+(t2-u13.d*a))==a+(t2+u13.d*a)) { retval = y; goto ret; } }

        /* Second stage */
        MUL2(a,da,a,da,x2,xx2,t1,t2,t3,t4,t5,t6,t7,t8)
@@ -325,11 +337,12 @@ tan(double x) {
        if (n) {
  	/* Second stage -cot */
  	DIV2(one.d,zero.d,c1,cc1,c2,cc2,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-	if ((y=c2+(cc2-u16.d*c2)) == c2+(cc2+u16.d*c2))  return (-y); }
+	if ((y=c2+(cc2-u16.d*c2)) == c2+(cc2+u16.d*c2)) { retval = (-y); goto ret; } }
        else {
  	/* Second stage tan */
-	if ((y=c1+(cc1-u15.d*c1)) == c1+(cc1+u15.d*c1))  return (y); }
-      return tanMp(x);
+	if ((y=c1+(cc1-u15.d*c1)) == c1+(cc1+u15.d*c1)) { retval = (y); goto ret; } }
+      retval = tanMp(x);
+      goto ret;
      }

      /* (IX) The case 25<  abs(x)<= 1e8,    0.0608<  abs(y)<= 0.787 */
@@ -342,17 +355,17 @@ tan(double x) {
      if (n) {
        /* -cot */
        t2 = pz*(fi+gi)/(fi+pz);
-      if ((y=gi-(t2-gi*u18.d))==gi-(t2+gi*u18.d))  return (-sy*y);
+      if ((y=gi-(t2-gi*u18.d))==gi-(t2+gi*u18.d)) { retval = (-sy*y); goto ret; }
        t3 = (t2<ZERO) ? -t2 : t2;
        t4 = gi*ua18.d+t3*ub18.d;
-      if ((y=gi-(t2-t4))==gi-(t2+t4))  return (-sy*y); }
+      if ((y=gi-(t2-t4))==gi-(t2+t4)) { retval = (-sy*y); goto ret; } }
      else   {
        /* tan */
        t2 = pz*(gi+fi)/(gi-pz);
-      if ((y=fi+(t2-fi*u17.d))==fi+(t2+fi*u17.d))  return (sy*y);
+      if ((y=fi+(t2-fi*u17.d))==fi+(t2+fi*u17.d)) { retval = (sy*y); goto ret; }
        t3 = (t2<ZERO) ? -t2 : t2;
        t4 = fi*ua17.d+t3*ub17.d;
-      if ((y=fi+(t2-t4))==fi+(t2+t4))  return (sy*y); }
+      if ((y=fi+(t2-t4))==fi+(t2+t4)) { retval = (sy*y); goto ret; } }

      /* Second stage */
      ffi = xfg[i][3].d;
@@ -373,12 +386,13 @@ tan(double x) {
      if (n) {
        /* -cot */
        DIV2(c1,cc1,c2,cc2,c3,cc3,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-      if ((y=c3+(cc3-u20.d*c3))==c3+(cc3+u20.d*c3))  return (-sy*y); }
+      if ((y=c3+(cc3-u20.d*c3))==c3+(cc3+u20.d*c3)) { retval = (-sy*y); goto ret; } }
      else {
        /* tan */
        DIV2(c2,cc2,c1,cc1,c3,cc3,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-      if ((y=c3+(cc3-u19.d*c3))==c3+(cc3+u19.d*c3))  return (sy*y); }
-    return tanMp(x);
+      if ((y=c3+(cc3-u19.d*c3))==c3+(cc3+u19.d*c3)) { retval = (sy*y); goto ret; } }
+    retval = tanMp(x);
+    goto ret;
    }

    /* (---) The case 1e8<  abs(x)<  2**1024 */
@@ -389,7 +403,7 @@ tan(double x) {
    else         {ya= a;  yya= da;  sy= ONE;}

    /* (+++) The case 1e8<  abs(x)<  2**1024,    abs(y)<= 1e-7 */
-  if (ya<=gy1.d)  return tanMp(x);
+  if (ya<=gy1.d) { retval = tanMp(x); goto ret; }

    /* (X) The case 1e8<  abs(x)<  2**1024,    1e-7<  abs(y)<= 0.0608 */
    if (ya<=gy2.d) {
@@ -399,10 +413,10 @@ tan(double x) {
        /* First stage -cot */
        EADD(a,t2,b,db)
        DIV2(one.d,zero.d,b,db,c,dc,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-      if ((y=c+(dc-u22.d*c))==c+(dc+u22.d*c))  return (-y); }
+      if ((y=c+(dc-u22.d*c))==c+(dc+u22.d*c)) { retval = (-y); goto ret; } }
      else {
        /* First stage tan */
-      if ((y=a+(t2-u21.d*a))==a+(t2+u21.d*a))  return y; }
+      if ((y=a+(t2-u21.d*a))==a+(t2+u21.d*a)) { retval = y; goto ret; } }

      /* Second stage */
      /* Reduction by algorithm iv */
@@ -431,11 +445,12 @@ tan(double x) {
      if (n) {
        /* Second stage -cot */
        DIV2(one.d,zero.d,c1,cc1,c2,cc2,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-      if ((y=c2+(cc2-u24.d*c2)) == c2+(cc2+u24.d*c2))  return (-y); }
+      if ((y=c2+(cc2-u24.d*c2)) == c2+(cc2+u24.d*c2)) { retval = (-y); goto ret; } }
      else {
        /* Second stage tan */
-      if ((y=c1+(cc1-u23.d*c1)) == c1+(cc1+u23.d*c1))  return y; }
-    return tanMp(x);
+      if ((y=c1+(cc1-u23.d*c1)) == c1+(cc1+u23.d*c1)) { retval = y; goto ret; } }
+    retval = tanMp(x);
+    goto ret;
    }

    /* (XI) The case 1e8<  abs(x)<  2**1024,    0.0608<  abs(y)<= 0.787 */
@@ -448,17 +463,17 @@ tan(double x) {
    if (n) {
      /* -cot */
      t2 = pz*(fi+gi)/(fi+pz);
-    if ((y=gi-(t2-gi*u26.d))==gi-(t2+gi*u26.d))  return (-sy*y);
+    if ((y=gi-(t2-gi*u26.d))==gi-(t2+gi*u26.d)) { retval = (-sy*y); goto ret; }
      t3 = (t2<ZERO) ? -t2 : t2;
      t4 = gi*ua26.d+t3*ub26.d;
-    if ((y=gi-(t2-t4))==gi-(t2+t4))  return (-sy*y); }
+    if ((y=gi-(t2-t4))==gi-(t2+t4)) { retval = (-sy*y); goto ret; } }
    else   {
      /* tan */
      t2 = pz*(gi+fi)/(gi-pz);
-    if ((y=fi+(t2-fi*u25.d))==fi+(t2+fi*u25.d))  return (sy*y);
+    if ((y=fi+(t2-fi*u25.d))==fi+(t2+fi*u25.d)) { retval = (sy*y); goto ret; }
      t3 = (t2<ZERO) ? -t2 : t2;
      t4 = fi*ua25.d+t3*ub25.d;
-    if ((y=fi+(t2-t4))==fi+(t2+t4))  return (sy*y); }
+    if ((y=fi+(t2-t4))==fi+(t2+t4)) { retval = (sy*y); goto ret; } }

    /* Second stage */
    ffi = xfg[i][3].d;
@@ -479,14 +494,18 @@ tan(double x) {
    if (n) {
      /* -cot */
      DIV2(c1,cc1,c2,cc2,c3,cc3,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-    if ((y=c3+(cc3-u28.d*c3))==c3+(cc3+u28.d*c3))  return (-sy*y); }
+    if ((y=c3+(cc3-u28.d*c3))==c3+(cc3+u28.d*c3)) { retval = (-sy*y); goto ret; } }
    else {
      /* tan */
      DIV2(c2,cc2,c1,cc1,c3,cc3,t1,t2,t3,t4,t5,t6,t7,t8,t9,t10)
-    if ((y=c3+(cc3-u27.d*c3))==c3+(cc3+u27.d*c3))  return (sy*y); }
-  return tanMp(x);
-}
+    if ((y=c3+(cc3-u27.d*c3))==c3+(cc3+u27.d*c3)) { retval = (sy*y); goto ret; } }
+  retval = tanMp(x);
+  goto ret;

This could be if ... else without the goto ret. Or the last goto ret could be removed (it doesn't harm but ret is the next instruction)

+ ret:
+  libc_feupdateenv (&env);
+  return retval;
+}


Thanks for the patch, please commit after considering my suggestions,
Andreas
--
 Andreas Jaeger aj@{suse.com,opensuse.org} Twitter/Identica: jaegerandi
  SUSE LINUX Products GmbH, Maxfeldstr. 5, 90409 Nürnberg, Germany
   GF: Jeff Hawn,Jennifer Guild,Felix Imendörffer,HRB16746 (AG Nürnberg)
    GPG fingerprint = 93A3 365E CE47 B889 DF7F  FED1 389A 563C C272 A126

Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]