[PATCH] range-op-float: frange_arithmetic tweaks for MODE_COMPOSITE_P

[Jakub Jelinek <jakub@redhat.com>]

Hi!

As mentioned in PR107967, ibm-ldouble-format documents that
+- has 1ulp accuracy, * 2ulps and / 3ulps.
So, even if the result is exact, we need to widen the range a little bit.

The following patch does that.  I just wonder what it means for reverse
division (the op1_range case), which we implement through multiplication,
when division has 3ulps error and multiplication just 2ulps.  In any case,
this format is a mess and for non-default rounding modes can't be trusted
at all, instead of +inf or something close to it it happily computes -inf.

2022-12-07  Jakub Jelinek  <jakub@redhat.com>

	* range-op-float.cc (frange_arithmetic): For mode_composite,
	on top of rounding in the right direction accept extra 1ulp
	error for PLUS/MINUS_EXPR, extra 2ulps error for MULT_EXPR
	and extra 3ulps error for RDIV_EXPR.

--- gcc/range-op-float.cc.jj	2022-12-07 12:46:01.536123757 +0100
+++ gcc/range-op-float.cc	2022-12-07 12:50:40.812085139 +0100
@@ -344,22 +344,70 @@ frange_arithmetic (enum tree_code code,
 	    }
 	}
     }
-  if (round && (inexact || !real_identical (&result, &value)))
+  if (!inexact && !real_identical (&result, &value))
+    inexact = true;
+  if (round && (inexact || mode_composite))
     {
       if (mode_composite)
 	{
-	  if (real_isdenormal (&result, mode)
-	      || real_iszero (&result))
+	  if (real_isdenormal (&result, mode) || real_iszero (&result))
 	    {
 	      // IBM extended denormals only have DFmode precision.
 	      REAL_VALUE_TYPE tmp;
 	      real_convert (&tmp, DFmode, &value);
-	      frange_nextafter (DFmode, tmp, inf);
+	      if (inexact)
+		frange_nextafter (DFmode, tmp, inf);
+	      switch (code)
+		{
+		case PLUS_EXPR:
+		case MINUS_EXPR:
+		  // ibm-ldouble-format documents 1ulp for + and -.
+		  frange_nextafter (DFmode, tmp, inf);
+		  break;
+		case MULT_EXPR:
+		  // ibm-ldouble-format documents 2ulps for *.
+		  frange_nextafter (DFmode, tmp, inf);
+		  frange_nextafter (DFmode, tmp, inf);
+		  break;
+		case RDIV_EXPR:
+		  // ibm-ldouble-format documents 3ulps for /.
+		  frange_nextafter (DFmode, tmp, inf);
+		  frange_nextafter (DFmode, tmp, inf);
+		  frange_nextafter (DFmode, tmp, inf);
+		  break;
+		default:
+		  if (!inexact)
+		    return;
+		  break;
+		}
 	      real_convert (&result, mode, &tmp);
 	      return;
 	    }
 	}
-      frange_nextafter (mode, result, inf);
+      if (inexact)
+	frange_nextafter (mode, result, inf);
+      if (mode_composite)
+	switch (code)
+	  {
+	  case PLUS_EXPR:
+	  case MINUS_EXPR:
+	    // ibm-ldouble-format documents 1ulp for + and -.
+	    frange_nextafter (mode, result, inf);
+	    break;
+	  case MULT_EXPR:
+	    // ibm-ldouble-format documents 2ulps for *.
+	    frange_nextafter (mode, result, inf);
+	    frange_nextafter (mode, result, inf);
+	    break;
+	  case RDIV_EXPR:
+	    // ibm-ldouble-format documents 3ulps for /.
+	    frange_nextafter (mode, result, inf);
+	    frange_nextafter (mode, result, inf);
+	    frange_nextafter (mode, result, inf);
+	    break;
+	  default:
+	    break;
+	  }
     }
 }
 

	Jakub