diff --git a/docs/calculation-performance/hex8-performance.md b/docs/calculation-performance/hex8-performance.md index 2cd86c1..4d86043 100644 --- a/docs/calculation-performance/hex8-performance.md +++ b/docs/calculation-performance/hex8-performance.md @@ -656,3 +656,61 @@ Dot3-term과 geometric9-term의 독립 b4 lane을 기존 zero-start·component/j 각 기존 diagonal accumulator Add에 쓰는 작은 대안도 명시했다. High/low 정밀도는 유지하며 일반 Float64 MKL 경로를 바꾸지 않는다. 이는 아직 구현·검증 전인 제안이고 현재 source를 변경하지 않았다. Whole-model 추가 반복이나 profile 추가 수집은 없다. + +## 기하 수축 packet 구현과 Release 상태 차이 수정 + +위 제안을 `9677265`에서 구현했다. 독립 b4 lane은 기존 `dots[b]` 및 material의 +`col=3*b+k`와 연결하며, gradient-stress의 원래 준비 제품과 Dot의 component 순서, +geometric의 j,l 순서 및 zero-start를 보존한다. 동일 `geo*weight`는 공유하지만 세 +diagonal accumulator의 Add는 각각 수행한다. Unsupported 환경과 linear 응답은 +원래 scalar loop를 사용한다. Private POD·baseline bridge·AVX2 경계를 유지하고 +일반 Float64 MKL, CMake, 공통 math 및 public ABI를 바꾸지 않았다. + +이 후보의 Debug 관련 81개 검사는 통과했으나 Release는 80/81 통과였다. +Ordered overflow fixture에서 raw bit는 같고 selected FE21/baseline5 및 MXCSR41/40, +prior 상태에서는43/42로 INVALID 차이가 있었다. 당시 wrapper가 stderr를 먼저 +예외 처리하여 CTest native는 unknown, wrapper native는1이었다. 이 실패를 현재 +통과로 소급하지 않는다. [실패 XML](../../.harness/performance/t3/packet-geometric-release-gtest/fesa_unit_tests_53.xml), +[실패 binary 보존](../../.harness/performance/t3/packet-geometric-9677265/failure-release/manifest.json)을 유지한다. +후속 Release codegen·snapshot·전체 모델·profile은 당시 수행하지 않았다. + +최소 수정 `0c29a9f`는 private AVX2 cpp와 기존 Hex8 검사 두 파일에 한정한다. +Vector TwoSum의 첫 sum이 nonfinite이면 error 계산 전에 false를 반환한다. 각 Add의 +head/tail/middle/joined/final 위치와 Product의 마지막 TwoSum 실패는 원래 입력과 +mode로 기존 noinline scalar bridge를 호출한다. Partial 결과를 재사용하지 않으며 +조건부 Debug 분기, flags clear, 기대값 완화, 산술식 및 진단 gate 변경은 없다. +[고정 manifest](../../.harness/performance/t3/packet-geometric-fix-0c29a9f/manifest.json)와 +[독립 source 검토](../../.harness/performance/t3/source-numerical-review-packet-geometric-fix.md)는 +10개 hash의 전후 일치와 source-only PASS를 기록한다. + +새 검사는 normalized near-max joined overflow, N=2^-1022와 d=2^-1074에 의한 +subnormal error, 다른 lane의 성공 prefix 이후 overflow, 각 lane 위치 및 prior 상태를 +포함한다. 원래 세 weight 제품과 서로 다른 세 accumulator의 raw bit·FE_ALL_EXCEPT· +MXCSR status를 기준으로 삼았다. 첫 narrow 실행은 Debug·Release 각각2/2, +actual AVX2의 기존 실패6개와 새 prefix34개가 통과했다. 이후 실제 discovery는 +82개였고 Debug 82/82(10.43초), Release 82/82(4.46초), 모든 단계 native0이었다. +[Debug 명령/native](../../.harness/performance/t3/packet-geometric-fix-full-debug-native.json), +[Release 명령/native](../../.harness/performance/t3/packet-geometric-fix-full-release-native.json), +[Release XML](../../.harness/performance/t3/packet-geometric-fix-full-release-gtest/)을 보존한다. +Configure/build는 corrected MSVC dependency prefix, 병렬4 이하로 순차 수행했다. + +두 구성의 전체 객체 symbol/disassembly는 AVX2 외부 실행 정의2개, 공통 math 정의0개, +baseline AVX·YMM·ANDN0개를 확인했다. Kernel의 fused 명령은 명시 residual VFMSUB +1개뿐이고 Product의 충분조건 true 분기4개는 각각 두 BSF를 건너뛰며 false는 원래 +두 BSF를 유지한다. [코드 생성 assertion](../../.harness/performance/t3/packet-geometric-fix-codegen-assertion.json)을 보존한다. +실제 Release dispatch 객체의 scalar Add는 first-sum0x16e 이후 finite 분기0x18a에서 +error 계산0x18f를 건너뛴다. TwoSum도0x25에서0xb1로 같은 순서의 return을 한다. +Linked unit exe에서 해당 Add byte signature가 한 곳에 존재함을 확인했다. +[Hoist 근거](../../.harness/performance/t3/packet-geometric-fix-release-scalar-hoist.json)는 +PDB를 통한 개별 callsite나 COMDAT winner 식별까지 주장하지 않는다. 실제 bridge의 +raw/status 동등성은 위 Debug·Release oracle로 별도 확인했다. + +기존 actual native0 원본 oracle의 정합 header/LIB/입력 연결을 검증한 뒤 새 candidate만 +compile/run했다. 두 native와 엄격한 비교 wrapper는0이며9,804행과 전체 byte가 동일하다. +양쪽 SHA는 `DED8D2C380D15E297A73BDD018B1B840C1ADC0D8372934D9FC2EF2B628A36F63`이다. +[Snapshot 명령/native](../../.harness/performance/t3/packet-geometric-fix-snapshot-native.json), +[엄격한 비교](../../.harness/performance/t3/packet-geometric-fix-snapshot-comparison.json)를 보존한다. +[종료 확인](../../.harness/performance/t3/packet-geometric-fix-full-cleanup.json)은 source/hash +불변, tracked clean 및 관련 actor0을 기록한다. CPU를 반환했고 이번 슬롯에서 전체 모델과 +profile은 실행하지 않았다. 이 후보의 전체 경과 시간 및 최종101-frame 비교는 미실행이며, +기존135.9302297초 개발 결과와 목표100초 미달, 최종 전체·공식 검증 미완료는 유지한다. diff --git a/src/fesa/elements/hex8_high_low.cpp b/src/fesa/elements/hex8_high_low.cpp index ccb873f..ed1f9dc 100644 --- a/src/fesa/elements/hex8_high_low.cpp +++ b/src/fesa/elements/hex8_high_low.cpp @@ -246,11 +246,38 @@ Result Hex8Element::EvaluateHighLow( gradient_stress[a][j][l] = gradients[a][j] * s[l][j]; for (std::size_t a = 0; a < 8; ++a) { std::array dots{}; + std::array weighted_geo{}; + if (use_packets) + for (std::size_t first_b = 0; first_b < 8; first_b += 4) { + hex8_internal::GeometryPacketInput input{}; + const auto raw_pair = [](HighLow pair) { + return hex8_internal::PacketPair{pair.high, pair.low}; + }; + input.weight = raw_pair(HighLow{point.weight}); + for (std::size_t j = 0; j < 3; ++j) { + input.gradient_a[j] = raw_pair(gradients[a][j]); + for (std::size_t l = 0; l < 3; ++l) + input.gradient_stress[j][l] = raw_pair(gradient_stress[a][j][l]); + for (std::size_t lane = 0; lane < 4; ++lane) + input.gradient_b[j].values[lane] = + raw_pair(gradients[first_b + lane][j]); + } + hex8_internal::GeometryPacketOutput output; + // Same selector contract as material; unsupported/linear stays below. + hex8_internal::EvaluateGeometryPacketAvx2(input, &output); + for (std::size_t lane = 0; lane < 4; ++lane) { + dots[first_b + lane] = {output.gradient_dot.values[lane].high, + output.gradient_dot.values[lane].low}; + weighted_geo[first_b + lane] = + {output.weighted_geometric.values[lane].high, + output.weighted_geometric.values[lane].low}; + } + } for (std::size_t b = 0; b < 8; ++b) { - const auto dot = Dot(gradients[a], gradients[b]); + const auto dot = use_packets ? dots[b] : Dot(gradients[a], gradients[b]); dots[b] = dot; HighLow geo; - if (nonlinear) + if (nonlinear && !use_packets) for (std::size_t j = 0; j < 3; ++j) for (std::size_t l = 0; l < 3; ++l) geo += gradient_stress[a][j][l] * gradients[b][l]; @@ -265,7 +292,8 @@ Result Hex8Element::EvaluateHighLow( material[row][column] += mat * HighLow{point.weight}; } if (i == k) - geometric[row][column] += geo * HighLow{point.weight}; + geometric[row][column] += + use_packets ? weighted_geo[b] : geo * HighLow{point.weight}; } } if (use_packets) diff --git a/src/fesa/elements/hex8_pair_simd.h b/src/fesa/elements/hex8_pair_simd.h index 6a0fbbd..4dffa07 100644 --- a/src/fesa/elements/hex8_pair_simd.h +++ b/src/fesa/elements/hex8_pair_simd.h @@ -23,6 +23,30 @@ struct MaterialPacketInput { PairPacket accumulated; PacketPair weight; }; +/// @brief One owned trial point/a and four independent b gradients, raw only. +struct GeometryPacketInput { + PacketPair gradient_a[3]; + PacketPair gradient_stress[3][3]; + PairPacket gradient_b[3]; + PacketPair weight; +}; +/// @brief Independent raw contractions; accumulator entries remain caller-owned. +struct GeometryPacketOutput { + PairPacket gradient_dot; + PairPacket weighted_geometric; +}; +using GeometryPacketEvaluator = void (*)(const GeometryPacketInput&, + GeometryPacketOutput*); +/// @brief Baseline checked zero-start Dot3/ordered j,l Geo9 and weight products. +/// @note Lane order is 0..3. Invalid pairs remain invalid for original gates. +void EvaluateGeometryPacketScalar(const GeometryPacketInput& input, + GeometryPacketOutput* output); +/// @brief Evaluates independent b lanes without horizontal reduction. +/// @pre The material selector has admitted CPU/OS/current masked FP state. +/// @note Keeps original zero-start, product association and j,l Add order. +void EvaluateGeometryPacketAvx2(const GeometryPacketInput& input, + GeometryPacketOutput* output); + using MaterialPacketEvaluator = void (*)(const MaterialPacketInput&, PairPacket*); /// @brief Evaluates each lane's original checked material graph in order 0..3. diff --git a/src/fesa/elements/hex8_pair_simd_avx2.cpp b/src/fesa/elements/hex8_pair_simd_avx2.cpp index 46e105d..570b2a5 100644 --- a/src/fesa/elements/hex8_pair_simd_avx2.cpp +++ b/src/fesa/elements/hex8_pair_simd_avx2.cpp @@ -111,6 +111,10 @@ VectorPair Load(const PairPacket& value) { _mm256_setr_pd(value.values[0].low, value.values[1].low, value.values[2].low, value.values[3].low)}; } +/// @brief Replicates one raw first product without changing either component. +VectorPair Splat(PacketPair value) { + return {_mm256_set1_pd(value.high), _mm256_set1_pd(value.low)}; +} /// @brief Stores raw high/low fields without rounding or lane reduction. PairPacket Store(VectorPair value) { double high[4]; @@ -120,9 +124,14 @@ PairPacket Store(VectorPair value) { 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) { +/// @brief Retains finite-sum error association; false requests original graph. +/// @note False leaves output unreadable and preserves the necessary sum flags. +/// Prefix flags still require optimized scalar-oracle/codegen validation. +bool TwoSum(__m256d left, __m256d right, VectorPair* output) { const auto sum = _mm256_add_pd(left, right); + // Do not speculate Inf-Inf error arithmetic beyond a failed first sum. + // The caller retains original inputs for its existing baseline bridge. + if (_mm256_movemask_pd(Finite(sum)) != 15) return false; const auto bridge = _mm256_sub_pd(sum, left); const auto restored = _mm256_sub_pd(sum, bridge); const auto error = _mm256_add_pd(_mm256_sub_pd(left, restored), @@ -130,8 +139,9 @@ VectorPair TwoSum(__m256d left, __m256d right) { const auto valid = _mm256_and_pd(_mm256_and_pd(Finite(sum), Finite(bridge)), _mm256_and_pd(Finite(restored), Finite(error))); const auto exact_zero = _mm256_and_pd(Zero(sum), Zero(error)); - return InvalidOutside({Select(exact_zero, _mm256_setzero_pd(), sum), - Select(exact_zero, _mm256_setzero_pd(), error)}, valid); + *output = InvalidOutside({Select(exact_zero, _mm256_setzero_pd(), sum), + Select(exact_zero, _mm256_setzero_pd(), error)}, valid); + return true; } /// @brief Routes exceptional packets to baseline checked Add in lane order. VectorPair ScalarAdd(VectorPair left, VectorPair right, bool pair_shortcuts) { @@ -163,10 +173,15 @@ VectorPair Add(VectorPair left, VectorPair right, bool pair_shortcuts) { if (left_zero != 0 || right_zero != 0) return ScalarAdd(left, right, true); } - const auto head = TwoSum(left.high, right.high); - const auto tail = TwoSum(left.low, right.low); - const auto middle = TwoSum(head.low, tail.high); - const auto joined = TwoSum(head.high, middle.high); + VectorPair head; + VectorPair tail; + VectorPair middle; + VectorPair joined; + if (!TwoSum(left.high, right.high, &head) || + !TwoSum(left.low, right.low, &tail) || + !TwoSum(head.low, tail.high, &middle) || + !TwoSum(head.high, middle.high, &joined)) + return ScalarAdd(left, right, pair_shortcuts); const auto remainder = _mm256_add_pd( _mm256_add_pd(joined.low, middle.low), tail.low); const auto valid = _mm256_and_pd(Finite(remainder), _mm256_and_pd( @@ -175,9 +190,11 @@ VectorPair Add(VectorPair left, VectorPair right, bool pair_shortcuts) { _mm256_and_pd(Finite(middle.high), Finite(middle.low)), _mm256_and_pd(Finite(joined.high), Finite(joined.low)))))); // Original Add returns before its last TwoSum when a prior finite gate fails. - return InvalidOutside(TwoSum( - Select(valid, joined.high, _mm256_setzero_pd()), - Select(valid, remainder, _mm256_setzero_pd())), valid); + VectorPair result; + if (!TwoSum(Select(valid, joined.high, _mm256_setzero_pd()), + Select(valid, remainder, _mm256_setzero_pd()), &result)) + return ScalarAdd(left, right, pair_shortcuts); + return InvalidOutside(result, valid); } /// @brief Checks finite/zero/tail premises before products; neutralizes rejected /// lanes before the explicit residual FMA and final TwoSum. @@ -207,9 +224,14 @@ VectorPair Product(__m256d left, __m256d right) { const auto error = _mm256_fmsub_pd(Select(product_valid, left, one), Select(product_valid, right, one), Select(product_valid, product, one)); const auto valid = _mm256_and_pd(product_valid, Finite(error)); - return InvalidOutside(TwoSum( - Select(valid, product, _mm256_setzero_pd()), - Select(valid, error, _mm256_setzero_pd())), valid); + VectorPair result; + if (!TwoSum(Select(valid, product, _mm256_setzero_pd()), + Select(valid, error, _mm256_setzero_pd()), &result)) { + PairPacket scalar; + ScalarPacketTwoProduct(a, b, &scalar); + return Load(scalar); + } + return InvalidOutside(result, valid); } /// @brief Uses uniform low-zero branches with original four-product association; /// mixed or invalid packets preserve the baseline checked path. @@ -251,6 +273,20 @@ void EvaluateMaterialPacketAvx2(const MaterialPacketInput& input, *output = Store(Add(Load(input.accumulated), Multiply(material, Load(weight)), true)); } +void EvaluateGeometryPacketAvx2(const GeometryPacketInput& input, + GeometryPacketOutput* output) { + VectorPair dot{_mm256_setzero_pd(), _mm256_setzero_pd()}; + for (std::size_t i = 0; i < 3; ++i) + dot = Add(dot, Multiply(Splat(input.gradient_a[i]), + Load(input.gradient_b[i])), true); + VectorPair geo{_mm256_setzero_pd(), _mm256_setzero_pd()}; + for (std::size_t j = 0; j < 3; ++j) + for (std::size_t l = 0; l < 3; ++l) + geo = Add(geo, Multiply(Splat(input.gradient_stress[j][l]), + Load(input.gradient_b[l])), true); + output->gradient_dot = Store(dot); + output->weighted_geometric = Store(Multiply(geo, Splat(input.weight))); +} } // namespace fesa::hex8_internal #else namespace fesa::hex8_internal { @@ -258,5 +294,9 @@ void EvaluateMaterialPacketAvx2(const MaterialPacketInput& input, PairPacket* output) { EvaluateMaterialPacketScalar(input, output); } +void EvaluateGeometryPacketAvx2(const GeometryPacketInput& input, + GeometryPacketOutput* output) { + EvaluateGeometryPacketScalar(input, output); +} } // namespace fesa::hex8_internal #endif diff --git a/src/fesa/elements/hex8_pair_simd_dispatch.cpp b/src/fesa/elements/hex8_pair_simd_dispatch.cpp index a4838d7..97e0118 100644 --- a/src/fesa/elements/hex8_pair_simd_dispatch.cpp +++ b/src/fesa/elements/hex8_pair_simd_dispatch.cpp @@ -42,6 +42,27 @@ void ScalarPacketTwoProduct(const double left[4], const double right[4], output->values[lane] = {value.high, value.low}; } } +FESA_HEX8_BASELINE_NOINLINE +void EvaluateGeometryPacketScalar(const GeometryPacketInput& input, + GeometryPacketOutput* output) { + const auto value = [](PacketPair pair) { + return math_internal::HighLow{pair.high, pair.low}; + }; + for (std::size_t lane = 0; lane < 4; ++lane) { + math_internal::HighLow dot; + for (std::size_t i = 0; i < 3; ++i) + dot = AddPair(dot, MultiplyPair(value(input.gradient_a[i]), + value(input.gradient_b[i].values[lane]))); + math_internal::HighLow geo; + for (std::size_t j = 0; j < 3; ++j) + for (std::size_t l = 0; l < 3; ++l) + geo = AddPair(geo, MultiplyPair(value(input.gradient_stress[j][l]), + value(input.gradient_b[l].values[lane]))); + const auto weighted = MultiplyPair(geo, value(input.weight)); + output->gradient_dot.values[lane] = {dot.high, dot.low}; + output->weighted_geometric.values[lane] = {weighted.high, weighted.low}; + } +} #undef FESA_HEX8_BASELINE_NOINLINE bool SupportsMaterialPacketInstructions(std::uint32_t leaf1_ecx, diff --git a/tests/unit/elements/hex8_high_low_test.cpp b/tests/unit/elements/hex8_high_low_test.cpp index ad01adf..d1df9e7 100644 --- a/tests/unit/elements/hex8_high_low_test.cpp +++ b/tests/unit/elements/hex8_high_low_test.cpp @@ -449,6 +449,296 @@ TEST(Hex8HighLow, PacketTailScreenRetainsExactBoundaryAndFallback) { } #endif +#if defined(_M_X64) || defined(__SSE2__) +fesa::hex8_internal::GeometryPacketOutput FrozenGeometryPacket( + const fesa::hex8_internal::GeometryPacketInput& input) { + const auto pair = [](fesa::hex8_internal::PacketPair value) { + return fesa::high_low_original::HighLow{value.high, value.low}; + }; + fesa::hex8_internal::GeometryPacketOutput output{}; + for (std::size_t lane = 0; lane < 4; ++lane) { + fesa::high_low_original::HighLow dot; + for (std::size_t i = 0; i < 3; ++i) + dot = fesa::high_low_original::Add(dot, + fesa::high_low_original::Multiply(pair(input.gradient_a[i]), + pair(input.gradient_b[i].values[lane]))); + fesa::high_low_original::HighLow geo; + for (std::size_t j = 0; j < 3; ++j) + for (std::size_t l = 0; l < 3; ++l) + geo = fesa::high_low_original::Add(geo, + fesa::high_low_original::Multiply(pair(input.gradient_stress[j][l]), + pair(input.gradient_b[l].values[lane]))); + const auto weighted = fesa::high_low_original::Multiply(geo, pair(input.weight)); + output.gradient_dot.values[lane] = {dot.high, dot.low}; + output.weighted_geometric.values[lane] = {weighted.high, weighted.low}; + } + return output; +} + +fesa::hex8_internal::GeometryPacketInput GeometryDyadicFixture() { + fesa::hex8_internal::GeometryPacketInput input{}; + input.weight = {0.125, 0}; + for (std::size_t j = 0; j < 3; ++j) { + input.gradient_a[j] = {static_cast(j + 1) * 0.25, 0}; + for (std::size_t l = 0; l < 3; ++l) { + // Deliberately asymmetric first-products and distinct b/component mapping. + input.gradient_stress[j][l] = + {static_cast(1 + 3 * j + l) * 0.125, 0}; + for (std::size_t lane = 0; lane < 4; ++lane) + input.gradient_b[l].values[lane] = + {static_cast(1 + 4 * l + lane) * 0.0625, 0}; + } + } + return input; +} + +// Original observation repeats geo*weight for three independent accumulators. +// Selected observation reuses the weighted pair but retains each checked Add. +void ExpectGeometryPacketBitsAndFlags( + const fesa::hex8_internal::GeometryPacketInput& input, bool prior_status) { + SCOPED_TRACE(prior_status); + std::fenv_t original; + ASSERT_EQ(std::fegetenv(&original), 0); + const auto original_csr = _mm_getcsr(); + struct Observation { + fesa::hex8_internal::GeometryPacketOutput output; + fesa::hex8_internal::PairPacket accumulators[3]; + int flags; + unsigned int status; + }; + const auto current_pair = [](fesa::hex8_internal::PacketPair value) { + return HighLow{value.high, value.low}; + }; + const auto observe = [&](bool original_order) { + EXPECT_EQ(std::fesetenv(&original), 0); + _mm_setcsr(original_csr & ~0x3fU); + EXPECT_EQ(std::feclearexcept(FE_ALL_EXCEPT), 0); + if (prior_status) { + EXPECT_EQ(std::feraiseexcept(FE_INEXACT), 0); + _mm_setcsr(_mm_getcsr() | 0x02U); + } + Observation value{}; + if (!original_order) + fesa::hex8_internal::EvaluateGeometryPacketAvx2(input, &value.output); + for (std::size_t lane = 0; lane < 4; ++lane) { + HighLow geo; + if (original_order) { + HighLow dot; + for (std::size_t i = 0; i < 3; ++i) + dot = fesa::hex8_internal::AddPair(dot, + fesa::hex8_internal::MultiplyPair(current_pair(input.gradient_a[i]), + current_pair(input.gradient_b[i].values[lane]))); + for (std::size_t j = 0; j < 3; ++j) + for (std::size_t l = 0; l < 3; ++l) + geo = fesa::hex8_internal::AddPair(geo, + fesa::hex8_internal::MultiplyPair( + current_pair(input.gradient_stress[j][l]), + current_pair(input.gradient_b[l].values[lane]))); + value.output.gradient_dot.values[lane] = {dot.high, dot.low}; + } + const HighLow accumulators[3]{{}, + {static_cast(lane + 1), std::ldexp(1.0, -54)}, {-1, -0.0}}; + for (std::size_t component = 0; component < 3; ++component) { + const auto weighted = original_order + ? fesa::hex8_internal::MultiplyPair(geo, current_pair(input.weight)) + : current_pair(value.output.weighted_geometric.values[lane]); + if (original_order) + value.output.weighted_geometric.values[lane] = + {weighted.high, weighted.low}; + const auto accumulated = fesa::hex8_internal::AddPair( + accumulators[component], weighted); + value.accumulators[component].values[lane] = + {accumulated.high, accumulated.low}; + } + } + value.flags = std::fetestexcept(FE_ALL_EXCEPT); + value.status = _mm_getcsr() & 0x3fU; + return value; + }; + const auto baseline = observe(true); + const auto selected = observe(false); + fesa::hex8_internal::GeometryPacketOutput scalar; + fesa::hex8_internal::EvaluateGeometryPacketScalar(input, &scalar); + const auto expected = FrozenGeometryPacket(input); + EXPECT_EQ(std::fesetenv(&original), 0); + _mm_setcsr(original_csr); + EXPECT_EQ(selected.flags, baseline.flags); + EXPECT_EQ(selected.status, baseline.status); + const auto same_bits = [](const fesa::hex8_internal::PairPacket& a, + const fesa::hex8_internal::PairPacket& b) { + for (std::size_t lane = 0; lane < 4; ++lane) { + EXPECT_EQ(PairBits(a.values[lane].high), PairBits(b.values[lane].high)); + EXPECT_EQ(PairBits(a.values[lane].low), PairBits(b.values[lane].low)); + } + }; + same_bits(selected.output.gradient_dot, expected.gradient_dot); + same_bits(selected.output.weighted_geometric, expected.weighted_geometric); + same_bits(scalar.gradient_dot, expected.gradient_dot); + same_bits(scalar.weighted_geometric, expected.weighted_geometric); + for (std::size_t component = 0; component < 3; ++component) + same_bits(selected.accumulators[component], baseline.accumulators[component]); +} + +TEST(Hex8HighLow, GeometryPacketPreservesAsymmetricLanesAndUniformShapes) { + if (fesa::hex8_internal::SelectMaterialPacketEvaluator() == + fesa::hex8_internal::EvaluateMaterialPacketScalar) + GTEST_SKIP() << "Geometry AVX2 masked environment is not selected."; + RecordProperty("private_geometry_packet_backend", "avx2"); + static_assert(std::is_trivial_v); + static_assert(std::is_trivial_v); + std::size_t checks = 0; + for (int shape = 0; shape < 4; ++shape) { + auto input = GeometryDyadicFixture(); + for (std::size_t j = 0; j < 3; ++j) { + if ((shape & 1) != 0) { + input.gradient_a[j].low = std::ldexp(input.gradient_a[j].high, -55); + for (std::size_t l = 0; l < 3; ++l) + input.gradient_stress[j][l].low = + std::ldexp(input.gradient_stress[j][l].high, -55); + } + if ((shape & 2) != 0) + for (std::size_t lane = 0; lane < 4; ++lane) + input.gradient_b[j].values[lane].low = + std::ldexp(input.gradient_b[j].values[lane].high, -55); + } + ExpectGeometryPacketBitsAndFlags(input, false); + ExpectGeometryPacketBitsAndFlags(input, true); + checks += 2; + } + RecordProperty("avx2_geometry_uniform_bit_flag_checks", static_cast(checks)); +} + +TEST(Hex8HighLow, GeometryPacketRetainsSpecialLanesAndExactTailBoundary) { + if (fesa::hex8_internal::SelectMaterialPacketEvaluator() == + fesa::hex8_internal::EvaluateMaterialPacketScalar) + GTEST_SKIP() << "Geometry AVX2 masked environment is not selected."; + using fesa::hex8_internal::PacketPair; + const double tiny = std::numeric_limits::denorm_min(); + const double maximum = std::numeric_limits::max(); + const std::vector special{{0, 0}, {-0.0, 0}, {1, -0.0}, {1, 1}, + {tiny, 0}, {maximum, maximum}, + {std::numeric_limits::signaling_NaN(), 0}, + {0, std::numeric_limits::signaling_NaN()}}; + std::size_t checks = 0; + for (const auto pair : special) + for (std::size_t lane = 0; lane < 4; ++lane) { + auto input = GeometryDyadicFixture(); + input.gradient_b[1].values[lane] = pair; + ExpectGeometryPacketBitsAndFlags(input, false); + ExpectGeometryPacketBitsAndFlags(input, true); + checks += 2; + } + const auto odd = [](int exponent) { + return std::nextafter(std::ldexp(1.0, exponent), + std::numeric_limits::infinity()); + }; + for (int exponent : {-470, -471}) + for (std::size_t lane = 0; lane < 4; ++lane) { + auto input = GeometryDyadicFixture(); + input.gradient_stress[0][0] = {odd(-500), 0}; + input.gradient_b[0].values[lane] = {odd(exponent), 0}; + ExpectGeometryPacketBitsAndFlags(input, false); + checks++; + } + RecordProperty("avx2_geometry_special_tail_checks", static_cast(checks)); +} + +TEST(Hex8HighLow, GeometryPacketPreservesOrderedOverflowAndPreparedInvalid) { + if (fesa::hex8_internal::SelectMaterialPacketEvaluator() == + fesa::hex8_internal::EvaluateMaterialPacketScalar) + GTEST_SKIP() << "Geometry AVX2 masked environment is not selected."; + const double maximum = std::numeric_limits::max(); + std::size_t checks = 0; + for (bool overflow_order : {false, true}) { + SCOPED_TRACE(overflow_order); + fesa::hex8_internal::GeometryPacketInput input{}; + input.weight = {0.5, 0}; + for (std::size_t l = 0; l < 3; ++l) + for (std::size_t lane = 0; lane < 4; ++lane) + input.gradient_b[l].values[lane] = {1, 0}; + input.gradient_stress[0][0] = {maximum, 0}; + input.gradient_stress[0][1] = {overflow_order ? maximum : -maximum, 0}; + input.gradient_stress[0][2] = {overflow_order ? -maximum : maximum, 0}; + ExpectGeometryPacketBitsAndFlags(input, false); + ExpectGeometryPacketBitsAndFlags(input, true); + checks += 2; + } + auto input = GeometryDyadicFixture(); + const auto invalid_first = fesa::high_low_original::Multiply( + {std::ldexp(1.0, -600), 0}, {std::ldexp(1.0, -500), 0}); + input.gradient_stress[0][0] = {invalid_first.high, invalid_first.low}; + for (std::size_t lane = 0; lane < 4; ++lane) + input.gradient_b[0].values[lane] = {std::ldexp(1.0, 600), 0}; + ExpectGeometryPacketBitsAndFlags(input, false); + ExpectGeometryPacketBitsAndFlags(input, true); + checks += 2; + RecordProperty("avx2_geometry_ordered_invalid_checks", static_cast(checks)); +} + +TEST(Hex8HighLow, GeometryPacketSumFallbackPreservesNearMaxAndSubnormalStatus) { + if (fesa::hex8_internal::SelectMaterialPacketEvaluator() == + fesa::hex8_internal::EvaluateMaterialPacketScalar) + GTEST_SKIP() << "Geometry AVX2 masked environment is not selected."; + const double maximum = std::numeric_limits::max(); + const double normal = std::ldexp(1.0, -1022); + const double tiny = std::ldexp(1.0, -1074); + const double quarter = std::ldexp(1.0, 969); + std::size_t checks = 0; + for (double tail : {0.0, normal, tiny}) + for (std::size_t bad = 0; bad < 4; ++bad) { + SCOPED_TRACE(tail); + SCOPED_TRACE(bad); + fesa::hex8_internal::GeometryPacketInput input{}; + input.weight = {0.5, 0}; + // Both pairs are normalized. Head sum is finite; middle contributes + // 2^970 and the joined sum overflows. Other lanes retain finite values. + input.gradient_stress[0][0] = {maximum, quarter}; + input.gradient_stress[0][1] = {quarter, tail}; + for (std::size_t lane = 0; lane < 4; ++lane) { + const double scale = lane == bad ? 1.0 : 0.5; + input.gradient_b[0].values[lane] = {scale, 0}; + input.gradient_b[1].values[lane] = {scale, 0}; + input.gradient_b[2].values[lane] = {1, 0}; + } + ExpectGeometryPacketBitsAndFlags(input, false); + ExpectGeometryPacketBitsAndFlags(input, true); + checks += 2; + } + fesa::hex8_internal::GeometryPacketInput finite{}; + finite.weight = {1, 0}; + finite.gradient_stress[0][0] = {normal, 0}; + finite.gradient_stress[0][1] = {normal + tiny, 0}; + for (std::size_t lane = 0; lane < 4; ++lane) { + finite.gradient_b[0].values[lane] = {1, 0}; + finite.gradient_b[1].values[lane] = {1, 0}; + } + // TwoSum(N,N+d) has finite normal inputs and a subnormal error d. + ExpectGeometryPacketBitsAndFlags(finite, false); + ExpectGeometryPacketBitsAndFlags(finite, true); + checks += 2; + // Mixed prefix: finite subnormal-error lanes run before another lane's head + // overflow. The bridge must preserve the whole packet's status union. + for (std::size_t bad = 0; bad < 4; ++bad) { + SCOPED_TRACE(bad); + fesa::hex8_internal::GeometryPacketInput mixed{}; + mixed.weight = {1, 0}; + mixed.gradient_stress[0][2] = {normal, 0}; + mixed.gradient_stress[1][2] = {normal + tiny, 0}; + mixed.gradient_stress[2][0] = {maximum, 0}; + mixed.gradient_stress[2][1] = {maximum, 0}; + for (std::size_t lane = 0; lane < 4; ++lane) { + mixed.gradient_b[0].values[lane] = {lane == bad ? 1.0 : 0.0, 0}; + mixed.gradient_b[1].values[lane] = {lane == bad ? 1.0 : 0.0, 0}; + mixed.gradient_b[2].values[lane] = {lane == bad ? 0.0 : 1.0, 0}; + } + ExpectGeometryPacketBitsAndFlags(mixed, false); + ExpectGeometryPacketBitsAndFlags(mixed, true); + checks += 2; + } + RecordProperty("avx2_geometry_sum_fallback_checks", static_cast(checks)); +} +#endif + TEST(Hex8HighLow, PacketBackendSelectionRecoversAfterEnvironmentRestoration) { const auto initial = fesa::hex8_internal::SelectMaterialPacketEvaluator(); const int rounding = std::fegetround();