Partially Addresses #100 - #111
Open
ryanmrichard wants to merge 19 commits into
Open
Conversation
Adds util/fp_traits.hpp as the single customization point for the floating-point operations used by the quadrature generators. The primary template dispatches via the two-step using-declaration idiom, so types whose math is reachable by ADL need no specialization. Separates double -> T conversion into from_exact (for values with an exact binary representation, e.g. 1.0, 0.5, 2.0) and from_inexact (for tabulated data and irrational constants, e.g. M_PI, ln 2, 16.0/3.0). Types that bound their own error must widen the latter but must not widen the former: widening an exact literal injects spurious width into expressions such as 1.0 - x, which sit on the cancellation-sensitive path of several radial transformations. Removes the std::is_floating_point_v gate on the Mura-Knowles quadrature_traits specialization. The restriction was spurious -- the arithmetic is closed-form -- and because it was expressed as SFINAE, a non-floating-point type fell through to the undefined primary template and produced an incomplete-type diagnostic rather than an explanatory one. Conversely, adds static_assert guards to GaussLegendre and GaussLobatto, whose nodes come from a scalar Newton iteration testing against numeric_limits<double>::epsilon(). That test cannot succeed for interval types, and naive interval Newton widens rather than contracts. A static_assert is used rather than SFINAE so the failure explains itself. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The tabulated solid-angle grids were declared
template <typename T> struct lebedev_laikov_N {
static constexpr std::array<cartesian_pt_t<T>,N> points = {...};
static constexpr std::array<T,N> weights = {...};
};
constexpr requires a literal type, so instantiating any of these over a type
that is not literal -- sigma::Interval (wraps boost::numeric::interval),
sigma::Uncertain (holds an unordered_map of dependencies) -- fails to compile.
This was the mechanical barrier to using the library with such types.
Retypes the storage to double across all 231 table files and converts in
copy_grid, which is their only consumer. The surrounding template parameter is
left in place though now unused by the table itself, so that none of the
dispatch branches in the four family headers has to change.
Conversion goes through fp_traits::from_inexact: the tabulated values are not
exact, carrying both the representation error of their decimal literals and
the residual error of the nonlinear solve that produced the abscissae.
scripts/retype_s2_tables.py performs the rewrite and is idempotent, so the
transformation can be replayed after syncing tables from upstream rather than
reapplied by hand across 231 files.
Also fixes three mixed-type expressions in UniformTrapezoid -- dividing by a
size_t, scaling by a size_t, and *= a double literal -- which compile silently
for the built-in types but have no viable overload for interval types. This
file has no transcendental calls, so an audit that looks only for std:: math
does not find it.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Adds include/integratorxx/config.hpp defining `ixx_int` and `ixx_real`, the `IXX_INT` / `IXX_REAL` macros that denote which kind a literal is, and the `ENABLE_STRING_REALS` build flag. `types.hpp` includes it, so the macros reach every table file through the existing `*_grids.hpp` umbrellas with no new include lines. Replaces fp_traits' from_exact/from_inexact with from_integer, from_real, and divide_integer. The previous split was between values exactly representable as a double and values that are not; measured across all 231 tables, of 1,919,804 literals 1,913,026 are non-integral and *none* of those is exactly representable. The only exactly-representable values in the corpus are the 6,778 occurrences of 0 and +/-1 at the axis points. So the distinction that carries information is integer vs. non-integer, which is what the new names capture. Reclassifies every literal accordingly -- integral by value, not by spelling, so a `2.0` in the source becomes `from_integer(2)` -- and restructures expressions to stay in exact integer arithmetic until a single conversion at the end. Where a constant is the ratio of two integers it now goes through divide_integer, so the division is the only rounding rather than being inherited from a constant the compiler already rounded. Several expressions collapse as a result, e.g. the Gauss-Legendre initial guess `pi*(i - 1/4)/(npts + 1/2)` becomes `pi * (4i-1)/(2(2*npts+1))`, and Gauss-Chebyshev-2-modified's `1 - 2i/(n+1)` becomes the single ratio `(n+1-2i)/(n+1)`. This perturbs double results by about an ulp in places -- it is strictly more accurate -- so the tests are now gated on their tolerances rather than on bitwise equality with the previous commit. `M_PI` cannot be used with IXX_REAL, since stringifying it yields the text "M_PI" rather than digits, so config.hpp spells pi out and the uses refer to that. `M_2_PI` is formed as `from_integer(2) / from_real(ixx_pi)` rather than being spelled as its own constant. The 231 s2 tables switch to `ixx_real` elements with every entry wrapped in IXX_REAL. Entries whose value is integral are wrapped too, despite the rule, because a std::array is homogeneous and they must share their neighbours' element type; they convert exactly either way. The 125 Womersley tables did not tabulate weights at all -- they built them with `create_array<N, T>(4.0 * M_PI / N.0)`, which is templated on T (so it breaks for non-literal types), pre-divides in double, and materializes N identical values. They now record only the point count, and copy_grid forms the weight as pi times the exact rational 4/N: one rounding instead of two, correct under ENABLE_STRING_REALS, and a true enclosure for types that bound their own error. Verified in both modes: the 7 buildable upstream tests pass, and a MuraKnowles x LebedevLaikov grid over sigma::Uncertain<double> still reproduces the double grid bitwise. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
The tabulated solid-angle grids are declared
template <typename T> struct lebedev_laikov_N {
static constexpr std::array<cartesian_pt_t<T>,N> points = {...};
static constexpr std::array<T,N> weights = {...};
};
`constexpr` requires a literal type, so instantiating any of these over a type
that is not literal -- an interval type wrapping boost::numeric::interval, or an
uncertainty-propagating scalar holding a map of dependencies -- fails to
compile. That is the mechanical barrier to using the library with such types.
This commit adds the machinery to store the tables in a fixed type and convert
them to the quadrature's value type on read; the tables themselves are rewritten
in a following commit.
config.hpp defines `ixx_int` and `ixx_real`, the `IXX_INT` / `IXX_REAL` macros
that denote which kind a literal is, and the `INTEGRATORXX_ENABLE_STRING_REALS`
build option. types.hpp includes it, so the macros reach every table file
through the existing `*_grids.hpp` umbrellas with no new include lines.
The split between the two kinds is integer vs. non-integer rather than exactly-
representable vs. not. Measured across all 231 tables, of 1,919,804 literals
1,913,026 are non-integral and *none* of those is exactly representable as a
double; the only exactly-representable values in the corpus are the 6,778
occurrences of 0 and +/-1 at the axis points. So integrality is the distinction
that carries information, which is what `from_integer` / `from_real` capture.
`divide_integer` covers the case of a constant that is the ratio of two
integers, so the division is the only rounding rather than being inherited from
a constant the compiler already rounded.
fp_traits is the single seam through which every operation whose meaning depends
on the value type is routed, so that supporting a new type is a matter of
specializing one class rather than auditing every quadrature. It carries only
the conversions here; the math functions arrive with the quadratures that call
them.
Under ENABLE_STRING_REALS `ixx_real` becomes a `std::string_view` and IXX_REAL
captures its argument's source text rather than a double. A double literal has
already lost precision by the time any code can inspect it, so a type more
precise than double -- or one that must bound its own error -- cannot recover
the intended value from it; carrying the decimal text lets such a type parse it
directly. The default `from_real` still parses via double, so the option is
inert until a type specializes it. It is applied to the `integratorxx` target
rather than to a directory because it changes the type of a table entry: in the
non-header-only build the library's own TUs in src/ and its consumers must agree
on what `ixx_real` is.
copy_grid, the tables' only consumer, does the conversion. It also grows a
branch for equal-weight grids, which record a point count in place of a weight
table; the Womersley tables take that form in the following commit.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Performs the rewrite of the 231 tabulated solid-angle grids to the ixx_real literal policy: element type `T` -> `ixx_real`, every table entry wrapped in `IXX_REAL(...)`, and the Womersley grids' computed weights replaced by the point count they are derived from. A script rather than a hand edit because the change touches ~1.9 million literals, and because it is idempotent: the transformation can be replayed after syncing tables from upstream instead of being reapplied by hand across 231 files. `--check` reports what would change without writing, so CI or a reviewer can confirm the committed tables are what the script produces. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Generated by
python3 scripts/retype_s2_tables.py
and reproducible from this commit's parent by re-running it; the script is
idempotent, so `--check` reports no changes against the result.
The 231 tables switch to `ixx_real` elements with every entry wrapped in
IXX_REAL, which is what lifts the literal-type requirement that `constexpr`
imposes on them. The surrounding `template <typename T>` is left in place though
T is no longer used by the table itself, so that none of the ~30 dispatch
branches in each of the four family headers has to change.
Entries whose value is integral (0 and +/-1 at the axis points) are wrapped in
IXX_REAL too, despite the rule that integral literals are spelled IXX_INT: a
std::array is homogeneous, so they must share their neighbours' element type.
They convert exactly either way.
The 125 Womersley tables did not tabulate weights at all -- they built them with
`create_array<N, T>(4.0 * M_PI / N.0)`, which is templated on T (so it breaks
for non-literal types), pre-divides in double, and materializes N identical
values. They now record only the point count, and copy_grid forms the weight as
pi times the exact rational 4/N: one rounding instead of two, correct under
ENABLE_STRING_REALS, and a true enclosure for types that bound their own error.
That reassociation moves the weight by one ulp for 29 of the 125 sizes -- it is
strictly the more accurate of the two -- and by nothing at all for the other 96.
The Lebedev-Laikov, Delley, and Ahrens-Beylkin grids are unchanged bitwise.
Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Use traits/macros to facilitate generic floating-point types.
Adds Molecular Grids
makes molecular grids type agnostic
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Accomplishes:
N.B. Conflicts with the approach in #110.