From 92dceb39aafaa41bf807f10a08f603218f6b2634 Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 09:03:48 +0200 Subject: [PATCH 01/15] refactor(cpp): make file_magic, scalar_size() and fixed_header_size() constexpr, add constexpr current_format_version --- include/ndtbl/detail/binary_io.hpp | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/include/ndtbl/detail/binary_io.hpp b/include/ndtbl/detail/binary_io.hpp index 0ca1040..6ca8031 100644 --- a/include/ndtbl/detail/binary_io.hpp +++ b/include/ndtbl/detail/binary_io.hpp @@ -25,7 +25,10 @@ namespace detail { * This constant is part of the binary file format implementation and is not * intended to be consumed directly by library users. */ -static const char file_magic[8] = { 'N', 'D', 'T', 'B', 'L', '\0', '\0', '\0' }; +static constexpr char file_magic[8] = { 'N', 'D', 'T', 'B', + 'L', '\0', '\0', '\0' }; + +static constexpr std::uint8_t current_format_version = 1u; /** * @brief Write one exact byte sequence to a binary stream. @@ -213,7 +216,7 @@ require_zero(std::uint64_t value, const std::string& what) } } -inline std::size_t +inline constexpr std::size_t fixed_header_size() { return sizeof(file_magic) + sizeof(std::uint8_t) + sizeof(std::uint8_t) + @@ -284,7 +287,7 @@ write_group_stream_impl(std::ostream& os, const std::size_t payload_offset = metadata_size(metadata); write_bytes(os, file_magic, sizeof(file_magic)); - write_uint_le(os, 1u); + write_uint_le(os, current_format_version); write_uint_le(os, static_cast(metadata.value_type)); write_uint_le(os, 0u); @@ -381,7 +384,7 @@ verify_magic(std::istream& is) } } -inline std::size_t +inline constexpr std::size_t scalar_size(scalar_type type) { if (type == scalar_type::float32) { @@ -413,7 +416,7 @@ read_group_layout_impl(std::istream& is) verify_magic(is); const std::uint8_t version = read_uint_le(is); - if (version != 1u) { + if (version != current_format_version) { throw std::runtime_error("unsupported ndtbl version"); } From 829b6582ae53ecab6bde83a14a346b5048d8aa33 Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 09:13:56 +0200 Subject: [PATCH 02/15] refactor(cpp): moved current_format_version to types.hpp, remove duplicate --- include/ndtbl/detail/binary_io.hpp | 2 -- include/ndtbl/metadata.hpp | 2 +- include/ndtbl/types.hpp | 5 +++++ 3 files changed, 6 insertions(+), 3 deletions(-) diff --git a/include/ndtbl/detail/binary_io.hpp b/include/ndtbl/detail/binary_io.hpp index 6ca8031..9fbb1bf 100644 --- a/include/ndtbl/detail/binary_io.hpp +++ b/include/ndtbl/detail/binary_io.hpp @@ -28,8 +28,6 @@ namespace detail { static constexpr char file_magic[8] = { 'N', 'D', 'T', 'B', 'L', '\0', '\0', '\0' }; -static constexpr std::uint8_t current_format_version = 1u; - /** * @brief Write one exact byte sequence to a binary stream. * diff --git a/include/ndtbl/metadata.hpp b/include/ndtbl/metadata.hpp index 4b64b7c..b410bbd 100644 --- a/include/ndtbl/metadata.hpp +++ b/include/ndtbl/metadata.hpp @@ -31,7 +31,7 @@ struct GroupMetadata /// Field names in payload storage order. std::vector field_names; /// ndtbl binary file format version. - std::uint8_t format_version = 1u; + std::uint8_t format_version = current_format_version; }; } // namespace ndtbl diff --git a/include/ndtbl/types.hpp b/include/ndtbl/types.hpp index e277a5c..78ce50e 100644 --- a/include/ndtbl/types.hpp +++ b/include/ndtbl/types.hpp @@ -32,6 +32,11 @@ enum class bounds_policy : std::uint8_t throw_error = 2 }; +/** + * @brief Current ndtbl binary file format version. + */ +static constexpr std::uint8_t current_format_version = 1u; + /** * @brief Map a supported C++ scalar type to the ndtbl on-disk type tag. * From db0543dd2d5f0bcdd37823a37be1c8c36b5e38c7 Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 12:29:38 +0200 Subject: [PATCH 03/15] chore: TensorStencil and Grid must have at least a dimension of 1 --- include/ndtbl/grid.hpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/include/ndtbl/grid.hpp b/include/ndtbl/grid.hpp index de46e21..0cc930a 100644 --- a/include/ndtbl/grid.hpp +++ b/include/ndtbl/grid.hpp @@ -34,6 +34,11 @@ pow_size(std::size_t base, std::size_t exponent) template class TensorStencil { + static_assert(Dim > 0, "ndtbl grids must have at least one dimension"); + static_assert( + PointsPerAxis > 0, + "ndtbl tensor stencils must contain at least one point per axis"); + public: /// Number of dimensions. static constexpr std::size_t dimensions = Dim; @@ -105,6 +110,8 @@ using CubicStencil = TensorStencil; template class Grid { + static_assert(Dim > 0, "ndtbl grids must have at least one dimension"); + public: /// Number of dimensions. static constexpr std::size_t dimensions = Dim; From 0961950cf02e3df149897fd591050fc9f0cc24e3 Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 13:49:33 +0200 Subject: [PATCH 04/15] refactor(cpp): template RuntimeFieldGroup on output type, use per-thread scratch, write without static_cast if payload type matches --- include/ndtbl/io.hpp | 23 ++--- include/ndtbl/runtime_field_group.hpp | 119 +++++++++++++++++++++----- tests/CMakeLists.txt | 2 + tests/ndtbl_io_t.cpp | 99 +++++++++++++++++++++ 4 files changed, 210 insertions(+), 33 deletions(-) diff --git a/include/ndtbl/io.hpp b/include/ndtbl/io.hpp index 0aef5d1..fe76b99 100644 --- a/include/ndtbl/io.hpp +++ b/include/ndtbl/io.hpp @@ -104,13 +104,15 @@ write_group(const std::string& path, * @brief Write a runtime-erased field group to a binary ndtbl file. * * @tparam Dim Grid dimensionality of the group. + * @tparam Output Scalar output type used by runtime-erased interpolation. * @param path Output file path. * @param group Runtime-erased field group to serialize. * @see RuntimeFieldGroup */ -template +template inline void -write_group(const std::string& path, const RuntimeFieldGroup& group) +write_group(const std::string& path, + const RuntimeFieldGroup& group) { std::ofstream os(path.c_str(), std::ios::binary); if (!os.is_open()) { @@ -146,14 +148,15 @@ read_group_metadata(const std::string& path) * The file dimension must match the compile-time `Dim` argument. * * @tparam Dim Expected grid dimensionality of the file. + * @tparam Output Output type used by runtime-erased interpolation. * @param path Input file path. - * @return Runtime-erased field group with either float or double - * payload storage. + * @return Runtime-erased field group with either float or double payload + * storage and `Output` interpolation results. * @see read_group_metadata * @see RuntimeFieldGroup */ -template -inline RuntimeFieldGroup +template +inline RuntimeFieldGroup read_group(const std::string& path) { std::ifstream is(path.c_str(), std::ios::binary); @@ -175,7 +178,7 @@ read_group(const std::string& path) const std::shared_ptr payload_owner = detail::map_payload_bytes( path, layout.payload_offset, layout.payload_size); - return RuntimeFieldGroup(FieldGroup( + return RuntimeFieldGroup(FieldGroup( grid, metadata.field_names, PayloadView(payload_owner.get(), layout.value_count), @@ -185,7 +188,7 @@ read_group(const std::string& path) // read-only FieldGroup storage instead of copied during load. std::vector values = detail::read_payload(is, layout.value_count); - return RuntimeFieldGroup( + return RuntimeFieldGroup( FieldGroup(grid, metadata.field_names, std::move(values))); #endif } @@ -195,7 +198,7 @@ read_group(const std::string& path) const std::shared_ptr payload_owner = detail::map_payload_bytes( path, layout.payload_offset, layout.payload_size); - return RuntimeFieldGroup(FieldGroup( + return RuntimeFieldGroup(FieldGroup( grid, metadata.field_names, PayloadView(payload_owner.get(), layout.value_count), @@ -205,7 +208,7 @@ read_group(const std::string& path) // read-only FieldGroup storage instead of copied during load. std::vector values = detail::read_payload(is, layout.value_count); - return RuntimeFieldGroup( + return RuntimeFieldGroup( FieldGroup(grid, metadata.field_names, std::move(values))); #endif } diff --git a/include/ndtbl/runtime_field_group.hpp b/include/ndtbl/runtime_field_group.hpp index 605fc9e..d0abab0 100644 --- a/include/ndtbl/runtime_field_group.hpp +++ b/include/ndtbl/runtime_field_group.hpp @@ -9,6 +9,7 @@ #include #include #include +#include #include namespace ndtbl { @@ -22,10 +23,16 @@ write_group_stream(std::ostream& os, const FieldGroup& group); * * The dimensionality remains part of the type. Only the stored scalar payload * type is selected from file metadata at runtime. + * + * @tparam Dim Grid dimensionality of the group. + * @tparam Output Floating-point type produced by interpolation calls. */ -template +template class RuntimeFieldGroup { + static_assert(std::is_floating_point::value, + "RuntimeFieldGroup output type must be floating point"); + public: /** * @brief Construct an empty runtime-erased group handle. @@ -124,14 +131,14 @@ class RuntimeFieldGroup * * @param coordinates Query coordinates in grid axis order. * @param policy Bounds handling behavior for out-of-domain coordinates. - * @return Interpolated field values converted to `double`. + * @return Interpolated field values converted to `Output`. * @see evaluate_all_linear_into */ - std::vector evaluate_all_linear( + std::vector evaluate_all_linear( const std::array& coordinates, bounds_policy policy = bounds_policy::clamp) const { - std::vector values(field_count(), 0.0); + std::vector values(field_count(), Output(0)); evaluate_all_linear_into(coordinates, values.data(), policy); return values; } @@ -147,7 +154,7 @@ class RuntimeFieldGroup */ void evaluate_all_linear_into( const std::array& coordinates, - double* values, + Output* values, bounds_policy policy = bounds_policy::clamp) const { if (!impl_) { @@ -165,14 +172,14 @@ class RuntimeFieldGroup * * @param coordinates Query coordinates in grid axis order. * @param policy Bounds handling behavior for out-of-domain coordinates. - * @return Cubically interpolated field values converted to `double`. + * @return Cubically interpolated field values converted to `Output`. * @see evaluate_all_cubic_into */ - std::vector evaluate_all_cubic( + std::vector evaluate_all_cubic( const std::array& coordinates, bounds_policy policy = bounds_policy::clamp) const { - std::vector values(field_count(), 0.0); + std::vector values(field_count(), Output(0)); evaluate_all_cubic_into(coordinates, values.data(), policy); return values; } @@ -188,7 +195,7 @@ class RuntimeFieldGroup */ void evaluate_all_cubic_into( const std::array& coordinates, - double* values, + Output* values, bounds_policy policy = bounds_policy::clamp) const { if (!impl_) { @@ -222,11 +229,11 @@ class RuntimeFieldGroup virtual std::size_t field_index(const std::string& field_name) const = 0; virtual void evaluate_all_linear_into( const std::array& coordinates, - double* values, + Output* values, bounds_policy policy) const = 0; virtual void evaluate_all_cubic_into( const std::array& coordinates, - double* values, + Output* values, bounds_policy policy) const = 0; virtual void write(std::ostream& os) const = 0; }; @@ -236,7 +243,6 @@ class RuntimeFieldGroup { explicit Model(const FieldGroup& group) : group_(group) - , scratch_(group.field_count(), Value(0)) { } @@ -257,23 +263,25 @@ class RuntimeFieldGroup } void evaluate_all_linear_into(const std::array& coordinates, - double* values, + Output* values, bounds_policy policy) const override { - group_.evaluate_all_linear_into(coordinates, scratch_.data(), policy); - for (std::size_t field = 0; field < scratch_.size(); ++field) { - values[field] = static_cast(scratch_[field]); - } + evaluate_all_linear_into_impl( + coordinates, + values, + policy, + typename std::is_same::type()); } void evaluate_all_cubic_into(const std::array& coordinates, - double* values, + Output* values, bounds_policy policy) const override { - group_.evaluate_all_cubic_into(coordinates, scratch_.data(), policy); - for (std::size_t field = 0; field < scratch_.size(); ++field) { - values[field] = static_cast(scratch_[field]); - } + evaluate_all_cubic_into_impl( + coordinates, + values, + policy, + typename std::is_same::type()); } void write(std::ostream& os) const override @@ -282,7 +290,72 @@ class RuntimeFieldGroup } FieldGroup group_; - mutable std::vector scratch_; + + private: + void evaluate_all_linear_into_impl( + const std::array& coordinates, + Output* values, + bounds_policy policy, + std::true_type) const + { + group_.evaluate_all_linear_into(coordinates, values, policy); + } + + void evaluate_all_linear_into_impl( + const std::array& coordinates, + Output* values, + bounds_policy policy, + std::false_type) const + { + std::vector& scratch = thread_scratch(); + const std::size_t count = group_.field_count(); + if (scratch.size() < count) { + scratch.resize(count); + } + + group_.evaluate_all_linear_into(coordinates, scratch.data(), policy); + copy_scratch_to_output(scratch, count, values); + } + + void evaluate_all_cubic_into_impl( + const std::array& coordinates, + Output* values, + bounds_policy policy, + std::true_type) const + { + group_.evaluate_all_cubic_into(coordinates, values, policy); + } + + void evaluate_all_cubic_into_impl( + const std::array& coordinates, + Output* values, + bounds_policy policy, + std::false_type) const + { + std::vector& scratch = thread_scratch(); + const std::size_t count = group_.field_count(); + if (scratch.size() < count) { + scratch.resize(count); + } + + group_.evaluate_all_cubic_into(coordinates, scratch.data(), policy); + copy_scratch_to_output(scratch, count, values); + } + + static std::vector& thread_scratch() + { + static thread_local std::vector scratch; + return scratch; + } + + static void copy_scratch_to_output(const std::vector& scratch, + std::size_t count, + Output* values) + { + for (std::size_t field = 0; field < count; ++field) { + values[field] = static_cast(scratch[field]); + } + } }; std::shared_ptr impl_; diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index c5c8398..506ee1d 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -15,6 +15,7 @@ endif() # Include Catch2's CMake helper functions (e.g., catch_discover_tests). include(Catch) +find_package(Threads REQUIRED) # Create a single test binary that compiles all test files. add_executable(ndtbl_tests) @@ -31,6 +32,7 @@ target_link_libraries(ndtbl_tests PUBLIC ndtbl Catch2::Catch2WithMain + Threads::Threads ) # allow user to run tests with `make test` or `ctest` diff --git a/tests/ndtbl_io_t.cpp b/tests/ndtbl_io_t.cpp index 2f24571..3c0384d 100644 --- a/tests/ndtbl_io_t.cpp +++ b/tests/ndtbl_io_t.cpp @@ -9,6 +9,7 @@ #include #include #include +#include #include TEST_CASE("typed loader round-trips metadata and float payloads", "[io]") @@ -71,6 +72,104 @@ TEST_CASE("runtime field group can be rewritten after reading", "[io]") std::remove(output_path.c_str()); } +TEST_CASE("runtime field group supports caller-selected output precision", + "[io]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + }; + const ndtbl::Grid<1> grid(axes); + + const ndtbl::FieldGroup float_group( + grid, { "A", "B" }, { 0.0f, 10.0f, 1.0f, 12.0f }); + const ndtbl::RuntimeFieldGroup<1, float> float_runtime(float_group); + std::array float_values = { 0.0f, 0.0f }; + float_runtime.evaluate_all_linear_into({ 0.5 }, float_values.data()); + REQUIRE(float_values[0] == Catch::Approx(0.5f)); + REQUIRE(float_values[1] == Catch::Approx(11.0f)); + + const ndtbl::FieldGroup double_group( + grid, { "A", "B" }, { 0.0, 10.0, 1.0, 12.0 }); + const ndtbl::RuntimeFieldGroup<1, float> converted_runtime(double_group); + converted_runtime.evaluate_all_linear_into({ 0.25 }, float_values.data()); + REQUIRE(float_values[0] == Catch::Approx(0.25f)); + REQUIRE(float_values[1] == Catch::Approx(10.5f)); + + const ndtbl::RuntimeFieldGroup<1, long double> long_double_runtime( + double_group); + const std::vector long_double_values = + long_double_runtime.evaluate_all_linear({ 0.75 }); + REQUIRE(static_cast(long_double_values[0]) == Catch::Approx(0.75)); + REQUIRE(static_cast(long_double_values[1]) == Catch::Approx(11.5)); +} + +TEST_CASE("typed loader supports caller-selected runtime output precision", + "[io]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + }; + const ndtbl::FieldGroup group( + ndtbl::Grid<1>(axes), { "A" }, { 2.0, 4.0 }); + + const std::string path = ndtbl_test::temporary_path(); + ndtbl::write_group(path, group); + + const ndtbl::RuntimeFieldGroup<1, float> loaded = + ndtbl::read_group<1, float>(path); + std::array values = { 0.0f }; + loaded.evaluate_all_linear_into({ 0.5 }, values.data()); + REQUIRE(values[0] == Catch::Approx(3.0f)); + + std::remove(path.c_str()); +} + +TEST_CASE("runtime field group can be evaluated concurrently", "[io]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 15.0, 16), + }; + std::vector payload; + for (std::size_t point = 0; point < axes[0].size(); ++point) { + const float coordinate = static_cast(axes[0].coordinate(point)); + payload.push_back(coordinate); + payload.push_back(10.0f + 2.0f * coordinate); + } + + const ndtbl::FieldGroup group( + ndtbl::Grid<1>(axes), { "A", "B" }, payload); + const ndtbl::RuntimeFieldGroup<1> runtime(group); + + const std::size_t thread_count = 8; + const std::size_t iterations = 1000; + std::vector failures(thread_count, 0); + std::vector threads; + + for (std::size_t thread = 0; thread < thread_count; ++thread) { + threads.push_back(std::thread([&, thread]() { + std::array values = { 0.0, 0.0 }; + for (std::size_t iteration = 0; iteration < iterations; ++iteration) { + const double coordinate = + static_cast((thread + iteration) % 16); + runtime.evaluate_all_linear_into({ coordinate }, values.data()); + if (values[0] != Catch::Approx(coordinate) || + values[1] != Catch::Approx(10.0 + 2.0 * coordinate)) { + failures[thread] = 1; + return; + } + } + })); + } + + for (std::size_t thread = 0; thread < threads.size(); ++thread) { + threads[thread].join(); + } + + for (std::size_t thread = 0; thread < failures.size(); ++thread) { + REQUIRE(failures[thread] == 0); + } +} + TEST_CASE("typed loader rejects mismatched dimensions", "[io]") { const std::array axes = { From 292ff1b36e9a0a65804ab6ed01856a963be5e99f Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 13:51:37 +0200 Subject: [PATCH 05/15] docs(cpp): suggest release build for benchmarks --- README.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.md b/README.md index 40f1719..33f69c7 100644 --- a/README.md +++ b/README.md @@ -159,7 +159,7 @@ runtime-erased combined lookup for representative 2D, 4D, and 6D tables. See Build the benchmark target: ```bash -cmake -B build -Dndtbl_BUILD_BENCHMARKS=ON +cmake -B build -DCMAKE_BUILD_TYPE=Release -Dndtbl_BUILD_BENCHMARKS=ON cmake --build build --target ndtbl_lookup_benchmarks ``` From 7f4f90ec3ef4efbbb6706a7bd68e291e9781a20f Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 14:05:58 +0200 Subject: [PATCH 06/15] refactor(cpp): rename template parameter Value to Stored --- include/ndtbl/detail/binary_io.hpp | 38 ++++++++++++------------ include/ndtbl/field_group.hpp | 42 +++++++++++++-------------- include/ndtbl/io.hpp | 32 ++++++++++---------- include/ndtbl/payload.hpp | 28 +++++++++--------- include/ndtbl/runtime_field_group.hpp | 36 +++++++++++------------ include/ndtbl/types.hpp | 12 ++++---- 6 files changed, 94 insertions(+), 94 deletions(-) diff --git a/include/ndtbl/detail/binary_io.hpp b/include/ndtbl/detail/binary_io.hpp index 9fbb1bf..ea54c5a 100644 --- a/include/ndtbl/detail/binary_io.hpp +++ b/include/ndtbl/detail/binary_io.hpp @@ -133,10 +133,10 @@ read_float_le(std::istream& is) return value; } -template +template struct payload_uint { - typedef typename std::conditional::type type; }; @@ -255,16 +255,16 @@ metadata_size(const GroupMetadata& metadata) /** * @brief Internal implementation for writing raw metadata and payload data. * - * @tparam Value Scalar payload type stored in the payload vector. + * @tparam Stored Scalar payload type stored in the payload vector. * @param os Destination stream in binary mode. * @param metadata File metadata describing the payload layout. * @param payload Point-major interleaved field payload. */ -template +template inline void write_group_stream_impl(std::ostream& os, const GroupMetadata& metadata, - const PayloadView& payload) + const PayloadView& payload) { if (metadata.axes.size() != metadata.dimension) { throw std::invalid_argument( @@ -327,7 +327,7 @@ write_group_stream_impl(std::ostream& os, payload.byte_size()); } else { for (std::size_t index = 0; index < payload.size(); ++index) { - write_float_le::type>( + write_float_le::type>( os, payload[index]); } } @@ -337,18 +337,18 @@ write_group_stream_impl(std::ostream& os, /** * @brief Internal implementation for writing a typed field group. * - * @tparam Value Scalar payload type stored in the group. + * @tparam Stored Scalar payload type stored in the group. * @tparam Dim Grid dimensionality of the group. * @param os Destination stream in binary mode. * @param group Typed field group to serialize. * @see write_group_stream_impl(std::ostream&, const GroupMetadata&, - * const std::vector&) + * const std::vector&) */ -template +template inline void -write_group_stream_impl(std::ostream& os, const FieldGroup& group) +write_group_stream_impl(std::ostream& os, const FieldGroup& group) { - GroupMetadata metadata = { scalar_type_of(), + GroupMetadata metadata = { scalar_type_of(), Dim, group.field_count(), group.point_count(), @@ -358,11 +358,11 @@ write_group_stream_impl(std::ostream& os, const FieldGroup& group) write_group_stream_impl(os, metadata, group.interleaved_values()); } -template +template inline void write_group_stream_impl(std::ostream& os, const GroupMetadata& metadata, - const std::vector& payload) + const std::vector& payload) { write_group_stream_impl(os, metadata, payload_view(payload)); } @@ -500,16 +500,16 @@ read_group_metadata_impl(std::istream& is) /** * @brief Read a contiguous payload block from a binary stream. * - * @tparam Value Scalar payload type to deserialize. + * @tparam Stored Scalar payload type to deserialize. * @param is Source stream positioned at the start of the payload. * @param value_count Number of scalar values to read. * @return Payload vector with `value_count` entries. */ -template -inline std::vector +template +inline std::vector read_payload(std::istream& is, std::size_t value_count) { - std::vector values(value_count); + std::vector values(value_count); if (value_count == 0) { return values; } @@ -517,13 +517,13 @@ read_payload(std::istream& is, std::size_t value_count) if (host_is_little_endian()) { read_bytes(is, reinterpret_cast(values.data()), - values.size() * sizeof(Value)); + values.size() * sizeof(Stored)); return values; } for (std::size_t index = 0; index < value_count; ++index) { values[index] = - read_float_le::type>(is); + read_float_le::type>(is); } return values; } diff --git a/include/ndtbl/field_group.hpp b/include/ndtbl/field_group.hpp index d6e6e98..1515c26 100644 --- a/include/ndtbl/field_group.hpp +++ b/include/ndtbl/field_group.hpp @@ -25,7 +25,7 @@ namespace ndtbl { * axis varies fastest before stepping to the next field tuple. One prepared * interpolation stencil can accumulate all fields together. */ -template +template class FieldGroup { public: @@ -38,8 +38,8 @@ class FieldGroup */ FieldGroup(const Grid& grid, const std::vector& field_names, - const std::vector& interleaved_values) - : FieldGroup(grid, field_names, std::vector(interleaved_values)) + const std::vector& interleaved_values) + : FieldGroup(grid, field_names, std::vector(interleaved_values)) { } @@ -53,7 +53,7 @@ class FieldGroup */ FieldGroup(const Grid& grid, const std::vector& field_names, - std::vector&& interleaved_values) + std::vector&& interleaved_values) : grid_(grid) , field_names_(field_names) { @@ -71,7 +71,7 @@ class FieldGroup */ FieldGroup(const Grid& grid, const std::vector& field_names, - const PayloadView& interleaved_values, + const PayloadView& interleaved_values, std::shared_ptr payload_owner) : grid_(grid) , field_names_(field_names) @@ -117,7 +117,7 @@ class FieldGroup * * @return Point-major interleaved payload view. */ - PayloadView interleaved_values() const { return interleaved_values_; } + PayloadView interleaved_values() const { return interleaved_values_; } /** * @brief Resolve a field name to its local field index. @@ -147,9 +147,9 @@ class FieldGroup * @see evaluate_all_linear(const std::array&) */ template - std::vector evaluate_all(const Stencil& stencil) const + std::vector evaluate_all(const Stencil& stencil) const { - std::vector results(field_count(), Value(0)); + std::vector results(field_count(), Stored(0)); evaluate_all_into(stencil, results.data()); return results; } @@ -164,15 +164,15 @@ class FieldGroup * @see evaluate_all(const Stencil&) */ template - void evaluate_all_into(const Stencil& stencil, Value* results) const + void evaluate_all_into(const Stencil& stencil, Stored* results) const { - std::fill(results, results + field_count(), Value(0)); + std::fill(results, results + field_count(), Stored(0)); for (std::size_t point = 0; point < Stencil::points; ++point) { const double weight = stencil.weight(point); const std::size_t base = stencil.point_index(point) * field_count(); for (std::size_t field = 0; field < field_count(); ++field) { results[field] += - static_cast(weight * interleaved_values_[base + field]); + static_cast(weight * interleaved_values_[base + field]); } } } @@ -186,7 +186,7 @@ class FieldGroup * @param policy Bounds handling behavior for out-of-domain coordinates. * @see evaluate_all(const Stencil&) */ - std::vector evaluate_all_linear( + std::vector evaluate_all_linear( const std::array& coordinates, bounds_policy policy = bounds_policy::clamp) const { @@ -200,11 +200,11 @@ class FieldGroup * @param coordinates Query coordinates in grid axis order. * @param results Output buffer with space for `field_count()` values. * @param policy Bounds handling behavior for out-of-domain coordinates. - * @see evaluate_all_into(const Stencil&, Value*) + * @see evaluate_all_into(const Stencil&, Stored*) */ void evaluate_all_linear_into( const std::array& coordinates, - Value* results, + Stored* results, bounds_policy policy = bounds_policy::clamp) const { evaluate_all_into(grid_.prepare_linear(coordinates, policy), results); @@ -223,7 +223,7 @@ class FieldGroup * @return Cubically interpolated field values in storage order. * @see Grid::prepare_cubic */ - std::vector evaluate_all_cubic( + std::vector evaluate_all_cubic( const std::array& coordinates, bounds_policy policy = bounds_policy::clamp) const { @@ -241,21 +241,21 @@ class FieldGroup */ void evaluate_all_cubic_into( const std::array& coordinates, - Value* results, + Stored* results, bounds_policy policy = bounds_policy::clamp) const { evaluate_all_into(grid_.prepare_cubic(coordinates, policy), results); } private: - void adopt_owned_payload(std::vector&& interleaved_values) + void adopt_owned_payload(std::vector&& interleaved_values) { - const std::shared_ptr> storage = - std::make_shared>(std::move(interleaved_values)); + const std::shared_ptr> storage = + std::make_shared>(std::move(interleaved_values)); const std::uint8_t* const data = storage->empty() ? nullptr : reinterpret_cast(storage->data()); - interleaved_values_ = PayloadView(data, storage->size()); + interleaved_values_ = PayloadView(data, storage->size()); payload_owner_ = std::shared_ptr(storage, data); } @@ -276,7 +276,7 @@ class FieldGroup private: Grid grid_; std::vector field_names_; - PayloadView interleaved_values_; + PayloadView interleaved_values_; std::shared_ptr payload_owner_; }; diff --git a/include/ndtbl/io.hpp b/include/ndtbl/io.hpp index fe76b99..b1a65c1 100644 --- a/include/ndtbl/io.hpp +++ b/include/ndtbl/io.hpp @@ -20,15 +20,15 @@ namespace ndtbl { /** * @brief Write a typed field group to an already opened binary stream. * - * @tparam Value Scalar payload type stored in the group. + * @tparam Stored Scalar payload type stored in the group. * @tparam Dim Grid dimensionality of the group. * @param os Destination stream in binary mode. * @param group Field group to serialize. - * @see write_group(const std::string&, const FieldGroup&) + * @see write_group(const std::string&, const FieldGroup&) */ -template +template inline void -write_group_stream(std::ostream& os, const FieldGroup& group) +write_group_stream(std::ostream& os, const FieldGroup& group) { detail::write_group_stream_impl(os, group); } @@ -40,18 +40,18 @@ write_group_stream(std::ostream& os, const FieldGroup& group) * interleaved point-major payload in row-major axis order, but not a typed * `FieldGroup`. * - * @tparam Value Scalar payload type stored in the payload vector. + * @tparam Stored Scalar payload type stored in the payload vector. * @param os Destination stream in binary mode. * @param metadata Group metadata to encode into the file header. * @param interleaved_values Point-major field payload to serialize in row-major * axis order. - * @see write_group_stream(std::ostream&, const FieldGroup&) + * @see write_group_stream(std::ostream&, const FieldGroup&) */ -template +template inline void write_group_stream(std::ostream& os, const GroupMetadata& metadata, - const std::vector& interleaved_values) + const std::vector& interleaved_values) { detail::write_group_stream_impl(os, metadata, interleaved_values); } @@ -59,15 +59,15 @@ write_group_stream(std::ostream& os, /** * @brief Write a typed field group to a binary ndtbl file. * - * @tparam Value Scalar payload type stored in the group. + * @tparam Stored Scalar payload type stored in the group. * @tparam Dim Grid dimensionality of the group. * @param path Output file path. * @param group Field group to serialize. - * @see write_group_stream(std::ostream&, const FieldGroup&) + * @see write_group_stream(std::ostream&, const FieldGroup&) */ -template +template inline void -write_group(const std::string& path, const FieldGroup& group) +write_group(const std::string& path, const FieldGroup& group) { std::ofstream os(path.c_str(), std::ios::binary); if (!os.is_open()) { @@ -79,19 +79,19 @@ write_group(const std::string& path, const FieldGroup& group) /** * @brief Write a raw ndtbl payload with explicit metadata to a file. * - * @tparam Value Scalar payload type stored in the payload vector. + * @tparam Stored Scalar payload type stored in the payload vector. * @param path Output file path. * @param metadata Group metadata to encode into the file header. * @param interleaved_values Point-major field payload to serialize in row-major * axis order. * @see write_group_stream(std::ostream&, const GroupMetadata&, - * const std::vector&) + * const std::vector&) */ -template +template inline void write_group(const std::string& path, const GroupMetadata& metadata, - const std::vector& interleaved_values) + const std::vector& interleaved_values) { std::ofstream os(path.c_str(), std::ios::binary); if (!os.is_open()) { diff --git a/include/ndtbl/payload.hpp b/include/ndtbl/payload.hpp index e96d8fc..0b286a3 100644 --- a/include/ndtbl/payload.hpp +++ b/include/ndtbl/payload.hpp @@ -11,13 +11,13 @@ namespace ndtbl { /** * @brief Lightweight read-only view over a contiguous ndtbl payload. * - * The payload is addressed in logical `Value` elements, but is stored as bytes + * The payload is addressed in logical `Stored` elements, but is stored as bytes * internally so the view remains valid for both aligned heap-backed storage and * potentially unaligned file-backed memory mappings. * - * @tparam Value Scalar payload type addressed through the view. + * @tparam Stored Scalar payload type addressed through the view. */ -template +template class PayloadView { public: @@ -34,7 +34,7 @@ class PayloadView * @brief Construct a payload view from raw bytes. * * @param data Pointer to the first payload byte. - * @param size Number of logical `Value` entries in the payload. + * @param size Number of logical `Stored` entries in the payload. */ PayloadView(const std::uint8_t* data, std::size_t size) : data_(data) @@ -54,7 +54,7 @@ class PayloadView * * @return Number of occupied bytes. */ - std::size_t byte_size() const { return size_ * sizeof(Value); } + std::size_t byte_size() const { return size_ * sizeof(Stored); } /** * @brief Return the underlying payload bytes. @@ -68,19 +68,19 @@ class PayloadView * * This uses `memcpy` rather than typed pointer dereferencing so the same * implementation stays valid for memory-mapped payloads whose file offset is - * not guaranteed to satisfy `alignof(Value)`. + * not guaranteed to satisfy `alignof(Stored)`. * * @param index Zero-based payload index. * @return Deserialized payload value. */ - Value operator[](std::size_t index) const + Stored operator[](std::size_t index) const { if (index >= size_) { throw std::out_of_range("ndtbl payload index out of range"); } - Value value; - std::memcpy(&value, data_ + index * sizeof(Value), sizeof(Value)); + Stored value; + std::memcpy(&value, data_ + index * sizeof(Stored), sizeof(Stored)); return value; } @@ -92,18 +92,18 @@ class PayloadView /** * @brief Build a read-only payload view over an existing vector. * - * @tparam Value Scalar payload type stored in the vector. + * @tparam Stored Scalar payload type stored in the vector. * @param values Contiguous payload storage to view. * @return Read-only view into `values`. */ -template -inline PayloadView -payload_view(const std::vector& values) +template +inline PayloadView +payload_view(const std::vector& values) { const std::uint8_t* data = values.empty() ? nullptr : reinterpret_cast(values.data()); - return PayloadView(data, values.size()); + return PayloadView(data, values.size()); } } // namespace ndtbl diff --git a/include/ndtbl/runtime_field_group.hpp b/include/ndtbl/runtime_field_group.hpp index d0abab0..51c5dba 100644 --- a/include/ndtbl/runtime_field_group.hpp +++ b/include/ndtbl/runtime_field_group.hpp @@ -14,12 +14,12 @@ namespace ndtbl { -template +template void -write_group_stream(std::ostream& os, const FieldGroup& group); +write_group_stream(std::ostream& os, const FieldGroup& group); /** - * @brief Runtime-erased wrapper around a typed `FieldGroup`. + * @brief Runtime-erased wrapper around a typed `FieldGroup`. * * The dimensionality remains part of the type. Only the stored scalar payload * type is selected from file metadata at runtime. @@ -42,13 +42,13 @@ class RuntimeFieldGroup /** * @brief Construct a runtime-erased wrapper from a typed field group. * - * @tparam Value Scalar payload type stored in the source group. + * @tparam Stored Scalar payload type stored in the source group. * @param group Typed field group to wrap. * @see FieldGroup */ - template - explicit RuntimeFieldGroup(const FieldGroup& group) - : impl_(std::make_shared>(group)) + template + explicit RuntimeFieldGroup(const FieldGroup& group) + : impl_(std::make_shared>(group)) { } @@ -238,17 +238,17 @@ class RuntimeFieldGroup virtual void write(std::ostream& os) const = 0; }; - template + template struct Model : Concept { - explicit Model(const FieldGroup& group) + explicit Model(const FieldGroup& group) : group_(group) { } std::size_t field_count() const override { return group_.field_count(); } - scalar_type value_type() const override { return scalar_type_of(); } + scalar_type value_type() const override { return scalar_type_of(); } std::vector field_names() const override { @@ -270,7 +270,7 @@ class RuntimeFieldGroup coordinates, values, policy, - typename std::is_same::type()); + typename std::is_same::type()); } void evaluate_all_cubic_into(const std::array& coordinates, @@ -281,7 +281,7 @@ class RuntimeFieldGroup coordinates, values, policy, - typename std::is_same::type()); + typename std::is_same::type()); } void write(std::ostream& os) const override @@ -289,7 +289,7 @@ class RuntimeFieldGroup write_group_stream(os, group_); } - FieldGroup group_; + FieldGroup group_; private: void evaluate_all_linear_into_impl( @@ -307,7 +307,7 @@ class RuntimeFieldGroup bounds_policy policy, std::false_type) const { - std::vector& scratch = thread_scratch(); + std::vector& scratch = thread_scratch(); const std::size_t count = group_.field_count(); if (scratch.size() < count) { scratch.resize(count); @@ -332,7 +332,7 @@ class RuntimeFieldGroup bounds_policy policy, std::false_type) const { - std::vector& scratch = thread_scratch(); + std::vector& scratch = thread_scratch(); const std::size_t count = group_.field_count(); if (scratch.size() < count) { scratch.resize(count); @@ -342,13 +342,13 @@ class RuntimeFieldGroup copy_scratch_to_output(scratch, count, values); } - static std::vector& thread_scratch() + static std::vector& thread_scratch() { - static thread_local std::vector scratch; + static thread_local std::vector scratch; return scratch; } - static void copy_scratch_to_output(const std::vector& scratch, + static void copy_scratch_to_output(const std::vector& scratch, std::size_t count, Output* values) { diff --git a/include/ndtbl/types.hpp b/include/ndtbl/types.hpp index 78ce50e..d9ec141 100644 --- a/include/ndtbl/types.hpp +++ b/include/ndtbl/types.hpp @@ -42,19 +42,19 @@ static constexpr std::uint8_t current_format_version = 1u; * * Only `float` and `double` are supported in the current implementation. * - * @tparam Value Supported scalar payload type. + * @tparam Stored Supported scalar payload type. * @return Corresponding ndtbl scalar type tag. */ -template +template constexpr scalar_type scalar_type_of() noexcept { - static_assert(std::is_same::value || - std::is_same::value, + static_assert(std::is_same::value || + std::is_same::value, "ndtbl supports only float and double payloads"); - return std::is_same::value ? scalar_type::float32 - : scalar_type::float64; + return std::is_same::value ? scalar_type::float32 + : scalar_type::float64; } } // namespace ndtbl From af26cea38f78304a7bfec3a7460db03037a9c73f Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 14:09:29 +0200 Subject: [PATCH 07/15] chore: bump clang-format version to v22.1.5 --- .pre-commit-config.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 3e67034..82d1a33 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -1,7 +1,7 @@ repos: # Format C++ code with Clang-Format - automatically applying the changes - repo: https://github.com/pre-commit/mirrors-clang-format - rev: v22.1.4 + rev: v22.1.5 hooks: - id: clang-format args: From d8e63c4f8575cfe9d414e165c1d4c1040813b1a9 Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 14:36:12 +0200 Subject: [PATCH 08/15] refactor(cpp): reject grids whose total point count or strides or grid.point_count() * field_count(), cannot fit into std::size_t --- CMakeLists.txt | 1 + include/ndtbl/detail/binary_io.hpp | 29 +-------------------- include/ndtbl/detail/size_math.hpp | 41 ++++++++++++++++++++++++++++++ include/ndtbl/field_group.hpp | 4 ++- include/ndtbl/grid.hpp | 7 +++-- tests/ndtbl_interpolation_t.cpp | 29 +++++++++++++++++++++ 6 files changed, 80 insertions(+), 31 deletions(-) create mode 100644 include/ndtbl/detail/size_math.hpp diff --git a/CMakeLists.txt b/CMakeLists.txt index b8dbbf7..d542f12 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -39,6 +39,7 @@ target_sources(ndtbl include/ndtbl/payload.hpp include/ndtbl/io.hpp include/ndtbl/ndtbl.hpp + include/ndtbl/detail/size_math.hpp include/ndtbl/detail/binary_io.hpp include/ndtbl/detail/mapped_payload.hpp ) diff --git a/include/ndtbl/detail/binary_io.hpp b/include/ndtbl/detail/binary_io.hpp index ea54c5a..b009076 100644 --- a/include/ndtbl/detail/binary_io.hpp +++ b/include/ndtbl/detail/binary_io.hpp @@ -1,5 +1,6 @@ #pragma once +#include "ndtbl/detail/size_math.hpp" #include "ndtbl/field_group.hpp" #include "ndtbl/metadata.hpp" @@ -178,34 +179,6 @@ read_string(std::istream& is) return value; } -inline std::size_t -checked_multiply_size(std::size_t lhs, std::size_t rhs, const std::string& what) -{ - if (lhs != 0 && rhs > std::numeric_limits::max() / lhs) { - throw std::runtime_error("ndtbl " + what + " exceeds supported size"); - } - return lhs * rhs; -} - -inline std::size_t -checked_add_size(std::size_t lhs, std::size_t rhs, const std::string& what) -{ - if (rhs > std::numeric_limits::max() - lhs) { - throw std::runtime_error("ndtbl " + what + " exceeds supported size"); - } - return lhs + rhs; -} - -inline std::size_t -narrow_u64_to_size(std::uint64_t value, const std::string& what) -{ - if (value > - static_cast(std::numeric_limits::max())) { - throw std::runtime_error("ndtbl " + what + " exceeds supported size"); - } - return static_cast(value); -} - inline void require_zero(std::uint64_t value, const std::string& what) { diff --git a/include/ndtbl/detail/size_math.hpp b/include/ndtbl/detail/size_math.hpp new file mode 100644 index 0000000..88a1cfd --- /dev/null +++ b/include/ndtbl/detail/size_math.hpp @@ -0,0 +1,41 @@ +#pragma once + +#include +#include +#include +#include +#include + +namespace ndtbl { +namespace detail { + +inline std::size_t +checked_multiply_size(std::size_t lhs, std::size_t rhs, const std::string& what) +{ + if (lhs != 0 && rhs > std::numeric_limits::max() / lhs) { + throw std::runtime_error("ndtbl " + what + " exceeds supported size"); + } + return lhs * rhs; +} + +inline std::size_t +checked_add_size(std::size_t lhs, std::size_t rhs, const std::string& what) +{ + if (rhs > std::numeric_limits::max() - lhs) { + throw std::runtime_error("ndtbl " + what + " exceeds supported size"); + } + return lhs + rhs; +} + +inline std::size_t +narrow_u64_to_size(std::uint64_t value, const std::string& what) +{ + if (value > + static_cast(std::numeric_limits::max())) { + throw std::runtime_error("ndtbl " + what + " exceeds supported size"); + } + return static_cast(value); +} + +} // namespace detail +} // namespace ndtbl diff --git a/include/ndtbl/field_group.hpp b/include/ndtbl/field_group.hpp index 1515c26..54e2089 100644 --- a/include/ndtbl/field_group.hpp +++ b/include/ndtbl/field_group.hpp @@ -1,5 +1,6 @@ #pragma once +#include "ndtbl/detail/size_math.hpp" #include "ndtbl/grid.hpp" #include "ndtbl/payload.hpp" @@ -266,7 +267,8 @@ class FieldGroup "field group must contain at least one field"); } - const std::size_t expected_size = grid_.point_count() * field_names_.size(); + const std::size_t expected_size = detail::checked_multiply_size( + grid_.point_count(), field_names_.size(), "field payload value count"); if (interleaved_values_.size() != expected_size) { throw std::invalid_argument( "field payload size does not match grid and field count"); diff --git a/include/ndtbl/grid.hpp b/include/ndtbl/grid.hpp index 0cc930a..852ca75 100644 --- a/include/ndtbl/grid.hpp +++ b/include/ndtbl/grid.hpp @@ -1,6 +1,7 @@ #pragma once #include "ndtbl/axis.hpp" +#include "ndtbl/detail/size_math.hpp" #include #include @@ -135,9 +136,11 @@ class Grid strides_[Dim - 1] = 1; for (std::size_t axis = Dim; axis-- > 0;) { extents_[axis] = axes_[axis].size(); - point_count_ *= extents_[axis]; + point_count_ = detail::checked_multiply_size( + point_count_, extents_[axis], "point count"); if (axis + 1 < Dim) { - strides_[axis] = strides_[axis + 1] * extents_[axis + 1]; + strides_[axis] = detail::checked_multiply_size( + strides_[axis + 1], extents_[axis + 1], "grid stride"); } } } diff --git a/tests/ndtbl_interpolation_t.cpp b/tests/ndtbl_interpolation_t.cpp index 4b4a5c7..7ac791e 100644 --- a/tests/ndtbl_interpolation_t.cpp +++ b/tests/ndtbl_interpolation_t.cpp @@ -7,6 +7,7 @@ #include #include +#include #include #include #include @@ -170,6 +171,34 @@ static_assert(ndtbl::LinearStencil<4>::points == 16, static_assert(ndtbl::CubicStencil<4>::points == 256, "4D cubic interpolation should use 256 table points"); +TEST_CASE("grid rejects shape products exceeding supported size", + "[grid][validation]") +{ + const std::size_t huge_extent = + std::numeric_limits::max() / 2u + 1u; + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + ndtbl::Axis::uniform(0.0, 1.0, huge_extent), + }; + + REQUIRE_THROWS_AS(ndtbl::Grid<2>(axes), std::runtime_error); +} + +TEST_CASE("field group rejects payload shape products exceeding supported size", + "[field_group][validation]") +{ + const std::size_t huge_extent = + std::numeric_limits::max() / 2u + 1u; + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, huge_extent), + }; + const ndtbl::Grid<1> grid(axes); + + REQUIRE_THROWS_AS( + (ndtbl::FieldGroup(grid, { "A", "B" }, std::vector())), + std::runtime_error); +} + TEST_CASE("field group exactly recovers linear fields on uniform axes in 2D", "[field_group][interpolation]") { From 9dedffa289f6b1cde16af5ac8c1c2b1844bb54fd Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 14:51:15 +0200 Subject: [PATCH 09/15] refactor(cpp): Axis::bracket() now rejects NaN, +inf, -inf regardless of bounds_policy --- include/ndtbl/axis.hpp | 5 +++++ tests/ndtbl_interpolation_t.cpp | 20 ++++++++++++++++++++ 2 files changed, 25 insertions(+) diff --git a/include/ndtbl/axis.hpp b/include/ndtbl/axis.hpp index 13b00cc..f836dbf 100644 --- a/include/ndtbl/axis.hpp +++ b/include/ndtbl/axis.hpp @@ -3,6 +3,7 @@ #include "ndtbl/types.hpp" #include +#include #include #include #include @@ -201,6 +202,10 @@ class Axis double value, bounds_policy policy = bounds_policy::clamp) const { + if (!std::isfinite(value)) { + throw std::out_of_range("axis query coordinate must be finite"); + } + if (policy == bounds_policy::throw_error && (value < min_ || value > max_)) { throw std::out_of_range("axis query coordinate outside bounds"); diff --git a/tests/ndtbl_interpolation_t.cpp b/tests/ndtbl_interpolation_t.cpp index 7ac791e..9c9c221 100644 --- a/tests/ndtbl_interpolation_t.cpp +++ b/tests/ndtbl_interpolation_t.cpp @@ -298,6 +298,26 @@ TEST_CASE("single-point axes reject out-of-domain coordinates in throw mode", std::out_of_range); } +TEST_CASE("axes reject non-finite query coordinates", "[axis][interpolation]") +{ + const ndtbl::Axis uniform = ndtbl::Axis::uniform(0.0, 1.0, 3); + const ndtbl::Axis explicit_axis = + ndtbl::Axis::from_coordinates({ 0.0, 0.5, 1.0 }); + const double nan = std::numeric_limits::quiet_NaN(); + const double infinity = std::numeric_limits::infinity(); + + REQUIRE_THROWS_AS(uniform.bracket(nan), std::out_of_range); + REQUIRE_THROWS_AS(uniform.bracket(nan, ndtbl::bounds_policy::throw_error), + std::out_of_range); + REQUIRE_THROWS_AS(explicit_axis.bracket(nan), std::out_of_range); + REQUIRE_THROWS_AS( + explicit_axis.bracket(nan, ndtbl::bounds_policy::throw_error), + std::out_of_range); + + REQUIRE_THROWS_AS(uniform.bracket(infinity), std::out_of_range); + REQUIRE_THROWS_AS(explicit_axis.bracket(-infinity), std::out_of_range); +} + TEST_CASE("field group rejects uniform-axis queries outside the domain", "[field_group][interpolation]") { From ecccefdebdf7ef9212e8407cb3d70cc4af3cc45f Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 15:01:31 +0200 Subject: [PATCH 10/15] refactor(cpp): reject inconsistent metadata before writing --- include/ndtbl/detail/binary_io.hpp | 28 +++++++++++++---- tests/ndtbl_io_t.cpp | 50 ++++++++++++++++++++++++++++++ 2 files changed, 72 insertions(+), 6 deletions(-) diff --git a/include/ndtbl/detail/binary_io.hpp b/include/ndtbl/detail/binary_io.hpp index b009076..b812ab3 100644 --- a/include/ndtbl/detail/binary_io.hpp +++ b/include/ndtbl/detail/binary_io.hpp @@ -225,6 +225,17 @@ metadata_size(const GroupMetadata& metadata) return total; } +inline std::size_t +axis_point_count(const std::vector& axes) +{ + std::size_t point_count = 1; + for (std::size_t axis = 0; axis < axes.size(); ++axis) { + point_count = + checked_multiply_size(point_count, axes[axis].size(), "point count"); + } + return point_count; +} + /** * @brief Internal implementation for writing raw metadata and payload data. * @@ -249,6 +260,16 @@ write_group_stream_impl(std::ostream& os, "ndtbl metadata field count does not match field names"); } + if (metadata.value_type != scalar_type_of()) { + throw std::invalid_argument( + "ndtbl metadata scalar type does not match payload type"); + } + + if (metadata.point_count != axis_point_count(metadata.axes)) { + throw std::invalid_argument( + "ndtbl point count does not match axis extents"); + } + const std::size_t expected_values = checked_multiply_size( metadata.point_count, metadata.field_count, "payload value count"); if (payload.size() != expected_values) { @@ -435,12 +456,7 @@ read_group_layout_impl(std::istream& is) metadata.field_names.push_back(read_string(is)); } - std::size_t expected_point_count = 1; - for (std::size_t axis = 0; axis < metadata.axes.size(); ++axis) { - expected_point_count = checked_multiply_size( - expected_point_count, metadata.axes[axis].size(), "point count"); - } - if (expected_point_count != metadata.point_count) { + if (axis_point_count(metadata.axes) != metadata.point_count) { throw std::runtime_error("ndtbl point count does not match axis extents"); } diff --git a/tests/ndtbl_io_t.cpp b/tests/ndtbl_io_t.cpp index 3c0384d..17b1b55 100644 --- a/tests/ndtbl_io_t.cpp +++ b/tests/ndtbl_io_t.cpp @@ -8,6 +8,7 @@ #include #include #include +#include #include #include #include @@ -124,6 +125,55 @@ TEST_CASE("typed loader supports caller-selected runtime output precision", std::remove(path.c_str()); } +TEST_CASE("raw writer round-trips explicit metadata", "[io]") +{ + const ndtbl::GroupMetadata metadata = { + ndtbl::scalar_type::float32, 1, 2, 3, + { ndtbl::Axis::uniform(0.0, 2.0, 3) }, { "A", "B" } + }; + const std::vector payload = { 0.0f, 10.0f, 1.0f, 11.0f, 2.0f, 12.0f }; + + const std::string path = ndtbl_test::temporary_path(); + ndtbl::write_group(path, metadata, payload); + + const ndtbl::RuntimeFieldGroup<1> loaded = ndtbl::read_group<1>(path); + std::array values = { 0.0, 0.0 }; + loaded.evaluate_all_linear_into({ 0.5 }, values.data()); + REQUIRE(values[0] == Catch::Approx(0.5)); + REQUIRE(values[1] == Catch::Approx(10.5)); + + std::remove(path.c_str()); +} + +TEST_CASE("raw writer rejects metadata scalar type mismatches", "[io]") +{ + const ndtbl::GroupMetadata metadata = { + ndtbl::scalar_type::float64, 1, 1, 2, + { ndtbl::Axis::uniform(0.0, 1.0, 2) }, { "A" } + }; + std::ostringstream os; + + REQUIRE_THROWS_AS( + ndtbl::write_group_stream(os, metadata, std::vector({ 0.0f, 1.0f })), + std::invalid_argument); +} + +TEST_CASE("raw writer rejects point counts mismatching axis extents", "[io]") +{ + const ndtbl::GroupMetadata metadata = { ndtbl::scalar_type::float32, + 2, + 1, + 3, + { ndtbl::Axis::uniform(0.0, 1.0, 2), + ndtbl::Axis::uniform(0.0, 1.0, 2) }, + { "A" } }; + std::ostringstream os; + + REQUIRE_THROWS_AS(ndtbl::write_group_stream( + os, metadata, std::vector({ 0.0f, 1.0f, 2.0f })), + std::invalid_argument); +} + TEST_CASE("runtime field group can be evaluated concurrently", "[io]") { const std::array axes = { From c30ca4db9a367d020b85492439cf042915e70df8 Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Fri, 8 May 2026 17:43:09 +0200 Subject: [PATCH 11/15] refactor(cpp): optimize PayloadView access paths --- include/ndtbl/detail/binary_io.hpp | 2 +- include/ndtbl/field_group.hpp | 56 +++++++++++++++++++++++------- include/ndtbl/payload.hpp | 44 ++++++++++++++++++++--- tests/ndtbl_interpolation_t.cpp | 42 ++++++++++++++++++++++ 4 files changed, 126 insertions(+), 18 deletions(-) diff --git a/include/ndtbl/detail/binary_io.hpp b/include/ndtbl/detail/binary_io.hpp index b812ab3..f51e14e 100644 --- a/include/ndtbl/detail/binary_io.hpp +++ b/include/ndtbl/detail/binary_io.hpp @@ -322,7 +322,7 @@ write_group_stream_impl(std::ostream& os, } else { for (std::size_t index = 0; index < payload.size(); ++index) { write_float_le::type>( - os, payload[index]); + os, payload.unchecked(index)); } } } diff --git a/include/ndtbl/field_group.hpp b/include/ndtbl/field_group.hpp index 54e2089..e1dfd5d 100644 --- a/include/ndtbl/field_group.hpp +++ b/include/ndtbl/field_group.hpp @@ -167,15 +167,16 @@ class FieldGroup template void evaluate_all_into(const Stencil& stencil, Stored* results) const { - std::fill(results, results + field_count(), Stored(0)); - for (std::size_t point = 0; point < Stencil::points; ++point) { - const double weight = stencil.weight(point); - const std::size_t base = stencil.point_index(point) * field_count(); - for (std::size_t field = 0; field < field_count(); ++field) { - results[field] += - static_cast(weight * interleaved_values_[base + field]); - } + const std::size_t fields = field_count(); + const Stored* const typed_values = interleaved_values_.typed_data(); + + std::fill(results, results + fields, Stored(0)); + if (typed_values != nullptr) { + evaluate_all_into_typed(stencil, typed_values, fields, results); + return; } + + evaluate_all_into_bytes(stencil, fields, results); } /** @@ -249,15 +250,46 @@ class FieldGroup } private: + template + void evaluate_all_into_typed(const Stencil& stencil, + const Stored* values, + std::size_t fields, + Stored* results) const + { + for (std::size_t point = 0; point < Stencil::points; ++point) { + const double weight = stencil.weight(point); + const std::size_t base = stencil.point_index(point) * fields; + for (std::size_t field = 0; field < fields; ++field) { + const Stored value = values[base + field]; + results[field] += static_cast(weight * value); + } + } + } + + template + void evaluate_all_into_bytes(const Stencil& stencil, + std::size_t fields, + Stored* results) const + { + for (std::size_t point = 0; point < Stencil::points; ++point) { + const double weight = stencil.weight(point); + const std::size_t base = stencil.point_index(point) * fields; + for (std::size_t field = 0; field < fields; ++field) { + const Stored value = interleaved_values_.unchecked(base + field); + results[field] += static_cast(weight * value); + } + } + } + void adopt_owned_payload(std::vector&& interleaved_values) { const std::shared_ptr> storage = std::make_shared>(std::move(interleaved_values)); - const std::uint8_t* const data = - storage->empty() ? nullptr - : reinterpret_cast(storage->data()); + const Stored* const data = storage->empty() ? nullptr : storage->data(); interleaved_values_ = PayloadView(data, storage->size()); - payload_owner_ = std::shared_ptr(storage, data); + payload_owner_ = std::shared_ptr( + storage, + data == nullptr ? nullptr : reinterpret_cast(data)); } void validate_payload_shape() const diff --git a/include/ndtbl/payload.hpp b/include/ndtbl/payload.hpp index 0b286a3..8dc5268 100644 --- a/include/ndtbl/payload.hpp +++ b/include/ndtbl/payload.hpp @@ -26,6 +26,7 @@ class PayloadView */ PayloadView() : data_(nullptr) + , typed_data_(nullptr) , size_(0) { } @@ -38,6 +39,21 @@ class PayloadView */ PayloadView(const std::uint8_t* data, std::size_t size) : data_(data) + , typed_data_(nullptr) + , size_(size) + { + } + + /** + * @brief Construct a payload view from typed scalar storage. + * + * @param data Pointer to the first typed payload value. + * @param size Number of logical `Stored` entries in the payload. + */ + PayloadView(const Stored* data, std::size_t size) + : data_(data == nullptr ? nullptr + : reinterpret_cast(data)) + , typed_data_(data) , size_(size) { } @@ -64,7 +80,15 @@ class PayloadView const std::uint8_t* byte_data() const { return data_; } /** - * @brief Read one payload value by index. + * @brief Return typed payload storage when the view was built from it. + * + * @return Typed pointer, or `nullptr` when direct typed loads would not be + * well-defined. + */ + const Stored* typed_data() const { return typed_data_; } + + /** + * @brief Read one payload value by index with bounds checking. * * This uses `memcpy` rather than typed pointer dereferencing so the same * implementation stays valid for memory-mapped payloads whose file offset is @@ -73,12 +97,23 @@ class PayloadView * @param index Zero-based payload index. * @return Deserialized payload value. */ - Stored operator[](std::size_t index) const + Stored at(std::size_t index) const { if (index >= size_) { throw std::out_of_range("ndtbl payload index out of range"); } + return unchecked(index); + } + + /** + * @brief Read one payload value without bounds checking. + * + * @param index Zero-based payload index known to be valid. + * @return Deserialized payload value. + */ + Stored unchecked(std::size_t index) const + { Stored value; std::memcpy(&value, data_ + index * sizeof(Stored), sizeof(Stored)); return value; @@ -86,6 +121,7 @@ class PayloadView private: const std::uint8_t* data_; + const Stored* typed_data_; std::size_t size_; }; @@ -100,9 +136,7 @@ template inline PayloadView payload_view(const std::vector& values) { - const std::uint8_t* data = - values.empty() ? nullptr - : reinterpret_cast(values.data()); + const Stored* data = values.empty() ? nullptr : values.data(); return PayloadView(data, values.size()); } diff --git a/tests/ndtbl_interpolation_t.cpp b/tests/ndtbl_interpolation_t.cpp index 9c9c221..840badc 100644 --- a/tests/ndtbl_interpolation_t.cpp +++ b/tests/ndtbl_interpolation_t.cpp @@ -6,8 +6,11 @@ #include #include +#include #include +#include #include +#include #include #include #include @@ -199,6 +202,45 @@ TEST_CASE("field group rejects payload shape products exceeding supported size", std::runtime_error); } +TEST_CASE("payload view keeps checked access separate from typed fast access", + "[payload]") +{ + const std::vector payload = { 1.0, 2.0, 3.0 }; + const ndtbl::PayloadView view = ndtbl::payload_view(payload); + + REQUIRE(view.typed_data() == payload.data()); + REQUIRE(view.at(1) == Catch::Approx(2.0)); + REQUIRE(view.unchecked(2) == Catch::Approx(3.0)); + REQUIRE_THROWS_AS(view.at(payload.size()), std::out_of_range); +} + +TEST_CASE("field group evaluates unaligned byte-backed payloads", + "[field_group][payload][interpolation]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + }; + const std::vector payload = { 1.0, 3.0 }; + + std::shared_ptr> storage = + std::make_shared>( + payload.size() * sizeof(double) + 1u); + std::uint8_t* const unaligned_data = storage->data() + 1u; + std::memcpy(unaligned_data, payload.data(), payload.size() * sizeof(double)); + + const ndtbl::PayloadView view(unaligned_data, payload.size()); + REQUIRE(view.typed_data() == nullptr); + + const std::shared_ptr owner(storage, unaligned_data); + const ndtbl::FieldGroup group( + ndtbl::Grid<1>(axes), { "A" }, view, owner); + + std::array values = { 0.0 }; + group.evaluate_all_linear_into({ 0.25 }, values.data()); + + REQUIRE(values[0] == Catch::Approx(1.5)); +} + TEST_CASE("field group exactly recovers linear fields on uniform axes in 2D", "[field_group][interpolation]") { From 0682830a847fe0aba47d5960dde1b2b01a80183e Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Tue, 12 May 2026 07:22:55 +0200 Subject: [PATCH 12/15] refactor(cpp): change order of template parameters, Dim first --- benchmarks/lookup_benchmarks.cpp | 6 +++--- include/ndtbl/detail/binary_io.hpp | 2 +- include/ndtbl/field_group.hpp | 5 ++++- include/ndtbl/grid.hpp | 9 +++++++++ include/ndtbl/io.hpp | 18 +++++++++--------- include/ndtbl/runtime_field_group.hpp | 10 +++++----- tests/ndtbl_interpolation_t.cpp | 26 +++++++++++++------------- tests/ndtbl_io_t.cpp | 22 +++++++++++----------- 8 files changed, 55 insertions(+), 43 deletions(-) diff --git a/benchmarks/lookup_benchmarks.cpp b/benchmarks/lookup_benchmarks.cpp index eaa35f6..39e0e48 100644 --- a/benchmarks/lookup_benchmarks.cpp +++ b/benchmarks/lookup_benchmarks.cpp @@ -50,13 +50,13 @@ template struct LookupContext { ndtbl::Grid grid; - ndtbl::FieldGroup group; + ndtbl::FieldGroup group; ndtbl::RuntimeFieldGroup runtime_group; std::vector> queries; ndtbl::LinearStencil prepared; LookupContext(const ndtbl::Grid& grid_in, - const ndtbl::FieldGroup& group_in, + const ndtbl::FieldGroup& group_in, const std::vector>& queries_in) : grid(grid_in) , group(group_in) @@ -229,7 +229,7 @@ make_context(std::size_t extent, const std::array shape = filled_shape(extent); const std::array axes = make_axes(shape, axis_kind); const ndtbl::Grid grid(axes); - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup group( grid, make_field_names(field_count), make_payload(grid, field_count)); return LookupContext(grid, group, make_queries(axes, query_count)); } diff --git a/include/ndtbl/detail/binary_io.hpp b/include/ndtbl/detail/binary_io.hpp index f51e14e..c7cda48 100644 --- a/include/ndtbl/detail/binary_io.hpp +++ b/include/ndtbl/detail/binary_io.hpp @@ -340,7 +340,7 @@ write_group_stream_impl(std::ostream& os, */ template inline void -write_group_stream_impl(std::ostream& os, const FieldGroup& group) +write_group_stream_impl(std::ostream& os, const FieldGroup& group) { GroupMetadata metadata = { scalar_type_of(), Dim, diff --git a/include/ndtbl/field_group.hpp b/include/ndtbl/field_group.hpp index e1dfd5d..3bbd188 100644 --- a/include/ndtbl/field_group.hpp +++ b/include/ndtbl/field_group.hpp @@ -25,8 +25,11 @@ namespace ndtbl { * `point0.field0, point0.field1, ..., point1.field0, ...` where the last grid * axis varies fastest before stepping to the next field tuple. One prepared * interpolation stencil can accumulate all fields together. + * + * @tparam Dim Grid dimensionality of the group. + * @tparam Stored Scalar payload type stored in the group. */ -template +template class FieldGroup { public: diff --git a/include/ndtbl/grid.hpp b/include/ndtbl/grid.hpp index 852ca75..441e275 100644 --- a/include/ndtbl/grid.hpp +++ b/include/ndtbl/grid.hpp @@ -31,6 +31,9 @@ pow_size(std::size_t base, std::size_t exponent) * A stencil stores linearized point indices and corresponding weights for one * interpolation point so that multiple fields on the same grid can reuse the * same lookup work. + * + * @tparam Dim Number of dimensions. + * @tparam PointsPerAxis Number of interpolation points along each axis. */ template class TensorStencil @@ -107,6 +110,12 @@ using CubicStencil = TensorStencil; /** * @brief N-dimensional grid metadata with stride and query preparation logic. + + * A grid is defined by one axis descriptor per dimension. It provides the + * logic to prepare interpolation stencils for query points and to validate + * compatibility with field groups. + * + * @tparam Dim Number of dimensions. */ template class Grid diff --git a/include/ndtbl/io.hpp b/include/ndtbl/io.hpp index b1a65c1..3dba326 100644 --- a/include/ndtbl/io.hpp +++ b/include/ndtbl/io.hpp @@ -24,11 +24,11 @@ namespace ndtbl { * @tparam Dim Grid dimensionality of the group. * @param os Destination stream in binary mode. * @param group Field group to serialize. - * @see write_group(const std::string&, const FieldGroup&) + * @see write_group(const std::string&, const FieldGroup&) */ template inline void -write_group_stream(std::ostream& os, const FieldGroup& group) +write_group_stream(std::ostream& os, const FieldGroup& group) { detail::write_group_stream_impl(os, group); } @@ -45,7 +45,7 @@ write_group_stream(std::ostream& os, const FieldGroup& group) * @param metadata Group metadata to encode into the file header. * @param interleaved_values Point-major field payload to serialize in row-major * axis order. - * @see write_group_stream(std::ostream&, const FieldGroup&) + * @see write_group_stream(std::ostream&, const FieldGroup&) */ template inline void @@ -63,11 +63,11 @@ write_group_stream(std::ostream& os, * @tparam Dim Grid dimensionality of the group. * @param path Output file path. * @param group Field group to serialize. - * @see write_group_stream(std::ostream&, const FieldGroup&) + * @see write_group_stream(std::ostream&, const FieldGroup&) */ template inline void -write_group(const std::string& path, const FieldGroup& group) +write_group(const std::string& path, const FieldGroup& group) { std::ofstream os(path.c_str(), std::ios::binary); if (!os.is_open()) { @@ -178,7 +178,7 @@ read_group(const std::string& path) const std::shared_ptr payload_owner = detail::map_payload_bytes( path, layout.payload_offset, layout.payload_size); - return RuntimeFieldGroup(FieldGroup( + return RuntimeFieldGroup(FieldGroup( grid, metadata.field_names, PayloadView(payload_owner.get(), layout.value_count), @@ -189,7 +189,7 @@ read_group(const std::string& path) std::vector values = detail::read_payload(is, layout.value_count); return RuntimeFieldGroup( - FieldGroup(grid, metadata.field_names, std::move(values))); + FieldGroup(grid, metadata.field_names, std::move(values))); #endif } @@ -198,7 +198,7 @@ read_group(const std::string& path) const std::shared_ptr payload_owner = detail::map_payload_bytes( path, layout.payload_offset, layout.payload_size); - return RuntimeFieldGroup(FieldGroup( + return RuntimeFieldGroup(FieldGroup( grid, metadata.field_names, PayloadView(payload_owner.get(), layout.value_count), @@ -209,7 +209,7 @@ read_group(const std::string& path) std::vector values = detail::read_payload(is, layout.value_count); return RuntimeFieldGroup( - FieldGroup(grid, metadata.field_names, std::move(values))); + FieldGroup(grid, metadata.field_names, std::move(values))); #endif } diff --git a/include/ndtbl/runtime_field_group.hpp b/include/ndtbl/runtime_field_group.hpp index 51c5dba..fb2785d 100644 --- a/include/ndtbl/runtime_field_group.hpp +++ b/include/ndtbl/runtime_field_group.hpp @@ -16,10 +16,10 @@ namespace ndtbl { template void -write_group_stream(std::ostream& os, const FieldGroup& group); +write_group_stream(std::ostream& os, const FieldGroup& group); /** - * @brief Runtime-erased wrapper around a typed `FieldGroup`. + * @brief Runtime-erased wrapper around a typed `FieldGroup`. * * The dimensionality remains part of the type. Only the stored scalar payload * type is selected from file metadata at runtime. @@ -47,7 +47,7 @@ class RuntimeFieldGroup * @see FieldGroup */ template - explicit RuntimeFieldGroup(const FieldGroup& group) + explicit RuntimeFieldGroup(const FieldGroup& group) : impl_(std::make_shared>(group)) { } @@ -241,7 +241,7 @@ class RuntimeFieldGroup template struct Model : Concept { - explicit Model(const FieldGroup& group) + explicit Model(const FieldGroup& group) : group_(group) { } @@ -289,7 +289,7 @@ class RuntimeFieldGroup write_group_stream(os, group_); } - FieldGroup group_; + FieldGroup group_; private: void evaluate_all_linear_into_impl( diff --git a/tests/ndtbl_interpolation_t.cpp b/tests/ndtbl_interpolation_t.cpp index 840badc..97420a3 100644 --- a/tests/ndtbl_interpolation_t.cpp +++ b/tests/ndtbl_interpolation_t.cpp @@ -49,7 +49,7 @@ require_linear_recovery(const std::array& axes, const double intercept0 = 1.25; const double intercept1 = -0.75; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup group( ndtbl::Grid(axes), { "A", "B" }, ndtbl_test::build_linear_payload( @@ -83,7 +83,7 @@ require_clamped_linear_recovery(const std::array& axes, const double intercept0 = 1.25; const double intercept1 = -0.75; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup group( ndtbl::Grid(axes), { "A", "B" }, ndtbl_test::build_linear_payload( @@ -108,7 +108,7 @@ require_loaded_linear_recovery(const std::array& axes, const double intercept0 = 1.25; const double intercept1 = -0.75; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup group( ndtbl::Grid(axes), { "A", "B" }, ndtbl_test::build_linear_payload( @@ -198,7 +198,7 @@ TEST_CASE("field group rejects payload shape products exceeding supported size", const ndtbl::Grid<1> grid(axes); REQUIRE_THROWS_AS( - (ndtbl::FieldGroup(grid, { "A", "B" }, std::vector())), + (ndtbl::FieldGroup<1, double>(grid, { "A", "B" }, std::vector())), std::runtime_error); } @@ -232,7 +232,7 @@ TEST_CASE("field group evaluates unaligned byte-backed payloads", REQUIRE(view.typed_data() == nullptr); const std::shared_ptr owner(storage, unaligned_data); - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "A" }, view, owner); std::array values = { 0.0 }; @@ -367,7 +367,7 @@ TEST_CASE("field group rejects uniform-axis queries outside the domain", ndtbl::Axis::uniform(-1.0, 2.0, 4), ndtbl::Axis::uniform(0.0, 6.0, 4), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<2, double> group( ndtbl::Grid<2>(axes), { "A", "B" }, ndtbl_test::build_linear_payload( @@ -392,7 +392,7 @@ TEST_CASE("field group rejects explicit-axis queries outside the domain", ndtbl::Axis::from_coordinates({ -2.0, -0.5, 1.0, 3.5 }), ndtbl::Axis::from_coordinates({ 1.0, 1.75, 4.0, 6.0 }), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<2, double> group( ndtbl::Grid<2>(axes), { "A", "B" }, ndtbl_test::build_linear_payload( @@ -417,7 +417,7 @@ TEST_CASE("runtime field group forwards throw bounds policy", ndtbl::Axis::uniform(-1.0, 2.0, 4), ndtbl::Axis::uniform(0.0, 6.0, 4), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<2, double> group( ndtbl::Grid<2>(axes), { "A", "B" }, ndtbl_test::build_linear_payload( @@ -466,7 +466,7 @@ TEST_CASE("field group cubic interpolation exactly recovers 1D cubic data", ndtbl::Axis::uniform(-2.0, 3.0, 6), }; const ndtbl::Grid<1> grid(axes); - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( grid, { "cubic" }, build_single_field_payload(axes, [](const std::array& coords) { @@ -492,7 +492,7 @@ TEST_CASE("field group cubic interpolation exactly recovers explicit-axis data", ndtbl::Axis::from_coordinates({ -1.0, -0.25, 0.5, 1.25, 2.5 }), }; const ndtbl::Grid<2> grid(axes); - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<2, double> group( grid, { "cubic" }, build_single_field_payload(axes, [](const std::array& coords) { @@ -511,7 +511,7 @@ TEST_CASE("field group cubic interpolation handles boundary windows", const std::array axes = { ndtbl::Axis::uniform(0.0, 4.0, 5), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "cubic" }, build_single_field_payload(axes, [](const std::array& coords) { @@ -543,7 +543,7 @@ TEST_CASE("cubic interpolation respects throw bounds policy", const std::array axes = { ndtbl::Axis::from_coordinates({ -1.0, 0.0, 1.0, 2.0 }), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "cubic" }, build_single_field_payload(axes, [](const std::array& coords) { @@ -568,7 +568,7 @@ TEST_CASE("runtime field groups expose explicit cubic interpolation", const std::array axes = { ndtbl::Axis::uniform(-2.0, 3.0, 6), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "cubic" }, build_single_field_payload(axes, [](const std::array& coords) { diff --git a/tests/ndtbl_io_t.cpp b/tests/ndtbl_io_t.cpp index 17b1b55..cb6d0a0 100644 --- a/tests/ndtbl_io_t.cpp +++ b/tests/ndtbl_io_t.cpp @@ -19,7 +19,7 @@ TEST_CASE("typed loader round-trips metadata and float payloads", "[io]") ndtbl::Axis::uniform(0.0, 1.0, 2), ndtbl::Axis::uniform(0.0, 1.0, 2), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<2, float> group( ndtbl::Grid<2>(axes), { "A", "B" }, { 0.0f, 10.0f, 1.0f, 11.0f, 2.0f, 12.0f, 3.0f, 13.0f }); @@ -50,7 +50,7 @@ TEST_CASE("runtime field group can be rewritten after reading", "[io]") ndtbl::Axis::uniform(0.0, 1.0, 2), ndtbl::Axis::uniform(0.0, 1.0, 2), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<2, double> group( ndtbl::Grid<2>(axes), { "A", "B" }, { 0.0, 10.0, 1.0, 11.0, 2.0, 12.0, 3.0, 13.0 }); @@ -81,7 +81,7 @@ TEST_CASE("runtime field group supports caller-selected output precision", }; const ndtbl::Grid<1> grid(axes); - const ndtbl::FieldGroup float_group( + const ndtbl::FieldGroup<1, float> float_group( grid, { "A", "B" }, { 0.0f, 10.0f, 1.0f, 12.0f }); const ndtbl::RuntimeFieldGroup<1, float> float_runtime(float_group); std::array float_values = { 0.0f, 0.0f }; @@ -89,7 +89,7 @@ TEST_CASE("runtime field group supports caller-selected output precision", REQUIRE(float_values[0] == Catch::Approx(0.5f)); REQUIRE(float_values[1] == Catch::Approx(11.0f)); - const ndtbl::FieldGroup double_group( + const ndtbl::FieldGroup<1, double> double_group( grid, { "A", "B" }, { 0.0, 10.0, 1.0, 12.0 }); const ndtbl::RuntimeFieldGroup<1, float> converted_runtime(double_group); converted_runtime.evaluate_all_linear_into({ 0.25 }, float_values.data()); @@ -110,7 +110,7 @@ TEST_CASE("typed loader supports caller-selected runtime output precision", const std::array axes = { ndtbl::Axis::uniform(0.0, 1.0, 2), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "A" }, { 2.0, 4.0 }); const std::string path = ndtbl_test::temporary_path(); @@ -186,7 +186,7 @@ TEST_CASE("runtime field group can be evaluated concurrently", "[io]") payload.push_back(10.0f + 2.0f * coordinate); } - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, float> group( ndtbl::Grid<1>(axes), { "A", "B" }, payload); const ndtbl::RuntimeFieldGroup<1> runtime(group); @@ -225,7 +225,7 @@ TEST_CASE("typed loader rejects mismatched dimensions", "[io]") const std::array axes = { ndtbl::Axis::uniform(0.0, 1.0, 2), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "A" }, { 0.0, 1.0 }); const std::string path = ndtbl_test::temporary_path(); @@ -241,7 +241,7 @@ TEST_CASE("typed loader rejects truncated payload files", "[io]") const std::array axes = { ndtbl::Axis::uniform(0.0, 1.0, 2), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "A" }, { 0.0, 1.0 }); const std::string path = ndtbl_test::temporary_path(); @@ -262,7 +262,7 @@ TEST_CASE("writer produces the documented little-endian layout", "[io]") const std::array axes = { ndtbl::Axis::uniform(2.0, 2.0, 1), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, float> group( ndtbl::Grid<1>(axes), { "A" }, { 1.5f }); const std::string path = ndtbl_test::temporary_path(); @@ -298,7 +298,7 @@ TEST_CASE("typed loader rejects nonzero reserved header fields", "[io]") const std::array axes = { ndtbl::Axis::uniform(0.0, 1.0, 2), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "A" }, { 0.0, 1.0 }); const std::string path = ndtbl_test::temporary_path(); @@ -319,7 +319,7 @@ TEST_CASE("typed loader rejects mismatched payload offsets", "[io]") const std::array axes = { ndtbl::Axis::uniform(0.0, 1.0, 2), }; - const ndtbl::FieldGroup group( + const ndtbl::FieldGroup<1, double> group( ndtbl::Grid<1>(axes), { "A" }, { 0.0, 1.0 }); const std::string path = ndtbl_test::temporary_path(); From 0428953a6c668ac7158618a64bb40963c84b6e8d Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Tue, 12 May 2026 12:55:49 +0200 Subject: [PATCH 13/15] ci: create GH release also for Python tag py-vX.X.X, including the source dist --- .github/workflows/release-pypi.yml | 21 +++++++++++++++++++++ 1 file changed, 21 insertions(+) diff --git a/.github/workflows/release-pypi.yml b/.github/workflows/release-pypi.yml index c0869c2..ab5e07a 100644 --- a/.github/workflows/release-pypi.yml +++ b/.github/workflows/release-pypi.yml @@ -138,3 +138,24 @@ jobs: - name: Publish distributions to PyPI run: | uv publish dist/* + + github-release: + name: Create GitHub Release + if: ${{ github.event_name == 'push' }} + needs: [build, publish-pypi] + runs-on: ubuntu-latest + permissions: + contents: write + + steps: + - name: Download built distributions + uses: actions/download-artifact@v8 + with: + name: python-package-distributions + path: dist/ + + - name: Create GitHub Release + uses: softprops/action-gh-release@v3 + with: + generate_release_notes: true + files: dist/* From 3ad70f32a1a5163c40cec51f875c8de8319eb497 Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Tue, 12 May 2026 13:14:20 +0200 Subject: [PATCH 14/15] docs: add Doxygen comments --- include/ndtbl/axis.hpp | 3 ++ include/ndtbl/grid.hpp | 67 ++++++++++++++++++++++++++++++++++++----- include/ndtbl/types.hpp | 6 ++++ 3 files changed, 68 insertions(+), 8 deletions(-) diff --git a/include/ndtbl/axis.hpp b/include/ndtbl/axis.hpp index f836dbf..ec7ecef 100644 --- a/include/ndtbl/axis.hpp +++ b/include/ndtbl/axis.hpp @@ -21,6 +21,9 @@ namespace ndtbl { class Axis { public: + /** + * @brief Construct a single-point uniform axis at coordinate zero. + */ Axis() : kind_(axis_kind::uniform) , size_(1) diff --git a/include/ndtbl/grid.hpp b/include/ndtbl/grid.hpp index 441e275..c952dcc 100644 --- a/include/ndtbl/grid.hpp +++ b/include/ndtbl/grid.hpp @@ -13,6 +13,16 @@ namespace ndtbl { namespace detail { +/** + * @brief Compute an integer power for `std::size_t` values. + * + * Computes `base` raised to the power `exponent` by repeated multiplication. + * Intended for small compile-time constants such as stencil sizes. + * + * @param base Base value. + * @param exponent Non-negative exponent. + * @return `base` raised to `exponent`. + */ constexpr std::size_t pow_size(std::size_t base, std::size_t exponent) { @@ -102,8 +112,8 @@ template using LinearStencil = TensorStencil; /** - * @brief Local tensor-product cubic interpolation stencil with four points per - * axis. + * @brief Local tensor-product cubic Lagrange interpolation stencil with four + * points per axis. */ template using CubicStencil = TensorStencil; @@ -126,6 +136,12 @@ class Grid /// Number of dimensions. static constexpr std::size_t dimensions = Dim; + /** + * @brief Construct an empty grid descriptor. + * + * The default-constructed grid has zero points and is mainly useful as a + * placeholder before assigning a fully constructed grid. + */ Grid() { point_count_ = 0; @@ -272,16 +288,34 @@ class Grid } /** - * @brief Precompute a local tensor-product cubic interpolation stencil for - * one query point. + * @brief Precompute a local tensor-product cubic Lagrange interpolation + * stencil for one query point. + * + * This method builds a local interpolation stencil using four support points + * per axis. Along each axis, the one-dimensional weights are the cubic + * Lagrange basis weights for the selected four axis coordinates: + * + * L_i(x) = prod_{j != i} (x - x_j) / (x_i - x_j) * - * Cubic interpolation uses four support points per axis and therefore `4^Dim` - * table values. It is intended as an experimental higher-order path. The - * linear path remains the default. + * The multidimensional stencil is formed as the tensor product of these + * one-dimensional Lagrange weights. Consequently, the stencil contains + * `4^Dim` table values. + * + * The basis weights are computed from the actual axis coordinates, so + * non-uniform axes are supported. Near domain boundaries, the four-point + * support window is shifted/clamped so that it remains inside the table. This + * makes the cubic stencil one-sided near boundaries. + * + * Cubic interpolation is intended as an experimental higher-order path. The + * linear path remains the recommended default for performance-critical + * high-dimensional table lookup. * * @param coordinates Query coordinates in axis order. * @param policy Bounds handling behavior for out-of-domain coordinates. - * @return Cubic stencil containing point indices and weights. + * @return Cubic stencil containing flattened table point indices and + * tensor-product weights. + * + * @throws std::invalid_argument if any axis contains fewer than four points. */ CubicStencil prepare_cubic( const std::array& coordinates, @@ -333,6 +367,23 @@ class Grid } private: + /** + * @brief Compute one-dimensional cubic Lagrange basis weights for one axis. + * + * Given four support point indices on the selected axis and a query + * coordinate, this function computes the four Lagrange basis weights + * + * w_i = prod_{j != i} (x - x_j) / (x_i - x_j) + * + * using the actual axis coordinates `x_i`. The resulting weights define the + * unique cubic polynomial, restricted to this local four-point stencil, that + * interpolates the four support values exactly. + * + * @param axis Axis for which the weights are computed. + * @param coordinate Query coordinate on that axis. + * @param indices Four support point indices on the selected axis. + * @param weights Output array receiving the four Lagrange basis weights. + */ void cubic_weights(std::size_t axis, double coordinate, const std::array& indices, diff --git a/include/ndtbl/types.hpp b/include/ndtbl/types.hpp index d9ec141..a58bd24 100644 --- a/include/ndtbl/types.hpp +++ b/include/ndtbl/types.hpp @@ -10,7 +10,9 @@ namespace ndtbl { */ enum class axis_kind : std::uint8_t { + /// Axis coordinates are generated from minimum, maximum, and size. uniform = 1, + /// Axis coordinates are stored explicitly in ascending order. explicit_coordinates = 2 }; @@ -19,7 +21,9 @@ enum class axis_kind : std::uint8_t */ enum class scalar_type : std::uint8_t { + /// IEEE-754 single-precision floating-point payload values. float32 = 1, + /// IEEE-754 double-precision floating-point payload values. float64 = 2 }; @@ -28,7 +32,9 @@ enum class scalar_type : std::uint8_t */ enum class bounds_policy : std::uint8_t { + /// Clamp out-of-domain interpolation queries to the table bounds. clamp = 1, + /// Reject out-of-domain interpolation queries with an exception. throw_error = 2 }; From 5a7fa10396352f00973f99587ff6cd66bcab2e5d Mon Sep 17 00:00:00 2001 From: Thomas Isensee <26852629+thomasisensee@users.noreply.github.com> Date: Tue, 12 May 2026 15:06:09 +0200 Subject: [PATCH 15/15] refactor(cpp): typed known-storage loader read_field_group, and runtime-erased loader read_runtime_field_group --- README.md | 9 +- doc/cppapi.rst | 4 + include/ndtbl/field_group.hpp | 147 +++++++++++++++++++++++--- include/ndtbl/io.hpp | 61 ++++++++++- include/ndtbl/runtime_field_group.hpp | 80 +------------- tests/ndtbl_interpolation_t.cpp | 47 +++++++- tests/ndtbl_io_t.cpp | 81 ++++++++++++-- 7 files changed, 321 insertions(+), 108 deletions(-) diff --git a/README.md b/README.md index 33f69c7..d1bd218 100644 --- a/README.md +++ b/README.md @@ -67,10 +67,11 @@ Relevant CMake options: - `ndtbl_BUILD_DOCS`: build the documentation, default `ON` for top-level builds - `ndtbl_ENABLE_MMAP`: enable POSIX-only `mmap`-backed payload reads, default `OFF` -When `ndtbl_ENABLE_MMAP=OFF` (the default), `read_group()` reads payload data -into owned heap storage. When `ndtbl_ENABLE_MMAP=ON`, supported POSIX builds -use read-only memory mapping instead, which can reduce heap usage for large -tables and enables shared memory access in multi-process environments. +When `ndtbl_ENABLE_MMAP=OFF` (the default), `read_field_group()` and +`read_runtime_field_group()` read payload data into owned heap storage. When +`ndtbl_ENABLE_MMAP=ON`, supported POSIX builds use read-only memory mapping +instead, which can reduce heap usage for large tables and enables shared memory +access in multi-process environments. If you want to install the C++ headers and CMake package metadata: diff --git a/doc/cppapi.rst b/doc/cppapi.rst index ff63548..4af8ff9 100644 --- a/doc/cppapi.rst +++ b/doc/cppapi.rst @@ -35,3 +35,7 @@ API reference .. doxygenenum:: ndtbl::bounds_policy .. doxygenfunction:: ndtbl::read_group_metadata(const std::string &path) + +.. doxygenfunction:: ndtbl::read_field_group(const std::string &path) + +.. doxygenfunction:: ndtbl::read_runtime_field_group(const std::string &path) diff --git a/include/ndtbl/field_group.hpp b/include/ndtbl/field_group.hpp index 3bbd188..76ea6de 100644 --- a/include/ndtbl/field_group.hpp +++ b/include/ndtbl/field_group.hpp @@ -12,6 +12,7 @@ #include #include #include +#include #include #include @@ -153,8 +154,27 @@ class FieldGroup template std::vector evaluate_all(const Stencil& stencil) const { - std::vector results(field_count(), Stored(0)); - evaluate_all_into(stencil, results.data()); + return evaluate_all_as(stencil); + } + + /** + * @brief Evaluate all fields using a previously prepared interpolation + * stencil and caller-selected output precision. + * + * @tparam Output Floating-point type produced by interpolation. + * @tparam Stencil Fixed-size interpolation stencil type. + * @param stencil Prepared stencil to reuse across fields. + * @return Interpolated field values in storage order. + * @see evaluate_all(const Stencil&) + */ + template + std::vector evaluate_all_as(const Stencil& stencil) const + { + static_assert(std::is_floating_point::value, + "FieldGroup output type must be floating point"); + + std::vector results(field_count(), Output(0)); + evaluate_all_into_as(stencil, results.data()); return results; } @@ -170,10 +190,29 @@ class FieldGroup template void evaluate_all_into(const Stencil& stencil, Stored* results) const { + evaluate_all_into_as(stencil, results); + } + + /** + * @brief Evaluate all fields using a previously prepared interpolation + * stencil into caller-provided storage with caller-selected output precision. + * + * @tparam Output Floating-point type produced by interpolation. + * @tparam Stencil Fixed-size interpolation stencil type. + * @param stencil Prepared stencil to reuse across fields. + * @param results Output buffer with space for `field_count()` values. + * @see evaluate_all_into(const Stencil&, Stored*) + */ + template + void evaluate_all_into_as(const Stencil& stencil, Output* results) const + { + static_assert(std::is_floating_point::value, + "FieldGroup output type must be floating point"); + const std::size_t fields = field_count(); const Stored* const typed_values = interleaved_values_.typed_data(); - std::fill(results, results + fields, Stored(0)); + std::fill(results, results + fields, Output(0)); if (typed_values != nullptr) { evaluate_all_into_typed(stencil, typed_values, fields, results); return; @@ -195,7 +234,25 @@ class FieldGroup const std::array& coordinates, bounds_policy policy = bounds_policy::clamp) const { - return evaluate_all(grid_.prepare_linear(coordinates, policy)); + return evaluate_all_linear_as(coordinates, policy); + } + + /** + * @brief Evaluate all fields directly from query coordinates using + * multilinear interpolation and caller-selected output precision. + * + * @tparam Output Floating-point type produced by interpolation. + * @param coordinates Query coordinates in grid axis order. + * @param policy Bounds handling behavior for out-of-domain coordinates. + * @return Interpolated field values in storage order. + * @see evaluate_all_linear(const std::array&, bounds_policy) + */ + template + std::vector evaluate_all_linear_as( + const std::array& coordinates, + bounds_policy policy = bounds_policy::clamp) const + { + return evaluate_all_as(grid_.prepare_linear(coordinates, policy)); } /** @@ -212,7 +269,28 @@ class FieldGroup Stored* results, bounds_policy policy = bounds_policy::clamp) const { - evaluate_all_into(grid_.prepare_linear(coordinates, policy), results); + evaluate_all_linear_into_as(coordinates, results, policy); + } + + /** + * @brief Evaluate all fields directly from query coordinates using + * multilinear interpolation into caller-provided storage with + * caller-selected output precision. + * + * @tparam Output Floating-point type produced by interpolation. + * @param coordinates Query coordinates in grid axis order. + * @param results Output buffer with space for `field_count()` values. + * @param policy Bounds handling behavior for out-of-domain coordinates. + * @see evaluate_all_linear_into(const std::array&, Stored*, + * bounds_policy) + */ + template + void evaluate_all_linear_into_as( + const std::array& coordinates, + Output* results, + bounds_policy policy = bounds_policy::clamp) const + { + evaluate_all_into_as(grid_.prepare_linear(coordinates, policy), results); } /** @@ -232,7 +310,25 @@ class FieldGroup const std::array& coordinates, bounds_policy policy = bounds_policy::clamp) const { - return evaluate_all(grid_.prepare_cubic(coordinates, policy)); + return evaluate_all_cubic_as(coordinates, policy); + } + + /** + * @brief Evaluate all fields directly from query coordinates using local + * cubic interpolation and caller-selected output precision. + * + * @tparam Output Floating-point type produced by interpolation. + * @param coordinates Query coordinates in grid axis order. + * @param policy Bounds handling behavior for out-of-domain coordinates. + * @return Cubically interpolated field values in storage order. + * @see evaluate_all_cubic(const std::array&, bounds_policy) + */ + template + std::vector evaluate_all_cubic_as( + const std::array& coordinates, + bounds_policy policy = bounds_policy::clamp) const + { + return evaluate_all_as(grid_.prepare_cubic(coordinates, policy)); } /** @@ -249,37 +345,58 @@ class FieldGroup Stored* results, bounds_policy policy = bounds_policy::clamp) const { - evaluate_all_into(grid_.prepare_cubic(coordinates, policy), results); + evaluate_all_cubic_into_as(coordinates, results, policy); + } + + /** + * @brief Evaluate all fields directly from query coordinates using local + * cubic interpolation into caller-provided storage with caller-selected + * output precision. + * + * @tparam Output Floating-point type produced by interpolation. + * @param coordinates Query coordinates in grid axis order. + * @param results Output buffer with space for `field_count()` values. + * @param policy Bounds handling behavior for out-of-domain coordinates. + * @see evaluate_all_cubic_into(const std::array&, Stored*, + * bounds_policy) + */ + template + void evaluate_all_cubic_into_as( + const std::array& coordinates, + Output* results, + bounds_policy policy = bounds_policy::clamp) const + { + evaluate_all_into_as(grid_.prepare_cubic(coordinates, policy), results); } private: - template + template void evaluate_all_into_typed(const Stencil& stencil, const Stored* values, std::size_t fields, - Stored* results) const + Output* results) const { for (std::size_t point = 0; point < Stencil::points; ++point) { - const double weight = stencil.weight(point); + const Output weight = static_cast(stencil.weight(point)); const std::size_t base = stencil.point_index(point) * fields; for (std::size_t field = 0; field < fields; ++field) { const Stored value = values[base + field]; - results[field] += static_cast(weight * value); + results[field] += weight * static_cast(value); } } } - template + template void evaluate_all_into_bytes(const Stencil& stencil, std::size_t fields, - Stored* results) const + Output* results) const { for (std::size_t point = 0; point < Stencil::points; ++point) { - const double weight = stencil.weight(point); + const Output weight = static_cast(stencil.weight(point)); const std::size_t base = stencil.point_index(point) * fields; for (std::size_t field = 0; field < fields; ++field) { const Stored value = interleaved_values_.unchecked(base + field); - results[field] += static_cast(weight * value); + results[field] += weight * static_cast(value); } } } diff --git a/include/ndtbl/io.hpp b/include/ndtbl/io.hpp index 3dba326..3cd2a88 100644 --- a/include/ndtbl/io.hpp +++ b/include/ndtbl/io.hpp @@ -129,7 +129,8 @@ write_group(const std::string& path, * * @param path Input file path. * @return Parsed metadata describing the stored group. - * @see read_group + * @see read_field_group + * @see read_runtime_field_group */ inline GroupMetadata read_group_metadata(const std::string& path) @@ -142,22 +143,76 @@ read_group_metadata(const std::string& path) return detail::read_group_metadata_impl(is); } +/** + * @brief Read an ndtbl file into a typed field group. + * + * The file dimension and scalar payload type must match the compile-time + * `Dim` and `Stored` arguments. + * + * @tparam Dim Expected grid dimensionality of the file. + * @tparam Stored Expected scalar payload type stored in the file. + * @param path Input file path. + * @return Typed field group with payload storage of type `Stored`. + * @see read_group_metadata + * @see FieldGroup + */ +template +inline FieldGroup +read_field_group(const std::string& path) +{ + std::ifstream is(path.c_str(), std::ios::binary); + if (!is.is_open()) { + throw std::runtime_error("failed to open ndtbl input file: " + path); + } + + const detail::parsed_group_layout layout = detail::read_group_layout_impl(is); + const GroupMetadata& metadata = layout.metadata; + if (metadata.dimension != Dim) { + throw std::runtime_error( + "ndtbl file dimension does not match typed loader"); + } + if (metadata.value_type != scalar_type_of()) { + throw std::runtime_error( + "ndtbl file scalar type does not match typed loader"); + } + + const std::array axes = detail::fixed_axes(metadata.axes); + const Grid grid(axes); +#if NDTBL_ENABLE_MMAP + const std::shared_ptr payload_owner = + detail::map_payload_bytes(path, layout.payload_offset, layout.payload_size); + return FieldGroup( + grid, + metadata.field_names, + PayloadView(payload_owner.get(), layout.value_count), + payload_owner); +#else + // Keep this non-const so the payload buffer can be moved into the read-only + // FieldGroup storage instead of copied during load. + std::vector values = + detail::read_payload(is, layout.value_count); + return FieldGroup(grid, metadata.field_names, std::move(values)); +#endif +} + /** * @brief Read an ndtbl file into a runtime-erased field group. * - * The file dimension must match the compile-time `Dim` argument. + * The file dimension must match the compile-time `Dim` argument. The scalar + * payload type is selected from file metadata at runtime. * * @tparam Dim Expected grid dimensionality of the file. * @tparam Output Output type used by runtime-erased interpolation. * @param path Input file path. * @return Runtime-erased field group with either float or double payload * storage and `Output` interpolation results. + * @see read_field_group * @see read_group_metadata * @see RuntimeFieldGroup */ template inline RuntimeFieldGroup -read_group(const std::string& path) +read_runtime_field_group(const std::string& path) { std::ifstream is(path.c_str(), std::ios::binary); if (!is.is_open()) { diff --git a/include/ndtbl/runtime_field_group.hpp b/include/ndtbl/runtime_field_group.hpp index fb2785d..02e167a 100644 --- a/include/ndtbl/runtime_field_group.hpp +++ b/include/ndtbl/runtime_field_group.hpp @@ -266,22 +266,16 @@ class RuntimeFieldGroup Output* values, bounds_policy policy) const override { - evaluate_all_linear_into_impl( - coordinates, - values, - policy, - typename std::is_same::type()); + group_.template evaluate_all_linear_into_as( + coordinates, values, policy); } void evaluate_all_cubic_into(const std::array& coordinates, Output* values, bounds_policy policy) const override { - evaluate_all_cubic_into_impl( - coordinates, - values, - policy, - typename std::is_same::type()); + group_.template evaluate_all_cubic_into_as( + coordinates, values, policy); } void write(std::ostream& os) const override @@ -290,72 +284,6 @@ class RuntimeFieldGroup } FieldGroup group_; - - private: - void evaluate_all_linear_into_impl( - const std::array& coordinates, - Output* values, - bounds_policy policy, - std::true_type) const - { - group_.evaluate_all_linear_into(coordinates, values, policy); - } - - void evaluate_all_linear_into_impl( - const std::array& coordinates, - Output* values, - bounds_policy policy, - std::false_type) const - { - std::vector& scratch = thread_scratch(); - const std::size_t count = group_.field_count(); - if (scratch.size() < count) { - scratch.resize(count); - } - - group_.evaluate_all_linear_into(coordinates, scratch.data(), policy); - copy_scratch_to_output(scratch, count, values); - } - - void evaluate_all_cubic_into_impl( - const std::array& coordinates, - Output* values, - bounds_policy policy, - std::true_type) const - { - group_.evaluate_all_cubic_into(coordinates, values, policy); - } - - void evaluate_all_cubic_into_impl( - const std::array& coordinates, - Output* values, - bounds_policy policy, - std::false_type) const - { - std::vector& scratch = thread_scratch(); - const std::size_t count = group_.field_count(); - if (scratch.size() < count) { - scratch.resize(count); - } - - group_.evaluate_all_cubic_into(coordinates, scratch.data(), policy); - copy_scratch_to_output(scratch, count, values); - } - - static std::vector& thread_scratch() - { - static thread_local std::vector scratch; - return scratch; - } - - static void copy_scratch_to_output(const std::vector& scratch, - std::size_t count, - Output* values) - { - for (std::size_t field = 0; field < count; ++field) { - values[field] = static_cast(scratch[field]); - } - } }; std::shared_ptr impl_; diff --git a/tests/ndtbl_interpolation_t.cpp b/tests/ndtbl_interpolation_t.cpp index 97420a3..ad3002f 100644 --- a/tests/ndtbl_interpolation_t.cpp +++ b/tests/ndtbl_interpolation_t.cpp @@ -117,7 +117,8 @@ require_loaded_linear_recovery(const std::array& axes, const std::string path = ndtbl_test::temporary_path(); ndtbl::write_group(path, group); - const ndtbl::RuntimeFieldGroup loaded = ndtbl::read_group(path); + const ndtbl::RuntimeFieldGroup loaded = + ndtbl::read_runtime_field_group(path); std::array values = { 0.0, 0.0 }; loaded.evaluate_all_linear_into(query, values.data()); @@ -241,6 +242,24 @@ TEST_CASE("field group evaluates unaligned byte-backed payloads", REQUIRE(values[0] == Catch::Approx(1.5)); } +TEST_CASE("field group evaluates float payloads into double output", + "[field_group][interpolation]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + }; + const ndtbl::FieldGroup<1, float> group( + ndtbl::Grid<1>(axes), { "A" }, { 16777216.0f, 16777218.0f }); + + const std::vector values = + group.evaluate_all_linear_as({ 0.5 }); + std::array into_values = { 0.0 }; + group.evaluate_all_linear_into_as({ 0.5 }, into_values.data()); + + REQUIRE(values[0] == Catch::Approx(16777217.0)); + REQUIRE(into_values[0] == Catch::Approx(16777217.0)); +} + TEST_CASE("field group exactly recovers linear fields on uniform axes in 2D", "[field_group][interpolation]") { @@ -425,7 +444,8 @@ TEST_CASE("runtime field group forwards throw bounds policy", const std::string path = ndtbl_test::temporary_path(); ndtbl::write_group(path, group); - const ndtbl::RuntimeFieldGroup<2> loaded = ndtbl::read_group<2>(path); + const ndtbl::RuntimeFieldGroup<2> loaded = + ndtbl::read_runtime_field_group<2>(path); std::array values = { 0.0, 0.0 }; REQUIRE_NOTHROW( loaded.evaluate_all_linear_into({ -2.0, 7.0 }, values.data())); @@ -484,6 +504,26 @@ TEST_CASE("field group cubic interpolation exactly recovers 1D cubic data", REQUIRE(into_values[0] == Catch::Approx(cubic_1d(0.4))); } +TEST_CASE("field group cubic interpolation supports double output from floats", + "[field_group][interpolation][cubic]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 3.0, 4), + }; + const ndtbl::FieldGroup<1, float> group( + ndtbl::Grid<1>(axes), + { "linear" }, + { 16777216.0f, 16777218.0f, 16777220.0f, 16777222.0f }); + + const std::vector values = + group.evaluate_all_cubic_as({ 0.5 }); + std::array into_values = { 0.0 }; + group.evaluate_all_cubic_into_as({ 0.5 }, into_values.data()); + + REQUIRE(values[0] == Catch::Approx(16777217.0)); + REQUIRE(into_values[0] == Catch::Approx(16777217.0)); +} + TEST_CASE("field group cubic interpolation exactly recovers explicit-axis data", "[field_group][interpolation][cubic]") { @@ -577,7 +617,8 @@ TEST_CASE("runtime field groups expose explicit cubic interpolation", const std::string path = ndtbl_test::temporary_path(); ndtbl::write_group(path, group); - const ndtbl::RuntimeFieldGroup<1> loaded = ndtbl::read_group<1>(path); + const ndtbl::RuntimeFieldGroup<1> loaded = + ndtbl::read_runtime_field_group<1>(path); const std::vector values = loaded.evaluate_all_cubic({ 0.4 }); REQUIRE(values[0] == Catch::Approx(cubic_1d(0.4))); diff --git a/tests/ndtbl_io_t.cpp b/tests/ndtbl_io_t.cpp index cb6d0a0..b623346 100644 --- a/tests/ndtbl_io_t.cpp +++ b/tests/ndtbl_io_t.cpp @@ -35,7 +35,8 @@ TEST_CASE("typed loader round-trips metadata and float payloads", "[io]") REQUIRE(metadata.value_type == ndtbl::scalar_type::float32); REQUIRE(metadata.field_names == std::vector({ "A", "B" })); - const ndtbl::RuntimeFieldGroup<2> loaded = ndtbl::read_group<2>(path); + const ndtbl::RuntimeFieldGroup<2> loaded = + ndtbl::read_runtime_field_group<2>(path); std::array values = { 0.0, 0.0 }; loaded.evaluate_all_linear_into({ 0.5, 0.5 }, values.data()); REQUIRE(values[0] == Catch::Approx(1.5)); @@ -44,6 +45,66 @@ TEST_CASE("typed loader round-trips metadata and float payloads", "[io]") std::remove(path.c_str()); } +TEST_CASE("typed field group loader reads known double payloads", "[io]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + }; + const ndtbl::FieldGroup<1, double> group( + ndtbl::Grid<1>(axes), { "A", "B" }, { 2.0, 10.0, 4.0, 14.0 }); + + const std::string path = ndtbl_test::temporary_path(); + ndtbl::write_group(path, group); + + const ndtbl::FieldGroup<1, double> loaded = + ndtbl::read_field_group<1, double>(path); + std::array values = { 0.0, 0.0 }; + loaded.evaluate_all_linear_into({ 0.25 }, values.data()); + + REQUIRE(values[0] == Catch::Approx(2.5)); + REQUIRE(values[1] == Catch::Approx(11.0)); + + std::remove(path.c_str()); +} + +TEST_CASE("typed field group loader reads known float payloads", "[io]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + }; + const ndtbl::FieldGroup<1, float> group( + ndtbl::Grid<1>(axes), { "A" }, { 1.0f, 3.0f }); + + const std::string path = ndtbl_test::temporary_path(); + ndtbl::write_group(path, group); + + const ndtbl::FieldGroup<1, float> loaded = + ndtbl::read_field_group<1, float>(path); + std::array values = { 0.0f }; + loaded.evaluate_all_linear_into({ 0.5 }, values.data()); + + REQUIRE(values[0] == Catch::Approx(2.0f)); + + std::remove(path.c_str()); +} + +TEST_CASE("typed field group loader rejects wrong scalar type", "[io]") +{ + const std::array axes = { + ndtbl::Axis::uniform(0.0, 1.0, 2), + }; + const ndtbl::FieldGroup<1, double> group( + ndtbl::Grid<1>(axes), { "A" }, { 0.0, 1.0 }); + + const std::string path = ndtbl_test::temporary_path(); + ndtbl::write_group(path, group); + + REQUIRE_THROWS_AS((ndtbl::read_field_group<1, float>(path)), + std::runtime_error); + + std::remove(path.c_str()); +} + TEST_CASE("runtime field group can be rewritten after reading", "[io]") { const std::array axes = { @@ -59,11 +120,12 @@ TEST_CASE("runtime field group can be rewritten after reading", "[io]") const std::string output_path = ndtbl_test::temporary_path(); ndtbl::write_group(input_path, group); - const ndtbl::RuntimeFieldGroup<2> loaded = ndtbl::read_group<2>(input_path); + const ndtbl::RuntimeFieldGroup<2> loaded = + ndtbl::read_runtime_field_group<2>(input_path); ndtbl::write_group(output_path, loaded); const ndtbl::RuntimeFieldGroup<2> rewritten = - ndtbl::read_group<2>(output_path); + ndtbl::read_runtime_field_group<2>(output_path); std::array values = { 0.0, 0.0 }; rewritten.evaluate_all_linear_into({ 0.5, 0.5 }, values.data()); REQUIRE(values[0] == Catch::Approx(1.5)); @@ -117,7 +179,7 @@ TEST_CASE("typed loader supports caller-selected runtime output precision", ndtbl::write_group(path, group); const ndtbl::RuntimeFieldGroup<1, float> loaded = - ndtbl::read_group<1, float>(path); + ndtbl::read_runtime_field_group<1, float>(path); std::array values = { 0.0f }; loaded.evaluate_all_linear_into({ 0.5 }, values.data()); REQUIRE(values[0] == Catch::Approx(3.0f)); @@ -136,7 +198,8 @@ TEST_CASE("raw writer round-trips explicit metadata", "[io]") const std::string path = ndtbl_test::temporary_path(); ndtbl::write_group(path, metadata, payload); - const ndtbl::RuntimeFieldGroup<1> loaded = ndtbl::read_group<1>(path); + const ndtbl::RuntimeFieldGroup<1> loaded = + ndtbl::read_runtime_field_group<1>(path); std::array values = { 0.0, 0.0 }; loaded.evaluate_all_linear_into({ 0.5 }, values.data()); REQUIRE(values[0] == Catch::Approx(0.5)); @@ -231,7 +294,10 @@ TEST_CASE("typed loader rejects mismatched dimensions", "[io]") const std::string path = ndtbl_test::temporary_path(); ndtbl::write_group(path, group); - REQUIRE_THROWS_AS(ndtbl::read_group<2>(path), std::runtime_error); + REQUIRE_THROWS_AS(ndtbl::read_runtime_field_group<2>(path), + std::runtime_error); + REQUIRE_THROWS_AS((ndtbl::read_field_group<2, double>(path)), + std::runtime_error); std::remove(path.c_str()); } @@ -252,7 +318,8 @@ TEST_CASE("typed loader rejects truncated payload files", "[io]") bytes.pop_back(); ndtbl_test::write_file_bytes(path, bytes); - REQUIRE_THROWS_AS(ndtbl::read_group<1>(path), std::runtime_error); + REQUIRE_THROWS_AS(ndtbl::read_runtime_field_group<1>(path), + std::runtime_error); std::remove(path.c_str()); }