From 93aa959b8adb53352282174da9d8f43818764229 Mon Sep 17 00:00:00 2001 From: Rasmus Munk Larsen <4643818-rmlarsen1@users.noreply.gitlab.com> Date: Thu, 12 Mar 2026 13:47:15 -0700 Subject: [PATCH] Add vectorization boundary tests for redux and visitor libeigen/eigen!2287 Co-authored-by: Rasmus Munk Larsen --- test/redux.cpp | 74 ++++++++++++++++++++++++++++++++++++++++++++++++ test/visitor.cpp | 55 +++++++++++++++++++++++++++++++++++ 2 files changed, 129 insertions(+) diff --git a/test/redux.cpp b/test/redux.cpp index af5c11d69..70b7261b1 100644 --- a/test/redux.cpp +++ b/test/redux.cpp @@ -193,6 +193,69 @@ void boolRedux(Index rows, Index cols) { } } +// Test reductions at sizes that hit vectorization boundaries in Redux.h: +// LinearVectorizedTraversal with 2-way unrolled packet loop, scalar pre/post loops. +template +void redux_vec_boundary() { + const Index PS = internal::packet_traits::size; + // Critical sizes: around packet multiples and at 2-way unroll boundaries + const Index sizes[] = {1, PS - 1, PS, PS + 1, 2 * PS - 1, 2 * PS, 2 * PS + 1, + 3 * PS, 3 * PS + 1, 4 * PS - 1, 4 * PS, 4 * PS + 1, 8 * PS, 8 * PS + 1}; + for (int si = 0; si < 14; ++si) { + const Index n = sizes[si]; + if (n <= 0) continue; + typedef Matrix Vec; + Vec v = Vec::Random(n); + // For prod, use values near 1 to avoid underflow (float) or overflow (int). + Vec v_for_prod = Vec::Ones(n) + Scalar(typename NumTraits::Real(0.2)) * v; + // Reference: scalar loops + Scalar ref_sum(0), ref_prod(1); + typename NumTraits::Real ref_min = numext::real(v(0)), ref_max = numext::real(v(0)); + for (Index k = 0; k < n; ++k) { + ref_sum += v(k); + ref_prod *= v_for_prod(k); + ref_min = (std::min)(ref_min, numext::real(v(k))); + ref_max = (std::max)(ref_max, numext::real(v(k))); + } + VERIFY_IS_APPROX(v.sum(), ref_sum); + VERIFY_IS_APPROX(v_for_prod.prod(), ref_prod); + VERIFY_IS_APPROX(v.real().minCoeff(), ref_min); + VERIFY_IS_APPROX(v.real().maxCoeff(), ref_max); + } +} + +// Test reductions on strided (non-contiguous) mapped data. +// This exercises SliceVectorizedTraversal or DefaultTraversal in Redux.h +// depending on stride and packet size. +template +void redux_strided() { + const Index n = 64; + typedef Matrix Vec; + Vec data = Vec::Random(2 * n); + // Map with inner stride of 2 — every other element + Map> strided(data.data(), n); + Scalar ref_sum(0); + typename NumTraits::Real ref_min = numext::real(strided(0)), ref_max = numext::real(strided(0)); + for (Index k = 0; k < n; ++k) { + ref_sum += strided(k); + ref_min = (std::min)(ref_min, numext::real(strided(k))); + ref_max = (std::max)(ref_max, numext::real(strided(k))); + } + VERIFY_IS_APPROX(strided.sum(), ref_sum); + VERIFY_IS_APPROX(strided.real().minCoeff(), ref_min); + VERIFY_IS_APPROX(strided.real().maxCoeff(), ref_max); + + // Also test reduction on a non-contiguous matrix block (SliceVectorizedTraversal) + typedef Matrix Mat; + Mat m = Mat::Random(16, 16); + for (Index bsz = 1; bsz <= 8; bsz *= 2) { + Scalar block_sum(0); + for (Index j = 0; j < bsz; ++j) + for (Index i = 0; i < bsz; ++i) block_sum += m(1 + i, 1 + j); + VERIFY_IS_APPROX(m.block(1, 1, bsz, bsz).sum(), block_sum); + } +} + EIGEN_DECLARE_TEST(redux) { // the max size cannot be too large, otherwise reduxion operations obviously generate large errors. int maxsize = (std::min)(100, EIGEN_TEST_MAX_SIZE); @@ -248,4 +311,15 @@ EIGEN_DECLARE_TEST(redux) { CALL_SUBTEST_11(boolRedux(4, 4)); CALL_SUBTEST_11(boolRedux(7, 13)); CALL_SUBTEST_11(boolRedux(63, 63)); + + // Vectorization boundary sizes — deterministic, run once. + // Integer types are excluded: full-range random ints overflow in sum/prod (UB). + // Integer reductions are already tested by matrixRedux/vectorRedux with clamped values. + CALL_SUBTEST_12(redux_vec_boundary()); + CALL_SUBTEST_12(redux_vec_boundary()); + + // Strided (non-contiguous) reductions. + CALL_SUBTEST_13(redux_strided()); + CALL_SUBTEST_13(redux_strided()); + CALL_SUBTEST_13(redux_strided>()); } diff --git a/test/visitor.cpp b/test/visitor.cpp index 27ca18415..0aa1f6520 100644 --- a/test/visitor.cpp +++ b/test/visitor.cpp @@ -315,6 +315,56 @@ void checkOptimalTraversal() { VERIFY(j == 0); } +// Test minCoeff/maxCoeff at vectorization boundary sizes. +// Visitor uses LinearVectorizedTraversal with packet-based min/max, +// so we test at sizes around packet multiples. +template +void visitor_vec_boundary() { + const Index PS = internal::packet_traits::size; + const Index sizes[] = {1, 2, 3, PS - 1, PS, PS + 1, 2 * PS - 1, 2 * PS, 2 * PS + 1, 4 * PS, 4 * PS + 1}; + for (int si = 0; si < 11; ++si) { + const Index n = sizes[si]; + if (n <= 0) continue; + typedef Matrix Vec; + Vec v = Vec::Random(n); + // Ensure all elements are distinct. + for (Index i = 0; i < n; ++i) + for (Index j = 0; j < i; ++j) + while (numext::equal_strict(v(i), v(j))) v(i) = internal::random(); + // Reference + Scalar ref_min = v(0), ref_max = v(0); + Index ref_minidx = 0, ref_maxidx = 0; + for (Index k = 0; k < n; ++k) { + if (v(k) < ref_min) { + ref_min = v(k); + ref_minidx = k; + } + if (v(k) > ref_max) { + ref_max = v(k); + ref_maxidx = k; + } + } + Index eigen_minidx, eigen_maxidx; + VERIFY_IS_APPROX(v.minCoeff(&eigen_minidx), ref_min); + VERIFY_IS_APPROX(v.maxCoeff(&eigen_maxidx), ref_max); + VERIFY(eigen_minidx == ref_minidx); + VERIFY(eigen_maxidx == ref_maxidx); + + // Also test matrix form at this size (exercises different inner/outer sizes). + if (n >= 2) { + typedef Matrix Mat; + // Test as n×1 and 1×n (different inner sizes for visitor traversal). + Mat mc = v; + Mat mr = v.transpose(); + Index ri, ci; + VERIFY_IS_APPROX(mc.minCoeff(&ri, &ci), ref_min); + VERIFY(ri == ref_minidx && ci == 0); + VERIFY_IS_APPROX(mr.minCoeff(&ri, &ci), ref_min); + VERIFY(ri == 0 && ci == ref_minidx); + } + } +} + EIGEN_DECLARE_TEST(visitor) { for (int i = 0; i < g_repeat; i++) { CALL_SUBTEST_1(matrixVisitor(Matrix())); @@ -332,4 +382,9 @@ EIGEN_DECLARE_TEST(visitor) { CALL_SUBTEST_10(vectorVisitor(VectorXf(33))); } CALL_SUBTEST_11(checkOptimalTraversal()); + + // Vectorization boundary sizes — deterministic, run once. + CALL_SUBTEST_12(visitor_vec_boundary()); + CALL_SUBTEST_12(visitor_vec_boundary()); + CALL_SUBTEST_12(visitor_vec_boundary()); }