diff --git a/cmake/FesaDependencies.cmake b/cmake/FesaDependencies.cmake index 0944d96..0e025e3 100644 --- a/cmake/FesaDependencies.cmake +++ b/cmake/FesaDependencies.cmake @@ -44,10 +44,13 @@ if(NOT TARGET TBB::tbb) message(FATAL_ERROR "The oneTBB CONFIG package did not define TBB::tbb") endif() -if(TARGET hdf5::hdf5-static) - set(_fesa_hdf5_target hdf5::hdf5-static) -elseif(TARGET hdf5::hdf5-shared) +# Prefer the packaged shared C runtime when both variants exist. The Windows +# static archive may carry compiler-runtime requirements that are not part of +# the HDF5 imported target's transitive link interface. +if(TARGET hdf5::hdf5-shared) set(_fesa_hdf5_target hdf5::hdf5-shared) +elseif(TARGET hdf5::hdf5-static) + set(_fesa_hdf5_target hdf5::hdf5-static) elseif(TARGET HDF5::HDF5) set(_fesa_hdf5_target HDF5::HDF5) else() diff --git a/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md b/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md index b7285e3..766e442 100644 --- a/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md +++ b/docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md @@ -1130,3 +1130,82 @@ the subsequent uncontended serial build, targeted, discovery, full sequence followed by a separate successful audit command. - concerns: none; no upstream-contract conflict or critical blocker remains. + +## Step 23 — hdf5-results-writer + +- task_id: `TASK-23` +- status: `completed` +- changed_files: `include/fesa/results/results_writer.hpp`, + `include/fesa/io/hdf5/hdf5_results_writer.hpp`, + `src/fesa/io/hdf5/hdf5_results_writer.cpp`, + `tests/unit/results/results_writer_test.cpp`, + `tests/unit/io/hdf5/hdf5_results_writer_test.cpp`, + `src/fesa/CMakeLists.txt`, `tests/CMakeLists.txt`, + `cmake/FesaDependencies.cmake`, + `docs/implementation-plans/linear-static-3d-euler-beam-implementation-report.md`, + `phases/linear-static-3d-euler-beam/index.json` +- requirement_ids: `FESA-REQ-LS3DEB-003`, `FESA-REQ-LS3DEB-015`, + `FESA-REQ-LS3DEB-019`, `FESA-REQ-LS3DEB-020`, + `FESA-REQ-LS3DEB-023`, `FESA-REQ-LS3DEB-025`, + `FESA-REQ-LS3DEB-028`, `FESA-REQ-LS3DEB-029`, + `FESA-REQ-LS3DEB-030`, `FESA-REQ-LS3DEB-031`, + `FESA-REQ-LS3DEB-032`, `FESA-REQ-LS3DEB-034` +- test_ids: `T23-HDF5-001`, `T23-HDF5-002`, `T23-HDF5-003`, + `T23-HDF5-004`, `T23-HDF5-005` + +| stage | exact command | exit_code | expected_or_observed_result | evidence_tail | +| --- | --- | ---: | --- | --- | +| RED-build | `cmake --build .harness/build --config Debug --target fesa_tests` | 1 | All five ledger tests were registered before either public production header existed | MSVC C1083 reported missing `fesa/results/results_writer.hpp` and `fesa/io/hdf5/hdf5_results_writer.hpp` | +| INFRA-link | `cmake --build .harness/build --config Debug --target fesa_tests` | 1 | Production and tests compiled, but the packaged static HDF5 target could not link in this environment | LNK1104 reported missing `libircmt.lib`; this is infrastructure evidence, not the feature RED | +| INFRA-runtime-discovery | `cmake --build .harness/build --config Debug --target fesa_tests` | 1 | Shared HDF5 linked after normalized-target repair, but post-build GoogleTest discovery could not load a transitive runtime | discovery exited `0xc0000135`; `dumpbin /dependents hdf5.dll` identified `libmmd.dll` | +| GREEN-build | `cmake --build .harness/build --config Debug --target fesa_tests` | 0 | Minimum schema-v0 writer, tests, solver library, and unit executable compile/link after the approved runtime repair | `hdf5_results_writer.cpp` and both tests built under `/W4 /WX`; required runtime DLLs were staged before discovery | +| GREEN-test | `ctest --test-dir .harness/build -C Debug -R Hdf5ResultsWriter --output-on-failure` | 0 | Exact schema, mandatory outputs, diagnostics/centroid, failure preservation, and successful replacement pass | 5/5 exact `Hdf5ResultsWriter` tests passed | +| VERIFY-configure | `cmake -S . -B .harness/build -A x64 -DFESA_GTEST_SOURCE_DIR=C:/git/googletest "-DMKL_DIR=C:/Program Files (x86)/Intel/oneAPI/mkl/2026.1/lib/cmake/mkl" "-DTBB_DIR=C:/Program Files (x86)/Intel/oneAPI/tbb/2023.1/lib/cmake/tbb" "-DHDF5_DIR=C:/Program Files/HDF_Group/HDF5/2.1.1/cmake"` | 0 | Approved explicit-dependency MSVC x64 build tree generates | Visual Studio 18/MSVC x64, Windows SDK 10.0.26100.0, oneMKL, oneTBB, and shared HDF5 resolved; configure/generate completed | +| VERIFY-build | `cmake --build .harness/build --config Debug` | 0 | Full Debug build passes without a new FESA warning | `fesa_solver.lib` and `fesa_unit_tests.exe` built under `/W4 /WX` | +| VERIFY-targeted | `ctest --test-dir .harness/build -C Debug -R Hdf5ResultsWriter --output-on-failure` | 0 | Focused Step 23 suite remains green | 5/5 exact `Hdf5ResultsWriter` tests passed | +| VERIFY-discovery | `ctest --test-dir .harness/build -C Debug --show-only=json-v1` | 0 | CTest discovers the accumulated suite and all five exact ledger names | 73 tests discovered, including exactly 5 `Hdf5ResultsWriter` tests | +| VERIFY-full | `ctest --test-dir .harness/build -C Debug --output-on-failure` | 0 | Full accumulated C++ suite has zero failures | 73/73 tests passed in 2.75 seconds | +| VERIFY-runtime | `& 'C:\Program Files\Microsoft Visual Studio\18\Community\VC\Tools\MSVC\14.51.36231\bin\Hostx64\x64\dumpbin.exe' /dependents 'C:\Program Files\HDF_Group\HDF5\2.1.1\bin\hdf5.dll'` | 0 | The package DLL's undeclared Intel runtime requirement is traceable | dependencies include `libmmd.dll`; the WIN32 staging guard copies `${OMP_DLL_DIR}/libmmd.dll` | +| VERIFY-contract-scans | Public HDF5/Win32 type scan, exact-test count, finalizer/RAII scope scan, forbidden shear/averaging/filter scan, `git diff --check`, and reference diff/status | 0 | Backend boundary, exact mandatory-output scope, atomic finalization, and reference immutability remain isolated | public backend leaks 0; exact tests 5; checked file closes and exact Replace/Move split present; forbidden scope 0; whitespace clean; reference unchanged | + +- contract_checks: `ResultsWriter` and `Hdf5ResultsWriter` expose only the + approved backend-free signatures. HDF5, Win32 handles, macros, RAII wrappers, + schema construction, self-check, and atomic finalization stay private to the + implementation; `Fesa::HDF5` is linked privately. +- contract_checks: schema v0 writes every literal metadata, model, and + `/steps/Step-1/frames/0` path, including `generalized_resultant`. Numeric + fields use IEEE little-endian float64, stable identities use little-endian + uint64 compound members, text uses UTF-8 variable-length strings, and exact + component/unit/coordinate/location/step/frame attributes are pinned. +- contract_checks: non-axis-aligned element axes, node/element identity, + connectivity, and representative first/last nonzero values for displacement, + reaction, end action, section resultant, generalized strain, and generalized + resultant prove solved state is connected to output. Output requests do not + filter any mandatory field; shear/averaged stress was not added. +- contract_checks: diagnostics are stable-sorted and preserved, an empty input + produces an exact zero-row compound dataset, and the default centroid stress + row serializes axial `S11` only. Dimensions, row identity/order, source + identity, and all numeric values are validated before finalization. +- atomicity: a unique same-directory temporary file is fully written, flushed, + checked-closed, reopened read-only for exact schema self-check, and + checked-closed before finalization. Existing finals use only + `ReplaceFileW(final,temp,nullptr,...)`; new finals use only + `MoveFileExW(temp,final,MOVEFILE_WRITE_THROUGH)`. Failure removes the temp, + leaves no partial new final, and preserves pre-existing final bytes. +- infrastructure_resolution: the approved normalized HDF5 dependency now + prefers the package shared target, with static and `HDF5::HDF5` fallbacks. + Tests stage `$` plus the package DLL's + verified `${OMP_DLL_DIR}/libmmd.dll` dependency without directly naming an + HDF5 package target or hardcoding an Intel installation path. +- generated_evidence: `.harness/build/src/fesa/Debug/fesa_solver.lib`, + `.harness/build/tests/Debug/fesa_unit_tests.exe`; all test HDF5/temp artifacts + are confined to build-tree temporary directories and cleaned by their owners. +- reference_diff: unchanged; `git diff --exit-code -- reference/` and + `git status --short -- reference/` both report no change. +- handoff: Step 24 can construct the backend-free writer and rely on mandatory + schema-v0 output and fail-atomic replacement. Its production executable must + stage the selected HDF5 target runtime DLLs and `${OMP_DLL_DIR}/libmmd.dll` + with the same portable target-expression/runtime-guard policy used by tests. +- concerns: no remaining critical implementation or upstream-contract blocker. + The HDF5 static-link and shared-runtime environment blockers were resolved + and are recorded separately from the valid missing-API RED. diff --git a/include/fesa/io/hdf5/hdf5_results_writer.hpp b/include/fesa/io/hdf5/hdf5_results_writer.hpp new file mode 100644 index 0000000..89fa817 --- /dev/null +++ b/include/fesa/io/hdf5/hdf5_results_writer.hpp @@ -0,0 +1,17 @@ +#pragma once + +#include "fesa/results/results_writer.hpp" + +namespace fesa { + +// Writes schema-v0 output while keeping backend and platform types private. +class Hdf5ResultsWriter final : public ResultsWriter { +public: + Status write( + const std::filesystem::path& outputPath, + const Domain& domain, + const AnalysisState& state, + const std::vector& diagnostics) override; +}; + +} // namespace fesa diff --git a/include/fesa/results/results_writer.hpp b/include/fesa/results/results_writer.hpp new file mode 100644 index 0000000..3352976 --- /dev/null +++ b/include/fesa/results/results_writer.hpp @@ -0,0 +1,25 @@ +#pragma once + +#include "fesa/analysis/analysis_state.hpp" +#include "fesa/core/diagnostic.hpp" +#include "fesa/core/status.hpp" +#include "fesa/model/domain.hpp" + +#include +#include + +namespace fesa { + +// Keeps the solver core independent of the authoritative result-storage backend. +class ResultsWriter { +public: + virtual ~ResultsWriter() = default; + + virtual Status write( + const std::filesystem::path& outputPath, + const Domain& domain, + const AnalysisState& state, + const std::vector& diagnostics) = 0; +}; + +} // namespace fesa diff --git a/src/fesa/CMakeLists.txt b/src/fesa/CMakeLists.txt index 0fe18f8..c966e6a 100644 --- a/src/fesa/CMakeLists.txt +++ b/src/fesa/CMakeLists.txt @@ -14,6 +14,7 @@ add_library( fem/dof_manager.cpp io/abaqus/domain_mapper.cpp io/abaqus/input_reader.cpp + io/hdf5/hdf5_results_writer.cpp math/matrix.cpp math/sparse_matrix.cpp math/vector.cpp @@ -33,6 +34,7 @@ target_link_libraries( PRIVATE Fesa::MKL Fesa::TBB + Fesa::HDF5 ) # Product warnings are strict without imposing FESA policy on external targets. diff --git a/src/fesa/io/hdf5/hdf5_results_writer.cpp b/src/fesa/io/hdf5/hdf5_results_writer.cpp new file mode 100644 index 0000000..b3b58b5 --- /dev/null +++ b/src/fesa/io/hdf5/hdf5_results_writer.cpp @@ -0,0 +1,1552 @@ +#define NOMINMAX +#include + +#include "fesa/io/hdf5/hdf5_results_writer.hpp" + +#include "fesa/build_info.hpp" + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace fesa { +namespace { + +constexpr std::size_t kDofsPerNode = 6U; +constexpr std::size_t kEndpointCount = 2U; +constexpr std::size_t kGaussPointCount = 2U; +constexpr std::size_t kEndActionComponentCount = 6U; +constexpr std::size_t kGeneralizedComponentCount = 4U; +constexpr const char* kStepName = "Step-1"; +constexpr std::size_t kFrameIndex = 0U; +constexpr const char* kStepRoot = "/steps/Step-1/frames/0"; + +class Hdf5Handle { +public: + using Closer = herr_t (*)(hid_t); + + Hdf5Handle() = default; + Hdf5Handle(const hid_t value, Closer closer) + : value_{value}, closer_{closer} {} + Hdf5Handle(const Hdf5Handle&) = delete; + Hdf5Handle& operator=(const Hdf5Handle&) = delete; + Hdf5Handle(Hdf5Handle&& other) noexcept + : value_{other.value_}, closer_{other.closer_} { + other.value_ = -1; + other.closer_ = nullptr; + } + Hdf5Handle& operator=(Hdf5Handle&& other) noexcept { + if (this != &other) { + reset(); + value_ = other.value_; + closer_ = other.closer_; + other.value_ = -1; + other.closer_ = nullptr; + } + return *this; + } + ~Hdf5Handle() { reset(); } + + hid_t get() const noexcept { return value_; } + herr_t closeChecked() noexcept { + if (value_ < 0 || closer_ == nullptr) { + return 0; + } + const hid_t value = value_; + const Closer closer = closer_; + value_ = -1; + closer_ = nullptr; + return closer(value); + } + +private: + void reset() noexcept { + if (value_ >= 0 && closer_ != nullptr) { + (void)closer_(value_); + } + value_ = -1; + closer_ = nullptr; + } + + hid_t value_{-1}; + Closer closer_{nullptr}; +}; + +// Expected backend failures become one structured FESA diagnostic, not an +// HDF5 error-stack dump mixed into deterministic CLI output. +class Hdf5ErrorSilencer { +public: + Hdf5ErrorSilencer() { + if (H5Eget_auto2(H5E_DEFAULT, &callback_, &clientData_) >= 0 && + H5Eset_auto2(H5E_DEFAULT, nullptr, nullptr) >= 0) { + active_ = true; + } + } + Hdf5ErrorSilencer(const Hdf5ErrorSilencer&) = delete; + Hdf5ErrorSilencer& operator=(const Hdf5ErrorSilencer&) = delete; + ~Hdf5ErrorSilencer() { + if (active_) { + (void)H5Eset_auto2(H5E_DEFAULT, callback_, clientData_); + } + } + +private: + H5E_auto2_t callback_{nullptr}; + void* clientData_{nullptr}; + bool active_{false}; +}; + +class TemporaryFileGuard { +public: + explicit TemporaryFileGuard(std::filesystem::path path) + : path_{std::move(path)} {} + TemporaryFileGuard(const TemporaryFileGuard&) = delete; + TemporaryFileGuard& operator=(const TemporaryFileGuard&) = delete; + ~TemporaryFileGuard() { + if (active_) { + std::error_code ignored; + (void)std::filesystem::remove(path_, ignored); + } + } + void release() noexcept { active_ = false; } + +private: + std::filesystem::path path_; + bool active_{true}; +}; + +class Hdf5Failure final : public std::runtime_error { +public: + explicit Hdf5Failure(const std::string& message) + : std::runtime_error{message} {} +}; + +void requireHdf5(const herr_t result, const char* message) { + if (result < 0) { + throw Hdf5Failure{message}; + } +} + +hid_t requireHdf5Id(const hid_t result, const char* message) { + if (result < 0) { + throw Hdf5Failure{message}; + } + return result; +} + +Status outputFailure(const std::string& code, const std::string& message) { + return Status::failure( + FailureCategory::output, + {{Severity::error, code, {}, "", "", message}}); +} + +bool isFinite(const std::array& values) { + return std::all_of(values.begin(), values.end(), [](const double value) { + return std::isfinite(value); + }); +} + +template +bool isFinite(const std::array& values) { + return std::all_of(values.begin(), values.end(), [](const double value) { + return std::isfinite(value); + }); +} + +bool isValidUtf8(const std::string& value) { + const auto* bytes = reinterpret_cast(value.data()); + std::size_t index = 0U; + while (index < value.size()) { + const unsigned char first = bytes[index]; + if (first <= 0x7fU) { + ++index; + continue; + } + + std::size_t continuationCount = 0U; + std::uint32_t codePoint = 0U; + if (first >= 0xc2U && first <= 0xdfU) { + continuationCount = 1U; + codePoint = first & 0x1fU; + } else if (first >= 0xe0U && first <= 0xefU) { + continuationCount = 2U; + codePoint = first & 0x0fU; + } else if (first >= 0xf0U && first <= 0xf4U) { + continuationCount = 3U; + codePoint = first & 0x07U; + } else { + return false; + } + if (index + continuationCount >= value.size()) { + return false; + } + for (std::size_t offset = 1U; offset <= continuationCount; ++offset) { + const unsigned char next = bytes[index + offset]; + if ((next & 0xc0U) != 0x80U) { + return false; + } + codePoint = (codePoint << 6U) | (next & 0x3fU); + } + if ((continuationCount == 2U && codePoint < 0x800U) || + (continuationCount == 3U && codePoint < 0x10000U) || + (codePoint >= 0xd800U && codePoint <= 0xdfffU) || + codePoint > 0x10ffffU) { + return false; + } + index += continuationCount + 1U; + } + return true; +} + +bool sameIdentity(const SourceEntityId& left, const SourceEntityId& right) { + return left.instanceName == right.instanceName && + left.sourceLabel == right.sourceLabel && + left.sourceLabelText == right.sourceLabelText; +} + +using AxisSet = std::array; + +bool computeLocalAxes( + const Domain& domain, const EulerBeam3DDefinition& element, AxisSet& axes) { + if (element.nodeIndices[0U] >= domain.nodes().size() || + element.nodeIndices[1U] >= domain.nodes().size() || + element.sectionIndex >= domain.sections().size()) { + return false; + } + const auto& first = domain.nodes()[element.nodeIndices[0U]].coordinates; + const auto& second = domain.nodes()[element.nodeIndices[1U]].coordinates; + const auto& guide = domain.sections()[element.sectionIndex].firstAxis; + const std::array delta = { + second[0U] - first[0U], + second[1U] - first[1U], + second[2U] - first[2U]}; + const double length = std::hypot(delta[0U], delta[1U], delta[2U]); + if (!isFinite(first) || !isFinite(second) || !isFinite(guide) || + !isFinite(delta) || !std::isfinite(length) || !(length > 0.0)) { + return false; + } + const std::array x = { + delta[0U] / length, delta[1U] / length, delta[2U] / length}; + const double projection = + guide[0U] * x[0U] + guide[1U] * x[1U] + guide[2U] * x[2U]; + const std::array yTrial = { + guide[0U] - projection * x[0U], + guide[1U] - projection * x[1U], + guide[2U] - projection * x[2U]}; + const double yNorm = std::hypot(yTrial[0U], yTrial[1U], yTrial[2U]); + if (!isFinite(yTrial) || !std::isfinite(yNorm) || !(yNorm > 0.0)) { + return false; + } + const std::array y = { + yTrial[0U] / yNorm, yTrial[1U] / yNorm, yTrial[2U] / yNorm}; + const std::array z = { + x[1U] * y[2U] - x[2U] * y[1U], + x[2U] * y[0U] - x[0U] * y[2U], + x[0U] * y[1U] - x[1U] * y[0U]}; + axes = { + x[0U], x[1U], x[2U], + y[0U], y[1U], y[2U], + z[0U], z[1U], z[2U]}; + return isFinite(axes); +} + +bool sizeProductFits( + const std::size_t left, const std::size_t right, std::size_t& product) { + if (right != 0U && left > (std::numeric_limits::max)() / right) { + return false; + } + product = left * right; + return true; +} + +Status validateWriterInput( + const std::filesystem::path& outputPath, + const Domain& domain, + const AnalysisState& state, + const std::vector& diagnostics, + std::vector& localAxes) { + if (outputPath.empty() || outputPath.filename().empty()) { + return outputFailure( + "invalid-output-path", "The HDF5 output path must name a file."); + } + if (state.identity().stepName != kStepName || + state.identity().frameIndex != kFrameIndex) { + return outputFailure( + "invalid-result-state", + "Schema v0 requires literal Step-1 and frame index 0."); + } + + std::size_t fullDofCount = 0U; + if (!sizeProductFits(domain.nodes().size(), kDofsPerNode, fullDofCount)) { + return outputFailure( + "invalid-result-state", "The nodal result shape overflows size_t."); + } + const std::array vectors = { + &state.displacement(), + &state.externalForce(), + &state.internalForce(), + &state.residual(), + &state.reaction()}; + for (const Vector* vector : vectors) { + if (vector->size() != fullDofCount) { + return outputFailure( + "invalid-result-state", + "Every V0 analysis vector must have node_count*6 values."); + } + for (std::size_t index = 0U; index < vector->size(); ++index) { + if (!std::isfinite((*vector)[index])) { + return outputFailure( + "invalid-result-state", + "Every V0 analysis vector value must be finite."); + } + } + } + + if (domain.sourcePath().empty() || domain.sourceContentIdentity().empty() || + !isValidUtf8(domain.sourceContentIdentity())) { + return outputFailure( + "invalid-result-identity", + "Source path and UTF-8 content identity are required."); + } + for (const Node& node : domain.nodes()) { + if (node.sourceId.sourceLabel <= 0 || + node.sourceId.sourceLabelText.empty() || + !isValidUtf8(node.sourceId.instanceName) || + !isValidUtf8(node.sourceId.sourceLabelText) || + !isFinite(node.coordinates)) { + return outputFailure( + "invalid-result-identity", + "Every node requires finite coordinates and UTF-8 source identity."); + } + } + + localAxes.clear(); + localAxes.reserve(domain.elements().size()); + for (const EulerBeam3DDefinition& element : domain.elements()) { + AxisSet axes{}; + if (element.sourceId.sourceLabel <= 0 || + element.sourceId.sourceLabelText.empty() || + !isValidUtf8(element.sourceId.instanceName) || + !isValidUtf8(element.sourceId.sourceLabelText) || + element.nodeIndices[0U] == element.nodeIndices[1U] || + element.materialIndex >= domain.materials().size() || + !computeLocalAxes(domain, element, axes)) { + return outputFailure( + "invalid-result-identity", + "Every element requires valid source, connectivity, property, and local-axis identity."); + } + localAxes.push_back(axes); + } + + std::size_t endpointCount = 0U; + std::size_t gaussCount = 0U; + if (!sizeProductFits(domain.elements().size(), kEndpointCount, endpointCount) || + !sizeProductFits(domain.elements().size(), kGaussPointCount, gaussCount) || + state.endpointResults().size() != endpointCount || + state.gaussResults().size() != gaussCount) { + return outputFailure( + "invalid-result-rows", + "Endpoint and Gauss row counts must match every element and location."); + } + for (std::size_t rowIndex = 0U; + rowIndex < state.endpointResults().size(); + ++rowIndex) { + const EntityIndex expectedElement = + static_cast(rowIndex / kEndpointCount); + const int expectedEndpoint = static_cast(rowIndex % kEndpointCount); + const EndpointResultRow& row = state.endpointResults()[rowIndex]; + const auto& element = domain.elements()[expectedElement]; + const auto& expectedNode = + domain.nodes()[element.nodeIndices[static_cast(expectedEndpoint)]]; + if (row.element != expectedElement || row.endpoint != expectedEndpoint || + !sameIdentity(row.node, expectedNode.sourceId) || + !isFinite(row.endAction) || !isFinite(row.sectionResultant)) { + return outputFailure( + "invalid-result-rows", + "Endpoint result rows must follow element/endpoint order and identity."); + } + } + for (std::size_t rowIndex = 0U; + rowIndex < state.gaussResults().size(); + ++rowIndex) { + const EntityIndex expectedElement = + static_cast(rowIndex / kGaussPointCount); + const int expectedGaussPoint = + static_cast(rowIndex % kGaussPointCount) + 1; + const GaussResultRow& row = state.gaussResults()[rowIndex]; + if (row.element != expectedElement || + row.gaussPoint != expectedGaussPoint || + !isFinite(row.generalizedStrain) || + !isFinite(row.generalizedResultant)) { + return outputFailure( + "invalid-result-rows", + "Gauss result rows must follow element/Gauss order and identity."); + } + } + + std::size_t stressIndex = 0U; + for (std::size_t elementIndex = 0U; + elementIndex < domain.elements().size(); + ++elementIndex) { + const auto& element = domain.elements()[elementIndex]; + const auto& sectionPoints = + domain.sections()[element.sectionIndex].sectionPoints; + for (std::size_t gauss = 0U; gauss < kGaussPointCount; ++gauss) { + const std::size_t count = sectionPoints.empty() ? 1U : sectionPoints.size(); + for (std::size_t point = 0U; point < count; ++point) { + if (stressIndex >= state.stressResults().size()) { + return outputFailure( + "invalid-result-rows", + "Axial stress rows are missing required element/Gauss/section locations."); + } + const StressS11Row& row = state.stressResults()[stressIndex++]; + const std::size_t expectedPoint = sectionPoints.empty() ? 0U : point + 1U; + const double expectedX1 = sectionPoints.empty() ? 0.0 : sectionPoints[point][0U]; + const double expectedX2 = sectionPoints.empty() ? 0.0 : sectionPoints[point][1U]; + const char* expectedSource = sectionPoints.empty() ? "fesa-default" : "input"; + if (row.element != static_cast(elementIndex) || + row.gaussPoint != static_cast(gauss + 1U) || + row.sectionPoint != expectedPoint || + row.x1 != expectedX1 || row.x2 != expectedX2 || + row.source != expectedSource || !isValidUtf8(row.source) || + !std::isfinite(row.x1) || !std::isfinite(row.x2) || + !std::isfinite(row.s11)) { + return outputFailure( + "invalid-result-rows", + "Axial stress rows must follow exact element/Gauss/section identity."); + } + } + } + } + if (stressIndex != state.stressResults().size()) { + return outputFailure( + "invalid-result-rows", "Axial stress output contains extra rows."); + } + + for (const Diagnostic& diagnostic : diagnostics) { + if (!isValidUtf8(diagnostic.code) || + !isValidUtf8(diagnostic.keyword) || + !isValidUtf8(diagnostic.entityIdentity) || + !isValidUtf8(diagnostic.message)) { + return outputFailure( + "invalid-result-diagnostic", + "Diagnostic text must be valid UTF-8."); + } + } + return Status::ok(); +} + +std::string normalizedPathString(const std::filesystem::path& path) { + if (path.empty()) { + return {}; + } + return std::filesystem::absolute(path).lexically_normal().generic_u8string(); +} + +std::string sourceInputIdentity(const Domain& domain) { + return "path=" + normalizedPathString(domain.sourcePath()) + + ";content_identity=" + domain.sourceContentIdentity(); +} + +Hdf5Handle makeUtf8StringType() { + Hdf5Handle type{ + requireHdf5Id(H5Tcopy(H5T_C_S1), "Unable to copy an HDF5 string type."), + H5Tclose}; + requireHdf5( + H5Tset_size(type.get(), H5T_VARIABLE), + "Unable to configure a variable-length HDF5 string."); + requireHdf5( + H5Tset_cset(type.get(), H5T_CSET_UTF8), + "Unable to configure UTF-8 HDF5 strings."); + requireHdf5( + H5Tset_strpad(type.get(), H5T_STR_NULLTERM), + "Unable to configure HDF5 string padding."); + return type; +} + +void writeStringAttribute( + const hid_t object, const char* name, const std::string& value) { + auto type = makeUtf8StringType(); + Hdf5Handle space{ + requireHdf5Id(H5Screate(H5S_SCALAR), "Unable to create an attribute space."), + H5Sclose}; + Hdf5Handle attribute{ + requireHdf5Id( + H5Acreate2(object, name, type.get(), space.get(), H5P_DEFAULT, H5P_DEFAULT), + "Unable to create an HDF5 string attribute."), + H5Aclose}; + const char* raw = value.c_str(); + requireHdf5( + H5Awrite(attribute.get(), type.get(), &raw), + "Unable to write an HDF5 string attribute."); +} + +void writeUint64Attribute( + const hid_t object, const char* name, const std::uint64_t value) { + Hdf5Handle space{ + requireHdf5Id(H5Screate(H5S_SCALAR), "Unable to create an attribute space."), + H5Sclose}; + Hdf5Handle attribute{ + requireHdf5Id( + H5Acreate2(object, name, H5T_STD_U64LE, space.get(), H5P_DEFAULT, H5P_DEFAULT), + "Unable to create an HDF5 integer attribute."), + H5Aclose}; + requireHdf5( + H5Awrite(attribute.get(), H5T_NATIVE_UINT64, &value), + "Unable to write an HDF5 integer attribute."); +} + +Hdf5Handle createGroup(const hid_t file, const char* path) { + Hdf5Handle linkProperties{ + requireHdf5Id(H5Pcreate(H5P_LINK_CREATE), "Unable to create link properties."), + H5Pclose}; + requireHdf5( + H5Pset_create_intermediate_group(linkProperties.get(), 1U), + "Unable to enable intermediate HDF5 groups."); + return Hdf5Handle{ + requireHdf5Id( + H5Gcreate2(file, path, linkProperties.get(), H5P_DEFAULT, H5P_DEFAULT), + "Unable to create an HDF5 group."), + H5Gclose}; +} + +Hdf5Handle createDatasetSpace(const std::vector& dimensions) { + return Hdf5Handle{ + requireHdf5Id( + H5Screate_simple( + static_cast(dimensions.size()), dimensions.data(), nullptr), + "Unable to create an HDF5 dataset space."), + H5Sclose}; +} + +void writeResultAttributes( + const hid_t dataset, + const std::string& componentNames, + const std::string& componentUnits, + const std::string& coordinateSystem, + const std::string& location) { + writeStringAttribute(dataset, "component_names", componentNames); + writeStringAttribute( + dataset, "component_unit_dimensions", componentUnits); + writeStringAttribute(dataset, "coordinate_system", coordinateSystem); + writeStringAttribute(dataset, "location", location); + writeStringAttribute(dataset, "step_name", kStepName); + writeUint64Attribute(dataset, "frame_index", 0U); +} + +void writeDoubleDataset( + const hid_t file, + const std::string& path, + const std::vector& dimensions, + const double* values, + const std::size_t valueCount, + const std::string& componentNames, + const std::string& componentUnits, + const std::string& coordinateSystem, + const std::string& location) { + auto space = createDatasetSpace(dimensions); + Hdf5Handle dataset{ + requireHdf5Id( + H5Dcreate2( + file, path.c_str(), H5T_IEEE_F64LE, space.get(), + H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT), + "Unable to create a floating-point HDF5 dataset."), + H5Dclose}; + if (valueCount != 0U) { + requireHdf5( + H5Dwrite( + dataset.get(), H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, + H5P_DEFAULT, values), + "Unable to write a floating-point HDF5 dataset."); + } + writeResultAttributes( + dataset.get(), componentNames, componentUnits, coordinateSystem, location); +} + +struct NodeWriteRow { + std::uint64_t internalNodeId; + const char* instanceName; + const char* sourceLabel; + double coordinates[3]; +}; + +struct ElementWriteRow { + std::uint64_t internalElementId; + const char* instanceName; + const char* sourceLabel; + std::uint64_t nodeInternalIds[2]; + double localAxes[9]; +}; + +struct StressWriteRow { + std::uint64_t internalElementId; + std::uint64_t gaussPointIndex; + std::uint64_t sectionPointIndex; + double x1; + double x2; + const char* source; + double s11; +}; + +struct DiagnosticWriteRow { + const char* severity; + const char* code; + const char* file; + std::uint64_t line; + const char* keyword; + const char* entityIdentity; + const char* message; +}; + +Hdf5Handle writeCompoundDataset( + const hid_t file, + const char* path, + const std::size_t rowCount, + const hid_t fileType, + const hid_t memoryType, + const void* rows) { + const std::vector dimensions = { + static_cast(rowCount)}; + auto space = createDatasetSpace(dimensions); + Hdf5Handle dataset{ + requireHdf5Id( + H5Dcreate2( + file, path, fileType, space.get(), H5P_DEFAULT, + H5P_DEFAULT, H5P_DEFAULT), + "Unable to create a compound HDF5 dataset."), + H5Dclose}; + if (rowCount != 0U) { + requireHdf5( + H5Dwrite( + dataset.get(), memoryType, H5S_ALL, H5S_ALL, + H5P_DEFAULT, rows), + "Unable to write a compound HDF5 dataset."); + } + return dataset; +} + +void writeMetadata(const hid_t file, const Domain& domain) { + auto metadata = createGroup(file, "/metadata"); + writeUint64Attribute(metadata.get(), "schema_version", 0U); + writeStringAttribute( + metadata.get(), "feature_id", "linear-static-3d-euler-beam"); + writeStringAttribute( + metadata.get(), "solver_version", std::string{solverVersion()}); + writeStringAttribute( + metadata.get(), "source_input_identity", sourceInputIdentity(domain)); + writeStringAttribute( + metadata.get(), "unit_system_label", "user-consistent-unspecified"); + writeStringAttribute( + metadata.get(), + "coordinate_convention", + "global-cartesian; beam-local=(t,n1,t-cross-n1)"); + writeStringAttribute( + metadata.get(), "element_formulation", "B33-3D-Euler-Bernoulli"); + writeStringAttribute(metadata.get(), "step_name", kStepName); + writeUint64Attribute(metadata.get(), "frame_index", 0U); +} + +void writeNodes(const hid_t file, const Domain& domain) { + std::vector rows; + rows.reserve(domain.nodes().size()); + for (std::size_t index = 0U; index < domain.nodes().size(); ++index) { + const Node& node = domain.nodes()[index]; + rows.push_back({ + static_cast(index), + node.sourceId.instanceName.c_str(), + node.sourceId.sourceLabelText.c_str(), + {node.coordinates[0U], node.coordinates[1U], node.coordinates[2U]}}); + } + + auto stringType = makeUtf8StringType(); + const hsize_t coordinateDimensions[] = {3U}; + Hdf5Handle fileCoordinates{ + requireHdf5Id( + H5Tarray_create2(H5T_IEEE_F64LE, 1, coordinateDimensions), + "Unable to create the node coordinate file type."), + H5Tclose}; + Hdf5Handle memoryCoordinates{ + requireHdf5Id( + H5Tarray_create2(H5T_NATIVE_DOUBLE, 1, coordinateDimensions), + "Unable to create the node coordinate memory type."), + H5Tclose}; + Hdf5Handle fileType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(NodeWriteRow)), + "Unable to create the node file type."), + H5Tclose}; + Hdf5Handle memoryType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(NodeWriteRow)), + "Unable to create the node memory type."), + H5Tclose}; + requireHdf5( + H5Tinsert(fileType.get(), "internal_node_id", + HOFFSET(NodeWriteRow, internalNodeId), H5T_STD_U64LE), + "Unable to define the node internal ID field."); + requireHdf5( + H5Tinsert(fileType.get(), "instance_name", + HOFFSET(NodeWriteRow, instanceName), stringType.get()), + "Unable to define the node instance field."); + requireHdf5( + H5Tinsert(fileType.get(), "source_label", + HOFFSET(NodeWriteRow, sourceLabel), stringType.get()), + "Unable to define the node source-label field."); + requireHdf5( + H5Tinsert(fileType.get(), "coordinates", + HOFFSET(NodeWriteRow, coordinates), fileCoordinates.get()), + "Unable to define the node coordinate field."); + requireHdf5( + H5Tinsert(memoryType.get(), "internal_node_id", + HOFFSET(NodeWriteRow, internalNodeId), H5T_NATIVE_UINT64), + "Unable to define the node internal ID memory field."); + requireHdf5( + H5Tinsert(memoryType.get(), "instance_name", + HOFFSET(NodeWriteRow, instanceName), stringType.get()), + "Unable to define the node instance memory field."); + requireHdf5( + H5Tinsert(memoryType.get(), "source_label", + HOFFSET(NodeWriteRow, sourceLabel), stringType.get()), + "Unable to define the node source-label memory field."); + requireHdf5( + H5Tinsert(memoryType.get(), "coordinates", + HOFFSET(NodeWriteRow, coordinates), memoryCoordinates.get()), + "Unable to define the node coordinate memory field."); + + auto dataset = writeCompoundDataset( + file, "/model/nodes", rows.size(), fileType.get(), memoryType.get(), + rows.data()); + writeStringAttribute(dataset.get(), "coordinate_system", "global-cartesian"); + writeStringAttribute(dataset.get(), "units_label", "length"); +} + +void writeElements( + const hid_t file, const Domain& domain, const std::vector& axes) { + std::vector rows; + rows.reserve(domain.elements().size()); + for (std::size_t index = 0U; index < domain.elements().size(); ++index) { + const auto& element = domain.elements()[index]; + ElementWriteRow row{ + static_cast(index), + element.sourceId.instanceName.c_str(), + element.sourceId.sourceLabelText.c_str(), + {static_cast(element.nodeIndices[0U]), + static_cast(element.nodeIndices[1U])}, + {}}; + std::copy(axes[index].begin(), axes[index].end(), row.localAxes); + rows.push_back(row); + } + + auto stringType = makeUtf8StringType(); + const hsize_t nodeDimensions[] = {2U}; + const hsize_t axisDimensions[] = {3U, 3U}; + Hdf5Handle fileNodes{ + requireHdf5Id( + H5Tarray_create2(H5T_STD_U64LE, 1, nodeDimensions), + "Unable to create the element connectivity file type."), + H5Tclose}; + Hdf5Handle memoryNodes{ + requireHdf5Id( + H5Tarray_create2(H5T_NATIVE_UINT64, 1, nodeDimensions), + "Unable to create the element connectivity memory type."), + H5Tclose}; + Hdf5Handle fileAxes{ + requireHdf5Id( + H5Tarray_create2(H5T_IEEE_F64LE, 2, axisDimensions), + "Unable to create the local-axis file type."), + H5Tclose}; + Hdf5Handle memoryAxes{ + requireHdf5Id( + H5Tarray_create2(H5T_NATIVE_DOUBLE, 2, axisDimensions), + "Unable to create the local-axis memory type."), + H5Tclose}; + Hdf5Handle fileType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(ElementWriteRow)), + "Unable to create the element file type."), + H5Tclose}; + Hdf5Handle memoryType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(ElementWriteRow)), + "Unable to create the element memory type."), + H5Tclose}; + const auto insertFields = [&](const hid_t type, + const hid_t integerType, + const hid_t nodeType, + const hid_t axesType) { + requireHdf5( + H5Tinsert(type, "internal_element_id", + HOFFSET(ElementWriteRow, internalElementId), integerType), + "Unable to define the element internal ID field."); + requireHdf5( + H5Tinsert(type, "instance_name", + HOFFSET(ElementWriteRow, instanceName), stringType.get()), + "Unable to define the element instance field."); + requireHdf5( + H5Tinsert(type, "source_label", + HOFFSET(ElementWriteRow, sourceLabel), stringType.get()), + "Unable to define the element source-label field."); + requireHdf5( + H5Tinsert(type, "node_internal_ids", + HOFFSET(ElementWriteRow, nodeInternalIds), nodeType), + "Unable to define the element connectivity field."); + requireHdf5( + H5Tinsert(type, "local_axes", + HOFFSET(ElementWriteRow, localAxes), axesType), + "Unable to define the element local-axis field."); + }; + insertFields(fileType.get(), H5T_STD_U64LE, fileNodes.get(), fileAxes.get()); + insertFields( + memoryType.get(), H5T_NATIVE_UINT64, memoryNodes.get(), memoryAxes.get()); + auto dataset = writeCompoundDataset( + file, "/model/elements", rows.size(), fileType.get(), memoryType.get(), + rows.data()); + writeStringAttribute( + dataset.get(), "formulation", "B33-3D-Euler-Bernoulli"); +} + +std::vector flattenEndpointValues( + const std::vector& rows, + const bool sectionResultants) { + const std::size_t components = + sectionResultants ? kGeneralizedComponentCount : kEndActionComponentCount; + std::vector values; + values.reserve(rows.size() * components); + for (const auto& row : rows) { + if (sectionResultants) { + values.insert( + values.end(), row.sectionResultant.begin(), row.sectionResultant.end()); + } else { + values.insert(values.end(), row.endAction.begin(), row.endAction.end()); + } + } + return values; +} + +std::vector flattenGaussValues( + const std::vector& rows, + const bool resultants) { + std::vector values; + values.reserve(rows.size() * kGeneralizedComponentCount); + for (const auto& row : rows) { + const auto& rowValues = + resultants ? row.generalizedResultant : row.generalizedStrain; + values.insert(values.end(), rowValues.begin(), rowValues.end()); + } + return values; +} + +void writeStress(const hid_t file, const AnalysisState& state) { + std::vector rows; + rows.reserve(state.stressResults().size()); + for (const auto& row : state.stressResults()) { + rows.push_back({ + static_cast(row.element), + static_cast(row.gaussPoint), + static_cast(row.sectionPoint), + row.x1, + row.x2, + row.source.c_str(), + row.s11}); + } + auto stringType = makeUtf8StringType(); + Hdf5Handle fileType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(StressWriteRow)), + "Unable to create the stress file type."), + H5Tclose}; + Hdf5Handle memoryType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(StressWriteRow)), + "Unable to create the stress memory type."), + H5Tclose}; + const auto insertFields = [&](const hid_t type, + const hid_t integerType, + const hid_t doubleType) { + requireHdf5( + H5Tinsert(type, "internal_element_id", + HOFFSET(StressWriteRow, internalElementId), integerType), + "Unable to define the stress element field."); + requireHdf5( + H5Tinsert(type, "gauss_point_index", + HOFFSET(StressWriteRow, gaussPointIndex), integerType), + "Unable to define the stress Gauss field."); + requireHdf5( + H5Tinsert(type, "section_point_index", + HOFFSET(StressWriteRow, sectionPointIndex), integerType), + "Unable to define the stress section-point field."); + requireHdf5( + H5Tinsert(type, "x1", HOFFSET(StressWriteRow, x1), doubleType), + "Unable to define the stress x1 field."); + requireHdf5( + H5Tinsert(type, "x2", HOFFSET(StressWriteRow, x2), doubleType), + "Unable to define the stress x2 field."); + requireHdf5( + H5Tinsert(type, "source", HOFFSET(StressWriteRow, source), stringType.get()), + "Unable to define the stress source field."); + requireHdf5( + H5Tinsert(type, "S11", HOFFSET(StressWriteRow, s11), doubleType), + "Unable to define the stress S11 field."); + }; + insertFields(fileType.get(), H5T_STD_U64LE, H5T_IEEE_F64LE); + insertFields(memoryType.get(), H5T_NATIVE_UINT64, H5T_NATIVE_DOUBLE); + auto dataset = writeCompoundDataset( + file, + "/steps/Step-1/frames/0/element/stress_s11", + rows.size(), + fileType.get(), + memoryType.get(), + rows.data()); + writeResultAttributes( + dataset.get(), "S11", "force/length^2", "beam-local", "section-point"); +} + +void writeDiagnostics( + const hid_t file, const std::vector& inputDiagnostics) { + std::vector diagnostics = inputDiagnostics; + sortDiagnostics(diagnostics); + std::vector files; + files.reserve(diagnostics.size()); + for (const auto& diagnostic : diagnostics) { + files.push_back(normalizedPathString(diagnostic.location.file)); + } + std::vector rows; + rows.reserve(diagnostics.size()); + for (std::size_t index = 0U; index < diagnostics.size(); ++index) { + const auto& diagnostic = diagnostics[index]; + rows.push_back({ + diagnostic.severity == Severity::warning ? "warning" : "error", + diagnostic.code.c_str(), + files[index].c_str(), + static_cast(diagnostic.location.line), + diagnostic.keyword.c_str(), + diagnostic.entityIdentity.c_str(), + diagnostic.message.c_str()}); + } + + auto stringType = makeUtf8StringType(); + Hdf5Handle fileType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(DiagnosticWriteRow)), + "Unable to create the diagnostic file type."), + H5Tclose}; + Hdf5Handle memoryType{ + requireHdf5Id( + H5Tcreate(H5T_COMPOUND, sizeof(DiagnosticWriteRow)), + "Unable to create the diagnostic memory type."), + H5Tclose}; + const auto insertFields = [&](const hid_t type, const hid_t integerType) { + requireHdf5( + H5Tinsert(type, "severity", HOFFSET(DiagnosticWriteRow, severity), + stringType.get()), + "Unable to define diagnostic severity."); + requireHdf5( + H5Tinsert(type, "code", HOFFSET(DiagnosticWriteRow, code), + stringType.get()), + "Unable to define diagnostic code."); + requireHdf5( + H5Tinsert(type, "file", HOFFSET(DiagnosticWriteRow, file), + stringType.get()), + "Unable to define diagnostic file."); + requireHdf5( + H5Tinsert(type, "line", HOFFSET(DiagnosticWriteRow, line), integerType), + "Unable to define diagnostic line."); + requireHdf5( + H5Tinsert(type, "keyword", HOFFSET(DiagnosticWriteRow, keyword), + stringType.get()), + "Unable to define diagnostic keyword."); + requireHdf5( + H5Tinsert(type, "entity_identity", + HOFFSET(DiagnosticWriteRow, entityIdentity), stringType.get()), + "Unable to define diagnostic entity identity."); + requireHdf5( + H5Tinsert(type, "message", HOFFSET(DiagnosticWriteRow, message), + stringType.get()), + "Unable to define diagnostic message."); + }; + insertFields(fileType.get(), H5T_STD_U64LE); + insertFields(memoryType.get(), H5T_NATIVE_UINT64); + (void)writeCompoundDataset( + file, "/diagnostics", rows.size(), fileType.get(), memoryType.get(), + rows.data()); +} + +void writeResultDatasets( + const hid_t file, const Domain& domain, const AnalysisState& state) { + const std::vector nodalDimensions = { + static_cast(domain.nodes().size()), kDofsPerNode}; + writeDoubleDataset( + file, + std::string{kStepRoot} + "/nodal/displacement", + nodalDimensions, + state.displacement().data(), + state.displacement().size(), + "UX,UY,UZ,URX,URY,URZ", + "length,length,length,radian,radian,radian", + "global-cartesian", + "nodal"); + writeDoubleDataset( + file, + std::string{kStepRoot} + "/nodal/reaction", + nodalDimensions, + state.reaction().data(), + state.reaction().size(), + "RF1,RF2,RF3,RM1,RM2,RM3", + "force,force,force,force*length,force*length,force*length", + "global-cartesian", + "nodal"); + + const std::vector endpointActionDimensions = { + static_cast(domain.elements().size()), + kEndpointCount, + kEndActionComponentCount}; + const std::vector generalizedDimensions = { + static_cast(domain.elements().size()), + kEndpointCount, + kGeneralizedComponentCount}; + const auto endActions = flattenEndpointValues(state.endpointResults(), false); + writeDoubleDataset( + file, + std::string{kStepRoot} + "/element/end_force_local", + endpointActionDimensions, + endActions.data(), + endActions.size(), + "FX,FY,FZ,MX,MY,MZ", + "force,force,force,force*length,force*length,force*length", + "beam-local", + "endpoint-outward-action"); + const auto sectionResultants = + flattenEndpointValues(state.endpointResults(), true); + writeDoubleDataset( + file, + std::string{kStepRoot} + "/element/section_resultant", + generalizedDimensions, + sectionResultants.data(), + sectionResultants.size(), + "N,T,My,Mz", + "force,force*length,force*length,force*length", + "beam-local", + "endpoint-positive-local-x-section-cut"); + const auto generalizedStrains = + flattenGaussValues(state.gaussResults(), false); + writeDoubleDataset( + file, + std::string{kStepRoot} + "/element/generalized_strain", + generalizedDimensions, + generalizedStrains.data(), + generalizedStrains.size(), + "epsilon0,kappa_x,kappa_y,kappa_z", + "1,1/length,1/length,1/length", + "beam-local", + "integration-point"); + const auto generalizedResultants = + flattenGaussValues(state.gaussResults(), true); + writeDoubleDataset( + file, + std::string{kStepRoot} + "/element/generalized_resultant", + generalizedDimensions, + generalizedResultants.data(), + generalizedResultants.size(), + "N,T,My,Mz", + "force,force*length,force*length,force*length", + "beam-local", + "integration-point"); + writeStress(file, state); +} + +void writeFile( + const std::filesystem::path& temporaryPath, + const Domain& domain, + const AnalysisState& state, + const std::vector& diagnostics, + const std::vector& localAxes) { + Hdf5Handle file{ + requireHdf5Id( + H5Fcreate( + temporaryPath.string().c_str(), H5F_ACC_EXCL, + H5P_DEFAULT, H5P_DEFAULT), + "Unable to create the temporary HDF5 file."), + H5Fclose}; + writeMetadata(file.get(), domain); + (void)createGroup(file.get(), "/model"); + (void)createGroup(file.get(), "/steps/Step-1/frames/0/nodal"); + (void)createGroup(file.get(), "/steps/Step-1/frames/0/element"); + writeNodes(file.get(), domain); + writeElements(file.get(), domain, localAxes); + writeResultDatasets(file.get(), domain, state); + writeDiagnostics(file.get(), diagnostics); + requireHdf5( + H5Fflush(file.get(), H5F_SCOPE_GLOBAL), + "Unable to flush the temporary HDF5 file."); + requireHdf5( + file.closeChecked(), + "Unable to close the temporary HDF5 file after writing."); +} + +Hdf5Handle openDatasetForCheck(const hid_t file, const char* path) { + return Hdf5Handle{ + requireHdf5Id( + H5Dopen2(file, path, H5P_DEFAULT), + "A mandatory HDF5 dataset is missing."), + H5Dclose}; +} + +std::vector checkedDimensions(const hid_t dataset) { + Hdf5Handle space{ + requireHdf5Id( + H5Dget_space(dataset), "Unable to inspect an HDF5 dataspace."), + H5Sclose}; + const int rank = H5Sget_simple_extent_ndims(space.get()); + if (rank < 0) { + throw Hdf5Failure{"Unable to inspect an HDF5 dataset rank."}; + } + std::vector dimensions(static_cast(rank)); + if (rank > 0) { + requireHdf5( + H5Sget_simple_extent_dims(space.get(), dimensions.data(), nullptr), + "Unable to inspect HDF5 dataset dimensions."); + } + return dimensions; +} + +void requireStringAttribute( + const hid_t object, const char* name, const std::string& expected) { + Hdf5Handle attribute{ + requireHdf5Id( + H5Aopen(object, name, H5P_DEFAULT), + "A mandatory HDF5 string attribute is missing."), + H5Aclose}; + Hdf5Handle type{ + requireHdf5Id( + H5Aget_type(attribute.get()), "Unable to inspect an HDF5 attribute type."), + H5Tclose}; + if (H5Tget_class(type.get()) != H5T_STRING || + H5Tis_variable_str(type.get()) <= 0 || + H5Tget_cset(type.get()) != H5T_CSET_UTF8) { + throw Hdf5Failure{"An HDF5 string attribute is not variable-length UTF-8."}; + } + char* raw = nullptr; + requireHdf5( + H5Aread(attribute.get(), type.get(), &raw), + "Unable to read an HDF5 string attribute."); + const std::string actual = raw == nullptr ? std::string{} : std::string{raw}; + if (raw != nullptr) { + requireHdf5(H5free_memory(raw), "Unable to release an HDF5 string value."); + } + if (actual != expected) { + throw Hdf5Failure{"An HDF5 string attribute has the wrong value."}; + } +} + +void requireUint64Attribute( + const hid_t object, const char* name, const std::uint64_t expected) { + Hdf5Handle attribute{ + requireHdf5Id( + H5Aopen(object, name, H5P_DEFAULT), + "A mandatory HDF5 integer attribute is missing."), + H5Aclose}; + Hdf5Handle type{ + requireHdf5Id( + H5Aget_type(attribute.get()), + "Unable to inspect an HDF5 integer attribute type."), + H5Tclose}; + if (H5Tget_class(type.get()) != H5T_INTEGER || + H5Tget_size(type.get()) != sizeof(std::uint64_t) || + H5Tget_sign(type.get()) != H5T_SGN_NONE || + H5Tequal(type.get(), H5T_STD_U64LE) <= 0) { + throw Hdf5Failure{"An HDF5 integer attribute is not portable uint64."}; + } + std::uint64_t actual = 0U; + requireHdf5( + H5Aread(attribute.get(), H5T_NATIVE_UINT64, &actual), + "Unable to read an HDF5 integer attribute."); + if (actual != expected) { + throw Hdf5Failure{"An HDF5 integer attribute has the wrong value."}; + } +} + +void requirePortableCompoundMember( + const hid_t compoundType, const unsigned index, const std::string& name) { + Hdf5Handle memberType{ + requireHdf5Id( + H5Tget_member_type(compoundType, index), + "Unable to inspect a compound HDF5 field type."), + H5Tclose}; + const bool isUint64 = + name == "internal_node_id" || name == "internal_element_id" || + name == "gauss_point_index" || name == "section_point_index" || + name == "line"; + const bool isFloat64 = name == "x1" || name == "x2" || name == "S11"; + const bool isString = + name == "instance_name" || name == "source_label" || + name == "source" || name == "severity" || name == "code" || + name == "file" || name == "keyword" || + name == "entity_identity" || name == "message"; + if (isUint64) { + if (H5Tget_class(memberType.get()) != H5T_INTEGER || + H5Tget_size(memberType.get()) != sizeof(std::uint64_t) || + H5Tget_sign(memberType.get()) != H5T_SGN_NONE || + H5Tequal(memberType.get(), H5T_STD_U64LE) <= 0) { + throw Hdf5Failure{"A compound HDF5 field is not portable uint64."}; + } + return; + } + if (isFloat64) { + if (H5Tget_class(memberType.get()) != H5T_FLOAT || + H5Tget_size(memberType.get()) != sizeof(double) || + H5Tequal(memberType.get(), H5T_IEEE_F64LE) <= 0) { + throw Hdf5Failure{"A compound HDF5 field is not float64."}; + } + return; + } + if (isString) { + if (H5Tget_class(memberType.get()) != H5T_STRING || + H5Tis_variable_str(memberType.get()) <= 0 || + H5Tget_cset(memberType.get()) != H5T_CSET_UTF8) { + throw Hdf5Failure{"A compound HDF5 field is not UTF-8."}; + } + return; + } + + if (H5Tget_class(memberType.get()) != H5T_ARRAY) { + throw Hdf5Failure{"A compound HDF5 array field has the wrong type."}; + } + const int rank = H5Tget_array_ndims(memberType.get()); + if (rank <= 0) { + throw Hdf5Failure{"Unable to inspect a compound HDF5 array rank."}; + } + std::vector dimensions(static_cast(rank)); + requireHdf5( + H5Tget_array_dims2(memberType.get(), dimensions.data()), + "Unable to inspect compound HDF5 array dimensions."); + Hdf5Handle baseType{ + requireHdf5Id( + H5Tget_super(memberType.get()), + "Unable to inspect a compound HDF5 array base type."), + H5Tclose}; + if (name == "node_internal_ids") { + if (dimensions != std::vector{2U} || + H5Tequal(baseType.get(), H5T_STD_U64LE) <= 0) { + throw Hdf5Failure{"Element connectivity is not uint64[2]."}; + } + } else if (name == "coordinates") { + if (dimensions != std::vector{3U} || + H5Tequal(baseType.get(), H5T_IEEE_F64LE) <= 0) { + throw Hdf5Failure{"Node coordinates are not float64[3]."}; + } + } else if (name == "local_axes") { + if (dimensions != std::vector{3U, 3U} || + H5Tequal(baseType.get(), H5T_IEEE_F64LE) <= 0) { + throw Hdf5Failure{"Element local axes are not float64[3,3]."}; + } + } else { + throw Hdf5Failure{"A compound HDF5 field has an unknown type contract."}; + } +} + +void requireResultAttributes( + const hid_t dataset, + const char* componentNames, + const char* componentUnits, + const char* coordinateSystem, + const char* location) { + requireStringAttribute(dataset, "component_names", componentNames); + requireStringAttribute( + dataset, "component_unit_dimensions", componentUnits); + requireStringAttribute(dataset, "coordinate_system", coordinateSystem); + requireStringAttribute(dataset, "location", location); + requireStringAttribute(dataset, "step_name", kStepName); + requireUint64Attribute(dataset, "frame_index", 0U); +} + +void requireDoubleDataset( + const hid_t file, + const char* path, + const std::vector& expectedDimensions, + const char* componentNames, + const char* componentUnits, + const char* coordinateSystem, + const char* location) { + auto dataset = openDatasetForCheck(file, path); + if (checkedDimensions(dataset.get()) != expectedDimensions) { + throw Hdf5Failure{"A floating-point HDF5 dataset has the wrong shape."}; + } + Hdf5Handle type{ + requireHdf5Id( + H5Dget_type(dataset.get()), "Unable to inspect an HDF5 dataset type."), + H5Tclose}; + if (H5Tget_class(type.get()) != H5T_FLOAT || + H5Tget_size(type.get()) != sizeof(double) || + H5Tequal(type.get(), H5T_IEEE_F64LE) <= 0) { + throw Hdf5Failure{"A result dataset is not IEEE-754 float64 little-endian."}; + } + requireResultAttributes( + dataset.get(), componentNames, componentUnits, coordinateSystem, location); +} + +void requireCompoundDataset( + const hid_t file, + const char* path, + const hsize_t expectedRows, + const std::vector& expectedMembers) { + auto dataset = openDatasetForCheck(file, path); + if (checkedDimensions(dataset.get()) != std::vector{expectedRows}) { + throw Hdf5Failure{"A compound HDF5 dataset has the wrong shape."}; + } + Hdf5Handle type{ + requireHdf5Id( + H5Dget_type(dataset.get()), "Unable to inspect a compound type."), + H5Tclose}; + if (H5Tget_class(type.get()) != H5T_COMPOUND || + H5Tget_nmembers(type.get()) != static_cast(expectedMembers.size())) { + throw Hdf5Failure{"A compound HDF5 dataset has the wrong field count."}; + } + for (std::size_t index = 0U; index < expectedMembers.size(); ++index) { + char* rawName = H5Tget_member_name( + type.get(), static_cast(index)); + if (rawName == nullptr) { + throw Hdf5Failure{"Unable to inspect a compound HDF5 field name."}; + } + const std::string actualName{rawName}; + requireHdf5( + H5free_memory(rawName), "Unable to release a compound field name."); + if (actualName != expectedMembers[index]) { + throw Hdf5Failure{"A compound HDF5 field has the wrong identity."}; + } + requirePortableCompoundMember( + type.get(), static_cast(index), expectedMembers[index]); + } +} + +void selfCheckFile( + const std::filesystem::path& path, + const Domain& domain, + const AnalysisState& state, + const std::size_t diagnosticCount) { + 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}; + { + Hdf5Handle metadata{ + requireHdf5Id( + H5Gopen2(file.get(), "/metadata", H5P_DEFAULT), + "The HDF5 metadata group is missing."), + H5Gclose}; + requireUint64Attribute(metadata.get(), "schema_version", 0U); + requireStringAttribute( + metadata.get(), "feature_id", "linear-static-3d-euler-beam"); + requireStringAttribute( + metadata.get(), "solver_version", std::string{solverVersion()}); + requireStringAttribute( + metadata.get(), "source_input_identity", sourceInputIdentity(domain)); + requireStringAttribute( + metadata.get(), "unit_system_label", "user-consistent-unspecified"); + requireStringAttribute( + metadata.get(), + "coordinate_convention", + "global-cartesian; beam-local=(t,n1,t-cross-n1)"); + requireStringAttribute( + metadata.get(), "element_formulation", "B33-3D-Euler-Bernoulli"); + requireStringAttribute(metadata.get(), "step_name", kStepName); + requireUint64Attribute(metadata.get(), "frame_index", 0U); + + requireCompoundDataset( + file.get(), + "/model/nodes", + static_cast(domain.nodes().size()), + {"internal_node_id", "instance_name", "source_label", "coordinates"}); + auto nodes = openDatasetForCheck(file.get(), "/model/nodes"); + requireStringAttribute(nodes.get(), "coordinate_system", "global-cartesian"); + requireStringAttribute(nodes.get(), "units_label", "length"); + requireCompoundDataset( + file.get(), + "/model/elements", + static_cast(domain.elements().size()), + {"internal_element_id", "instance_name", "source_label", + "node_internal_ids", "local_axes"}); + auto elements = openDatasetForCheck(file.get(), "/model/elements"); + requireStringAttribute( + elements.get(), "formulation", "B33-3D-Euler-Bernoulli"); + + const std::vector nodalDimensions = { + static_cast(domain.nodes().size()), kDofsPerNode}; + const std::vector endDimensions = { + static_cast(domain.elements().size()), + kEndpointCount, + kEndActionComponentCount}; + const std::vector generalizedDimensions = { + static_cast(domain.elements().size()), + kGaussPointCount, + kGeneralizedComponentCount}; + requireDoubleDataset( + file.get(), + "/steps/Step-1/frames/0/nodal/displacement", + nodalDimensions, + "UX,UY,UZ,URX,URY,URZ", + "length,length,length,radian,radian,radian", + "global-cartesian", + "nodal"); + requireDoubleDataset( + file.get(), + "/steps/Step-1/frames/0/nodal/reaction", + nodalDimensions, + "RF1,RF2,RF3,RM1,RM2,RM3", + "force,force,force,force*length,force*length,force*length", + "global-cartesian", + "nodal"); + requireDoubleDataset( + file.get(), + "/steps/Step-1/frames/0/element/end_force_local", + endDimensions, + "FX,FY,FZ,MX,MY,MZ", + "force,force,force,force*length,force*length,force*length", + "beam-local", + "endpoint-outward-action"); + requireDoubleDataset( + file.get(), + "/steps/Step-1/frames/0/element/section_resultant", + generalizedDimensions, + "N,T,My,Mz", + "force,force*length,force*length,force*length", + "beam-local", + "endpoint-positive-local-x-section-cut"); + requireDoubleDataset( + file.get(), + "/steps/Step-1/frames/0/element/generalized_strain", + generalizedDimensions, + "epsilon0,kappa_x,kappa_y,kappa_z", + "1,1/length,1/length,1/length", + "beam-local", + "integration-point"); + requireDoubleDataset( + file.get(), + "/steps/Step-1/frames/0/element/generalized_resultant", + generalizedDimensions, + "N,T,My,Mz", + "force,force*length,force*length,force*length", + "beam-local", + "integration-point"); + requireCompoundDataset( + file.get(), + "/steps/Step-1/frames/0/element/stress_s11", + static_cast(state.stressResults().size()), + {"internal_element_id", "gauss_point_index", "section_point_index", + "x1", "x2", "source", "S11"}); + auto stress = openDatasetForCheck( + file.get(), "/steps/Step-1/frames/0/element/stress_s11"); + requireResultAttributes( + stress.get(), "S11", "force/length^2", "beam-local", "section-point"); + requireCompoundDataset( + file.get(), + "/diagnostics", + static_cast(diagnosticCount), + {"severity", "code", "file", "line", "keyword", + "entity_identity", "message"}); + } + requireHdf5( + file.closeChecked(), + "Unable to close the read-only HDF5 schema self-check handle."); +} + +std::filesystem::path makeTemporaryPath( + const std::filesystem::path& outputPath) { + static std::atomic sequence{0U}; + const std::filesystem::path parent = outputPath.parent_path(); + for (std::size_t attempt = 0U; attempt < 1024U; ++attempt) { + const std::wstring filename = + L"." + outputPath.filename().wstring() + L".tmp." + + std::to_wstring(GetCurrentProcessId()) + L"." + + std::to_wstring(sequence.fetch_add(1U)); + const std::filesystem::path candidate = parent / filename; + std::error_code error; + const bool exists = std::filesystem::exists(candidate, error); + if (error) { + throw Hdf5Failure{"Unable to inspect the temporary output path."}; + } + if (!exists) { + return candidate; + } + } + throw Hdf5Failure{"Unable to allocate a unique temporary output path."}; +} + +bool finalizeFile( + const std::filesystem::path& temporaryPath, + const std::filesystem::path& outputPath) { + const DWORD attributes = GetFileAttributesW(outputPath.c_str()); + if (attributes != INVALID_FILE_ATTRIBUTES) { + // ReplaceFileW is the one-operation existing-final policy: failure + // leaves the old authoritative bytes in place. + return ReplaceFileW( + outputPath.c_str(), + temporaryPath.c_str(), + nullptr, + 0U, + nullptr, + nullptr) != FALSE; + } + const DWORD error = GetLastError(); + if (error != ERROR_FILE_NOT_FOUND && error != ERROR_PATH_NOT_FOUND) { + return false; + } + // A new final uses a same-directory write-through move after close/check. + return MoveFileExW( + temporaryPath.c_str(), + outputPath.c_str(), + MOVEFILE_WRITE_THROUGH) != FALSE; +} + +} // namespace + +Status Hdf5ResultsWriter::write( + const std::filesystem::path& outputPath, + const Domain& domain, + const AnalysisState& state, + const std::vector& diagnostics) { + std::vector localAxes; + try { + const Status validation = validateWriterInput( + outputPath, domain, state, diagnostics, localAxes); + if (!validation.isOk()) { + return validation; + } + + Hdf5ErrorSilencer silenceBackendErrors; + const std::filesystem::path temporaryPath = makeTemporaryPath(outputPath); + TemporaryFileGuard cleanup{temporaryPath}; + writeFile(temporaryPath, domain, state, diagnostics, localAxes); + // Close occurs when writeFile returns; only the reopened read-only file + // can approve this temp artifact for authoritative replacement. + selfCheckFile( + temporaryPath, domain, state, diagnostics.size()); + if (!finalizeFile(temporaryPath, outputPath)) { + return outputFailure( + "hdf5-finalization-failure", + "The checked temporary HDF5 file could not replace the final output."); + } + cleanup.release(); + return Status::ok(); + } catch (const Hdf5Failure& failure) { + return outputFailure("hdf5-write-failure", failure.what()); + } catch (const std::exception& failure) { + return outputFailure("hdf5-write-failure", failure.what()); + } +} + +} // namespace fesa diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 0e8b38b..7bfd53e 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -20,10 +20,12 @@ add_executable( unit/io/abaqus/domain_mapper_test.cpp unit/io/abaqus/input_reader_test.cpp unit/io/abaqus/input_syntax_test.cpp + unit/io/hdf5/hdf5_results_writer_test.cpp unit/model/domain_test.cpp unit/model/model_types_test.cpp unit/results/result_records_test.cpp unit/results/result_recovery_test.cpp + unit/results/results_writer_test.cpp unit/solvers/linear/linear_solver_test.cpp unit/solvers/linear/mkl_pardiso_solver_test.cpp ) @@ -32,6 +34,7 @@ target_link_libraries( fesa_unit_tests PRIVATE fesa_solver + Fesa::HDF5 GTest::gtest_main ) @@ -40,6 +43,9 @@ if(WIN32) NOT OMP_DLL_DIR OR NOT OMP_DLLNAME) message(FATAL_ERROR "The oneMKL runtime files required by tests were not resolved") endif() + if(NOT EXISTS "${OMP_DLL_DIR}/libmmd.dll") + message(FATAL_ERROR "The Intel runtime required by the selected HDF5 DLL was not resolved") + endif() file(GLOB _fesa_mkl_dispatch_dlls "${MKL_ROOT}/bin/mkl_def.*.dll") if(NOT _fesa_mkl_dispatch_dlls) @@ -65,9 +71,18 @@ if(WIN32) COMMAND "${CMAKE_COMMAND}" -E copy_if_different "${OMP_DLL_DIR}/${OMP_DLLNAME}" "$" + COMMAND "${CMAKE_COMMAND}" -E copy_if_different + "${OMP_DLL_DIR}/libmmd.dll" + "$" COMMAND "${CMAKE_COMMAND}" -E copy_if_different "${_fesa_mkl_dispatch_dll}" "$" + # Stage runtime DLLs selected through imported targets, including the + # normalized HDF5 backend, without naming a package-specific target. + COMMAND "${CMAKE_COMMAND}" -E copy_if_different + "$" + "$" + COMMAND_EXPAND_LISTS ) unset(_fesa_mkl_dispatch_dll) diff --git a/tests/unit/io/hdf5/hdf5_results_writer_test.cpp b/tests/unit/io/hdf5/hdf5_results_writer_test.cpp new file mode 100644 index 0000000..35a7f9a --- /dev/null +++ b/tests/unit/io/hdf5/hdf5_results_writer_test.cpp @@ -0,0 +1,1005 @@ +#define NOMINMAX +#include + +#include "fesa/io/hdf5/hdf5_results_writer.hpp" + +#include "fesa/analysis/analysis_model.hpp" +#include "fesa/analysis/analysis_state.hpp" +#include "fesa/build_info.hpp" +#include "fesa/fem/dof_manager.hpp" +#include "fesa/model/domain.hpp" + +#include + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace { + +constexpr const char* kStepRoot = "/steps/Step-1/frames/0"; + +class Hdf5Handle { +public: + using Closer = herr_t (*)(hid_t); + + Hdf5Handle() = default; + Hdf5Handle(const hid_t value, Closer closer) : value_{value}, closer_{closer} {} + Hdf5Handle(const Hdf5Handle&) = delete; + Hdf5Handle& operator=(const Hdf5Handle&) = delete; + Hdf5Handle(Hdf5Handle&& other) noexcept + : value_{other.value_}, closer_{other.closer_} { + other.value_ = -1; + other.closer_ = nullptr; + } + Hdf5Handle& operator=(Hdf5Handle&& other) noexcept { + if (this != &other) { + reset(); + value_ = other.value_; + closer_ = other.closer_; + other.value_ = -1; + other.closer_ = nullptr; + } + return *this; + } + ~Hdf5Handle() { reset(); } + + hid_t get() const noexcept { return value_; } + +private: + void reset() noexcept { + if (value_ >= 0 && closer_ != nullptr) { + (void)closer_(value_); + } + value_ = -1; + closer_ = nullptr; + } + + hid_t value_{-1}; + Closer closer_{nullptr}; +}; + +class WinHandle { +public: + explicit WinHandle(HANDLE value) : value_{value} {} + WinHandle(const WinHandle&) = delete; + WinHandle& operator=(const WinHandle&) = delete; + ~WinHandle() { + if (value_ != INVALID_HANDLE_VALUE) { + (void)CloseHandle(value_); + } + } + HANDLE get() const noexcept { return value_; } + +private: + HANDLE value_{INVALID_HANDLE_VALUE}; +}; + +class TempDirectory { +public: + explicit TempDirectory(const std::string& label) { + static std::atomic sequence{0U}; + path_ = std::filesystem::temp_directory_path() / + ("fesa-step23-" + label + "-" + + std::to_string(GetCurrentProcessId()) + "-" + + std::to_string(sequence.fetch_add(1U))); + std::error_code error; + if (!std::filesystem::create_directory(path_, error) || error) { + throw std::runtime_error{"Unable to create the Step 23 test directory."}; + } + } + TempDirectory(const TempDirectory&) = delete; + TempDirectory& operator=(const TempDirectory&) = delete; + ~TempDirectory() { + std::error_code ignored; + std::filesystem::remove_all(path_, ignored); + } + const std::filesystem::path& path() const noexcept { return path_; } + +private: + std::filesystem::path path_; +}; + +struct WriterFixture { + std::unique_ptr domain; + std::unique_ptr dofs; + std::unique_ptr state; +}; + +fesa::ModelDefinition makeDefinition( + const std::filesystem::path& source, + const bool useDefaultCentroid) { + fesa::ModelDefinition definition{}; + definition.sourcePath = source; + definition.sourceContentIdentity = "fnv1a64:0123456789abcdef"; + definition.nodes = { + {{u8"Beam-\u03b1", 101, "101"}, {0.0, 0.0, 0.0}, {source, 10U}}, + {{u8"Beam-\u03b1", 202, "202"}, {3.0, 4.0, 0.0}, {source, 11U}}}; + definition.materials = { + {"Steel", 210.0e9, 0.3, {source, 20U}}}; + definition.sections = {{ + "General", + 0.02, + 3.0e-5, + 0.0, + 4.0e-5, + 5.0e-5, + {0.0, 0.0, 1.0}, + useDefaultCentroid + ? std::vector>{} + : std::vector>{{{-0.1, 0.2}, {0.3, -0.4}}}, + {source, 30U}}}; + definition.elements = {{ + {u8"Beam-\u03b1", 303, "303"}, + {0U, 1U}, + 0U, + 0U, + {source, 40U}}}; + definition.steps = {{ + "Step-1", {}, {}, 0.1, 1.0, 0.01, 1.0, {source, 50U}}}; + return definition; +} + +WriterFixture makeFixture( + const std::filesystem::path& source, + const bool useDefaultCentroid = false) { + auto domainResult = fesa::Domain::create( + makeDefinition(source, useDefaultCentroid)); + if (!domainResult.hasValue()) { + throw std::runtime_error{"Writer fixture Domain construction failed."}; + } + auto domain = std::make_unique( + std::move(domainResult.value())); + + auto modelResult = fesa::AnalysisModel::create(*domain); + if (!modelResult.hasValue()) { + throw std::runtime_error{"Writer fixture AnalysisModel construction failed."}; + } + const fesa::AnalysisModel model = std::move(modelResult.value()); + auto dofsResult = fesa::DofManager::create(model); + if (!dofsResult.hasValue()) { + throw std::runtime_error{"Writer fixture DofManager construction failed."}; + } + auto dofs = std::make_unique( + std::move(dofsResult.value())); + auto state = std::make_unique( + fesa::AnalysisState::create(*dofs, {"Step-1", 0U})); + + for (std::size_t index = 0U; index < state->displacement().size(); ++index) { + state->displacement()[index] = 0.25 + static_cast(index); + state->externalForce()[index] = 100.0 + static_cast(index); + state->internalForce()[index] = 200.0 + 2.0 * static_cast(index); + state->residual()[index] = 100.0 + static_cast(index); + state->reaction()[index] = 100.0 + static_cast(index); + } + + const auto& nodes = domain->nodes(); + state->endpointResults() = { + {0U, + 0, + nodes[0U].sourceId, + {1.0, 2.0, 3.0, 4.0, 5.0, 6.0}, + {11.0, 12.0, 13.0, 14.0}}, + {0U, + 1, + nodes[1U].sourceId, + {7.0, 8.0, 9.0, 10.0, 11.0, 12.0}, + {15.0, 16.0, 17.0, 18.0}}}; + state->gaussResults() = { + {0U, 1, {0.01, 0.02, 0.03, 0.04}, {21.0, 22.0, 23.0, 24.0}}, + {0U, 2, {0.05, 0.06, 0.07, 0.08}, {25.0, 26.0, 27.0, 28.0}}}; + if (useDefaultCentroid) { + state->stressResults() = { + {0U, 1, 0U, 0.0, 0.0, 31.0, "fesa-default"}, + {0U, 2, 0U, 0.0, 0.0, 32.0, "fesa-default"}}; + } else { + state->stressResults() = { + {0U, 1, 1U, -0.1, 0.2, 31.0, "input"}, + {0U, 1, 2U, 0.3, -0.4, 32.0, "input"}, + {0U, 2, 1U, -0.1, 0.2, 33.0, "input"}, + {0U, 2, 2U, 0.3, -0.4, 34.0, "input"}}; + } + return {std::move(domain), std::move(dofs), std::move(state)}; +} + +Hdf5Handle openFile(const std::filesystem::path& path) { + const hid_t file = H5Fopen(path.string().c_str(), H5F_ACC_RDONLY, H5P_DEFAULT); + if (file < 0) { + throw std::runtime_error{"Unable to open test HDF5 output."}; + } + return Hdf5Handle{file, H5Fclose}; +} + +Hdf5Handle openDataset(const hid_t file, const std::string& path) { + const hid_t dataset = H5Dopen2(file, path.c_str(), H5P_DEFAULT); + if (dataset < 0) { + throw std::runtime_error{"Unable to open expected HDF5 dataset: " + path}; + } + return Hdf5Handle{dataset, H5Dclose}; +} + +std::vector datasetDimensions( + const hid_t file, const std::string& path) { + const auto dataset = openDataset(file, path); + Hdf5Handle space{H5Dget_space(dataset.get()), H5Sclose}; + if (space.get() < 0) { + throw std::runtime_error{"Unable to inspect HDF5 dataspace."}; + } + const int rank = H5Sget_simple_extent_ndims(space.get()); + if (rank < 0) { + throw std::runtime_error{"Unable to inspect HDF5 rank."}; + } + std::vector dimensions(static_cast(rank)); + if (rank > 0 && + H5Sget_simple_extent_dims(space.get(), dimensions.data(), nullptr) < 0) { + throw std::runtime_error{"Unable to inspect HDF5 dimensions."}; + } + return dimensions; +} + +std::vector readDoubleDataset( + const hid_t file, const std::string& path) { + const auto dimensions = datasetDimensions(file, path); + std::size_t valueCount = 1U; + for (const hsize_t dimension : dimensions) { + valueCount *= static_cast(dimension); + } + const auto dataset = openDataset(file, path); + std::vector values(valueCount); + if (!values.empty() && + H5Dread(dataset.get(), H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL, + H5P_DEFAULT, values.data()) < 0) { + throw std::runtime_error{"Unable to read numeric HDF5 dataset."}; + } + return values; +} + +std::string readStringAttribute(const hid_t object, const char* name) { + Hdf5Handle attribute{H5Aopen(object, name, H5P_DEFAULT), H5Aclose}; + Hdf5Handle type{H5Aget_type(attribute.get()), H5Tclose}; + if (attribute.get() < 0 || type.get() < 0 || + H5Tget_class(type.get()) != H5T_STRING || + H5Tis_variable_str(type.get()) <= 0 || + H5Tget_cset(type.get()) != H5T_CSET_UTF8) { + throw std::runtime_error{"Expected a variable-length UTF-8 attribute."}; + } + char* raw = nullptr; + if (H5Aread(attribute.get(), type.get(), &raw) < 0 || raw == nullptr) { + throw std::runtime_error{"Unable to read UTF-8 HDF5 attribute."}; + } + const std::string value{raw}; + (void)H5free_memory(raw); + return value; +} + +std::uint64_t readUint64Attribute(const hid_t object, const char* name) { + Hdf5Handle attribute{H5Aopen(object, name, H5P_DEFAULT), H5Aclose}; + Hdf5Handle type{H5Aget_type(attribute.get()), H5Tclose}; + if (attribute.get() < 0 || type.get() < 0 || + H5Tget_class(type.get()) != H5T_INTEGER || + H5Tget_size(type.get()) != sizeof(std::uint64_t) || + H5Tget_sign(type.get()) != H5T_SGN_NONE || + H5Tequal(type.get(), H5T_STD_U64LE) <= 0) { + throw std::runtime_error{"Expected a portable uint64 HDF5 attribute."}; + } + std::uint64_t value = 0U; + if (H5Aread(attribute.get(), H5T_NATIVE_UINT64, &value) < 0) { + throw std::runtime_error{"Unable to read uint64 HDF5 attribute."}; + } + return value; +} + +void expectPortableCompoundMember( + const hid_t compoundType, + const unsigned index, + const std::string& name) { + Hdf5Handle memberType{ + H5Tget_member_type(compoundType, index), H5Tclose}; + ASSERT_GE(memberType.get(), 0); + const bool isUint64 = + name == "internal_node_id" || name == "internal_element_id" || + name == "gauss_point_index" || name == "section_point_index" || + name == "line"; + const bool isFloat64 = name == "x1" || name == "x2" || name == "S11"; + const bool isString = + name == "instance_name" || name == "source_label" || + name == "source" || name == "severity" || name == "code" || + name == "file" || name == "keyword" || + name == "entity_identity" || name == "message"; + if (isUint64) { + EXPECT_EQ(H5Tget_class(memberType.get()), H5T_INTEGER); + EXPECT_EQ(H5Tget_size(memberType.get()), sizeof(std::uint64_t)); + EXPECT_EQ(H5Tget_sign(memberType.get()), H5T_SGN_NONE); + EXPECT_GT(H5Tequal(memberType.get(), H5T_STD_U64LE), 0); + return; + } + if (isFloat64) { + EXPECT_EQ(H5Tget_class(memberType.get()), H5T_FLOAT); + EXPECT_EQ(H5Tget_size(memberType.get()), sizeof(double)); + EXPECT_GT(H5Tequal(memberType.get(), H5T_IEEE_F64LE), 0); + return; + } + if (isString) { + EXPECT_EQ(H5Tget_class(memberType.get()), H5T_STRING); + EXPECT_GT(H5Tis_variable_str(memberType.get()), 0); + EXPECT_EQ(H5Tget_cset(memberType.get()), H5T_CSET_UTF8); + return; + } + + ASSERT_EQ(H5Tget_class(memberType.get()), H5T_ARRAY); + const int rank = H5Tget_array_ndims(memberType.get()); + ASSERT_GT(rank, 0); + std::vector dimensions(static_cast(rank)); + ASSERT_GE(H5Tget_array_dims2(memberType.get(), dimensions.data()), 0); + Hdf5Handle baseType{H5Tget_super(memberType.get()), H5Tclose}; + ASSERT_GE(baseType.get(), 0); + if (name == "node_internal_ids") { + EXPECT_EQ(dimensions, std::vector({2U})); + EXPECT_GT(H5Tequal(baseType.get(), H5T_STD_U64LE), 0); + } else if (name == "coordinates") { + EXPECT_EQ(dimensions, std::vector({3U})); + EXPECT_GT(H5Tequal(baseType.get(), H5T_IEEE_F64LE), 0); + } else { + EXPECT_EQ(name, "local_axes"); + EXPECT_EQ(dimensions, std::vector({3U, 3U})); + EXPECT_GT(H5Tequal(baseType.get(), H5T_IEEE_F64LE), 0); + } +} + +void expectCompoundMembers( + const hid_t file, + const std::string& path, + const std::vector& expectedNames) { + const auto dataset = openDataset(file, path); + Hdf5Handle type{H5Dget_type(dataset.get()), H5Tclose}; + ASSERT_EQ(H5Tget_class(type.get()), H5T_COMPOUND); + ASSERT_EQ( + H5Tget_nmembers(type.get()), static_cast(expectedNames.size())); + for (std::size_t index = 0U; index < expectedNames.size(); ++index) { + char* rawName = H5Tget_member_name(type.get(), + static_cast(index)); + ASSERT_NE(rawName, nullptr); + const std::string actualName{rawName}; + (void)H5free_memory(rawName); + EXPECT_EQ(actualName, expectedNames[index]); + expectPortableCompoundMember( + type.get(), static_cast(index), expectedNames[index]); + } +} + +void expectNumericDataset( + const hid_t file, + const std::string& path, + const std::vector& dimensions, + const std::string& components, + const std::string& units, + const std::string& coordinateSystem, + const std::string& location) { + EXPECT_EQ(datasetDimensions(file, path), dimensions); + const auto dataset = openDataset(file, path); + Hdf5Handle type{H5Dget_type(dataset.get()), H5Tclose}; + ASSERT_EQ(H5Tget_class(type.get()), H5T_FLOAT); + EXPECT_EQ(H5Tget_size(type.get()), 8U); + EXPECT_GT(H5Tequal(type.get(), H5T_IEEE_F64LE), 0); + EXPECT_EQ(readStringAttribute(dataset.get(), "component_names"), components); + EXPECT_EQ( + readStringAttribute(dataset.get(), "component_unit_dimensions"), units); + EXPECT_EQ( + readStringAttribute(dataset.get(), "coordinate_system"), coordinateSystem); + EXPECT_EQ(readStringAttribute(dataset.get(), "location"), location); + EXPECT_EQ(readStringAttribute(dataset.get(), "step_name"), "Step-1"); + EXPECT_EQ(readUint64Attribute(dataset.get(), "frame_index"), 0U); +} + +Hdf5Handle makeUtf8StringType() { + Hdf5Handle type{H5Tcopy(H5T_C_S1), H5Tclose}; + if (type.get() < 0 || H5Tset_size(type.get(), H5T_VARIABLE) < 0 || + H5Tset_cset(type.get(), H5T_CSET_UTF8) < 0) { + throw std::runtime_error{"Unable to create a test UTF-8 memory type."}; + } + return type; +} + +struct NodeReadRow { + std::uint64_t internalNodeId; + char* instanceName; + char* sourceLabel; + double coordinates[3]; +}; + +std::vector readNodeRows(const hid_t file) { + const auto dataset = openDataset(file, "/model/nodes"); + Hdf5Handle space{H5Dget_space(dataset.get()), H5Sclose}; + auto stringType = makeUtf8StringType(); + const hsize_t coordinateDimensions[] = {3U}; + Hdf5Handle coordinatesType{ + H5Tarray_create2(H5T_NATIVE_DOUBLE, 1, coordinateDimensions), H5Tclose}; + Hdf5Handle memoryType{H5Tcreate(H5T_COMPOUND, sizeof(NodeReadRow)), H5Tclose}; + if (H5Tinsert(memoryType.get(), "internal_node_id", + HOFFSET(NodeReadRow, internalNodeId), H5T_NATIVE_UINT64) < 0 || + H5Tinsert(memoryType.get(), "instance_name", + HOFFSET(NodeReadRow, instanceName), stringType.get()) < 0 || + H5Tinsert(memoryType.get(), "source_label", + HOFFSET(NodeReadRow, sourceLabel), stringType.get()) < 0 || + H5Tinsert(memoryType.get(), "coordinates", + HOFFSET(NodeReadRow, coordinates), coordinatesType.get()) < 0) { + throw std::runtime_error{"Unable to create the node memory type."}; + } + std::vector rows(2U); + if (H5Dread(dataset.get(), memoryType.get(), H5S_ALL, H5S_ALL, + H5P_DEFAULT, rows.data()) < 0) { + throw std::runtime_error{"Unable to read node rows."}; + } + return rows; +} + +void reclaimNodeRows(const hid_t file, std::vector& rows) { + const auto dataset = openDataset(file, "/model/nodes"); + Hdf5Handle space{H5Dget_space(dataset.get()), H5Sclose}; + auto stringType = makeUtf8StringType(); + const hsize_t coordinateDimensions[] = {3U}; + Hdf5Handle coordinatesType{ + H5Tarray_create2(H5T_NATIVE_DOUBLE, 1, coordinateDimensions), H5Tclose}; + Hdf5Handle memoryType{H5Tcreate(H5T_COMPOUND, sizeof(NodeReadRow)), H5Tclose}; + (void)H5Tinsert(memoryType.get(), "internal_node_id", + HOFFSET(NodeReadRow, internalNodeId), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "instance_name", + HOFFSET(NodeReadRow, instanceName), stringType.get()); + (void)H5Tinsert(memoryType.get(), "source_label", + HOFFSET(NodeReadRow, sourceLabel), stringType.get()); + (void)H5Tinsert(memoryType.get(), "coordinates", + HOFFSET(NodeReadRow, coordinates), coordinatesType.get()); + (void)H5Dvlen_reclaim( + memoryType.get(), space.get(), H5P_DEFAULT, rows.data()); +} + +struct ElementReadRow { + std::uint64_t internalElementId; + char* instanceName; + char* sourceLabel; + std::uint64_t nodeInternalIds[2]; + double localAxes[9]; +}; + +std::vector readElementRows(const hid_t file) { + const auto dataset = openDataset(file, "/model/elements"); + auto stringType = makeUtf8StringType(); + const hsize_t nodeDimensions[] = {2U}; + const hsize_t axesDimensions[] = {3U, 3U}; + Hdf5Handle nodeType{ + H5Tarray_create2(H5T_NATIVE_UINT64, 1, nodeDimensions), H5Tclose}; + Hdf5Handle axesType{ + H5Tarray_create2(H5T_NATIVE_DOUBLE, 2, axesDimensions), H5Tclose}; + Hdf5Handle memoryType{ + H5Tcreate(H5T_COMPOUND, sizeof(ElementReadRow)), H5Tclose}; + (void)H5Tinsert(memoryType.get(), "internal_element_id", + HOFFSET(ElementReadRow, internalElementId), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "instance_name", + HOFFSET(ElementReadRow, instanceName), stringType.get()); + (void)H5Tinsert(memoryType.get(), "source_label", + HOFFSET(ElementReadRow, sourceLabel), stringType.get()); + (void)H5Tinsert(memoryType.get(), "node_internal_ids", + HOFFSET(ElementReadRow, nodeInternalIds), nodeType.get()); + (void)H5Tinsert(memoryType.get(), "local_axes", + HOFFSET(ElementReadRow, localAxes), axesType.get()); + std::vector rows(1U); + if (H5Dread(dataset.get(), memoryType.get(), H5S_ALL, H5S_ALL, + H5P_DEFAULT, rows.data()) < 0) { + throw std::runtime_error{"Unable to read element rows."}; + } + Hdf5Handle space{H5Dget_space(dataset.get()), H5Sclose}; + EXPECT_EQ(rows[0U].internalElementId, 0U); + EXPECT_STREQ(rows[0U].instanceName, u8"Beam-\u03b1"); + EXPECT_STREQ(rows[0U].sourceLabel, "303"); + EXPECT_EQ(rows[0U].nodeInternalIds[0U], 0U); + EXPECT_EQ(rows[0U].nodeInternalIds[1U], 1U); + const std::array expectedAxes = { + 0.6, 0.8, 0.0, + 0.0, 0.0, 1.0, + 0.8, -0.6, 0.0}; + for (std::size_t index = 0U; index < expectedAxes.size(); ++index) { + EXPECT_NEAR(rows[0U].localAxes[index], expectedAxes[index], 1.0e-15); + } + (void)H5Dvlen_reclaim( + memoryType.get(), space.get(), H5P_DEFAULT, rows.data()); + return rows; +} + +struct StressReadRow { + std::uint64_t internalElementId; + std::uint64_t gaussPointIndex; + std::uint64_t sectionPointIndex; + double x1; + double x2; + char* source; + double s11; +}; + +std::vector readStressRows(const hid_t file) { + const std::string path = std::string{kStepRoot} + "/element/stress_s11"; + const auto dataset = openDataset(file, path); + const auto dimensions = datasetDimensions(file, path); + auto stringType = makeUtf8StringType(); + Hdf5Handle memoryType{ + H5Tcreate(H5T_COMPOUND, sizeof(StressReadRow)), H5Tclose}; + (void)H5Tinsert(memoryType.get(), "internal_element_id", + HOFFSET(StressReadRow, internalElementId), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "gauss_point_index", + HOFFSET(StressReadRow, gaussPointIndex), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "section_point_index", + HOFFSET(StressReadRow, sectionPointIndex), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "x1", + HOFFSET(StressReadRow, x1), H5T_NATIVE_DOUBLE); + (void)H5Tinsert(memoryType.get(), "x2", + HOFFSET(StressReadRow, x2), H5T_NATIVE_DOUBLE); + (void)H5Tinsert(memoryType.get(), "source", + HOFFSET(StressReadRow, source), stringType.get()); + (void)H5Tinsert(memoryType.get(), "S11", + HOFFSET(StressReadRow, s11), H5T_NATIVE_DOUBLE); + std::vector rows( + dimensions.empty() ? 0U : static_cast(dimensions[0U])); + if (!rows.empty() && + H5Dread(dataset.get(), memoryType.get(), H5S_ALL, H5S_ALL, + H5P_DEFAULT, rows.data()) < 0) { + throw std::runtime_error{"Unable to read stress rows."}; + } + return rows; +} + +void reclaimStressRows(const hid_t file, std::vector& rows) { + const std::string path = std::string{kStepRoot} + "/element/stress_s11"; + const auto dataset = openDataset(file, path); + Hdf5Handle space{H5Dget_space(dataset.get()), H5Sclose}; + auto stringType = makeUtf8StringType(); + Hdf5Handle memoryType{ + H5Tcreate(H5T_COMPOUND, sizeof(StressReadRow)), H5Tclose}; + (void)H5Tinsert(memoryType.get(), "internal_element_id", + HOFFSET(StressReadRow, internalElementId), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "gauss_point_index", + HOFFSET(StressReadRow, gaussPointIndex), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "section_point_index", + HOFFSET(StressReadRow, sectionPointIndex), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "x1", HOFFSET(StressReadRow, x1), H5T_NATIVE_DOUBLE); + (void)H5Tinsert(memoryType.get(), "x2", HOFFSET(StressReadRow, x2), H5T_NATIVE_DOUBLE); + (void)H5Tinsert(memoryType.get(), "source", + HOFFSET(StressReadRow, source), stringType.get()); + (void)H5Tinsert(memoryType.get(), "S11", HOFFSET(StressReadRow, s11), H5T_NATIVE_DOUBLE); + if (!rows.empty()) { + (void)H5Dvlen_reclaim( + memoryType.get(), space.get(), H5P_DEFAULT, rows.data()); + } +} + +struct DiagnosticReadRow { + char* severity; + char* code; + char* file; + std::uint64_t line; + char* keyword; + char* entityIdentity; + char* message; +}; + +std::vector readDiagnosticRows(const hid_t file) { + const auto dataset = openDataset(file, "/diagnostics"); + const auto dimensions = datasetDimensions(file, "/diagnostics"); + auto stringType = makeUtf8StringType(); + Hdf5Handle memoryType{ + H5Tcreate(H5T_COMPOUND, sizeof(DiagnosticReadRow)), H5Tclose}; + (void)H5Tinsert(memoryType.get(), "severity", + HOFFSET(DiagnosticReadRow, severity), stringType.get()); + (void)H5Tinsert(memoryType.get(), "code", + HOFFSET(DiagnosticReadRow, code), stringType.get()); + (void)H5Tinsert(memoryType.get(), "file", + HOFFSET(DiagnosticReadRow, file), stringType.get()); + (void)H5Tinsert(memoryType.get(), "line", + HOFFSET(DiagnosticReadRow, line), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "keyword", + HOFFSET(DiagnosticReadRow, keyword), stringType.get()); + (void)H5Tinsert(memoryType.get(), "entity_identity", + HOFFSET(DiagnosticReadRow, entityIdentity), stringType.get()); + (void)H5Tinsert(memoryType.get(), "message", + HOFFSET(DiagnosticReadRow, message), stringType.get()); + std::vector rows( + dimensions.empty() ? 0U : static_cast(dimensions[0U])); + if (!rows.empty() && + H5Dread(dataset.get(), memoryType.get(), H5S_ALL, H5S_ALL, + H5P_DEFAULT, rows.data()) < 0) { + throw std::runtime_error{"Unable to read diagnostic rows."}; + } + return rows; +} + +void reclaimDiagnosticRows( + const hid_t file, std::vector& rows) { + const auto dataset = openDataset(file, "/diagnostics"); + Hdf5Handle space{H5Dget_space(dataset.get()), H5Sclose}; + auto stringType = makeUtf8StringType(); + Hdf5Handle memoryType{ + H5Tcreate(H5T_COMPOUND, sizeof(DiagnosticReadRow)), H5Tclose}; + (void)H5Tinsert(memoryType.get(), "severity", + HOFFSET(DiagnosticReadRow, severity), stringType.get()); + (void)H5Tinsert(memoryType.get(), "code", + HOFFSET(DiagnosticReadRow, code), stringType.get()); + (void)H5Tinsert(memoryType.get(), "file", + HOFFSET(DiagnosticReadRow, file), stringType.get()); + (void)H5Tinsert(memoryType.get(), "line", + HOFFSET(DiagnosticReadRow, line), H5T_NATIVE_UINT64); + (void)H5Tinsert(memoryType.get(), "keyword", + HOFFSET(DiagnosticReadRow, keyword), stringType.get()); + (void)H5Tinsert(memoryType.get(), "entity_identity", + HOFFSET(DiagnosticReadRow, entityIdentity), stringType.get()); + (void)H5Tinsert(memoryType.get(), "message", + HOFFSET(DiagnosticReadRow, message), stringType.get()); + if (!rows.empty()) { + (void)H5Dvlen_reclaim( + memoryType.get(), space.get(), H5P_DEFAULT, rows.data()); + } +} + +std::vector readBytes(const std::filesystem::path& path) { + std::ifstream input{path, std::ios::binary}; + return {std::istreambuf_iterator{input}, std::istreambuf_iterator{}}; +} + +void writeBytes(const std::filesystem::path& path, const std::vector& bytes) { + std::ofstream output{path, std::ios::binary | std::ios::trunc}; + output.write(bytes.data(), static_cast(bytes.size())); + if (!output) { + throw std::runtime_error{"Unable to write atomicity sentinel bytes."}; + } +} + +std::size_t entryCount(const std::filesystem::path& directory) { + return static_cast( + std::distance(std::filesystem::directory_iterator{directory}, + std::filesystem::directory_iterator{})); +} + +void expectOutputFailure( + const fesa::Status& status, const std::string& expectedCode) { + ASSERT_FALSE(status.isOk()); + EXPECT_EQ(status.failureCategory(), fesa::FailureCategory::output); + ASSERT_EQ(status.diagnostics().size(), 1U); + EXPECT_EQ(status.diagnostics()[0U].severity, fesa::Severity::error); + EXPECT_EQ(status.diagnostics()[0U].code, expectedCode); +} + +} // namespace + +TEST(Hdf5ResultsWriter, WritesExactSchemaShapesAttributesAndIdentity) { + TempDirectory directory{"schema"}; + const auto source = directory.path() / "model.inp"; + auto fixture = makeFixture(source); + const auto output = directory.path() / "results.h5"; + + fesa::Hdf5ResultsWriter writer; + ASSERT_TRUE(writer.write(output, *fixture.domain, *fixture.state, {}).isOk()); + ASSERT_GT(H5Fis_hdf5(output.string().c_str()), 0); + + const auto file = openFile(output); + for (const char* path : { + "/metadata", + "/model/nodes", + "/model/elements", + "/steps/Step-1/frames/0/nodal/displacement", + "/steps/Step-1/frames/0/nodal/reaction", + "/steps/Step-1/frames/0/element/end_force_local", + "/steps/Step-1/frames/0/element/section_resultant", + "/steps/Step-1/frames/0/element/generalized_strain", + "/steps/Step-1/frames/0/element/generalized_resultant", + "/steps/Step-1/frames/0/element/stress_s11", + "/diagnostics"}) { + EXPECT_GT(H5Lexists(file.get(), path, H5P_DEFAULT), 0) << path; + } + + Hdf5Handle metadata{ + H5Gopen2(file.get(), "/metadata", H5P_DEFAULT), H5Gclose}; + ASSERT_GE(metadata.get(), 0); + EXPECT_EQ(readUint64Attribute(metadata.get(), "schema_version"), 0U); + EXPECT_EQ( + readStringAttribute(metadata.get(), "feature_id"), + "linear-static-3d-euler-beam"); + EXPECT_EQ( + readStringAttribute(metadata.get(), "solver_version"), + std::string{fesa::solverVersion()}); + const std::string normalizedSource = + std::filesystem::absolute(source).lexically_normal().generic_u8string(); + EXPECT_EQ( + readStringAttribute(metadata.get(), "source_input_identity"), + "path=" + normalizedSource + + ";content_identity=fnv1a64:0123456789abcdef"); + EXPECT_EQ( + readStringAttribute(metadata.get(), "unit_system_label"), + "user-consistent-unspecified"); + EXPECT_EQ( + readStringAttribute(metadata.get(), "coordinate_convention"), + "global-cartesian; beam-local=(t,n1,t-cross-n1)"); + EXPECT_EQ( + readStringAttribute(metadata.get(), "element_formulation"), + "B33-3D-Euler-Bernoulli"); + EXPECT_EQ(readStringAttribute(metadata.get(), "step_name"), "Step-1"); + EXPECT_EQ(readUint64Attribute(metadata.get(), "frame_index"), 0U); + + EXPECT_EQ(datasetDimensions(file.get(), "/model/nodes"), + std::vector({2U})); + expectCompoundMembers( + file.get(), "/model/nodes", + {"internal_node_id", "instance_name", "source_label", "coordinates"}); + const auto nodesDataset = openDataset(file.get(), "/model/nodes"); + EXPECT_EQ( + readStringAttribute(nodesDataset.get(), "coordinate_system"), + "global-cartesian"); + EXPECT_EQ(readStringAttribute(nodesDataset.get(), "units_label"), "length"); + auto nodes = readNodeRows(file.get()); + ASSERT_EQ(nodes.size(), 2U); + EXPECT_EQ(nodes[0U].internalNodeId, 0U); + EXPECT_STREQ(nodes[0U].instanceName, u8"Beam-\u03b1"); + EXPECT_STREQ(nodes[0U].sourceLabel, "101"); + EXPECT_DOUBLE_EQ(nodes[1U].coordinates[0U], 3.0); + EXPECT_DOUBLE_EQ(nodes[1U].coordinates[1U], 4.0); + reclaimNodeRows(file.get(), nodes); + + EXPECT_EQ(datasetDimensions(file.get(), "/model/elements"), + std::vector({1U})); + expectCompoundMembers( + file.get(), "/model/elements", + {"internal_element_id", "instance_name", "source_label", + "node_internal_ids", "local_axes"}); + const auto elementsDataset = openDataset(file.get(), "/model/elements"); + EXPECT_EQ( + readStringAttribute(elementsDataset.get(), "formulation"), + "B33-3D-Euler-Bernoulli"); + (void)readElementRows(file.get()); + + expectNumericDataset( + file.get(), std::string{kStepRoot} + "/nodal/displacement", {2U, 6U}, + "UX,UY,UZ,URX,URY,URZ", + "length,length,length,radian,radian,radian", + "global-cartesian", "nodal"); + const auto displacement = readDoubleDataset( + file.get(), std::string{kStepRoot} + "/nodal/displacement"); + ASSERT_EQ(displacement.size(), 12U); + EXPECT_DOUBLE_EQ(displacement.front(), 0.25); + EXPECT_DOUBLE_EQ(displacement.back(), 11.25); + expectNumericDataset( + file.get(), std::string{kStepRoot} + "/nodal/reaction", {2U, 6U}, + "RF1,RF2,RF3,RM1,RM2,RM3", + "force,force,force,force*length,force*length,force*length", + "global-cartesian", "nodal"); + const auto reaction = readDoubleDataset( + file.get(), std::string{kStepRoot} + "/nodal/reaction"); + ASSERT_EQ(reaction.size(), 12U); + EXPECT_DOUBLE_EQ(reaction.front(), 100.0); + EXPECT_DOUBLE_EQ(reaction.back(), 111.0); + expectNumericDataset( + file.get(), std::string{kStepRoot} + "/element/end_force_local", + {1U, 2U, 6U}, "FX,FY,FZ,MX,MY,MZ", + "force,force,force,force*length,force*length,force*length", + "beam-local", "endpoint-outward-action"); + const auto endForce = readDoubleDataset( + file.get(), std::string{kStepRoot} + "/element/end_force_local"); + ASSERT_EQ(endForce.size(), 12U); + EXPECT_DOUBLE_EQ(endForce.front(), 1.0); + EXPECT_DOUBLE_EQ(endForce.back(), 12.0); + expectNumericDataset( + file.get(), std::string{kStepRoot} + "/element/section_resultant", + {1U, 2U, 4U}, "N,T,My,Mz", + "force,force*length,force*length,force*length", + "beam-local", "endpoint-positive-local-x-section-cut"); + const auto sectionResultant = readDoubleDataset( + file.get(), std::string{kStepRoot} + "/element/section_resultant"); + ASSERT_EQ(sectionResultant.size(), 8U); + EXPECT_DOUBLE_EQ(sectionResultant.front(), 11.0); + EXPECT_DOUBLE_EQ(sectionResultant.back(), 18.0); + expectNumericDataset( + file.get(), std::string{kStepRoot} + "/element/generalized_strain", + {1U, 2U, 4U}, "epsilon0,kappa_x,kappa_y,kappa_z", + "1,1/length,1/length,1/length", + "beam-local", "integration-point"); + const auto generalizedStrain = readDoubleDataset( + file.get(), std::string{kStepRoot} + "/element/generalized_strain"); + ASSERT_EQ(generalizedStrain.size(), 8U); + EXPECT_DOUBLE_EQ(generalizedStrain.front(), 0.01); + EXPECT_DOUBLE_EQ(generalizedStrain.back(), 0.08); + expectNumericDataset( + file.get(), std::string{kStepRoot} + "/element/generalized_resultant", + {1U, 2U, 4U}, "N,T,My,Mz", + "force,force*length,force*length,force*length", + "beam-local", "integration-point"); + const auto generalizedResultant = readDoubleDataset( + file.get(), std::string{kStepRoot} + "/element/generalized_resultant"); + ASSERT_EQ(generalizedResultant.size(), 8U); + EXPECT_DOUBLE_EQ(generalizedResultant.front(), 21.0); + EXPECT_DOUBLE_EQ(generalizedResultant.back(), 28.0); + + const std::string stressPath = std::string{kStepRoot} + "/element/stress_s11"; + EXPECT_EQ(datasetDimensions(file.get(), stressPath), + std::vector({4U})); + expectCompoundMembers( + file.get(), stressPath, + {"internal_element_id", "gauss_point_index", "section_point_index", + "x1", "x2", "source", "S11"}); + const auto stressDataset = openDataset(file.get(), stressPath); + EXPECT_EQ(readStringAttribute(stressDataset.get(), "component_names"), "S11"); + EXPECT_EQ( + readStringAttribute(stressDataset.get(), "component_unit_dimensions"), + "force/length^2"); + EXPECT_EQ( + readStringAttribute(stressDataset.get(), "coordinate_system"), + "beam-local"); + EXPECT_EQ(readStringAttribute(stressDataset.get(), "location"), "section-point"); + EXPECT_EQ(readStringAttribute(stressDataset.get(), "step_name"), "Step-1"); + EXPECT_EQ(readUint64Attribute(stressDataset.get(), "frame_index"), 0U); + auto stressRows = readStressRows(file.get()); + ASSERT_EQ(stressRows.size(), 4U); + EXPECT_EQ(stressRows[0U].internalElementId, 0U); + EXPECT_EQ(stressRows[0U].gaussPointIndex, 1U); + EXPECT_EQ(stressRows[0U].sectionPointIndex, 1U); + EXPECT_DOUBLE_EQ(stressRows[0U].x1, -0.1); + EXPECT_DOUBLE_EQ(stressRows[0U].x2, 0.2); + EXPECT_STREQ(stressRows[0U].source, "input"); + EXPECT_DOUBLE_EQ(stressRows[3U].s11, 34.0); + reclaimStressRows(file.get(), stressRows); + + EXPECT_EQ(datasetDimensions(file.get(), "/diagnostics"), + std::vector({0U})); + expectCompoundMembers( + file.get(), "/diagnostics", + {"severity", "code", "file", "line", "keyword", + "entity_identity", "message"}); + EXPECT_EQ( + H5Lexists(file.get(), + "/steps/Step-1/frames/0/element/transverse_shear_stress", + H5P_DEFAULT), + 0); +} + +TEST(Hdf5ResultsWriter, WritesMandatoryOutputsDespiteOutputRequests) { + TempDirectory directory{"mandatory"}; + auto fixture = makeFixture(directory.path() / "request-model.inp"); + const fesa::Diagnostic ignoredRequest{ + fesa::Severity::warning, + "ignored-output-request", + {fixture.domain->sourcePath(), 70U}, + "*OUTPUT", + "FIELD", + "Abaqus output requests do not filter FESA mandatory results."}; + const auto output = directory.path() / "results.h5"; + + fesa::Hdf5ResultsWriter writer; + ASSERT_TRUE( + writer.write(output, *fixture.domain, *fixture.state, {ignoredRequest}) + .isOk()); + const auto file = openFile(output); + for (const char* suffix : { + "/nodal/displacement", + "/nodal/reaction", + "/element/end_force_local", + "/element/section_resultant", + "/element/generalized_strain", + "/element/generalized_resultant", + "/element/stress_s11"}) { + const std::string path = std::string{kStepRoot} + suffix; + EXPECT_GT(H5Lexists(file.get(), path.c_str(), H5P_DEFAULT), 0) << path; + } + EXPECT_EQ(datasetDimensions(file.get(), "/diagnostics"), + std::vector({1U})); +} + +TEST(Hdf5ResultsWriter, WritesWarningsAndDefaultCentroid) { + TempDirectory directory{"warnings"}; + auto fixture = makeFixture(directory.path() / "centroid.inp", true); + std::vector diagnostics = { + {fesa::Severity::warning, + "ignored-output-request", + {fixture.domain->sourcePath(), 80U}, + "*OUTPUT", + "FIELD", + "Ignored output request."}, + {fesa::Severity::warning, + "ignored-keyword", + {fixture.domain->sourcePath(), 20U}, + "*PREPRINT", + "", + "Ignored generator control."}}; + const auto output = directory.path() / "results.h5"; + + fesa::Hdf5ResultsWriter writer; + ASSERT_TRUE(writer.write(output, *fixture.domain, *fixture.state, diagnostics).isOk()); + const auto file = openFile(output); + auto stressRows = readStressRows(file.get()); + ASSERT_EQ(stressRows.size(), 2U); + for (std::size_t index = 0U; index < stressRows.size(); ++index) { + EXPECT_EQ(stressRows[index].internalElementId, 0U); + EXPECT_EQ(stressRows[index].gaussPointIndex, index + 1U); + EXPECT_EQ(stressRows[index].sectionPointIndex, 0U); + EXPECT_DOUBLE_EQ(stressRows[index].x1, 0.0); + EXPECT_DOUBLE_EQ(stressRows[index].x2, 0.0); + EXPECT_STREQ(stressRows[index].source, "fesa-default"); + } + reclaimStressRows(file.get(), stressRows); + + auto rows = readDiagnosticRows(file.get()); + ASSERT_EQ(rows.size(), 2U); + EXPECT_STREQ(rows[0U].severity, "warning"); + EXPECT_STREQ(rows[0U].code, "ignored-keyword"); + EXPECT_STREQ( + rows[0U].file, + std::filesystem::absolute(fixture.domain->sourcePath()) + .lexically_normal() + .generic_u8string() + .c_str()); + EXPECT_EQ(rows[0U].line, 20U); + EXPECT_STREQ(rows[0U].keyword, "*PREPRINT"); + EXPECT_STREQ(rows[0U].entityIdentity, ""); + EXPECT_STREQ(rows[0U].message, "Ignored generator control."); + EXPECT_STREQ(rows[1U].code, "ignored-output-request"); + EXPECT_EQ(rows[1U].line, 80U); + reclaimDiagnosticRows(file.get(), rows); +} + +TEST(Hdf5ResultsWriter, FailureLeavesNoPartialAndPreservesExistingFinal) { + TempDirectory directory{"failure"}; + auto fixture = makeFixture(directory.path() / "failure.inp"); + fesa::Hdf5ResultsWriter writer; + + fixture.state->displacement()[0U] = + std::numeric_limits::quiet_NaN(); + const auto invalidOutput = directory.path() / "invalid-results.h5"; + expectOutputFailure( + writer.write(invalidOutput, *fixture.domain, *fixture.state, {}), + "invalid-result-state"); + EXPECT_FALSE(std::filesystem::exists(invalidOutput)); + EXPECT_EQ(entryCount(directory.path()), 0U); + fixture.state->displacement()[0U] = 0.25; + + const auto final = directory.path() / "results.h5"; + const std::vector sentinel = {'p', 'r', 'e', 'v', 'i', 'o', 'u', 's'}; + writeBytes(final, sentinel); + WinHandle lock{CreateFileW( + final.c_str(), + GENERIC_READ, + FILE_SHARE_READ | FILE_SHARE_WRITE, + nullptr, + OPEN_EXISTING, + FILE_ATTRIBUTE_NORMAL, + nullptr)}; + ASSERT_NE(lock.get(), INVALID_HANDLE_VALUE); + + expectOutputFailure( + writer.write(final, *fixture.domain, *fixture.state, {}), + "hdf5-finalization-failure"); + EXPECT_EQ(readBytes(final), sentinel); + EXPECT_EQ(entryCount(directory.path()), 1U); +} + +TEST(Hdf5ResultsWriter, SuccessfullyReplacesExistingFinal) { + TempDirectory directory{"replace"}; + auto fixture = makeFixture(directory.path() / "replace.inp"); + const auto final = directory.path() / "results.h5"; + writeBytes(final, {'o', 'l', 'd'}); + + fesa::Hdf5ResultsWriter writer; + ASSERT_TRUE(writer.write(final, *fixture.domain, *fixture.state, {}).isOk()); + EXPECT_GT(H5Fis_hdf5(final.string().c_str()), 0); + EXPECT_EQ(entryCount(directory.path()), 1U); + const auto file = openFile(final); + Hdf5Handle metadata{ + H5Gopen2(file.get(), "/metadata", H5P_DEFAULT), H5Gclose}; + ASSERT_GE(metadata.get(), 0); + EXPECT_EQ(readUint64Attribute(metadata.get(), "schema_version"), 0U); +} diff --git a/tests/unit/results/results_writer_test.cpp b/tests/unit/results/results_writer_test.cpp new file mode 100644 index 0000000..4ced0b2 --- /dev/null +++ b/tests/unit/results/results_writer_test.cpp @@ -0,0 +1,24 @@ +#include "fesa/results/results_writer.hpp" + +#include "fesa/analysis/analysis_state.hpp" +#include "fesa/core/diagnostic.hpp" +#include "fesa/core/status.hpp" +#include "fesa/model/domain.hpp" + +#include +#include +#include + +namespace { + +using WriteSignature = fesa::Status (fesa::ResultsWriter::*)( + const std::filesystem::path&, + const fesa::Domain&, + const fesa::AnalysisState&, + const std::vector&); + +static_assert(std::has_virtual_destructor_v); +static_assert(std::is_abstract_v); +static_assert(std::is_same_v); + +} // namespace