From a73d3fc9ca4f95914fdc45169a56b9ef979bf4dd Mon Sep 17 00:00:00 2001 From: "KOKO\\Mimi" Date: Tue, 20 Oct 2026 09:45:58 +0900 Subject: [PATCH] Add MITC4 drilling diagnostics output and comparison contract --- msvc/Fesa.Common.props | 2 +- scripts/verification/compare_hht_dynamic.py | 37 +++- .../io/hdf5/dynamic_hdf5_results_writer.cpp | 38 +++- .../io/hdf5/eigen_hdf5_results_writer.cpp | 21 ++- src/fesa/io/hdf5/eigen_hdf5_self_check.cpp | 14 +- src/fesa/io/hdf5/eigen_hdf5_self_check.h | 5 +- src/fesa/io/hdf5/hdf5_results_writer.cpp | 25 ++- src/fesa/io/hdf5/hdf5_self_check.cpp | 50 ++++- src/fesa/io/hdf5/hdf5_self_check.h | 6 +- src/fesa/io/hdf5/mitc4_drilling_diagnostics.h | 172 ++++++++++++++++++ src/fesa/io/hdf5/mitc4_file_hooks.h | 42 +++++ .../io/hdf5/nonlinear_hdf5_results_writer.cpp | 27 ++- tests/CMakeLists.txt | 2 + .../integration/app/fesa_application_test.cpp | 38 ++++ .../reference/eigen_reference_comparison.cpp | 5 + tests/reference/eigen_reference_comparison.h | 2 + .../eigen_reference_comparison_test.cpp | 19 ++ tests/reference/mitc4_drilling_reader.h | 21 +++ tests/reference/mitc4_drilling_report.h | 58 ++++++ .../reference/mitc4_reference_comparison.cpp | 10 +- tests/reference/mitc4_reference_comparison.h | 2 + .../mitc4_reference_comparison_test.cpp | 55 ++++++ .../nonlinear_static_reference_cases_test.cpp | 74 ++++++++ .../nonlinear_static_reference_comparison.cpp | 25 ++- .../nonlinear_static_reference_comparison.h | 2 + tests/test_compare_hht_dynamic.py | 85 +++++++++ .../hdf5/dynamic_hdf5_results_writer_test.cpp | 61 +++++++ .../hdf5/eigen_hdf5_results_writer_test.cpp | 51 +++++- .../unit/io/hdf5/hdf5_results_writer_test.cpp | 159 +++++++++++++++- tests/unit/io/hdf5/mitc4_fault_test_helpers.h | 42 +++++ .../nonlinear_hdf5_results_writer_test.cpp | 70 ++++++- 31 files changed, 1175 insertions(+), 45 deletions(-) create mode 100644 src/fesa/io/hdf5/mitc4_drilling_diagnostics.h create mode 100644 src/fesa/io/hdf5/mitc4_file_hooks.h create mode 100644 tests/reference/mitc4_drilling_reader.h create mode 100644 tests/reference/mitc4_drilling_report.h create mode 100644 tests/unit/io/hdf5/mitc4_fault_test_helpers.h diff --git a/msvc/Fesa.Common.props b/msvc/Fesa.Common.props index 34bf8c1..f8a2f99 100644 --- a/msvc/Fesa.Common.props +++ b/msvc/Fesa.Common.props @@ -45,7 +45,7 @@ WIN32;_WINDOWS;%(PreprocessorDefinitions) NDEBUG;WIN32;_WINDOWS;%(PreprocessorDefinitions) $(FesaRoot)include;%(AdditionalIncludeDirectories) - $(FesaRoot)include;$(FesaRoot)src\fesa;%(AdditionalIncludeDirectories) + $(FesaRoot)include;$(FesaRoot)src\fesa;%(AdditionalIncludeDirectories) $(FesaExternalOptions) /MP2 %(AdditionalOptions) TurnOffAllWarnings diff --git a/scripts/verification/compare_hht_dynamic.py b/scripts/verification/compare_hht_dynamic.py index 074830a..98114f0 100644 --- a/scripts/verification/compare_hht_dynamic.py +++ b/scripts/verification/compare_hht_dynamic.py @@ -197,6 +197,38 @@ def dataset(file: h5py.File, path: str, shape: tuple[int, ...]) -> h5py.Dataset: return result +def drilling_diagnostics(file: h5py.File) -> dict[str, Any] | None: + """Validate additive diagnostics without introducing a balance verdict.""" + metadata = file["metadata"] + mode_present = "mitc4_drilling_formulation" in metadata.attrs + version_present = "drilling_diagnostics_version" in metadata.attrs + require(mode_present == version_present, "partial-drilling-markers", "Both drilling markers are required.") + if not mode_present: + return None + mode = string_attribute(metadata, "mitc4_drilling_formulation") + require(mode in ("direct-initial-director-spring", "legacy-corotated-twist"), + "unknown-drilling-formulation", mode) + numeric_attribute(metadata, "drilling_diagnostics_version", " tuple[list[np.ndarray], dict[str, Any]]: """Precheck every result node/frame/component, then project by source key.""" @@ -319,9 +351,10 @@ def read_results(path: Path, nodes: dict[int, tuple[float, ...]], if frame == 0: require((values == 0.0).all(), "nonzero-initial-displacement", "Frame zero must have q0=0.") projected.append(values[selection]) + drilling = drilling_diagnostics(file) return projected, {"source_input_identity": source_identity, "source_node_count": NODE_COUNT, "accepted_frame_count": FRAME_COUNT, "checked_result_scalars": NODE_COUNT * FRAME_COUNT * 6, - "projected_node_labels": list(NODE_LABELS)} + "projected_node_labels": list(NODE_LABELS), "drilling_diagnostics": drilling} def evaluate_family(increment: int, quantity: str, result: list[float], @@ -374,6 +407,7 @@ def evaluate_family(increment: int, quantity: str, result: list[float], def compare(input_path: Path, reference_path: Path, results_path: Path) -> dict[str, Any]: report: dict[str, Any] = { + "drilling_diagnostics": None, "schema": "hht-reference-comparison-v1", "case": CASE, "policy": dict(POLICY), "identities": {}, "precheck": {"passed": False}, "families": [], "blocking_failures": [], "warnings": [], "blocking_passed": False, "verdict": "failed", @@ -386,6 +420,7 @@ def compare(input_path: Path, reference_path: Path, results_path: Path) -> dict[ nodes, content_identity = input_inventory(input_path) reference = read_reference(reference_path) result, evidence = read_results(results_path, nodes, content_identity) + report["drilling_diagnostics"] = evidence.pop("drilling_diagnostics") families = [] for frame in range(FRAME_COUNT): for quantity, first in (("U", 0), ("UR", 3)): diff --git a/src/fesa/io/hdf5/dynamic_hdf5_results_writer.cpp b/src/fesa/io/hdf5/dynamic_hdf5_results_writer.cpp index db70de4..67bc522 100644 --- a/src/fesa/io/hdf5/dynamic_hdf5_results_writer.cpp +++ b/src/fesa/io/hdf5/dynamic_hdf5_results_writer.cpp @@ -18,6 +18,10 @@ #include "io/hdf5/hdf5_primitives.h" #include "io/hdf5/hdf5_raii.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +#include "io/hdf5/mitc4_file_hooks.h" + namespace fesa { namespace { @@ -849,13 +853,13 @@ void CheckTables(const hid_t file, const DynamicAnalysisState& state, Status SelfCheck(const std::filesystem::path& path, const Domain& domain, const DynamicAnalysisState& state, const hdf5_internal::WriterModelData& model, - const std::vector& diagnostics) { + const std::vector& diagnostics, herr_t (*close_reader)(hid_t)) { try { Hdf5Handle file{ hdf5_internal::RequireHdf5Id( H5Fopen(path.string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT), "Unable to reopen dynamic candidate."), - H5Fclose}; + close_reader}; Hdf5Handle metadata{hdf5_internal::RequireHdf5Id( H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), "The dynamic metadata group is missing."), @@ -902,9 +906,16 @@ Status SelfCheck(const std::filesystem::path& path, const Domain& domain, "Unable to inspect accepted frame group."); Check(information.nlinks == state.AcceptedFrames().size(), "Unexpected frame group inventory."); + const auto drilling = hdf5_internal::RequireDrillingDiagnostics( + file.Get(), state.AcceptedFrames().front().mitc4_drilling->mode); + hdf5_internal::RequireDrilling(drilling.frames.size() == state.AcceptedFrames().size()); + std::size_t drilling_index = 0U; for (const auto& frame : state.AcceptedFrames()) { + hdf5_internal::CheckDrillingFrame(drilling.frames[drilling_index++], frame.frame_index, + frame.global_equilibrium, *frame.mitc4_drilling); WriteFrame(file.Get(), domain, frame, true); } + hdf5_internal::RequireHdf5(file.CloseChecked(), "Unable to close dynamic readback."); return Status::Ok(); } catch (const std::exception& error) { return Failure("dynamic-output-schema", error.what()); @@ -913,10 +924,10 @@ Status SelfCheck(const std::filesystem::path& path, const Domain& domain, } // namespace -Status DynamicHdf5ResultsWriter::Write( +Status hdf5_internal::WriteDynamicFile( const std::filesystem::path& output_path, const Domain& domain, const DynamicAnalysisState& state, - const std::vector& diagnostics) { + const std::vector& diagnostics, const Mitc4FileHooks& hooks) { if (output_path.empty() || output_path.filename().empty() || !state.IsInitialized() || state.AcceptedFrames().empty() || domain.Dynamic() == nullptr) @@ -924,6 +935,11 @@ Status DynamicHdf5ResultsWriter::Write( "Dynamic output requires a named path, initialized state, " "and accepted frames."); try { + for (const auto& frame : state.AcceptedFrames()) { + hdf5_internal::ValidateDrillingEvidence(frame.global_equilibrium, frame.mitc4_drilling); + hdf5_internal::RequireDrilling(frame.mitc4_drilling->mode == + state.AcceptedFrames().front().mitc4_drilling->mode); + } hdf5_internal::Hdf5ErrorSilencer silence; const auto& step = *domain.Dynamic(); Check(state.AcceptedFrames().size() == @@ -965,24 +981,28 @@ Status DynamicHdf5ResultsWriter::Write( H5Fcreate(candidate_path.string().c_str(), H5F_ACC_EXCL, H5P_DEFAULT, H5P_DEFAULT), "Unable to create dynamic candidate."), - H5Fclose}; + hooks.candidate_close}; (void)hdf5_internal::CreateGroup(file.Get(), "/model/shell"); hdf5_internal::WriteModel(file.Get(), domain, model); WriteMetadata(file.Get(), domain, state.AcceptedFrames().front()); + hdf5_internal::WriteDrillingMetadata(file.Get(), state.AcceptedFrames().front().mitc4_drilling->mode); WriteDensity(file.Get(), domain); (void)hdf5_internal::CreateGroup(file.Get(), "/steps"); (void)hdf5_internal::CreateGroup(file.Get(), "/steps/Step-1"); WriteAcceptedFrames(file.Get(), state.AcceptedFrames()); for (const auto& frame : state.AcceptedFrames()) { WriteFrame(file.Get(), domain, frame); + hdf5_internal::WriteDrillingFrame(file.Get(), + "/steps/Step-1/frames/" + std::to_string(frame.frame_index), *frame.mitc4_drilling); } WriteDiagnostics(file.Get(), diagnostics); hdf5_internal::RequireHdf5(H5Fflush(file.Get(), H5F_SCOPE_GLOBAL), "Unable to flush dynamic candidate."); hdf5_internal::RequireHdf5(file.CloseChecked(), "Unable to close dynamic candidate."); + if (hooks.before_self_check) hooks.before_self_check(candidate_path); const Status check = - SelfCheck(candidate_path, domain, state, model, diagnostics); + SelfCheck(candidate_path, domain, state, model, diagnostics, hooks.reader_close); const Status finalization = hdf5_internal::AtomicFinalizeValidatedCandidate( candidate_path, output_path, check); if (finalization.IsOk()) cleanup.Release(); @@ -994,4 +1014,10 @@ Status DynamicHdf5ResultsWriter::Write( } } +Status DynamicHdf5ResultsWriter::Write( + const std::filesystem::path& output_path, const Domain& domain, + const DynamicAnalysisState& state, const std::vector& diagnostics) { + return hdf5_internal::WriteDynamicFile(output_path, domain, state, diagnostics); +} + } // namespace fesa diff --git a/src/fesa/io/hdf5/eigen_hdf5_results_writer.cpp b/src/fesa/io/hdf5/eigen_hdf5_results_writer.cpp index 60688de..dc6ee25 100644 --- a/src/fesa/io/hdf5/eigen_hdf5_results_writer.cpp +++ b/src/fesa/io/hdf5/eigen_hdf5_results_writer.cpp @@ -22,6 +22,10 @@ #include "io/hdf5/hdf5_raii.h" #include "math/finite_values.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +#include "io/hdf5/mitc4_file_hooks.h" + namespace fesa { namespace { @@ -875,13 +879,14 @@ void WriteDiagnostics(hid_t file, const std::vector& diagnostics) { } // namespace -Status EigenHdf5ResultsWriter::Write( +Status hdf5_internal::WriteEigenFile( const std::filesystem::path& output_path, const Domain& domain, const EigenAnalysisState& state, - const std::vector& diagnostics) { + const std::vector& diagnostics, const Mitc4FileHooks& hooks) { const Status valid = Validate(domain, state, diagnostics); if (!valid.IsOk()) return valid; try { + hdf5_internal::RequireDrilling(state.drilling_mode.has_value()); hdf5_internal::Hdf5ErrorSilencer silence; const auto candidate_path = hdf5_internal::MakeCandidatePath(output_path); hdf5_internal::CandidateFileGuard cleanup{candidate_path}; @@ -890,7 +895,7 @@ Status EigenHdf5ResultsWriter::Write( H5Fcreate(candidate_path.string().c_str(), H5F_ACC_EXCL, H5P_DEFAULT, H5P_DEFAULT), "Unable to create eigen candidate."), - H5Fclose}; + hooks.candidate_close}; hdf5_internal::WriterModelData model{ {}, state.constraint_mask, state.prescribed_displacement}; (void)hdf5_internal::CreateGroup(file.Get(), "/model"); @@ -899,6 +904,7 @@ Status EigenHdf5ResultsWriter::Write( WriteDensities(file.Get(), domain); WriteDiagnostics(file.Get(), diagnostics); WriteEigenMetadata(file.Get()); + hdf5_internal::WriteDrillingMetadata(file.Get(), *state.drilling_mode); WriteStep(file.Get(), domain, state); (void)hdf5_internal::CreateGroup(file.Get(), kEigenRoot); WriteModes(file.Get(), state, domain.Frequency()->source_step_name_raw); @@ -940,8 +946,9 @@ Status EigenHdf5ResultsWriter::Write( "Unable to flush eigen candidate."); hdf5_internal::RequireHdf5(file.CloseChecked(), "Unable to close eigen candidate."); + if (hooks.before_self_check) hooks.before_self_check(candidate_path); const Status checked = hdf5_internal::SelfCheckEigenFile( - candidate_path, domain, state, diagnostics); + candidate_path, domain, state, diagnostics, hooks.reader_close); const Status result = hdf5_internal::AtomicFinalizeValidatedCandidate( candidate_path, output_path, checked); if (result.IsOk()) cleanup.Release(); @@ -953,4 +960,10 @@ Status EigenHdf5ResultsWriter::Write( } } +Status EigenHdf5ResultsWriter::Write( + const std::filesystem::path& output_path, const Domain& domain, + const EigenAnalysisState& state, const std::vector& diagnostics) { + return hdf5_internal::WriteEigenFile(output_path, domain, state, diagnostics); +} + } // namespace fesa diff --git a/src/fesa/io/hdf5/eigen_hdf5_self_check.cpp b/src/fesa/io/hdf5/eigen_hdf5_self_check.cpp index b41907e..9a923a2 100644 --- a/src/fesa/io/hdf5/eigen_hdf5_self_check.cpp +++ b/src/fesa/io/hdf5/eigen_hdf5_self_check.cpp @@ -15,6 +15,8 @@ #include "io/hdf5/hdf5_primitives.h" #include "io/hdf5/hdf5_raii.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + namespace fesa::hdf5_internal { namespace { @@ -726,13 +728,21 @@ bool CheckModalRows(hid_t file, const Domain& domain, Status SelfCheckEigenFile(const std::filesystem::path& path, const Domain& domain, const EigenAnalysisState& state, - const std::vector& diagnostics) { + const std::vector& diagnostics, + herr_t (*close_reader)(hid_t)) { try { Hdf5ErrorSilencer silence; Hdf5Handle file{RequireHdf5Id(H5Fopen(path.string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT), "Unable to reopen eigen candidate."), - H5Fclose}; + close_reader}; + const auto drilling = ReadDrillingDiagnostics(file.Get()); + if (state.drilling_mode) { + RequireDrilling(drilling.has_value() && drilling->formulation == + Mitc4DrillingModeName(*state.drilling_mode) && drilling->frames.empty()); + } else { + RequireDrilling(!drilling.has_value()); + } for (const char* required : {"/metadata", "/model/nodes", "/model/elements", "/model/shell/densities", "/diagnostics", "/steps/Step-1/eigen/modes", diff --git a/src/fesa/io/hdf5/eigen_hdf5_self_check.h b/src/fesa/io/hdf5/eigen_hdf5_self_check.h index da7458e..2cbd068 100644 --- a/src/fesa/io/hdf5/eigen_hdf5_self_check.h +++ b/src/fesa/io/hdf5/eigen_hdf5_self_check.h @@ -1,6 +1,8 @@ #ifndef FESA_IO_HDF5_EIGEN_HDF5_SELF_CHECK_H_ #define FESA_IO_HDF5_EIGEN_HDF5_SELF_CHECK_H_ +#include + #include #include @@ -17,7 +19,8 @@ namespace fesa::hdf5_internal { /// @brief Reopens and validates a closed authoritative eigen HDF5 candidate. Status SelfCheckEigenFile(const std::filesystem::path& path, const Domain& domain, const EigenAnalysisState& state, - const std::vector& diagnostics); + const std::vector& diagnostics, + herr_t (*close_reader)(hid_t) = H5Fclose); } // namespace fesa::hdf5_internal diff --git a/src/fesa/io/hdf5/hdf5_results_writer.cpp b/src/fesa/io/hdf5/hdf5_results_writer.cpp index 077ca81..4640766 100644 --- a/src/fesa/io/hdf5/hdf5_results_writer.cpp +++ b/src/fesa/io/hdf5/hdf5_results_writer.cpp @@ -12,12 +12,16 @@ #include "io/hdf5/hdf5_result_writer.h" #include "io/hdf5/hdf5_self_check.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +#include "io/hdf5/mitc4_file_hooks.h" + namespace fesa { -Status Hdf5ResultsWriter::Write(const std::filesystem::path& output_path, +Status hdf5_internal::WriteStaticFile(const std::filesystem::path& output_path, const Domain& domain, const AnalysisState& state, - const std::vector& diagnostics) { + const std::vector& diagnostics, const Mitc4FileHooks& hooks) { auto plan_result = hdf5_internal::BuildWriterPlan(output_path, domain, state, diagnostics); if (!plan_result.HasValue()) { @@ -26,6 +30,8 @@ Status Hdf5ResultsWriter::Write(const std::filesystem::path& output_path, hdf5_internal::WriterPlan plan = std::move(plan_result.Value()); try { + if (hdf5_internal::IsShellDomain(domain)) + hdf5_internal::ValidateDrillingEvidence(state.Equilibrium(), state.Mitc4Drilling()); hdf5_internal::Hdf5ErrorSilencer silence_backend_errors; const std::filesystem::path candidate_path = hdf5_internal::MakeCandidatePath(output_path); @@ -35,7 +41,7 @@ Status Hdf5ResultsWriter::Write(const std::filesystem::path& output_path, H5Fcreate(candidate_path.string().c_str(), H5F_ACC_EXCL, H5P_DEFAULT, H5P_DEFAULT), "Unable to create the temporary HDF5 file."), - H5Fclose}; + hooks.candidate_close}; (void)hdf5_internal::CreateGroup(file.Get(), "/model"); (void)hdf5_internal::CreateGroup(file.Get(), "/steps/Step-1/frames/0/nodal"); @@ -50,14 +56,19 @@ Status Hdf5ResultsWriter::Write(const std::filesystem::path& output_path, } hdf5_internal::WriteModel(file.Get(), domain, plan.model_data); hdf5_internal::WriteResults(file.Get(), domain, state, diagnostics, plan); + if (hdf5_internal::IsShellDomain(domain)) { + hdf5_internal::WriteDrillingMetadata(file.Get(), state.Mitc4Drilling()->mode); + hdf5_internal::WriteDrillingFrame(file.Get(), hdf5_internal::kStepRoot, *state.Mitc4Drilling()); + } hdf5_internal::RequireHdf5(H5Fflush(file.Get(), H5F_SCOPE_GLOBAL), "Unable to flush the temporary HDF5 file."); hdf5_internal::RequireHdf5( file.CloseChecked(), "Unable to close the temporary HDF5 file after writing."); + if (hooks.before_self_check) hooks.before_self_check(candidate_path); const Status self_check = hdf5_internal::SelfCheck( - candidate_path, domain, state, diagnostics.size(), plan.model_data); + candidate_path, domain, state, diagnostics.size(), plan.model_data, hooks.reader_close); const Status atomic = hdf5_internal::AtomicFinalizeValidatedCandidate( candidate_path, output_path, self_check); if (atomic.IsOk()) { @@ -71,4 +82,10 @@ Status Hdf5ResultsWriter::Write(const std::filesystem::path& output_path, } } +Status Hdf5ResultsWriter::Write( + const std::filesystem::path& output_path, const Domain& domain, + const AnalysisState& state, const std::vector& diagnostics) { + return hdf5_internal::WriteStaticFile(output_path, domain, state, diagnostics); +} + } // namespace fesa diff --git a/src/fesa/io/hdf5/hdf5_self_check.cpp b/src/fesa/io/hdf5/hdf5_self_check.cpp index 2779079..95e8973 100644 --- a/src/fesa/io/hdf5/hdf5_self_check.cpp +++ b/src/fesa/io/hdf5/hdf5_self_check.cpp @@ -13,6 +13,8 @@ #include "io/hdf5/hdf5_raii.h" #include "io/hdf5/hdf5_primitives.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + namespace fesa::hdf5_internal { Hdf5Handle OpenDatasetForCheck(const hid_t file, const char* path) { @@ -516,18 +518,30 @@ void RequireSchemaV1ShellIdentity(const hid_t file, const char* path, void SelfCheckFile(const std::filesystem::path& path, const Domain& domain, const AnalysisState& state, const std::size_t diagnostic_count, - const WriterModelData& model_data) { + const WriterModelData& model_data, + herr_t (*close_reader)(hid_t)) { if (H5Fis_hdf5(path.string().c_str()) <= 0) { throw Hdf5Failure{"The temporary output is not an HDF5 file."}; } Hdf5Handle file{ RequireHdf5Id(H5Fopen(path.string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT), "Unable to reopen the temporary HDF5 file read-only."), - H5Fclose}; + close_reader}; Hdf5Handle metadata{ RequireHdf5Id(H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), "The HDF5 metadata group is missing."), H5Gclose}; + if (IsShellDomain(domain)) { + const auto drilling = ReadDrillingDiagnostics(file.Get()); + if (state.Mitc4Drilling()) { + ValidateDrillingEvidence(state.Equilibrium(), state.Mitc4Drilling()); + RequireDrilling(drilling.has_value() && drilling->formulation == + Mitc4DrillingModeName(state.Mitc4Drilling()->mode) && drilling->frames.size() == 1U); + CheckDrillingFrame(drilling->frames.front(), 0U, state.Equilibrium(), *state.Mitc4Drilling()); + } else { + RequireDrilling(!drilling.has_value()); + } + } RequireUint64Attribute(metadata.Get(), "schema_version", 0U); RequireStringAttribute(metadata.Get(), "solver_version", std::string{SolverVersion()}); @@ -751,18 +765,36 @@ void SelfCheckFile(const std::filesystem::path& path, const Domain& domain, void NonlinearSelfCheckFile(const std::filesystem::path& path, const Domain& domain, const NonlinearAnalysisState& state, - const std::size_t diagnostic_count) { + const std::size_t diagnostic_count, + herr_t (*close_reader)(hid_t)) { if (H5Fis_hdf5(path.string().c_str()) <= 0) { throw Hdf5Failure{"The nonlinear candidate is not an HDF5 file."}; } Hdf5Handle file{ RequireHdf5Id(H5Fopen(path.string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT), "Unable to reopen nonlinear candidate read-only."), - H5Fclose}; + close_reader}; Hdf5Handle metadata{ RequireHdf5Id(H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), "The nonlinear metadata group is missing."), H5Gclose}; + RequireDrilling(!state.AcceptedFrames().empty()); + const auto drilling = ReadDrillingDiagnostics(file.Get()); + const auto& expected_drilling = state.AcceptedFrames().front().mitc4_drilling; + if (expected_drilling) { + RequireDrilling(drilling.has_value() && drilling->formulation == + Mitc4DrillingModeName(expected_drilling->mode) && + drilling->frames.size() == state.AcceptedFrames().size()); + for (std::size_t i = 0; i < drilling->frames.size(); ++i) { + const auto& frame = state.AcceptedFrames()[i]; + ValidateDrillingEvidence(frame.global_equilibrium, frame.mitc4_drilling); + RequireDrilling(frame.mitc4_drilling->mode == expected_drilling->mode); + CheckDrillingFrame(drilling->frames[i], frame.frame_index, frame.global_equilibrium, *frame.mitc4_drilling); + } + } else { + RequireDrilling(!drilling.has_value()); + for (const auto& frame : state.AcceptedFrames()) RequireDrilling(!frame.mitc4_drilling); + } RequireUint64Attribute(metadata.Get(), "schema_version", 1U); RequireStringAttribute(metadata.Get(), "feature_id", "nonlinear-static-newton-raphson"); @@ -919,9 +951,10 @@ void NonlinearSelfCheckFile(const std::filesystem::path& path, Status SelfCheck(const std::filesystem::path& path, const Domain& domain, const AnalysisState& state, const std::size_t diagnostic_count, - const WriterModelData& model_data) { + const WriterModelData& model_data, + herr_t (*close_reader)(hid_t)) { try { - SelfCheckFile(path, domain, state, diagnostic_count, model_data); + SelfCheckFile(path, domain, state, diagnostic_count, model_data, close_reader); return Status::Ok(); } catch (const Hdf5Failure& failure) { return OutputFailure("hdf5-write-failure", failure.what()); @@ -933,9 +966,10 @@ Status SelfCheck(const std::filesystem::path& path, const Domain& domain, Status NonlinearSelfCheck(const std::filesystem::path& path, const Domain& domain, const NonlinearAnalysisState& state, - const std::size_t diagnostic_count) { + const std::size_t diagnostic_count, + herr_t (*close_reader)(hid_t)) { try { - NonlinearSelfCheckFile(path, domain, state, diagnostic_count); + NonlinearSelfCheckFile(path, domain, state, diagnostic_count, close_reader); return Status::Ok(); } catch (const Hdf5Failure& failure) { return OutputFailure("hdf5-write-failure", failure.what()); diff --git a/src/fesa/io/hdf5/hdf5_self_check.h b/src/fesa/io/hdf5/hdf5_self_check.h index be6fb94..7369eda 100644 --- a/src/fesa/io/hdf5/hdf5_self_check.h +++ b/src/fesa/io/hdf5/hdf5_self_check.h @@ -15,13 +15,15 @@ namespace fesa::hdf5_internal { /// @brief Reopens and validates the complete closed schema-v0 candidate. Status SelfCheck(const std::filesystem::path& path, const Domain& domain, const AnalysisState& state, std::size_t diagnostic_count, - const WriterModelData& model_data); + const WriterModelData& model_data, + herr_t (*close_reader)(hid_t) = H5Fclose); /// @brief Reopens and validates the complete closed schema-v1 candidate. Status NonlinearSelfCheck(const std::filesystem::path& path, const Domain& domain, const NonlinearAnalysisState& state, - std::size_t diagnostic_count); + std::size_t diagnostic_count, + herr_t (*close_reader)(hid_t) = H5Fclose); } // namespace fesa::hdf5_internal diff --git a/src/fesa/io/hdf5/mitc4_drilling_diagnostics.h b/src/fesa/io/hdf5/mitc4_drilling_diagnostics.h new file mode 100644 index 0000000..5a0f826 --- /dev/null +++ b/src/fesa/io/hdf5/mitc4_drilling_diagnostics.h @@ -0,0 +1,172 @@ +#ifndef FESA_IO_HDF5_MITC4_DRILLING_DIAGNOSTICS_H_ +#define FESA_IO_HDF5_MITC4_DRILLING_DIAGNOSTICS_H_ + +#include +#include +#include +#include +#include +#include + +#include "fesa/results/result_records.h" +#include "io/hdf5/hdf5_primitives.h" + +namespace fesa::hdf5_internal { + +struct StoredDrillingFrame { + std::string storage_step; + std::uint64_t frame_index; + std::array raw_equilibrium; + std::array drilling_moment; + std::array equilibrium_with_drilling; +}; + +struct StoredDrillingDiagnostics { + std::string formulation; + std::uint64_t version; + std::vector frames; +}; + +inline void RequireDrilling(bool valid) { + if (!valid) throw Hdf5Failure{"Invalid MITC4 drilling diagnostics contract."}; +} + +inline void ValidateDrillingEvidence( + const std::array& raw, + const std::optional& evidence) { + RequireDrilling(evidence.has_value()); + RequireDrilling(Mitc4DrillingModeName(evidence->mode)[0] != '\0' && + IsFinite(raw) && IsFinite(evidence->drilling_moment) && + IsFinite(evidence->equilibrium_with_drilling)); + for (std::size_t i = 0; i < 6U; ++i) { + RequireDrilling(evidence->equilibrium_with_drilling[i] == + (i < 3U ? raw[i] : raw[i] - evidence->drilling_moment[i - 3U])); + } +} + +inline void WriteDrillingMetadata(hid_t file, Mitc4DrillingMode mode) { + RequireDrilling(Mitc4DrillingModeName(mode)[0] != '\0'); + Hdf5Handle metadata{RequireHdf5Id(H5Gopen2(file, "/metadata", H5P_DEFAULT), + "Missing drilling metadata group."), H5Gclose}; + WriteStringAttribute(metadata.Get(), "mitc4_drilling_formulation", + Mitc4DrillingModeName(mode)); + WriteUint64Attribute(metadata.Get(), "drilling_diagnostics_version", 1U); +} + +inline void WriteDrillingFrame(hid_t file, const std::string& root, + const Mitc4DrillingEvidence& evidence) { + WriteModelDoubleDataset(file, root + "/global/drilling_moment", {3U}, + evidence.drilling_moment.data(), 3U, "MDRILL_X,MDRILL_Y,MDRILL_Z", + "force*length,force*length,force*length", "global-cartesian", "global-origin"); + WriteModelDoubleDataset(file, root + "/global/equilibrium_with_drilling", {6U}, + evidence.equilibrium_with_drilling.data(), 6U, + "FORCE_1,FORCE_2,FORCE_3,MOMENT_1,MOMENT_2,MOMENT_3", + "force,force,force,force*length,force*length,force*length", + "global-cartesian", "global-origin"); +} + +inline std::string ReadDrillingString(hid_t object, const char* name) { + Hdf5Handle attribute{RequireHdf5Id(H5Aopen(object, name, H5P_DEFAULT), + "Missing drilling string."), H5Aclose}; + Hdf5Handle type{H5Aget_type(attribute.Get()), H5Tclose}; + Hdf5Handle space{H5Aget_space(attribute.Get()), H5Sclose}; + RequireDrilling(H5Sget_simple_extent_type(space.Get()) == H5S_SCALAR && + H5Tget_class(type.Get()) == H5T_STRING && + H5Tis_variable_str(type.Get()) > 0 && + H5Tget_cset(type.Get()) == H5T_CSET_UTF8); + char* value = nullptr; + RequireHdf5(H5Aread(attribute.Get(), type.Get(), &value), "Cannot read drilling string."); + const std::string result = value == nullptr ? "" : value; + if (value != nullptr) H5free_memory(value); + RequireDrilling(IsValidUtf8(result)); + return result; +} + +template +inline std::array ReadDrillingValues(hid_t file, const std::string& path) { + Hdf5Handle dataset{RequireHdf5Id(H5Dopen2(file, path.c_str(), H5P_DEFAULT), + "Missing drilling dataset."), H5Dclose}; + Hdf5Handle type{H5Dget_type(dataset.Get()), H5Tclose}; + Hdf5Handle space{H5Dget_space(dataset.Get()), H5Sclose}; + hsize_t dimension = 0U; + RequireDrilling(H5Tequal(type.Get(), H5T_IEEE_F64LE) > 0 && + H5Sget_simple_extent_ndims(space.Get()) == 1 && + H5Sget_simple_extent_dims(space.Get(), &dimension, nullptr) >= 0 && + dimension == N); + std::array values{}; + RequireHdf5(H5Dread(dataset.Get(), H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, + H5P_DEFAULT, values.data()), "Cannot read drilling dataset."); + RequireDrilling(IsFinite(values)); + return values; +} + +// Markerless historical files retain their old contract; no mode is inferred. +inline std::optional ReadDrillingDiagnostics(hid_t file) { + Hdf5Handle metadata{RequireHdf5Id(H5Gopen2(file, "/metadata", H5P_DEFAULT), + "Missing drilling metadata group."), H5Gclose}; + const htri_t mode = H5Aexists(metadata.Get(), "mitc4_drilling_formulation"); + const htri_t version = H5Aexists(metadata.Get(), "drilling_diagnostics_version"); + RequireDrilling(mode >= 0 && version >= 0 && (mode > 0) == (version > 0)); + if (mode == 0) return std::nullopt; + StoredDrillingDiagnostics result{}; + result.formulation = ReadDrillingString(metadata.Get(), "mitc4_drilling_formulation"); + RequireDrilling(result.formulation == Mitc4DrillingModeName(Mitc4DrillingMode::kDirectSpring) || + result.formulation == Mitc4DrillingModeName(Mitc4DrillingMode::kLegacyCorotatedTwist)); + Hdf5Handle attribute{H5Aopen(metadata.Get(), "drilling_diagnostics_version", H5P_DEFAULT), H5Aclose}; + Hdf5Handle type{H5Aget_type(attribute.Get()), H5Tclose}; + Hdf5Handle space{H5Aget_space(attribute.Get()), H5Sclose}; + RequireDrilling(H5Tequal(type.Get(), H5T_STD_U64LE) > 0 && + H5Sget_simple_extent_type(space.Get()) == H5S_SCALAR); + RequireHdf5(H5Aread(attribute.Get(), H5T_NATIVE_UINT64, &result.version), "Cannot read drilling version."); + RequireDrilling(result.version == 1U); + if (H5Lexists(file, "/steps/Step-1/frames", H5P_DEFAULT) > 0) { + Hdf5Handle frames{H5Gopen2(file, "/steps/Step-1/frames", H5P_DEFAULT), H5Gclose}; + H5G_info_t info{}; + RequireHdf5(H5Gget_info(frames.Get(), &info), "Cannot inspect drilling frames."); + for (hsize_t i = 0; i < info.nlinks; ++i) { + const auto size = H5Lget_name_by_idx(frames.Get(), ".", H5_INDEX_NAME, H5_ITER_INC, + i, nullptr, 0U, H5P_DEFAULT); + RequireDrilling(size > 0); + std::vector name(static_cast(size) + 1U); + RequireDrilling(H5Lget_name_by_idx(frames.Get(), ".", H5_INDEX_NAME, H5_ITER_INC, + i, name.data(), name.size(), H5P_DEFAULT) == size); + const std::string token{name.data()}; + const std::string root = "/steps/Step-1/frames/" + token + "/global/"; + // Each new MITC4 static/dynamic frame has an existing raw equilibrium. + StoredDrillingFrame frame{}; + frame.storage_step = "Step-1"; + RequireDrilling(!token.empty() && token.find_first_not_of("0123456789") == std::string::npos); + frame.frame_index = std::stoull(token); + RequireDrilling(std::to_string(frame.frame_index) == token); + frame.raw_equilibrium = ReadDrillingValues<6>(file, root + "equilibrium"); + frame.drilling_moment = ReadDrillingValues<3>(file, root + "drilling_moment"); + frame.equilibrium_with_drilling = ReadDrillingValues<6>(file, root + "equilibrium_with_drilling"); + for (std::size_t j = 0; j < 6U; ++j) + RequireDrilling(frame.equilibrium_with_drilling[j] == (j < 3U ? frame.raw_equilibrium[j] : + frame.raw_equilibrium[j] - frame.drilling_moment[j - 3U])); + result.frames.push_back(frame); + } + std::sort(result.frames.begin(), result.frames.end(), + [](const auto& left, const auto& right) { return left.frame_index < right.frame_index; }); + } + return result; +} + +inline StoredDrillingDiagnostics RequireDrillingDiagnostics(hid_t file, Mitc4DrillingMode mode) { + const auto result = ReadDrillingDiagnostics(file); + RequireDrilling(result.has_value() && result->formulation == Mitc4DrillingModeName(mode)); + return *result; +} + +inline void CheckDrillingFrame(const StoredDrillingFrame& stored, + std::uint64_t index, + const std::array& raw, + const Mitc4DrillingEvidence& evidence) { + RequireDrilling(stored.frame_index == index && stored.raw_equilibrium == raw && + stored.drilling_moment == evidence.drilling_moment && + stored.equilibrium_with_drilling == evidence.equilibrium_with_drilling); +} + +} // namespace fesa::hdf5_internal + +#endif // FESA_IO_HDF5_MITC4_DRILLING_DIAGNOSTICS_H_ diff --git a/src/fesa/io/hdf5/mitc4_file_hooks.h b/src/fesa/io/hdf5/mitc4_file_hooks.h new file mode 100644 index 0000000..4e629b7 --- /dev/null +++ b/src/fesa/io/hdf5/mitc4_file_hooks.h @@ -0,0 +1,42 @@ +#ifndef FESA_IO_HDF5_MITC4_FILE_HOOKS_H_ +#define FESA_IO_HDF5_MITC4_FILE_HOOKS_H_ + +#include + +#include +#include +#include + +#include "fesa/core/status.h" + +namespace fesa { +class Domain; +class AnalysisState; +class NonlinearAnalysisState; +class DynamicAnalysisState; +struct EigenAnalysisState; +struct Diagnostic; + +namespace hdf5_internal { + +// Private lifecycle fault seams. Production uses real closers and no callback. +// Injected closers must release the handle even when reporting failure. +struct Mitc4FileHooks { + herr_t (*candidate_close)(hid_t) = H5Fclose; + herr_t (*reader_close)(hid_t) = H5Fclose; + std::function before_self_check; +}; + +Status WriteStaticFile(const std::filesystem::path&, const Domain&, + const AnalysisState&, const std::vector&, const Mitc4FileHooks& = {}); +Status WriteNonlinearFile(const std::filesystem::path&, const Domain&, + const NonlinearAnalysisState&, const std::vector&, const Mitc4FileHooks& = {}); +Status WriteDynamicFile(const std::filesystem::path&, const Domain&, + const DynamicAnalysisState&, const std::vector&, const Mitc4FileHooks& = {}); +Status WriteEigenFile(const std::filesystem::path&, const Domain&, + const EigenAnalysisState&, const std::vector&, const Mitc4FileHooks& = {}); + +} // namespace hdf5_internal +} // namespace fesa + +#endif // FESA_IO_HDF5_MITC4_FILE_HOOKS_H_ diff --git a/src/fesa/io/hdf5/nonlinear_hdf5_results_writer.cpp b/src/fesa/io/hdf5/nonlinear_hdf5_results_writer.cpp index cadd201..481dcb2 100644 --- a/src/fesa/io/hdf5/nonlinear_hdf5_results_writer.cpp +++ b/src/fesa/io/hdf5/nonlinear_hdf5_results_writer.cpp @@ -22,6 +22,10 @@ #include "math/finite_values.h" #include "results/shell_result_recovery.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +#include "io/hdf5/mitc4_file_hooks.h" + namespace fesa { namespace { @@ -955,10 +959,10 @@ Status ValidateDiagnostics(const std::vector& diagnostics) { } // namespace -Status NonlinearHdf5ResultsWriter::Write( +Status hdf5_internal::WriteNonlinearFile( const std::filesystem::path& output_path, const Domain& domain, const NonlinearAnalysisState& state, - const std::vector& diagnostics) { + const std::vector& diagnostics, const Mitc4FileHooks& hooks) { if (output_path.empty() || output_path.filename().empty()) { return Failure("invalid-output-path", "The HDF5 output path must name a file."); } @@ -981,6 +985,11 @@ Status NonlinearHdf5ResultsWriter::Write( } try { + for (const auto& frame : state.AcceptedFrames()) { + hdf5_internal::ValidateDrillingEvidence(frame.global_equilibrium, frame.mitc4_drilling); + hdf5_internal::RequireDrilling(frame.mitc4_drilling->mode == + state.AcceptedFrames().front().mitc4_drilling->mode); + } hdf5_internal::Hdf5ErrorSilencer silence_backend_errors; const std::filesystem::path candidate_path = hdf5_internal::MakeCandidatePath(output_path); @@ -989,25 +998,29 @@ Status NonlinearHdf5ResultsWriter::Write( H5Fcreate(candidate_path.string().c_str(), H5F_ACC_EXCL, H5P_DEFAULT, H5P_DEFAULT), "Unable to create nonlinear temporary HDF5 file."), - H5Fclose}; + hooks.candidate_close}; (void)hdf5_internal::CreateGroup(file.Get(), "/model"); (void)hdf5_internal::CreateGroup(file.Get(), "/model/shell"); (void)hdf5_internal::CreateGroup(file.Get(), "/steps/Step-1"); WriteMetadata(file.Get(), domain); + hdf5_internal::WriteDrillingMetadata(file.Get(), state.AcceptedFrames().front().mitc4_drilling->mode); WriteNodes(file.Get(), domain); WriteElements(file.Get(), domain); WriteModelData(file.Get(), domain, model_data); WriteAcceptedFrameTable(file.Get(), state.AcceptedFrames()); for (const auto& frame : state.AcceptedFrames()) { WriteFrame(file.Get(), domain, frame); + hdf5_internal::WriteDrillingFrame(file.Get(), + "/steps/Step-1/frames/" + std::to_string(frame.frame_index), *frame.mitc4_drilling); } WriteDiagnostics(file.Get(), diagnostics); hdf5_internal::RequireHdf5(H5Fflush(file.Get(), H5F_SCOPE_GLOBAL), "Unable to flush nonlinear temporary HDF5 file."); hdf5_internal::RequireHdf5( file.CloseChecked(), "Unable to close nonlinear temporary HDF5 file."); + if (hooks.before_self_check) hooks.before_self_check(candidate_path); const Status self_check = hdf5_internal::NonlinearSelfCheck( - candidate_path, domain, state, diagnostics.size()); + candidate_path, domain, state, diagnostics.size(), hooks.reader_close); const Status atomic = hdf5_internal::AtomicFinalizeValidatedCandidate( candidate_path, output_path, self_check); if (atomic.IsOk()) { @@ -1021,4 +1034,10 @@ Status NonlinearHdf5ResultsWriter::Write( } } +Status NonlinearHdf5ResultsWriter::Write( + const std::filesystem::path& output_path, const Domain& domain, + const NonlinearAnalysisState& state, const std::vector& diagnostics) { + return hdf5_internal::WriteNonlinearFile(output_path, domain, state, diagnostics); +} + } // namespace fesa diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 9c4a638..1c81290 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -162,6 +162,8 @@ add_executable( reference/thermal_reference_cases_test.cpp ) +target_include_directories(fesa_reference_tests PRIVATE ${PROJECT_SOURCE_DIR}/src/fesa) + target_link_libraries( fesa_reference_tests PRIVATE diff --git a/tests/integration/app/fesa_application_test.cpp b/tests/integration/app/fesa_application_test.cpp index efdbb4f..04942be 100644 --- a/tests/integration/app/fesa_application_test.cpp +++ b/tests/integration/app/fesa_application_test.cpp @@ -17,6 +17,8 @@ #include #include +#include "io/hdf5/mitc4_drilling_diagnostics.h" + namespace { constexpr const char* kStepRoot = "/steps/Step-1/frames/0"; @@ -525,4 +527,40 @@ TEST(Mitc4ShellCli, WritesShellHdf5ThroughTheExistingApplicationRoute) { const auto result = RunApplication({input.string()}); ASSERT_EQ(result.exit_code, 0) << result.standard_error; ExpectShellHdf5Identity(output, input); + const auto file = OpenFile(output); + const auto drilling = fesa::hdf5_internal::ReadDrillingDiagnostics(file.Get()); + ASSERT_TRUE(drilling); + EXPECT_EQ(drilling->formulation, "direct-initial-director-spring"); + EXPECT_EQ(drilling->version, 1U); + ASSERT_EQ(drilling->frames.size(), 1U); + EXPECT_EQ(drilling->frames.front().frame_index, 0U); +} + + +TEST(Mitc4ShellCli, DrillingActualHhtPipelineWritesInitialAndAcceptedDiagnostics) { + TempDirectory directory{"shell-drilling-hht"}; + const auto input = directory.Path() / "shell.inp"; + const auto output = directory.Path() / "shell_results.h5"; + std::string deck = AllConstrainedShellDeck(); + const auto material = deck.find("*Material, name=Steel\n"); + ASSERT_NE(material, std::string::npos); + deck.insert(material + std::string{"*Material, name=Steel\n"}.size(), "*Density\n2.\n"); + const auto step = deck.find("*Step,"); + ASSERT_NE(step, std::string::npos); + deck.replace(step, std::string::npos, + "*Step, name=Load, nlgeom=YES\n*Dynamic, Direct\n0.25, 0.5\n*Cload\nAll, 3, 1.\n*End Step\n"); + WriteText(input, deck); + const auto result = RunApplication({input.string()}); + ASSERT_EQ(result.exit_code, 0) << result.standard_error; + const auto file = OpenFile(output); + const auto drilling = fesa::hdf5_internal::ReadDrillingDiagnostics(file.Get()); + ASSERT_TRUE(drilling); + EXPECT_EQ(drilling->formulation, "direct-initial-director-spring"); + ASSERT_EQ(drilling->frames.size(), 3U); + for (std::size_t i = 0; i < 3U; ++i) { + EXPECT_EQ(drilling->frames[i].frame_index, i); + EXPECT_EQ(drilling->frames[i].raw_equilibrium, (std::array{})); + EXPECT_EQ(drilling->frames[i].drilling_moment, (std::array{})); + EXPECT_EQ(drilling->frames[i].equilibrium_with_drilling, (std::array{})); + } } diff --git a/tests/reference/eigen_reference_comparison.cpp b/tests/reference/eigen_reference_comparison.cpp index 6a80dd5..d0920dd 100644 --- a/tests/reference/eigen_reference_comparison.cpp +++ b/tests/reference/eigen_reference_comparison.cpp @@ -20,6 +20,8 @@ #include "fesa/io/abaqus/domain_mapper.h" #include "fesa/io/abaqus/input_reader.h" +#include "mitc4_drilling_reader.h" + namespace fesa::test { namespace { @@ -385,6 +387,9 @@ Result EigenReferenceComparison::Compare( } EigenReferenceComparisonReport report{}; + report.drilling_diagnostics = ReadMitc4DrillingReport(file.Get()); + if (report.drilling_diagnostics && !report.drilling_diagnostics->frames.empty()) + Fail("schema-mismatch", "Eigen drilling diagnostics must have no equilibrium frames."); report.case_id = reference_case.case_id; report.source_step_name = frequency.source_step_name_raw; report.requested_mode_count = requested; diff --git a/tests/reference/eigen_reference_comparison.h b/tests/reference/eigen_reference_comparison.h index cb6c002..35e53b4 100644 --- a/tests/reference/eigen_reference_comparison.h +++ b/tests/reference/eigen_reference_comparison.h @@ -7,6 +7,7 @@ #include #include "fesa/core/status.h" +#include "mitc4_drilling_report.h" #include "reference_tolerance_policy.h" namespace fesa::test { @@ -36,6 +37,7 @@ struct EigenReferenceComparisonReport { bool passed; std::string policy_id; double relative_rms_tolerance; + std::optional drilling_diagnostics; }; /// @brief Compares modal eigenvalues without changing the approved input or diff --git a/tests/reference/eigen_reference_comparison_test.cpp b/tests/reference/eigen_reference_comparison_test.cpp index 8b3dc6b..3e0a5ff 100644 --- a/tests/reference/eigen_reference_comparison_test.cpp +++ b/tests/reference/eigen_reference_comparison_test.cpp @@ -123,6 +123,25 @@ TEST(EigenReferenceComparison, EXPECT_EQ(ok.Value().policy_id, "mitc4-eigen-first-implementation"); EXPECT_DOUBLE_EQ(ok.Value().relative_rms_tolerance, 0.05); EXPECT_EQ(ok.Value().families.size(), 2U); + EXPECT_FALSE(ok.Value().drilling_diagnostics); + for (const std::string mode : {"direct-initial-director-spring", "legacy-corotated-twist"}) { + f = H5Fopen(h5.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT); + m = H5Gopen2(f, "/metadata", H5P_DEFAULT); + if (H5Aexists(m, "mitc4_drilling_formulation") > 0) H5Adelete(m, "mitc4_drilling_formulation"); + if (H5Aexists(m, "drilling_diagnostics_version") > 0) H5Adelete(m, "drilling_diagnostics_version"); + Str(m, "mitc4_drilling_formulation", mode); + U(m, "drilling_diagnostics_version", 1U); + H5Gclose(m); + H5Fclose(f); + const auto compared = fesa::test::EigenReferenceComparison::Compare(reference_case); + ASSERT_TRUE(compared.HasValue()); + ASSERT_TRUE(compared.Value().drilling_diagnostics); + EXPECT_EQ(compared.Value().drilling_diagnostics->formulation, mode); + EXPECT_EQ(compared.Value().drilling_diagnostics->version, 1U); + EXPECT_TRUE(compared.Value().drilling_diagnostics->frames.empty()); + EXPECT_EQ(compared.Value().families.size(), 2U); + } + // Only synthetic FESA output is changed; each TXT mode remains an independent // scalar family, including at both inclusive five-percent boundaries. for (const auto& sample : diff --git a/tests/reference/mitc4_drilling_reader.h b/tests/reference/mitc4_drilling_reader.h new file mode 100644 index 0000000..878f9ac --- /dev/null +++ b/tests/reference/mitc4_drilling_reader.h @@ -0,0 +1,21 @@ +#ifndef FESA_TESTS_REFERENCE_MITC4_DRILLING_READER_H_ +#define FESA_TESTS_REFERENCE_MITC4_DRILLING_READER_H_ + +#include "mitc4_drilling_report.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +namespace fesa::test { + +inline std::optional ReadMitc4DrillingReport(hid_t file) { + const auto stored = hdf5_internal::ReadDrillingDiagnostics(file); + if (!stored) return std::nullopt; + Mitc4DrillingDiagnosticsReport result{stored->formulation, stored->version, {}}; + for (const auto& frame : stored->frames) + result.frames.push_back({frame.storage_step, frame.frame_index, frame.raw_equilibrium, + frame.drilling_moment, frame.equilibrium_with_drilling}); + return result; +} + +} // namespace fesa::test + +#endif // FESA_TESTS_REFERENCE_MITC4_DRILLING_READER_H_ diff --git a/tests/reference/mitc4_drilling_report.h b/tests/reference/mitc4_drilling_report.h new file mode 100644 index 0000000..660a518 --- /dev/null +++ b/tests/reference/mitc4_drilling_report.h @@ -0,0 +1,58 @@ +#ifndef FESA_TESTS_REFERENCE_MITC4_DRILLING_REPORT_H_ +#define FESA_TESTS_REFERENCE_MITC4_DRILLING_REPORT_H_ + +#include +#include +#include +#include +#include +#include +#include + +namespace fesa::test { + +struct Mitc4DrillingFrameReport { + std::string storage_step; + std::uint64_t frame_index; + std::array raw_equilibrium; + std::array drilling_moment; + std::array equilibrium_with_drilling; +}; + +struct Mitc4DrillingDiagnosticsReport { + std::string formulation; + std::uint64_t version; + std::vector frames; +}; + +inline void WriteDrillingDiagnosticsJson( + std::ostream& stream, const std::optional& report) { + if (!report) { stream << "null"; return; } + stream << "{\"formulation\":\"" << report->formulation + << "\",\"version\":" << report->version << ",\"frames\":["; + const auto numbers = [&](const auto& values) { + stream << '['; + for (std::size_t i = 0; i < values.size(); ++i) { + if (i != 0U) stream << ','; + stream << values[i]; + } + stream << ']'; + }; + for (std::size_t i = 0; i < report->frames.size(); ++i) { + if (i != 0U) stream << ','; + const auto& frame = report->frames[i]; + stream << "{\"storage_step\":\"Step-1\",\"frame_index\":" << frame.frame_index + << ",\"raw_equilibrium\":"; + numbers(frame.raw_equilibrium); + stream << ",\"drilling_moment\":"; + numbers(frame.drilling_moment); + stream << ",\"equilibrium_with_drilling\":"; + numbers(frame.equilibrium_with_drilling); + stream << '}'; + } + stream << "]}"; +} + +} // namespace fesa::test + +#endif // FESA_TESTS_REFERENCE_MITC4_DRILLING_REPORT_H_ diff --git a/tests/reference/mitc4_reference_comparison.cpp b/tests/reference/mitc4_reference_comparison.cpp index bc8100c..f90181a 100644 --- a/tests/reference/mitc4_reference_comparison.cpp +++ b/tests/reference/mitc4_reference_comparison.cpp @@ -26,6 +26,8 @@ #include "reference_tolerance_policy.h" +#include "mitc4_drilling_reader.h" + namespace fesa::test { namespace { @@ -573,6 +575,7 @@ std::vector ReadDisplacement(const hid_t file, } struct Hdf5Projection { + std::optional drilling_diagnostics; std::vector rows; std::string source_element_type; std::string internal_formulation; @@ -619,7 +622,7 @@ Hdf5Projection ReadHdf5(const Mitc4ReferenceCase& reference_case) { values.begin() + static_cast(row * kComponents.size()), kComponents.size(), rows[row].values.begin()); } - return {std::move(rows), reference_case.expected_source_element_type, + return {ReadMitc4DrillingReport(file.Get()), std::move(rows), reference_case.expected_source_element_type, kInternalFormulation, kIntegrationRule}; } @@ -739,6 +742,7 @@ Result Mitc4ReferenceComparison::Compare( } Mitc4ComparisonReport report{}; + report.drilling_diagnostics = std::move(hdf5.drilling_diagnostics); report.case_id = reference_case.case_id; report.source_element_type = std::move(hdf5.source_element_type); report.internal_formulation = std::move(hdf5.internal_formulation); @@ -889,7 +893,9 @@ Status Mitc4ReferenceComparison::WriteDeterministicJson( } stream.imbue(std::locale::classic()); stream << std::setprecision(std::numeric_limits::max_digits10); - stream << "{\"case_id\":\"" << JsonEscape(report.case_id) + stream << "{\"drilling_diagnostics\":"; + WriteDrillingDiagnosticsJson(stream, report.drilling_diagnostics); + stream << ",\"case_id\":\"" << JsonEscape(report.case_id) << "\",\"source_element_type\":\"" << JsonEscape(report.source_element_type) << "\",\"internal_formulation\":\"" diff --git a/tests/reference/mitc4_reference_comparison.h b/tests/reference/mitc4_reference_comparison.h index 2fb07b8..e3753e8 100644 --- a/tests/reference/mitc4_reference_comparison.h +++ b/tests/reference/mitc4_reference_comparison.h @@ -8,6 +8,7 @@ #include #include "fesa/core/status.h" +#include "mitc4_drilling_report.h" namespace fesa::test { @@ -77,6 +78,7 @@ struct Mitc4ComparisonReport { std::vector warnings; std::size_t worst_row; bool passed; + std::optional drilling_diagnostics; }; class Mitc4ReferenceComparison { diff --git a/tests/reference/mitc4_reference_comparison_test.cpp b/tests/reference/mitc4_reference_comparison_test.cpp index a65397f..0a35450 100644 --- a/tests/reference/mitc4_reference_comparison_test.cpp +++ b/tests/reference/mitc4_reference_comparison_test.cpp @@ -22,6 +22,8 @@ #error FESA_TEST_BINARY_DIR must identify the CMake binary root. #endif +#include "io/hdf5/mitc4_drilling_diagnostics.h" + namespace { constexpr const char* kInstanceName = "PART-1-1"; @@ -771,3 +773,56 @@ TEST(Mitc4ReferenceComparison, RequiresOnlyDeclaredInputCsvAndHdf5) { } } // namespace + + +TEST(Mitc4ReferenceComparison, DrillingOldMarkerlessReportIsNull) { + ContractFixture fixture{"drilling-old"}; + const auto compared = fesa::test::Mitc4ReferenceComparison::Compare(fixture.ReferenceCase()); + ASSERT_TRUE(compared.HasValue()); + EXPECT_FALSE(compared.Value().drilling_diagnostics); + const auto json = fixture.Root() / "drilling-old.json"; + ASSERT_TRUE(fesa::test::Mitc4ReferenceComparison::WriteDeterministicJson(compared.Value(), json).IsOk()); + EXPECT_NE(ReadBytes(json).find("\"drilling_diagnostics\":null"), std::string::npos); +} + +TEST(Mitc4ReferenceComparison, DrillingNewMarkersPreserveAllFramesWithoutNewVerdict) { + for (const auto mode : {fesa::Mitc4DrillingMode::kDirectSpring, + fesa::Mitc4DrillingMode::kLegacyCorotatedTwist}) { + ContractFixture fixture{"drilling-new"}; + { + Hdf5Handle file{H5Fopen(fixture.Results().string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + fesa::hdf5_internal::WriteDrillingMetadata(file.Get(), mode); + // Nonzero raw balance is diagnostic only for the linear displacement comparator. + for (const int index : {0, 2, 10}) { + const std::string root = "/steps/Step-1/frames/" + std::to_string(index); + (void)fesa::hdf5_internal::CreateGroup(file.Get(), (root + "/global").c_str()); + const std::array raw{1., 2., 3., 4., 5., 6.}; + const fesa::Mitc4DrillingEvidence evidence{mode, {4., 5., 6.}, {1., 2., 3., 0., 0., 0.}}; + fesa::hdf5_internal::WriteModelDoubleDataset(file.Get(), root + "/global/equilibrium", {6U}, + raw.data(), 6U, "", "", "", ""); + fesa::hdf5_internal::WriteDrillingFrame(file.Get(), root, evidence); + } + } + const auto compared = fesa::test::Mitc4ReferenceComparison::Compare(fixture.ReferenceCase()); + ASSERT_TRUE(compared.HasValue()); + EXPECT_TRUE(compared.Value().passed); + EXPECT_EQ(compared.Value().metrics.size(), 2U); + const auto& diagnostics = compared.Value().drilling_diagnostics; + ASSERT_TRUE(diagnostics); + EXPECT_EQ(diagnostics->formulation, fesa::Mitc4DrillingModeName(mode)); + ASSERT_EQ(diagnostics->frames.size(), 3U); + EXPECT_EQ(diagnostics->frames[0].frame_index, 0U); + EXPECT_EQ(diagnostics->frames[1].frame_index, 2U); + EXPECT_EQ(diagnostics->frames[2].frame_index, 10U); + const auto json = fixture.Root() / "drilling-new.json"; + ASSERT_TRUE(fesa::test::Mitc4ReferenceComparison::WriteDeterministicJson(compared.Value(), json).IsOk()); + const auto bytes = ReadBytes(json); + EXPECT_NE(bytes.find("\"drilling_diagnostics\":{\"formulation\":\"" + + diagnostics->formulation + "\",\"version\":1,\"frames\":["), std::string::npos); + { + Hdf5Handle file{H5Fopen(fixture.Results().string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + ASSERT_GE(H5Ldelete(file.Get(), "/steps/Step-1/frames/2/global/drilling_moment", H5P_DEFAULT), 0); + } + EXPECT_FALSE(fesa::test::Mitc4ReferenceComparison::Compare(fixture.ReferenceCase()).HasValue()); + } +} diff --git a/tests/reference/nonlinear_static_reference_cases_test.cpp b/tests/reference/nonlinear_static_reference_cases_test.cpp index 0bf9c0d..04e7668 100644 --- a/tests/reference/nonlinear_static_reference_cases_test.cpp +++ b/tests/reference/nonlinear_static_reference_cases_test.cpp @@ -21,6 +21,8 @@ #error FESA_TEST_BINARY_DIR must identify the CMake binary root. #endif +#include "io/hdf5/mitc4_drilling_diagnostics.h" + namespace { std::string ReadBytes(const std::filesystem::path& path) { @@ -409,3 +411,75 @@ TEST(NonlinearStaticReference, EXPECT_FALSE(nonzero.normalized_metric.has_value()) << identity; } } + + +TEST(NonlinearStaticReference, DrillingSyntheticMomentSelectsOnlyExistingInternalFamily) { + const auto reference_case = MakeCase("shellnl-s4-final"); + const auto copied = fesa::test::CopyInputToExecutionDirectory(reference_case); + ASSERT_TRUE(copied.HasValue()); + const auto generated = fesa::test::RunFesaApplicationOnCopiedInput(reference_case, copied.Value()); + ASSERT_TRUE(generated.HasValue()); + const auto& path = generated.Value().generated_results_hdf5; + const auto baseline = fesa::test::CompareHdf5ToReferenceCsv(reference_case, path); + ASSERT_TRUE(baseline.HasValue()); + ASSERT_TRUE(baseline.Value().drilling_diagnostics); + const auto final_index = baseline.Value().drilling_diagnostics->frames.back().frame_index; + const std::string root = "/steps/Step-1/frames/" + std::to_string(final_index) + "/global/"; + const auto moment_gate = [](const auto& report) -> const auto& { + const auto found = std::find_if(report.internal_equilibrium.begin(), report.internal_equilibrium.end(), + [](const auto& family) { return family.family_identity == "coarse.global_moment_balance"; }); + if (found == report.internal_equilibrium.end()) throw std::runtime_error{"Missing existing moment family."}; + return *found; + }; + // This mutates only a generated working output into an algebra fixture. + // Its comparator verdict is not physics reference evidence. + for (const std::string mode : {"direct-initial-director-spring", "legacy-corotated-twist", "markerless"}) { + { + fesa::hdf5_internal::Hdf5Handle file{H5Fopen(path.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + fesa::hdf5_internal::Hdf5Handle metadata{H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), H5Gclose}; + ASSERT_GE(H5Adelete(metadata.Get(), "mitc4_drilling_formulation"), 0); + ASSERT_GE(H5Adelete(metadata.Get(), "drilling_diagnostics_version"), 0); + if (mode != "markerless") { + fesa::hdf5_internal::WriteStringAttribute(metadata.Get(), "mitc4_drilling_formulation", mode); + fesa::hdf5_internal::WriteUint64Attribute(metadata.Get(), "drilling_diagnostics_version", 1U); + } + for (const auto& entry : std::vector>>{ + {"equilibrium", {0.,0.,0.,1e6,0.,0.}}, {"drilling_moment", {1e6,0.,0.}}, + {"equilibrium_with_drilling", {0.,0.,0.,0.,0.,0.}}}) { + fesa::hdf5_internal::Hdf5Handle dataset{H5Dopen2(file.Get(), (root + entry.first).c_str(), H5P_DEFAULT), H5Dclose}; + ASSERT_GE(H5Dwrite(dataset.Get(), H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, + H5P_DEFAULT, entry.second.data()), 0); + } + } + const auto compared = fesa::test::CompareHdf5ToReferenceCsv(reference_case, path); + ASSERT_TRUE(compared.HasValue()); + EXPECT_TRUE(compared.Value().external_blocking_passed); + EXPECT_EQ(compared.Value().families.size(), baseline.Value().families.size()); + EXPECT_EQ(compared.Value().internal_equilibrium.size(), baseline.Value().internal_equilibrium.size()); + const auto& gate = moment_gate(compared.Value()); + EXPECT_EQ(gate.scale, moment_gate(baseline.Value()).scale); + EXPECT_EQ(gate.numerator, mode == "direct-initial-director-spring" ? 0.0 : 1e6); + EXPECT_EQ(gate.passed, mode == "direct-initial-director-spring"); + EXPECT_EQ(compared.Value().drilling_diagnostics.has_value(), mode != "markerless"); + } + // Restore new metadata, retain exact raw-minus-M algebra, and expose a real imbalance. + { + fesa::hdf5_internal::Hdf5Handle file{H5Fopen(path.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + fesa::hdf5_internal::WriteDrillingMetadata(file.Get(), fesa::Mitc4DrillingMode::kDirectSpring); + for (const auto& entry : std::vector>>{ + {"equilibrium", {0.,0.,0.,2e6,0.,0.}}, {"equilibrium_with_drilling", {0.,0.,0.,1e6,0.,0.}}}) { + fesa::hdf5_internal::Hdf5Handle dataset{H5Dopen2(file.Get(), (root + entry.first).c_str(), H5P_DEFAULT), H5Dclose}; + ASSERT_GE(H5Dwrite(dataset.Get(), H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, + H5P_DEFAULT, entry.second.data()), 0); + } + } + const auto imbalanced = fesa::test::CompareHdf5ToReferenceCsv(reference_case, path); + ASSERT_TRUE(imbalanced.HasValue()); + EXPECT_FALSE(moment_gate(imbalanced.Value()).passed); + EXPECT_EQ(moment_gate(imbalanced.Value()).scale, moment_gate(baseline.Value()).scale); + { + fesa::hdf5_internal::Hdf5Handle file{H5Fopen(path.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + ASSERT_GE(H5Ldelete(file.Get(), "/steps/Step-1/frames/0/global/drilling_moment", H5P_DEFAULT), 0); + } + EXPECT_FALSE(fesa::test::CompareHdf5ToReferenceCsv(reference_case, path).HasValue()); +} diff --git a/tests/reference/nonlinear_static_reference_comparison.cpp b/tests/reference/nonlinear_static_reference_comparison.cpp index 874e701..b55dae6 100644 --- a/tests/reference/nonlinear_static_reference_comparison.cpp +++ b/tests/reference/nonlinear_static_reference_comparison.cpp @@ -37,6 +37,8 @@ #error FESA_TEST_BINARY_DIR must identify the CMake binary root. #endif +#include "mitc4_drilling_reader.h" + namespace fesa::test { namespace { @@ -896,6 +898,8 @@ std::vector ReadConstraintMask(const hid_t file, } struct Hdf5Projection { + std::uint64_t selected_frame_index; + std::optional drilling_diagnostics; std::filesystem::path copied_input_path; std::vector nodes; std::vector displacement; @@ -972,7 +976,7 @@ Hdf5Projection ReadNonlinearHdf5( const std::uint64_t frame_index = ReadFactorOneFrame(file.Get()); const std::string root = "/steps/Step-1/frames/" + std::to_string(frame_index) + "/nodal/"; - return {copied_input, + return {frame_index, ReadMitc4DrillingReport(file.Get()), copied_input, nodes, ReadNodalDataset(file.Get(), root + "displacement", "UX,UY,UZ,URX,URY,URZ", contract.node_count), @@ -1320,6 +1324,16 @@ void AppendInternalEquilibrium(NonlinearReferenceComparisonReport* report, applied_moment_scale += applied_resultant.Norm(); reaction_moment_scale += constrained_resultant.Norm(); } + const std::array* moment_balance = &projection.global_equilibrium; + if (projection.drilling_diagnostics && + projection.drilling_diagnostics->formulation == "direct-initial-director-spring") { + const auto& frames = projection.drilling_diagnostics->frames; + const auto selected = std::find_if(frames.begin(), frames.end(), [&](const auto& frame) { + return frame.frame_index == projection.selected_frame_index; + }); + if (selected == frames.end()) Fail("schema-mismatch", "Missing factor-one drilling diagnostics."); + moment_balance = &selected->equilibrium_with_drilling; + } report->constrained_reaction_resultant = { reaction_force.Components(), reaction_moment.Components()}; const std::string prefix = ContractFor(reference_case).family_prefix; @@ -1329,9 +1343,9 @@ void AppendInternalEquilibrium(NonlinearReferenceComparisonReport* report, projection.global_equilibrium[2U]}, applied_force_scale, reaction_force_scale); AppendGlobalBalance(report, prefix + "global_moment_balance", - {projection.global_equilibrium[3U], - projection.global_equilibrium[4U], - projection.global_equilibrium[5U]}, + {(*moment_balance)[3U], + (*moment_balance)[4U], + (*moment_balance)[5U]}, applied_moment_scale, reaction_moment_scale); AppendFreeResidual(report, prefix + "accepted_free_force_residual", projection, false); @@ -1385,6 +1399,7 @@ void CompareLinearCase(NonlinearReferenceComparisonReport& report, : diagnostics.front().message; Fail("schema-mismatch", message); } + report.drilling_diagnostics = legacy.Value().drilling_diagnostics; for (const auto& metric : legacy.Value().metrics) { ReferenceToleranceFamilyDecision tolerance{}; tolerance.identity = metric.blocking ? "linear.translation" : "linear.rotation"; @@ -1418,6 +1433,7 @@ void CompareNonlinearCase(NonlinearReferenceComparisonReport& report, report.generated_results_hdf5); report.copied_input_path = hdf5.copied_input_path; report.generated_log_path = DerivedLogPath(hdf5.copied_input_path); + report.drilling_diagnostics = hdf5.drilling_diagnostics; RequireExactNodeLabels(hdf5.nodes, contract.node_count); const auto displacement_rows = ReadCsvWideRows( reference_case.displacement_csv_path, true, kDisplacementComponents, @@ -1475,6 +1491,7 @@ void ComparePressureCase(NonlinearReferenceComparisonReport& report, ReadNonlinearHdf5(reference_case, report.generated_results_hdf5); report.copied_input_path = hdf5.copied_input_path; report.generated_log_path = DerivedLogPath(hdf5.copied_input_path); + report.drilling_diagnostics = hdf5.drilling_diagnostics; RequireExactNodeLabels(hdf5.nodes, contract.node_count); RequirePressureInstance(hdf5.nodes); const auto displacement_rows = diff --git a/tests/reference/nonlinear_static_reference_comparison.h b/tests/reference/nonlinear_static_reference_comparison.h index a0ac0aa..f5de95d 100644 --- a/tests/reference/nonlinear_static_reference_comparison.h +++ b/tests/reference/nonlinear_static_reference_comparison.h @@ -10,6 +10,7 @@ #include #include "fesa/core/status.h" +#include "mitc4_drilling_report.h" #include "reference_tolerance_policy.h" namespace fesa::test { @@ -114,6 +115,7 @@ struct NonlinearReferenceComparisonReport { ConstrainedReactionResultantReport constrained_reaction_resultant; bool external_blocking_passed; bool internal_blocking_passed; + std::optional drilling_diagnostics; }; enum class NonlinearReferenceClassification { diff --git a/tests/test_compare_hht_dynamic.py b/tests/test_compare_hht_dynamic.py index 3bbe9a7..6328712 100644 --- a/tests/test_compare_hht_dynamic.py +++ b/tests/test_compare_hht_dynamic.py @@ -349,3 +349,88 @@ def test_cli_writes_json_and_exit_two_for_missing_results(tmp_path, artifacts): evidence = json.loads(report.read_text()) assert not evidence["precheck"]["passed"] assert evidence["families"] == [] + + +DRILLING_MODES = ("direct-initial-director-spring", "legacy-corotated-twist") + + +def drilling_result(tmp_path, artifacts, mode=DRILLING_MODES[0]): + path = changed_result(tmp_path, artifacts) + with h5py.File(path, "r+") as file: + file["metadata"].attrs["mitc4_drilling_formulation"] = mode + file["metadata"].attrs["drilling_diagnostics_version"] = np.uint64(1) + for frame in range(101): + root = file.require_group(f"steps/Step-1/frames/{frame}/global") + raw = np.array([1., 2., 3., 4., 5., 6.], dtype="mode = mode; + ASSERT_TRUE(fixture.state.CommitAcceptedFrame(frame).IsOk()); + } + TempFile path{"drilling-roundtrip"}; + ASSERT_TRUE(fesa::DynamicHdf5ResultsWriter{}.Write(path.Path(), *fixture.domain, fixture.state, {}).IsOk()); + Hdf5Handle file{H5Fopen(path.Path().string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT), H5Fclose}; + const auto stored = fesa::hdf5_internal::ReadDrillingDiagnostics(file.Get()); + ASSERT_TRUE(stored); + EXPECT_EQ(stored->formulation, fesa::Mitc4DrillingModeName(mode)); + ASSERT_EQ(stored->frames.size(), 3U); + for (std::size_t i = 0; i < 3U; ++i) { + EXPECT_EQ(stored->frames[i].frame_index, i); + EXPECT_EQ(stored->frames[i].raw_equilibrium, fixture.state.AcceptedFrames()[i].global_equilibrium); + EXPECT_EQ(stored->frames[i].equilibrium_with_drilling, + fixture.state.AcceptedFrames()[i].mitc4_drilling->equilibrium_with_drilling); + } + } +} + +TEST(DynamicHdf5ResultsWriter, DrillingMissingActualStatePreservesPriorFinal) { + auto fixture = MakeFixture(); + for (std::uint32_t i = 0; i < 3U; ++i) { + auto frame = MakeFrame(*fixture.domain, i); + frame.mitc4_drilling.reset(); + ASSERT_TRUE(fixture.state.CommitAcceptedFrame(frame).IsOk()); + } + TempFile path{"drilling-missing"}; + std::ofstream(path.Path(), std::ios::binary) << "prior"; + EXPECT_FALSE(fesa::DynamicHdf5ResultsWriter{}.Write(path.Path(), *fixture.domain, fixture.state, {}).IsOk()); + std::ifstream stream{path.Path(), std::ios::binary}; + EXPECT_EQ((std::string{std::istreambuf_iterator{stream}, {}}), "prior"); +} + + +TEST(DynamicHdf5ResultsWriter, DrillingLifecycleFaultsPreserveFinalBytes) { + auto fixture = MakeFixture(); + CommitHistory(fixture); + TempFile path{"drilling-faults"}; + for (const std::string fault : {"candidate-close", "reader-close", "markers", "datasets"}) { + SCOPED_TRACE(fault); + std::ofstream(path.Path(), std::ios::binary | std::ios::trunc) << "prior"; + EXPECT_FALSE(fesa::hdf5_internal::WriteDynamicFile(path.Path(), *fixture.domain, fixture.state, {}, + fesa::test::Mitc4LifecycleFault(fault)).IsOk()); + std::ifstream stream{path.Path(), std::ios::binary}; + EXPECT_EQ((std::string{std::istreambuf_iterator{stream}, {}}), "prior"); + } +} diff --git a/tests/unit/io/hdf5/eigen_hdf5_results_writer_test.cpp b/tests/unit/io/hdf5/eigen_hdf5_results_writer_test.cpp index 926d667..e3510d1 100644 --- a/tests/unit/io/hdf5/eigen_hdf5_results_writer_test.cpp +++ b/tests/unit/io/hdf5/eigen_hdf5_results_writer_test.cpp @@ -19,6 +19,10 @@ #include "fesa/model/domain.h" #include "io/hdf5/eigen_hdf5_self_check.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +#include "mitc4_fault_test_helpers.h" + namespace { fesa::ModelDefinition Definition(const std::filesystem::path& source) { @@ -81,6 +85,10 @@ TEST(EigenHdf5ResultsWriter, WritesModalSchemaAndRejectsTamperedState) { const hid_t file = H5Fopen(path.string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT); ASSERT_GE(file, 0); + const auto drilling = fesa::hdf5_internal::ReadDrillingDiagnostics(file); + ASSERT_TRUE(drilling); + EXPECT_EQ(drilling->formulation, fesa::Mitc4DrillingModeName(*state.drilling_mode)); + EXPECT_TRUE(drilling->frames.empty()); EXPECT_GT(H5Lexists(file, "/steps/Step-1/eigen/modes", H5P_DEFAULT), 0); const hid_t controls = H5Dopen2(file, "/steps/Step-1/eigen/backend/base_fpm", H5P_DEFAULT); @@ -105,7 +113,8 @@ TEST(EigenHdf5ResultsWriter, WritesModalSchemaAndRejectsTamperedState) { [](auto& candidate) { candidate.constraint_mask[0] = 1U; }, [](auto& candidate) { candidate.nodal_reduction[0].rcond_inf = 1e-15; }, [](auto& candidate) { candidate.backend.base_fpm[5U] = 1; }, - [](auto& candidate) { candidate.eta_gram = 1e-5; }}; + [](auto& candidate) { candidate.eta_gram = 1e-5; }, + [](auto& candidate) { candidate.drilling_mode.reset(); }}; for (const auto& mutate : mutations) { auto candidate = state; mutate(candidate); @@ -133,6 +142,27 @@ TEST(EigenHdf5ResultsWriter, WritesModalSchemaAndRejectsTamperedState) { locked_final.close(); ASSERT_TRUE( fesa::hdf5_internal::SelfCheckEigenFile(path, domain, state, {}).IsOk()); + for (const auto mode : {fesa::Mitc4DrillingMode::kDirectSpring, + fesa::Mitc4DrillingMode::kLegacyCorotatedTwist}) { + state.drilling_mode = mode; + ASSERT_TRUE(writer.Write(path, domain, state, {}).IsOk()); + fesa::hdf5_internal::Hdf5Handle stored_file{ + H5Fopen(path.string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT), H5Fclose}; + const auto stored = fesa::hdf5_internal::ReadDrillingDiagnostics(stored_file.Get()); + ASSERT_TRUE(stored); + EXPECT_EQ(stored->formulation, fesa::Mitc4DrillingModeName(mode)); + EXPECT_TRUE(stored->frames.empty()); + } + { + fesa::hdf5_internal::Hdf5Handle stored_file{ + H5Fopen(path.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + fesa::hdf5_internal::Hdf5Handle metadata{ + H5Gopen2(stored_file.Get(), "/metadata", H5P_DEFAULT), H5Gclose}; + ASSERT_GE(H5Adelete(metadata.Get(), "mitc4_drilling_formulation"), 0); + ASSERT_GE(H5Adelete(metadata.Get(), "drilling_diagnostics_version"), 0); + } + EXPECT_FALSE(fesa::hdf5_internal::SelfCheckEigenFile(path, domain, state, {}).IsOk()); + ASSERT_TRUE(writer.Write(path, domain, state, {}).IsOk()); const hid_t corrupt_file = H5Fopen(path.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT); ASSERT_GE(corrupt_file, 0); @@ -153,3 +183,22 @@ TEST(EigenHdf5ResultsWriter, WritesModalSchemaAndRejectsTamperedState) { } } // namespace + + +TEST(EigenHdf5ResultsWriter, DrillingLifecycleFaultsPreserveFinalBytes) { + const auto path = std::filesystem::temp_directory_path() / "fesa-eigen-drilling-fault.h5"; + auto domain = fesa::Domain::Create(Definition(path.parent_path() / "shell.inp")); + ASSERT_TRUE(domain.HasValue()); + auto state = fesa::LinearFrequencyProcedureExecutor{}.Execute(domain.Value()); + ASSERT_TRUE(state.HasValue()); + for (const std::string fault : {"candidate-close", "reader-close", "markers"}) { + SCOPED_TRACE(fault); + std::ofstream(path, std::ios::binary | std::ios::trunc) << "prior"; + EXPECT_FALSE(fesa::hdf5_internal::WriteEigenFile(path, domain.Value(), state.Value(), {}, + fesa::test::Mitc4LifecycleFault(fault)).IsOk()); + std::ifstream stream{path, std::ios::binary}; + EXPECT_EQ((std::string{std::istreambuf_iterator{stream}, {}}), "prior"); + } + std::error_code ignored; + std::filesystem::remove(path, ignored); +} diff --git a/tests/unit/io/hdf5/hdf5_results_writer_test.cpp b/tests/unit/io/hdf5/hdf5_results_writer_test.cpp index 5250cd4..aeb0ab9 100644 --- a/tests/unit/io/hdf5/hdf5_results_writer_test.cpp +++ b/tests/unit/io/hdf5/hdf5_results_writer_test.cpp @@ -15,6 +15,7 @@ #include #include #include +#include #include #include #include @@ -27,6 +28,13 @@ #include "fesa/math/vector3.h" #include "fesa/model/domain.h" +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +#include "mitc4_fault_test_helpers.h" + +#include "io/hdf5/hdf5_result_writer.h" +#include "io/hdf5/hdf5_self_check.h" + namespace { constexpr const char* kStepRoot = "/steps/Step-1/frames/0"; @@ -241,7 +249,8 @@ fesa::ModelDefinition MakeShellDefinition(const std::filesystem::path& source) { return definition; } -WriterFixture MakeShellFixture(const std::filesystem::path& source) { +WriterFixture MakeShellFixture(const std::filesystem::path& source, + std::optional mode = fesa::Mitc4DrillingMode::kDirectSpring) { auto domain_result = fesa::Domain::Create(MakeShellDefinition(source)); if (!domain_result.HasValue()) { throw std::runtime_error{ @@ -303,6 +312,8 @@ WriterFixture MakeShellFixture(const std::filesystem::path& source) { candidate.physical_strain_energy = 123.5; candidate.equilibrium = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0}; candidate.verification_metrics = {1.0e-13, 2.0e-13, 3.0e-13}; + if (mode) candidate.mitc4_drilling = fesa::Mitc4DrillingEvidence{ + *mode, {0.5, 1.0, 1.5}, {1.0, 2.0, 3.0, 3.5, 4.0, 4.5}}; const fesa::Status commit = state->CommitShellResults({0U}, std::move(candidate)); if (!commit.IsOk()) { @@ -1343,3 +1354,149 @@ TEST(Hdf5ResultsWriter, EXPECT_EQ(ReadBytes(final), valid_bytes); EXPECT_EQ(EntryCount(directory.Path()), 1U); } + + +TEST(Hdf5ResultsWriter, DrillingActualModesRoundTripWithoutDefaultInference) { + for (const auto mode : {fesa::Mitc4DrillingMode::kDirectSpring, + fesa::Mitc4DrillingMode::kLegacyCorotatedTwist}) { + TempDirectory directory{"drilling-roundtrip"}; + auto fixture = MakeShellFixture(directory.Path() / "shell.inp", mode); + const auto output = directory.Path() / "results.h5"; + ASSERT_TRUE(fesa::Hdf5ResultsWriter{}.Write(output, *fixture.domain, *fixture.state, {}).IsOk()); + const auto file = OpenFile(output); + const auto stored = fesa::hdf5_internal::ReadDrillingDiagnostics(file.Get()); + ASSERT_TRUE(stored); + EXPECT_EQ(stored->formulation, fesa::Mitc4DrillingModeName(mode)); + ASSERT_EQ(stored->frames.size(), 1U); + EXPECT_EQ(stored->frames[0].raw_equilibrium, fixture.state->Equilibrium()); + EXPECT_EQ(stored->frames[0].drilling_moment, fixture.state->Mitc4Drilling()->drilling_moment); + EXPECT_EQ(stored->frames[0].equilibrium_with_drilling, + fixture.state->Mitc4Drilling()->equilibrium_with_drilling); + } +} + +TEST(Hdf5ResultsWriter, DrillingMissingActualStatePreservesFinal) { + TempDirectory directory{"drilling-missing"}; + auto fixture = MakeShellFixture(directory.Path() / "shell.inp", std::nullopt); + const auto output = directory.Path() / "results.h5"; + const std::vector sentinel{'p','r','i','o','r'}; + WriteBytes(output, sentinel); + EXPECT_FALSE(fesa::Hdf5ResultsWriter{}.Write(output, *fixture.domain, *fixture.state, {}).IsOk()); + EXPECT_EQ(ReadBytes(output), sentinel); + EXPECT_EQ(EntryCount(directory.Path()), 1U); +} + +TEST(Hdf5ResultsWriter, DrillingReadbackCannotTreatRemovedMarkersAsHistorical) { + for (const bool remove_markers : {true, false}) { + TempDirectory directory{"drilling-readback"}; + auto fixture = MakeShellFixture(directory.Path() / "shell.inp"); + const auto output = directory.Path() / "results.h5"; + ASSERT_TRUE(fesa::Hdf5ResultsWriter{}.Write(output, *fixture.domain, *fixture.state, {}).IsOk()); + { + fesa::hdf5_internal::Hdf5Handle file{ + H5Fopen(output.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + if (remove_markers) { + fesa::hdf5_internal::Hdf5Handle metadata{ + H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), H5Gclose}; + ASSERT_GE(H5Adelete(metadata.Get(), "mitc4_drilling_formulation"), 0); + ASSERT_GE(H5Adelete(metadata.Get(), "drilling_diagnostics_version"), 0); + EXPECT_FALSE(fesa::hdf5_internal::ReadDrillingDiagnostics(file.Get())); + } else { + ASSERT_GE(H5Ldelete(file.Get(), "/steps/Step-1/frames/0/global/drilling_moment", H5P_DEFAULT), 0); + ASSERT_GE(H5Ldelete(file.Get(), "/steps/Step-1/frames/0/global/equilibrium_with_drilling", H5P_DEFAULT), 0); + } + } + const auto plan = fesa::hdf5_internal::BuildWriterPlan(output, *fixture.domain, *fixture.state, {}); + ASSERT_TRUE(plan.HasValue()); + EXPECT_FALSE(fesa::hdf5_internal::SelfCheck(output, *fixture.domain, *fixture.state, 0U, + plan.Value().model_data).IsOk()); + } +} + + +TEST(Hdf5ResultsWriter, DrillingReaderRejectsMalformedAndInconsistentEvidence) { + for (const std::string mutation : {"partial-mode", "partial-version", "mode", "version", "version-type", + "mode-type", "shape", "type", "endian", "nan", "inf", "algebra", "force"}) { + SCOPED_TRACE(mutation); + TempDirectory directory{"drilling-invalid"}; + auto fixture = MakeShellFixture(directory.Path() / "shell.inp"); + const auto output = directory.Path() / "results.h5"; + ASSERT_TRUE(fesa::Hdf5ResultsWriter{}.Write(output, *fixture.domain, *fixture.state, {}).IsOk()); + fesa::hdf5_internal::Hdf5ErrorSilencer silence; + fesa::hdf5_internal::Hdf5Handle file{ + H5Fopen(output.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + fesa::hdf5_internal::Hdf5Handle metadata{ + H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), H5Gclose}; + const std::string path = "/steps/Step-1/frames/0/global/drilling_moment"; + if (mutation == "partial-mode" || mutation == "partial-version") { + ASSERT_GE(H5Adelete(metadata.Get(), mutation == "partial-mode" ? + "mitc4_drilling_formulation" : "drilling_diagnostics_version"), 0); + } else if (mutation == "mode" || mutation == "mode-type") { + ASSERT_GE(H5Adelete(metadata.Get(), "mitc4_drilling_formulation"), 0); + if (mutation == "mode") fesa::hdf5_internal::WriteStringAttribute(metadata.Get(), "mitc4_drilling_formulation", "unknown"); + else fesa::hdf5_internal::WriteUint64Attribute(metadata.Get(), "mitc4_drilling_formulation", 1U); + } else if (mutation == "version" || mutation == "version-type") { + ASSERT_GE(H5Adelete(metadata.Get(), "drilling_diagnostics_version"), 0); + if (mutation == "version") fesa::hdf5_internal::WriteUint64Attribute(metadata.Get(), "drilling_diagnostics_version", 2U); + else fesa::hdf5_internal::WriteDoubleAttribute(metadata.Get(), "drilling_diagnostics_version", 1.0); + } else if (mutation == "shape" || mutation == "type" || mutation == "endian") { + ASSERT_GE(H5Ldelete(file.Get(), path.c_str(), H5P_DEFAULT), 0); + const hsize_t size = mutation == "shape" ? 2U : 3U; + fesa::hdf5_internal::Hdf5Handle space{H5Screate_simple(1, &size, nullptr), H5Sclose}; + fesa::hdf5_internal::Hdf5Handle dataset{H5Dcreate2(file.Get(), path.c_str(), + mutation == "type" ? H5T_IEEE_F32LE : (mutation == "endian" ? H5T_IEEE_F64BE : H5T_IEEE_F64LE), + space.Get(), H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT), H5Dclose}; + ASSERT_GE(dataset.Get(), 0); + } else { + const std::string target = mutation == "inf" ? "/steps/Step-1/frames/0/global/equilibrium" : + (mutation == "nan" ? path : "/steps/Step-1/frames/0/global/equilibrium_with_drilling"); + fesa::hdf5_internal::Hdf5Handle dataset{H5Dopen2(file.Get(), target.c_str(), H5P_DEFAULT), H5Dclose}; + std::array values{1., 2., 3., 3.5, 4., 4.5}; + if (mutation == "nan") values[0] = std::numeric_limits::quiet_NaN(); + else if (mutation == "inf") values[0] = std::numeric_limits::infinity(); + else values[mutation == "force" ? 0U : 3U] = 1e-100; + ASSERT_GE(H5Dwrite(dataset.Get(), H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, H5P_DEFAULT, values.data()), 0); + } + EXPECT_THROW((void)fesa::hdf5_internal::ReadDrillingDiagnostics(file.Get()), + fesa::hdf5_internal::Hdf5Failure); + } +} + + +TEST(Hdf5ResultsWriter, DrillingHistoricalSelfCheckDoesNotInferMode) { + TempDirectory directory{"drilling-historical"}; + auto fixture = MakeShellFixture(directory.Path() / "shell.inp"); + const auto output = directory.Path() / "results.h5"; + ASSERT_TRUE(fesa::Hdf5ResultsWriter{}.Write(output, *fixture.domain, *fixture.state, {}).IsOk()); + { + fesa::hdf5_internal::Hdf5Handle file{ + H5Fopen(output.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + fesa::hdf5_internal::Hdf5Handle metadata{H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), H5Gclose}; + ASSERT_GE(H5Adelete(metadata.Get(), "mitc4_drilling_formulation"), 0); + ASSERT_GE(H5Adelete(metadata.Get(), "drilling_diagnostics_version"), 0); + ASSERT_GE(H5Ldelete(file.Get(), "/steps/Step-1/frames/0/global/drilling_moment", H5P_DEFAULT), 0); + ASSERT_GE(H5Ldelete(file.Get(), "/steps/Step-1/frames/0/global/equilibrium_with_drilling", H5P_DEFAULT), 0); + } + auto historical = MakeShellFixture(directory.Path() / "shell.inp", std::nullopt); + const auto plan = fesa::hdf5_internal::BuildWriterPlan(output, *historical.domain, *historical.state, {}); + ASSERT_TRUE(plan.HasValue()); + EXPECT_TRUE(fesa::hdf5_internal::SelfCheck(output, *historical.domain, *historical.state, 0U, + plan.Value().model_data).IsOk()); + EXPECT_FALSE(historical.state->Mitc4Drilling()); +} + + +TEST(Hdf5ResultsWriter, DrillingLifecycleFaultsPreserveFinalBytes) { + TempDirectory directory{"drilling-faults"}; + auto fixture = MakeShellFixture(directory.Path() / "shell.inp"); + const auto path = directory.Path() / "results.h5"; + const std::vector prior{'p','r','i','o','r'}; + for (const std::string fault : {"candidate-close", "reader-close", "markers", "datasets"}) { + SCOPED_TRACE(fault); + WriteBytes(path, prior); + EXPECT_FALSE(fesa::hdf5_internal::WriteStaticFile(path, *fixture.domain, *fixture.state, {}, + fesa::test::Mitc4LifecycleFault(fault)).IsOk()); + EXPECT_EQ(ReadBytes(path), prior); + EXPECT_EQ(EntryCount(directory.Path()), 1U); + } +} diff --git a/tests/unit/io/hdf5/mitc4_fault_test_helpers.h b/tests/unit/io/hdf5/mitc4_fault_test_helpers.h new file mode 100644 index 0000000..6479fd3 --- /dev/null +++ b/tests/unit/io/hdf5/mitc4_fault_test_helpers.h @@ -0,0 +1,42 @@ +#ifndef FESA_TESTS_UNIT_IO_HDF5_MITC4_FAULT_TEST_HELPERS_H_ +#define FESA_TESTS_UNIT_IO_HDF5_MITC4_FAULT_TEST_HELPERS_H_ + +#include + +#include "io/hdf5/mitc4_file_hooks.h" +#include "io/hdf5/hdf5_raii.h" + +namespace fesa::test { + +inline herr_t CloseMitc4FileAndReportFailure(hid_t file) { + (void)H5Fclose(file); + return -1; +} + +inline hdf5_internal::Mitc4FileHooks Mitc4LifecycleFault(const std::string& fault) { + hdf5_internal::Mitc4FileHooks hooks; + if (fault == "candidate-close") hooks.candidate_close = CloseMitc4FileAndReportFailure; + if (fault == "reader-close") hooks.reader_close = CloseMitc4FileAndReportFailure; + if (fault == "markers" || fault == "datasets") { + hooks.before_self_check = [fault](const std::filesystem::path& path) { + hdf5_internal::Hdf5Handle file{hdf5_internal::RequireHdf5Id( + H5Fopen(path.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), "Cannot corrupt test candidate."), H5Fclose}; + if (fault == "markers") { + hdf5_internal::Hdf5Handle metadata{H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), H5Gclose}; + hdf5_internal::RequireHdf5(H5Adelete(metadata.Get(), "mitc4_drilling_formulation"), "Cannot remove mode."); + hdf5_internal::RequireHdf5(H5Adelete(metadata.Get(), "drilling_diagnostics_version"), "Cannot remove version."); + } else { + hdf5_internal::RequireHdf5(H5Ldelete(file.Get(), + "/steps/Step-1/frames/0/global/drilling_moment", H5P_DEFAULT), "Cannot remove moment."); + hdf5_internal::RequireHdf5(H5Ldelete(file.Get(), + "/steps/Step-1/frames/0/global/equilibrium_with_drilling", H5P_DEFAULT), "Cannot remove corrected balance."); + } + hdf5_internal::RequireHdf5(file.CloseChecked(), "Cannot close corrupted candidate."); + }; + } + return hooks; +} + +} // namespace fesa::test + +#endif // FESA_TESTS_UNIT_IO_HDF5_MITC4_FAULT_TEST_HELPERS_H_ diff --git a/tests/unit/io/hdf5/nonlinear_hdf5_results_writer_test.cpp b/tests/unit/io/hdf5/nonlinear_hdf5_results_writer_test.cpp index ee8d73b..c357821 100644 --- a/tests/unit/io/hdf5/nonlinear_hdf5_results_writer_test.cpp +++ b/tests/unit/io/hdf5/nonlinear_hdf5_results_writer_test.cpp @@ -13,6 +13,7 @@ #include #include #include +#include #include #include #include @@ -25,6 +26,11 @@ #include "fesa/model/domain.h" #include "fesa/results/nonlinear_results_writer.h" #include "io/hdf5/hdf5_raii.h" +#include "io/hdf5/hdf5_self_check.h" + +#include "io/hdf5/mitc4_drilling_diagnostics.h" + +#include "mitc4_fault_test_helpers.h" namespace { @@ -184,7 +190,8 @@ std::vector MakeShellRows(const double seed) { fesa::AcceptedFrameState MakeFrame(const fesa::NonlinearTrialState& committed, const std::uint32_t frame_index, const std::uint32_t increment_index, - const double load_factor) { + const double load_factor, + std::optional mode = fesa::Mitc4DrillingMode::kDirectSpring) { fesa::AcceptedFrameState frame{ frame_index, increment_index, @@ -223,10 +230,13 @@ fesa::AcceptedFrameState MakeFrame(const fesa::NonlinearTrialState& committed, frame.reaction[offset + component] = frame.residual[offset + component]; } } + if (mode) frame.mitc4_drilling = fesa::Mitc4DrillingEvidence{ + *mode, {}, frame.global_equilibrium}; return frame; } -Fixture MakeFixture(const std::filesystem::path& source) { +Fixture MakeFixture(const std::filesystem::path& source, + std::optional mode = fesa::Mitc4DrillingMode::kDirectSpring) { auto domain_result = fesa::Domain::Create(MakeDefinition(source)); if (!domain_result.HasValue()) { throw std::runtime_error{"Unable to construct nonlinear HDF5 Domain."}; @@ -248,10 +258,10 @@ Fixture MakeFixture(const std::filesystem::path& source) { } fesa::NonlinearAnalysisState state = std::move(state_result.Value()); if (!state.CommitAcceptedFrame( - MakeFrame(state.CommittedState(), 0U, 1U, 0.5)) + MakeFrame(state.CommittedState(), 0U, 1U, 0.5, mode)) .IsOk() || !state.CommitAcceptedFrame( - MakeFrame(state.CommittedState(), 1U, 2U, 1.0)) + MakeFrame(state.CommittedState(), 1U, 2U, 1.0, mode)) .IsOk()) { throw std::runtime_error{"Unable to commit nonlinear HDF5 test frames."}; } @@ -615,3 +625,55 @@ TEST(NonlinearHdf5ResultsWriter, CandidateFinalizationFailurePreservesPriorFinal EXPECT_EQ(status.Diagnostics()[0U].code, "hdf5-finalization-failure"); EXPECT_EQ(ReadBytes(output), sentinel); } + + +TEST(NonlinearHdf5ResultsWriter, DrillingModesAndEveryAcceptedFrameRoundTrip) { + for (const auto mode : {fesa::Mitc4DrillingMode::kDirectSpring, + fesa::Mitc4DrillingMode::kLegacyCorotatedTwist}) { + TempDirectory directory{"drilling"}; + auto fixture = MakeFixture(directory.Path() / "shell.inp", mode); + const auto path = directory.Path() / "results.h5"; + ASSERT_TRUE(fesa::NonlinearHdf5ResultsWriter{}.Write(path, *fixture.domain, fixture.state, {}).IsOk()); + { + const auto file = OpenFile(path); + const auto stored = fesa::hdf5_internal::ReadDrillingDiagnostics(file.Get()); + ASSERT_TRUE(stored); + EXPECT_EQ(stored->formulation, fesa::Mitc4DrillingModeName(mode)); + ASSERT_EQ(stored->frames.size(), fixture.state.AcceptedFrames().size()); + EXPECT_EQ(stored->frames.front().frame_index, 0U); + EXPECT_EQ(stored->frames.back().raw_equilibrium, fixture.state.AcceptedFrames().back().global_equilibrium); + } + { + Hdf5Handle file{H5Fopen(path.string().c_str(), H5F_ACC_RDWR, H5P_DEFAULT), H5Fclose}; + Hdf5Handle metadata{H5Gopen2(file.Get(), "/metadata", H5P_DEFAULT), H5Gclose}; + ASSERT_GE(H5Adelete(metadata.Get(), "mitc4_drilling_formulation"), 0); + ASSERT_GE(H5Adelete(metadata.Get(), "drilling_diagnostics_version"), 0); + } + EXPECT_FALSE(fesa::hdf5_internal::NonlinearSelfCheck(path, *fixture.domain, fixture.state, 0U).IsOk()); + } +} + +TEST(NonlinearHdf5ResultsWriter, DrillingMissingActualStatePreservesPriorFinal) { + TempDirectory directory{"drilling-missing"}; + auto fixture = MakeFixture(directory.Path() / "shell.inp", std::nullopt); + const auto path = directory.Path() / "results.h5"; + std::ofstream(path, std::ios::binary) << "prior"; + EXPECT_FALSE(fesa::NonlinearHdf5ResultsWriter{}.Write(path, *fixture.domain, fixture.state, {}).IsOk()); + const auto bytes = ReadBytes(path); + EXPECT_EQ(std::string(bytes.begin(), bytes.end()), "prior"); +} + + +TEST(NonlinearHdf5ResultsWriter, DrillingLifecycleFaultsPreserveFinalBytes) { + TempDirectory directory{"drilling-faults"}; + auto fixture = MakeFixture(directory.Path() / "shell.inp"); + const auto path = directory.Path() / "results.h5"; + for (const std::string fault : {"candidate-close", "reader-close", "markers", "datasets"}) { + SCOPED_TRACE(fault); + std::ofstream(path, std::ios::binary | std::ios::trunc) << "prior"; + EXPECT_FALSE(fesa::hdf5_internal::WriteNonlinearFile(path, *fixture.domain, fixture.state, {}, + fesa::test::Mitc4LifecycleFault(fault)).IsOk()); + const auto bytes = ReadBytes(path); + EXPECT_EQ(std::string(bytes.begin(), bytes.end()), "prior"); + } +}