Skip to content

Commit f428d4d

Browse files
committed
SIMD, TST: Apply Marten’s suggestions
- Extend the C++ doc scope to better explain precision control, which is chosen based on the data type. - Add test cases for sine/cosine facts—for example, `cos(0) == 1`.
1 parent 09d5154 commit f428d4d

2 files changed

Lines changed: 64 additions & 5 deletions

File tree

numpy/_core/src/common/simd/simd.hpp

Lines changed: 16 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -76,17 +76,25 @@ namespace sr = npsr::HWY_NAMESPACE;
7676
/// Default precision configrations for NumPy SIMD Routines
7777
// TODO: Allow user to configure the default precision at build time.
7878
namespace detail {
79-
using PreciseHigh = decltype(npsr::Precise{});
80-
using PreciseLow = decltype(npsr::Precise{npsr::kLowAccuracy});
79+
// `npsr::Precise` is a meta-RAII class responsible for managing the floating-point environment at runtime
80+
// and holding compile-time configurations for intrinsic floating-point behavior.
81+
// Default configurations are set to the highest IEEE compliance intent.
82+
// For more information, see:
83+
// https://github.com/numpy/numpy-simd-routines/blob/1acd14571b9f58072d1f7f17886945000d78b56a/npsr/precise.h#L93
84+
// However, when float64 SIMD is not available or disabled, we define a dummy type to avoid
85+
// clearing floating-point exceptions within its destructor.
8186
struct PresiceDummy {};
87+
8288
template <typename T>
8389
struct PreciseByType {};
8490
template <>
85-
struct PreciseByType<float> { using Type = PreciseLow; };
91+
struct PreciseByType<float> {
92+
using Type = decltype(npsr::Precise{npsr::kLowAccuracy});
93+
};
8694
template <>
8795
struct PreciseByType<double> {
8896
#if NPY_HWY_F64
89-
using Type = PreciseHigh;
97+
using Type = decltype(npsr::Precise{});
9098
#else
9199
// If float64 SIMD isn’t available, use a dummy type.
92100
// The scalar path will run, but `Type` must still be defined.
@@ -96,9 +104,12 @@ struct PreciseByType<double> {
96104
#endif
97105
};
98106
} // namespace detail
107+
108+
/// `npsr::Precise<...>` for the given type `T`.
109+
/// If `T` is `float`, it uses the low-accuracy configuration by default.
110+
/// If `T` is `double`, it uses the high-accuracy configuration by default(if float64 SIMD is available).
99111
template <typename T>
100112
using Precise = typename detail::PreciseByType<T>::Type;
101-
102113
#endif
103114
#include "simd.inc.hpp"
104115
} // namespace simd

numpy/_core/tests/test_umath.py

Lines changed: 48 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1649,6 +1649,54 @@ def test_sincos_overlaps(self, callable, dtype, stride):
16491649
callable(x[::stride], out=x[:M])
16501650
assert_equal(x[:M], y)
16511651

1652+
@pytest.mark.parametrize("func", [np.sin, np.cos])
1653+
@pytest.mark.parametrize("dtype", ["f", "d"])
1654+
def test_trig_facts(self, dtype, func):
1655+
atol, max_ulp = (1e-6, 3) if dtype == "f" else (1e-15, 1)
1656+
test_data = {
1657+
np.sin: [
1658+
(0, 0),
1659+
(np.pi, 0),
1660+
(-np.pi, 0),
1661+
(2 * np.pi, 0),
1662+
(-2 * np.pi, 0),
1663+
(np.pi / 2, 1),
1664+
(-np.pi / 2, -1),
1665+
(3 * np.pi / 2, -1),
1666+
(-3 * np.pi / 2, 1),
1667+
(np.pi / 6, 0.5),
1668+
(5 * np.pi / 6, 0.5),
1669+
(-np.pi / 6, -0.5),
1670+
(-5 * np.pi / 6, -0.5),
1671+
],
1672+
np.cos: [
1673+
(0, 1),
1674+
(np.pi, -1),
1675+
(-np.pi, -1),
1676+
(2 * np.pi, 1),
1677+
(-2 * np.pi, 1),
1678+
(np.pi / 2, 0),
1679+
(-np.pi / 2, 0),
1680+
(3 * np.pi / 2, 0),
1681+
(-3 * np.pi / 2, 0),
1682+
(np.pi / 3, 0.5),
1683+
(5 * np.pi / 3, 0.5),
1684+
(-np.pi / 3, 0.5),
1685+
(-5 * np.pi / 3, 0.5),
1686+
],
1687+
}[func]
1688+
for s_x, s_expected in test_data:
1689+
x = np.array([s_x] * (128 + 1), dtype=dtype)
1690+
result = func(x)
1691+
# Check that all results are identical, to catch issues with
1692+
# vectorization
1693+
assert np.all(result == result[0])
1694+
expected = np.array([s_expected] * (128 + 1), dtype=dtype)
1695+
if s_expected == 0:
1696+
assert_allclose(result, expected, atol=atol, rtol=0)
1697+
else:
1698+
assert_array_max_ulp(result, expected, maxulp=max_ulp)
1699+
16521700
@pytest.mark.parametrize('dt', ['e', 'f', 'd', 'g'])
16531701
def test_sqrt_values(self, dt):
16541702
with np.errstate(all='ignore'):

0 commit comments

Comments
 (0)