Fix convert_to_real bugs (PR 36578 etc.)

["Joseph S. Myers" <joseph at codesourcery dot com>]

convert_to_real handles conversion to TYPE of arithmetic on ITYPE.
ITYPE may in turn result from extensions to ARG0 and ARG1 being in
narrower types.  I.e., (type)((itype)arg0 OP (itype)arg1) where the
inner conversions are extensions (or null).

There are several wrong-code bugs in this function, one of which is PR
36578; this patch fixes these bugs (for conversions of binary
arithmetic operations +-*/, not any other bugs for other conversions
it transforms).  The various cases are as follows:

(a) Args are decimal, TYPE is decimal.

Arithmetic is done in the widest decimal arg type (which may be
narrower than both ITYPE, the type originally used for the arithmetic,
and TYPE, the result type of the conversion) and converted to TYPE.

This can cause a lot of loss of precision, not just a change in the
last place; see the test convert-dfp-fold-2.c this patch adds for
example.  This patch changes this case to use the widest of the args
and TYPE; that may still suffer from a last-bit double-rounding
problem like that in PR 36578 (see (d) below), but I propose to leave
working out the correct logic for this to the decimal floating-point
experts.

(b) Args are binary, TYPE is decimal.

Arithmetic is done in TYPE if TYPE actually equals one of the decimal
float nodes (otherwise, e.g. qualification or typedef, falls through
to a binary case, and may use a binary or a decimal type depending on
the precisions involved).  See the test convert-bfp-13.c this patch
adds for example.

I believe this function should not try to optimize cases involving
binary-to-decimal conversions at all.  This patch disables such
optimization.

(c) Args are decimal, TYPE is binary.

Arithmetic is done in the widest decimal arg type (which may be
narrower than ITYPE) and converted to TYPE.  See the test
convert-bfp-14.c for example.

I believe this function should not try to optimize decimal-to-binary
conversions at all.  This patch disables such optimization.

(d) Args and TYPE are binary.

If the widest of TYPE and the arg types (NEWTYPE) is narrower than
ITYPE then the arithmetic is done in that NEWTYPE, and the result is
converted to TYPE.  (Otherwise, the arithmetic is done as expected in
ITYPE and the result converted to TYPE.)

The args are no wider than NEWTYPE, which is narrower than ITYPE.  So
the question is whether converting (type)((itype)arg0 OP (itype)arg1)
to (type)(arg0 OP arg1) is valid for ARG0 and ARG1 of type NEWTYPE.
NEWTYPE may be the same as TYPE, or wider than TYPE.

If NEWTYPE is the same as TYPE, e.g. (float)((double)float +
(double)float) coverted to (float + float), then the transformations
are correct if both types ITYPE and NEWTYPE follow IEEE rounding
semantics with the same rounding mode and the wider type has strictly
more than twice as many mantissa bits as the smaller one, can
represent infinities and NaNs if the smaller one can, and has
sufficient exponent range for the product or ratio of two values
representable in the smaller type to be within the range of normal
values of the wider type; this patch implements this logic.  (Proof
omitted.)  When ITYPE is not wide enough, this is PR 36578.  See
pr36578-1.c.  Note that the case of float and double above is safe and
seems to be a case where the transformation is worthwhile to preserve
(considering code such as "f = 2.0 * f;" as in a comment in
strip_float_extensions).

If NEWTYPE is wider than TYPE, e.g. (float)((long double)double +
(long double)double) converted to (float)(double + double), the
transformation is unsafe regardless of the details of the types
involved; double rounding can arise if the result of double arithmetic
is a double value half way between two representable float values but
the exact value is sufficiently different (in the right direction) for
this difference to be visible in long double arithmetic.  See
pr36578-2.c.

Bootstrapped with no regressions on i686-pc-linux-gnu.  OK to commit?

2008-10-29  Joseph Myers  <joseph@codesourcery.com>

	PR middle-end/36578
	* convert.c (convert_to_real): Do not optimize conversions of
	binary arithmetic operations between binary and decimal
	floating-point types.  Consider mode of target type in determining
	decimal type for arithmetic.  Unless
	flag_unsafe_math_optimizations, do not optimize binary conversions
	where this may change rounding behavior.

testsuite:
2008-10-29  Joseph Myers  <joseph@codesourcery.com>

	PR middle-end/36578
	* gcc.dg/dfp/convert-bfp-13.c, gcc.dg/dfp/convert-bfp-14.c,
	gcc.dg/dfp/convert-dfp-fold-2.c, gcc.target/i386/pr36578-1.c,
	gcc.target/i386/pr36578-2.c: New tests.

Index: convert.c
===================================================================
--- convert.c	(revision 141407)
+++ convert.c	(working copy)
@@ -263,18 +263,23 @@
 	     tree arg1 = strip_float_extensions (TREE_OPERAND (expr, 1));
 
 	     if (FLOAT_TYPE_P (TREE_TYPE (arg0))
-		 && FLOAT_TYPE_P (TREE_TYPE (arg1)))
+		 && FLOAT_TYPE_P (TREE_TYPE (arg1))
+		 && DECIMAL_FLOAT_TYPE_P (itype) == DECIMAL_FLOAT_TYPE_P (type))
 	       {
 		  tree newtype = type;
+		  const struct real_format *tfmt, *ifmt;
 
 		  if (TYPE_MODE (TREE_TYPE (arg0)) == SDmode
-		      || TYPE_MODE (TREE_TYPE (arg1)) == SDmode)
+		      || TYPE_MODE (TREE_TYPE (arg1)) == SDmode
+		      || TYPE_MODE (type) == SDmode)
 		    newtype = dfloat32_type_node;
 		  if (TYPE_MODE (TREE_TYPE (arg0)) == DDmode
-		      || TYPE_MODE (TREE_TYPE (arg1)) == DDmode)
+		      || TYPE_MODE (TREE_TYPE (arg1)) == DDmode
+		      || TYPE_MODE (type) == DDmode)
 		    newtype = dfloat64_type_node;
 		  if (TYPE_MODE (TREE_TYPE (arg0)) == TDmode
-		      || TYPE_MODE (TREE_TYPE (arg1)) == TDmode)
+		      || TYPE_MODE (TREE_TYPE (arg1)) == TDmode
+		      || TYPE_MODE (type) == TDmode)
                     newtype = dfloat128_type_node;
 		  if (newtype == dfloat32_type_node
 		      || newtype == dfloat64_type_node
@@ -292,7 +297,54 @@
 		    newtype = TREE_TYPE (arg0);
 		  if (TYPE_PRECISION (TREE_TYPE (arg1)) > TYPE_PRECISION (newtype))
 		    newtype = TREE_TYPE (arg1);
-		  if (TYPE_PRECISION (newtype) < TYPE_PRECISION (itype))
+		  tfmt = REAL_MODE_FORMAT (TYPE_MODE (type));
+		  ifmt = REAL_MODE_FORMAT (TYPE_MODE (itype));
+		  /* Sometimes this transformation is safe (cannot
+		     change results through affecting double rounding
+		     cases) and sometimes it is not.  If NEWTYPE is
+		     wider than TYPE, e.g. (float)((long double)double
+		     + (long double)double) converted to
+		     (float)(double + double), the transformation is
+		     unsafe regardless of the details of the types
+		     involved; double rounding can arise if the result
+		     of NEWTYPE arithmetic is a NEWTYPE value half way
+		     between two representable TYPE values but the
+		     exact value is sufficiently different (in the
+		     right direction) for this difference to be
+		     visible in ITYPE arithmetic.  If NEWTYPE is the
+		     same as TYPE, however, the transformation may be
+		     safe depending on the types involved: it is safe
+		     if the ITYPE has strictly more than twice as many
+		     mantissa bits as TYPE, can represent infinities
+		     and NaNs if the TYPE can, and has sufficient
+		     exponent range for the product or ratio of two
+		     values representable in the TYPE to be within the
+		     range of normal values of ITYPE.  */
+		  if (TYPE_PRECISION (newtype) < TYPE_PRECISION (itype)
+		      && (flag_unsafe_math_optimizations
+			  || (TYPE_PRECISION (newtype) == TYPE_PRECISION (type)
+			      /* These conditions are conservative
+				 rather than trying to catch the exact
+				 boundary conditions; the main case to
+				 allow is IEEE float and double.  */
+			      && ifmt->b == tfmt->b
+			      && ifmt->p > 2 * tfmt->p
+			      && ifmt->emin < 2 * tfmt->emin - tfmt->p - 2
+			      && ifmt->emin < (tfmt->emin - tfmt->emax
+					       - tfmt->p - 2)
+			      && ifmt->emax > 2 * tfmt->emax + 2
+			      && ifmt->emax > (tfmt->emax - tfmt->emin
+					       + tfmt->p + 2)
+			      && (ifmt->round_towards_zero
+				  == tfmt->round_towards_zero)
+			      && (ifmt->has_sign_dependent_rounding
+				  == tfmt->has_sign_dependent_rounding)
+			      && ifmt->has_nans >= tfmt->has_nans
+			      && ifmt->has_inf >= tfmt->has_inf
+			      && (ifmt->has_signed_zero
+				  >= tfmt->has_signed_zero)
+			      && !MODE_COMPOSITE_P (TYPE_MODE (type))
+			      && !MODE_COMPOSITE_P (TYPE_MODE (itype)))))
 		    {
 		      expr = build2 (TREE_CODE (expr), newtype,
 				     fold (convert_to_real (newtype, arg0)),
Index: testsuite/gcc.target/i386/pr36578-2.c
===================================================================
--- testsuite/gcc.target/i386/pr36578-2.c	(revision 0)
+++ testsuite/gcc.target/i386/pr36578-2.c	(revision 0)
@@ -0,0 +1,23 @@
+/* Test for unsafe floating-point conversions.  */
+/* { dg-do run } */
+/* { dg-options "-msse2 -mfpmath=sse" } */
+
+#include "sse2-check.h"
+
+extern void abort (void);
+extern void exit (int);
+extern int printf(const char *, ...);
+
+volatile double d1 = 0x1.000001p0;
+volatile double d2 = 0x1p-54;
+volatile float f = 0x1.000002p0f;
+volatile float f2;
+
+static void
+sse2_test (void)
+{
+  f2 = (float)((long double)d1 + (long double)d2);
+  if (f != f2)
+    abort ();
+  exit (0);
+}
Index: testsuite/gcc.target/i386/pr36578-1.c
===================================================================
--- testsuite/gcc.target/i386/pr36578-1.c	(revision 0)
+++ testsuite/gcc.target/i386/pr36578-1.c	(revision 0)
@@ -0,0 +1,22 @@
+/* Test for unsafe floating-point conversions.  PR 36578.  */
+/* { dg-do run } */
+/* { dg-options "-msse2 -mfpmath=sse" } */
+
+#include "sse2-check.h"
+
+extern void abort (void);
+extern void exit (int);
+extern int printf(const char *, ...);
+
+volatile double d1 = 1.0;
+volatile double d2 = 0x1.00001p-53;
+volatile double d3;
+
+static void
+sse2_test (void)
+{
+  d3 = (double)((long double)d1 + (long double)d2);
+  if (d3 != d1)
+    abort ();
+  exit (0);
+}
Index: testsuite/gcc.dg/dfp/convert-bfp-13.c
===================================================================
--- testsuite/gcc.dg/dfp/convert-bfp-13.c	(revision 0)
+++ testsuite/gcc.dg/dfp/convert-bfp-13.c	(revision 0)
@@ -0,0 +1,20 @@
+/* Test for bug where fold changed binary operation to decimal
+   depending on typedefs.  */
+/* { dg-options "-std=gnu99" } */
+
+extern void abort (void);
+extern void exit (int);
+
+volatile double d = 1.2345675;
+
+typedef const volatile _Decimal32 d32;
+
+int
+main (void)
+{
+  _Decimal32 a = (d * d);
+  d32 b = (d * d);
+  if (a != b)
+    abort ();
+  exit (0);
+}
Index: testsuite/gcc.dg/dfp/convert-dfp-fold-2.c
===================================================================
--- testsuite/gcc.dg/dfp/convert-dfp-fold-2.c	(revision 0)
+++ testsuite/gcc.dg/dfp/convert-dfp-fold-2.c	(revision 0)
@@ -0,0 +1,17 @@
+/* Test for bug where fold narrowed decimal floating-point
+   operations.  */
+/* { dg-options "-std=gnu99" } */
+
+extern void abort (void);
+extern void exit (int);
+
+volatile _Decimal32 f = 1.23456DF;
+volatile _Decimal64 d = 1.23456DD;
+
+int
+main (void)
+{
+  if ((_Decimal128)((_Decimal64)f * (_Decimal64)f) != (_Decimal128)(d * d))
+    abort ();
+  exit (0);
+}
Index: testsuite/gcc.dg/dfp/convert-bfp-14.c
===================================================================
--- testsuite/gcc.dg/dfp/convert-bfp-14.c	(revision 0)
+++ testsuite/gcc.dg/dfp/convert-bfp-14.c	(revision 0)
@@ -0,0 +1,17 @@
+/* Test for bug where fold narrowed decimal floating-point
+   operations.  */
+/* { dg-options "-std=gnu99" } */
+
+extern void abort (void);
+extern void exit (int);
+
+volatile _Decimal32 f = 1.23456DF;
+volatile _Decimal64 d = 1.23456DD;
+
+int
+main (void)
+{
+  if ((double)((_Decimal64)f * (_Decimal64)f) != (double)(d * d))
+    abort ();
+  exit (0);
+}

-- 
Joseph S. Myers
joseph@codesourcery.com