about summary refs log tree commit diff
diff options
context:
space:
mode:
authorBenjamin Saunders <ben.e.saunders@gmail.com>2020-05-22 20:49:20 -0700
committerBenjamin Saunders <ben.e.saunders@gmail.com>2020-05-24 00:03:35 -0700
commit730f7366bba5abdf5ae0c2f1222795e40d48f90c (patch)
tree84536380e5fc3796d268a5f3e7ed9c0e770bd9de
parent215f2d3294b08dbdcf8f7d40de21ef1e7eae0a2d (diff)
downloadrust-730f7366bba5abdf5ae0c2f1222795e40d48f90c.tar.gz
rust-730f7366bba5abdf5ae0c2f1222795e40d48f90c.zip
Fix asinh of negative values
When `x` has large magnitude, `x + ((x * x) + 1.0).sqrt()` approaches
`x + x.abs()`. For negative values of `x`, this leads to catastrophic
cancellation, resulting in large errors or even 0 being passed to
`ln`, producing incorrect results including `-inf`.

Becuase asinh is an odd function, i.e. -asinh(x) = asinh(-x) for all
x, we can avoid the catastrophic cancellation and obtain correct
results by taking the absolute value of `self` for the first
term. `self * self` is always positive, so in effect this gives us
`x.abs().asinh().copysign(x)` which as discussed above is
algebraically equivalent, but is much more accurate.
-rw-r--r--src/libstd/f32.rs4
-rw-r--r--src/libstd/f64.rs4
2 files changed, 6 insertions, 2 deletions
diff --git a/src/libstd/f32.rs b/src/libstd/f32.rs
index 8e743ace99b..ab468f42b3b 100644
--- a/src/libstd/f32.rs
+++ b/src/libstd/f32.rs
@@ -835,7 +835,7 @@ impl f32 {
         if self == Self::NEG_INFINITY {
             Self::NEG_INFINITY
         } else {
-            (self + ((self * self) + 1.0).sqrt()).ln().copysign(self)
+            (self.abs() + ((self * self) + 1.0).sqrt()).ln().copysign(self)
         }
     }
 
@@ -1414,6 +1414,8 @@ mod tests {
         assert!((-0.0f32).asinh().is_sign_negative()); // issue 63271
         assert_approx_eq!(2.0f32.asinh(), 1.443635475178810342493276740273105f32);
         assert_approx_eq!((-2.0f32).asinh(), -1.443635475178810342493276740273105f32);
+        // regression test for the catastrophic cancellation fixed in 72486
+        assert_approx_eq!((-3000.0f32).asinh(), -8.699514775987968673236893537700647f32);
     }
 
     #[test]
diff --git a/src/libstd/f64.rs b/src/libstd/f64.rs
index fe64d27b1ef..c033198f021 100644
--- a/src/libstd/f64.rs
+++ b/src/libstd/f64.rs
@@ -837,7 +837,7 @@ impl f64 {
         if self == Self::NEG_INFINITY {
             Self::NEG_INFINITY
         } else {
-            (self + ((self * self) + 1.0).sqrt()).ln().copysign(self)
+            (self.abs() + ((self * self) + 1.0).sqrt()).ln().copysign(self)
         }
     }
 
@@ -1443,6 +1443,8 @@ mod tests {
         // issue 63271
         assert_approx_eq!(2.0f64.asinh(), 1.443635475178810342493276740273105f64);
         assert_approx_eq!((-2.0f64).asinh(), -1.443635475178810342493276740273105f64);
+        // regression test for the catastrophic cancellation fixed in 72486
+        assert_approx_eq!((-67452098.07139316f64).asinh(), -18.72007542627454439398548429400083);
     }
 
     #[test]