Implement and use libc_feholdexcept_setround_53bit and libc_feupdateenv_53bit

so that double arithmetic in s_sin is done in 53 bit (without extend i386 double precision)

3 files changed, 39 insertions, 5 deletions

diff --git a/sysdeps/generic/math_private.h b/sysdeps/generic/math_private.h
index 777762dd33..be1e4d2314 100644
--- a/sysdeps/generic/math_private.h
+++ b/sysdeps/generic/math_private.h
@@ -370,6 +370,9 @@ extern void __docos (double __x, double __dx, double __v[]);
 #define libc_feholdexcept_setroundl(e, r) \
   do { feholdexcept (e); fesetround (r); } while (0)
 
+#define libc_feholdexcept_setround_53bit(e, r) \
+  libc_feholdexcept_setround (e, r)
+
 #define libc_fetestexcept(e) fetestexcept (e)
 #define libc_fetestexceptf(e) fetestexcept (e)
 #define libc_fetestexceptl(e) fetestexcept (e)
@@ -382,6 +385,8 @@ extern void __docos (double __x, double __dx, double __v[]);
 #define libc_feupdateenvf(e) (void) feupdateenv (e)
 #define libc_feupdateenvl(e) (void) feupdateenv (e)
 
+#define libc_feupdateenv_53bit(e) libc_feupdateenv (e)
+
 #define __nan(str) \
   (__builtin_constant_p (str) && str[0] == '\0' ? NAN : __nan (str))
 #define __nanf(str) \
diff --git a/sysdeps/i386/fpu/math_private.h b/sysdeps/i386/fpu/math_private.h
index 5253998a57..d96996fadf 100644
--- a/sysdeps/i386/fpu/math_private.h
+++ b/sysdeps/i386/fpu/math_private.h
@@ -16,4 +16,33 @@ do							\
 while (0)
 
 #include_next <math_private.h>
+
+# include <fpu_control.h>
+
+# undef libc_feholdexcept_setround_53bit
+# define libc_feholdexcept_setround_53bit(e, r)	\
+  do						\
+    {						\
+      fpu_control_t cw;				\
+      libc_feholdexcept_setround (e, r);	\
+      _FPU_GETCW (cw);				\
+      cw &= ~(fpu_control_t) _FPU_EXTENDED;	\
+      cw |= _FPU_DOUBLE;			\
+      _FPU_SETCW (cw);				\
+    }						\
+  while (0)
+
+# undef libc_feupdateenv_53bit
+# define libc_feupdateenv_53bit(e)		\
+  do						\
+    {						\
+      fpu_control_t cw;				\
+      libc_feupdateenv (e);			\
+      _FPU_GETCW (cw);				\
+      cw &= ~(fpu_control_t) _FPU_EXTENDED;	\
+      cw |= _FPU_EXTENDED;			\
+      _FPU_SETCW (cw);				\
+    }						\
+  while (0)
+
 #endif
diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index e3e3a2a965..4b4b67573d 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -1,7 +1,7 @@
 /*
  * IBM Accurate Mathematical Library
  * written by International Business Machines Corp.
- * Copyright (C) 2001, 2009, 2011 Free Software Foundation
+ * Copyright (C) 2001-2012 Free Software Foundation
  *
  * 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
@@ -111,7 +111,7 @@ __sin(double x){
 	fenv_t env;
 	double retval = 0;
 
-	libc_feholdexcept_setround (&env, FE_TONEAREST);
+	libc_feholdexcept_setround_53bit (&env, FE_TONEAREST);
 
 	u.x = x;
 	m = u.i[HIGH_HALF];
@@ -365,7 +365,7 @@ __sin(double x){
 	}
 
  ret:
-	libc_feupdateenv (&env);
+	libc_feupdateenv_53bit (&env);
 	return retval;
 }
 
@@ -386,7 +386,7 @@ __cos(double x)
   fenv_t env;
   double retval = 0;
 
-  libc_feholdexcept_setround (&env, FE_TONEAREST);
+  libc_feholdexcept_setround_53bit (&env, FE_TONEAREST);
 
   u.x = x;
   m = u.i[HIGH_HALF];
@@ -635,7 +635,7 @@ __cos(double x)
   }
 
  ret:
-  libc_feupdateenv (&env);
+  libc_feupdateenv_53bit (&env);
   return retval;
 }