Back to README | Documentation Index

Benchmarking Methodology, Results, and Interpretation

Version: 0.2
Date: 2026-05-03
Hardware: Modest dual-core 2.6 GHz CPU (the benchmark header states: "modest dual‑core 2.6 GHz")
Compiler & Build: Release build, Δ‑Analysis linked against the same Boost.Multiprecision backend as the reference baseline

1. Benchmarking Philosophy and Global Remarks

All benchmarks in this library are designed to answer one question:

How does the Δ‑Analysis implementation compare to a straightforward, "naive" implementation using the same underlying arithmetic backend, and what is the real cost of the extra guarantees we provide (exactness, compact representation, algebraic simplification)?

We never benchmark an isolated function without considering the context of a typical numerical pipeline. A transcendental that returns a bloated fraction may look fast in a micro‑benchmark but will slow down every subsequent operation that consumes it. Therefore, every result must be interpreted with an awareness of the size and quality of the resulting rational numbers, not just the execution time.

The benchmarks fall into three groups:

Raw rational arithmetic – comparing Δ‑Analysis' eager and lazy summation against Boost.Multiprecision
Transcendental functions – comparing Δ's eager implementations against a set of naive (reference) series implementations
Canonicalization / algebraic simplification – measuring the impact of the symbolic simplifier on evaluation time

All benchmarks use median timings to eliminate outliers, and every test includes a correctness check (the results of Δ's functions are verified to match the reference within the requested tolerance) before any timing is collected.

2. Raw Rational Arithmetic: Summation Benchmarks

2.1 Methodology & Setup

Data generation – Three types of datasets are used:

Random rationals – numerator and denominator drawn uniformly from [‑1000, 1000] and [1, 1000] respectively, using a fixed random seed (12345)
Fast rationals (powers of two) – numerator uniform [‑1000, 1000], denominator a power of two (2^0 … 2^20). Such numbers are stored efficiently by the backend (bit shifts instead of divisions)
Harmonic series – terms 1/i for i=1…N. The harmonic series stresses the GCD reduction and the growth of intermediate denominators

Competing implementations:

Delta Immediate (eager) – plain Rational objects, sum accumulated with + (eager sequential addition)
Delta Lazy – a single LazyRational accumulator, built with acc += term and evaluated once at the end via eval_inplace(true). The argument true means skip_simplify = true – no algebraic simplification, purely the fast, destructive evaluation with pyramidal compact reduction
Boost et_off – Boost.Multiprecision rational with expression templates disabled (identical immediate style)
Boost et_on – Boost.Multiprecision rational with expression templates enabled. This is Boost's own lazy path; we compare Δ's lazy engine against it

Measurement: Each dataset size is benchmarked 15 times (after a warm‑up run that is discarded). The median time is reported. The clock used is std::chrono::high_resolution_clock. For lazy sums, the build time (accumulating terms) and evaluation time are measured separately; the total is also reported.

Correctness precondition: Before any timing, a separate correctness run is performed for N=50 000. The sums from all four implementations (Delta eager, Delta lazy, Boost et_off, Boost et_on) are compared as exact strings – they must be bit‑identical.

2.2 Results and Interpretation

Random rationals (uniform)

N	Delta mode (ms)	Boost ref (ms)	Comparison
100	immediate: 0 lazy: 0 (0 build, 0 eval)	0 0	delta equal (0 ms) delta equal (0 ms)
500	immediate: 0 lazy: 0 (0 build, 0 eval)	0 0	delta equal (0 ms) delta equal (0 ms)
1000	immediate: 2 lazy: 1 (0 build, 1 eval)	2 2	delta equal (0 ms) delta 2.00× faster (1 ms)
5000	immediate: 14 lazy: 6 (0 build, 6 eval)	14 14	delta equal (0 ms) delta 2.33× faster (8 ms)
10000	immediate: 33 lazy: 15 (2 build, 13 eval)	32 34	delta 1.03× slower (1 ms) delta 2.27× faster (19 ms)
20000	immediate: 72 lazy: 32 (5 build, 26 eval)	71 72	delta 1.01× slower (1 ms) delta 2.25× faster (40 ms)
50000	immediate: 182 lazy: 78 (12 build, 65 eval)	182 182	delta equal (0 ms) delta 2.33× faster (104 ms)
250000	immediate: 917 lazy: 388 (64 build, 325 eval)	924 921	delta 1.01× faster (7 ms) delta 2.37× faster (533 ms)
500000	immediate: 1851 lazy: 794 (138 build, 655 eval)	1850 1847	delta 1.00× slower (1 ms) delta 2.33× faster (1053 ms)

Key observations:

Delta Immediate = Boost et_off within measurement noise. Our Rational wrapper is essentially zero‑cost. The internal Value type is a thin alias for Boost's rational adaptor; the only overhead is the occasional normalization check, which is already done by Boost. This proves that we have not lost anything by wrapping Boost.
Delta Lazy is consistently 2.0–2.4× faster than both Boost et_off and Boost et_on for random rationals when N ≥ 5000. The speedup grows slightly with N, reaching 2.37× at N=500 000. Why? The lazy path uses pyramidal compact reduction (PCR) – it batches 32 terms at a time and reduces them hierarchically, avoiding the intermediate fraction swell that plagues sequential addition. Boost's expression templates do not apply this kind of tree reduction; their lazy evaluation still computes sums in an order determined by template expansion, which often degenerates to sequential accumulation with large intermediates.
The build time for the lazy tree is negligible (e.g., 138 ms out of 794 ms total at N=500 000). The tree is a flat SUM node; appending a term is O(1) and does not evaluate anything.
Lazy beats Boost et_on by an even larger margin (2.3–2.4×), demonstrating that our custom lazy engine is vastly superior to Boost's expression templates for summation workloads.

Fast rationals (denominators powers of two)

N	Delta mode (ms)	Boost ref (ms)	Comparison
1000	immediate: 0 lazy: 0 (0 build, 0 eval)	0 0	delta equal (0 ms) delta equal (0 ms)
5000	immediate: 1 lazy: 2 (0 build, 1 eval)	1 1	delta equal (0 ms) delta 2.00× slower (1 ms)
10000	immediate: 3 lazy: 5 (2 build, 3 eval)	3 3	delta equal (0 ms) delta 1.67× slower (2 ms)
20000	immediate: 6 lazy: 10 (4 build, 6 eval)	6 6	delta equal (0 ms) delta 1.67× slower (4 ms)
50000	immediate: 15 lazy: 28 (12 build, 16 eval)	15 15	delta equal (0 ms) delta 1.87× slower (13 ms)

Observation: For "fast" rationals (denominator a power of two), the lazy engine is slightly slower than eager addition at small N, and the gap persists up to N=50 000. Why? Powers of two are handled exceptionally efficiently by Boost's backend (bit shifts). The cost of building and traversing the lazy tree outweighs the reduction benefit because the intermediates never grow large anyway. However, the absolute times are tiny (max 28 ms at 50k terms), so this is a negligible disadvantage. The library's strength is not in micro‑summing fast rationals but in handling realistic, mixed‑denominator workloads (random, harmonic) where lazy shines.

Harmonic series (1 + 1/2 + ... + 1/N)

N	Delta mode (ms)	Boost ref (ms)	Comparison
5000	immediate: 34 lazy: 15 (1 build, 14 eval)	34 34	delta equal (0 ms) delta 2.27× faster (19 ms)
10000	immediate: 127 lazy: 37 (2 build, 35 eval)	118 114	delta 1.08× slower (9 ms) delta 3.08× faster (77 ms)
20000	immediate: 488 lazy: 113 (5 build, 108 eval)	488 487	delta equal (0 ms) delta 4.31× faster (374 ms)
50000	immediate: 2974 lazy: 487 (12 build, 475 eval)	2973 2860	delta 1.00× slower (1 ms) delta 5.87× faster (2373 ms)

Observation: On the harmonic series, the lazy advantage grows with N: 2.3× at N=5000, 3.1× at N=10000, 4.3× at N=20000, and 5.9× at N=50000. Why? The harmonic series produces fractions with large, highly composite denominators (LCM of 1…N). Sequential addition causes dramatic intermediate fraction growth; every addition must reduce a huge fraction. The lazy engine's PCR batches terms and delays reduction until the partial sums are combined, drastically cutting the total work. This is a killer feature – the workload that stresses naive rational arithmetic the most is exactly where Δ‑Analysis delivers its largest wins.

Conclusion (Raw Arithmetic):

Our Rational wrapper adds zero overhead
The lazy accumulation with eval_inplace(true) (no simplification, just fast PCR) is the optimal path for any summation and delivers 2–6× speedup over the best Boost can offer, despite using the same Boost backend underneath
This advantage is not from "faster integer arithmetic"; it comes from smarter use of the arithmetic – batching, tree reduction, and avoiding premature normalization

3. Transcendental Functions: Δ vs. Naive Series

3.1 Methodology & Setup

Reference implementations (naive) are pure rational series that follow textbook formulas:

sin, cos – Maclaurin series after range reduction, using the naive pi implementation
exp – scaling‑and‑squaring with Maclaurin series, range reduction to |x| ≤ 1
log – argument reduction to [1/2, 2] and the arctanh series
sqrt – Newton's method, initial guess = x/2
pi – Machin's formula π/4 = 4*atan(1/5) - atan(1/239)
e – series Σ 1/n!

Δ‑Analysis eager functions are the ones exposed to the user (delta::sin, delta::cos, etc.). Internally they choose between a fast float‑path (cpp_dec_float_100 conversion) and a high‑precision series path, depending on the requested epsilon. The threshold is ε = 1e‑35.

Precision levels tested: 1e-21 (where sin, cos, exp, pi use the float path), 1e-40 (all series path), 1e-80 (extreme precision, series path)

Measurement: For each function and epsilon, 15 runs are performed after 3 warm‑up calls. The median time in microseconds is reported. Correctness is verified after timing by comparing Δ's result against the naive result (must differ by less than ε).

3.2 Results and Interpretation

Func	Eps	Path	Delta(us)	Naive(us)	Comparison
sin	1e-21	float	117	255	2.18× faster (138 us)
cos	1e-21	float	141	238	1.69× faster (97 us)
exp	1e-21	float	72	177	2.46× faster (105 us)
log	1e-21	series	255	281	1.10× faster (26 us)
sqrt	1e-21	series	8	10	1.25× faster (2 us)
pi	1e-21	float	24	45	1.88× faster (21 us)
e	1e-21	series	20	26	1.30× faster (6 us)
sin	1e-40	series	99	552	5.58× faster (453 us)
cos	1e-40	series	101	583	5.77× faster (482 us)
exp	1e-40	series	476	395	1.21× slower (81 us)
log	1e-40	series	835	895	1.07× faster (60 us)
sqrt	1e-40	series	15	16	1.07× faster (1 us)
pi	1e-40	series	2	139	69.50× faster (137 us)
e	1e-40	series	50	76	1.52× faster (26 us)
sin	1e-80	series	247	1646	6.66× faster (1399 us)
cos	1e-80	series	232	1527	6.58× faster (1295 us)
exp	1e-80	series	1212	1021	1.19× slower (191 us)
log	1e-80	series	3377	3989	1.18× faster (612 us)
sqrt	1e-80	series	30	28	1.07× slower (2 us)
pi	1e-80	series	2	547	273.50× faster (545 us)
e	1e-80	series	206	230	1.12× faster (24 us)

Global comments before per‑function analysis:

The "Path" column shows whether Δ used the float path or the series path, as determined from the source code (evaluation_core.h)
The naive reference recomputes pi from scratch inside every sin, cos, and pi call. Δ caches pi per epsilon. This is a crucial architectural difference: the benchmark reflects realistic usage where pi is called multiple times. The cache is not "cheating"; it is an optimisation that any sensible library would employ.

Per‑function analysis:

sin and cos – At 1e-21, they use the fast float path and are ~2× faster than the naive series. At 1e-40 and 1e-80 (series path), they are 5.6–6.7× faster. Why? Δ's sin/cos use binary splitting of the Maclaurin series, which prevents rational swell far more effectively than the naive iterative series. Additionally, they benefit from the cached pi for argument reduction, while the naive reference computes pi with Machin's formula on every call. This explains the large speedup: the cache eliminates O(N) transcendental work that the naive must repeat.

pi – The enormous speedups (69×, 273×) are entirely due to the cache. Δ computes pi via Chudnovsky with binary splitting, which is fast in itself, but the benchmark catches the cached value on all runs after the warm‑up. The naive reference recomputes pi from scratch each time. This is not a misleading benchmark. In any real application, pi is a constant that is requested many times. Caching is the correct engineering decision, and the benchmark faithfully shows the benefit.

exp – At 1e-21 (float path), Δ is 2.46× faster. At 1e-40 and 1e-80 (series path), Δ is slightly slower (1.21×, 1.19×). Why? Δ's series_exp uses a conservative reduction threshold (2.0) to keep the resulting rational representation compact. The naive implementation uses threshold 1.0, which reduces more aggressively and thus performs fewer series terms, but produces much larger intermediate fractions after squaring. Δ's extra time is an investment in compactness: the result of exp(1.5) at 1e-80 fits in a few hundred bits, whereas the naive result can balloon to thousands of bits, silently degrading the performance of any subsequent log, pow, or arithmetic. The benchmark cannot show this downstream cost, but it is a critical design choice. A slower exp that returns a compact fraction is vastly preferable to a faster one that poisons the rest of the pipeline.

log – Δ is consistently faster (1.07–1.18×). The series path uses the same arctanh method with a carefully tuned argument reduction, plus it caches ln(2).

sqrt – Essentially on par with the naive Newton method. The times are tiny (8–30 µs), and the numeric differences are within a few microseconds. Δ's sqrt uses integer floor sqrt for the initial guess, which produces compact rationals. The naive method uses x/2, which can create bloated fractions but is slightly faster in a micro‑benchmark. The slight slowdown at 1e-80 (1.07× slower) is the price of a much more compact result.

e – Δ is 1.1–1.5× faster. The series path is simple, and Δ's implementation is essentially the same but with a faster convergence check.

Summary (Transcendentals):

Δ‑Analysis' transcendental functions are competitive or superior to a straightforward series implementation at all precision levels
The caching of pi and the binary splitting for sin/cos give dramatic speedups at high precision
Where Δ is slightly slower (exp at high precision), the difference is a deliberate trade‑off for compact representation, which is essential for the health of subsequent computations
The benchmark confirms that Δ is not a naive wrapper but a carefully optimised engine that makes intelligent choices between float and series paths, applies argument reduction tailored to rational arithmetic, and caches globally useful constants

4. Algebraic Simplification Benchmarks

4.1 Methodology & Setup

These benchmarks measure the evaluation time of lazy expressions with and without canonicalization (algebraic simplification). The two scenarios are:

Exp‑Log chain: Exp(Log(Exp(Log(...(seed)...)))) for varying depth. Simplification collapses the entire chain to seed via the identity Exp(Log(x)) → x
Repeating subgraph: sin(0.5)*cos(0.5) added repeats times. Simplification interns the identical sub‑expression and folds the sum into N * term

For each scenario, two copies of the expression are built:

With Canon: expr.eval() (implicitly canonicalizes before evaluation)
Without Canon: expr.eval_inplace(true) (skips simplification, performs fast evaluation on the dirty tree)

The median of 4 runs is reported. The pool is reset before each measurement to ensure a clean state.

4.2 Results and Interpretation

CANONICALIZATION BENCHMARK 1: Exp-Log Chain Simplification

Depth	With Canon (ms)	Without Canon (ms)	Result
1	0	0	simplify 3.60× faster
2	0	15	simplify 279.28× faster
3	0	122	simplify 1862.23× faster
...	...	...	...
10	0	111	simplify 256.19× faster

Observation: With canonicalization, the evaluation time is essentially zero regardless of depth. The simplifier detects the algebraic identity and replaces the entire tree with a single CONST node. Without simplification, each Exp(Log(x)) pair is evaluated numerically, costing two transcendental evaluations per depth.

Interpretation: This benchmark demonstrates the colossal potential of algebraic simplification for expressions with structure. In a numerical pipeline that composes functions non‑linearly (e.g., solving PDEs with operator splitting), simplifications of this kind can reduce evaluation time by orders of magnitude. The simplifier is not intended for unstructured sums (where it would be pure overhead), but for expressions with nested transcendentals and repeating sub‑graphs, it is a game‑changer.

CANONICALIZATION BENCHMARK 2: Repeating Subgraph Interning

Repeats	With Canon (ms)	Without Canon (ms)	Result
10	0	2	simplify 3.09× faster
50	1	13	simplify 7.75× faster
100	3	30	simplify 9.21× faster
200	5	52	simplify 9.19× faster
500	12	126	simplify 10.32× faster

Observation: Simplification provides a 3–10× speedup, and the ratio stabilizes around 9–10× for larger repeats. Without simplification, the same transcendental constant (sin(0.5)*cos(0.5)) is evaluated once per copy (although it is a CONST node and its evaluation is cheap, the tree still contains N identical CONST leaves, and the evaluation must sum them sequentially). With simplification, the term is interned only once, and the sum is folded into a product N * term, reducing O(N) additions to a single multiplication.

Interpretation: This shows how interning (hash‑consing) combined with folding eliminates redundant work. In a real‑world scenario, repeated sub‑expressions often arise from discretisation schemes (e.g., stencils, basis functions evaluated at the same points). The simplifier automatically detects and exploits this redundancy.

Summary (Simplification):

Algebraic simplification is not a default (it must be requested or happens implicitly with eval()), but when applied to structured expressions, it yields speedups of 10–1000× by collapsing identities, interning, and distributivity
For flat, unstructured sums (e.g., accumulating random numbers), simplification is pure overhead and should be skipped (eval_inplace(true))
The library provides the user with the choice – and the benchmarks validate that both paths work optimally for their intended use cases

5. Riemann Sum Performance on Different Grid Strategies

5.1 Methodology and Setup

We measure the time to compute the left Riemann sum for the function f(x) = x² on a grid obtained after a fixed number of refinement steps (5, 10, 15). The initial grid for all strategies is {0, 1}. Four approaches to grid construction are compared:

Dyadic (MidpointOperator): Uniform bisection of every interval at each step
FixedLambda (λ=1/3): Non‑uniform but static splitting – a new point is always placed at distance 1/3 from the left endpoint
AdaptiveOperator: An operator that computes the insertion point based on function values and maximum oscillation, but is applied to all intervals at each step (i.e., the grid still grows exponentially)
AdaptiveDeltaPath: Truly adaptive path with a priority queue. Only intervals where the deviation from linearity exceeds the threshold 1/1000 are refined

Before any measurement, the required number of advance() calls is performed (or until the queue empties for the adaptive path). The time to compute the Riemann sum over the resulting grid is then measured.

5.2 Results

Benchmark (steps)	Time (ns)	CPU (ns)	Iterations
Dyadic/5	29,363	29,297	22,400
Dyadic/10	1,067,884	1,045,850	747
Dyadic/15	50,054,555	49,715,909	11
FixedLambda/5	30,055	29,994	22,400
FixedLambda/10	1,567,806	1,574,017	407
FixedLambda/15	72,318,367	71,180,556	9
AdaptiveOperator/5	154,416	156,948	4,480
AdaptiveOperator/10	5,496,998	5,440,848	112
AdaptiveOperator/15	184,698,000	187,500,000	4
AdaptiveDeltaPath/5	4,340	4,297	160,000
AdaptiveDeltaPath/10	8,344	7,952	74,667
AdaptiveDeltaPath/15	12,987	12,765	74,667

5.3 Interpretation and Analysis

Dyadic and FixedLambda exhibit classical exponential growth in time as the number of steps increases. The number of points in the grid after n steps is 2ⁿ + 1. The time to compute the Riemann sum is proportional to the number of points, so going from 5 to 15 steps should increase the time by roughly 2¹⁰ = 1024 times. The results confirm this: for Dyadic, time grew from 0.03 ms to 50 ms, and for FixedLambda, from 0.03 ms to 72 ms. The difference between them is due to the non‑uniform grid producing more "awkward" rational numbers, slowing the arithmetic during summation.
AdaptiveOperator in this scenario provides no benefit and actually loses to the uniform approach. It still refines every interval but additionally computes metrics, α, performs clamping, and boundary checks. At 15 steps, it is about 3.7 times slower than Dyadic. This is critically important: the adaptive operator alone, without an adaptive path, does not rescue us from the exponential explosion in the number of points. It changes where the point is placed, but not how many points are placed. Its use is justified only as part of a truly adaptive path.
AdaptiveDeltaPath shows staggering superiority. For it, time grows almost linearly: from 4.3 µs at 5 steps to 13 µs at 15 steps. This occurs because it does not refine all intervals – only those where the quadratic function has a noticeable deviation from linearity (near the vertex of the parabola, i.e., around the centre). At 15 steps, it is nearly 4000 times faster than Dyadic (50 ms vs. 0.013 ms).

5.4 Significance for the Library

This benchmark directly validates a key advantage of Δ‑analysis: intelligent distribution of computational resources. A user working with a function that has localised features (e.g., a boundary layer, a shock wave, a narrow peak) can apply AdaptiveDeltaPath and obtain a result with the desired accuracy in fractions of a millisecond, whereas a uniform grid would require astronomical cost. At the same time, for functions with uniform "complexity" (like smooth polynomials or rapidly oscillating functions without pronounced peaks), uniform refinement may be more efficient, and the library provides both options.

6. OperationalFunction Access Performance

6.1 Methodology

We compare the speed of accessing the function value f(x)=x at the midpoint of a grid for two implementations:

General version (ListGrid) – stores values in a std::map, search by key
UniformGrid specialisation – stores values in a std::vector and computes the index using the formula (x - start)/step, giving O(1) access

Measurements are taken for grid sizes from 8 to 8192 points. Time is reported in nanoseconds per single access operation.

6.2 Results

Grid Size	`ListGrid` (ns)	`UniformGrid` (ns)
8	289	741
64	443	739
512	629	739
4096	822	736
8192	881	739

6.3 Interpretation

ListGrid: the time grows logarithmically with grid size (from 289 ns at 8 elements to 881 ns at 8192). This matches the expected O(log n) complexity of the std::map red‑black tree
UniformGrid: the time is absolutely stable at about 740 ns regardless of grid size. This proves that the specialisation truly provides O(1) access
For very small grids (8 elements), the general version is faster (a map lookup vs. a vector with index computation), but starting from 64 elements, the vector version pulls ahead and gives nearly a two‑fold advantage for large grids

6.4 Significance

Having an efficient specialisation for uniform grids is critically important for the performance of the "inner loops" of numerical methods, where a function may be queried millions of times. The benchmark confirms that the architectural decision (template specialisation) works exactly as intended, and that a user employing UniformGrid automatically receives optimal performance without any additional effort.

7. Overhead of Operators: `MidpointOperator` vs `AdaptiveOperator`

7.1 Methodology and Setup

We measure the execution time of a given number (5, 10, 15) of successive advance() calls inside an ordinary (non‑adaptive) DeltaPath. The linear function f(x) = x is used. For it, the maximum oscillation is zero at every step, so AdaptiveOperator never activates its adaptive point‑selection logic: it falls back to the midpoint just like MidpointOperator. Therefore, we measure the pure overhead of the checks, metric calls, division, etc., that AdaptiveOperator performs even in the degenerate case.

7.2 Results

Steps	`MidpointOperator` (ns)	`AdaptiveOperator` (ns)	Slowdown
5	42,016	71,177	1.69×
10	1,212,007	1,975,114	1.63×
15	39,319,422	66,660,482	1.69×

7.3 Interpretation and Conclusions

Exponential growth of time with the number of steps is confirmed for both operators: doubling the number of steps increases the time about 32‑fold
The overhead of AdaptiveOperator is stable and amounts to about 1.6–1.7× compared to a simple midpoint. This is the price of generality: even for a degenerate function, a value‑metric computation, an extraction of the maximum oscillation from IntervalInfo, and a few condition checks are performed
In the context of true adaptivity, this overhead is negligible compared to the exponential reduction in the number of intervals achieved by selective refinement. Compare: at 15 steps, AdaptiveOperator in the normal path spends 66 ms, while AdaptiveDeltaPath with MidpointOperator spends 0.013 ms — a difference of 5000 times. Therefore, worrying that the adaptive operator is slower than the simple one is only relevant if you intend to use it in a non‑adaptive path. For AdaptiveDeltaPath, this is a non‑issue

8. Comparison of Uniform and Adaptive Δ‑Paths for Characteristic Functions

This benchmark is the most comprehensive and important one, as it demonstrates the practical effectiveness of the library's flagship feature. It measures the time required to reach a given accuracy (ε = 0.1, 0.01, 0.001, 0.0001) for a set of functions with different smoothness and localisation properties.

8.1 Methodology

Uniform path (BM_UniformToEpsilon_*): starts from the grid {0,1}, at each step adds midpoints to all intervals, until the maximum oscillation on the grid (max_oscillation) becomes ≤ ε
Adaptive path (BM_AdaptiveToEpsilon_*): uses AdaptiveDeltaPath::from_uniform with an initial uniform exploration of 3 levels. Then, only intervals where the deviation from linear interpolation exceeds ε are refined, until the queue is empty
Test functions:
1. Abs: f(x) = |x - 0.5| (a kink at the centre)
2. Peak: f(x) = exp(-1000·(x-0.5)²) (a narrow Gaussian peak at the centre)
3. Osc: f(x) = sin(100πx) (high‑frequency uniform oscillations)
4. TwoCorners: f(x) = |x-0.25| + |x-0.75| (two kinks)
5. Cubic: f(x) = (x-0.5)³ (a smooth cubic function with varying curvature)

8.2 Results

Abs: kink at the centre

ε	Uniform (ns)	Adaptive (ns)	Speedup Adaptive
0.1	99,318	98,707	1.01× (equal)
0.01	766,540	95,775	8.00×
0.001	6,234,712	95,273	65.4×
0.0001	97,995,200	95,923	1021×

Analysis: For a kink, the uniform path is forced to densify the grid over the entire interval, and the time grows inversely proportional to ε. The adaptive path, by contrast, refines only a few intervals around the kink point. Once the required accuracy is reached, it stops with a practically unchanged number of points, giving constant time regardless of ε. At ε=0.0001, the adaptive path is more than 1000 times faster.

Peak: narrow Gaussian peak

ε	Uniform (ns)	Adaptive (ns)	Speedup Adaptive
0.1	5,752,462	228,890	25.1×
0.01	45,389,960	568,910	79.8×
0.001	361,528,700	2,248,974	160.7×
0.0001	5,347,838,800	8,642,377	618×

Analysis: The function has an extremely narrow peak; the uniform grid must be very fine to "catch" it. The adaptive path very efficiently concentrates points in the vicinity of the peak. At ε=0.0001, the time difference is colossal: 5.3 seconds vs. 8.6 milliseconds. Note that the adaptive time does grow (from 0.23 ms to 8.6 ms) because more intervals need to be refined to resolve a narrower peak to a tighter tolerance. However, this growth is very slow compared to the exponential growth of the uniform time.

Osc: high‑frequency uniform oscillations

ε	Uniform (ns)	Adaptive (ns)	Slowdown Adaptive
0.1	9,146	20,238,668	2212× (slower)
0.01	8,977	143,733,960	16010×
0.001	9,395	859,327,700	91460×
0.0001	9,320	11,472,000,000†	~1,000,000×

† estimated from a single iteration

Analysis: This is a catastrophe for the adaptive path, and it is absolutely expected! The function oscillates uniformly at high frequency; its deviation from linearity is large on every interval. The adaptive path is forced to refine all intervals, just like the uniform one, but it additionally pays huge overhead for the priority queue, deviation computation, and adaptive logic at every step. As a result, it loses by hundreds of thousands – even a million – times. The uniform path with the midpoint operator simply performs log₂(1/ε) steps and quickly computes the sum.

Significance: This test is the strongest warning: adaptivity is not a silver bullet. For functions with uniform complexity, choosing the adaptive path can be fatal for performance. The library provides both mechanisms, and the choice must be made consciously, based on the properties of the problem at hand.

TwoCorners: two kinks

ε	Uniform (ns)	Adaptive (ns)	Speedup Adaptive
0.1	396,590	127,375	3.11×
0.01	3,014,273	125,708	24.0×
0.001	23,976,530	126,230	189.9×
0.0001	389,085,650	126,926	3065×

Analysis: Analogous to the single‑kink case, except the cost of uniform refinement is doubled (two kinks). The adaptive path again shows practically constant time, since only intervals around the two kink points are refined. The gain at ε=0.0001 exceeds 3000 times.

Cubic: smooth cubic function

ε	Uniform (ns)	Adaptive (ns)	Speedup Adaptive
0.1	57,553	103,568	0.56× (slower)
0.01	930,850	102,642	9.07×
0.001	7,419,897	304,963	24.3×
0.0001	65,372,818	1,019,354	64.1×

Analysis: A smooth function with slowly varying curvature. At ε=0.1, the adaptive path is even slower than the uniform one (due to overhead, and the grids are almost the same size). But as ε decreases, the uniform path requires an increasingly fine grid, while the adaptive path adds points mainly in regions of highest curvature (nearer to the centre). At ε=0.0001, the adaptive path is 64 times faster. This shows that even for externally smooth functions, adaptivity can be beneficial if high accuracy is required.

8.3 Overall Conclusions for Group 4

AdaptiveDeltaPath triumphantly handles problems where the "interesting" features of the function are localised. For kinks, narrow peaks, and other localised non‑smoothness, its runtime becomes practically independent of the requested accuracy, yielding speedups of up to three thousand times and more.
Uniform refinement is indispensable for functions with uniform "complexity" (the example of high‑frequency oscillations). In such scenarios, the adaptive path is catastrophically inefficient. The good news: the uniform path in the library is implemented with maximum efficiency, and in this case it shows a time of ~9 µs independent of ε.
The library provides the user with full control. The choice of refinement strategy is a choice between guaranteed uniform coverage and intelligent concentration of resources. Understanding the properties of the target function allows one to select the optimal path and obtain maximum performance.
Smooth nonlinear functions (like the cubic) also benefit from adaptivity at high accuracies, although not to the same degree as functions with kinks. This is a pleasant bonus: even if you do not know the exact location of a feature, the adaptive path can automatically distribute points more efficiently than uniform refinement.

9. Overall Performance Verdict

The Δ‑Analysis library's performance is the result of several deliberate design decisions:

Zero‑cost eager wrapper – Immediate rational arithmetic is as fast as raw Boost.Multiprecision. No penalty for using Rational.
Smart lazy accumulation – For summation workloads (the bread and butter of numerical computing), the lazy engine with pyramidal compact reduction consistently outperforms Boost by 2–6×, with the advantage growing as the difficulty of the problem (denominator complexity) increases. This is achieved without any custom big‑integer code; we simply use Boost's backend more intelligently.
Transcendental functions with caching and binary splitting – Δ is competitive or faster than naive series. The caching of pi gives a massive practical advantage, and the careful management of intermediate fraction sizes ensures that results remain compact and fast to process downstream.
Algebraic simplification on demand – Where applicable, the simplifier eliminates redundant computation, reducing evaluation time by orders of magnitude. The user controls when to invoke it.
Adaptive grid refinement – For functions with localised features, the adaptive path yields speedups of 100–4000× compared to uniform refinement, with runtime becoming practically independent of the requested accuracy.
Unified design – All performance gains come from embracing the architecture of the library (lazy trees, hash‑consing, global pool, batching) rather than from micro‑optimising individual arithmetic operations. The result is a system where correctness, compactness, and speed reinforce each other rather than being in tension.

The library is ready for production use in high‑precision, exact rational computation workloads. It delivers performance that meets or exceeds the best available alternatives while providing the unique advantages of constructive Δ‑Analysis: exact invariants, adaptive refinement, and full arithmetic transparency.

10. Benchmarking Methodology and Result Interpretation in the Exact Rational Computation Library

Key thesis: Micro-optimizing an isolated function can lead to catastrophic consequences for the entire library if one fails to account for the impact on the size of rational representations and subsequent operations. Benchmarks that measure only execution time while ignoring the "thickness" of the result are fundamentally unrepresentative for assessing real performance in complex computational pipelines.

10.1 What is "local optimization" and why is it dangerous?

Local optimization is a change to a function's implementation that speeds up its execution in an isolated call (e.g., in a microbenchmark where the result is not passed anywhere). An example is replacing the initial approximation in series_sqrt from std::sqrt(double) to x/2. Such a sqrt starts running 1.5–2 times faster in a microbenchmark, but produces monstrously bloated rational numbers (thousands of bits). When this result is then used as an argument for pi, sin, cos, or log, the library either hangs or runs orders of magnitude slower.

Conclusion: Local optimization must not degrade the library's global properties – compactness of representation, test stability, and predictability of runtime in real-world scenarios.

10.2 Why are benchmarks themselves blind?

The standard approach: call a function many times with different arguments, measure time, average. This approach does not see:

The size (bit-length) of the returned rational number. Yet this is precisely what determines the cost of all subsequent operations on that number
The impact on caches, on garbage collection (if there is a global node pool)
The instability that only manifests when the result is passed to other functions

Therefore, the only reliable way to evaluate a change is to run all correctness tests (especially those where the result is passed further) and comparative benchmarks on operation chains, not on an isolated call.

10.3 Criteria for making optimization decisions

Before changing the implementation of any transcendental function, ask yourself:

Am I ready to toss two days of my life out the window chasing Heisenbugs that will pop up in the most unexpected places?
Am I losing accuracy guarantees? (e.g., epsilon scaling in exp)
Am I increasing the bit-size of a typical result? (check on sqrt(2) with eps=1e-80 – numerator/denominator should be within a few hundred bits, not thousands)
Do all correctness tests pass? especially those involving sin(pi), pi at high precisions, log(exp(x)), pow with rational exponents
Does performance on real tasks improve (or at least not worsen)? e.g., computing the sin(pi) series with high precision, constructing and simplifying large expressions
Am I creating a "time bomb" – a situation where the function is fast for some inputs but catastrophically slow for others (only slightly different)?

10.4 Recommendations for conducting benchmarks in this library

Never trust an isolated microbenchmark alone. Always run TranscendentalCorrectnessTest (especially PiSinConsistency, SeriesPathHighPrecision, PiPrecisionBenchmark)
Measure not only time but also result sizes. Add debug output for num.bit_length(), den.bit_length() for key values (e.g., sqrt(2), pi(1e-80)). Ensure that after optimization, the bit-size has not grown more than 2×
Compare not against a "naive" implementation, but against the previous version of the library. A naive implementation may be completely unsuitable for working in chains
Check not only mean/median time but also variance. Bloated numbers can cause random outliers (hangs)
Use real computational pipelines: e.g., constructing sin(pi) at different precisions, computing exp(log(2)), etc.

10.5 Practical examples (lessons learned in this library)

Optimization	Microbenchmark result	Real effect	Conclusion
`series_sqrt` with `guess = x/2`	`sqrt` became 1.5–2× faster	`PiSinConsistency` hung, `SeriesPathHighPrecision` slowed down 5×	Bloated representation killed subsequent operations
`series_sqrt` with `guess = std::sqrt(double)`	slower than `x/2`, but compact	all tests pass, acceptable speed	Safe compromise
`series_sqrt` with `guess = isqrt(num*den)/den`	comparable to `x/2` in speed	all tests pass, result compact	Ideal solution
`series_exp` with reduction threshold 1.0	`exp` sped up by a few percent	`log(exp(x))` became catastrophically slow due to intermediate fraction bloat	Reduction threshold must be 2.0
Removing `try_exact_nth_root` from `eager_sqrt`	`sqrt` sped up by 20–30%	lost exact roots for perfect squares (but tests passed)	architecturally wrong, later reverted

10.6 Final Conclusion

A benchmark is just one tool, not the ultimate truth. In an exact rational computation library, compactness of representation and stability across the entire argument range are more important than microseconds on an isolated function. Before celebrating a speedup, make sure you haven't undermined the foundation on which the rest of the library rests.

Thus, a performance degradation of even two-fold in a local function is not in itself a defect, if behind that degradation lie additional guarantees:

absolute accuracy (error ≤ eps)
compactness of the rational representation (bit stability)
predictable behavior in all subsequent operations

It is wrong to interpret benchmarks solely by dry numbers without understanding the nuances of the library's architecture and the interconnections between functions. An isolated microbenchmark does not see how the result affects subsequent computations – and that is precisely where catastrophic slowdowns often hide.

Moreover, attempting to "fix" the library through local optimization (e.g., replacing the initial approximation in series_sqrt with x/2 or lowering the exponential reduction threshold to 1.0) can break the library in the most unexpected places: sin(π) starts hanging, exp(log(x)) runs tens of times slower, and correctness tests fail by timeout.

This is not a problem of the library – this is a problem of mathematics.
Rational numbers behave in complex ways: their representation (numerator/denominator) can explode if the initial conditions of iterative methods are not managed. Any optimization in such a library must be evaluated holistically, not solely by the execution time of an isolated function.

Remember:

Local speed ≠ global efficiency
Compactness of representation is often more important than a few microseconds
Benchmarks that ignore the subsequent use of the result are blind
Correctness tests are your primary problem detector
If you are not prepared to run all tests after an optimization – do not optimize

"If you optimize a function that returns a monster fraction, you haven't optimized anything – you've just moved the monster elsewhere."

"You are not optimizing a function in a vacuum; you are optimizing the entire library, where every result will sooner or later become someone else's argument."

FilesExpand file tree

benchmark_results.md

Latest commit

History