std.math.complex: fix cosh/tanh

igor84 · Aug 11, 2024 · 659056a · 659056a
1 parent c1c1b2c
commit 659056a
Show file tree

Hide file tree

Showing 2 changed files with 20 additions and 2 deletions.
diff --git a/lib/std/math/complex/cosh.zig b/lib/std/math/complex/cosh.zig
@@ -123,7 +123,7 @@ fn cosh64(z: Complex(f64)) Complex(f64) {
         }
         // x >= 1455: result always overflows
         else {
-            const h = 0x1p1023;
+            const h = 0x1p1023 * x;
             return Complex(f64).init(h * h * @cos(y), h * @sin(y));
         }
     }
@@ -170,3 +170,12 @@ test cosh64 {
     try testing.expectApproxEqAbs(-73.46729221264526, c.re, epsilon);
     try testing.expectApproxEqAbs(10.471557674805572, c.im, epsilon);
 }
+
+test "cosh64 musl" {
+    const epsilon = math.floatEps(f64);
+    const a = Complex(f64).init(7.44648873421389e17, 1.6008058402057622e19);
+    const c = cosh(a);
+
+    try testing.expectApproxEqAbs(std.math.inf(f64), c.re, epsilon);
+    try testing.expectApproxEqAbs(std.math.inf(f64), c.im, epsilon);
+}
diff --git a/lib/std/math/complex/tanh.zig b/lib/std/math/complex/tanh.zig
@@ -70,7 +70,7 @@ fn tanh64(z: Complex(f64)) Complex(f64) {
     const ix = hx & 0x7fffffff;
 
     if (ix >= 0x7ff00000) {
-        if ((ix & 0x7fffff) | lx != 0) {
+        if ((ix & 0xfffff) | lx != 0) {
             const r = if (y == 0) y else x * y;
             return Complex(f64).init(x, r);
         }
@@ -118,3 +118,12 @@ test tanh64 {
     try testing.expectApproxEqAbs(0.9999128201513536, c.re, epsilon);
     try testing.expectApproxEqAbs(-0.00002536867620767604, c.im, epsilon);
 }
+
+test "tanh64 musl" {
+    const epsilon = math.floatEps(f64);
+    const a = Complex(f64).init(std.math.inf(f64), std.math.inf(f64));
+    const c = tanh(a);
+
+    try testing.expectApproxEqAbs(1, c.re, epsilon);
+    try testing.expectApproxEqAbs(0, c.im, epsilon);
+}