From 0ac5d47fa6cd4a283bd80d4551d5c581eb794bb0 Mon Sep 17 00:00:00 2001 From: Martino Ferrari Date: Mon, 28 Sep 2026 00:05:15 +0200 Subject: [PATCH] perf(data): absolute-space pyramid so acquire() is O(produced), not O(window) Frame cost was O(capacity x levels) once the ring saturated: acquire() saw the window base move (it always moves once the ring wraps), cleared the pyramid and re-appended all 1,048,576 retained samples, and append() walked every live level with no early stop -- 13,647,936 bucket merges and 73.8 ms per frame at a 1M-sample window, scaling with the window, which is what capped the 1 MS/s demo at 32-58 ms/frame while it was still filling. Index the pyramid in absolute bucket space instead. Level l bucket k covers the samples [k * (kBucket << l), (k + 1) * (kBucket << l)), so a bucket's span never depends on where the window starts and the window can slide over the aggregate instead of forcing it to be rebuilt. Each level keeps its complete buckets in a deque plus the absolute base of its first entry; append() folds a sample into the single in-progress level-0 bucket and, only when that bucket completes, cascades upward, stopping at the first level that gains no new entry; a new drop_before() slides the window by discarding leading buckets and remembers where the window begins. acquire() is now: drop to the new base, then append exactly what was produced -- 256 merges and 0.034 ms per frame for 8192 new samples, independent of the window, and memory bounded by the window (32,767 entries after 50M samples through a 1M ring, the same as after 10M). A producer that laps the consumer by more than capacity() between two acquires restarts the window rather than splicing a gap. Interface: `level(l) -> std::span` cannot survive a ring and is replaced by level_base/level_size/bucket in absolute bucket coordinates, plus window_base/window_end/drop_before. The *_buckets queries stay in window coordinates (index 0 = oldest retained sample), so Series keeps its window semantics for callers; its exactness split now aligns in absolute space (window_base() + the ceil/floor) because the bucket grid is absolute. The sums move from per-level prefix sums to per-bucket ones -- the only representation a sliding front can keep exact -- and StreamSeries::y_bounds() becomes reduce(0, size()), exact in O(log n + kBucket) per frame instead of a scan of the window. The pyramid owns the published window now, so there is no second copy of base/end to keep in step, and gather() emits only buckets that are both complete and wholly retained. Tests: test-only merge counters (mirroring PLOT_PYRAMID_COUNT_COVERS, absent from every non-test build) pin an acquire proportional to production and independent of capacity, an append cost amortised over 16 live levels, and bounded memory; brute-force checks pin absolute bucket aggregates across slides, sliding-window reduce/sum/sumsq over unaligned ranges, and the level-1 gather mapping. The demo now draws 946 samples/px at 1.03 ms/frame (967 fps) where the saturated window cost 54.96 ms/frame (18 fps). --- CMakeLists.txt | 10 +- .../plans/2026-09-27-plot-library-scaffold.md | 16 +- .../specs/2026-09-27-plot-library-design.md | 60 +++- include/plot/data/detail/merge_count.hpp | 31 +++ include/plot/data/minmax_pyramid.hpp | 107 +++++-- include/plot/data/series.hpp | 4 +- include/plot/data/stream_series.hpp | 27 +- src/data/minmax_pyramid.cpp | 263 ++++++++++-------- src/data/series.cpp | 53 ++-- src/data/stream_series.cpp | 92 +++--- tests/CMakeLists.txt | 7 +- tests/test_minmax_pyramid.cpp | 157 ++++++++++- tests/test_series.cpp | 186 +++++++++++++ 13 files changed, 769 insertions(+), 244 deletions(-) create mode 100644 include/plot/data/detail/merge_count.hpp diff --git a/CMakeLists.txt b/CMakeLists.txt index 0685e31..e2f7a9e 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -35,10 +35,12 @@ target_link_libraries(plot_core PRIVATE ${CMAKE_DL_LIBS}) plot_apply_warnings(plot_core) if(PLOT_BUILD_TESTS) - # Test-only: counts the pyramid's covering-decomposition entries so a test can - # pin its O(log n) cost, which no query answer reveals. Off in every non-test - # build (see include/plot/data/detail/cover_count.hpp). - target_compile_definitions(plot_core PRIVATE PLOT_PYRAMID_COUNT_COVERS) + # Test-only: counts the pyramid's covering-decomposition entries and its bucket + # merges so a test can pin the query's O(log n) cost and the append path's + # amortised O(1) cost, neither of which any answer reveals. Off in every + # non-test build (see include/plot/data/detail/{cover,merge}_count.hpp). + target_compile_definitions(plot_core PRIVATE PLOT_PYRAMID_COUNT_COVERS + PLOT_PYRAMID_COUNT_MERGES) enable_testing() add_subdirectory(tests) diff --git a/docs/superpowers/plans/2026-09-27-plot-library-scaffold.md b/docs/superpowers/plans/2026-09-27-plot-library-scaffold.md index 0982616..1ab752e 100644 --- a/docs/superpowers/plans/2026-09-27-plot-library-scaffold.md +++ b/docs/superpowers/plans/2026-09-27-plot-library-scaffold.md @@ -1063,12 +1063,20 @@ Take every quantity from `Series` (Task 3), never from the pyramid directly: `s. ```cpp struct PlotLayout { - Rect axes; // plot area (excludes tick labels) - Rect y_label_area, x_label_area; - f32 y_axis_width = 56.0f; - f32 x_axis_height = 22.0f; + Rect axes; // plot area: excludes tick labels AND titles + Rect y_label_area, x_label_area; // one text line each, for tick labels + Rect y_title_area, x_title_area; // one text line each, for axis titles + f32 y_axis_width = 56.0f; + f32 x_axis_height = 22.0f; + f32 y_title_width = 0.0f; // set from the title's measured extent + f32 x_title_height = 0.0f; }; ``` +`layout()` reserves one text line for tick labels and a second for the axis title on +each occupied edge, and insets `axes` by all of them. A title drawn inside `axes` +sits on top of the data and is a defect — the reserved strips exist to make that +impossible rather than merely unlikely. An axis with no title reserves nothing. + `Plot` owns vectors of axes/traces/cursors/measures, computes `PlotLayout` from `viewport()` and the text backend metrics, and recomputes measurement results for every measure whose cursors moved (dirty flag, not every frame). `Grid` produces tick-aligned line segments in data space; `Trace` binds a `Series` to a color/width/axis pair; `Cursor` carries `f64 x`, `SnapMode { Free, NearestSample }`, color and an optional label; `Annotation` is a pinned point with text. diff --git a/docs/superpowers/specs/2026-09-27-plot-library-design.md b/docs/superpowers/specs/2026-09-27-plot-library-design.md index 7fbd50a..c96e3bd 100644 --- a/docs/superpowers/specs/2026-09-27-plot-library-design.md +++ b/docs/superpowers/specs/2026-09-27-plot-library-design.md @@ -82,8 +82,19 @@ approach that fails at high sample rates. Vertex buffers hold raw *data* coordinates. Screen mapping is `x' = x * sx + tx`, `y' = y * sy + ty` evaluated in the vertex shader from a -`Transform` uniform. Panning and zooming therefore touch one uniform and zero -vertex data. Only new samples or an LOD-level change cause an upload. +`Transform` uniform. The CPU never transforms a vertex, and a pan or zoom changes +one uniform and no vertex *values*. + +The staged window is view-dependent: each frame the renderer gathers just the +visible sample range (or the visible buckets of the chosen LOD level) and uploads +it. That per-frame cost is **O(visible pixels)** — about two segments per pixel +column, flat in series length — so the property this section exists for holds. +A staging cache keyed on `(series, level, first, last, viewport)` would let a +static plot skip the redundant gather and upload entirely; it is not implemented, +because it changes no asymptotics and would add an invalidation surface (new +samples, LOD changes, series mutation, streaming window slide) whose bugs are +worse than the redundant work it removes. Recorded here as deliberate future +work, not as a missing feature. Log axes are handled inside the same transform (`log10` before scaling) so no separate code path or pre-transformed buffer is needed. @@ -94,8 +105,18 @@ Every series owns a mip chain: level 0 buckets `kBucket = 64` consecutive sample and stores one `(min, max)` pair per bucket; each subsequent level halves the bucket count. Memory cost is `n/8` bytes at level 0 (2 floats per 64 samples) times the geometric sum `1 + ½ + ¼ + …`, giving `n/4` bytes total — ≈6% of the 4n bytes of raw samples, not 2×. -- Built incrementally: `append()` updates only the last bucket and the O(log n) - chain above it — amortised O(1) per sample. +- Built incrementally, in **absolute bucket space**: level `l` bucket `k` covers + samples `[k · (kBucket << l), (k + 1) · (kBucket << l))`. A bucket's span + therefore never depends on where the current window starts, which is what lets a + window slide over a pyramid instead of forcing it to be rebuilt. `append()` + folds a sample into the single in-progress level-0 bucket and, exactly when that + bucket completes, pushes it and merges upward, stopping at the first level that + gains no new entry: level `l` gains one entry per `kBucket << l` samples, so an + append costs amortised O(1) (`1/kBucket + 1/(2·kBucket) + … = 2/kBucket` merges) + and never walks the live levels. +- A sliding window (a `StreamSeries` ring) drops the leading buckets it no longer + holds and appends only the newly produced samples, so nothing is recomputed and + memory stays bounded by the retained window rather than the total ever pushed. - The renderer picks the level where samples-per-pixel ≈ 2, so a 100M-sample series rasterises ≈ `2 · viewport_width` vertices instead of 100M. - Buckets preserve extremes, so spikes are never dropped by decimation. @@ -127,6 +148,12 @@ otherwise silently include samples outside the range. Exactness therefore lives on `Series` (which owns the raw samples), and `MinMaxPyramid` exposes only bucket-exact queries, named `*_buckets` so they cannot be mistaken for exact ones. +The bucket grid is absolute while a sliding window's base is not, so that split is +computed in absolute coordinates — `Series` asks the pyramid for `window_base()` +and aligns `ceil`/`floor` there before mapping back to window indices. That is +also how a stream's `y_bounds()` stays exact in O(log n + kBucket) per frame: it +is `reduce(0, size())` over the window rather than a scan of it. + ### 4.5 Decimation approach — decision | | Approach | Verdict | @@ -196,9 +223,17 @@ namespace data { class MinMaxPyramid { public: void reserve(u64 n); void append(f32 v); void clear(); - u64 size() const; + // Absolute bucket space: level l entry i covers the absolute samples + // [i * (kBucket << l), (i + 1) * (kBucket << l)). The window slides over the + // same buckets; queries are in window coordinates (index 0 = oldest retained). + void drop_before(u64 absolute_sample); + u64 size() const; // retained window size + u64 window_base() const; // abs index of index 0 + u64 window_end() const; // abs index past the last + u32 level_count() const; + u64 level_base(u32 l) const, level_size(u32 l) const; + MinMax bucket(u32 l, u64 absolute_bucket) const; u32 level_for(u64 n, f64 samples_per_pixel) const; - std::span level(u32 l) const; u64 snap(f64 x, f64 x0, f64 dx) const; // nearest sample index // Whole-bucket only: aggregates over the level-0 buckets fully contained in // [begin, end). Partial edges are excluded, not approximated. Callers wanting @@ -308,8 +343,17 @@ namespace input { - Producers write to `data::StreamSeries` from an acquisition thread through a lock-free SPSC ring. No mutex on the sample path. - A frame begins with `StreamSeries::acquire()`, which publishes a consistent - `[0, produced)` snapshot and extends the pyramid by exactly the newly produced - samples. Readers use the snapshot; the producer never blocks. + `[0, produced)` snapshot: it slides the pyramid's retained window to the new + base, dropping the leading buckets the ring has overwritten, then extends the + pyramid by exactly the newly produced samples. Its cost is **O(samples produced + since the last acquire)** — independent of both the window size and the series + length — and so is its memory, because the drops bound the pyramid to the + retained window instead of the total ever pushed. A frame's total cost is + therefore O(new samples + visible pixels), and no part of the frame path is + proportional to the window. The one exception is a producer that laps the + consumer by more than `capacity()` samples between two acquires: the samples in + between were never appended, so that acquire re-reads the retained window from + scratch. Readers use the snapshot; the producer never blocks. - `Renderer`, `Plot`, and `Controller` are confined to the render thread. ## 8. Error handling diff --git a/include/plot/data/detail/merge_count.hpp b/include/plot/data/detail/merge_count.hpp new file mode 100644 index 0000000..0612f72 --- /dev/null +++ b/include/plot/data/detail/merge_count.hpp @@ -0,0 +1,31 @@ +#pragma once + +#ifdef PLOT_PYRAMID_COUNT_MERGES + +#include + +#include + +namespace plot::data::detail { + +/// Test-only instrumentation. Counts the bucket merges the append path has +/// performed since the counter was last reset: one for the bucket pushed at +/// level 0 when it completes, plus one for every parent bucket written above it. +/// +/// The append path's *cost* is invisible in every answer -- a pyramid that walks +/// every live level per sample returns exactly the same min/max/sum as one that +/// stops at the first level that gained no entry -- so a test that wants to pin +/// "amortised O(1) per sample, independent of the number of live levels" has to +/// count merges. The counter exists only when PLOT_PYRAMID_COUNT_MERGES is +/// defined (the test build defines it), so every normal build carries no counter +/// and no increment on the append path. +/// +/// Not part of the query surface and not meant for production use. +inline std::atomic& merge_count() { + static std::atomic count{0}; + return count; +} + +} // namespace plot::data::detail + +#endif // PLOT_PYRAMID_COUNT_MERGES diff --git a/include/plot/data/minmax_pyramid.hpp b/include/plot/data/minmax_pyramid.hpp index 6b65439..9c60ddb 100644 --- a/include/plot/data/minmax_pyramid.hpp +++ b/include/plot/data/minmax_pyramid.hpp @@ -1,5 +1,5 @@ #pragma once -#include +#include #include #include @@ -16,13 +16,33 @@ struct MinMax { /// series and the level maths is a shift. inline constexpr u32 kBucket = 64; -/// Incremental min/max LOD pyramid. +/// Incremental min/max LOD pyramid, indexed in absolute bucket space. /// -/// Level 0 buckets kBucket consecutive samples and keeps one (min, max) per -/// bucket; level l buckets 2^l level-0 buckets. Memory is n/4 bytes -- about 6% -/// of 4-byte raw samples. append() merges into the last bucket at each level and -/// walks up the O(log n) chain above it, so building is amortised O(1) per sample -/// and nothing is ever recomputed from raw samples. +/// Level `l` bucket `k` covers the absolute samples +/// `[k * (kBucket << l), (k + 1) * (kBucket << l))`. Bucket boundaries are +/// therefore fixed for the pyramid's whole life: appending never re-buckets an +/// existing aggregate, and the window a caller reads can slide over the same +/// buckets instead of being rebuilt from the raw samples. +/// +/// Each level holds its complete buckets in a deque (O(1) push_back *and* +/// pop_front) plus a `base`, the absolute level-`l` bucket index of its first +/// entry. `append()` folds a sample into the in-progress level-0 bucket and, when +/// that bucket completes, pushes it and walks up -- stopping at the first level +/// that gains no new entry. Level `l` gains one entry per `kBucket << l` samples, +/// so an append costs amortised O(1) and never walks every live level. +/// +/// `drop_before(u)` slides the window's leading edge to absolute sample `u`: each +/// level discards the buckets whose span starts before `u`, so the memory is +/// bounded by the retained window (about n/4 bytes, ~6% of 4-byte raw samples) +/// and not by the number of samples ever appended. Every stored entry is the +/// complete aggregate of its own span and every stored span lies inside the +/// retained window, so no aggregate is ever recomputed or invalidated. +/// +/// Queries are in *window* coordinates: index 0 is the oldest retained sample, +/// `size()` is the retained sample count, and `window_base()` is the absolute +/// index window index 0 refers to. That is the same window `Series` publishes, +/// which is why a `StreamSeries` can hand the pyramid a moving base and keep +/// answering in its own coordinates. /// /// A pyramid stores *aggregates*, so it can only answer exactly for whole /// buckets; a range with partial edge buckets would silently pull in samples @@ -30,23 +50,49 @@ inline constexpr u32 kBucket = 64; /// *_buckets, and exact range queries live on Series, which owns the raw samples. class MinMaxPyramid { public: - MinMaxPyramid() { add_level(); } + MinMaxPyramid() = default; - /// Sizes every level (buckets and prefix sums) for `n` samples so a following - /// run of appends performs no reallocation. The pyramid itself stays empty. + /// Sizes the level table for a window of `n` samples so a following run of + /// appends never grows the level vector. The bucket deques themselves grow in + /// chunks: a sliding window cannot know their final length up front. void reserve(u64 n); void append(f32 v); + /// Slides the retained window's leading edge to absolute sample index + /// `absolute_sample`, dropping every aggregate that covered a sample before it. + /// Monotonic for arguments ahead of the base; a no-op for arguments at or + /// behind it, and for an empty pyramid. + /// + /// An argument past `window_end()` (a producer that outran the consumer by more + /// than the ring held) drops everything and restarts the window there: the + /// samples in between were never appended, so they are not representable. + void drop_before(u64 absolute_sample); void clear(); - u64 size() const noexcept { return size_; } + /// Samples in the retained window. + u64 size() const noexcept { return count_ - base_; } + /// Absolute index of the sample window index 0 refers to. + u64 window_base() const noexcept { return base_; } + /// Absolute index one past the last appended sample: the first sample index a + /// following append() will occupy. + u64 window_end() const noexcept { return count_; } /// At least 1, even when empty, so level_count() - 1 is always valid. u32 level_count() const noexcept { return live_levels_; } - /// Level `l`'s buckets; empty span when `l` is beyond level_count(). - std::span level(u32 l) const; - /// Min/max over the whole level-0 buckets fully contained in [begin, end). - /// Partial edge buckets are excluded, never approximated; an empty result is a - /// degenerate Range with span() == 0. + /// Absolute level-`l` bucket index of level `l`'s first entry, and how many + /// entries it holds. They cover `[level_base(l), level_base(l) + level_size(l))`: + /// entry `i` is the aggregate of the absolute samples + /// `[i * (kBucket << l), (i + 1) * (kBucket << l))`. A level that exists but + /// holds nothing, or a level past level_count(), reports size 0 (its base is + /// then meaningless). + u64 level_base(u32 l) const noexcept; + u64 level_size(u32 l) const noexcept; + /// Aggregate of the absolute level-`l` bucket `i`. Precondition: + /// level_base(l) <= i < level_base(l) + level_size(l). + MinMax bucket(u32 l, u64 i) const; + + /// Min/max over the whole level-0 buckets fully contained in the window range + /// [begin, end). Partial edge buckets are excluded, never approximated; an + /// empty result is a degenerate Range with span() == 0. Range reduce_buckets(u64 begin, u64 end) const; /// Sum of the samples in those same whole buckets, accumulated in double. f64 sum_buckets(u64 begin, u64 end) const; @@ -62,13 +108,32 @@ class MinMaxPyramid { u64 snap(f64 x, f64 x0, f64 dx) const; private: - void add_level(); + /// One bucket: its extremes plus the raw-sum aggregates Series folds over. + struct Entry { + MinMax mm; + f64 sum = 0.0; + f64 sumsq = 0.0; + }; + /// One level's complete buckets, in absolute bucket-index order. + struct Level { + std::deque entries; + u64 base = 0; ///< absolute level-l bucket index of entries.front() + }; - u64 size_ = 0; - u64 reserved_n_ = 0; + void add_level(); + /// Appends `e` at absolute level-`l` bucket `k`, which must be one past the + /// level's last entry. + void push(u32 l, u64 k, const Entry& e); + /// Entry for absolute level-`l` bucket `i`; the caller has checked the range. + const Entry& at(u32 l, u64 i) const { + return levels_[l].entries[static_cast(i - levels_[l].base)]; + } + + u64 count_ = 0; ///< absolute index one past the last appended sample + u64 base_ = 0; ///< absolute index of the window's oldest sample u32 live_levels_ = 1; - std::vector> levels_; // levels_[0] = kBucket granularity - std::vector> sums_, sumsqs_; // inclusive prefix sums per level + std::vector levels_{1}; + Entry pending_; ///< the in-progress level-0 bucket, not yet an entry }; } // namespace plot::data diff --git a/include/plot/data/series.hpp b/include/plot/data/series.hpp index 43fe9e8..3b4b5f9 100644 --- a/include/plot/data/series.hpp +++ b/include/plot/data/series.hpp @@ -65,7 +65,9 @@ class Series { u64 left_begin = 0, left_end = 0; // partial head bucket u64 right_begin = 0, right_end = 0; // partial tail bucket }; - static QueryParts split_query(u64 begin, u64 end); + /// Needs the pyramid's window base: the bucket grid it aligns to is absolute, + /// while begin/end are window coordinates. + QueryParts split_query(u64 begin, u64 end) const; /// Reads both edge spans through copy_raw() into `out` (which must hold at /// least 2*kBucket samples) and returns how many were written. u64 copy_edges(const QueryParts& parts, f32* out) const; diff --git a/include/plot/data/stream_series.hpp b/include/plot/data/stream_series.hpp index 817ee8d..e80b2c9 100644 --- a/include/plot/data/stream_series.hpp +++ b/include/plot/data/stream_series.hpp @@ -17,13 +17,16 @@ namespace plot::data { /// /// Indexing is the *retained window published by the last acquire()*: size(), /// reduce(), sum(), sumsq(), gather() and pyramid() all use it, and index 0 is -/// the oldest sample still held. Once the producer has written more than -/// capacity() samples the oldest ones are gone, so acquire() notices that the -/// window's base moved and rebuilds the pyramid from the retained samples -- an -/// incremental append would splice a dropped prefix onto a retained suffix and -/// mis-align every bucket (bucket boundaries are absolute, not relative to the -/// window). Until the ring saturates, acquire() is O(new samples); once it wraps, -/// an acquire whose base moved costs O(capacity()). +/// the oldest sample still held. +/// +/// The pyramid behind those queries is indexed in absolute sample space, so a +/// window that slides only makes it drop the leading buckets the ring has +/// overwritten and append the samples produced since the last acquire(): acquire() +/// costs O(samples produced), not O(capacity()), and the aggregates that stay in +/// the window are never recomputed. The one exception is a producer that laps the +/// consumer -- more than capacity() samples between two acquires -- where the +/// samples in between were never appended and the pyramid restarts at the new +/// base, so that acquire re-reads the whole window. /// /// x0() is therefore the x of the oldest retained sample, x_bounds().min == x0(), /// and y_bounds() covers exactly the published samples. @@ -42,8 +45,9 @@ namespace plot::data { /// one of those k newer samples (or a still-newer sample, once the producer has /// lapped again) instead of the sample the snapshot names. Nothing else moves: the /// count, the absolute index -> x mapping and the newest end of the window are -/// derived from the one published count, and the next acquire() rebuilds from the -/// base it observes, so the pyramid never mixes two windows. Ingest with no reuse +/// derived from one published count, and the next acquire() slides the pyramid's +/// window to the base it observes (dropping the leading buckets, appending what is +/// new), so the pyramid never mixes two windows. Ingest with no reuse /// underneath a read at all needs the producer to stop before unread data /// (drop-on-full) or a claim protocol, i.e. a different contract from the one this /// class specifies. @@ -76,8 +80,9 @@ class StreamSeries final : public Series { RingBuffer ring_; f64 x0_ = 0.0; f64 dx_ = 1.0; - u64 base_ = 0; // absolute index of retained index 0 in the published snapshot - u64 end_ = 0; // absolute index one past the published snapshot + // The pyramid owns the published window: its window_base() is the absolute + // index of retained index 0 and its size()/window_end() are the snapshot's + // extent, so there is no second copy of that state to keep in step. MinMaxPyramid pyramid_; Range y_bounds_{}; }; diff --git a/src/data/minmax_pyramid.cpp b/src/data/minmax_pyramid.cpp index 573e13e..8300d12 100644 --- a/src/data/minmax_pyramid.cpp +++ b/src/data/minmax_pyramid.cpp @@ -1,8 +1,10 @@ #include +#include #include #include #include +#include namespace plot::data { namespace { @@ -16,18 +18,19 @@ constexpr u64 buckets_for(u64 n, u32 l) { return n / span + (n % span == 0 ? u64{0} : u64{1}); } +/// Ceiling division: first level-l bucket whose span starts at or after v. +constexpr u64 ceil_div(u64 v, u64 d) { return v / d + (v % d == 0 ? u64{0} : u64{1}); } + /// First level-0 bucket boundary at or after v. -constexpr u64 ceil_to_bucket(u64 v) { - return (v / kBucket + (v % kBucket == 0 ? u64{0} : u64{1})) * kBucket; -} +constexpr u64 ceil_to_bucket(u64 v) { return ceil_div(v, kBucket) * kBucket; } /// Canonical, coarsest-first decomposition of the level-0 bucket span [b0, b1) -/// into pyramid entries, calling fn(level, index) for each. +/// into pyramid entries, calling fn(level, absolute_bucket_index) for each. /// /// A level-l entry is emitted only when the whole level-l bucket exists (all /// kBucket << l samples are present, so its aggregate is final). Callers -/// guarantee b1 <= ceil(size() / kBucket), so the only possible remainder is the -/// single trailing partial level-0 bucket, which is emitted at level 0. +/// guarantee b0 >= the first whole bucket of the window and b1 <= the number of +/// complete buckets, so every emission both exists and lies inside the window. /// /// The level is chosen per *emission*, not once per pass: at each cursor the /// largest level is taken whose bucket is aligned with the cursor, lies wholly @@ -40,16 +43,18 @@ constexpr u64 ceil_to_bucket(u64 v) { /// no coarser level is ever reachable from an odd cursor. template void cover_buckets(const MinMaxPyramid& p, u64 b0, u64 b1, Fn&& fn) { - const u64 complete = p.size() / kBucket; - const u64 limit = std::min(b1, complete); + // One past the last complete level-0 bucket: levels sit on top of those, so + // this also bounds every coarser bucket the cursor may take. + const u64 limit = std::min(b1, p.level_base(0) + p.level_size(0)); const u32 top = p.level_count() - 1; // level_count() is at least 1 u64 b = b0; while (b < b1) { u32 li = top; for (;;) { const u64 lb = u64{1} << li; - if ((b % lb) == 0 && b + lb <= limit && - (b / lb) < static_cast(p.level(li).size())) { + const u64 i = b >> li; + if ((b % lb) == 0 && b + lb <= limit && i >= p.level_base(li) && + i - p.level_base(li) < p.level_size(li)) { break; } if (li == 0) break; @@ -69,7 +74,7 @@ Range fold_range(const MinMaxPyramid& p, u64 b0, u64 b1) { Range r{}; bool first = true; cover_buckets(p, b0, b1, [&](u32 l, u64 i) { - const MinMax m = p.level(l)[static_cast(i)]; + const MinMax m = p.bucket(l, i); if (first) { r.min = f64(m.min); r.max = f64(m.max); @@ -84,140 +89,170 @@ Range fold_range(const MinMaxPyramid& p, u64 b0, u64 b1) { } // namespace -void MinMaxPyramid::add_level() { - levels_.emplace_back(); - sums_.emplace_back(); - sumsqs_.emplace_back(); -} +void MinMaxPyramid::add_level() { levels_.emplace_back(); } void MinMaxPyramid::reserve(u64 n) { - if (n > reserved_n_) reserved_n_ = n; - const u64 want = std::max(std::max(n, size_), u64{2}); + // Materialise every level n samples can need -- level 0 plus one per doubling + // of the bucket -- so a run of appends never grows the level vector. See the + // header for why the bucket deques themselves cannot be pre-sized. + u32 want = 1; + for (u64 span = level_span(0); span <= n && want < 63; span <<= 1) ++want; + if (levels_.size() < want) levels_.resize(want); +} - // Materialise every level n samples can need -- buckets *and* prefix sums -- - // and size each one up front, so a run of appends allocates nothing at all. - u32 max_levels = 1; - for (u64 span = level_span(0); span <= want && max_levels < 63; span <<= 1) ++max_levels; - while (levels_.size() < max_levels) add_level(); - for (u32 l = 0; l < levels_.size(); ++l) { - const usize want_l = static_cast(buckets_for(want, l)); - levels_[l].reserve(want_l); - sums_[l].reserve(want_l); - sumsqs_[l].reserve(want_l); - } +void MinMaxPyramid::push(u32 l, u64 k, const Entry& e) { + Level& lv = levels_[l]; + // The first entry of a level fixes its base; every later entry is the next + // absolute bucket, which is what keeps a level one contiguous run. + if (lv.entries.empty()) lv.base = k; + lv.entries.push_back(e); +#ifdef PLOT_PYRAMID_COUNT_MERGES + detail::merge_count().fetch_add(1, std::memory_order_relaxed); +#endif } void MinMaxPyramid::append(f32 v) { - const u64 prev = size_; const f64 vd = f64(v); - // sum_[l].back() is the sum over every sample appended so far, and the newest - // sample always lands in the last bucket of every level, so the last prefix - // entry of every level is exactly the running total. - const f64 total_before = sums_[0].empty() ? 0.0 : sums_[0].back(); - const f64 total = total_before + vd; - const f64 sq_before = sumsqs_[0].empty() ? 0.0 : sumsqs_[0].back(); - const f64 total_sq = sq_before + vd * vd; - - const u64 idx = prev / kBucket; - if (idx == levels_[0].size()) { - levels_[0].push_back(MinMax{v, v}); - sums_[0].push_back(total); - sumsqs_[0].push_back(total_sq); + // The new sample always lands in the in-progress level-0 bucket, which is the + // only bucket any append ever touches: a bucket becomes an entry exactly once, + // complete, and is never rewritten. + if (count_ % kBucket == 0) { + pending_.mm.min = v; + pending_.mm.max = v; + pending_.sum = vd; + pending_.sumsq = vd * vd; } else { - MinMax& cur = levels_[0][static_cast(idx)]; - cur.min = std::min(cur.min, v); - cur.max = std::max(cur.max, v); - sums_[0][static_cast(idx)] = total; - sumsqs_[0][static_cast(idx)] = total_sq; + pending_.mm.min = std::min(pending_.mm.min, v); + pending_.mm.max = std::max(pending_.mm.max, v); + pending_.sum += vd; + pending_.sumsq += vd * vd; + } + ++count_; + if (count_ % kBucket != 0) return; // the bucket is still partial: nothing to push + + const u64 k0 = count_ / kBucket - 1; // the level-0 bucket that just completed + if (k0 < level_base(0)) return; // its span starts before the window: never read + const Entry done = pending_; + push(0, k0, done); + + // Cascade: level l gains an entry exactly when kBucket << l divides the sample + // count, so the walk stops at the first level that gained nothing -- amortised + // O(1) per sample (1/kBucket + 1/(2*kBucket) + ... = 1/(kBucket/2) merges). + for (u32 l = 1; l < 63; ++l) { + const u64 span = level_span(l); + if (count_ % span != 0) break; + const u64 k = count_ / span - 1; // the level-l bucket that just completed + const u64 child = 2 * k; // its left child, at level l-1 + if (child < level_base(l - 1)) break; // a child is gone: this parent is not retained + if (l == static_cast(levels_.size())) add_level(); + if (l >= live_levels_) live_levels_ = l + 1; + const Entry& a = at(l - 1, child); + const Entry& b = at(l - 1, child + 1); + Entry m; + m.mm.min = std::min(a.mm.min, b.mm.min); + m.mm.max = std::max(a.mm.max, b.mm.max); + m.sum = a.sum + b.sum; + m.sumsq = a.sumsq + b.sumsq; + push(l, k, m); + } +} + +void MinMaxPyramid::drop_before(u64 absolute_sample) { + if (absolute_sample <= base_) return; + if (absolute_sample > count_) { + // The producer produced more than the ring held between two acquires: every + // appended sample is gone and the ones in between were never appended, so the + // pyramid restarts empty at the new base rather than splicing a gap. + base_ = absolute_sample; + count_ = absolute_sample; + pending_ = Entry{}; + live_levels_ = 1; + for (u32 l = 0; l < levels_.size(); ++l) { + levels_[l].entries.clear(); + levels_[l].base = ceil_div(absolute_sample, level_span(l)); + } + return; } - for (u32 l = 1; l < 63; ++l) { - // A level only carries data once the level below has more than one bucket: a - // one-child parent would just duplicate it. - if (levels_[l - 1].size() < 2) break; - if (l == live_levels_) { - ++live_levels_; - if (l == levels_.size()) { // no reserve() sized this level up front - add_level(); - if (reserved_n_ != 0) { - const usize want = static_cast(buckets_for(reserved_n_, l)); - levels_[l].reserve(want); - sums_[l].reserve(want); - sumsqs_[l].reserve(want); - } - } - } - const u64 ci = prev / level_span(l - 1); // current bucket at level l-1 - const u64 parent = ci / 2; // its parent at level l - MinMax m; - if ((ci & 1u) != 0) { - const MinMax& a = levels_[l - 1][static_cast(ci - 1)]; - const MinMax& b = levels_[l - 1][static_cast(ci)]; - m.min = std::min(a.min, b.min); - m.max = std::max(a.max, b.max); + base_ = absolute_sample; + for (u32 l = 0; l < levels_.size(); ++l) { + Level& lv = levels_[l]; + // A bucket is dropped once its span starts before the window: the first + // retained bucket of every level is the first one wholly inside the window, + // which is exactly what the queries assume. + const u64 keep = ceil_div(base_, level_span(l)); + if (keep <= lv.base) continue; + const u64 have = static_cast(lv.entries.size()); + const u64 drop = std::min(keep - lv.base, have); + lv.entries.erase(lv.entries.begin(), + lv.entries.begin() + static_cast(drop)); + // The base is where the next entry goes: past the last one if the level was + // emptied, which is how the cursor catches up with the sliding window. + if (drop == have) { + lv.base = keep; } else { - m = levels_[l - 1][static_cast(ci)]; - } - if (parent == levels_[l].size()) { - levels_[l].push_back(m); - sums_[l].push_back(total); - sumsqs_[l].push_back(total_sq); - } else { - levels_[l][static_cast(parent)] = m; - sums_[l][static_cast(parent)] = total; - sumsqs_[l][static_cast(parent)] = total_sq; + lv.base += drop; } } - size_ = prev + 1; + // The in-progress level-0 bucket may now start before the window's base: its + // aggregate covers samples that are gone. It belongs to the bucket the window + // starts inside, which is never read (the queries take whole buckets and the + // Series scans the two partial edges raw), so discarding it is lossless. + if (ceil_div(base_, kBucket) > count_ / kBucket) pending_ = Entry{}; + while (live_levels_ > 1 && level_size(live_levels_ - 1) == 0) --live_levels_; } void MinMaxPyramid::clear() { - size_ = 0; + count_ = 0; + base_ = 0; live_levels_ = 1; - for (auto& lv : levels_) lv.clear(); - for (auto& s : sums_) s.clear(); - for (auto& s : sumsqs_) s.clear(); + levels_.resize(1); + levels_[0].entries.clear(); + levels_[0].base = 0; + pending_ = Entry{}; } -std::span MinMaxPyramid::level(u32 l) const { - if (l >= level_count()) return {}; - return levels_[l]; +u64 MinMaxPyramid::level_base(u32 l) const noexcept { + return l < levels_.size() ? levels_[l].base : u64{0}; } +u64 MinMaxPyramid::level_size(u32 l) const noexcept { + return l < levels_.size() ? static_cast(levels_[l].entries.size()) : u64{0}; +} + +MinMax MinMaxPyramid::bucket(u32 l, u64 i) const { return at(l, i).mm; } + Range MinMaxPyramid::reduce_buckets(u64 begin, u64 end) const { - const u64 e = std::min(end, size_); - const u64 b = std::min(begin, e); - const u64 b0 = ceil_to_bucket(b) / kBucket; // first whole bucket inside - const u64 b1 = e / kBucket; // one past the last whole bucket + const u64 e = std::min(end, size()); + if (begin >= e) return Range{0.0, 0.0}; + // The bucket grid is absolute while these are window coordinates: window index + // i is absolute window_base() + i, so the first whole bucket inside the range is + // the first absolute bucket boundary at or after that. + const u64 b0 = ceil_to_bucket(base_ + begin) / kBucket; + const u64 b1 = (base_ + e) / kBucket; if (b1 <= b0) return Range{0.0, 0.0}; return fold_range(*this, b0, b1); } f64 MinMaxPyramid::sum_buckets(u64 begin, u64 end) const { - const u64 e = std::min(end, size_); - const u64 b = std::min(begin, e); - const u64 b0 = ceil_to_bucket(b) / kBucket; - const u64 b1 = e / kBucket; + const u64 e = std::min(end, size()); + if (begin >= e) return 0.0; + const u64 b0 = ceil_to_bucket(base_ + begin) / kBucket; + const u64 b1 = (base_ + e) / kBucket; if (b1 <= b0) return 0.0; f64 total = 0.0; - cover_buckets(*this, b0, b1, [&](u32 l, u64 i) { - const f64 prev = i == 0 ? 0.0 : sums_[l][static_cast(i - 1)]; - total += sums_[l][static_cast(i)] - prev; - }); + cover_buckets(*this, b0, b1, [&](u32 l, u64 i) { total += at(l, i).sum; }); return total; } f64 MinMaxPyramid::sumsq_buckets(u64 begin, u64 end) const { - const u64 e = std::min(end, size_); - const u64 b = std::min(begin, e); - const u64 b0 = ceil_to_bucket(b) / kBucket; - const u64 b1 = e / kBucket; + const u64 e = std::min(end, size()); + if (begin >= e) return 0.0; + const u64 b0 = ceil_to_bucket(base_ + begin) / kBucket; + const u64 b1 = (base_ + e) / kBucket; if (b1 <= b0) return 0.0; f64 total = 0.0; - cover_buckets(*this, b0, b1, [&](u32 l, u64 i) { - const f64 prev = i == 0 ? 0.0 : sumsqs_[l][static_cast(i - 1)]; - total += sumsqs_[l][static_cast(i)] - prev; - }); + cover_buckets(*this, b0, b1, [&](u32 l, u64 i) { total += at(l, i).sumsq; }); return total; } @@ -229,7 +264,7 @@ u32 MinMaxPyramid::level_for(u64 n, f64 samples_per_pixel) const { const f64 lg = std::floor(std::log2(ratio)); u32 l = lg >= f64(62) ? 62u : static_cast(lg); while (l > 0) { - const bool exists = l < level_count() && level(l).size() >= 2; + const bool exists = l < level_count() && level_size(l) >= 2; if (exists && buckets_for(n, l) >= 2) break; --l; } @@ -237,14 +272,14 @@ u32 MinMaxPyramid::level_for(u64 n, f64 samples_per_pixel) const { } u64 MinMaxPyramid::snap(f64 x, f64 x0, f64 dx) const { - if (dx <= 0.0 || size_ == 0) return 0; + if (dx <= 0.0 || size() == 0) return 0; const f64 i = (x - x0) / dx; // Bound the round before it: llround() is undefined for values a long long // cannot hold, and this also maps NaN to index 0. if (!(i > 0.0)) return 0; - if (i >= f64(size_)) return size_ - 1; + if (i >= f64(size())) return size() - 1; const u64 idx = static_cast(std::llround(i)); - return idx >= size_ ? size_ - 1 : idx; + return idx >= size() ? size() - 1 : idx; } } // namespace plot::data diff --git a/src/data/series.cpp b/src/data/series.cpp index cf4c893..73ad2aa 100644 --- a/src/data/series.cpp +++ b/src/data/series.cpp @@ -12,23 +12,30 @@ constexpr u64 kGatherChunk = 256; } // namespace -Series::QueryParts Series::split_query(u64 begin, u64 end) { +Series::QueryParts Series::split_query(u64 begin, u64 end) const { QueryParts p; - // First level-0 bucket boundary at or after begin, and the last at or before - // end: [interior_begin, interior_end) is exactly the part the pyramid can - // answer with whole buckets. - const u64 up = (begin / kBucket + (begin % kBucket == 0 ? u64{0} : u64{1})) * kBucket; - const u64 down = (end / kBucket) * kBucket; + // The pyramid's buckets are absolute while these are window coordinates, so its + // level-0 boundaries sit at window index `window_base() rounded up to the next + // multiple of kBucket`, and every kBucket after that. Align the split in + // absolute coordinates and come back to the window: aligning in window + // coordinates would treat index 0 as a boundary and, for a window that slid + // into the middle of a bucket, route samples through the pyramid that it does + // not cover. + const u64 origin = pyramid().window_base(); + const u64 a_begin = origin + begin; + const u64 a_end = origin + end; + const u64 up = (a_begin / kBucket + (a_begin % kBucket == 0 ? u64{0} : u64{1})) * kBucket; + const u64 down = (a_end / kBucket) * kBucket; if (up < down) { - p.interior_begin = up; - p.interior_end = down; + p.interior_begin = up - origin; // up >= a_begin >= origin + p.interior_end = down - origin; } // The two edge spans never overlap and never exceed one bucket each: when the // interior is empty they tile [begin, end) together, otherwise they are the // head and tail buckets half inside the range. p.left_begin = begin; - p.left_end = std::min(end, up); - p.right_begin = std::max(down, p.left_end); + p.left_end = std::min(end, up - origin); + p.right_begin = std::max(down > origin ? down - origin : u64{0}, p.left_end); p.right_end = end; return p; } @@ -132,17 +139,25 @@ void Series::gather_impl(u32 level, u64 begin, u64 end, std::vector& out_xy return; } - // Level >= 1: one envelope pair per bucket, x at the bucket's first sample - // (min) and at its last sample (max), which keeps x monotonic while covering - // the bucket's span. + // Level >= 1: one envelope pair per *retained* bucket the range overlaps, x at + // the bucket's first sample (min) and at its last sample (max), which keeps x + // monotonic while covering the bucket's span. The pyramid's buckets are + // absolute, so map the window range onto its grid and clamp to the buckets that + // are both complete and wholly inside the window: a bucket that is only partly + // retained aggregates samples the window no longer has, and a partial trailing + // bucket is not stored at all. const u64 span = u64{kBucket} << level; - const std::span buckets = pyramid().level(level); - const u64 first = b / span; - const u64 last = std::min((e - 1) / span + 1, static_cast(buckets.size())); + const u64 origin = pyramid().window_base(); + const u64 lbase = pyramid().level_base(level); + const u64 lend = lbase + pyramid().level_size(level); + if (lbase >= lend) return; + const u64 a_end = origin + e; + const u64 first = std::max((origin + b) / span, lbase); + const u64 last = std::min((a_end - 1) / span + 1, lend); for (u64 bi = first; bi < last; ++bi) { - const MinMax m = buckets[static_cast(bi)]; - const u64 first_sample = bi * span; - const u64 last_sample = std::min(first_sample + span, n) - 1; + const MinMax m = pyramid().bucket(level, bi); + const u64 first_sample = bi * span - origin; // >= 0: stored buckets start inside the window + const u64 last_sample = std::min((bi + 1) * span, a_end) - 1 - origin; out_xy.push_back(f32(x_0 + f64(first_sample) * step)); out_xy.push_back(m.min); out_xy.push_back(f32(x_0 + f64(last_sample) * step)); diff --git a/src/data/stream_series.cpp b/src/data/stream_series.cpp index 788dd55..42add8e 100644 --- a/src/data/stream_series.cpp +++ b/src/data/stream_series.cpp @@ -6,14 +6,18 @@ namespace plot::data { namespace { -/// Samples copied per ring read on the ingest path, so rebuilding a saturated -/// window costs no allocation. +/// Samples copied per ring read on the ingest path, so appending the samples +/// produced between two frames costs no allocation. constexpr u64 kReadChunk = 256; } // namespace StreamSeries::StreamSeries(u64 capacity, f64 x0, f64 dx) - : ring_(capacity), x0_(x0), dx_(dx) {} + : ring_(capacity), x0_(x0), dx_(dx) { + // Size the pyramid's level table for a full window up front: the ring never + // holds more than capacity() samples, so appending never grows it. + pyramid_.reserve(capacity); +} void StreamSeries::push(f32 y) { ring_.push(y); } @@ -27,62 +31,35 @@ void StreamSeries::acquire() { const u64 p = ring_.produced(); const u64 cap = ring_.capacity(); const u64 base = p > cap ? p - cap : 0; // oldest index the ring still holds - const u64 n = p - base; // retained, published sample count - if (base != base_) { - // The window slid: samples the pyramid covered are gone, and because bucket - // boundaries are absolute, extending it incrementally would mis-align every - // bucket. Rebuild from the retained samples instead. - pyramid_.clear(); - pyramid_.reserve(n); - f32 lo = 0.0f, hi = 0.0f; - bool first = true; - f32 chunk[kReadChunk]; - for (u64 i = 0; i < n;) { - const u64 take = std::min(kReadChunk, n - i); - ring_.copy_from(base + i, take, chunk); - for (u64 k = 0; k < take; ++k) { - const f32 v = chunk[k]; - pyramid_.append(v); - if (first) { - lo = v; - hi = v; - first = false; - } else { - lo = std::min(lo, v); - hi = std::max(hi, v); - } - } - i += take; - } - y_bounds_ = first ? Range{0.0, 0.0} : Range{f64(lo), f64(hi)}; - base_ = base; - } else { - // Contiguous: only the samples produced since the last acquire() are new. - if (n == 0) y_bounds_ = Range{0.0, 0.0}; - f32 chunk[kReadChunk]; - for (u64 i = end_ - base_; i < n;) { - const u64 take = std::min(kReadChunk, n - i); - ring_.copy_from(base + i, take, chunk); - for (u64 k = 0; k < take; ++k) { - const f32 v = chunk[k]; - pyramid_.append(v); - if (i + k == 0) { - y_bounds_ = Range{f64(v), f64(v)}; - } else { - y_bounds_.min = std::min(y_bounds_.min, f64(v)); - y_bounds_.max = std::max(y_bounds_.max, f64(v)); - } - } - i += take; - } + // Slide the pyramid's window: the buckets that covered samples the ring has + // overwritten are dropped, and the rest stay exactly as they were built -- + // their boundaries are absolute, so nothing has to be recomputed. This is also + // where a producer that lapped the consumer by more than the ring holds is + // handled: the pyramid restarts at the new base. + pyramid_.drop_before(base); + + // Append exactly the samples produced since the last acquire(). Normally that + // is one frame's worth; after a lapped producer it is the whole window. + f32 chunk[kReadChunk]; + for (u64 a = pyramid_.window_end(); a < p;) { + const u64 take = std::min(kReadChunk, p - a); + ring_.copy_from(a, take, chunk); + for (u64 k = 0; k < take; ++k) pyramid_.append(chunk[k]); + a += take; } - end_ = p; + + // Exact bounds for the published window, in the same O(log n + kBucket) as any + // other exact query over it: the bucket-aligned interior comes from the + // pyramid and the two partial edge spans are scanned through copy_raw(). An + // incrementally extended bound cannot work here -- it would have to know + // whether the extreme it is holding has just slid out of the window. + y_bounds_ = reduce(0, size()); } -u64 StreamSeries::size() const { return end_ - base_; } +u64 StreamSeries::size() const { return pyramid_.size(); } -f64 StreamSeries::x0() const { return x0_ + f64(base_) * dx_; } +f64 StreamSeries::x0() const { return x0_ + f64(pyramid_.window_base()) * dx_; } f64 StreamSeries::dx() const { return dx_; } @@ -105,9 +82,10 @@ void StreamSeries::copy_raw(u64 begin, u64 end, f32* out) const { const u64 n = size(); const u64 b = std::min(begin, n); const u64 e = std::min(std::max(end, b), n); - // Retained index b is absolute base_ + b: fixed for the whole copy, so a - // pushing producer cannot shift the window midway. - ring_.copy_from(base_ + b, e - b, out); + // Retained index b is absolute window_base() + b: fixed for the whole copy, so + // a pushing producer cannot shift the window midway (the pyramid is written by + // this thread alone, so its window cannot move under the read either). + ring_.copy_from(pyramid_.window_base() + b, e - b, out); std::fill(out + static_cast(e - b), out + static_cast(want), 0.0f); } diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 78741b3..a7267f5 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -40,9 +40,10 @@ add_executable(plot_tests egl_headless.cpp ) target_link_libraries(plot_tests PRIVATE plot::core ${plot_egl_library} Threads::Threads) -# Same switch the library is built with: the decomposition-cost test reads the -# counter behind it (see include/plot/data/detail/cover_count.hpp). -target_compile_definitions(plot_tests PRIVATE PLOT_PYRAMID_COUNT_COVERS) +# Same switches the library is built with: the cost tests read the counters +# behind them (see include/plot/data/detail/{cover,merge}_count.hpp). +target_compile_definitions(plot_tests PRIVATE PLOT_PYRAMID_COUNT_COVERS + PLOT_PYRAMID_COUNT_MERGES) plot_apply_warnings(plot_tests) add_test(NAME plot_tests COMMAND plot_tests) diff --git a/tests/test_minmax_pyramid.cpp b/tests/test_minmax_pyramid.cpp index 46a3af1..fd5cf2f 100644 --- a/tests/test_minmax_pyramid.cpp +++ b/tests/test_minmax_pyramid.cpp @@ -1,12 +1,36 @@ #include "harness.hpp" #include +#include #include +#include #include using plot::data::kBucket; +using plot::data::MinMax; using plot::data::MinMaxPyramid; using plot::Range; +namespace { + +/// Brute-force aggregate of the absolute span of level-`l` bucket `b`, from the +/// same value function the pyramid was fed. This is the oracle the absolute +/// bucket indexing is judged against: a bucket that has slid inside the window +/// must still aggregate exactly its own absolute samples. +template +MinMax brute_bucket(u32 l, u64 b, const ValueFn& value) { + const u64 span = u64{kBucket} << l; + const u64 first = b * span; + MinMax m{value(first), value(first)}; + for (u64 i = first + 1; i < first + span; ++i) { + const f32 v = value(i); + m.min = std::min(m.min, v); + m.max = std::max(m.max, v); + } + return m; +} + +} // namespace + PLOT_TEST(pyramid_bucket_extremes_preserved) { MinMaxPyramid p; p.reserve(1000); @@ -107,7 +131,8 @@ PLOT_TEST(pyramid_clear_and_rebuild) { // bucket and ask the bucket-grained query for its extrema. for (u32 i = 1; i < kBucket; ++i) p.append(1.0f); PLOT_CHECK(p.size() == kBucket); - PLOT_CHECK(p.level(0).size() == 1); + PLOT_CHECK(p.level_base(0) == 0); + PLOT_CHECK(p.level_size(0) == 1); PLOT_CHECK_NEAR(p.reduce_buckets(0, kBucket).max, 7.0, 1e-6); PLOT_CHECK_NEAR(p.reduce_buckets(0, kBucket).min, 1.0, 1e-6); } @@ -121,7 +146,7 @@ PLOT_TEST(pyramid_cover_cost_is_logarithmic_for_unaligned_starts) { MinMaxPyramid p; p.reserve(1'000'000); for (u64 i = 0; i < 1'000'000; ++i) p.append(0.0f); - PLOT_CHECK(p.level(0).size() == 1'000'000 / kBucket); // 15625, no partial bucket + PLOT_CHECK(p.level_size(0) == 1'000'000 / kBucket); // 15625, no partial bucket auto folds_for = [&](u64 begin, u64 end) { plot::data::detail::cover_entry_count().store(0, std::memory_order_relaxed); @@ -134,3 +159,131 @@ PLOT_TEST(pyramid_cover_cost_is_logarithmic_for_unaligned_starts) { PLOT_CHECK(unaligned < 100); PLOT_CHECK(unaligned <= u64{4} * aligned); // misalignment costs entries, not a scan } + +PLOT_TEST(pyramid_drop_before_slides_the_window) { + MinMaxPyramid p; + for (u64 i = 0; i < 4096; ++i) p.append(f32(i % 13)); + PLOT_CHECK(p.window_base() == 0); + PLOT_CHECK(p.window_end() == 4096); + PLOT_CHECK(p.size() == 4096); + + p.drop_before(1000); + PLOT_CHECK(p.window_base() == 1000); + PLOT_CHECK(p.size() == 3096); + // The window's first *whole* bucket starts at 1024: the bucket the drop cut + // into is not representable as a whole bucket and is never queried. + PLOT_CHECK(p.level_base(0) == 1024 / kBucket); + PLOT_CHECK(p.level_base(0) + p.level_size(0) == 4096 / kBucket); + + f32 lo = 1.0f, hi = -1.0f; + for (u64 i = 1024; i < 4096; ++i) { + lo = std::min(lo, f32(i % 13)); + hi = std::max(hi, f32(i % 13)); + } + const Range r = p.reduce_buckets(0, p.size()); + PLOT_CHECK_NEAR(r.min, lo, 1e-6); + PLOT_CHECK_NEAR(r.max, hi, 1e-6); + // Window [1000, 1024) holds no whole bucket at all: an empty answer, not one + // that quietly reaches into the dropped samples. + PLOT_CHECK(p.reduce_buckets(0, 24).span() == 0.0); + PLOT_CHECK(p.sum_buckets(0, 24) == 0.0); + PLOT_CHECK(p.sumsq_buckets(0, 24) == 0.0); + + p.drop_before(1000); // at the base: no-op + PLOT_CHECK(p.window_base() == 1000); + p.drop_before(4096); // at the end: the window is empty + PLOT_CHECK(p.size() == 0); + p.append(5.0f); // and appending starts a fresh window + PLOT_CHECK(p.size() == 1); + PLOT_CHECK(p.window_base() == 4096); +} + +PLOT_TEST(pyramid_bucket_aggregates_stay_absolute_across_slides) { + // The invariant the sliding window rests on: a bucket that has slid inside the + // window is never recomputed, so it must already hold the aggregate of its own + // absolute span -- at every level, after any number of drops. + MinMaxPyramid p; + const auto value = [](u64 i) { return f32((i * 37u) % 101u) - 50.0f; }; + constexpr u64 kTotal = 20'000; + constexpr u64 kWindow = 4'000; + for (u64 i = 0; i < kTotal; ++i) { + p.append(value(i)); + if ((i + 1) % 997 == 0 && i + 1 > kWindow) p.drop_before(i + 1 - kWindow); + } + // The last drop is the last multiple of 997 at or below kTotal. + PLOT_CHECK(p.window_base() == (kTotal / 997) * 997 - kWindow); + for (u32 l = 0; l < p.level_count(); ++l) { + PLOT_CHECK(p.level_size(l) > 0); + const u64 span = u64{kBucket} << l; + for (u64 b = p.level_base(l); b < p.level_base(l) + p.level_size(l); ++b) { + const MinMax want = brute_bucket(l, b, value); + const MinMax got = p.bucket(l, b); + PLOT_CHECK_NEAR(got.min, want.min, 1e-6); + PLOT_CHECK_NEAR(got.max, want.max, 1e-6); + PLOT_CHECK(b * span >= p.window_base()); // wholly inside the window + PLOT_CHECK((b + 1) * span <= p.window_end()); // and complete + } + } +} + +PLOT_TEST(pyramid_sliding_window_bucket_queries_match_bruteforce) { + MinMaxPyramid p; + const auto value = [](u64 i) { return f32(i % 89) * 0.25f - 5.0f; }; + constexpr u64 kTotal = 12'000; + constexpr u64 kWindow = 3'000; + for (u64 i = 0; i < kTotal; ++i) { + p.append(value(i)); + if ((i + 1) % 613 != 0) continue; + const u64 abs0 = i + 1 > kWindow ? i + 1 - kWindow : 0; + p.drop_before(abs0); + const u64 n = p.size(); + // Window coordinates: window index j is absolute abs0 + j. + const u64 first = abs0 / kBucket + (abs0 % kBucket == 0 ? u64{0} : u64{1}); + const u64 last = (abs0 + n) / kBucket; + if (last <= first) continue; + f32 lo = 1e30f, hi = -1e30f; + f64 sum = 0.0, sq = 0.0; + for (u64 a = first * kBucket; a < last * kBucket; ++a) { + const f32 v = value(a); + lo = std::min(lo, v); + hi = std::max(hi, v); + sum += f64(v); + sq += f64(v) * f64(v); + } + const u64 wb = first * kBucket - abs0; + const u64 we = last * kBucket - abs0; + const Range r = p.reduce_buckets(wb, we); + PLOT_CHECK_NEAR(r.min, lo, 1e-6); + PLOT_CHECK_NEAR(r.max, hi, 1e-6); + PLOT_CHECK_NEAR(p.sum_buckets(wb, we), sum, 1e-3); + PLOT_CHECK_NEAR(p.sumsq_buckets(wb, we), sq, 1e-1); + } +} + +PLOT_TEST(pyramid_append_cost_is_amortised_constant) { + // A sample must not walk the live levels: level l gains an entry once per + // kBucket << l samples, so N appends perform ~N/32 merges however many levels + // the pyramid holds. Walking every level per append costs ~N * levels. + MinMaxPyramid p; + p.reserve(1u << 20); + constexpr u64 kAppends = 1u << 20; + for (u64 i = 0; i < kAppends; ++i) p.append(0.5f); + PLOT_CHECK(p.level_count() >= 10); // 2^20 samples: 16384 level-0 buckets + + plot::data::detail::merge_count().store(0, std::memory_order_relaxed); + for (u64 i = 0; i < kAppends; ++i) p.append(0.5f); + const u64 merges = plot::data::detail::merge_count().load(std::memory_order_relaxed); + // ~N/64 + N/128 + ... = ~N/32; the bound leaves 2x room, and is four orders of + // magnitude below the N * live_levels an unbounded walk would perform. + PLOT_CHECK(merges <= kAppends / 16); + + // Sliding must not add per-sample work either: the appends below run against a + // window whose base moves with every sample. + plot::data::detail::merge_count().store(0, std::memory_order_relaxed); + for (u64 i = 0; i < kAppends; ++i) { + p.drop_before(i + 1); + p.append(0.5f); + } + const u64 slid = plot::data::detail::merge_count().load(std::memory_order_relaxed); + PLOT_CHECK(slid <= kAppends / 16); +} diff --git a/tests/test_series.cpp b/tests/test_series.cpp index e3ccf3f..3316ea5 100644 --- a/tests/test_series.cpp +++ b/tests/test_series.cpp @@ -1,13 +1,18 @@ #include "harness.hpp" #include +#include #include #include +#include #include +#include #include #include #include #include +using plot::data::kBucket; + PLOT_TEST(array_series_bounds_and_gather) { std::vector y{1, 5, -2, 3}; plot::data::ArraySeries s(std::move(y), /*x0=*/10.0, /*dx=*/0.5); @@ -159,3 +164,184 @@ PLOT_TEST(stream_series_push_and_read_are_race_free) { PLOT_CHECK(reads >= 1); PLOT_CHECK(sizes_ok); } + +namespace { + +/// Buckets the pyramid holds across every level: its footprint, not its window. +u64 pyramid_entries(const plot::data::MinMaxPyramid& p) { + u64 total = 0; + for (u32 l = 0; l < p.level_count(); ++l) total += p.level_size(l); + return total; +} + +} // namespace + +PLOT_TEST(stream_series_sliding_window_queries_match_bruteforce) { + // After the base has moved many times, an exact query must still match a scan + // of exactly the retained samples. The pyramid's buckets live in absolute + // sample space while the window's base is arbitrary, and that meeting point is + // where an off-by-one hides -- so test the sliding case, not just the static + // ArraySeries one. + constexpr u64 kCapacity = 1u << 14; + constexpr u64 kPushes = 60'000; + const auto value = [](u64 i) { return f32((i * 29u) % 173u) - 86.0f; }; + plot::data::StreamSeries s(kCapacity, 0.0, 1e-6); + for (u64 i = 0; i < kPushes; i += 337) { + const u64 n = std::min(337, kPushes - i); + for (u64 k = 0; k < n; ++k) s.push(value(i + k)); + s.acquire(); // the window's base moves on every one of these + } + + const u64 n = s.size(); + const u64 base = s.produced() - n; + PLOT_CHECK(n == kCapacity); + PLOT_CHECK(base > 2 * kCapacity); // the window really did slide + PLOT_CHECK(s.pyramid().size() == n); // pyramid() indexes the window + PLOT_CHECK(s.pyramid().window_base() == base); + PLOT_CHECK(s.pyramid().window_end() == s.produced()); + PLOT_CHECK_NEAR(s.x0(), f64(base) * 1e-6, 1e-12); + PLOT_CHECK(s.x_bounds().min == s.x0()); + PLOT_CHECK_NEAR(s.x_bounds().max, f64(base + n - 1) * 1e-6, 1e-9); + + // Brute force in window coordinates over the retained samples only. + const auto brute = [&](u64 a, u64 b) { + f64 lo = f64(value(base + a)), hi = lo; + f64 sum = 0.0, sq = 0.0; + for (u64 i = a; i < b; ++i) { + const f64 v = f64(value(base + i)); + lo = std::min(lo, v); + hi = std::max(hi, v); + sum += v; + sq += v * v; + } + return std::array{lo, hi, sum, sq}; + }; + + const u64 starts[] = {0, 1, 63, 64, 65, 1000, 7919, n - 200}; + for (const u64 a : starts) { + for (u64 b = a + 1; b <= n; b += 1301) { + const std::array want = brute(a, b); + const plot::Range r = s.reduce(a, b); + PLOT_CHECK_NEAR(r.min, want[0], 1e-4); + PLOT_CHECK_NEAR(r.max, want[1], 1e-4); + PLOT_CHECK_NEAR(s.sum(a, b), want[2], 1e-2); + PLOT_CHECK_NEAR(s.sumsq(a, b), want[3], 1e-1); + } + } + const std::array whole = brute(0, n); + PLOT_CHECK_NEAR(s.y_bounds().min, whole[0], 1e-6); + PLOT_CHECK_NEAR(s.y_bounds().max, whole[1], 1e-6); + + // The level accessors are in absolute bucket coordinates: a stored bucket must + // aggregate its own absolute span, the sliding window included. + const u64 span = u64{kBucket} << 1; + const u64 lbase = s.pyramid().level_base(1); + const u64 lsize = s.pyramid().level_size(1); + PLOT_CHECK(lsize >= 2); + for (u64 b = lbase; b < lbase + std::min(lsize, 4); ++b) { + const std::array want = brute(b * span - base, (b + 1) * span - base); + const plot::data::MinMax m = s.pyramid().bucket(1, b); + PLOT_CHECK_NEAR(m.min, want[0], 1e-6); + PLOT_CHECK_NEAR(m.max, want[1], 1e-6); + } + + // level-1 gather: one envelope pair per retained bucket the range overlaps, x + // at the bucket's own first and last sample. This is the mapping the renderer + // draws from. + std::vector out; + s.gather(1, 0, n, out); + const u64 first = std::max(base / span, lbase); + const u64 last = std::min((base + n - 1) / span + 1, lbase + lsize); + PLOT_CHECK(out.size() == 4 * (last - first)); + for (u64 bi = first; bi < last; ++bi) { + const std::array want = brute(bi * span - base, (bi + 1) * span - base); + const plot::usize k = static_cast(bi - first); + const u64 first_sample = bi * span - base; + const u64 last_sample = std::min((bi + 1) * span, base + n) - 1 - base; + // x is in data coordinates: x0() is the retained window's first sample, so + // window index j plots at (base + j) * dx. + PLOT_CHECK_NEAR(out[4 * k], f32(f64(base + first_sample) * 1e-6), 1e-9); + PLOT_CHECK_NEAR(out[4 * k + 1], want[0], 1e-6); + PLOT_CHECK_NEAR(out[4 * k + 2], f32(f64(base + last_sample) * 1e-6), 1e-9); + PLOT_CHECK_NEAR(out[4 * k + 3], want[1], 1e-6); + } +} + +PLOT_TEST(stream_series_steady_state_acquire_is_proportional_to_produced) { + // The defect this pins: a saturated ring whose base moved every frame used to + // rebuild the whole window, costing O(capacity * live_levels) merges per + // acquire. It must cost O(samples produced since the last acquire) instead, so + // the same production through a 16x larger window must cost the same. + const auto merges_for = [](u64 capacity) { + plot::data::StreamSeries s(capacity, 0.0, 1.0); + constexpr u64 kSaturate = 1u << 20; + for (u64 i = 0; i < kSaturate; ++i) s.push(f32(i % 977)); + s.acquire(); + PLOT_CHECK(s.size() == capacity); + + constexpr u64 kFrames = 64; + constexpr u64 kFrame = 1u << 13; // 8192 samples: 1 MS/s at ~122 fps + plot::data::detail::merge_count().store(0, std::memory_order_relaxed); + for (u64 f = 0; f < kFrames; ++f) { + for (u64 k = 0; k < kFrame; ++k) s.push(f32((kSaturate + f * kFrame + k) % 977)); + s.acquire(); + } + return plot::data::detail::merge_count().load(std::memory_order_relaxed); + }; + + constexpr u64 kProduced = 64 * (1u << 13); + const u64 small = merges_for(1u << 16); // 64k-sample window + const u64 big = merges_for(1u << 20); // 1M-sample window: 16x to rebuild + // ~kProduced/32 in practice; rebuilding the 1M window alone would merge + // kCapacity/64 * live_levels ~ 245'000 buckets per acquire. + PLOT_CHECK(small <= kProduced / 16); + PLOT_CHECK(big <= kProduced / 16); + PLOT_CHECK(big <= 2 * small + kBucket); +} + +PLOT_TEST(stream_series_steady_state_acquire_is_fast) { + // Wall-clock backstop for the merge count above: rebuilding the 1M-sample + // window measured ~44 ms, so 128 such frames would need ~5.6 s. The sliding + // path does the same work in a few milliseconds; the bound is two orders of + // magnitude clear of both, so it is not a thin timing margin. + plot::data::StreamSeries s(1u << 20, 0.0, 1.0); + constexpr u64 kSaturate = 1u << 20; + for (u64 i = 0; i < kSaturate; ++i) s.push(f32(i % 977)); + s.acquire(); + + constexpr u64 kFrames = 128; + constexpr u64 kFrame = 1u << 13; + const auto t0 = std::chrono::steady_clock::now(); + for (u64 f = 0; f < kFrames; ++f) { + for (u64 k = 0; k < kFrame; ++k) s.push(f32((kSaturate + f * kFrame + k) % 977)); + s.acquire(); + } + const f64 ms = + std::chrono::duration(std::chrono::steady_clock::now() - t0).count() * 1.0e3; + PLOT_CHECK(s.size() == 1u << 20); + PLOT_CHECK(ms < 250.0); +} + +PLOT_TEST(stream_series_pyramid_memory_is_bounded_by_the_window) { + // 50M samples through a 1M-sample ring: the pyramid's footprint must track the + // retained window (n/64 level-0 buckets plus the levels over them), not the + // samples ever pushed. + constexpr u64 kCapacity = 1u << 20; + constexpr u64 kTotal = 50u * (1u << 20); + constexpr u64 kChunk = 1u << 19; + plot::data::StreamSeries s(kCapacity, 0.0, 1.0); + u64 entries_after_10m = 0; + for (u64 p = 0; p < kTotal; p += kChunk) { + for (u64 i = 0; i < kChunk; ++i) s.push(f32((p + i) % 251)); + s.acquire(); + if (p + kChunk == 10u * (1u << 20)) entries_after_10m = pyramid_entries(s.pyramid()); + } + const u64 entries = pyramid_entries(s.pyramid()); + const u64 bound = 2 * (kCapacity / kBucket) + kBucket; + PLOT_CHECK(s.produced() == kTotal); + PLOT_CHECK(s.size() == kCapacity); + PLOT_CHECK(entries <= bound); + // 5x the samples pushed, same footprint: levels deeper than one chunk can + // differ by a bucket each, hence the slack. + PLOT_CHECK(entries <= entries_after_10m + kBucket); +}