From 127286cf2b61fbd61116769ee059d569fb94ebd8 Mon Sep 17 00:00:00 2001 From: "KOKO\\Mimi" Date: Sat, 1 Aug 2026 23:36:17 +0900 Subject: [PATCH] =?UTF-8?q?feat(deterministic-parallel-assembly):=20step?= =?UTF-8?q?=202=20=E2=80=94=20thread-count-determinism?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- CMakeLists.txt | 10 + tests/CMakeLists.txt | 36 +++ .../thread_count_determinism_test.cpp | 211 ++++++++++++++++++ tests/performance/assembly_benchmark.cpp | 142 ++++++++++++ 4 files changed, 399 insertions(+) create mode 100644 tests/integration/assembly/thread_count_determinism_test.cpp create mode 100644 tests/performance/assembly_benchmark.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 388f957..1a6d26b 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -67,5 +67,15 @@ target_compile_features(fesa PRIVATE cxx_std_20) target_compile_options(fesa PRIVATE /W4 /permissive- /EHsc) if(BUILD_TESTING) + add_executable(fesa_assembly_benchmark + tests/performance/assembly_benchmark.cpp + ) + + target_link_libraries(fesa_assembly_benchmark PRIVATE fesa_core) + target_compile_features(fesa_assembly_benchmark PRIVATE cxx_std_20) + target_compile_options( + fesa_assembly_benchmark PRIVATE /W4 /permissive- /EHsc + ) + add_subdirectory(tests) endif() diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 1aad74c..ca3f6b2 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -446,6 +446,42 @@ set_property( ${FESA_DEPENDENCY_RUNTIME_MODIFICATIONS} ) +add_executable(fesa_thread_count_determinism_tests + integration/assembly/thread_count_determinism_test.cpp +) + +target_compile_features( + fesa_thread_count_determinism_tests PRIVATE cxx_std_20 +) +target_compile_options( + fesa_thread_count_determinism_tests + PRIVATE + /W4 + /permissive- + /EHsc +) + +target_link_libraries(fesa_thread_count_determinism_tests + PRIVATE + fesa_core + GTest::gtest_main +) + +add_test( + NAME ThreadCountDeterminism + COMMAND "$" +) + +set_property( + TEST ThreadCountDeterminism + PROPERTY ENVIRONMENT "MKL_NUM_THREADS=1" +) +set_property( + TEST ThreadCountDeterminism + PROPERTY ENVIRONMENT_MODIFICATION + ${FESA_DEPENDENCY_RUNTIME_MODIFICATIONS} +) + add_executable(fesa_constraint_tests unit/constraints/essential_bc_test.cpp ) diff --git a/tests/integration/assembly/thread_count_determinism_test.cpp b/tests/integration/assembly/thread_count_determinism_test.cpp new file mode 100644 index 0000000..d0c9f47 --- /dev/null +++ b/tests/integration/assembly/thread_count_determinism_test.cpp @@ -0,0 +1,211 @@ +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +#include + +namespace { + +constexpr std::size_t branch_count = 24; +constexpr std::size_t repeat_count = 10; + +fesa::Domain build_branched_cantilever() { + fesa::DomainBuilder builder; + builder.add_node({ + fesa::NodeId{0}, + fesa::EntityOrigin{"BeamPart", "Beam-1", 1}, + fesa::Vec3{0.0, 0.0, 0.0}, + }); + builder.add_node({ + fesa::NodeId{1}, + fesa::EntityOrigin{"BeamPart", "Beam-1", 2}, + fesa::Vec3{1.0, 0.0, 0.0}, + }); + for (std::size_t index = 0; index < branch_count; ++index) { + const double y = static_cast( + static_cast(index % 7) - 3) * 0.35; + const double z = static_cast( + static_cast((index * 3) % 11) - 5) * 0.22 + 0.1; + builder.add_node({ + fesa::NodeId{static_cast(index + 2)}, + fesa::EntityOrigin{ + "BeamPart", + "Beam-1", + static_cast(index + 3), + }, + fesa::Vec3{ + 2.0 + static_cast(index) * 0.05, + y, + z, + }, + }); + } + builder.add_material({ + fesa::MaterialId{0}, + "Steel", + 210.0e9, + 0.3, + }); + builder.add_section({ + fesa::SectionId{0}, + "General", + 0.02, + 3.0e-5, + 4.0e-5, + 2.0e-5, + 0.015, + 0.016, + fesa::ShearPropertySource::input, + fesa::Vec3{0.0, 1.0, 0.0}, + {}, + }); + + for (std::size_t offset = 0; offset < branch_count; ++offset) { + const std::size_t index = branch_count - offset - 1; + builder.add_beam_element({ + fesa::ElementId{static_cast(index + 1)}, + fesa::EntityOrigin{ + "BeamPart", + "Beam-1", + static_cast(index + 2), + }, + { + fesa::NodeId{1}, + fesa::NodeId{static_cast(index + 2)}, + }, + fesa::MaterialId{0}, + fesa::SectionId{0}, + }); + } + builder.add_beam_element({ + fesa::ElementId{0}, + fesa::EntityOrigin{"BeamPart", "Beam-1", 1}, + {fesa::NodeId{0}, fesa::NodeId{1}}, + fesa::MaterialId{0}, + fesa::SectionId{0}, + }); + + std::vector fixed; + fixed.reserve(6); + for (std::uint8_t dof = 1; dof <= 6; ++dof) { + fixed.push_back({fesa::NodeId{0}, dof, 0.0}); + } + builder.set_step({ + "Load", + std::move(fixed), + {{ + fesa::NodeId{static_cast(branch_count + 1)}, + {1000.0, -2500.0, 1750.0, 20.0, -30.0, 40.0}, + }}, + }); + + auto result = std::move(builder).build(); + if (!result.domain.has_value()) { + throw std::runtime_error{"Test Domain failed validation."}; + } + return std::move(*result.domain); +} + +std::vector bits(const std::vector& values) { + std::vector result; + result.reserve(values.size()); + for (const double value : values) { + result.push_back(std::bit_cast(value)); + } + return result; +} + +void expect_bitwise_equal( + const fesa::EquationSystem& expected, + const fesa::EquationSystem& actual) { + EXPECT_EQ(actual.stiffness.order, expected.stiffness.order); + EXPECT_EQ(actual.stiffness.row_offsets, expected.stiffness.row_offsets); + EXPECT_EQ( + actual.stiffness.column_indices, + expected.stiffness.column_indices); + EXPECT_EQ(bits(actual.stiffness.values), bits(expected.stiffness.values)); + EXPECT_EQ(bits(actual.force), bits(expected.force)); +} + +struct LinearState final { + std::vector displacement; + std::vector reaction; +}; + +LinearState solve_linear_system( + const fesa::EquationSystem& system, + const fesa::DofManager& dofs) { + fesa::ConstraintResult constrained = + fesa::eliminate_essential_bcs(system, dofs); + if (!constrained.reduced_system.has_value()) { + throw std::runtime_error{"Constraint elimination failed."}; + } + + fesa::PardisoLinearSolver solver; + fesa::LinearSolveResult solved = solver.solve( + constrained.reduced_system->stiffness, + constrained.reduced_system->force); + if (!solved.diagnostics.empty()) { + throw std::runtime_error{"Linear solve failed."}; + } + + std::vector displacement = + dofs.reconstruct_full(solved.solution); + return { + displacement, + fesa::recover_reaction(system, displacement), + }; +} + +std::size_t available_concurrency() { + return std::max( + std::size_t{1}, + static_cast(std::thread::hardware_concurrency())); +} + +TEST(ThreadCountDeterminism, AssemblyAndLinearStateAreBitwiseStable) { + const fesa::Domain domain = build_branched_cantilever(); + const fesa::DofManager dofs = fesa::DofManager::build(domain); + const fesa::EquationSystem serial = fesa::assemble_serial(domain, dofs); + const LinearState serial_state = solve_linear_system(serial, dofs); + const std::array thread_counts{ + 1, + 2, + available_concurrency(), + }; + + for (std::size_t repeat = 0; repeat < repeat_count; ++repeat) { + for (const std::size_t threads : thread_counts) { + SCOPED_TRACE(::testing::Message{} + << "repeat=" << repeat << ", threads=" << threads); + const fesa::EquationSystem parallel = + fesa::assemble_parallel(domain, dofs, {threads, 1}); + expect_bitwise_equal(serial, parallel); + + const LinearState parallel_state = + solve_linear_system(parallel, dofs); + EXPECT_EQ( + bits(parallel_state.displacement), + bits(serial_state.displacement)); + EXPECT_EQ( + bits(parallel_state.reaction), + bits(serial_state.reaction)); + } + } +} + +} // namespace diff --git a/tests/performance/assembly_benchmark.cpp b/tests/performance/assembly_benchmark.cpp new file mode 100644 index 0000000..85ac38f --- /dev/null +++ b/tests/performance/assembly_benchmark.cpp @@ -0,0 +1,142 @@ +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +namespace { + +constexpr std::size_t benchmark_element_count = 1000; + +fesa::Domain build_beam_chain() { + fesa::DomainBuilder builder; + for (std::size_t index = 0; index <= benchmark_element_count; ++index) { + builder.add_node({ + fesa::NodeId{static_cast(index)}, + fesa::EntityOrigin{ + "BenchmarkPart", + "Benchmark-1", + static_cast(index + 1), + }, + fesa::Vec3{static_cast(index), 0.0, 0.0}, + }); + } + builder.add_material({ + fesa::MaterialId{0}, + "Steel", + 210.0e9, + 0.3, + }); + builder.add_section({ + fesa::SectionId{0}, + "General", + 0.02, + 3.0e-5, + 4.0e-5, + 2.0e-5, + 0.015, + 0.016, + fesa::ShearPropertySource::input, + fesa::Vec3{0.0, 1.0, 0.0}, + {}, + }); + for (std::size_t index = 0; index < benchmark_element_count; ++index) { + builder.add_beam_element({ + fesa::ElementId{static_cast(index)}, + fesa::EntityOrigin{ + "BenchmarkPart", + "Benchmark-1", + static_cast(index + 1), + }, + { + fesa::NodeId{static_cast(index)}, + fesa::NodeId{static_cast(index + 1)}, + }, + fesa::MaterialId{0}, + fesa::SectionId{0}, + }); + } + builder.set_step({"Benchmark", {}, {}}); + + auto result = std::move(builder).build(); + if (!result.domain.has_value()) { + throw std::runtime_error{"Benchmark Domain failed validation."}; + } + return std::move(*result.domain); +} + +template +std::pair measure(Operation&& operation) { + const auto start = std::chrono::steady_clock::now(); + fesa::EquationSystem result = std::forward(operation)(); + const auto finish = std::chrono::steady_clock::now(); + const double milliseconds = + std::chrono::duration(finish - start).count(); + return {std::move(result), milliseconds}; +} + +std::vector bits(const std::vector& values) { + std::vector result; + result.reserve(values.size()); + for (const double value : values) { + result.push_back(std::bit_cast(value)); + } + return result; +} + +bool same_system( + const fesa::EquationSystem& left, + const fesa::EquationSystem& right) { + return left.stiffness.order == right.stiffness.order && + left.stiffness.row_offsets == right.stiffness.row_offsets && + left.stiffness.column_indices == right.stiffness.column_indices && + bits(left.stiffness.values) == bits(right.stiffness.values) && + bits(left.force) == bits(right.force); +} + +std::size_t available_concurrency() { + const std::size_t available = + static_cast(std::thread::hardware_concurrency()); + return available == 0 ? 1 : available; +} + +} // namespace + +int main() { + try { + const fesa::Domain domain = build_beam_chain(); + const fesa::DofManager dofs = fesa::DofManager::build(domain); + const std::size_t threads = available_concurrency(); + + auto [serial, serial_ms] = measure( + [&] { return fesa::assemble_serial(domain, dofs); }); + auto [parallel, parallel_ms] = measure([&] { + return fesa::assemble_parallel(domain, dofs, {threads, 1}); + }); + + if (!same_system(serial, parallel)) { + std::cerr << "Serial and parallel assembly results differ.\n"; + return 1; + } + + std::cout << "elements=" << domain.beam_elements().size() + << " serial_ms=" << serial_ms + << " parallel_ms=" << parallel_ms + << " parallel_threads=" << threads << '\n'; + return 0; + } catch (const std::exception& error) { + std::cerr << "Assembly benchmark failed: " << error.what() << '\n'; + return 1; + } +}