From dbbf796c1c15c19691e50d6f310324297be37cb8 Mon Sep 17 00:00:00 2001 From: Susi Lehtola Date: Sat, 29 Aug 2026 18:40:10 +0300 Subject: [PATCH 1/2] Fix container type and empty-vector defects in the composite quadratures Three unrelated defects in the same two headers. SubQuadrature's generator declared its weights in the *point* container type and returned it as the weight_container half of the tuple: point_container points ( ... ); point_container weights( ... ); // wrong type This compiles today only because every radial quadrature happens to use std::vector for both. It breaks as soon as the two differ, which is exactly what mixed precision needs: SubQuadrature> does not compile before this change. Two of the surrounding assertions also tested `>= 0` on size_t, which is vacuous; drop them and keep the bounds that say something. RadialGridPartition::finalize() read partition_idx_.back() without checking that anything had been added, which is undefined behaviour on an empty vector. make_pruned_grid() guards against an empty region list one layer up, but finalize() is public and the class is default-constructible, so the guard is not structural. Throw instead. Six adaptor overrides in the two spherical quadrature headers were written as `override {;` -- a stray empty statement after the opening brace. Harmless, but remove them. Co-Authored-By: Claude Opus 5 (1M context) Claude-Session: https://claude.ai/code/session_01FDTFYJMQ76iujDFNHzZyXF --- .../pruned_spherical_quadrature.hpp | 10 +++++++--- .../spherical_quadrature.hpp | 6 +++--- .../composite_quadratures/sub_quadrature.hpp | 18 ++++++++---------- 3 files changed, 18 insertions(+), 16 deletions(-) diff --git a/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp b/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp index 02e729a..bad0a62 100644 --- a/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp +++ b/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp @@ -5,6 +5,7 @@ #include #include +#include #include #include #include @@ -48,6 +49,9 @@ class RadialGridPartition { template void finalize( const RadialQuad& rq ) { + if( partition_idx_.empty() ) + throw std::runtime_error( + "RadialGridPartition: cannot finalize without any angular quadrature"); if( partition_idx_.back() != rq.npts() ) partition_idx_.emplace_back( rq.npts() ); } @@ -183,13 +187,13 @@ class PrunedSphericalQuadrature : const point_container& sph_points_adaptor() const override { return quad_base_type::points(); } - point_container& sph_points_adaptor() override {; + point_container& sph_points_adaptor() override { return quad_base_type::points(); } - const weight_container& sph_weights_adaptor() const override {; + const weight_container& sph_weights_adaptor() const override { return quad_base_type::weights(); } - weight_container& sph_weights_adaptor() override {; + weight_container& sph_weights_adaptor() override { return quad_base_type::weights(); } diff --git a/include/integratorxx/composite_quadratures/spherical_quadrature.hpp b/include/integratorxx/composite_quadratures/spherical_quadrature.hpp index 09ef168..e7c0c0d 100644 --- a/include/integratorxx/composite_quadratures/spherical_quadrature.hpp +++ b/include/integratorxx/composite_quadratures/spherical_quadrature.hpp @@ -133,13 +133,13 @@ class SphericalQuadrature : const point_container& sph_points_adaptor() const override { return quad_base_type::points(); } - point_container& sph_points_adaptor() override {; + point_container& sph_points_adaptor() override { return quad_base_type::points(); } - const weight_container& sph_weights_adaptor() const override {; + const weight_container& sph_weights_adaptor() const override { return quad_base_type::weights(); } - weight_container& sph_weights_adaptor() override {; + weight_container& sph_weights_adaptor() override { return quad_base_type::weights(); } diff --git a/include/integratorxx/composite_quadratures/sub_quadrature.hpp b/include/integratorxx/composite_quadratures/sub_quadrature.hpp index 0f35375..176afe9 100644 --- a/include/integratorxx/composite_quadratures/sub_quadrature.hpp +++ b/include/integratorxx/composite_quadratures/sub_quadrature.hpp @@ -49,16 +49,14 @@ struct quadrature_traits< inline static std::tuple generate( range_idx_type idx_range, const Quadrature& quad ) { const auto [begin_idx, end_idx] = idx_range; - assert( begin_idx >= 0 ); - assert( begin_idx <= quad.npts() ); - assert( end_idx >= 0 ); - assert( end_idx <= quad.npts() ); - assert( begin_idx < end_idx ); - - point_container points( quad.points().begin() + begin_idx, - quad.points().begin() + end_idx ); - point_container weights( quad.weights().begin() + begin_idx, - quad.weights().begin() + end_idx ); + // begin_idx and end_idx are size_t, so the lower bounds are vacuous + assert( begin_idx < end_idx ); + assert( end_idx <= quad.npts() ); + + point_container points ( quad.points().begin() + begin_idx, + quad.points().begin() + end_idx ); + weight_container weights( quad.weights().begin() + begin_idx, + quad.weights().begin() + end_idx ); return std::make_tuple( points, weights ); } From 7ab221b260a4d6928bf1cdc34e6bbbda3653402f Mon Sep 17 00:00:00 2001 From: Susi Lehtola Date: Tue, 6 Oct 2026 09:40:17 +0300 Subject: [PATCH 2/2] Close the remaining empty-partition UB and cover both fixes finalize() now refuses an empty partition, but the same default-constructed object still reached undefined behaviour through end(), cend() and end() const: each computes partition_idx_.end()-1 to skip the npts sentinel, which steps before begin() when there are no entries at all. Under -D_GLIBCXX_DEBUG, RadialGridPartition>{}.end() aborts with "attempt to retreat a past-the-end iterator". Present an empty range in that case instead. Both fixes were untested. The new "Unfinalized Partition" section pins the finalize() throw and the empty range; the new "Sub Quadratures" case covers the container type, which nothing in the suite instantiated. The latter is not a silent regression but a hard one: with master's sub_quadrature.hpp, SubQuadrature> fails to compile, since the weights are built as a vector that make_tuple cannot convert to vector. Also drop a comment that described the removed asserts rather than the code left behind. Co-Authored-By: Claude Opus 5 Claude-Session: https://claude.ai/code/session_01FDTFYJMQ76iujDFNHzZyXF --- .../pruned_spherical_quadrature.hpp | 7 +++ .../composite_quadratures/sub_quadrature.hpp | 1 - test/composite_quadratures.cxx | 45 +++++++++++++++++++ 3 files changed, 52 insertions(+), 1 deletion(-) diff --git a/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp b/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp index bad0a62..87aa160 100644 --- a/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp +++ b/include/integratorxx/composite_quadratures/pruned_spherical_quadrature.hpp @@ -138,6 +138,11 @@ class RadialGridPartition { return iterator( partition_idx_.begin(), quads_.begin() ); } iterator end() { + // partition_idx_ carries one entry past the last partition (the npts + // sentinel added by finalize), so the end index iterator is one before its + // end. An unfinalized partition holds no entries at all, and end()-1 would + // step before begin(); the range is simply empty in that case. + if( partition_idx_.empty() ) return begin(); return iterator( partition_idx_.end()-1, quads_.end() ); } @@ -145,6 +150,7 @@ class RadialGridPartition { return const_iterator( partition_idx_.cbegin(), quads_.cbegin() ); } const_iterator cend() const { + if( partition_idx_.empty() ) return cbegin(); return const_iterator( partition_idx_.cend()-1, quads_.cend() ); } @@ -152,6 +158,7 @@ class RadialGridPartition { return const_iterator( partition_idx_.cbegin(), quads_.cbegin() ); } const_iterator end() const { + if( partition_idx_.empty() ) return begin(); return const_iterator( partition_idx_.cend()-1, quads_.cend() ); } }; diff --git a/include/integratorxx/composite_quadratures/sub_quadrature.hpp b/include/integratorxx/composite_quadratures/sub_quadrature.hpp index 176afe9..9ff10f4 100644 --- a/include/integratorxx/composite_quadratures/sub_quadrature.hpp +++ b/include/integratorxx/composite_quadratures/sub_quadrature.hpp @@ -49,7 +49,6 @@ struct quadrature_traits< inline static std::tuple generate( range_idx_type idx_range, const Quadrature& quad ) { const auto [begin_idx, end_idx] = idx_range; - // begin_idx and end_idx are size_t, so the lower bounds are vacuous assert( begin_idx < end_idx ); assert( end_idx <= quad.npts() ); diff --git a/test/composite_quadratures.cxx b/test/composite_quadratures.cxx index 02fdf79..acba24e 100644 --- a/test/composite_quadratures.cxx +++ b/test/composite_quadratures.cxx @@ -4,6 +4,7 @@ #include #include #include +#include #include #include @@ -239,6 +240,29 @@ TEST_CASE( "Pruned Spherical Quadratures", "[sph-quad]" ) { } + SECTION("Unfinalized Partition") { + + size_t nrad = 4; + IntegratorXX::MuraKnowles r(nrad); + IntegratorXX::LebedevLaikov q(6); + IntegratorXX::RadialGridPartition rgp; + + // Nothing has been added, so there is no last index to compare against + // rq.npts(): finalize must say so rather than read back() off an empty + // vector. + REQUIRE_THROWS( rgp.finalize(r) ); + + // The iterators must present an empty range over it. partition_idx_ is + // empty, so the end()-1 that skips the npts sentinel would step before + // begin(). + CHECK( rgp.begin() == rgp.end() ); + CHECK( rgp.cbegin() == rgp.cend() ); + const auto& crgp = rgp; + CHECK( crgp.begin() == crgp.end() ); + + } + + SECTION("Validity") { size_t nrad = 12; @@ -300,9 +324,30 @@ TEST_CASE( "Pruned Spherical Quadratures", "[sph-quad]" ) { } +TEST_CASE( "Sub Quadratures", "[sph-quad]" ) { + // The slice's weights must be built in the parent's weight_container, not in + // its point_container. The two coincide for the usual case, + // which is why nothing caught it; under mixed precision the weights would be + // narrowed to the point type. + IntegratorXX::Becke r(20); + const std::pair range{ 5ul, 12ul }; + IntegratorXX::SubQuadrature sq( range, r ); + static_assert( std::is_same_v< decltype(sq)::weight_container, + decltype(r)::weight_container >, + "sub-quadrature weights must keep the parent's weight type" ); + static_assert( std::is_same_v< decltype(sq)::point_container, + decltype(r)::point_container >, + "sub-quadrature points must keep the parent's point type" ); + REQUIRE( sq.npts() == range.second - range.first ); + // The slice is a copy, so it must agree with the parent exactly + for( size_t i = 0; i < sq.npts(); ++i ) { + CHECK( sq.points()[i] == r.points()[i + range.first] ); + CHECK( sq.weights()[i] == r.weights()[i + range.first] ); + } +}