From 446d3fa059fce8141ed30d985dc8a4fa7e26a719 Mon Sep 17 00:00:00 2001 From: "KOKO\\Mimi" Date: Wed, 14 Oct 2026 21:10:27 +0900 Subject: [PATCH] fix(hex8): isolate AVX2 packet code from baseline COMDATs --- src/fesa/elements/hex8_high_low.cpp | 24 ++- src/fesa/elements/hex8_pair_simd.h | 48 ++++- src/fesa/elements/hex8_pair_simd_avx2.cpp | 130 ++++++++---- src/fesa/elements/hex8_pair_simd_dispatch.cpp | 76 ++++++- tests/unit/elements/hex8_high_low_test.cpp | 195 ++++++++++++------ 5 files changed, 344 insertions(+), 129 deletions(-) diff --git a/src/fesa/elements/hex8_high_low.cpp b/src/fesa/elements/hex8_high_low.cpp index 6373dcd..ccb873f 100644 --- a/src/fesa/elements/hex8_high_low.cpp +++ b/src/fesa/elements/hex8_high_low.cpp @@ -271,25 +271,29 @@ Result Hex8Element::EvaluateHighLow( if (use_packets) for (std::size_t i = 0; i < 3; ++i) for (std::size_t first = 0; first < 24; first += 4) { - hex8_internal::MaterialPacketInput input; - input.weight = HighLow{point.weight}; + hex8_internal::MaterialPacketInput input{}; + const auto raw_pair = [](HighLow pair) { + return hex8_internal::PacketPair{pair.high, pair.low}; + }; + input.weight = raw_pair(HighLow{point.weight}); const auto row = 3U * a + i; for (std::size_t lane = 0; lane < 4; ++lane) { const auto column = first + lane; const auto b = column / 3U; const auto k = column % 3U; - input.deviatoric_left[lane] = deviatoric_fg[a][i]; - input.deviatoric_right[lane] = fg[b][k]; - input.shear_left[lane] = shear_fg[b][i]; - input.shear_right[lane] = fg[a][k]; - input.tensor_left[lane] = shear_fft[i][k]; - input.gradient_dot[lane] = dots[b]; - input.accumulated[lane] = material[row][column]; + input.deviatoric_left.values[lane] = raw_pair(deviatoric_fg[a][i]); + input.deviatoric_right.values[lane] = raw_pair(fg[b][k]); + input.shear_left.values[lane] = raw_pair(shear_fg[b][i]); + input.shear_right.values[lane] = raw_pair(fg[a][k]); + input.tensor_left.values[lane] = raw_pair(shear_fft[i][k]); + input.gradient_dot.values[lane] = raw_pair(dots[b]); + input.accumulated.values[lane] = raw_pair(material[row][column]); } hex8_internal::PairPacket output; packet_evaluator(input, &output); for (std::size_t lane = 0; lane < 4; ++lane) - material[row][first + lane] = output[lane]; + material[row][first + lane] = + {output.values[lane].high, output.values[lane].low}; } } } diff --git a/src/fesa/elements/hex8_pair_simd.h b/src/fesa/elements/hex8_pair_simd.h index 73c2f53..6a0fbbd 100644 --- a/src/fesa/elements/hex8_pair_simd.h +++ b/src/fesa/elements/hex8_pair_simd.h @@ -1,28 +1,60 @@ #ifndef FESA_SRC_FESA_ELEMENTS_HEX8_PAIR_SIMD_H_ #define FESA_SRC_FESA_ELEMENTS_HEX8_PAIR_SIMD_H_ -#include - -#include "../math/high_low_arithmetic.h" +#include namespace fesa::hex8_internal { -using PairPacket = std::array; +/// @brief Trivial raw carrier; neither component is normalized or rounded here. +struct PacketPair { + double high; + double low; +}; +/// @brief Four independent lanes with no executable constructor or operator. +struct PairPacket { + PacketPair values[4]; +}; struct MaterialPacketInput { - PairPacket deviatoric_left, deviatoric_right; - PairPacket shear_left, shear_right; - PairPacket tensor_left, gradient_dot; + PairPacket deviatoric_left; + PairPacket deviatoric_right; + PairPacket shear_left; + PairPacket shear_right; + PairPacket tensor_left; + PairPacket gradient_dot; PairPacket accumulated; - math_internal::HighLow weight; + PacketPair weight; }; using MaterialPacketEvaluator = void (*)(const MaterialPacketInput&, PairPacket*); +/// @brief Evaluates each lane's original checked material graph in order 0..3. +/// @note Invalid arithmetic remains a raw invalid pair for the caller's gates. void EvaluateMaterialPacketScalar(const MaterialPacketInput& input, PairPacket* output); +/// @brief Evaluates independent material lanes with the checked vector graph. +/// @pre The selector has admitted the CPU/OS ISA and current masked FP state. +/// @note No lane reduction or model failure is introduced; special packets bridge +/// to original scalar arithmetic, preserving components and operation order. void EvaluateMaterialPacketAvx2(const MaterialPacketInput& input, PairPacket* output); +/// @brief Baseline-ISA checked addition bridge in lane order 0..3. +/// @note pair_shortcuts chooses the original local AddPair or general Add graph; +/// invalid pairs are preserved, with no FP-state mutation or normalization here. +void ScalarPacketAdd(const PairPacket& left, const PairPacket& right, + bool pair_shortcuts, PairPacket* output); +/// @brief Baseline-ISA original checked MultiplyPair bridge in lane order 0..3. +/// @note Raw invalid outputs remain invalid; no component is rounded away. +void ScalarPacketMultiply(const PairPacket& left, const PairPacket& right, + PairPacket* output); +/// @brief Baseline-ISA original checked TwoProduct bridge in lane order 0..3. +/// @note Keeps finite/zero/underflow/FMA checks and the caller's FP state. +void ScalarPacketTwoProduct(const double left[4], const double right[4], + PairPacket* output); /// @brief Selects a private backend from CPU/OS and this thread's FP state. /// @note Never mutates FP state or caches a thread's environment decision. MaterialPacketEvaluator SelectMaterialPacketEvaluator(); +/// @brief Pure instruction/OS-state predicate; does not execute ISA probes. +bool SupportsMaterialPacketInstructions(std::uint32_t leaf1_ecx, + std::uint64_t xcr0, + std::uint32_t leaf7_ebx); } // namespace fesa::hex8_internal #endif // FESA_SRC_FESA_ELEMENTS_HEX8_PAIR_SIMD_H_ diff --git a/src/fesa/elements/hex8_pair_simd_avx2.cpp b/src/fesa/elements/hex8_pair_simd_avx2.cpp index bc4bef9..ffbebef 100644 --- a/src/fesa/elements/hex8_pair_simd_avx2.cpp +++ b/src/fesa/elements/hex8_pair_simd_avx2.cpp @@ -1,17 +1,61 @@ #include "hex8_pair_simd.h" -#include "hex8_pair_arithmetic.h" - #if defined(FESA_HEX8_PAIR_AVX2) && defined(_MSC_VER) && defined(_M_X64) +#include #include +#include +#include namespace fesa::hex8_internal { namespace { -struct VectorPair { __m256d high, low; }; +struct VectorPair { + __m256d high; + __m256d low; +}; +// TU-local predicates mirror the checked scalar premises without publishing +// AVX2-compiled common math or HighLow COMDAT definitions to baseline callers. +/// @brief Reads binary64 fields without FP classification or comparison. +std::uint64_t ScalarBits(double value) { + static_assert(sizeof(double) == sizeof(std::uint64_t) && + std::numeric_limits::is_iec559 && + std::numeric_limits::digits == 53); + std::uint64_t bits = 0; + std::memcpy(&bits, &value, sizeof(bits)); + return bits; +} +/// @brief Bitwise finite gate before any scalar FP comparison. +bool FiniteScalar(double value) { + return (ScalarBits(value) & 0x7ff0000000000000ULL) != + 0x7ff0000000000000ULL; +} +/// @brief Preserves finite, canonical-zero and RN64 normalization gates. +bool NormalizedScalar(PacketPair value) { + if (!FiniteScalar(value.high) || !FiniteScalar(value.low)) return false; + if (value.high == 0 && value.low == 0) + return (ScalarBits(value.high) >> 63U) == 0 && + (ScalarBits(value.low) >> 63U) == 0; + return value.high + value.low == value.high; +} +/// @brief Exact least-bit exponent; caller has checked finite/nonzero first. +int LowestBitExponentScalar(double value) { + const auto bits = ScalarBits(value); + const auto encoded_exponent = (bits >> 52U) & 0x7ffU; + auto significand = bits & 0xfffffffffffffULL; + const int exponent = encoded_exponent == 0U + ? -1074 + : static_cast(encoded_exponent) - 1075; + if (encoded_exponent != 0U) significand |= 1ULL << 52U; + unsigned long trailing = 0; + _BitScanForward64(&trailing, significand); + return exponent + static_cast(trailing); +} + +/// @brief Selects each lane without evaluating additional FP arithmetic. __m256d Select(__m256d mask, __m256d selected, __m256d other) { return _mm256_blendv_pd(other, selected, mask); } +/// @brief Produces independent finite masks without signaling NaN compares. __m256d Finite(__m256d value) { const auto exponent = _mm256_set1_epi64x(0x7ff0000000000000LL); const auto nonfinite = _mm256_cmpeq_epi64( @@ -19,30 +63,36 @@ __m256d Finite(__m256d value) { return _mm256_castsi256_pd(_mm256_xor_si256( nonfinite, _mm256_set1_epi64x(-1))); } +/// @brief Recognizes both signed zero encodings using integer fields. __m256d Zero(__m256d value) { return _mm256_castsi256_pd(_mm256_cmpeq_epi64( _mm256_and_si256(_mm256_castpd_si256(value), _mm256_set1_epi64x(0x7fffffffffffffffLL)), _mm256_setzero_si256())); } +/// @brief Preserves valid lanes and encodes invalid lanes as quiet NaN/+0. VectorPair InvalidOutside(VectorPair value, __m256d valid) { - return {Select(valid, value.high, _mm256_set1_pd( - std::numeric_limits::quiet_NaN())), + return {Select(valid, value.high, _mm256_castsi256_pd( + _mm256_set1_epi64x(0x7ff8000000000000LL))), Select(valid, value.low, _mm256_setzero_pd())}; } +/// @brief Loads raw high/low fields in lane order 0..3, with no normalization. VectorPair Load(const PairPacket& value) { - return {_mm256_setr_pd(value[0].high, value[1].high, - value[2].high, value[3].high), - _mm256_setr_pd(value[0].low, value[1].low, - value[2].low, value[3].low)}; + return {_mm256_setr_pd(value.values[0].high, value.values[1].high, + value.values[2].high, value.values[3].high), + _mm256_setr_pd(value.values[0].low, value.values[1].low, + value.values[2].low, value.values[3].low)}; } +/// @brief Stores raw high/low fields without rounding or lane reduction. PairPacket Store(VectorPair value) { - double high[4], low[4]; + double high[4]; + double low[4]; _mm256_storeu_pd(high, value.high); _mm256_storeu_pd(low, value.low); return {{{high[0], low[0]}, {high[1], low[1]}, {high[2], low[2]}, {high[3], low[3]}}}; } +/// @brief Retains the original per-lane subtraction/addition association and gates. VectorPair TwoSum(__m256d left, __m256d right) { const auto sum = _mm256_add_pd(left, right); const auto bridge = _mm256_sub_pd(sum, left); @@ -55,29 +105,32 @@ VectorPair TwoSum(__m256d left, __m256d right) { return InvalidOutside({Select(exact_zero, _mm256_setzero_pd(), sum), Select(exact_zero, _mm256_setzero_pd(), error)}, valid); } +/// @brief Routes exceptional packets to baseline checked Add in lane order. VectorPair ScalarAdd(VectorPair left, VectorPair right, bool pair_shortcuts) { - const auto a = Store(left), b = Store(right); + const auto a = Store(left); + const auto b = Store(right); PairPacket result; - for (std::size_t lane = 0; lane < result.size(); ++lane) - result[lane] = pair_shortcuts ? AddPair(a[lane], b[lane]) - : math_internal::Add(a[lane], b[lane]); + ScalarPacketAdd(a, b, pair_shortcuts, &result); return Load(result); } +/// @brief Routes exceptional packets to baseline checked MultiplyPair. VectorPair ScalarMultiply(const PairPacket& left, const PairPacket& right) { PairPacket result; - for (std::size_t lane = 0; lane < result.size(); ++lane) - result[lane] = MultiplyPair(left[lane], right[lane]); + ScalarPacketMultiply(left, right, &result); return Load(result); } +/// @brief Retains normalized/zero branches and the checked per-lane Add graph. VectorPair Add(VectorPair left, VectorPair right, bool pair_shortcuts) { - const auto a = Store(left), b = Store(right); - int left_zero = 0, right_zero = 0; - for (std::size_t lane = 0; lane < a.size(); ++lane) { - if (!math_internal::IsNormalized(a[lane]) || - !math_internal::IsNormalized(b[lane])) + const auto a = Store(left); + const auto b = Store(right); + int left_zero = 0; + int right_zero = 0; + for (std::size_t lane = 0; lane < 4; ++lane) { + if (!NormalizedScalar(a.values[lane]) || + !NormalizedScalar(b.values[lane])) return ScalarAdd(left, right, pair_shortcuts); - if (a[lane].high == 0 && a[lane].low == 0) left_zero |= 1 << lane; - if (b[lane].high == 0 && b[lane].low == 0) right_zero |= 1 << lane; + if (a.values[lane].high == 0 && a.values[lane].low == 0) left_zero |= 1 << lane; + if (b.values[lane].high == 0 && b.values[lane].low == 0) right_zero |= 1 << lane; } if (pair_shortcuts) { if (left_zero == 15) @@ -103,22 +156,24 @@ VectorPair Add(VectorPair left, VectorPair right, bool pair_shortcuts) { Select(valid, joined.high, _mm256_setzero_pd()), Select(valid, remainder, _mm256_setzero_pd())), valid); } +/// @brief Checks finite/zero/tail premises before products; neutralizes rejected +/// lanes before the explicit residual FMA and final TwoSum. VectorPair Product(__m256d left, __m256d right) { - double a[4], b[4]; + double a[4]; + double b[4]; _mm256_storeu_pd(a, left); _mm256_storeu_pd(b, right); bool special = false; for (std::size_t lane = 0; lane < 4; ++lane) - if (!math_internal::IsFiniteBinary64(a[lane]) || - !math_internal::IsFiniteBinary64(b[lane]) || + if (!FiniteScalar(a[lane]) || + !FiniteScalar(b[lane]) || a[lane] == 0 || b[lane] == 0 || - math_internal::LowestBitExponent(a[lane]) + - math_internal::LowestBitExponent(b[lane]) < -1074) + LowestBitExponentScalar(a[lane]) + + LowestBitExponentScalar(b[lane]) < -1074) special = true; if (special) { PairPacket result; - for (std::size_t lane = 0; lane < 4; ++lane) - result[lane] = math_internal::TwoProduct(a[lane], b[lane]); + ScalarPacketTwoProduct(a, b, &result); return Load(result); } const auto product = _mm256_mul_pd(left, right); @@ -132,15 +187,18 @@ VectorPair Product(__m256d left, __m256d right) { Select(valid, product, _mm256_setzero_pd()), Select(valid, error, _mm256_setzero_pd())), valid); } +/// @brief Uses uniform low-zero branches with original four-product association; +/// mixed or invalid packets preserve the baseline checked path. VectorPair Multiply(VectorPair left, VectorPair right) { - const auto a = Store(left), b = Store(right); + const auto a = Store(left); + const auto b = Store(right); int branch = -1; - for (std::size_t lane = 0; lane < a.size(); ++lane) { - const bool normalized = math_internal::IsNormalized(a[lane]) && - math_internal::IsNormalized(b[lane]); + for (std::size_t lane = 0; lane < 4; ++lane) { + const bool normalized = NormalizedScalar(a.values[lane]) && + NormalizedScalar(b.values[lane]); if (!normalized) return ScalarMultiply(a, b); - const int shape = (a[lane].low != 0 ? 1 : 0) | - (b[lane].low != 0 ? 2 : 0); + const int shape = (a.values[lane].low != 0 ? 1 : 0) | + (b.values[lane].low != 0 ? 2 : 0); if (branch != -1 && shape != branch) return ScalarMultiply(a, b); branch = shape; } diff --git a/src/fesa/elements/hex8_pair_simd_dispatch.cpp b/src/fesa/elements/hex8_pair_simd_dispatch.cpp index facc870..a4838d7 100644 --- a/src/fesa/elements/hex8_pair_simd_dispatch.cpp +++ b/src/fesa/elements/hex8_pair_simd_dispatch.cpp @@ -9,18 +9,70 @@ #endif namespace fesa::hex8_internal { +#if defined(_MSC_VER) +#define FESA_HEX8_BASELINE_NOINLINE __declspec(noinline) +#else +#define FESA_HEX8_BASELINE_NOINLINE +#endif +FESA_HEX8_BASELINE_NOINLINE +void ScalarPacketAdd(const PairPacket& left, const PairPacket& right, + bool pair_shortcuts, PairPacket* output) { + for (std::size_t lane = 0; lane < 4; ++lane) { + const math_internal::HighLow a{left.values[lane].high, left.values[lane].low}; + const math_internal::HighLow b{right.values[lane].high, right.values[lane].low}; + const auto value = pair_shortcuts ? AddPair(a, b) : math_internal::Add(a, b); + output->values[lane] = {value.high, value.low}; + } +} +FESA_HEX8_BASELINE_NOINLINE +void ScalarPacketMultiply(const PairPacket& left, const PairPacket& right, + PairPacket* output) { + for (std::size_t lane = 0; lane < 4; ++lane) { + const auto value = MultiplyPair( + {left.values[lane].high, left.values[lane].low}, + {right.values[lane].high, right.values[lane].low}); + output->values[lane] = {value.high, value.low}; + } +} +FESA_HEX8_BASELINE_NOINLINE +void ScalarPacketTwoProduct(const double left[4], const double right[4], + PairPacket* output) { + for (std::size_t lane = 0; lane < 4; ++lane) { + const auto value = math_internal::TwoProduct(left[lane], right[lane]); + output->values[lane] = {value.high, value.low}; + } +} +#undef FESA_HEX8_BASELINE_NOINLINE + +bool SupportsMaterialPacketInstructions(std::uint32_t leaf1_ecx, + std::uint64_t xcr0, + std::uint32_t leaf7_ebx) { + constexpr std::uint32_t required_leaf1 = + (1U << 12) | (1U << 26) | (1U << 27) | (1U << 28); + // MSVC /arch:AVX2 permits BMI codegen; ANDN is present in this kernel. + constexpr std::uint32_t required_leaf7 = + (1U << 3) | (1U << 5) | (1U << 8); + return (leaf1_ecx & required_leaf1) == required_leaf1 && + (xcr0 & 6U) == 6U && + (leaf7_ebx & required_leaf7) == required_leaf7; +} + void EvaluateMaterialPacketScalar(const MaterialPacketInput& input, PairPacket* output) { - for (std::size_t lane = 0; lane < output->size(); ++lane) { - const auto first = MultiplyPair(input.deviatoric_left[lane], - input.deviatoric_right[lane]); - const auto second = MultiplyPair(input.shear_left[lane], - input.shear_right[lane]); - const auto third = MultiplyPair(input.tensor_left[lane], - input.gradient_dot[lane]); + const auto value = [](PacketPair pair) { + return math_internal::HighLow{pair.high, pair.low}; + }; + for (std::size_t lane = 0; lane < 4; ++lane) { + const auto first = MultiplyPair(value(input.deviatoric_left.values[lane]), + value(input.deviatoric_right.values[lane])); + const auto second = MultiplyPair(value(input.shear_left.values[lane]), + value(input.shear_right.values[lane])); + const auto third = MultiplyPair(value(input.tensor_left.values[lane]), + value(input.gradient_dot.values[lane])); const auto material = AddPair(AddPair(first, second), third); - (*output)[lane] = AddPair(input.accumulated[lane], - MultiplyPair(material, input.weight)); + const auto accumulated = AddPair(value(input.accumulated.values[lane]), + MultiplyPair(material, value(input.weight))); + output->values[lane] = {accumulated.high, accumulated.low}; } } @@ -33,9 +85,11 @@ MaterialPacketEvaluator SelectMaterialPacketEvaluator() { __cpuidex(info, 1, 0); constexpr int required = (1 << 12) | (1 << 26) | (1 << 27) | (1 << 28); if ((info[2] & required) != required) return false; - if ((_xgetbv(0) & 6U) != 6U) return false; + const auto leaf1 = static_cast(info[2]); + const auto xcr0 = _xgetbv(0); __cpuidex(info, 7, 0); - return (info[1] & (1 << 5)) != 0; + return SupportsMaterialPacketInstructions( + leaf1, xcr0, static_cast(info[1])); }(); if (!hardware_supported || !math_internal::ArithmeticEnvironmentSupported()) return EvaluateMaterialPacketScalar; diff --git a/tests/unit/elements/hex8_high_low_test.cpp b/tests/unit/elements/hex8_high_low_test.cpp index 5f55a4d..a377275 100644 --- a/tests/unit/elements/hex8_high_low_test.cpp +++ b/tests/unit/elements/hex8_high_low_test.cpp @@ -12,6 +12,7 @@ #include #include #include +#include #include #include @@ -85,27 +86,31 @@ TEST(Hex8HighLow, ScalarCrossProductShortcutMatchesOriginalBits) { } } +fesa::hex8_internal::PacketPair RawPacketPair(HighLow value) { + return {value.high, value.low}; +} + fesa::hex8_internal::PairPacket FrozenMaterialPacket( const fesa::hex8_internal::MaterialPacketInput& input) { using Original = fesa::high_low_original::HighLow; - const auto original = [](HighLow value) { + const auto original = [](fesa::hex8_internal::PacketPair value) { return Original{value.high, value.low}; }; fesa::hex8_internal::PairPacket result; for (std::size_t lane = 0; lane < 4; ++lane) { const auto first = fesa::high_low_original::Multiply( - original(input.deviatoric_left[lane]), - original(input.deviatoric_right[lane])); + original(input.deviatoric_left.values[lane]), + original(input.deviatoric_right.values[lane])); const auto second = fesa::high_low_original::Multiply( - original(input.shear_left[lane]), original(input.shear_right[lane])); + original(input.shear_left.values[lane]), original(input.shear_right.values[lane])); const auto third = fesa::high_low_original::Multiply( - original(input.tensor_left[lane]), original(input.gradient_dot[lane])); + original(input.tensor_left.values[lane]), original(input.gradient_dot.values[lane])); const auto value = fesa::high_low_original::Add( fesa::high_low_original::Add(first, second), third); const auto accumulated = fesa::high_low_original::Add( - original(input.accumulated[lane]), + original(input.accumulated.values[lane]), fesa::high_low_original::Multiply(value, original(input.weight))); - result[lane] = {accumulated.high, accumulated.low}; + result.values[lane] = {accumulated.high, accumulated.low}; } return result; } @@ -119,10 +124,10 @@ bool ExpectFrozenMaterialPacket( evaluator(input, &selected); for (std::size_t lane = 0; lane < 4; ++lane) { SCOPED_TRACE(lane); - EXPECT_EQ(PairBits(scalar[lane].high), PairBits(expected[lane].high)); - EXPECT_EQ(PairBits(scalar[lane].low), PairBits(expected[lane].low)); - EXPECT_EQ(PairBits(selected[lane].high), PairBits(expected[lane].high)); - EXPECT_EQ(PairBits(selected[lane].low), PairBits(expected[lane].low)); + EXPECT_EQ(PairBits(scalar.values[lane].high), PairBits(expected.values[lane].high)); + EXPECT_EQ(PairBits(scalar.values[lane].low), PairBits(expected.values[lane].low)); + EXPECT_EQ(PairBits(selected.values[lane].high), PairBits(expected.values[lane].high)); + EXPECT_EQ(PairBits(selected.values[lane].low), PairBits(expected.values[lane].low)); } return evaluator != fesa::hex8_internal::EvaluateMaterialPacketScalar; } @@ -147,51 +152,55 @@ TEST(Hex8HighLow, FourColumnMaterialPacketMatchesFrozenScalarBits) { {0, std::numeric_limits::signaling_NaN()}}; for (std::size_t first = 0; first < boundary.size(); ++first) for (std::size_t second = 0; second < boundary.size(); ++second) { - MaterialPacketInput input; + MaterialPacketInput input{}; input.weight = {0.125}; for (std::size_t lane = 0; lane < 4; ++lane) { - input.deviatoric_left[lane] = boundary[first]; - input.deviatoric_right[lane] = boundary[second]; - input.shear_left[lane] = {2, std::ldexp(1.0, -54)}; - input.shear_right[lane] = {-0.5, std::ldexp(1.0, -56)}; - input.tensor_left[lane] = {0.75}; - input.gradient_dot[lane] = {-0.25}; - input.accumulated[lane] = lane == 0 ? HighLow{} : HighLow{0.125}; + input.deviatoric_left.values[lane] = RawPacketPair(boundary[first]); + input.deviatoric_right.values[lane] = RawPacketPair(boundary[second]); + input.shear_left.values[lane] = {2, std::ldexp(1.0, -54)}; + input.shear_right.values[lane] = {-0.5, std::ldexp(1.0, -56)}; + input.tensor_left.values[lane] = {0.75}; + input.gradient_dot.values[lane] = {-0.25}; + input.accumulated.values[lane] = lane == 0 ? fesa::hex8_internal::PacketPair{0, 0} + : fesa::hex8_internal::PacketPair{0.125, 0}; } avx2_evaluator_checks += ExpectFrozenMaterialPacket(input) ? 1U : 0U; ++packet_checks; } for (int exponent = -1074; exponent <= 1023; ++exponent) { - MaterialPacketInput input; + MaterialPacketInput input{}; input.weight = {1}; for (std::size_t lane = 0; lane < 4; ++lane) { - input.deviatoric_left[lane] = { + input.deviatoric_left.values[lane] = { std::ldexp(lane % 2U == 0 ? 1.0 : -1.0, exponent)}; - input.deviatoric_right[lane] = {lane < 2 ? 0.5 : 2.0}; + input.deviatoric_right.values[lane] = {lane < 2 ? 0.5 : 2.0}; } avx2_evaluator_checks += ExpectFrozenMaterialPacket(input) ? 1U : 0U; ++packet_checks; } std::mt19937_64 random(0x4c0116U); for (std::size_t sample = 0; sample < 2048; ++sample) { - MaterialPacketInput input; + MaterialPacketInput input{}; const auto pair = [&]() { const double high = std::ldexp( 1.0 + static_cast(random() % 256U) / 256.0, static_cast(random() % 201U) - 100); const double signed_high = (random() & 1U) != 0 ? high : -high; const double low = sample % 3U == 0 ? 0.0 : std::ldexp(high, -55); - return HighLow{signed_high, (random() & 1U) != 0 ? low : -low}; + return fesa::hex8_internal::PacketPair{ + signed_high, (random() & 1U) != 0 ? low : -low}; }; - input.weight = sample % 2U == 0 ? HighLow{0.125} : pair(); + input.weight = sample % 2U == 0 + ? fesa::hex8_internal::PacketPair{0.125, 0} : pair(); for (std::size_t lane = 0; lane < 4; ++lane) { - input.deviatoric_left[lane] = pair(); - input.deviatoric_right[lane] = pair(); - input.shear_left[lane] = pair(); - input.shear_right[lane] = pair(); - input.tensor_left[lane] = pair(); - input.gradient_dot[lane] = pair(); - input.accumulated[lane] = sample % 4U == 0 ? HighLow{} : pair(); + input.deviatoric_left.values[lane] = pair(); + input.deviatoric_right.values[lane] = pair(); + input.shear_left.values[lane] = pair(); + input.shear_right.values[lane] = pair(); + input.tensor_left.values[lane] = pair(); + input.gradient_dot.values[lane] = pair(); + input.accumulated.values[lane] = sample % 4U == 0 + ? fesa::hex8_internal::PacketPair{0, 0} : pair(); } avx2_evaluator_checks += ExpectFrozenMaterialPacket(input) ? 1U : 0U; ++packet_checks; @@ -202,15 +211,15 @@ TEST(Hex8HighLow, FourColumnMaterialPacketMatchesFrozenScalarBits) { } TEST(Hex8HighLow, PacketSpecialLaneDoesNotChangeOtherColumns) { - fesa::hex8_internal::MaterialPacketInput input; + fesa::hex8_internal::MaterialPacketInput input{}; input.weight = {0.125}; for (std::size_t lane = 0; lane < 4; ++lane) { - input.deviatoric_left[lane] = {1, std::ldexp(1.0, -54)}; - input.deviatoric_right[lane] = {2, std::ldexp(1.0, -55)}; - input.shear_left[lane] = {0.5}; - input.shear_right[lane] = {-0.25}; - input.tensor_left[lane] = {0.75}; - input.gradient_dot[lane] = {0.125}; + input.deviatoric_left.values[lane] = {1, std::ldexp(1.0, -54)}; + input.deviatoric_right.values[lane] = {2, std::ldexp(1.0, -55)}; + input.shear_left.values[lane] = {0.5}; + input.shear_right.values[lane] = {-0.25}; + input.tensor_left.values[lane] = {0.75}; + input.gradient_dot.values[lane] = {0.125}; } const auto reference = FrozenMaterialPacket(input); const std::vector special{{}, {-0.0, 0}, {1}, @@ -220,14 +229,14 @@ TEST(Hex8HighLow, PacketSpecialLaneDoesNotChangeOtherColumns) { for (std::size_t lane = 0; lane < 4; ++lane) for (const auto value : special) { auto changed = input; - changed.deviatoric_left[lane] = value; + changed.deviatoric_left.values[lane] = RawPacketPair(value); ExpectFrozenMaterialPacket(changed); fesa::hex8_internal::PairPacket actual; fesa::hex8_internal::SelectMaterialPacketEvaluator()(changed, &actual); for (std::size_t other = 0; other < 4; ++other) if (other != lane) { - EXPECT_EQ(PairBits(actual[other].high), PairBits(reference[other].high)); - EXPECT_EQ(PairBits(actual[other].low), PairBits(reference[other].low)); + EXPECT_EQ(PairBits(actual.values[other].high), PairBits(reference.values[other].high)); + EXPECT_EQ(PairBits(actual.values[other].low), PairBits(reference.values[other].low)); } } } @@ -237,7 +246,7 @@ TEST(Hex8HighLow, UniformDyadicPacketsExerciseEachLowZeroShape) { for (int shape = 0; shape < 4; ++shape) for (int low_weight = 0; low_weight < 2; ++low_weight) for (int zero_accumulator = 0; zero_accumulator < 2; ++zero_accumulator) { - fesa::hex8_internal::MaterialPacketInput input; + fesa::hex8_internal::MaterialPacketInput input{}; input.weight = {0.5, low_weight != 0 ? std::ldexp(1.0, -57) : 0.0}; for (std::size_t lane = 0; lane < 4; ++lane) { const bool left_low = (shape & 1) != 0; @@ -247,14 +256,15 @@ TEST(Hex8HighLow, UniformDyadicPacketsExerciseEachLowZeroShape) { // Thus each Multiply low-shape is uniform and its component products // bypass zero/unsupported-tail fallback. No mixed-zero Add fallback // is possible; the initial zero accumulator is uniform across lanes. - input.deviatoric_left[lane] = {1, left_low ? std::ldexp(1.0, -54) : 0}; - input.deviatoric_right[lane] = {2, right_low ? std::ldexp(1.0, -54) : 0}; - input.shear_left[lane] = {0.5, left_low ? std::ldexp(1.0, -56) : 0}; - input.shear_right[lane] = {0.25, right_low ? std::ldexp(1.0, -57) : 0}; - input.tensor_left[lane] = {0.75, left_low ? std::ldexp(1.0, -56) : 0}; - input.gradient_dot[lane] = {0.125, right_low ? std::ldexp(1.0, -58) : 0}; - input.accumulated[lane] = zero_accumulator != 0 ? HighLow{} - : HighLow{0.0625}; + input.deviatoric_left.values[lane] = {1, left_low ? std::ldexp(1.0, -54) : 0}; + input.deviatoric_right.values[lane] = {2, right_low ? std::ldexp(1.0, -54) : 0}; + input.shear_left.values[lane] = {0.5, left_low ? std::ldexp(1.0, -56) : 0}; + input.shear_right.values[lane] = {0.25, right_low ? std::ldexp(1.0, -57) : 0}; + input.tensor_left.values[lane] = {0.75, left_low ? std::ldexp(1.0, -56) : 0}; + input.gradient_dot.values[lane] = {0.125, right_low ? std::ldexp(1.0, -58) : 0}; + input.accumulated.values[lane] = zero_accumulator != 0 + ? fesa::hex8_internal::PacketPair{0, 0} + : fesa::hex8_internal::PacketPair{0.0625, 0}; } avx2_checks += ExpectFrozenMaterialPacket(input) ? 1U : 0U; ++checks; @@ -270,17 +280,17 @@ TEST(Hex8HighLow, OverflowLaneIsNeutralizedBeforeResidualFma) { std::fenv_t original; ASSERT_EQ(std::fegetenv(&original), 0); for (std::size_t bad = 0; bad < 4; ++bad) { - fesa::hex8_internal::MaterialPacketInput input; + fesa::hex8_internal::MaterialPacketInput input{}; input.weight = {0.5}; for (std::size_t lane = 0; lane < 4; ++lane) { - input.deviatoric_left[lane] = { + input.deviatoric_left.values[lane] = { lane == bad ? std::numeric_limits::max() : 1}; - input.deviatoric_right[lane] = {2}; - input.shear_left[lane] = {0.5}; - input.shear_right[lane] = {0.25}; - input.tensor_left[lane] = {0.75}; - input.gradient_dot[lane] = {0.125}; - input.accumulated[lane] = {0.0625}; + input.deviatoric_right.values[lane] = {2}; + input.shear_left.values[lane] = {0.5}; + input.shear_right.values[lane] = {0.25}; + input.tensor_left.values[lane] = {0.75}; + input.gradient_dot.values[lane] = {0.125}; + input.accumulated.values[lane] = {0.0625}; } // All operands are finite/nonzero and max*2 has LBE sum 972, so Product // reaches its vector multiply; only the post-product finite mask rejects it. @@ -292,8 +302,8 @@ TEST(Hex8HighLow, OverflowLaneIsNeutralizedBeforeResidualFma) { EXPECT_EQ(std::fesetenv(&original), 0); EXPECT_FALSE(extra_invalid); for (std::size_t lane = 0; lane < 4; ++lane) { - EXPECT_EQ(PairBits(output[lane].high), PairBits(expected[lane].high)); - EXPECT_EQ(PairBits(output[lane].low), PairBits(expected[lane].low)); + EXPECT_EQ(PairBits(output.values[lane].high), PairBits(expected.values[lane].high)); + EXPECT_EQ(PairBits(output.values[lane].low), PairBits(expected.values[lane].low)); } } RecordProperty("avx2_product_overflow_positions", 4); @@ -320,16 +330,73 @@ TEST(Hex8HighLow, PacketBackendSelectionRecoversAfterEnvironmentRestoration) { EXPECT_EQ(fesa::hex8_internal::SelectMaterialPacketEvaluator(), initial); } +TEST(Hex8HighLow, PacketScalarBridgesPreserveOriginalRawBits) { + using fesa::hex8_internal::PacketPair; + using fesa::hex8_internal::PairPacket; + static_assert(std::is_trivial_v && + std::is_standard_layout_v); + static_assert(std::is_trivial_v && + std::is_standard_layout_v); + static_assert(std::is_trivial_v); + const double tiny = std::numeric_limits::denorm_min(); + const std::vector pairs{{0, 0}, {-0.0, 0}, {0, -0.0}, + {1, -0.0}, {1, std::ldexp(1.0, -54)}, {tiny, 0}, {1, tiny}, + {std::numeric_limits::max(), 0}, {1, 1}, + {std::numeric_limits::infinity(), 0}, + {std::numeric_limits::signaling_NaN(), 0}}; + for (const auto left : pairs) + for (const auto right : pairs) { + PairPacket a{{left, left, left, left}}; + PairPacket b{{right, right, right, right}}; + PairPacket result; + const fesa::high_low_original::HighLow x{left.high, left.low}; + const fesa::high_low_original::HighLow y{right.high, right.low}; + const auto expect = [&](fesa::high_low_original::HighLow expected) { + for (std::size_t lane = 0; lane < 4; ++lane) { + EXPECT_EQ(PairBits(result.values[lane].high), PairBits(expected.high)); + EXPECT_EQ(PairBits(result.values[lane].low), PairBits(expected.low)); + } + }; + fesa::hex8_internal::ScalarPacketAdd(a, b, false, &result); + expect(fesa::high_low_original::Add(x, y)); + fesa::hex8_internal::ScalarPacketAdd(a, b, true, &result); + expect(fesa::high_low_original::Add(x, y)); + fesa::hex8_internal::ScalarPacketMultiply(a, b, &result); + expect(fesa::high_low_original::Multiply(x, y)); + const double high_left[4]{left.high, left.high, left.high, left.high}; + const double high_right[4]{right.high, right.high, right.high, right.high}; + fesa::hex8_internal::ScalarPacketTwoProduct(high_left, high_right, &result); + expect(fesa::high_low_original::TwoProduct(left.high, right.high)); + } +} + +TEST(Hex8HighLow, PacketCapabilityRejectsMissingIntegerIsaAndOsState) { + constexpr std::uint32_t leaf1 = + (1U << 12) | (1U << 26) | (1U << 27) | (1U << 28); + constexpr std::uint32_t leaf7 = (1U << 3) | (1U << 5) | (1U << 8); + EXPECT_TRUE(fesa::hex8_internal::SupportsMaterialPacketInstructions( + leaf1, 6, leaf7)); + for (const auto missing : {1U << 3, 1U << 5, 1U << 8}) + EXPECT_FALSE(fesa::hex8_internal::SupportsMaterialPacketInstructions( + leaf1, 6, leaf7 & ~missing)); + for (const auto missing : {1U << 12, 1U << 26, 1U << 27, 1U << 28}) + EXPECT_FALSE(fesa::hex8_internal::SupportsMaterialPacketInstructions( + leaf1 & ~missing, 6, leaf7)); + for (std::uint64_t state : {0ULL, 2ULL, 4ULL}) + EXPECT_FALSE(fesa::hex8_internal::SupportsMaterialPacketInstructions( + leaf1, state, leaf7)); +} + TEST(Hex8HighLow, PacketDoesNotClearPriorMaskedStickyFlags) { std::fenv_t original; ASSERT_EQ(std::fegetenv(&original), 0); const auto evaluator = fesa::hex8_internal::SelectMaterialPacketEvaluator(); ASSERT_EQ(std::feraiseexcept(FE_INEXACT), 0); - fesa::hex8_internal::MaterialPacketInput input; + fesa::hex8_internal::MaterialPacketInput input{}; input.weight = {1}; for (std::size_t lane = 0; lane < 4; ++lane) { - input.deviatoric_left[lane] = {1, std::ldexp(1.0, -54)}; - input.deviatoric_right[lane] = {2, std::ldexp(1.0, -54)}; + input.deviatoric_left.values[lane] = {1, std::ldexp(1.0, -54)}; + input.deviatoric_right.values[lane] = {2, std::ldexp(1.0, -54)}; } fesa::hex8_internal::PairPacket output; evaluator(input, &output);