Forward-merge release/26.12 into main

cpp/include/cudf/detail/operators/arithmetic.cuh

@@ -126,9 +126,18 @@ __device__ auto pymod(A a, B b) -> decltype((a % b + b) % b)
 }
 
 template <cuda::std::floating_point A, cuda::std::floating_point B>
-__device__ auto pymod(A a, B b) -> decltype(cuda::std::fmod(cuda::std::fmod(a, b) + b, b))
-{
-  return cuda::std::fmod(cuda::std::fmod(a, b) + b, b);
+__device__ auto pymod(A a, B b) -> cuda::std::common_type_t<A, B>
+{
+  using common_t = cuda::std::common_type_t<A, B>;
+  auto const x   = static_cast<common_t>(a);
+  auto const y   = static_cast<common_t>(b);
+  auto mod       = cuda::std::fmod(x, y);
+  if (mod != common_t{0.0}) {
+    if ((y < common_t{0.0}) != (mod < common_t{0.0})) { mod += y; }
+  } else {
+    mod = cuda::std::copysign(common_t{0.0}, y);
+  }
+  return mod;
 }
 
 /**

cpp/src/binaryop/compiled/operation.cuh

@@ -6,6 +6,7 @@
 #pragma once
 
 #include <cudf/detail/integral_math.cuh>
+#include <cudf/detail/operators/arithmetic.cuh>
 #include <cudf/utilities/concepts.hpp>
 #include <cudf/utilities/traits.hpp>
 
@@ -192,9 +193,7 @@ struct PyMod {
   __device__ inline auto operator()(TypeLhs x, TypeRhs y) -> double
     requires(cuda::std::floating_point<std::common_type_t<TypeLhs, TypeRhs>>)
   {
-    auto x1 = static_cast<double>(x);
-    auto y1 = static_cast<double>(y);
-    return fmod(fmod(x1, y1) + y1, y1);
+    return cudf::detail::ops::pymod(static_cast<double>(x), static_cast<double>(y));
   }
 
   template <cudf::concepts::duration TypeLhs, typename TypeRhs>

cpp/tests/ast/transform_tests.cpp

@@ -973,6 +973,25 @@ TYPED_TEST(TransformTest, PyMod)
   CUDF_TEST_EXPECT_COLUMNS_EQUAL(expected, result->view(), verbosity);
 }
 
+TYPED_TEST(TransformTest, PyModFloatingEdgeCases)
+{
+  using Executor = TypeParam;
+
+  auto const inf = std::numeric_limits<double>::infinity();
+  auto c_0       = column_wrapper<double>{1.0, -1.0, 1.0, -1.0, -1e-20, 1e-20};
+  auto c_1       = column_wrapper<double>{inf, inf, -1e18, -1e18, 1.0, -1.0};
+  auto table     = cudf::table_view{{c_0, c_1}};
+
+  auto col_ref_0  = cudf::ast::column_reference(0);
+  auto col_ref_1  = cudf::ast::column_reference(1);
+  auto expression = cudf::ast::operation(cudf::ast::ast_operator::PYMOD, col_ref_0, col_ref_1);
+
+  auto result   = Executor::compute_column(table, expression);
+  auto expected = column_wrapper<double>{1.0, inf, -1e18, -1.0, 1.0, -1.0};
+
+  CUDF_TEST_EXPECT_COLUMNS_EQUAL(expected, result->view(), verbosity);
+}
+
 TYPED_TEST(TransformTest, BasicEqualityNullEqualNoNulls)
 {
   using Executor = TypeParam;

cpp/tests/binaryop/operator_parity_test.cpp

@@ -12,6 +12,7 @@
 #include <cudf/scalar/scalar.hpp>
 #include <cudf/utilities/error.hpp>
 
+#include <cmath>
 #include <limits>
 
 struct BinaryOperatorParityTest : public cudf::test::BaseFixture {};
@@ -210,3 +211,59 @@ TEST_F(BinaryOperatorParityTest, CheckedDecimalScaleRules)
       lhs, rhs, cudf::binary_operator::ADD_OVERFLOW, cudf::data_type{cudf::type_id::DECIMAL32, -1}),
     cudf::data_type_error);
 }
+
+TEST_F(BinaryOperatorParityTest, FloatingPyModEdgeCases)
+{
+  auto const inf  = std::numeric_limits<double>::infinity();
+  auto const type = cudf::data_type{cudf::type_id::FLOAT64};
+
+  // Vector - Vector
+  auto lhs = cudf::test::fixed_width_column_wrapper<double>{1.0, -1.0, 1.0, -1.0, 1.0, -1.0};
+  auto rhs = cudf::test::fixed_width_column_wrapper<double>{inf, inf, -inf, -inf, -1e18, -1e18};
+  auto expected = cudf::test::fixed_width_column_wrapper<double>{1.0, inf, -inf, -1.0, -1e18, -1.0};
+
+  auto result = cudf::binary_operation(lhs, rhs, cudf::binary_operator::PYMOD, type);
+  CUDF_TEST_EXPECT_COLUMNS_EQUAL(expected, result->view());
+
+  // Vector - Scalar with positive infinity divisor (#24898)
+  auto s_inf      = cudf::numeric_scalar<double>{inf};
+  auto v_lhs      = cudf::test::fixed_width_column_wrapper<double>{1.0, 2.5};
+  auto exp_inf    = cudf::test::fixed_width_column_wrapper<double>{1.0, 2.5};
+  auto result_inf = cudf::binary_operation(v_lhs, s_inf, cudf::binary_operator::PYMOD, type);
+  CUDF_TEST_EXPECT_COLUMNS_EQUAL(exp_inf, result_inf->view());
+
+  // Vector - Scalar with large negative divisor (#24895)
+  auto s_neg      = cudf::numeric_scalar<double>{-1e18};
+  auto exp_neg    = cudf::test::fixed_width_column_wrapper<double>{-1e18, -1e18};
+  auto result_neg = cudf::binary_operation(v_lhs, s_neg, cudf::binary_operator::PYMOD, type);
+  CUDF_TEST_EXPECT_COLUMNS_EQUAL(exp_neg, result_neg->view());
+
+  // mod += y rounds to y; Python returns y: -1e-20 % 1.0 == 1.0, 1e-20 % -1.0 == -1.0
+  auto lhs_tiny    = cudf::test::fixed_width_column_wrapper<double>{-1e-20, 1e-20};
+  auto rhs_tiny    = cudf::test::fixed_width_column_wrapper<double>{1.0, -1.0};
+  auto exp_tiny    = cudf::test::fixed_width_column_wrapper<double>{1.0, -1.0};
+  auto result_tiny = cudf::binary_operation(lhs_tiny, rhs_tiny, cudf::binary_operator::PYMOD, type);
+  CUDF_TEST_EXPECT_COLUMNS_EQUAL(exp_tiny, result_tiny->view());
+
+  // Float32
+  auto const inf_f  = std::numeric_limits<float>::infinity();
+  auto const type_f = cudf::data_type{cudf::type_id::FLOAT32};
+  auto lhs_f        = cudf::test::fixed_width_column_wrapper<float>{1.0f, -1.0f};
+  auto rhs_f        = cudf::test::fixed_width_column_wrapper<float>{inf_f, inf_f};
+  auto expected_f   = cudf::test::fixed_width_column_wrapper<float>{1.0f, inf_f};
+  auto result_f     = cudf::binary_operation(lhs_f, rhs_f, cudf::binary_operator::PYMOD, type_f);
+  CUDF_TEST_EXPECT_COLUMNS_EQUAL(expected_f, result_f->view());
+
+  // Exact-multiple signed zero tests (1.0 % -1.0 == -0.0, -1.0 % 1.0 == +0.0)
+  auto s_neg_one    = cudf::numeric_scalar<double>{-1.0};
+  auto s_pos_one    = cudf::numeric_scalar<double>{1.0};
+  auto v_one        = cudf::test::fixed_width_column_wrapper<double>{1.0};
+  auto v_neg_one    = cudf::test::fixed_width_column_wrapper<double>{-1.0};
+  auto res_neg_zero = cudf::binary_operation(v_one, s_neg_one, cudf::binary_operator::PYMOD, type);
+  auto res_pos_zero =
+    cudf::binary_operation(v_neg_one, s_pos_one, cudf::binary_operator::PYMOD, type);
+  auto const neg_data = cudf::test::to_host<double>(res_neg_zero->view()).first;
+  auto const pos_data = cudf::test::to_host<double>(res_pos_zero->view()).first;
+  EXPECT_TRUE(std::signbit(neg_data.front()));
+  EXPECT_FALSE(std::signbit(pos_data.front()));
+}

cpp/tests/binaryop/util/operation.h

@@ -323,9 +323,15 @@ struct PyMod {
   TypeOut operator()(TypeLhs x, TypeRhs y) const
   {
     if constexpr (std::is_floating_point_v<TypeLhs> or std::is_floating_point_v<TypeRhs>) {
-      auto x1 = static_cast<double>(x);
-      auto y1 = static_cast<double>(y);
-      return fmod(fmod(x1, y1) + y1, y1);
+      auto const x1 = static_cast<double>(x);
+      auto const y1 = static_cast<double>(y);
+      auto mod      = std::fmod(x1, y1);
+      if (mod != 0.0) {
+        if ((y1 < 0.0) != (mod < 0.0)) { mod += y1; }
+      } else {
+        mod = std::copysign(0.0, y1);
+      }
+      return static_cast<TypeOut>(mod);
     } else {
       return ((x % y) + y) % y;
     }

python/cudf/cudf/tests/series/test_binops.py

@@ -3338,3 +3338,38 @@ def test_arrow_backed_unsigned_subtract_overflow():
     gother = cudf.from_pandas(other)
     with pytest.raises(pa.ArrowInvalid):
         gs - gother
+
+
+@pytest.mark.parametrize("divisor", [float("inf"), float("-inf"), -1e18, 1e18])
+def test_float_pymod_edge_cases(divisor):
+    val = [1.0, -1.0, 2.5, 0.0, -0.0]
+    psr = pd.Series(val)
+    gsr = cudf.Series(val)
+    res_p = psr % divisor
+    res_g = gsr % divisor
+    assert_eq(res_p, res_g)
+    assert_eq(psr.mod(divisor), gsr.mod(divisor))
+    res_rp = divisor % psr
+    res_rg = divisor % gsr
+    assert_eq(res_rp, res_rg)
+    mask = ~np.isnan(res_rp.to_numpy())
+    np.testing.assert_equal(
+        np.signbit(res_rg.to_numpy()[mask]),
+        np.signbit(res_rp.to_numpy()[mask]),
+    )
+    np.testing.assert_equal(
+        np.signbit(res_g.to_numpy()), np.signbit(res_p.to_numpy())
+    )
+
+
+def test_float_pymod_series():
+    psr1 = pd.Series([1.0, -1.0, 1.0, -1.0])
+    psr2 = pd.Series([float("inf"), float("inf"), -1e18, -1e18])
+    gsr1 = cudf.from_pandas(psr1)
+    gsr2 = cudf.from_pandas(psr2)
+    res_p = psr1 % psr2
+    res_g = gsr1 % gsr2
+    assert_eq(res_p, res_g)
+    np.testing.assert_equal(
+        np.signbit(res_g.to_numpy()), np.signbit(res_p.to_numpy())
+    )