Fix catan, catanh inaccuracy through use of log (bug 15394).

diff --git a/ChangeLog b/ChangeLog
index 2d44b8e..c3404ac 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,19 @@
+2013-04-24  Joseph Myers  <joseph@codesourcery.com>
+
+       [BZ #15394]
+       * math/s_catan.c (__catan): Calculate imaginary part of result
+       with log1p not log unless computing log of number close to 0.
+       * math/s_catanf.c (__catanf): Likewise.
+       * math/s_catanl.c (__catanl): Likewise.
+       * math/s_catanh.c (__catanh): Calculate real part of result with
+       log1p not log unless computing log of number close to 0.
+       * math/s_catanhf.c (__catanhf): Likewise.
+       * math/s_catanhl.c (__catanhl): Likewise.
+       * math/libm-test.inc (catan_test): Add more tests.
+       (catanh_test): Likewise.
+       * sysdeps/i386/fpu/libm-test-ulps: Update.
+       * sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
+
 2013-04-24  Siddhesh Poyarekar  <siddhesh@redhat.com>
 
        * benchtests/Makefile: Mention files in which fast and slow

diff --git a/NEWS b/NEWS
index a1aeb42..e5c4f90 100644 (file)
--- a/NEWS
+++ b/NEWS
@@ -14,7 +14,7 @@ Version 2.18
   14888, 14920, 14964, 14981, 14982, 14985, 14994, 14996, 15003, 15006,
   15020, 15023, 15036, 15054, 15055, 15062, 15078, 15160, 15214, 15232,
   15234, 15283, 15285, 15287, 15304, 15305, 15307, 15309, 15327, 15330,
-  15335, 15336, 15337, 15342, 15346, 15361.
+  15335, 15336, 15337, 15342, 15346, 15361, 15394.
 
 * CVE-2013-0242 Buffer overrun in regexp matcher has been fixed (Bugzilla
   #15078).

diff --git a/math/libm-test.inc b/math/libm-test.inc
index 0049fcd..447b603 100644 (file)
--- a/math/libm-test.inc
+++ b/math/libm-test.inc
@@ -4299,6 +4299,35 @@ catan_test (void)
 
   TEST_c_c (catan, qnan_value, qnan_value, qnan_value, qnan_value);
 
+  TEST_c_c (catan, 0x1p50L, 0.0L, 1.570796326794895731052901991514519103193L, 0.0L);
+  TEST_c_c (catan, 0x1p50L, -0.0L, 1.570796326794895731052901991514519103193L, -0.0L);
+  TEST_c_c (catan, -0x1p50L, 0.0L, -1.570796326794895731052901991514519103193L, 0.0L);
+  TEST_c_c (catan, -0x1p50L, -0.0L, -1.570796326794895731052901991514519103193L, -0.0L);
+  TEST_c_c (catan, 0.0L, 0x1p50L, 1.570796326794896619231321691639751442099L, 8.881784197001252323389053344728897997441e-16L);
+  TEST_c_c (catan, -0.0L, 0x1p50L, -1.570796326794896619231321691639751442099L, 8.881784197001252323389053344728897997441e-16L);
+  TEST_c_c (catan, 0.0L, -0x1p50L, 1.570796326794896619231321691639751442099L, -8.881784197001252323389053344728897997441e-16L);
+  TEST_c_c (catan, -0.0L, -0x1p50L, -1.570796326794896619231321691639751442099L, -8.881784197001252323389053344728897997441e-16L);
+#ifndef TEST_FLOAT
+  TEST_c_c (catan, 0x1p500L, 0.0L, 1.570796326794896619231321691639751442099L, 0.0L);
+  TEST_c_c (catan, 0x1p500L, -0.0L, 1.570796326794896619231321691639751442099L, -0.0L);
+  TEST_c_c (catan, -0x1p500L, 0.0L, -1.570796326794896619231321691639751442099L, 0.0L);
+  TEST_c_c (catan, -0x1p500L, -0.0L, -1.570796326794896619231321691639751442099L, -0.0L);
+  TEST_c_c (catan, 0.0L, 0x1p500L, 1.570796326794896619231321691639751442099L, 3.054936363499604682051979393213617699789e-151L);
+  TEST_c_c (catan, -0.0L, 0x1p500L, -1.570796326794896619231321691639751442099L, 3.054936363499604682051979393213617699789e-151L);
+  TEST_c_c (catan, 0.0L, -0x1p500L, 1.570796326794896619231321691639751442099L, -3.054936363499604682051979393213617699789e-151L);
+  TEST_c_c (catan, -0.0L, -0x1p500L, -1.570796326794896619231321691639751442099L, -3.054936363499604682051979393213617699789e-151L);
+#endif
+#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384
+  TEST_c_c (catan, 0x1p5000L, 0.0L, 1.570796326794896619231321691639751442099L, 0.0L);
+  TEST_c_c (catan, 0x1p5000L, -0.0L, 1.570796326794896619231321691639751442099L, -0.0L);
+  TEST_c_c (catan, -0x1p5000L, 0.0L, -1.570796326794896619231321691639751442099L, 0.0L);
+  TEST_c_c (catan, -0x1p5000L, -0.0L, -1.570796326794896619231321691639751442099L, -0.0L);
+  TEST_c_c (catan, 0.0L, 0x1p5000L, 1.570796326794896619231321691639751442099L, 7.079811261048172892385615158694057552948e-1506L);
+  TEST_c_c (catan, -0.0L, 0x1p5000L, -1.570796326794896619231321691639751442099L, 7.079811261048172892385615158694057552948e-1506L);
+  TEST_c_c (catan, 0.0L, -0x1p5000L, 1.570796326794896619231321691639751442099L, -7.079811261048172892385615158694057552948e-1506L);
+  TEST_c_c (catan, -0.0L, -0x1p5000L, -1.570796326794896619231321691639751442099L, -7.079811261048172892385615158694057552948e-1506L);
+#endif
+
   TEST_c_c (catan, 0.75L, 1.25L, 1.10714871779409050301706546017853704L, 0.549306144334054845697622618461262852L);
   TEST_c_c (catan, -2, -3, -1.4099210495965755225306193844604208L, -0.22907268296853876629588180294200276L);
 
@@ -4365,6 +4394,35 @@ catanh_test (void)
 
   TEST_c_c (catanh, qnan_value, qnan_value, qnan_value, qnan_value);
 
+  TEST_c_c (catanh, 0x1p50L, 0.0L, 8.881784197001252323389053344728897997441e-16L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0x1p50L, -0.0L, 8.881784197001252323389053344728897997441e-16L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0x1p50L, 0.0L, -8.881784197001252323389053344728897997441e-16L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0x1p50L, -0.0L, -8.881784197001252323389053344728897997441e-16L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0.0L, 0x1p50L, 0.0L, 1.570796326794895731052901991514519103193L);
+  TEST_c_c (catanh, -0.0L, 0x1p50L, -0.0L, 1.570796326794895731052901991514519103193L);
+  TEST_c_c (catanh, 0.0L, -0x1p50L, 0.0L, -1.570796326794895731052901991514519103193L);
+  TEST_c_c (catanh, -0.0L, -0x1p50L, -0.0L, -1.570796326794895731052901991514519103193L);
+#ifndef TEST_FLOAT
+  TEST_c_c (catanh, 0x1p500L, 0.0L, 3.054936363499604682051979393213617699789e-151L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0x1p500L, -0.0L, 3.054936363499604682051979393213617699789e-151L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0x1p500L, 0.0L, -3.054936363499604682051979393213617699789e-151L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0x1p500L, -0.0L, -3.054936363499604682051979393213617699789e-151L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0.0L, 0x1p500L, 0.0L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0.0L, 0x1p500L, -0.0L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0.0L, -0x1p500L, 0.0L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0.0L, -0x1p500L, -0.0L, -1.570796326794896619231321691639751442099L);
+#endif
+#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384
+  TEST_c_c (catanh, 0x1p5000L, 0.0L, 7.079811261048172892385615158694057552948e-1506L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0x1p5000L, -0.0L, 7.079811261048172892385615158694057552948e-1506L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0x1p5000L, 0.0L, -7.079811261048172892385615158694057552948e-1506L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0x1p5000L, -0.0L, -7.079811261048172892385615158694057552948e-1506L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0.0L, 0x1p5000L, 0.0L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0.0L, 0x1p5000L, -0.0L, 1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, 0.0L, -0x1p5000L, 0.0L, -1.570796326794896619231321691639751442099L);
+  TEST_c_c (catanh, -0.0L, -0x1p5000L, -0.0L, -1.570796326794896619231321691639751442099L);
+#endif
+
   TEST_c_c (catanh, 0.75L, 1.25L, 0.261492138795671927078652057366532140L, 0.996825126463918666098902241310446708L);
   TEST_c_c (catanh, -2, -3, -0.14694666622552975204743278515471595L, -1.3389725222944935611241935759091443L);
 

diff --git a/math/s_catan.c b/math/s_catan.c
index 46c18bf..783941a 100644 (file)
--- a/math/s_catan.c
+++ b/math/s_catan.c
@@ -61,7 +61,7 @@ __catan (__complex__ double x)
     }
   else
     {
-      double r2, num, den;
+      double r2, num, den, f;
 
       r2 = __real__ x * __real__ x;
 
@@ -75,7 +75,14 @@ __catan (__complex__ double x)
       den = __imag__ x - 1.0;
       den = r2 + den * den;
 
-      __imag__ res = 0.25 * __ieee754_log (num / den);
+      f = num / den;
+      if (f < 0.5)
+       __imag__ res = 0.25 * __ieee754_log (f);
+      else
+       {
+         num = 4.0 * __imag__ x;
+         __imag__ res = 0.25 * __log1p (num / den);
+       }
     }
 
   return res;

diff --git a/math/s_catanf.c b/math/s_catanf.c
index 5a43247..0dc85ff 100644 (file)
--- a/math/s_catanf.c
+++ b/math/s_catanf.c
@@ -61,7 +61,7 @@ __catanf (__complex__ float x)
     }
   else
     {
-      float r2, num, den;
+      float r2, num, den, f;
 
       r2 = __real__ x * __real__ x;
 
@@ -75,7 +75,14 @@ __catanf (__complex__ float x)
       den = __imag__ x - 1.0;
       den = r2 + den * den;
 
-      __imag__ res = 0.25 * __ieee754_logf (num / den);
+      f = num / den;
+      if (f < 0.5)
+       __imag__ res = 0.25 * __ieee754_logf (f);
+      else
+       {
+         num = 4.0 * __imag__ x;
+         __imag__ res = 0.25 * __log1pf (num / den);
+       }
     }
 
   return res;

diff --git a/math/s_catanh.c b/math/s_catanh.c
index 5371f44..0ee8c64 100644 (file)
--- a/math/s_catanh.c
+++ b/math/s_catanh.c
@@ -64,7 +64,14 @@ __catanh (__complex__ double x)
       double den = 1.0 - __real__ x;
       den = i2 + den * den;
 
-      __real__ res = 0.25 * (__ieee754_log (num) - __ieee754_log (den));
+      double f = num / den;
+      if (f < 0.5)
+       __real__ res = 0.25 * __ieee754_log (f);
+      else
+       {
+         num = 4.0 * __real__ x;
+         __real__ res = 0.25 * __log1p (num / den);
+       }
 
       den = 1 - __real__ x * __real__ x - i2;
 

diff --git a/math/s_catanhf.c b/math/s_catanhf.c
index 8385af4..ca9a301 100644 (file)
--- a/math/s_catanhf.c
+++ b/math/s_catanhf.c
@@ -64,7 +64,14 @@ __catanhf (__complex__ float x)
       float den = 1.0 - __real__ x;
       den = i2 + den * den;
 
-      __real__ res = 0.25 * (__ieee754_logf (num) - __ieee754_logf (den));
+      float f = num / den;
+      if (f < 0.5)
+       __real__ res = 0.25 * __ieee754_logf (f);
+      else
+       {
+         num = 4.0 * __real__ x;
+         __real__ res = 0.25 * __log1pf (num / den);
+       }
 
       den = 1 - __real__ x * __real__ x - i2;
 

diff --git a/math/s_catanhl.c b/math/s_catanhl.c
index 6844f03..4897c0c 100644 (file)
--- a/math/s_catanhl.c
+++ b/math/s_catanhl.c
@@ -64,7 +64,14 @@ __catanhl (__complex__ long double x)
       long double den = 1.0 - __real__ x;
       den = i2 + den * den;
 
-      __real__ res = 0.25 * (__ieee754_logl (num) - __ieee754_logl (den));
+      long double f = num / den;
+      if (f < 0.5)
+       __real__ res = 0.25 * __ieee754_logl (f);
+      else
+       {
+         num = 4.0 * __real__ x;
+         __real__ res = 0.25 * __log1pl (num / den);
+       }
 
       den = 1 - __real__ x * __real__ x - i2;
 

diff --git a/math/s_catanl.c b/math/s_catanl.c
index 57d2e59..e04dba7 100644 (file)
--- a/math/s_catanl.c
+++ b/math/s_catanl.c
@@ -61,7 +61,7 @@ __catanl (__complex__ long double x)
     }
   else
     {
-      long double r2, num, den;
+      long double r2, num, den, f;
 
       r2 = __real__ x * __real__ x;
 
@@ -75,7 +75,14 @@ __catanl (__complex__ long double x)
       den = __imag__ x - 1.0;
       den = r2 + den * den;
 
-      __imag__ res = 0.25 * __ieee754_logl (num / den);
+      f = num / den;
+      if (f < 0.5)
+       __imag__ res = 0.25 * __ieee754_logl (f);
+      else
+       {
+         num = 4.0 * __imag__ x;
+         __imag__ res = 0.25 * __log1pl (num / den);
+       }
     }
 
   return res;

diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps
index fecaa94..a6d9b0f 100644 (file)
--- a/sysdeps/i386/fpu/libm-test-ulps
+++ b/sysdeps/i386/fpu/libm-test-ulps
@@ -3417,6 +3417,8 @@ ldouble: 1
 Test "Real part of: catanh (0.75 + 1.25 i) == 0.261492138795671927078652057366532140 + 0.996825126463918666098902241310446708 i":
 double: 1
 idouble: 1
+ildouble: 1
+ldouble: 1
 
 # cbrt
 Test "cbrt (-27.0) == -3.0":

diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps
index 9c9c473..2747fc5 100644 (file)
--- a/sysdeps/x86_64/fpu/libm-test-ulps
+++ b/sysdeps/x86_64/fpu/libm-test-ulps
@@ -3985,6 +3985,8 @@ ifloat: 4
 Test "Real part of: catanh (0.75 + 1.25 i) == 0.261492138795671927078652057366532140 + 0.996825126463918666098902241310446708 i":
 double: 1
 idouble: 1
+ildouble: 1
+ldouble: 1
 Test "Imaginary part of: catanh (0.75 + 1.25 i) == 0.261492138795671927078652057366532140 + 0.996825126463918666098902241310446708 i":
 float: 6
 ifloat: 6