From 76a0c8ecde5488bc0fe2904668396a85d2549324 Mon Sep 17 00:00:00 2001 From: "KOKO\\Mimi" Date: Sat, 1 Aug 2026 00:51:34 +0900 Subject: [PATCH] =?UTF-8?q?feat(results-and-pipeline):=20step=201=20?= =?UTF-8?q?=E2=80=94=20minimal-hdf5-schema?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- CMakeLists.txt | 3 +- docs/HDF5_SCHEMA.md | 90 ++ include/fesa/io/hdf5/writer.hpp | 60 + src/fesa/io/hdf5/writer.cpp | 1582 ++++++++++++++++++++ tests/CMakeLists.txt | 43 + tests/integration/io/hdf5_results_test.cpp | 566 +++++++ 6 files changed, 2343 insertions(+), 1 deletion(-) create mode 100644 docs/HDF5_SCHEMA.md create mode 100644 include/fesa/io/hdf5/writer.hpp create mode 100644 src/fesa/io/hdf5/writer.cpp create mode 100644 tests/integration/io/hdf5_results_test.cpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 2ba8079..dbd7088 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -36,6 +36,7 @@ add_library(fesa_core STATIC src/fesa/fem/line2_shape.cpp src/fesa/io/abaqus/parser.cpp src/fesa/io/abaqus/semantic_mapper.cpp + src/fesa/io/hdf5/writer.cpp src/fesa/model/domain.cpp src/fesa/model/domain_builder.cpp src/fesa/results/result_database.cpp @@ -49,7 +50,7 @@ target_include_directories(fesa_core target_compile_features(fesa_core PUBLIC cxx_std_20) target_compile_options(fesa_core PRIVATE /W4 /permissive- /EHsc) -target_link_libraries(fesa_core PRIVATE MKL::MKL) +target_link_libraries(fesa_core PRIVATE MKL::MKL HDF5::HDF5) add_executable(fesa src/fesa/cli/main.cpp diff --git a/docs/HDF5_SCHEMA.md b/docs/HDF5_SCHEMA.md new file mode 100644 index 0000000..90e2c6d --- /dev/null +++ b/docs/HDF5_SCHEMA.md @@ -0,0 +1,90 @@ +# FESA HDF5 Schema 1.0.0 + +## 1. 범위 + +Schema `1.0.0`은 `results-and-pipeline` Phase의 최소 수직 슬라이스를 정의한다. +파일은 활성 `Domain`의 절점, Beam 연결성, 적용된 전단면적과 그 출처, 그리고 전역 +좌표계 절점 변위·회전 및 반력·반력모멘트를 저장한다. 단위 변환은 수행하지 않는다. + +이 버전에는 재료 전체 속성, 집합, 하중·경계조건, solver 설정, 요소 결과, history, +reference CSV 및 진단 dataset을 저장하지 않는다. 이후 같은 major version에서 +dataset을 추가할 수 있지만 아래 required object의 의미, 형상 또는 datatype을 +변경해서는 안 된다. + +## 2. 공통 규칙 + +- root attribute `schema_version`은 UTF-8 문자열 `1.0.0`이다. +- 정수 dataset은 명시한 little-endian 고정폭 타입을 사용한다. +- 실수 dataset은 IEEE 754 little-endian 64-bit 타입을 사용한다. +- 문자열 dataset과 attribute는 UTF-8 variable-length string을 사용한다. +- `dense_index`는 해당 dataset 행의 0-based index이며 연속적이다. +- `internal_id`와 결과의 `node_ids`는 FESA semantic model의 nonnegative ID다. +- flat/orphan mesh의 `part_name`과 `instance_name`은 빈 문자열이다. +- 결과의 6개 component 순서는 + `(Ux, Uy, Uz, Rx, Ry, Rz)` 및 `(RFx, RFy, RFz, RMx, RMy, RMz)`다. +- Step과 frame group 이름은 각각 0부터 연속된 decimal index다. 원래 Step 이름은 + Step group의 `name` attribute에 저장한다. + +## 3. Required objects + +```text +/ +├── @schema_version UTF-8 = "1.0.0" +├── model +│ ├── nodes +│ │ ├── dense_index uint64 [node_count] +│ │ ├── internal_id int64 [node_count] +│ │ ├── part_name UTF-8 [node_count] +│ │ ├── instance_name UTF-8 [node_count] +│ │ ├── local_label int64 [node_count] +│ │ └── coordinates float64[node_count, 3] +│ ├── elements +│ │ ├── dense_index uint64 [element_count] +│ │ ├── internal_id int64 [element_count] +│ │ ├── connectivity uint64 [element_count, 2] +│ │ └── section_id int64 [element_count] +│ └── sections +│ ├── internal_id int64 [section_count] +│ ├── shear_area_y float64[section_count] +│ ├── shear_area_z float64[section_count] +│ └── shear_source uint8 [section_count] +└── results + └── steps + └── + ├── @name UTF-8 + └── frames + └── + ├── @step_time float64 + └── nodal + ├── node_ids int64 [result_node_count] + ├── displacement float64[result_node_count, 6] + └── reaction float64[result_node_count, 6] +``` + +`model/elements/connectivity`는 `model/nodes/dense_index`를 참조한다. 따라서 +`internal_id`가 연속적이거나 Domain 저장 순서와 같다고 가정하지 않는다. +`section_id`는 `model/sections/internal_id`를 참조한다. + +`shear_source` 값은 다음과 같다. + +| 값 | 의미 | +|---:|---| +| `0` | 입력에서 명시된 전단강성으로부터 구성한 값 (`input`) | +| `1` | Phase 1 기본값 `Asy=Asz=5A/6`, `SCF=0` (`phase1_default`) | + +## 4. Writer와 reader 계약 + +- writer는 쓰기 전에 `ResultDatabase` 유효성과 schema version을 검사한다. +- schema `1.0.0`이 표현하지 않는 non-empty frame diagnostics는 파일을 만들기 전에 + `hdf5.unsupported_result_diagnostics`로 거부한다. +- required object 생성·쓰기·flush·close 중 HDF5 오류가 발생하면 성공으로 반환하지 + 않고 `DiagnosticStage::results` 오류로 변환한다. +- reader는 schema version, required object, datatype, rank와 shape를 검사한다. +- reader는 model internal ID의 uniqueness, finite coordinates 및 finite positive + shear area를 검사한다. +- writer와 reader는 모든 nodal result ID가 `model/nodes/internal_id`에 존재하는지 + 검사하며, 없는 ID를 성공 결과로 반환하지 않는다. +- reader는 nodal result를 `ResultDatabase`로, model dataset을 HDF5 adapter 전용 + read-only inspection model로 반환한다. `Domain`과 `ResultDatabase`에는 HDF5 + 저장 계약을 추가하지 않는다. +- malformed 또는 지원하지 않는 파일은 부분 database를 반환하지 않는다. diff --git a/include/fesa/io/hdf5/writer.hpp b/include/fesa/io/hdf5/writer.hpp new file mode 100644 index 0000000..f7c0fd1 --- /dev/null +++ b/include/fesa/io/hdf5/writer.hpp @@ -0,0 +1,60 @@ +#pragma once + +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +namespace fesa { + +struct Hdf5NodeSnapshot final { + std::uint64_t dense_index; + NodeId id; + EntityOrigin origin; + Vec3 coordinates; +}; + +struct Hdf5ElementSnapshot final { + std::uint64_t dense_index; + ElementId id; + std::array connectivity; + SectionId section; +}; + +struct Hdf5SectionSnapshot final { + SectionId id; + double shear_area_y; + double shear_area_z; + ShearPropertySource shear_source; +}; + +struct Hdf5ModelSnapshot final { + std::vector nodes; + std::vector elements; + std::vector sections; +}; + +struct Hdf5ReadResult final { + std::optional database; + std::vector diagnostics; + std::optional model; +}; + +[[nodiscard]] std::vector write_hdf5( + const std::filesystem::path& path, + const Domain& domain, + const ResultDatabase& database); + +[[nodiscard]] Hdf5ReadResult read_hdf5_results( + const std::filesystem::path& path); + +} // namespace fesa diff --git a/src/fesa/io/hdf5/writer.cpp b/src/fesa/io/hdf5/writer.cpp new file mode 100644 index 0000000..bbeb485 --- /dev/null +++ b/src/fesa/io/hdf5/writer.cpp @@ -0,0 +1,1582 @@ +#include + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace fesa { +namespace { + +constexpr std::string_view schema_version = "1.0.0"; + +class Hdf5Error final : public std::runtime_error { +public: + Hdf5Error(std::string code, std::string message) + : std::runtime_error{std::move(message)}, code_{std::move(code)} {} + + [[nodiscard]] const std::string& code() const noexcept { + return code_; + } + +private: + std::string code_; +}; + +template +auto call_hdf5(Function&& function) { + using Return = std::invoke_result_t; + Return result{}; + H5E_BEGIN_TRY { + result = std::forward(function)(); + } + H5E_END_TRY; + return result; +} + +struct Hdf5Context final { + std::string close_error; + + void record_close_result( + const herr_t result, + const std::string_view resource) { + if (result < 0 && close_error.empty()) { + close_error = "Failed to release HDF5 " + + std::string{resource} + "."; + } + } +}; + +using CloseFunction = herr_t (*)(hid_t); + +class Hdf5Handle final { +public: + Hdf5Handle() = default; + + Hdf5Handle( + Hdf5Context& context, + const hid_t id, + CloseFunction close, + std::string resource) + : context_{&context}, + id_{id}, + close_{close}, + resource_{std::move(resource)} {} + + ~Hdf5Handle() { + reset(); + } + + Hdf5Handle(const Hdf5Handle&) = delete; + Hdf5Handle& operator=(const Hdf5Handle&) = delete; + + Hdf5Handle(Hdf5Handle&& other) noexcept + : context_{std::exchange(other.context_, nullptr)}, + id_{std::exchange(other.id_, H5I_INVALID_HID)}, + close_{std::exchange(other.close_, nullptr)}, + resource_{std::move(other.resource_)} {} + + Hdf5Handle& operator=(Hdf5Handle&& other) noexcept { + if (this != &other) { + reset(); + context_ = std::exchange(other.context_, nullptr); + id_ = std::exchange(other.id_, H5I_INVALID_HID); + close_ = std::exchange(other.close_, nullptr); + resource_ = std::move(other.resource_); + } + return *this; + } + + [[nodiscard]] hid_t get() const noexcept { + return id_; + } + + void reset() noexcept { + if (id_ >= 0 && close_ != nullptr && context_ != nullptr) { + const herr_t result = call_hdf5([this] { return close_(id_); }); + context_->record_close_result(result, resource_); + } + id_ = H5I_INVALID_HID; + close_ = nullptr; + context_ = nullptr; + } + +private: + Hdf5Context* context_{nullptr}; + hid_t id_{H5I_INVALID_HID}; + CloseFunction close_{nullptr}; + std::string resource_; +}; + +class Hdf5AllocatedString final { +public: + Hdf5AllocatedString(Hdf5Context& context, char* value) + : context_{context}, value_{value} {} + + ~Hdf5AllocatedString() { + if (value_ != nullptr) { + context_.record_close_result( + call_hdf5([this] { return H5free_memory(value_); }), + "variable-length string"); + } + } + + Hdf5AllocatedString(const Hdf5AllocatedString&) = delete; + Hdf5AllocatedString& operator=(const Hdf5AllocatedString&) = delete; + +private: + Hdf5Context& context_; + char* value_; +}; + +class Hdf5VlenBuffer final { +public: + Hdf5VlenBuffer( + Hdf5Context& context, + const hid_t type, + const hid_t space, + void* data) + : context_{context}, type_{type}, space_{space}, data_{data} {} + + ~Hdf5VlenBuffer() { + context_.record_close_result( + call_hdf5([this] { + return H5Treclaim(type_, space_, H5P_DEFAULT, data_); + }), + "variable-length dataset buffer"); + } + + Hdf5VlenBuffer(const Hdf5VlenBuffer&) = delete; + Hdf5VlenBuffer& operator=(const Hdf5VlenBuffer&) = delete; + +private: + Hdf5Context& context_; + hid_t type_; + hid_t space_; + void* data_; +}; + +[[noreturn]] void fail(std::string code, std::string message) { + throw Hdf5Error{std::move(code), std::move(message)}; +} + +void check_status( + const herr_t result, + const std::string_view operation, + const std::string_view code) { + if (result < 0) { + fail( + std::string{code}, + "HDF5 operation failed: " + std::string{operation} + "."); + } +} + +Hdf5Handle checked_handle( + Hdf5Context& context, + const hid_t id, + CloseFunction close, + std::string resource, + const std::string_view code) { + if (id < 0) { + fail( + std::string{code}, + "Failed to open or create HDF5 " + resource + "."); + } + return {context, id, close, std::move(resource)}; +} + +Hdf5Handle create_group( + Hdf5Context& context, + const hid_t parent, + const std::string_view name) { + const std::string owned_name{name}; + return checked_handle( + context, + call_hdf5([&] { + return H5Gcreate2( + parent, + owned_name.c_str(), + H5P_DEFAULT, + H5P_DEFAULT, + H5P_DEFAULT); + }), + &H5Gclose, + "group '" + owned_name + "'", + "hdf5.write_failed"); +} + +Hdf5Handle open_group( + Hdf5Context& context, + const hid_t parent, + const std::string_view name) { + const std::string owned_name{name}; + return checked_handle( + context, + call_hdf5([&] { + return H5Gopen2(parent, owned_name.c_str(), H5P_DEFAULT); + }), + &H5Gclose, + "group '" + owned_name + "'", + "hdf5.read_failed"); +} + +Hdf5Handle create_space( + Hdf5Context& context, + const std::span dimensions, + const std::string_view code) { + return checked_handle( + context, + call_hdf5([&] { + return H5Screate_simple( + static_cast(dimensions.size()), + dimensions.data(), + nullptr); + }), + &H5Sclose, + "dataspace", + code); +} + +Hdf5Handle create_scalar_space( + Hdf5Context& context, + const std::string_view code) { + return checked_handle( + context, + call_hdf5([] { return H5Screate(H5S_SCALAR); }), + &H5Sclose, + "scalar dataspace", + code); +} + +Hdf5Handle create_utf8_type( + Hdf5Context& context, + const std::string_view code) { + auto type = checked_handle( + context, + call_hdf5([] { return H5Tcopy(H5T_C_S1); }), + &H5Tclose, + "UTF-8 string datatype", + code); + check_status( + call_hdf5([&] { return H5Tset_size(type.get(), H5T_VARIABLE); }), + "set variable string size", + code); + check_status( + call_hdf5([&] { return H5Tset_cset(type.get(), H5T_CSET_UTF8); }), + "set UTF-8 character set", + code); + return type; +} + +std::string path_utf8(const std::filesystem::path& path) { + const std::u8string value = path.u8string(); + return {reinterpret_cast(value.data()), value.size()}; +} + +Diagnostic error_diagnostic(std::string code, std::string message) { + return { + DiagnosticStage::results, + Severity::error, + std::move(code), + std::move(message), + std::nullopt, + }; +} + +std::optional find_unknown_result_node( + const ResultDatabase& database, + const std::unordered_set& model_node_ids) { + for (const ResultStep& step : database.steps) { + for (const ResultFrame& frame : step.frames) { + for (const NodeId node_id : frame.nodal.node_ids) { + if (!model_node_ids.contains(node_id.value())) { + return node_id.value(); + } + } + } + } + return std::nullopt; +} + +void write_string_attribute( + Hdf5Context& context, + const hid_t parent, + const std::string_view name, + const std::string_view value) { + auto type = create_utf8_type(context, "hdf5.write_failed"); + auto space = create_scalar_space(context, "hdf5.write_failed"); + const std::string owned_name{name}; + auto attribute = checked_handle( + context, + call_hdf5([&] { + return H5Acreate2( + parent, + owned_name.c_str(), + type.get(), + space.get(), + H5P_DEFAULT, + H5P_DEFAULT); + }), + &H5Aclose, + "attribute '" + owned_name + "'", + "hdf5.write_failed"); + const std::string owned_value{value}; + const char* pointer = owned_value.c_str(); + check_status( + call_hdf5([&] { + return H5Awrite(attribute.get(), type.get(), &pointer); + }), + "write string attribute '" + owned_name + "'", + "hdf5.write_failed"); +} + +void write_double_attribute( + Hdf5Context& context, + const hid_t parent, + const std::string_view name, + const double value) { + auto space = create_scalar_space(context, "hdf5.write_failed"); + const std::string owned_name{name}; + auto attribute = checked_handle( + context, + call_hdf5([&] { + return H5Acreate2( + parent, + owned_name.c_str(), + H5T_IEEE_F64LE, + space.get(), + H5P_DEFAULT, + H5P_DEFAULT); + }), + &H5Aclose, + "attribute '" + owned_name + "'", + "hdf5.write_failed"); + check_status( + call_hdf5([&] { + return H5Awrite(attribute.get(), H5T_NATIVE_DOUBLE, &value); + }), + "write floating-point attribute '" + owned_name + "'", + "hdf5.write_failed"); +} + +template +void write_numeric_1d( + Hdf5Context& context, + const hid_t parent, + const std::string_view name, + const hid_t file_type, + const hid_t memory_type, + const std::span values) { + const std::array dimensions{ + static_cast(values.size())}; + auto space = create_space(context, dimensions, "hdf5.write_failed"); + const std::string owned_name{name}; + auto dataset = checked_handle( + context, + call_hdf5([&] { + return H5Dcreate2( + parent, + owned_name.c_str(), + file_type, + space.get(), + H5P_DEFAULT, + H5P_DEFAULT, + H5P_DEFAULT); + }), + &H5Dclose, + "dataset '" + owned_name + "'", + "hdf5.write_failed"); + if (!values.empty()) { + check_status( + call_hdf5([&] { + return H5Dwrite( + dataset.get(), + memory_type, + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + values.data()); + }), + "write dataset '" + owned_name + "'", + "hdf5.write_failed"); + } +} + +template +void write_numeric_2d( + Hdf5Context& context, + const hid_t parent, + const std::string_view name, + const hid_t file_type, + const hid_t memory_type, + const std::size_t row_count, + const std::size_t column_count, + const std::span values) { + if (row_count * column_count != values.size()) { + fail( + "hdf5.write_failed", + "Internal HDF5 dataset shape does not match its values."); + } + const std::array dimensions{ + static_cast(row_count), + static_cast(column_count), + }; + auto space = create_space(context, dimensions, "hdf5.write_failed"); + const std::string owned_name{name}; + auto dataset = checked_handle( + context, + call_hdf5([&] { + return H5Dcreate2( + parent, + owned_name.c_str(), + file_type, + space.get(), + H5P_DEFAULT, + H5P_DEFAULT, + H5P_DEFAULT); + }), + &H5Dclose, + "dataset '" + owned_name + "'", + "hdf5.write_failed"); + if (!values.empty()) { + check_status( + call_hdf5([&] { + return H5Dwrite( + dataset.get(), + memory_type, + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + values.data()); + }), + "write dataset '" + owned_name + "'", + "hdf5.write_failed"); + } +} + +void write_strings_1d( + Hdf5Context& context, + const hid_t parent, + const std::string_view name, + const std::span values) { + auto type = create_utf8_type(context, "hdf5.write_failed"); + const std::array dimensions{ + static_cast(values.size())}; + auto space = create_space(context, dimensions, "hdf5.write_failed"); + const std::string owned_name{name}; + auto dataset = checked_handle( + context, + call_hdf5([&] { + return H5Dcreate2( + parent, + owned_name.c_str(), + type.get(), + space.get(), + H5P_DEFAULT, + H5P_DEFAULT, + H5P_DEFAULT); + }), + &H5Dclose, + "dataset '" + owned_name + "'", + "hdf5.write_failed"); + if (!values.empty()) { + std::vector pointers; + pointers.reserve(values.size()); + for (const std::string& value : values) { + pointers.push_back(value.c_str()); + } + check_status( + call_hdf5([&] { + return H5Dwrite( + dataset.get(), + type.get(), + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + pointers.data()); + }), + "write string dataset '" + owned_name + "'", + "hdf5.write_failed"); + } +} + +void write_model( + Hdf5Context& context, + const hid_t file, + const Domain& domain) { + auto model_group = create_group(context, file, "model"); + + std::unordered_map node_dense_indices; + std::vector node_dense; + std::vector node_ids; + std::vector part_names; + std::vector instance_names; + std::vector local_labels; + std::vector coordinates; + node_dense.reserve(domain.nodes().size()); + node_ids.reserve(domain.nodes().size()); + part_names.reserve(domain.nodes().size()); + instance_names.reserve(domain.nodes().size()); + local_labels.reserve(domain.nodes().size()); + coordinates.reserve(domain.nodes().size() * 3U); + for (std::size_t index = 0; index < domain.nodes().size(); ++index) { + const Node& node = domain.nodes()[index]; + const auto dense_index = static_cast(index); + if (!node_dense_indices.emplace(node.id.value(), dense_index).second) { + fail("hdf5.write_failed", "Domain contains duplicate node IDs."); + } + node_dense.push_back(dense_index); + node_ids.push_back(node.id.value()); + part_names.push_back(node.origin.part_name); + instance_names.push_back(node.origin.instance_name); + local_labels.push_back(node.origin.local_label); + coordinates.push_back(node.position.x); + coordinates.push_back(node.position.y); + coordinates.push_back(node.position.z); + } + + { + auto nodes_group = create_group(context, model_group.get(), "nodes"); + write_numeric_1d( + context, + nodes_group.get(), + "dense_index", + H5T_STD_U64LE, + H5T_NATIVE_UINT64, + std::span{node_dense}); + write_numeric_1d( + context, + nodes_group.get(), + "internal_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64, + std::span{node_ids}); + write_strings_1d(context, nodes_group.get(), "part_name", part_names); + write_strings_1d( + context, nodes_group.get(), "instance_name", instance_names); + write_numeric_1d( + context, + nodes_group.get(), + "local_label", + H5T_STD_I64LE, + H5T_NATIVE_INT64, + std::span{local_labels}); + write_numeric_2d( + context, + nodes_group.get(), + "coordinates", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + domain.nodes().size(), + 3U, + std::span{coordinates}); + } + + std::vector element_dense; + std::vector element_ids; + std::vector connectivity; + std::vector element_section_ids; + element_dense.reserve(domain.beam_elements().size()); + element_ids.reserve(domain.beam_elements().size()); + connectivity.reserve(domain.beam_elements().size() * 2U); + element_section_ids.reserve(domain.beam_elements().size()); + for (std::size_t index = 0; index < domain.beam_elements().size(); ++index) { + const BeamElement& element = domain.beam_elements()[index]; + const auto first = node_dense_indices.find(element.nodes[0].value()); + const auto second = node_dense_indices.find(element.nodes[1].value()); + if (first == node_dense_indices.end() || + second == node_dense_indices.end()) { + fail( + "hdf5.write_failed", + "Beam connectivity references a missing node ID."); + } + element_dense.push_back(static_cast(index)); + element_ids.push_back(element.id.value()); + connectivity.push_back(first->second); + connectivity.push_back(second->second); + element_section_ids.push_back(element.section.value()); + } + + { + auto elements_group = + create_group(context, model_group.get(), "elements"); + write_numeric_1d( + context, + elements_group.get(), + "dense_index", + H5T_STD_U64LE, + H5T_NATIVE_UINT64, + std::span{element_dense}); + write_numeric_1d( + context, + elements_group.get(), + "internal_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64, + std::span{element_ids}); + write_numeric_2d( + context, + elements_group.get(), + "connectivity", + H5T_STD_U64LE, + H5T_NATIVE_UINT64, + domain.beam_elements().size(), + 2U, + std::span{connectivity}); + write_numeric_1d( + context, + elements_group.get(), + "section_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64, + std::span{element_section_ids}); + } + + std::vector section_ids; + std::vector shear_area_y; + std::vector shear_area_z; + std::vector shear_sources; + section_ids.reserve(domain.sections().size()); + shear_area_y.reserve(domain.sections().size()); + shear_area_z.reserve(domain.sections().size()); + shear_sources.reserve(domain.sections().size()); + for (const BeamSection& section : domain.sections()) { + section_ids.push_back(section.id.value()); + shear_area_y.push_back(section.shear_area_y); + shear_area_z.push_back(section.shear_area_z); + shear_sources.push_back( + section.shear_source == ShearPropertySource::input ? 0U : 1U); + } + + { + auto sections_group = + create_group(context, model_group.get(), "sections"); + write_numeric_1d( + context, + sections_group.get(), + "internal_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64, + std::span{section_ids}); + write_numeric_1d( + context, + sections_group.get(), + "shear_area_y", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + std::span{shear_area_y}); + write_numeric_1d( + context, + sections_group.get(), + "shear_area_z", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + std::span{shear_area_z}); + write_numeric_1d( + context, + sections_group.get(), + "shear_source", + H5T_STD_U8LE, + H5T_NATIVE_UINT8, + std::span{shear_sources}); + } +} + +void write_results( + Hdf5Context& context, + const hid_t file, + const ResultDatabase& database) { + auto results_group = create_group(context, file, "results"); + auto steps_group = create_group(context, results_group.get(), "steps"); + for (std::size_t step_index = 0; + step_index < database.steps.size(); + ++step_index) { + const ResultStep& step = database.steps[step_index]; + auto step_group = create_group( + context, steps_group.get(), std::to_string(step_index)); + write_string_attribute(context, step_group.get(), "name", step.name); + auto frames_group = + create_group(context, step_group.get(), "frames"); + for (std::size_t frame_index = 0; + frame_index < step.frames.size(); + ++frame_index) { + const ResultFrame& frame = step.frames[frame_index]; + auto frame_group = create_group( + context, frames_group.get(), std::to_string(frame_index)); + write_double_attribute( + context, frame_group.get(), "step_time", frame.step_time); + auto nodal_group = + create_group(context, frame_group.get(), "nodal"); + + std::vector node_ids; + std::vector displacement; + std::vector reaction; + node_ids.reserve(frame.nodal.node_ids.size()); + displacement.reserve(frame.nodal.displacement.size() * 6U); + reaction.reserve(frame.nodal.reaction.size() * 6U); + for (const NodeId node_id : frame.nodal.node_ids) { + node_ids.push_back(node_id.value()); + } + for (const auto& values : frame.nodal.displacement) { + displacement.insert( + displacement.end(), values.begin(), values.end()); + } + for (const auto& values : frame.nodal.reaction) { + reaction.insert(reaction.end(), values.begin(), values.end()); + } + write_numeric_1d( + context, + nodal_group.get(), + "node_ids", + H5T_STD_I64LE, + H5T_NATIVE_INT64, + std::span{node_ids}); + write_numeric_2d( + context, + nodal_group.get(), + "displacement", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + frame.nodal.node_ids.size(), + 6U, + std::span{displacement}); + write_numeric_2d( + context, + nodal_group.get(), + "reaction", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + frame.nodal.node_ids.size(), + 6U, + std::span{reaction}); + } + } +} + +void require_numeric_type( + Hdf5Context& context, + const hid_t dataset, + const hid_t expected, + const std::string_view name) { + auto type = checked_handle( + context, + call_hdf5([&] { return H5Dget_type(dataset); }), + &H5Tclose, + "datatype for dataset '" + std::string{name} + "'", + "hdf5.read_failed"); + const htri_t equal = + call_hdf5([&] { return H5Tequal(type.get(), expected); }); + if (equal <= 0) { + fail( + "hdf5.read_failed", + "Dataset '" + std::string{name} + + "' has an unexpected datatype."); + } +} + +std::size_t checked_size( + const hsize_t value, + const std::string_view description) { + if (value > static_cast(std::numeric_limits::max())) { + fail( + "hdf5.read_failed", + std::string{description} + " exceeds addressable memory."); + } + return static_cast(value); +} + +template +std::vector read_numeric_1d( + Hdf5Context& context, + const hid_t parent, + const std::string_view name, + const hid_t file_type, + const hid_t memory_type) { + const std::string owned_name{name}; + auto dataset = checked_handle( + context, + call_hdf5([&] { + return H5Dopen2(parent, owned_name.c_str(), H5P_DEFAULT); + }), + &H5Dclose, + "dataset '" + owned_name + "'", + "hdf5.read_failed"); + require_numeric_type(context, dataset.get(), file_type, owned_name); + auto space = checked_handle( + context, + call_hdf5([&] { return H5Dget_space(dataset.get()); }), + &H5Sclose, + "dataspace for dataset '" + owned_name + "'", + "hdf5.read_failed"); + const int rank = call_hdf5( + [&] { return H5Sget_simple_extent_ndims(space.get()); }); + if (rank != 1) { + fail( + "hdf5.read_failed", + "Dataset '" + owned_name + "' must have rank 1."); + } + std::array dimensions{}; + check_status( + call_hdf5([&] { + return H5Sget_simple_extent_dims( + space.get(), dimensions.data(), nullptr); + }), + "read shape for dataset '" + owned_name + "'", + "hdf5.read_failed"); + std::vector values( + checked_size(dimensions[0], "Dataset row count")); + if (!values.empty()) { + check_status( + call_hdf5([&] { + return H5Dread( + dataset.get(), + memory_type, + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + values.data()); + }), + "read dataset '" + owned_name + "'", + "hdf5.read_failed"); + } + return values; +} + +template +std::vector read_numeric_2d( + Hdf5Context& context, + const hid_t parent, + const std::string_view name, + const hid_t file_type, + const hid_t memory_type, + const std::size_t expected_columns) { + const std::string owned_name{name}; + auto dataset = checked_handle( + context, + call_hdf5([&] { + return H5Dopen2(parent, owned_name.c_str(), H5P_DEFAULT); + }), + &H5Dclose, + "dataset '" + owned_name + "'", + "hdf5.read_failed"); + require_numeric_type(context, dataset.get(), file_type, owned_name); + auto space = checked_handle( + context, + call_hdf5([&] { return H5Dget_space(dataset.get()); }), + &H5Sclose, + "dataspace for dataset '" + owned_name + "'", + "hdf5.read_failed"); + const int rank = call_hdf5( + [&] { return H5Sget_simple_extent_ndims(space.get()); }); + if (rank != 2) { + fail( + "hdf5.read_failed", + "Dataset '" + owned_name + "' must have rank 2."); + } + std::array dimensions{}; + check_status( + call_hdf5([&] { + return H5Sget_simple_extent_dims( + space.get(), dimensions.data(), nullptr); + }), + "read shape for dataset '" + owned_name + "'", + "hdf5.read_failed"); + if (dimensions[1] != static_cast(expected_columns)) { + fail( + "hdf5.read_failed", + "Dataset '" + owned_name + + "' has an unexpected column count."); + } + const std::size_t row_count = + checked_size(dimensions[0], "Dataset row count"); + if (row_count > + std::numeric_limits::max() / expected_columns) { + fail("hdf5.read_failed", "Dataset element count overflows."); + } + std::vector values(row_count * expected_columns); + if (!values.empty()) { + check_status( + call_hdf5([&] { + return H5Dread( + dataset.get(), + memory_type, + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + values.data()); + }), + "read dataset '" + owned_name + "'", + "hdf5.read_failed"); + } + return values; +} + +void require_utf8_string_type( + const hid_t type, + const std::string_view description) { + const H5T_class_t type_class = + call_hdf5([&] { return H5Tget_class(type); }); + const htri_t is_variable = + call_hdf5([&] { return H5Tis_variable_str(type); }); + const H5T_cset_t character_set = + call_hdf5([&] { return H5Tget_cset(type); }); + if ( + type_class != H5T_STRING || is_variable <= 0 || + character_set != H5T_CSET_UTF8) { + fail( + "hdf5.read_failed", + std::string{description} + + " must be a variable-length UTF-8 string."); + } +} + +std::vector read_strings_1d( + Hdf5Context& context, + const hid_t parent, + const std::string_view name) { + const std::string owned_name{name}; + auto dataset = checked_handle( + context, + call_hdf5([&] { + return H5Dopen2(parent, owned_name.c_str(), H5P_DEFAULT); + }), + &H5Dclose, + "dataset '" + owned_name + "'", + "hdf5.read_failed"); + auto type = checked_handle( + context, + call_hdf5([&] { return H5Dget_type(dataset.get()); }), + &H5Tclose, + "datatype for dataset '" + owned_name + "'", + "hdf5.read_failed"); + require_utf8_string_type(type.get(), "Dataset '" + owned_name + "'"); + auto space = checked_handle( + context, + call_hdf5([&] { return H5Dget_space(dataset.get()); }), + &H5Sclose, + "dataspace for dataset '" + owned_name + "'", + "hdf5.read_failed"); + const int rank = call_hdf5( + [&] { return H5Sget_simple_extent_ndims(space.get()); }); + if (rank != 1) { + fail( + "hdf5.read_failed", + "Dataset '" + owned_name + "' must have rank 1."); + } + std::array dimensions{}; + check_status( + call_hdf5([&] { + return H5Sget_simple_extent_dims( + space.get(), dimensions.data(), nullptr); + }), + "read shape for dataset '" + owned_name + "'", + "hdf5.read_failed"); + const std::size_t count = + checked_size(dimensions[0], "String dataset row count"); + std::vector raw_values(count, nullptr); + if (raw_values.empty()) { + return {}; + } + check_status( + call_hdf5([&] { + return H5Dread( + dataset.get(), + type.get(), + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + raw_values.data()); + }), + "read string dataset '" + owned_name + "'", + "hdf5.read_failed"); + Hdf5VlenBuffer reclaim{ + context, type.get(), space.get(), raw_values.data()}; + std::vector values; + values.reserve(count); + for (const char* value : raw_values) { + if (value == nullptr) { + fail( + "hdf5.read_failed", + "String dataset '" + owned_name + + "' contains a null value."); + } + values.emplace_back(value); + } + return values; +} + +std::string read_string_attribute( + Hdf5Context& context, + const hid_t parent, + const std::string_view name) { + const std::string owned_name{name}; + auto attribute = checked_handle( + context, + call_hdf5([&] { + return H5Aopen(parent, owned_name.c_str(), H5P_DEFAULT); + }), + &H5Aclose, + "attribute '" + owned_name + "'", + "hdf5.read_failed"); + auto type = checked_handle( + context, + call_hdf5([&] { return H5Aget_type(attribute.get()); }), + &H5Tclose, + "datatype for attribute '" + owned_name + "'", + "hdf5.read_failed"); + require_utf8_string_type(type.get(), "Attribute '" + owned_name + "'"); + auto space = checked_handle( + context, + call_hdf5([&] { return H5Aget_space(attribute.get()); }), + &H5Sclose, + "dataspace for attribute '" + owned_name + "'", + "hdf5.read_failed"); + const H5S_class_t space_class = + call_hdf5([&] { return H5Sget_simple_extent_type(space.get()); }); + if (space_class != H5S_SCALAR) { + fail( + "hdf5.read_failed", + "Attribute '" + owned_name + "' must be scalar."); + } + char* raw_value = nullptr; + check_status( + call_hdf5([&] { + return H5Aread(attribute.get(), type.get(), &raw_value); + }), + "read string attribute '" + owned_name + "'", + "hdf5.read_failed"); + Hdf5AllocatedString allocated{context, raw_value}; + if (raw_value == nullptr) { + fail( + "hdf5.read_failed", + "Attribute '" + owned_name + "' contains a null string."); + } + return std::string{raw_value}; +} + +double read_double_attribute( + Hdf5Context& context, + const hid_t parent, + const std::string_view name) { + const std::string owned_name{name}; + auto attribute = checked_handle( + context, + call_hdf5([&] { + return H5Aopen(parent, owned_name.c_str(), H5P_DEFAULT); + }), + &H5Aclose, + "attribute '" + owned_name + "'", + "hdf5.read_failed"); + auto type = checked_handle( + context, + call_hdf5([&] { return H5Aget_type(attribute.get()); }), + &H5Tclose, + "datatype for attribute '" + owned_name + "'", + "hdf5.read_failed"); + const htri_t equal = + call_hdf5([&] { return H5Tequal(type.get(), H5T_IEEE_F64LE); }); + if (equal <= 0) { + fail( + "hdf5.read_failed", + "Attribute '" + owned_name + + "' has an unexpected datatype."); + } + auto space = checked_handle( + context, + call_hdf5([&] { return H5Aget_space(attribute.get()); }), + &H5Sclose, + "dataspace for attribute '" + owned_name + "'", + "hdf5.read_failed"); + const H5S_class_t space_class = + call_hdf5([&] { return H5Sget_simple_extent_type(space.get()); }); + if (space_class != H5S_SCALAR) { + fail( + "hdf5.read_failed", + "Attribute '" + owned_name + "' must be scalar."); + } + double value{}; + check_status( + call_hdf5([&] { + return H5Aread(attribute.get(), H5T_NATIVE_DOUBLE, &value); + }), + "read floating-point attribute '" + owned_name + "'", + "hdf5.read_failed"); + return value; +} + +std::size_t group_link_count( + const hid_t group, + const std::string_view description) { + H5G_info_t information{}; + check_status( + call_hdf5([&] { return H5Gget_info(group, &information); }), + "inspect " + std::string{description}, + "hdf5.read_failed"); + return checked_size(information.nlinks, "Group link count"); +} + +Hdf5ModelSnapshot read_model( + Hdf5Context& context, + const hid_t file) { + auto model_group = open_group(context, file, "model"); + Hdf5ModelSnapshot model; + + { + auto nodes_group = open_group(context, model_group.get(), "nodes"); + const auto dense_indices = read_numeric_1d( + context, + nodes_group.get(), + "dense_index", + H5T_STD_U64LE, + H5T_NATIVE_UINT64); + const auto ids = read_numeric_1d( + context, + nodes_group.get(), + "internal_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64); + const auto part_names = + read_strings_1d(context, nodes_group.get(), "part_name"); + const auto instance_names = + read_strings_1d(context, nodes_group.get(), "instance_name"); + const auto local_labels = read_numeric_1d( + context, + nodes_group.get(), + "local_label", + H5T_STD_I64LE, + H5T_NATIVE_INT64); + const auto coordinates = read_numeric_2d( + context, + nodes_group.get(), + "coordinates", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + 3U); + const std::size_t count = dense_indices.size(); + if ( + ids.size() != count || part_names.size() != count || + instance_names.size() != count || local_labels.size() != count || + coordinates.size() != count * 3U) { + fail( + "hdf5.read_failed", + "Node datasets have inconsistent row counts."); + } + model.nodes.reserve(count); + std::unordered_set unique_ids; + for (std::size_t index = 0; index < count; ++index) { + if (dense_indices[index] != static_cast(index)) { + fail( + "hdf5.read_failed", + "Node dense_index must be contiguous and 0-based."); + } + if (!unique_ids.insert(ids[index]).second) { + fail( + "hdf5.invalid_model_data", + "Serialized model contains duplicate node internal IDs."); + } + const Vec3 position{ + coordinates[index * 3U], + coordinates[index * 3U + 1U], + coordinates[index * 3U + 2U], + }; + if (!is_finite(position)) { + fail( + "hdf5.invalid_model_data", + "Serialized model contains nonfinite node coordinates."); + } + model.nodes.push_back({ + dense_indices[index], + NodeId{ids[index]}, + EntityOrigin{ + part_names[index], + instance_names[index], + local_labels[index], + }, + position, + }); + } + } + + { + auto sections_group = + open_group(context, model_group.get(), "sections"); + const auto ids = read_numeric_1d( + context, + sections_group.get(), + "internal_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64); + const auto shear_area_y = read_numeric_1d( + context, + sections_group.get(), + "shear_area_y", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE); + const auto shear_area_z = read_numeric_1d( + context, + sections_group.get(), + "shear_area_z", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE); + const auto shear_sources = read_numeric_1d( + context, + sections_group.get(), + "shear_source", + H5T_STD_U8LE, + H5T_NATIVE_UINT8); + if ( + shear_area_y.size() != ids.size() || + shear_area_z.size() != ids.size() || + shear_sources.size() != ids.size()) { + fail( + "hdf5.read_failed", + "Section datasets have inconsistent row counts."); + } + model.sections.reserve(ids.size()); + std::unordered_set unique_ids; + for (std::size_t index = 0; index < ids.size(); ++index) { + if (!unique_ids.insert(ids[index]).second) { + fail( + "hdf5.invalid_model_data", + "Serialized model contains duplicate section internal " + "IDs."); + } + if ( + !std::isfinite(shear_area_y[index]) || + !std::isfinite(shear_area_z[index]) || + shear_area_y[index] <= 0.0 || shear_area_z[index] <= 0.0) { + fail( + "hdf5.invalid_model_data", + "Serialized model requires finite positive shear areas."); + } + ShearPropertySource source; + if (shear_sources[index] == 0U) { + source = ShearPropertySource::input; + } else if (shear_sources[index] == 1U) { + source = ShearPropertySource::phase1_default; + } else { + fail( + "hdf5.read_failed", + "Section shear_source contains an unknown value."); + } + model.sections.push_back({ + SectionId{ids[index]}, + shear_area_y[index], + shear_area_z[index], + source, + }); + } + } + + { + auto elements_group = + open_group(context, model_group.get(), "elements"); + const auto dense_indices = read_numeric_1d( + context, + elements_group.get(), + "dense_index", + H5T_STD_U64LE, + H5T_NATIVE_UINT64); + const auto ids = read_numeric_1d( + context, + elements_group.get(), + "internal_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64); + const auto connectivity = read_numeric_2d( + context, + elements_group.get(), + "connectivity", + H5T_STD_U64LE, + H5T_NATIVE_UINT64, + 2U); + const auto section_ids = read_numeric_1d( + context, + elements_group.get(), + "section_id", + H5T_STD_I64LE, + H5T_NATIVE_INT64); + const std::size_t count = dense_indices.size(); + if ( + ids.size() != count || connectivity.size() != count * 2U || + section_ids.size() != count) { + fail( + "hdf5.read_failed", + "Element datasets have inconsistent row counts."); + } + std::unordered_set available_sections; + for (const Hdf5SectionSnapshot& section : model.sections) { + available_sections.insert(section.id.value()); + } + model.elements.reserve(count); + std::unordered_set unique_ids; + for (std::size_t index = 0; index < count; ++index) { + if (dense_indices[index] != static_cast(index)) { + fail( + "hdf5.read_failed", + "Element dense_index must be contiguous and 0-based."); + } + if (!unique_ids.insert(ids[index]).second) { + fail( + "hdf5.invalid_model_data", + "Serialized model contains duplicate element internal " + "IDs."); + } + const std::uint64_t first = connectivity[index * 2U]; + const std::uint64_t second = connectivity[index * 2U + 1U]; + if ( + first >= static_cast(model.nodes.size()) || + second >= static_cast(model.nodes.size())) { + fail( + "hdf5.read_failed", + "Element connectivity references a missing dense node."); + } + if (!available_sections.contains(section_ids[index])) { + fail( + "hdf5.read_failed", + "Element section_id references a missing section."); + } + model.elements.push_back({ + dense_indices[index], + ElementId{ids[index]}, + {first, second}, + SectionId{section_ids[index]}, + }); + } + } + + return model; +} + +ResultDatabase read_database( + Hdf5Context& context, + const hid_t file, + const std::string& version) { + ResultDatabase database{version, {}}; + auto results_group = open_group(context, file, "results"); + auto steps_group = open_group(context, results_group.get(), "steps"); + const std::size_t step_count = + group_link_count(steps_group.get(), "results/steps group"); + database.steps.reserve(step_count); + for (std::size_t step_index = 0; step_index < step_count; ++step_index) { + auto step_group = open_group( + context, steps_group.get(), std::to_string(step_index)); + ResultStep step; + step.name = read_string_attribute(context, step_group.get(), "name"); + auto frames_group = + open_group(context, step_group.get(), "frames"); + const std::size_t frame_count = + group_link_count(frames_group.get(), "frames group"); + step.frames.reserve(frame_count); + for (std::size_t frame_index = 0; + frame_index < frame_count; + ++frame_index) { + auto frame_group = open_group( + context, frames_group.get(), std::to_string(frame_index)); + ResultFrame frame; + frame.step_time = read_double_attribute( + context, frame_group.get(), "step_time"); + auto nodal_group = + open_group(context, frame_group.get(), "nodal"); + const auto node_ids = read_numeric_1d( + context, + nodal_group.get(), + "node_ids", + H5T_STD_I64LE, + H5T_NATIVE_INT64); + const auto displacement = read_numeric_2d( + context, + nodal_group.get(), + "displacement", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + 6U); + const auto reaction = read_numeric_2d( + context, + nodal_group.get(), + "reaction", + H5T_IEEE_F64LE, + H5T_NATIVE_DOUBLE, + 6U); + if ( + displacement.size() != node_ids.size() * 6U || + reaction.size() != node_ids.size() * 6U) { + fail( + "hdf5.read_failed", + "Nodal result datasets have inconsistent row counts."); + } + frame.nodal.node_ids.reserve(node_ids.size()); + frame.nodal.displacement.reserve(node_ids.size()); + frame.nodal.reaction.reserve(node_ids.size()); + for (std::size_t index = 0; index < node_ids.size(); ++index) { + frame.nodal.node_ids.emplace_back(node_ids[index]); + std::array displacement_values{}; + std::array reaction_values{}; + std::copy_n( + displacement.begin() + + static_cast(index * 6U), + 6U, + displacement_values.begin()); + std::copy_n( + reaction.begin() + + static_cast(index * 6U), + 6U, + reaction_values.begin()); + frame.nodal.displacement.push_back(displacement_values); + frame.nodal.reaction.push_back(reaction_values); + } + step.frames.push_back(std::move(frame)); + } + database.steps.push_back(std::move(step)); + } + return database; +} + +} // namespace + +std::vector write_hdf5( + const std::filesystem::path& path, + const Domain& domain, + const ResultDatabase& database) { + const Status validation = validate_result_database(database); + if (!validation.succeeded) { + return validation.diagnostics; + } + if (database.schema_version != schema_version) { + return {error_diagnostic( + "hdf5.unsupported_schema", + "HDF5 writer supports schema version 1.0.0 only.")}; + } + for (const ResultStep& step : database.steps) { + for (const ResultFrame& frame : step.frames) { + if (!frame.diagnostics.empty()) { + return {error_diagnostic( + "hdf5.unsupported_result_diagnostics", + "HDF5 schema 1.0.0 does not store frame diagnostics.")}; + } + } + } + std::unordered_set model_node_ids; + for (const Node& node : domain.nodes()) { + model_node_ids.insert(node.id.value()); + } + if (const auto unknown = + find_unknown_result_node(database, model_node_ids)) { + return {error_diagnostic( + "hdf5.result_node_not_in_model", + "Result node " + std::to_string(*unknown) + + " is not present in the serialized model.")}; + } + + Hdf5Context context; + try { + const std::string encoded_path = path_utf8(path); + auto file = checked_handle( + context, + call_hdf5([&] { + return H5Fcreate( + encoded_path.c_str(), + H5F_ACC_TRUNC, + H5P_DEFAULT, + H5P_DEFAULT); + }), + &H5Fclose, + "file '" + encoded_path + "'", + "hdf5.create_failed"); + write_string_attribute( + context, file.get(), "schema_version", schema_version); + write_model(context, file.get(), domain); + write_results(context, file.get(), database); + check_status( + call_hdf5([&] { return H5Fflush(file.get(), H5F_SCOPE_GLOBAL); }), + "flush output file", + "hdf5.write_failed"); + file.reset(); + if (!context.close_error.empty()) { + fail("hdf5.close_failed", context.close_error); + } + return {}; + } catch (const Hdf5Error& error) { + return {error_diagnostic(error.code(), error.what())}; + } catch (const std::exception& error) { + return {error_diagnostic("hdf5.write_failed", error.what())}; + } +} + +Hdf5ReadResult read_hdf5_results(const std::filesystem::path& path) { + Hdf5Context context; + try { + const std::string encoded_path = path_utf8(path); + auto file = checked_handle( + context, + call_hdf5([&] { + return H5Fopen( + encoded_path.c_str(), H5F_ACC_RDONLY, H5P_DEFAULT); + }), + &H5Fclose, + "file '" + encoded_path + "'", + "hdf5.open_failed"); + const std::string version = + read_string_attribute(context, file.get(), "schema_version"); + if (version != schema_version) { + fail( + "hdf5.unsupported_schema", + "HDF5 reader supports schema version 1.0.0 only."); + } + Hdf5ModelSnapshot model = read_model(context, file.get()); + ResultDatabase database = + read_database(context, file.get(), version); + const Status validation = validate_result_database(database); + if (!validation.succeeded) { + file.reset(); + return {std::nullopt, validation.diagnostics, std::nullopt}; + } + std::unordered_set model_node_ids; + for (const Hdf5NodeSnapshot& node : model.nodes) { + model_node_ids.insert(node.id.value()); + } + if (const auto unknown = + find_unknown_result_node(database, model_node_ids)) { + fail( + "hdf5.result_node_not_in_model", + "Result node " + std::to_string(*unknown) + + " is not present in the serialized model."); + } + file.reset(); + if (!context.close_error.empty()) { + fail("hdf5.close_failed", context.close_error); + } + return { + std::move(database), + {}, + std::move(model), + }; + } catch (const Hdf5Error& error) { + return { + std::nullopt, + {error_diagnostic(error.code(), error.what())}, + std::nullopt, + }; + } catch (const std::exception& error) { + return { + std::nullopt, + {error_diagnostic("hdf5.read_failed", error.what())}, + std::nullopt, + }; + } +} + +} // namespace fesa diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index dd083b2..a7f4e48 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -406,3 +406,46 @@ add_test( COMMAND "$" --gtest_filter=ResultDatabase.* ) + +add_executable(fesa_hdf5_results_tests + integration/io/hdf5_results_test.cpp +) + +target_compile_features(fesa_hdf5_results_tests PRIVATE cxx_std_20) +target_compile_options( + fesa_hdf5_results_tests + PRIVATE + /W4 + /permissive- + /EHsc +) +target_compile_definitions( + fesa_hdf5_results_tests + PRIVATE + FESA_TEST_BINARY_DIR="${CMAKE_BINARY_DIR}" +) + +target_link_libraries(fesa_hdf5_results_tests + PRIVATE + fesa_core + HDF5::HDF5 + GTest::gtest_main +) + +add_test( + NAME Hdf5 + COMMAND "$" + --gtest_filter=Hdf5.* +) + +add_test( + NAME ResultRoundTrip + COMMAND "$" + --gtest_filter=ResultRoundTrip.* +) + +set_property( + TEST Hdf5 ResultRoundTrip + PROPERTY ENVIRONMENT_MODIFICATION + ${FESA_DEPENDENCY_RUNTIME_MODIFICATIONS} +) diff --git a/tests/integration/io/hdf5_results_test.cpp b/tests/integration/io/hdf5_results_test.cpp new file mode 100644 index 0000000..2d4221c --- /dev/null +++ b/tests/integration/io/hdf5_results_test.cpp @@ -0,0 +1,566 @@ +#include +#include + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +namespace { + +class TestHdf5Handle final { +public: + using CloseFunction = herr_t (*)(hid_t); + + TestHdf5Handle(hid_t id, CloseFunction close) + : id_{id}, close_{close} { + if (id_ < 0) { + throw std::runtime_error{"Failed to open test HDF5 resource."}; + } + } + + ~TestHdf5Handle() { + if (id_ >= 0) { + close_(id_); + } + } + + TestHdf5Handle(const TestHdf5Handle&) = delete; + TestHdf5Handle& operator=(const TestHdf5Handle&) = delete; + + TestHdf5Handle(TestHdf5Handle&& other) noexcept + : id_{std::exchange(other.id_, H5I_INVALID_HID)}, + close_{other.close_} {} + + TestHdf5Handle& operator=(TestHdf5Handle&&) = delete; + + [[nodiscard]] hid_t get() const noexcept { + return id_; + } + +private: + hid_t id_; + CloseFunction close_; +}; + +std::string hdf5_path(const std::filesystem::path& path) { + const std::u8string value = path.u8string(); + return {reinterpret_cast(value.data()), value.size()}; +} + +void require_hdf5_status(const herr_t status) { + if (status < 0) { + throw std::runtime_error{"Test HDF5 mutation failed."}; + } +} + +void replace_step_time_with_vector(const std::filesystem::path& path) { + const std::string encoded_path = hdf5_path(path); + TestHdf5Handle file{ + H5Fopen(encoded_path.c_str(), H5F_ACC_RDWR, H5P_DEFAULT), + &H5Fclose, + }; + TestHdf5Handle frame{ + H5Gopen2( + file.get(), "/results/steps/0/frames/0", H5P_DEFAULT), + &H5Gclose, + }; + require_hdf5_status(H5Adelete(frame.get(), "step_time")); + const std::array dimensions{2U}; + TestHdf5Handle space{ + H5Screate_simple(1, dimensions.data(), nullptr), + &H5Sclose, + }; + TestHdf5Handle attribute{ + H5Acreate2( + frame.get(), + "step_time", + H5T_IEEE_F64LE, + space.get(), + H5P_DEFAULT, + H5P_DEFAULT), + &H5Aclose, + }; + const std::array values{1.25, 2.5}; + require_hdf5_status( + H5Awrite(attribute.get(), H5T_NATIVE_DOUBLE, values.data())); +} + +void write_int64_dataset( + const std::filesystem::path& path, + const std::string_view dataset_path, + const std::span values) { + const std::string encoded_path = hdf5_path(path); + const std::string owned_dataset_path{dataset_path}; + TestHdf5Handle file{ + H5Fopen(encoded_path.c_str(), H5F_ACC_RDWR, H5P_DEFAULT), + &H5Fclose, + }; + TestHdf5Handle dataset{ + H5Dopen2( + file.get(), owned_dataset_path.c_str(), H5P_DEFAULT), + &H5Dclose, + }; + require_hdf5_status(H5Dwrite( + dataset.get(), + H5T_NATIVE_INT64, + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + values.data())); +} + +void write_double_dataset( + const std::filesystem::path& path, + const std::string_view dataset_path, + const std::span values) { + const std::string encoded_path = hdf5_path(path); + const std::string owned_dataset_path{dataset_path}; + TestHdf5Handle file{ + H5Fopen(encoded_path.c_str(), H5F_ACC_RDWR, H5P_DEFAULT), + &H5Fclose, + }; + TestHdf5Handle dataset{ + H5Dopen2( + file.get(), owned_dataset_path.c_str(), H5P_DEFAULT), + &H5Dclose, + }; + require_hdf5_status(H5Dwrite( + dataset.get(), + H5T_NATIVE_DOUBLE, + H5S_ALL, + H5S_ALL, + H5P_DEFAULT, + values.data())); +} + +std::filesystem::path round_trip_path() { + return std::filesystem::path{FESA_TEST_BINARY_DIR} / "Testing" / + "Temporary" / "fesa-round-trip.h5"; +} + +fesa::Domain make_domain() { + fesa::DomainBuilder builder; + builder.add_node({ + fesa::NodeId{42}, + fesa::EntityOrigin{"BeamPart", "Beam-1", 1001}, + fesa::Vec3{1.25, -2.5, 3.75}, + }); + builder.add_node({ + fesa::NodeId{7}, + fesa::EntityOrigin{"BeamPart", "Beam-1", 1002}, + fesa::Vec3{4.0, 5.5, -6.25}, + }); + builder.add_material({ + fesa::MaterialId{6}, + "Steel", + 210.0e9, + 0.3, + }); + builder.add_section({ + fesa::SectionId{4}, + "InputShear", + 0.04, + 1.2e-4, + 1.4e-4, + 2.0e-4, + 0.031, + 0.032, + fesa::ShearPropertySource::input, + fesa::Vec3{0.0, 1.0, 0.0}, + {}, + }); + builder.add_section({ + fesa::SectionId{12}, + "DefaultShear", + 0.06, + 1.5e-4, + 1.7e-4, + 2.2e-4, + 0.05, + 0.05, + fesa::ShearPropertySource::phase1_default, + fesa::Vec3{0.0, 1.0, 0.0}, + {}, + }); + builder.add_beam_element({ + fesa::ElementId{9}, + fesa::EntityOrigin{"BeamPart", "Beam-1", 2001}, + {fesa::NodeId{7}, fesa::NodeId{42}}, + fesa::MaterialId{6}, + fesa::SectionId{4}, + }); + builder.add_beam_element({ + fesa::ElementId{17}, + fesa::EntityOrigin{"BeamPart", "Beam-1", 2002}, + {fesa::NodeId{42}, fesa::NodeId{7}}, + fesa::MaterialId{6}, + fesa::SectionId{4}, + }); + builder.set_step({"Load/Case", {}, {}}); + + auto built = std::move(builder).build(); + EXPECT_TRUE(built.domain.has_value()); + EXPECT_TRUE(built.diagnostics.empty()); + return std::move(*built.domain); +} + +fesa::ResultDatabase make_database() { + return { + "1.0.0", + {{ + "Load/Case", + {{ + 1.25, + { + {fesa::NodeId{7}, fesa::NodeId{42}}, + { + {1.0, 2.0, 3.0, 4.0, 5.0, 6.0}, + {-1.0, -2.0, -3.0, -4.0, -5.0, -6.0}, + }, + { + {10.0, 20.0, 30.0, 40.0, 50.0, 60.0}, + {-10.0, -20.0, -30.0, -40.0, -50.0, -60.0}, + }, + }, + {}, + }}, + }}, + }; +} + +bool has_results_error( + const std::vector& diagnostics, + const std::string_view code) { + return std::ranges::any_of( + diagnostics, + [code](const fesa::Diagnostic& diagnostic) { + return diagnostic.stage == fesa::DiagnosticStage::results && + diagnostic.severity == fesa::Severity::error && + diagnostic.code == code; + }); +} + +TEST(ResultRoundTrip, PreservesMinimalSchemaModelAndNodalResults) { + const auto path = round_trip_path(); + std::filesystem::create_directories(path.parent_path()); + std::filesystem::remove(path); + const auto domain = make_domain(); + const auto database = make_database(); + + const auto write_diagnostics = + fesa::write_hdf5(path, domain, database); + ASSERT_TRUE(write_diagnostics.empty()); + + const auto read = fesa::read_hdf5_results(path); + ASSERT_TRUE(read.diagnostics.empty()); + ASSERT_TRUE(read.database.has_value()); + ASSERT_TRUE(read.model.has_value()); + + EXPECT_EQ(read.database->schema_version, "1.0.0"); + ASSERT_EQ(read.database->steps.size(), 1U); + const auto& step = read.database->steps[0]; + EXPECT_EQ(step.name, "Load/Case"); + ASSERT_EQ(step.frames.size(), 1U); + const auto& frame = step.frames[0]; + EXPECT_DOUBLE_EQ(frame.step_time, 1.25); + EXPECT_EQ( + frame.nodal.node_ids, + (std::vector{fesa::NodeId{7}, fesa::NodeId{42}})); + EXPECT_EQ( + frame.nodal.displacement[0], + (std::array{1.0, 2.0, 3.0, 4.0, 5.0, 6.0})); + EXPECT_EQ( + frame.nodal.displacement[1], + (std::array{-1.0, -2.0, -3.0, -4.0, -5.0, -6.0})); + EXPECT_EQ( + frame.nodal.reaction[0], + (std::array{10.0, 20.0, 30.0, 40.0, 50.0, 60.0})); + EXPECT_EQ( + frame.nodal.reaction[1], + (std::array{ + -10.0, -20.0, -30.0, -40.0, -50.0, -60.0})); + EXPECT_TRUE(frame.diagnostics.empty()); + + ASSERT_EQ(read.model->nodes.size(), 2U); + EXPECT_EQ(read.model->nodes[0].dense_index, 0U); + EXPECT_EQ(read.model->nodes[0].id, fesa::NodeId{42}); + EXPECT_EQ( + read.model->nodes[0].origin, + (fesa::EntityOrigin{"BeamPart", "Beam-1", 1001})); + EXPECT_DOUBLE_EQ(read.model->nodes[0].coordinates.x, 1.25); + EXPECT_DOUBLE_EQ(read.model->nodes[0].coordinates.y, -2.5); + EXPECT_DOUBLE_EQ(read.model->nodes[0].coordinates.z, 3.75); + EXPECT_EQ(read.model->nodes[1].dense_index, 1U); + EXPECT_EQ(read.model->nodes[1].id, fesa::NodeId{7}); + EXPECT_EQ( + read.model->nodes[1].origin, + (fesa::EntityOrigin{"BeamPart", "Beam-1", 1002})); + EXPECT_DOUBLE_EQ(read.model->nodes[1].coordinates.x, 4.0); + EXPECT_DOUBLE_EQ(read.model->nodes[1].coordinates.y, 5.5); + EXPECT_DOUBLE_EQ(read.model->nodes[1].coordinates.z, -6.25); + + ASSERT_EQ(read.model->elements.size(), 2U); + EXPECT_EQ(read.model->elements[0].dense_index, 0U); + EXPECT_EQ(read.model->elements[0].id, fesa::ElementId{9}); + EXPECT_EQ( + read.model->elements[0].connectivity, + (std::array{1U, 0U})); + EXPECT_EQ(read.model->elements[0].section, fesa::SectionId{4}); + EXPECT_EQ(read.model->elements[1].dense_index, 1U); + EXPECT_EQ(read.model->elements[1].id, fesa::ElementId{17}); + EXPECT_EQ( + read.model->elements[1].connectivity, + (std::array{0U, 1U})); + EXPECT_EQ(read.model->elements[1].section, fesa::SectionId{4}); + + ASSERT_EQ(read.model->sections.size(), 2U); + EXPECT_EQ(read.model->sections[0].id, fesa::SectionId{4}); + EXPECT_DOUBLE_EQ(read.model->sections[0].shear_area_y, 0.031); + EXPECT_DOUBLE_EQ(read.model->sections[0].shear_area_z, 0.032); + EXPECT_EQ( + read.model->sections[0].shear_source, + fesa::ShearPropertySource::input); + EXPECT_EQ(read.model->sections[1].id, fesa::SectionId{12}); + EXPECT_DOUBLE_EQ(read.model->sections[1].shear_area_y, 0.05); + EXPECT_DOUBLE_EQ(read.model->sections[1].shear_area_z, 0.05); + EXPECT_EQ( + read.model->sections[1].shear_source, + fesa::ShearPropertySource::phase1_default); +} + +TEST(Hdf5, RejectsInvalidResultDatabaseBeforeWriting) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / "fesa-invalid.h5"; + std::filesystem::remove(path); + const auto domain = make_domain(); + auto database = make_database(); + database.steps[0].frames[0].nodal.reaction.pop_back(); + + const auto diagnostics = fesa::write_hdf5(path, domain, database); + + EXPECT_TRUE( + has_results_error(diagnostics, "results.nodal_size_mismatch")); + EXPECT_FALSE(std::filesystem::exists(path)); +} + +TEST(Hdf5, RejectsFrameDiagnosticsThatSchemaCannotRepresent) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-unsupported-diagnostics.h5"; + std::filesystem::remove(path); + const auto domain = make_domain(); + auto database = make_database(); + database.steps[0].frames[0].diagnostics.push_back({ + fesa::DiagnosticStage::solver, + fesa::Severity::warning, + "solver.residual", + "Residual diagnostic", + std::nullopt, + }); + + const auto diagnostics = fesa::write_hdf5(path, domain, database); + + EXPECT_TRUE(has_results_error( + diagnostics, "hdf5.unsupported_result_diagnostics")); + EXPECT_FALSE(std::filesystem::exists(path)); +} + +TEST(Hdf5, RejectsResultNodeMissingFromDomainBeforeWriting) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-missing-result-node.h5"; + std::filesystem::remove(path); + const auto domain = make_domain(); + auto database = make_database(); + database.steps[0].frames[0].nodal.node_ids[0] = fesa::NodeId{999}; + + const auto diagnostics = fesa::write_hdf5(path, domain, database); + + EXPECT_TRUE(has_results_error( + diagnostics, "hdf5.result_node_not_in_model")); + EXPECT_FALSE(std::filesystem::exists(path)); +} + +TEST(Hdf5, ReportsMissingFileAtResultsStage) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / "fesa-missing.h5"; + std::filesystem::remove(path); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE(has_results_error(read.diagnostics, "hdf5.open_failed")); +} + +TEST(Hdf5, RejectsNonScalarStepTimeAttribute) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-nonscalar-step-time.h5"; + std::filesystem::remove(path); + const auto write_diagnostics = + fesa::write_hdf5(path, make_domain(), make_database()); + ASSERT_TRUE(write_diagnostics.empty()); + replace_step_time_with_vector(path); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE(has_results_error(read.diagnostics, "hdf5.read_failed")); +} + +TEST(Hdf5, RejectsResultNodeMissingFromSerializedModel) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-invalid-result-node.h5"; + std::filesystem::remove(path); + const auto write_diagnostics = + fesa::write_hdf5(path, make_domain(), make_database()); + ASSERT_TRUE(write_diagnostics.empty()); + const std::array node_ids{999, 42}; + write_int64_dataset( + path, + "/results/steps/0/frames/0/nodal/node_ids", + node_ids); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE(has_results_error( + read.diagnostics, "hdf5.result_node_not_in_model")); +} + +TEST(Hdf5, RejectsDuplicateSerializedNodeIds) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-duplicate-node-ids.h5"; + std::filesystem::remove(path); + ASSERT_TRUE( + fesa::write_hdf5(path, make_domain(), make_database()).empty()); + const std::array ids{42, 42}; + write_int64_dataset(path, "/model/nodes/internal_id", ids); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE( + has_results_error(read.diagnostics, "hdf5.invalid_model_data")); +} + +TEST(Hdf5, RejectsDuplicateSerializedElementIds) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-duplicate-element-ids.h5"; + std::filesystem::remove(path); + ASSERT_TRUE( + fesa::write_hdf5(path, make_domain(), make_database()).empty()); + const std::array ids{9, 9}; + write_int64_dataset(path, "/model/elements/internal_id", ids); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE( + has_results_error(read.diagnostics, "hdf5.invalid_model_data")); +} + +TEST(Hdf5, RejectsDuplicateSerializedSectionIds) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-duplicate-section-ids.h5"; + std::filesystem::remove(path); + ASSERT_TRUE( + fesa::write_hdf5(path, make_domain(), make_database()).empty()); + const std::array ids{4, 4}; + write_int64_dataset(path, "/model/sections/internal_id", ids); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE( + has_results_error(read.diagnostics, "hdf5.invalid_model_data")); +} + +TEST(Hdf5, RejectsNonfiniteSerializedCoordinates) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-nonfinite-coordinates.h5"; + std::filesystem::remove(path); + ASSERT_TRUE( + fesa::write_hdf5(path, make_domain(), make_database()).empty()); + const std::array coordinates{ + std::numeric_limits::quiet_NaN(), + -2.5, + 3.75, + 4.0, + 5.5, + -6.25, + }; + write_double_dataset(path, "/model/nodes/coordinates", coordinates); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE( + has_results_error(read.diagnostics, "hdf5.invalid_model_data")); +} + +TEST(Hdf5, RejectsNonfiniteSerializedShearArea) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-nonfinite-shear-area.h5"; + std::filesystem::remove(path); + ASSERT_TRUE( + fesa::write_hdf5(path, make_domain(), make_database()).empty()); + const std::array shear_areas{ + std::numeric_limits::infinity(), + 0.05, + }; + write_double_dataset(path, "/model/sections/shear_area_y", shear_areas); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE( + has_results_error(read.diagnostics, "hdf5.invalid_model_data")); +} + +TEST(Hdf5, RejectsNonpositiveSerializedShearArea) { + const auto path = std::filesystem::path{FESA_TEST_BINARY_DIR} / + "Testing" / "Temporary" / + "fesa-nonpositive-shear-area.h5"; + std::filesystem::remove(path); + ASSERT_TRUE( + fesa::write_hdf5(path, make_domain(), make_database()).empty()); + const std::array shear_areas{-0.031, 0.05}; + write_double_dataset(path, "/model/sections/shear_area_y", shear_areas); + + const auto read = fesa::read_hdf5_results(path); + + EXPECT_FALSE(read.database.has_value()); + EXPECT_FALSE(read.model.has_value()); + EXPECT_TRUE( + has_results_error(read.diagnostics, "hdf5.invalid_model_data")); +} + +} // namespace