From c0f06e9f64d40c5fae89eafe8bd5bb205ad77d2e Mon Sep 17 00:00:00 2001 From: "KOKO\\Mimi" Date: Wed, 21 Oct 2026 07:04:15 +0900 Subject: [PATCH] test: add immutable shell Jet baseline oracle --- FESA.sln | 6 + msvc/fesa_jet_baseline_tests.vcxproj | 21 + tests/CMakeLists.txt | 19 + .../nonlinear_mitc4_shell_jet_baseline.cpp | 561 ++++++++++++++++++ .../nonlinear_mitc4_shell_jet_baseline.md | 40 ++ 5 files changed, 647 insertions(+) create mode 100644 msvc/fesa_jet_baseline_tests.vcxproj create mode 100644 tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.cpp create mode 100644 tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.md diff --git a/FESA.sln b/FESA.sln index 64295f3..03a1a75 100644 --- a/FESA.sln +++ b/FESA.sln @@ -14,6 +14,8 @@ Project("{BC8A1FFA-BEE3-4634-8014-F334798102B3}") = "fesa_reference_tests", "msv EndProject Project("{BC8A1FFA-BEE3-4634-8014-F334798102B3}") = "gtest", "msvc\gtest.vcxproj", "{03194BA4-3F3B-5399-BCAF-FB4FA47300CE}" EndProject +Project("{BC8A1FFA-BEE3-4634-8014-F334798102B3}") = "fesa_jet_baseline_tests", "msvc\fesa_jet_baseline_tests.vcxproj", "{5FD398DA-6058-4CE3-9440-C9E0FDD180C3}" +EndProject Global GlobalSection(SolutionConfigurationPlatforms) = preSolution Debug|x64 = Debug|x64 @@ -44,5 +46,9 @@ Global {03194BA4-3F3B-5399-BCAF-FB4FA47300CE}.Debug|x64.Build.0 = Debug|x64 {03194BA4-3F3B-5399-BCAF-FB4FA47300CE}.Release|x64.ActiveCfg = Release|x64 {03194BA4-3F3B-5399-BCAF-FB4FA47300CE}.Release|x64.Build.0 = Release|x64 + {5FD398DA-6058-4CE3-9440-C9E0FDD180C3}.Debug|x64.ActiveCfg = Debug|x64 + {5FD398DA-6058-4CE3-9440-C9E0FDD180C3}.Debug|x64.Build.0 = Debug|x64 + {5FD398DA-6058-4CE3-9440-C9E0FDD180C3}.Release|x64.ActiveCfg = Release|x64 + {5FD398DA-6058-4CE3-9440-C9E0FDD180C3}.Release|x64.Build.0 = Release|x64 EndGlobalSection EndGlobal diff --git a/msvc/fesa_jet_baseline_tests.vcxproj b/msvc/fesa_jet_baseline_tests.vcxproj new file mode 100644 index 0000000..0ab92c0 --- /dev/null +++ b/msvc/fesa_jet_baseline_tests.vcxproj @@ -0,0 +1,21 @@ + + + Debugx64Releasex64 + {5FD398DA-6058-4CE3-9440-C9E0FDD180C3}Win32Proj10.0 + + Applicationv145true + + + + + $(FesaRoot)src\fesa;$(FESA_MKL_ROOT)\include;$(FESA_TBB_ROOT)\include;$(FESA_HDF5_ROOT)\include;%(AdditionalIncludeDirectories) + MKL_ILP64;__TBB_NO_IMPLICIT_LINKAGE=1;H5_BUILT_AS_DYNAMIC_LIB;%(PreprocessorDefinitions) + + + + + {3AEC39D9-4641-5161-9107-265EB50CA882}true + + + + diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index af52699..94358c3 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -179,6 +179,17 @@ target_compile_definitions( FESA_TEST_BINARY_DIR="${PROJECT_BINARY_DIR}" ) +add_executable(fesa_jet_baseline_tests + unit/elements/nonlinear_mitc4_shell_jet_baseline.cpp) +target_include_directories(fesa_jet_baseline_tests PRIVATE + "${PROJECT_SOURCE_DIR}/src/fesa") +target_compile_definitions(fesa_jet_baseline_tests PRIVATE + MKL_ILP64 __TBB_NO_IMPLICIT_LINKAGE=1 H5_BUILT_AS_DYNAMIC_LIB) +target_compile_options(fesa_jet_baseline_tests PRIVATE + $<$:/fp:precise /EHsc>) +target_link_libraries(fesa_jet_baseline_tests PRIVATE + fesa_solver Fesa::MKL Fesa::TBB Fesa::HDF5) + if(WIN32) if(NOT TARGET ${FESA_MKL_THREAD_TARGET} OR NOT TARGET MKL::mkl_core) message(FATAL_ERROR "The oneMKL runtime files required by tests were not resolved") @@ -226,6 +237,7 @@ if(WIN32) fesa_stage_test_runtime(fesa_unit_tests) fesa_stage_test_runtime(fesa_integration_tests) fesa_stage_test_runtime(fesa_reference_tests) + fesa_stage_test_runtime(fesa_jet_baseline_tests) unset(_fesa_mkl_dispatch_dll) unset(_fesa_mkl_dispatch_dlls) @@ -258,4 +270,11 @@ add_custom_target( fesa_unit_tests fesa_integration_tests fesa_reference_tests + fesa_jet_baseline_tests ) + +# Verification never creates or replaces an immutable baseline snapshot. +add_test(NAME NonlinearMitc4Jet.IndependentOracle + COMMAND fesa_jet_baseline_tests --verify) +set_tests_properties(NonlinearMitc4Jet.IndependentOracle PROPERTIES + LABELS "unit;jet-baseline") diff --git a/tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.cpp b/tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.cpp new file mode 100644 index 0000000..746c9b8 --- /dev/null +++ b/tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.cpp @@ -0,0 +1,561 @@ +// Standalone private baseline fixture: compile this translation unit instead of +// nonlinear_mitc4_shell.cpp, then link the ordinary solver archive. No +// production source, public interface, or official CLI instrumentation is +// needed. Build target: fesa_jet_baseline_tests (MSBuild/CMake). --verify runs +// independent checks without writing evidence. Compare immutable evidence with +// --compare PATH, Debug to Debug or Release to Release; never regenerate a +// snapshot during regression. See nonlinear_mitc4_shell_jet_baseline.md. +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include "../../../src/fesa/elements/nonlinear_mitc4_shell.cpp" + +namespace fesa { +namespace { +using Coefficients = std::array; +std::fenv_t initial_environment; +unsigned initial_mxcsr; + +std::string FpState() { + unsigned control = 0; + _controlfp_s(&control, 0, 0); + std::ostringstream out; + out << std::fetestexcept(FE_ALL_EXCEPT) << ':' << std::fegetround() << ':' + << _statusfp() << ':' << control << ':' << _mm_getcsr(); + return out.str(); +} + +std::string ResetFp() { + std::fesetenv(&initial_environment); + std::feclearexcept(FE_ALL_EXCEPT); + _clearfp(); + _mm_setcsr(initial_mxcsr & ~0x3fu); + return FpState(); +} + +struct Monomial { + double coefficient; + std::array powers{}; +}; + +// Independent oracle: differentiate monomial exponent lists, never Jet rules. +double Differentiate(const std::vector& polynomial, + const Coefficients& point, + const std::vector& derivative) { + long double sum = 0; + for (auto term : polynomial) { + long double value = term.coefficient; + for (auto index : derivative) { + value *= term.powers[index]; + if (term.powers[index] == 0) { + value = 0; + break; + } + --term.powers[index]; + } + for (std::size_t i = 0; i < 24; ++i) + for (unsigned p = 0; p < term.powers[i]; ++p) value *= point[i]; + sum += value; + } + return static_cast(sum); +} + +std::uint64_t Bits(double value) { + std::uint64_t result; + std::memcpy(&result, &value, sizeof(result)); + return result; +} +void Check(double actual, double expected, const char* channel) { + if (!std::isfinite(actual) || + std::abs(actual - expected) > 2e-12 * std::max(1.0, std::abs(expected))) + throw std::runtime_error(std::string("independent oracle: ") + channel); +} +Jet Polynomial(const std::vector& terms, + const std::array& variables) { + Jet result = Constant(0); + for (const auto& term : terms) { + Jet value = Constant(term.coefficient); + for (std::size_t i = 0; i < 24; ++i) + for (unsigned p = 0; p < term.powers[i]; ++p) + value = value * variables[i]; + result = result + value; + } + return result; +} +std::string Snapshot(const char* name, const Jet& jet, int exceptions, + unsigned status, unsigned mxcsr) { + std::ostringstream out; + out << name << ' ' << jet.constant_derivatives << ' ' + << jet.second.IsImplicitZero() << ' ' << IsFinite(jet) << ' ' + << exceptions << ' ' << status << ' ' << mxcsr << std::hex; + auto put = [&](double x) { + if (std::isnan(x)) + out << " nan"; + else + out << ' ' << Bits(x); + }; + put(jet.value); + for (auto x : jet.first) put(x); + for (auto x : jet.second) put(x); + put(jet.directional_first); + put(jet.directional_second); + for (auto x : jet.mixed_first) put(x); + for (auto x : jet.mixed_second) put(x); + return out.str() + "\n"; +} +template +std::string Observe(const char* name, F expression, Jet* returned = nullptr) { + const auto before = ResetFp(); + Jet value = expression(); + const int exceptions = std::fetestexcept(FE_ALL_EXCEPT); + const unsigned status = _statusfp(); + const unsigned mxcsr = _mm_getcsr(); + const auto after = FpState(); + if (returned != nullptr) *returned = value; + return std::string("initial-fp=") + before + " final-fp=" + after + " " + + Snapshot(name, value, exceptions, status, mxcsr); +} + +std::string Verify(bool corrupt_oracle) { + Coefficients x{}, direction{}; + std::array variables{}; + std::vector terms; + for (std::size_t i = 0; i < 24; ++i) { + x[i] = (static_cast(i % 7) - 3) / 8; + direction[i] = (static_cast(i % 5) - 2) / 16; + variables[i] = Variable(x[i], i, direction[i]); + Monomial term{(i % 2 ? -1.0 : 1.0) / 8}; + term.powers[i] = 3; + term.powers[(i + 1) % 24] += 2; + terms.push_back(term); + Monomial cross{0.5}; + cross.powers[i] = 1; + cross.powers[(i + 11) % 24] += 1; + terms.push_back(cross); + } + Jet polynomial = Polynomial(terms, variables); + if (corrupt_oracle) polynomial.value += 1; // Deliberate fixture-only mutant. + Check(polynomial.value, Differentiate(terms, x, {}), "value"); + double d1 = 0, d2 = 0; + for (std::size_t i = 0; i < 24; ++i) { + const double gradient = Differentiate(terms, x, {i}); + Check(polynomial.first[i], gradient, "gradient"); + d1 += gradient * direction[i]; + double mixed1 = 0, mixed2 = 0; + for (std::size_t j = 0; j < 24; ++j) { + const double hessian = Differentiate(terms, x, {i, j}); + Check(polynomial.second.ReadData()[24 * i + j], hessian, "Hessian"); + d2 += hessian * direction[i] * direction[j]; + mixed1 += hessian * direction[j]; + for (std::size_t k = 0; k < 24; ++k) + mixed2 += + Differentiate(terms, x, {i, j, k}) * direction[j] * direction[k]; + } + Check(polynomial.mixed_first[i], mixed1, "mixed_first"); + Check(polynomial.mixed_second[i], mixed2, "mixed_second"); + } + Check(polynomial.directional_first, d1, "directional_first"); + Check(polynomial.directional_second, d2, "directional_second"); + // Affine composition has a rank-one Hessian; its third derivative contracts + // independently to the mixed-second channel. + Jet affine = Constant(2); + Coefficients weights{}; + double slope = 0; + for (std::size_t i = 0; i < 24; ++i) { + weights[i] = (static_cast(i % 3) - 1) / 64; + affine = affine + weights[i] * variables[i]; + slope += weights[i] * direction[i]; + } + auto unary = [&](const Jet& jet, double value, double a, double b, double c) { + Check(jet.value, value, "unary value"); + Check(jet.directional_first, a * slope, "unary direction1"); + Check(jet.directional_second, b * slope * slope, "unary direction2"); + for (std::size_t i = 0; i < 24; ++i) { + Check(jet.first[i], a * weights[i], "unary gradient"); + Check(jet.mixed_first[i], b * weights[i] * slope, "unary mixed1"); + Check(jet.mixed_second[i], c * weights[i] * slope * slope, + "unary mixed2"); + for (std::size_t j = 0; j < 24; ++j) + Check(jet.second.ReadData()[24 * i + j], b * weights[i] * weights[j], + "unary Hessian"); + } + }; + const double t = affine.value, r = std::sqrt(t); + unary(Sin(affine), std::sin(t), std::cos(t), -std::sin(t), -std::cos(t)); + unary(Cos(affine), std::cos(t), -std::sin(t), -std::cos(t), std::sin(t)); + unary(Reciprocal(affine), 1 / t, -1 / (t * t), 2 / (t * t * t), + -6 / (t * t * t * t)); + unary(Sqrt(affine), r, 0.5 / r, -0.25 / (r * r * r), + 0.375 / (r * r * r * r * r)); + const double den = 1 + t * t; + unary(Atan2(affine, Constant(1)), std::atan(t), 1 / den, -2 * t / (den * den), + (6 * t * t - 2) / (den * den * den)); + const Jet original = polynomial; + const Jet alias = original; + const auto before = Snapshot("alias", alias, 0, 0, 0); + const Jet transformed = alias * alias + original; + if (before != Snapshot("alias", original, 0, 0, 0)) + throw std::runtime_error("alias modified"); + std::string out; + std::vector>> futures; + for (int i = 0; i < 4; ++i) + futures.push_back(std::async(std::launch::async, [&] { + Jet value; + const auto evidence = Observe( + "shared-full", [&] { return original * original + original; }, + &value); + return std::make_pair(std::move(value), evidence); + })); + for (std::size_t i = 0; i < futures.size(); ++i) { + const auto worker = futures[i].get(); + if (Snapshot("shared", worker.first, 0, 0, 0) != + Snapshot("shared", transformed, 0, 0, 0)) + throw std::runtime_error("shared complete-channel snapshot"); + out += "worker-" + std::to_string(i) + " " + worker.second; + } + auto add = [&](const char* name, auto f) { out += Observe(name, f); }; + add("constant", [] { return Constant(-0.0); }); + add("variable-zero-direction", [] { return Variable(-0.0, 23); }); + add("variable-direction", [] { return Variable(0.25, 0, -0.5); }); + add("polynomial", [&] { return Polynomial(terms, variables); }); + add("polynomial-zero-direction", [&] { + auto zero_direction = variables; + for (auto& variable : zero_direction) variable.directional_first = 0; + return Polynomial(terms, zero_direction); + }); + add("actual-zero-hessian", [&] { return original - original; }); + add("implicit-zero-nonconstant", + [] { return Variable(1, 0) - Variable(1, 0); }); + add("implicit-add", [] { return Variable(1, 0) + Constant(-1); }); + add("implicit-subtract", [] { return Variable(1, 0) - Constant(1); }); + for (int kind = 0; kind < 4; ++kind) { + const auto name = "surviving-copy-" + std::to_string(kind); + add(name.c_str(), [&] { + Jet survivor; + std::string expected; + { + Jet source = kind == 0 ? Constant(-0.0) + : kind == 1 ? Variable(0.5, 7, -0.25) + : Polynomial(terms, variables); + if (kind == 3) source = source - source; + const Jet another_alias = source; + expected = Snapshot("survivor", source, 0, 0, 0); + survivor = another_alias; + } + if (Snapshot("survivor", survivor, 0, 0, 0) != expected) + throw std::runtime_error("surviving copy complete-channel snapshot"); + return survivor; + }); + } + add("constant-multiply", [] { return Constant(-0.0) * Constant(3); }); + add("constant-compose", [] { return Sin(Constant(0)); }); + add("scale", [&] { return -0.5 * original; }); + add("multiply", [&] { return original * original; }); + add("compose-sin", [&] { return Sin(original); }); + add("compose-cos", [&] { return Cos(original); }); + add("divide", [&] { return original / affine; }); + add("sqrt", [&] { return Sqrt(affine); }); + add("atan2", [&] { return Atan2(original, affine); }); + add("dot", [&] { + return Dot(JetVector3{original, affine, variables[0]}, + JetVector3{affine, variables[1], original}); + }); + add("cancellation", + [] { return (Variable(1e16, 0, 1) + Constant(1)) - Constant(1e16); }); + add("overflow", [] { + return Variable(std::numeric_limits::max(), 0, 1) * Constant(2); + }); + add("subnormal", [] { + return std::numeric_limits::denorm_min() * Variable(0.5, 0, 0.5); + }); + add("nan", [] { + return Constant(std::numeric_limits::quiet_NaN()) * Variable(1, 0); + }); + add("inf", [] { + return Constant(std::numeric_limits::infinity()) * Constant(0); + }); + add("singular", [] { return Reciprocal(Variable(0, 0, 1)); }); + add("singular-zero-direction", [] { return Reciprocal(Variable(0, 0)); }); + add("constant-singular-reciprocal", [] { return Reciprocal(Constant(0)); }); + add("constant-singular-sqrt", [] { return Sqrt(Constant(0)); }); + add("negative-sqrt", [] { return Sqrt(Variable(-1, 0, 1)); }); + add("after-failure", [&] { return Polynomial(terms, variables); }); + // Exercise the production MKL k=1 outer path directly, including entries + // whose rounded bits agree while their IEEE exception state can differ. + auto outer = [](double a, double b) { + Jet jet{}; + Coefficients left{}, right{}; + for (std::size_t i = 0; i < 24; ++i) { + left[i] = i % 2 ? -a : a; + right[i] = i % 3 ? -b : b; + } + Outer(left, right, jet.second.OverwriteData()); + return jet; + }; + add("outer-signed-zero", [&] { return outer(-0.0, 1); }); + add("outer-underflow", + [&] { return outer(std::numeric_limits::denorm_min(), 0.5); }); + add("outer-overflow", + [&] { return outer(std::numeric_limits::max(), 2); }); + add("outer-invalid", + [&] { return outer(std::numeric_limits::infinity(), 0); }); + if (!Constant(1).constant_derivatives || + !Constant(1).second.IsImplicitZero() || + Variable(1, 0).constant_derivatives || + !Variable(1, 0).second.IsImplicitZero() || + polynomial.second.IsImplicitZero()) + throw std::runtime_error("provenance"); + return out; +} + +// Public request-mode evidence for later direction-channel separation. Existing +// tangent tests check nonzero-direction derivatives but do not freeze the full +// omitted-direction response inventory and invalid-direction diagnostics. +std::string PublicSnapshots() { + std::array nodes; + const std::array, 4> coordinates{ + {{{-1, -1, 0}}, {{1, -1, 0}}, {{1, 1, 0}}, {{-1, 1, 0}}}}; + std::array pointers{}; + std::array frames{}; + NonlinearMitc4TrialDofs trial{}; + Vector zero(24), direction(24), wrong(23), nan(24), inf(24); + for (std::size_t i = 0; i < 4; ++i) { + nodes[i] = {{"Jet-Baseline", static_cast(i + 1), + std::to_string(i + 1)}, + coordinates[i], + {"jet-baseline.inp", i + 1}}; + pointers[i] = &nodes[i]; + frames[i] = {static_cast(i), {0, 0, 1}, {1, 0, 0}, {0, 1, 0}}; + trial.displacement[i] = {0.01 * i, -0.02 * i, 0.015 * i}; + trial.rotation_vector[i] = {0.02 * i, -0.01 * i, 0.005 * i}; + } + for (std::size_t i = 0; i < 24; ++i) + direction[i] = (static_cast(i % 5) - 2) / 8; + nan[7] = std::numeric_limits::quiet_NaN(); + inf[7] = std::numeric_limits::infinity(); + const auto created = NonlinearMitc4Shell::Create( + pointers, frames, {"Section", 0.2, 0, {"jet-baseline.inp", 20}}, + {"Material", 120, 0.25, {"jet-baseline.inp", 21}}); + if (!created.HasValue()) throw std::runtime_error("public fixture creation"); + std::ostringstream out; + for (const auto mode : {Mitc4DrillingMode::kDirectSpring, + Mitc4DrillingMode::kLegacyCorotatedTwist}) { + const auto shell = + mitc4_drilling_internal::TestAccess::WithMode(created.Value(), mode); + auto observe = [&](const char* name, const NonlinearMitc4TrialDofs& input, + const Vector* request, + std::optional expected_success) { + const auto before = ResetFp(); + const auto result = + EvaluateInternalResidualAndTangent(shell, input, request); + const int exceptions = std::fetestexcept(FE_ALL_EXCEPT); + const unsigned status = _statusfp(), mxcsr = _mm_getcsr(); + if (expected_success && result.HasValue() != *expected_success) + throw std::runtime_error(std::string("public outcome: ") + name); + out << "public " << static_cast(mode) << ' ' << name << ' ' + << result.HasValue() << ' ' << exceptions << ' ' << status << ' ' + << mxcsr << " initial-fp=" << before << " final-fp=" << FpState(); + auto put = [&](double value) { + if (!std::isfinite(value)) + throw std::runtime_error("public nonfinite field"); + out << ' ' << std::hex << Bits(value) << std::dec; + }; + auto vector = [&](const Vector& value) { + out << " v" << value.Size(); + for (std::size_t i = 0; i < value.Size(); ++i) put(value[i]); + }; + auto matrix = [&](const Matrix& value) { + out << " m" << value.Rows() << 'x' << value.Columns(); + for (std::size_t i = 0; i < value.Rows() * value.Columns(); ++i) + put(value.Data()[i]); + }; + if (result.HasValue()) { + const auto& r = result.Value(); + vector(r.internal_residual); + matrix(r.internal_tangent); + vector(r.internal_force); + put(r.physical_strain_energy); + vector(r.physical_force); + vector(r.drilling_force); + matrix(r.drilling_tangent); + matrix(r.drilling_metric); + put(r.drilling_energy); + vector(r.force_directional_second); + vector(r.drilling_force_directional_second); + vector(r.drilling_constraint); + matrix(r.drilling_constraint_gradient); + for (const auto& hessian : r.drilling_constraint_hessian) + matrix(hessian); + for (const auto& third : + r.drilling_constraint_third_directional_gradient) + vector(third); + put(r.recovery.physical_strain_energy); + for (const auto& point : r.recovery.gauss_points) { + out << " location=" << static_cast(point.location); + for (double value : point.natural_coordinates) put(value); + for (const auto& axis : {point.local_frame.e1, point.local_frame.e2, + point.local_frame.e3}) + for (double value : axis) put(value); + for (double value : point.generalized_strain) put(value); + for (double value : point.section_resultant) put(value); + for (const auto& stress : point.stress) { + out << " section=" << static_cast(stress.position); + put(stress.zeta); + for (double value : stress.components) put(value); + } + } + } else { + const auto& failure = result.GetStatus(); + out << " category=" + << (failure.Category() ? static_cast(*failure.Category()) + : -1); + for (const auto& diagnostic : failure.Diagnostics()) + out << ' ' << static_cast(diagnostic.severity) << ' ' + << std::quoted(diagnostic.code) << ' ' + << std::quoted(diagnostic.location.file.generic_string()) << ' ' + << diagnostic.location.line << ' ' + << std::quoted(diagnostic.keyword) << ' ' + << std::quoted(diagnostic.entity_identity) << ' ' + << std::quoted(diagnostic.message); + } + out << '\n'; + }; + observe("none", trial, nullptr, true); + observe("zero", trial, &zero, true); + observe("nonzero", trial, &direction, true); + observe("invalid-dimension", trial, &wrong, false); + observe("nan-direction", trial, &nan, false); + observe("inf-direction", trial, &inf, false); + auto invalid = trial; + invalid.displacement[1] = {std::numeric_limits::quiet_NaN(), 0, 0}; + observe("nan-trial-none", invalid, nullptr, false); + observe("none-after-failure", trial, nullptr, true); + auto requests = [&](const std::string& name, + const NonlinearMitc4TrialDofs& input) { + observe((name + "-none").c_str(), input, nullptr, std::nullopt); + observe((name + "-zero").c_str(), input, &zero, std::nullopt); + observe((name + "-nonzero").c_str(), input, &direction, std::nullopt); + }; + requests("zero-rotation", NonlinearMitc4TrialDofs{}); + for (const double threshold : {1e-8, 0.25}) { + const double angle = std::sqrt(threshold); + const std::array angles{ + std::nextafter(angle, 0.0), angle, + std::nextafter(angle, std::numeric_limits::infinity())}; + for (std::size_t i = 0; i < angles.size(); ++i) { + NonlinearMitc4TrialDofs boundary{}; + for (auto& rotation : boundary.rotation_vector) + rotation = {angles[i], 0, 0}; + requests((threshold == 0.25 ? "rotation-squared-quarter-" + : "rotation-squared-small-") + + std::to_string(i), + boundary); + } + } + auto translated = trial; + for (auto& displacement : translated.displacement) + displacement = displacement + Vector3{1e9, -2e9, 3e9}; + requests("large-translation", translated); + } + return out.str(); +} + +__declspec(noinline) double Consume(const Jet& value) { + double result = + value.value + value.directional_first + value.directional_second; + for (double entry : value.first) result += entry; + for (double entry : value.second) result += entry; + for (double entry : value.mixed_first) result += entry; + for (double entry : value.mixed_second) result += entry; + return result + value.constant_derivatives + value.second.IsImplicitZero(); +} + +void Bench() { + volatile double consumed = 0; + for (int kind = 0; kind < 8; ++kind) { + Jet left = + kind == 0 ? Constant(0.5) : Variable(0.5, 0, kind % 2 ? 0.0 : 0.25); + Jet right = Variable(0.75, 1, kind == 0 || kind % 2 ? 0.0 : -0.5); + if (kind >= 3) { + left = left * left; + right = right * right; + } + if (kind == 5) left = left - left; + if (kind == 7) { + left = Constant(0.5); + right = Constant(0.75); + } + const auto start = std::chrono::steady_clock::now(); + for (int i = 0; i < 20000; ++i) { + left.value += (i % 2 ? -1.0 : 1.0) / 1048576; + const Jet value = kind == 6 ? Dot(JetVector3{left, right, left}, + JetVector3{right, left, right}) + : left * right; + consumed = consumed + Consume(value); + } + std::cout << "bench," << kind << ",20000," + << std::chrono::duration( + std::chrono::steady_clock::now() - start) + .count() + << ",allocation-included," << consumed << '\n'; + } +} +} // namespace +} // namespace fesa + +int main(int argc, char** argv) { + try { + std::fegetenv(&fesa::initial_environment); + fesa::initial_mxcsr = _mm_getcsr(); + if (argc == 2 && std::string(argv[1]) == "--bench") { + fesa::Bench(); + return 0; + } + const bool verify_only = argc == 2 && std::string(argv[1]) == "--verify"; + if (argc < 3 && !verify_only) + throw std::runtime_error( + "--verify OR --write|--compare PATH [--mutate|--corrupt-oracle]"); + std::string actual = + fesa::Verify(argc > 3 && std::string(argv[3]) == "--corrupt-oracle"); + actual += fesa::PublicSnapshots(); + if (argc > 3 && std::string(argv[3]) == "--mutate") + actual[actual.find('\n') - 1] ^= 1; + if (verify_only) { + // CTest and ordinary regression do not create baseline files. + } else if (std::string(argv[1]) == "--write") { + if (std::ifstream(argv[2], std::ios::binary)) + throw std::runtime_error("snapshot already exists"); + std::ofstream output(argv[2], std::ios::binary); + output << actual; + output.close(); + if (!output) throw std::runtime_error("snapshot write failure"); + } else if (std::string(argv[1]) == "--compare") { + std::ifstream input(argv[2], std::ios::binary); + const std::string expected((std::istreambuf_iterator(input)), {}); + if (expected != actual) + throw std::runtime_error("snapshot mismatch (bits/provenance/fenv)"); + } else { + throw std::runtime_error("unknown verification mode"); + } + std::cout + << argv[1] + << " passed; independent oracle and 42 private cases " + "(651 floating channels each), 64 public request/error cases\n"; + return 0; + } catch (const std::exception& error) { + std::cerr << error.what() << '\n'; + return 1; + } +} diff --git a/tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.md b/tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.md new file mode 100644 index 0000000..1b8d834 --- /dev/null +++ b/tests/unit/elements/nonlinear_mitc4_shell_jet_baseline.md @@ -0,0 +1,40 @@ +# Nonlinear MITC4 Jet 기준 검사 + +`fesa_jet_baseline_tests`는 제품의 private 구현을 별도 번역 단위에 포함하고 +일반 solver archive를 연결한다. 제품 source와 공개 API를 계측하지 않는다. +독립 다항식 지수 미분 및 미분 항등식과 원래 경로의 비트 snapshot은 별도 근거다. + +VS x64 환경과 저장소의 승인된 dependency 환경변수에서 실행한다. + +```bat +MSBuild msvc\fesa_jet_baseline_tests.vcxproj /m:2 /p:Configuration=Debug /p:Platform=x64 +.build\native\Debug\tests\fesa_jet_baseline_tests.exe --verify +cmake --preset windows-msvc-ninja-debug +cmake --build --preset windows-msvc-ninja-debug --target fesa_jet_baseline_tests --parallel 2 +ctest --preset windows-msvc-ninja-debug -R "^NonlinearMitc4Jet\.IndependentOracle$" --output-on-failure +``` + +Release는 Native의 `Configuration=Release` 또는 해당 CMake Release preset을 쓴다. +기존 unit/integration/reference 세 GTest 실행물의 발견 identity와 이 전용 검사의 +선택 수는 구분한다. CTest의 `--verify`는 baseline 파일을 생성하지 않는다. + +승인된 기준 revision에서만 새 경로에 `--write `을 한 번 실행한다. +기존 파일은 덮어쓰지 않는다. 이후에는 원래 파일과 hash를 보존하고 같은 빌드 경로와 +구성끼리 비교한다. 일반 회귀나 후보 코드에서 `--write`로 기대값을 만들지 않는다. + +```bat +.build\native\Debug\tests\fesa_jet_baseline_tests.exe --compare "" +.build\native\Release\tests\fesa_jet_baseline_tests.exe --compare "" +.build\cmake-debug\tests\fesa_jet_baseline_tests.exe --compare "" +``` + +`--compare --mutate`와 `--corrupt-oracle`는 각각 snapshot 비교와 독립 +oracle의 fixture 전용 오류 주입이며 exit 1이 기대값이다. `--bench`는 별도 micro +측정이다. 연산의 할당·해제와 교대 입력값 갱신을 포함하고 fixture 생성은 제외한다. +전체 채널을 읽는 noinline checksum과 volatile consume을 포함해 결과 제거를 막는다. +공식 모델 CLI의 시간과 혼합하지 않는다. + +snapshot은 finite 값·signed zero·Inf 비트, NaN 분류, provenance와 implicit Hessian, +복원한 초기 FP 환경 및 연산 뒤 FE/rounding/control word/full MXCSR를 기록한다. +공개 반환의 방향 요청 유무·실패 진단·회복량 identity도 포함한다. 수식의 정확성은 +독립 oracle로 검사하고, snapshot은 기준의 계산 순서와 관측 의미 보존을 검사한다.