Enforce zero-sensitive min/max

Author: juntyr

quaddtype/numpy_quaddtype/src/ops.hpp

@@ -514,33 +514,53 @@ quad_mod(const Sleef_quad *a, const Sleef_quad *b)
 static inline Sleef_quad
 quad_minimum(const Sleef_quad *in1, const Sleef_quad *in2)
 {
-    return Sleef_iunordq1(*in1, *in2)   ? (Sleef_iunordq1(*in1, *in1) ? *in1 : *in2)
-        : Sleef_icmpltq1(*in1, *in2) ? *in1
-                                    : *in2;
+    if (Sleef_iunordq1(*in1, *in2)) {
+        return Sleef_iunordq1(*in1, *in1) ? *in1 : *in2;
+    }
+    // minimum(-0.0, +0.0) = -0.0
+    if (Sleef_icmpeqq1(*in1, QUAD_ZERO) && Sleef_icmpeqq1(*in2, QUAD_ZERO)) {
+        return Sleef_icmpleq1(Sleef_copysignq1(QUAD_ONE, *in1), Sleef_copysignq1(QUAD_ONE, *in2)) ? *in1 : *in2;
+    }
+    return Sleef_fminq1(*in1, *in2);
 }
 
 static inline Sleef_quad
 quad_maximum(const Sleef_quad *in1, const Sleef_quad *in2)
 {
-    return Sleef_iunordq1(*in1, *in2)   ? (Sleef_iunordq1(*in1, *in1) ? *in1 : *in2)
-        : Sleef_icmpgtq1(*in1, *in2) ? *in1
-                                    : *in2;
+    if (Sleef_iunordq1(*in1, *in2)) {
+        return Sleef_iunordq1(*in1, *in1) ? *in1 : *in2;
+    }
+    // maximum(-0.0, +0.0) = +0.0
+    if (Sleef_icmpeqq1(*in1, QUAD_ZERO) && Sleef_icmpeqq1(*in2, QUAD_ZERO)) {
+        return Sleef_icmpgeq1(Sleef_copysignq1(QUAD_ONE, *in1), Sleef_copysignq1(QUAD_ONE, *in2)) ? *in1 : *in2;
+    }
+    return Sleef_fmaxq1(*in1, *in2);
 }
 
 static inline Sleef_quad
 quad_fmin(const Sleef_quad *in1, const Sleef_quad *in2)
 {
-    return Sleef_iunordq1(*in1, *in2)   ? (Sleef_iunordq1(*in2, *in2) ? *in1 : *in2)
-           : Sleef_icmpleq1(*in1, *in2) ? *in1
-                                        : *in2;
+    if (Sleef_iunordq1(*in1, *in2)) {
+        return Sleef_iunordq1(*in2, *in2) ? *in1 : *in2;
+    }
+    // fmin(-0.0, +0.0) = -0.0
+    if (Sleef_icmpeqq1(*in1, QUAD_ZERO) && Sleef_icmpeqq1(*in2, QUAD_ZERO)) {
+        return Sleef_icmpleq1(Sleef_copysignq1(QUAD_ONE, *in1), Sleef_copysignq1(QUAD_ONE, *in2)) ? *in1 : *in2;
+    }
+    return Sleef_fminq1(*in1, *in2);
 }
 
 static inline Sleef_quad
 quad_fmax(const Sleef_quad *in1, const Sleef_quad *in2)
 {
-    return Sleef_iunordq1(*in1, *in2)   ? (Sleef_iunordq1(*in2, *in2) ? *in1 : *in2)
-        : Sleef_icmpgeq1(*in1, *in2) ? *in1
-                                    : *in2;
+    if (Sleef_iunordq1(*in1, *in2)) {
+        return Sleef_iunordq1(*in2, *in2) ? *in1 : *in2;
+    }
+    // maximum(-0.0, +0.0) = +0.0
+    if (Sleef_icmpeqq1(*in1, QUAD_ZERO) && Sleef_icmpeqq1(*in2, QUAD_ZERO)) {
+        return Sleef_icmpgeq1(Sleef_copysignq1(QUAD_ONE, *in1), Sleef_copysignq1(QUAD_ONE, *in2)) ? *in1 : *in2;
+    }
+    return Sleef_fmaxq1(*in1, *in2);
 }
 
 static inline Sleef_quad

quaddtype/tests/test_quaddtype.py

@@ -97,6 +97,13 @@ def test_array_minmax(op, a, b):
     quad_res = op_func(quad_a, quad_b)
     float_res = op_func(float_a, float_b)
 
+    # native implementation may not be sensitive to zero signs
+    #  but we want to enforce it for the quad dtype
+    # e.g. min(+0.0, -0.0) = -0.0
+    if float_a == 0.0 and float_b == 0.0:
+        assert float_res == 0.0
+        float_res = np.copysign(0.0, op_fun(np.copysign(1.0, float_a), np.copysign(1.0, float_b)))
+
     np.testing.assert_array_equal(quad_res.astype(float), float_res)
 
     # Check sign for zero results
@@ -116,6 +123,13 @@ def test_array_aminmax(op, a, b):
     quad_res = op_func(quad_ab)
     float_res = op_func(float_ab)
 
+    # native implementation may not be sensitive to zero signs
+    #  but we want to enforce it for the quad dtype
+    # e.g. min(+0.0, -0.0) = -0.0
+    if float(a) == 0.0 and float(b) == 0.0:
+        assert float_res == 0.0
+        float_res = np.copysign(0.0, op_fun(np.array([np.copysign(1.0, float(a)), np.copysign(1.0, float(b))])))
+
     np.testing.assert_array_equal(np.array(quad_res).astype(float), float_res)
 
     # Check sign for zero results
@@ -490,7 +504,7 @@ def test_mod(a, b, backend, op):
     "0.9", "-0.9", "0.9999", "-0.9999",
     "1.1", "-1.1", "1.0001", "-1.0001",
     # Medium values
-    "10.0", "-10.0", "20.0", "-20.0",
+    "10.0", "-10.0", "20.0x", "-20.0",
     # Large values
     "100.0", "200.0", "700.0", "1000.0",
     "-100.0", "-200.0", "-700.0", "-1000.0",