Skip to content

Reduce_mul implementation #1132

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
wants to merge 2 commits into
base: master
Choose a base branch
from
Open

Reduce_mul implementation #1132

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
8cf4a1f
added reduce_mul implementation for common arch and specialized for x…
emrys53 Jul 1, 2025
1ebf5f6
Merge branch 'xtensor-stack:master' into master
emrys53 Jul 1, 2025
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
14 changes: 14 additions & 0 deletions 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 @@ -141,6 +141,20 @@ namespace xsimd
return res;
}

// hmul
template <class A, class T, class /*=typename std::enable_if<std::is_integral<T>::value, void>::type*/>
XSIMD_INLINE T hmul(batch<T, A> const& self, requires_arch<common>) noexcept
Copy link
Contributor

Choose a reason for hiding this comment

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

shouldn't this just be the generic implementation of reduce_mul, and if so, why not using the current generic reduction implementation?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

This function is only used for integral types with sizeof(T) < 4. The reason I implemented like this was because for example when testing reduce_mul for uint8_t with AVX machine. What happens is that it has 32 elements which means general reduce function multiplies and swizzle the batch 5 times (32->16->8->4->2->1) and in the end it calls self.get(0) (which is just too inefficient anyway issue #1133). So I thought it would be more efficient to implement it this way. What it does with uint8_t in my implementation is that it first splits the batch into 2 equal sse4_2 batches multiplies them together and then calls hmul.

Copy link
Contributor Author

@emrys53 emrys53 Jul 1, 2025 •
edited
Loading

Choose a reason for hiding this comment

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

Also the swizzle function for uint8_t is just loading the batch to a buffer and storing the result back. So using that 5 times (in AVX) is just too inefficient compared to split_avx function which splits it using SIMD intrinsics.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

@serge-sans-paille you haven't responded back. If you are open to the idea of reduce_mul implementation. I can fix the code so that it passes all the tests. If you have feedbacks or recommendations to the code I have written, you can write it here and I can fix it after we discuss.

{
alignas(A::alignment()) T buffer[batch<T, A>::size];
self.store_aligned(buffer);
T res = 1;
for (T val : buffer)
{
res *= val;
}
return res;
}

// incr
template <class A, class T>
XSIMD_INLINE batch<T, A> incr(batch<T, A> const& self, requires_arch<common>) noexcept
Expand Down
2 changes: 2 additions & 0 deletions include/xsimd/arch/common/xsimd_common_details.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -74,6 +74,8 @@ namespace xsimd
XSIMD_INLINE batch<as_integer_t<T>, A> nearbyint_as_int(const batch<T, A>& x) noexcept;
template <class T, class A>
XSIMD_INLINE T reduce_add(batch<T, A> const&) noexcept;
template<class T, class A>
XSIMD_INLINE T reduce_mul(batch<T, A> const&) noexcept;
template <class T, class A>
XSIMD_INLINE batch<T, A> select(batch_bool<T, A> const&, batch<T, A> const&, batch<T, A> const&) noexcept;
template <class T, class A>
Expand Down
34 changes: 34 additions & 0 deletions include/xsimd/arch/common/xsimd_common_math.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2114,6 +2114,24 @@ namespace xsimd
}
};

template <class T, T N>
struct reverse
{
static constexpr T get(T i, T)
{
return i >= N ? i : N - i - 1;
}
};

template <class A, class T, as_unsigned_integer_t<T> Lvl>
XSIMD_INLINE std::complex<T> reduce_mul(batch<std::complex<T>, A> const& self, std::integral_constant<as_unsigned_integer_t<T>, Lvl>) noexcept
{
XSIMD_IF_CONSTEXPR(Lvl == 1) { return self.get(0); }
using index_type = as_unsigned_integer_t<T>;
using batch_type = batch<std::complex<T>, A>;
batch_type rev = swizzle(self, make_batch_constant<index_type, reverse<index_type, Lvl>, A>());
return reduce_mul(self * rev, std::integral_constant<index_type, Lvl / 2>());
}
template <class Op, class A, class T>
XSIMD_INLINE T reduce(Op, batch<T, A> const& self, std::integral_constant<unsigned, 1>) noexcept
{
Expand All @@ -2129,6 +2147,13 @@ namespace xsimd
}
}

// reduce_mul
template <class A, class T>
XSIMD_INLINE std::complex<T> reduce_mul(batch<std::complex<T>, A> const& self, requires_arch<common>) noexcept
{
return detail::reduce_mul(self, std::integral_constant<as_unsigned_integer_t<T>, batch<std::complex<T>, A>::size>());
}

// reduce_max
template <class A, class T>
XSIMD_INLINE T reduce_max(batch<T, A> const& self, requires_arch<common>) noexcept
Expand All @@ -2147,6 +2172,15 @@ namespace xsimd
self, std::integral_constant<unsigned, batch<T, A>::size>());
}

// reduce_mul
template <class A, class T>
XSIMD_INLINE T reduce_mul(batch<T, A> const& self, requires_arch<common>) noexcept
{
return detail::reduce([](batch<T, A> const& x, batch<T, A> const& y)
{ return mul(x, y); },
self, std::integral_constant<unsigned, batch<T, A>::size>());
}

// remainder
template <class A>
XSIMD_INLINE batch<float, A> remainder(batch<float, A> const& self, batch<float, A> const& other, requires_arch<common>) noexcept
Expand Down
9 changes: 9 additions & 0 deletions include/xsimd/arch/xsimd_avx.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1076,6 +1076,15 @@ namespace xsimd
__m128i low = _mm256_castsi256_si128(acc);
return reduce_min(batch<T, sse4_2>(low));
}
// reduce_mul
template <class A, class T, class = typename std::enable_if<std::is_integral<T>::value || std::is_same<T, float>::value || std::is_same<T, double>::value, void>::type>
XSIMD_INLINE T reduce_mul(batch<T, A> const& self, requires_arch<avx>)
{
typename batch<T, sse4_2>::register_type low, high;
detail::split_avx(self, low, high);
batch<T, sse4_2> blow(low), bhigh(high);
return reduce_mul(blow * bhigh);
}

// rsqrt
template <class A>
Expand Down
33 changes: 33 additions & 0 deletions include/xsimd/arch/xsimd_avx512f.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,11 @@
#include <type_traits>

#include "../types/xsimd_avx512f_register.hpp"
#include "xsimd/config/xsimd_inline.hpp"
#include "xsimd/types/xsimd_register.hpp"
#include "xsimd/xsimd.hpp"

#include <avx512fintrin.h>

namespace xsimd
{
Expand Down Expand Up @@ -1545,6 +1550,34 @@ namespace xsimd
__m256i low = _mm512_castsi512_si256(acc);
return reduce_min(batch<T, avx2>(low));
}
// reduce_mul
template <class A, class T, class = typename std::enable_if<std::is_integral<T>::value || std::is_same<T, float>::value || std::is_same<T, double>::value, void>::type>
XSIMD_INLINE T reduce_mul(batch<T, A> const& self, requires_arch<avx512f>) noexcept
{
XSIMD_IF_CONSTEXPR(std::is_same<T,float>::value)
{
return _mm512_reduce_mul_ps(self);
}
else XSIMD_IF_CONSTEXPR(std::is_same<T,double>::value)
{
return _mm512_reduce_mul_pd(self);
}
else XSIMD_IF_CONSTEXPR(sizeof(T) == 4)
{
return _mm512_reduce_mul_epi32(self);
}
else XSIMD_IF_CONSTEXPR(sizeof(T) == 8)
{
return _mm512_reduce_mul_epi64(self);
}
else
{
typename batch<T, avx2>::register_type low, high;
detail::split_avx512(self, low, high);
batch<T, avx2> blow(low), bhigh(high);
return reduce_mul(blow * bhigh);
}
}

// rsqrt
template <class A>
Expand Down
2 changes: 2 additions & 0 deletions include/xsimd/arch/xsimd_common_fwd.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,8 @@ namespace xsimd
XSIMD_INLINE batch<T, A> ssub(batch<T, A> const& self, batch<T, A> const& other, requires_arch<common>) noexcept;
template <class A, class T, class = typename std::enable_if<std::is_integral<T>::value, void>::type>
XSIMD_INLINE T hadd(batch<T, A> const& self, requires_arch<common>) noexcept;
template <class A, class T, class = typename std::enable_if<std::is_integral<T>::value, void>::type>
XSIMD_INLINE T hmul(batch<T, A> const& self, requires_arch<common>) noexcept;

}
}
Expand Down
37 changes: 37 additions & 0 deletions include/xsimd/arch/xsimd_sse2.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1252,6 +1252,43 @@ namespace xsimd
{
return _mm_cvtsd_f64(_mm_add_sd(self, _mm_unpackhi_pd(self, self)));
}
// reduce_mul
template <class A>
XSIMD_INLINE float reduce_mul(batch<float, A> const& self, requires_arch<sse2>) noexcept
{
__m128 tmp1 = _mm_mul_ps(self, _mm_movehl_ps(self, self));
__m128 tmp2 = _mm_mul_ps(self, _mm_shuffle_ps(tmp1, tmp1, 0x1));
return _mm_cvtss_f32(tmp2);
}

template <class A, class T, class = typename std::enable_if<std::is_integral<T>::value, void>::type>
XSIMD_INLINE T reduce_mul(batch<T, A> const& self, requires_arch<sse2>) noexcept
{
XSIMD_IF_CONSTEXPR(sizeof(T) == 4)
{
batch<T, A> tmp1 = _mm_shuffle_epi32(self, _MM_SHUFFLE(0,1,2,3));
tmp1 = tmp1 * self;
batch<T, A> tmp2 = _mm_unpackhi_epi32(tmp1, tmp1);
tmp2 = tmp2 * tmp1;
return _mm_cvtsi128_si32(tmp2);
}
else XSIMD_IF_CONSTEXPR(sizeof(T) == 8)
{
batch<T, A> tmp1 = _mm_unpackhi_epi64(self, self);
tmp1 = tmp1 * self;
return _mm_cvtsi128_si64(tmp1);
}
else
{
return hmul(self, common {});
}
}

template <class A>
XSIMD_INLINE double reduce_mul(batch<double,A> const& self, requires_arch<sse2>) noexcept
{
return _mm_cvtsd_f64(_mm_mul_pd(self, _mm_unpackhi_pd(self, self)));
}

// reduce_max
template <class A, class T, class _ = typename std::enable_if<(sizeof(T) <= 2), void>::type>
Expand Down
9 changes: 9 additions & 0 deletions include/xsimd/arch/xsimd_sse3.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -51,6 +51,15 @@ namespace xsimd
return _mm_cvtss_f32(tmp1);
}

// reduce_mul
template <class A>
XSIMD_INLINE float reduce_mul(batch<float, A> const& self, requires_arch<sse3>) noexcept
{
__m128 tmp1 = _mm_mul_ps(self, _mm_movehl_ps(self, self));
__m128 tmp2 = _mm_mul_ps(tmp1, _mm_movehdup_ps(tmp1));
return _mm_cvtss_f32(tmp2);
}

}

}
Expand Down
14 changes: 14 additions & 0 deletions include/xsimd/types/xsimd_api.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1864,6 +1864,20 @@ namespace xsimd
return kernel::reduce_min<A>(x, A {});
}

/**
* @ingroup batch_reducers
*
* Multiplication of all the scalars of the batch \c x.
* @param x batch involved in the reduction
* @return the result of the reduction.
*/
template <class T, class A>
XSIMD_INLINE T reduce_mul(batch<T, A> const& x) noexcept
{
detail::static_check_supported_config<T, A>();
return kernel::reduce_mul<A>(x, A {});
}

/**
* @ingroup batch_math
*
Expand Down
8 changes: 8 additions & 0 deletions test/test_batch.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -800,6 +800,14 @@ struct batch_test
INFO("reduce_min");
CHECK_SCALAR_EQ(res, expected);
}

// reduce_mul
{
value_type expected = std::accumulate(lhs.cbegin(), lhs.cend(), value_type{1}, std::multiplies<value_type>());
value_type res = reduce_mul(batch_lhs());
INFO("reduce_mul");
CHECK_SCALAR_EQ(res, expected);
}
}

template <size_t N>
Expand Down
7 changes: 7 additions & 0 deletions test/test_batch_complex.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -572,6 +572,13 @@ struct batch_complex_test
value_type res = reduce_add(batch_lhs());
CHECK_SCALAR_EQ(res, expected);
}

// reduce_mul
{
value_type expected = std::accumulate(lhs.cbegin(), lhs.cend(), value_type{1}, std::multiplies<value_type>());
value_type res = reduce_mul(batch_lhs());
CHECK_SCALAR_EQ(res, expected);
}
}

void test_fused_operations() const
Expand Down
Loading