This is the mail archive of the gcc-patches@gcc.gnu.org mailing list for the GCC 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]

[BIB PATCH] Compile-time sqrt error analysis


To test the accuracy of GCC's real_sqrt routine, I've just run a
modified version of Gaston Gonnet's "sqrt" testsuite against the
gcc-3_4-basic-improvement-branch.

Details of the Gonnet FPAccuracy test set can be found at
http://www.inf.ethz.ch/personal/gonnet/FPAccuracy/Analysis.html
For sqrt, this test suite includes 1055 "difficult" arguments.
A minor modification to this software, allowed me to evaluate
the accuracy of GCC's optimization that evaluates sqrt of a
constant argument at compile-time.

The results of using real.c's real_sqrt routine:

20 largest ulp errors (stored in a double)
   0.50000 ulp for sqrt(0.999999999999999889) = 0.999999999999999889)
   0.50000 ulp for sqrt(15.9999999999999982) = 3.99999999999999956)
   0.50000 ulp for sqrt(0.0156249999999999983) = 0.124999999999999986)
   0.50000 ulp for sqrt(0.249999999999999972) = 0.499999999999999944)
   0.50000 ulp for sqrt(63.9999999999999929) = 7.99999999999999911)
   0.50000 ulp for sqrt(0.0624999999999999931) = 0.249999999999999972)
   0.50000 ulp for sqrt(3.99999999999999956) = 1.99999999999999978)
   0.50000 ulp for sqrt(16.0000000000000036) = 4)
   0.50000 ulp for sqrt(0.0156250000000000035) = 0.125)
   0.50000 ulp for sqrt(4.00000000000000089) = 2)
   0.50000 ulp for sqrt(64.0000000000000142) = 8)
   0.50000 ulp for sqrt(1.00000000000000022) = 1)
   0.50000 ulp for sqrt(0.250000000000000056) = 0.5)
   0.50000 ulp for sqrt(0.0625000000000000139) = 0.25)
   0.50000 ulp for sqrt(2.53884156217198101e-308) = 1.59337426933284621e-154)
   0.50000 ulp for sqrt(7.05037261855188913) = 2.65525377667594897)
   0.50000 ulp for sqrt(9.53144511356327229) = 3.08730385831444742)
   0.50000 ulp for sqrt(491.057774316958671) = 22.159823427025735)
   0.50000 ulp for sqrt(0.374608653801551006) = 0.6120528194539675)
   0.50000 ulp for sqrt(3.25760941642394334e-308) = 1.80488487622450423e-154)


The good news is that our current implementation achieves the perfect
score of 0.5 ulp for all 1055 test cases (for atleast IEEE double
precision).


This is particularly significant when compared against running
the same testsuite with mainline GCC on x86 (i.e. at run-time).

20 largest ulp errors (stored in a double)
   0.50024 ulp for sqrt(3.63861067050296029e-308) = 1.90751426482292942e-154)
   0.50024 ulp for sqrt(376.094201960404519) = 19.3931483251277328)
   0.50023 ulp for sqrt(0.598293376600028348) = 0.77349426410286215)
   0.50023 ulp for sqrt(3.38197028074290468e-308) = 1.83901339873936353e-154)
   0.50022 ulp for sqrt(0.111069327000007101) = 0.333270651273116281)
   0.50021 ulp for sqrt(0.403139530700020954) = 0.634932697771992904)
   0.50021 ulp for sqrt(382.516844720041036) = 19.5580378545507756)
   0.50019 ulp for sqrt(9.81613492599978166) = 3.13307116516682083)
   0.50019 ulp for sqrt(0.542851097300020236) = 0.736784294960214359)
   0.50019 ulp for sqrt(6.30543715450036402) = 2.51106295311375316)
   0.50019 ulp for sqrt(416.244840379996958) = 20.4020793151089634)
   0.50018 ulp for sqrt(491.057774320022588) = 22.1598234270948709)
   0.50018 ulp for sqrt(7.5635572969987992) = 2.75019222909941341)
   0.50018 ulp for sqrt(3.83071649195763422e-308) = 1.95722162566165073e-154)
   0.50017 ulp for sqrt(0.249590039300010919) = 0.499589871094291649)
   0.50017 ulp for sqrt(0.741966908100119227) = 0.861375010143734743)
   0.50017 ulp for sqrt(474.497615900024471) = 21.7829661869090856)
   0.50017 ulp for sqrt(0.42561435630000477) = 0.652391260134594919)
   0.50016 ulp for sqrt(0.374701946099967154) = 0.612129027329996189)
   0.50016 ulp for sqrt(3.38197028074333106e-308) = 1.83901339873947946e-154)


Not quite so good!  The problem here is not with the Pentium/Athlon
implementation of "fsqrt" [which also achieves a perfect 0.5 ulp].
Instead if the sqrt isn't evaluated at compile-time, the result
can't be constant folded with a multiplication, which then takes
place at run-time.  Unfortunately, Intel/AMD architectures have a
0.00024 ulp rounding error for multiplication in their default
rounding mode [caused by double rounding].


I'd therefore like to propose the following patch to enable constant
folding of sqrt without "-ffast-math" or "-funsafe-math-optimizations".
As demonstrated above not only does this produce identical results for
sqrt, but actually improves accuracy of this particular validation set.


The following patch has been tested by a complete bootstrap, all
languages except Ada and treelang, on i686-pc-linux-gnu, and with
a full "make -k check" with no new regressions.


Ok for the gcc-3_4-basic-improvements-branch?


2002-12-08  Roger Sayle  <roger@eyesopen.com>

	* builtins.c (fold_builtin): Remove -funsafe-math-optimizations
	check for evaluating sqrt of a constant at compile time.
	* simplify-rtx.c (simplify_unary_operation): Likewise.


Index: builtins.c
===================================================================
RCS file: /cvs/gcc/gcc/gcc/builtins.c,v
retrieving revision 1.159.4.7
diff -c -3 -p -r1.159.4.7 builtins.c
*** builtins.c	3 Dec 2002 17:34:43 -0000	1.159.4.7
--- builtins.c	8 Dec 2002 18:36:10 -0000
*************** fold_builtin (exp)
*** 4299,4311 ****
  	  enum built_in_function fcode;
  	  tree arg = TREE_VALUE (arglist);

- 	  /* Optimize sqrt(0.0) = 0.0 and sqrt(1.0) = 1.0.  */
- 	  if (real_zerop (arg) || real_onep (arg))
- 	    return arg;
-
  	  /* Optimize sqrt of constant value.  */
! 	  if (flag_unsafe_math_optimizations
! 	      && TREE_CODE (arg) == REAL_CST
  	      && ! TREE_CONSTANT_OVERFLOW (arg))
  	    {
  	      enum machine_mode mode;
--- 4299,4306 ----
  	  enum built_in_function fcode;
  	  tree arg = TREE_VALUE (arglist);

  	  /* Optimize sqrt of constant value.  */
! 	  if (TREE_CODE (arg) == REAL_CST
  	      && ! TREE_CONSTANT_OVERFLOW (arg))
  	    {
  	      enum machine_mode mode;
Index: simplify-rtx.c
===================================================================
RCS file: /cvs/gcc/gcc/gcc/simplify-rtx.c,v
retrieving revision 1.118.4.9
diff -c -3 -p -r1.118.4.9 simplify-rtx.c
*** simplify-rtx.c	1 Dec 2002 05:42:18 -0000	1.118.4.9
--- simplify-rtx.c	8 Dec 2002 18:36:13 -0000
*************** simplify_unary_operation (code, mode, op
*** 579,586 ****
        switch (code)
  	{
  	case SQRT:
- 	  if (! flag_unsafe_math_optimizations)
- 	    return 0;
  	  if (HONOR_SNANS (mode) && real_isnan (&d))
  	    return 0;
  	  real_sqrt (&t, mode, &d);
--- 579,584 ----


Roger
--
Roger Sayle,                         E-mail: roger@eyesopen.com
OpenEye Scientific Software,         WWW: http://www.eyesopen.com/
Suite 1107, 3600 Cerrillos Road,     Tel: (+1) 505-473-7385
Santa Fe, New Mexico, 87507.         Fax: (+1) 505-473-0833


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