Fixing improper conversion from sin() to sinf() in optimization mode.

When a sin() (cos(), log(), etc.) function is called on a value of
float type and the returned double value is converted to another value
of float type, GCC converts this function call into a float version
(sinf()) in the optimization mode. This avoids two type conversions
and the float version function call usually takes less time. However,
this can result in different result and therefore is unsafe.

For example, the following code produces different results using -O0
(correct), but the same results using -Ox other than -O0 (incorrect).


#include <stdio.h>
#include <math.h>

int main()
{
  float v = 1;
  printf("%.20f\n", (float)sin(v));
  printf("%.20f\n", sinf(v));
}


In this patch, we do this conversion only when the flag
-funsafe-math-optimizations is set. The patch is shown below.


thanks,

Cong
Index: gcc/testsuite/gcc.c-torture/execute/20030125-1.c
===================================================================
--- gcc/testsuite/gcc.c-torture/execute/20030125-1.c	(revision 201891)
+++ gcc/testsuite/gcc.c-torture/execute/20030125-1.c	(working copy)
@@ -44,11 +44,11 @@ __attribute__ ((noinline))
 double
 sin(double a)
 {
-	abort ();
+	return a;
 }
 __attribute__ ((noinline))
 float
 sinf(float a)
 {
-	return a;
+	abort ();
 }
Index: gcc/convert.c
===================================================================
--- gcc/convert.c	(revision 201891)
+++ gcc/convert.c	(working copy)
@@ -99,7 +99,7 @@ convert_to_real (tree type, tree expr)
   /* Disable until we figure out how to decide whether the functions are
      present in runtime.  */
   /* Convert (float)sqrt((double)x) where x is float into sqrtf(x) */
-  if (optimize
+  if (optimize && flag_unsafe_math_optimizations
       && (TYPE_MODE (type) == TYPE_MODE (double_type_node)
           || TYPE_MODE (type) == TYPE_MODE (float_type_node)))
     {

On Tue, 20 Aug 2013, Cong Hou wrote:

> When a sin() (cos(), log(), etc.) function is called on a value of
> float type and the returned double value is converted to another value
> of float type, GCC converts this function call into a float version
> (sinf()) in the optimization mode. This avoids two type conversions
> and the float version function call usually takes less time. However,
> this can result in different result and therefore is unsafe.

Whether it's safe depends on the existence of input values for which the 
double-rounded result is different from the single-rounded result.  It's 
certainly safe for some of the functions for which the optimization is 
enabled, regardless of the precisions involved (fabs, logb).  For sqrt, if 
the smaller precision is p then the larger precision being at least 2p+3 
(I think) makes things same.  For the other cases, exhaustive searches may 
be needed to determine when the conversion is safe.

(Actually, this code appears to deal with cases with up to three types 
involved - the operand, the result and the function that's called.  This 
complicates the analysis a bit, but the same principles apply.)

On Fri, 30 Aug 2013, Cong Hou wrote:

> I have done a simple experiment and found that it may only be safe to do
> this conversion for fabs() and sqrt(). For fabs() it is easy to see the
> safety. For sqrt(), I tried many random floating point values on two
> versions and did not see a difference. Although it cannot prove its safety,
> my proposed patch (shown below) does change the situation.
> 
> However, for other functions, my experiment shows it is unsafe to do this
> conversion. Those functions are:
> 
> sin, cos, sinh, cosh, asin, acos, asinh, acosh, tan, tanh, atan, atanh,
> log, log10, log1p, log2, logb, cbrt, erf, erfc, exp, exp2, expm1.

I don't see why it would be unsafe for logb - can you give an example 
(exact float input value as hex float, and the values you believe logb 
should return for float and double).

For sqrt you appear to have ignored my explanation of the necessary and 
sufficient condition on the precisions of the respective types.  The 
conversion is not safe for sqrt if the two types are double and long 
double and long double is x86 extended, for example.

On Sat, 31 Aug 2013, Cong Hou wrote:

> > I don't see why it would be unsafe for logb - can you give an example
> > (exact float input value as hex float, and the values you believe logb
> > should return for float and double).
> >
> 
> Please try the following code (you will get different results whether to
> use optimization mode):
> 
> #include <math.h>
> #include <stdio.h>
> 
> int main()
> {
>   int i = 0x3edc67d5;
>   float f = *((float*)&i);
>   float r1 = logb(f);
>   float r2 = logbf(f);
>   printf("%x %x\n", *((int*)&r1), *((int*)&r2));
> }

(a) Please stop sending HTML email, so your messages reach the mailing 
list, and resend your messages so far to the list.  The mailing list needs 
to see the whole of both sides of the discussion of any patch being 
proposed for GCC.

(b) I referred to the values *you believe logb should return*.  
Optimization is not meant to preserve library bugs; the comparison should 
be on the basis of correctly rounded results from both float and double 
functions.  The correct return from logb appears to be -2 here, and I get 
that from both logb and logbf with current git glibc.  The existence of a 
bug in some old library is not relevant here.

(c) I always advise writing such tests as *valid C code* using hex floats 
(or if really necessary, unions), not *invalid C code* with aliasing 
violations.