Skip to content

Fix saturating subtraction when the subtrahend is the type minimum #1380

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
vsaraikin wants to merge 1 commit into
base: master
Choose a base branch
from
Open

Fix saturating subtraction when the subtrahend is the type minimum #1380

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 6 additions & 1 deletion include/xsimd/arch/common/xsimd_common_arithmetic.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -398,7 +398,12 @@ namespace xsimd
{
if (std::is_signed<T>::value)
{
return sadd(self, -other);
// Saturating self - other, mirroring the signed sadd above.
// sadd(self, -other) is wrong when other == numeric_limits<T>::min(),
// since -other is not representable.
auto self_underflow_branch = max(std::numeric_limits<T>::min() + other, self);
auto self_overflow_branch = min(std::numeric_limits<T>::max() + other, self);
return select(other >= 0, self_underflow_branch, self_overflow_branch) - other;
}
else
{
Expand Down
5 changes: 4 additions & 1 deletion include/xsimd/arch/xsimd_avx.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1711,7 +1711,10 @@ namespace xsimd
{
if (std::is_signed<T>::value)
{
return sadd(self, -other);

@DiamonDinoia DiamonDinoia Jul 14, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

why directly to commmon here? are there no intrinsics for it?

@serge-sans-paille serge-sans-paille Jul 15, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

indeed there's

 _mm256_subs_epi16
 _mm256_subs_epu16
 _mm256_subs_epi8
 _mm256_subs_epu8

@vsaraikin vsaraikin Jul 15, 2026

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You're right, I've since dropped the common detour and implemented it inline here, mirroring the sadd right above

And yeah, as Serge noted, AVX2 has subs_epi8/epi16, those stay on the native path in avx2::ssub. This overload only handles the 32/64-bit case, which has no saturating-sub intrinsic

@DiamonDinoia DiamonDinoia Jul 15, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I mean you can also emulate saturating using sub + select no? Would that be faster that a scalar loop?

auto mask = (other >> (8 * sizeof(T) - 1));
auto self_overflow_branch = min(std::numeric_limits<T>::max() + other, self);
auto self_underflow_branch = max(std::numeric_limits<T>::min() + other, self);
return select(batch_bool<T, A>(mask.data), self_overflow_branch, self_underflow_branch) - other;
}
else
{
Expand Down
15 changes: 13 additions & 2 deletions include/xsimd/arch/xsimd_avx512f.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2419,9 +2419,20 @@ namespace xsimd
template <class A, class T, class = std::enable_if_t<std::is_integral<T>::value>>
XSIMD_INLINE batch<T, A> ssub(batch<T, A> const& self, batch<T, A> const& other, requires_arch<avx512f>) noexcept
{
if (std::is_signed<T>::value)
// Saturating sub for 8/16-bit integers needs AVX512BW; on AVX512F
// fall back to the AVX2 implementation on each 256-bit half.
XSIMD_IF_CONSTEXPR(sizeof(T) == 1 || sizeof(T) == 2)
{
return detail::fwd_to_avx([](__m256i s, __m256i o) noexcept
{ return ssub(batch<T, avx2>(s), batch<T, avx2>(o), avx2 {}); },
self, other);
}
else if (std::is_signed<T>::value)
{
return sadd(self, -other);

@DiamonDinoia DiamonDinoia Jul 14, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

same here.

@serge-sans-paille serge-sans-paille Jul 15, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

it comes with AVX512BW. Yet falling back to avx for [u]int[8[16]_t is probably better.

@vsaraikin vsaraikin Jul 15, 2026

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

happy to make the 8/16 path fall back to avx instead if you'd prefer

auto mask = other < 0;
auto self_overflow_branch = min(std::numeric_limits<T>::max() + other, self);
auto self_underflow_branch = max(std::numeric_limits<T>::min() + other, self);
return select(mask, self_overflow_branch, self_underflow_branch) - other;
}
else
{
Expand Down
16 changes: 15 additions & 1 deletion include/xsimd/arch/xsimd_scalar.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -720,7 +720,21 @@ namespace xsimd
{
if (std::numeric_limits<T>::is_signed)
{
return sadd(lhs, (T)-rhs);
// Compute the saturating difference directly. Using sadd(lhs, -rhs)
// is wrong when rhs == numeric_limits<T>::min(), because -rhs is not
// representable.
if ((rhs < 0) && (lhs > std::numeric_limits<T>::max() + rhs))
{
return std::numeric_limits<T>::max();
}
else if ((rhs > 0) && (lhs < std::numeric_limits<T>::lowest() + rhs))
{
return std::numeric_limits<T>::lowest();

@DiamonDinoia DiamonDinoia Jul 14, 2026

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

why max vs lowest?

@vsaraikin vsaraikin Jul 15, 2026

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

matching the sadd just above it, which uses max() / lowest()

lowest() (not min()) is deliberate so it's also correct for floating-point T, where min() is the smallest positive value

for integers they're identical

}
else
{
return lhs - rhs;
}
}
else
{
Expand Down
19 changes: 19 additions & 0 deletions test/test_xsimd_api.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -500,6 +500,25 @@ TEST_CASE_TEMPLATE("[xsimd api | integral types functions]", B, INTEGRAL_TYPES)
}
}

// Subtracting the type minimum must saturate, not wrap. Previously this was
// computed as sadd(x, -min), and -min is not representable.
TEST_CASE_TEMPLATE("[xsimd api | ssub at type minimum]", B, INTEGRAL_TYPES)
{
using value_type = typename scalar_type<B>::type;
value_type lo = std::numeric_limits<value_type>::min();
value_type hi = std::numeric_limits<value_type>::max();
if (std::numeric_limits<value_type>::is_signed)
{
CHECK_EQ(extract(xsimd::ssub(B(value_type(0)), B(lo))), hi);
CHECK_EQ(extract(xsimd::ssub(B(value_type(1)), B(lo))), hi);
}
else
{
// lo == 0 for unsigned types.
CHECK_EQ(extract(xsimd::ssub(B(value_type(5)), B(lo))), value_type(5));
}
}

/*
* Functions that apply on floating points types only
*/
Expand Down
Loading