chore(merge): pull request #309 from fiddlerwoaroof/docs/felo4-adr-incremental-evaluation

docs(adr): ADR 009 incremental evaluation (proposed) (improvise-felo.4)
This commit is contained in:
Ed Langley
2026-10-09 22:50:58 -07:00
committed by GitHub
3 changed files with 906 additions and 0 deletions
+447
View File
@@ -0,0 +1,447 @@
# ADR 009 — Incremental evaluation: an in-house DBSP core, plugged in behind a trait
**Status:** proposed. ·
**Ticket:** improvise-felo.4 (the record of this ADR; implementation goes to its own
epic) ·
**Relates to:** ADR 001 (the basis a write claims), ADR 002 (surrogates on the wire),
ADR 003, ADR 005, ADR 006, ADR 007, ADR 008
## Problem
Every change recomputes every formula cell (`Model::recompute_formulas`, called by
`rebuild_layout` after each batch), and raw aggregates are re-folded on every read. That
does not scale to a live source (ADR 003) streaming changes into a large table, nor to
ADR 008's virtualized grid, which needs occupancy, totals and widths maintained without
building the grid. A change should cost in proportion to what it touches.
## Decision
Formulas are evaluated incrementally by an in-house implementation of the DBSP core. A
cell write becomes a delta, aggregates are maintained, and only the formula cells whose
inputs moved are re-evaluated by the existing evaluator. The incremental results must
equal today's eager results exactly.
**Principle: everything is as incremental as possible**, so the work done per model
mutation is minimised. Where this ADR does work that grows with more than the
mutation, it says so and names the incremental path that normally applies.
### 1. The `zset` crate
A problem-agnostic workspace crate with no dependencies beyond `std` (so it is
wasm32-clean and serde-free). It provides:
- `ZSet<K>`, a map from `K: Ord` to a non-zero `i64` weight, which is an abelian group;
- a synchronous, single-threaded `Circuit` per process, with `delay`, `integrate`,
`differentiate`, `map`, `filter`, `flat_map`, `index_by`, the bilinear delta join
`join_indexed`, `distinct` and `aggregate`;
- resumable rebuilds: building an operator from a snapshot is a `step(deadline)`, so it
can run as a sliced job;
- operator state behind one trace interface (look up, read a range, apply a delta), so
later store-memory work can page or compress it.
**Operators are self-contained.** An operator's state and behaviour depend only on its
inputs and its own state, never on a shared in-process global or on another operator's
internals. That keeps the graph distributable: in principle the model could run across a
cluster of machines. Nothing more is designed for that here.
There are no nested streams, no recursion operator and no timestamps beyond one clock.
**The formula language is the query language.** What a request asks of the model beyond
its axes and window (filters, aggregation, derived values) is written in formula text
over ids (ADR 007). It goes through the same evaluator and circuit as a document
formula, so there is no second query language.
### 2. The plug
improvise-core gains a trait, `Derived`, for the two things the evaluator reads: the
formula result at a key, and the raw fold of the cells that match a partial key. The
incremental backend implements it in a seam crate, `improvise-incremental`, which
depends on `zset` and improvise-core. improvise-core never depends on `zset`.
A live source (ADR 003) that reports its changes feeds cell deltas like `DataStore`. A
source that answers aggregate queries but reports no changes implements the derived-value
operator directly, by query pushdown, with no change feed.
The circuit decides *which* formula cells to re-evaluate. The one existing evaluator
computes each value, so formulas are never compiled into operators.
**`Eager` is temporary.** Today's code moves behind the trait as `Eager`, unchanged. It
exists only to show discrepancies between the new evaluator and the old one: a
differential property test applies one random edit sequence to both backends and
compares every drawn cell, total and occupancy count. Eventually it is deleted, and
nothing depends on it surviving. Long-term tests compare against exact references, such
as `ExactSum` against rational sums, not against `Eager`.
### 3. The cell delta
`DataStore::set` and `DataStore::remove` emit a `ZSet<(CellKey, CellValue)>`:
- `+(k, v)` for a new cell;
- `−(k, v₀) + (k, v₁)` for an overwrite;
- `−(k, v₀)` for a removal.
Values compare by normal float rules, except that every NaN is one canonical value:
`-0.0` equals `+0.0`, and an equal-value overwrite consolidates to nothing. NaN payloads
are not preserved, so nothing can be smuggled through them. Today's record tier compares
bits, so overwriting `0.0` with `-0.0` stops counting as a change. Multisets order values
by `f64::total_cmp` after mapping every NaN to the canonical one and `-0.0` to `+0.0`. `CellKey` gains `Ord`, which is its canonical order: the
natural order of ids (ADR 007 amendment 1, clause 4), the order eager folds in. The
circuit steps once per model effect.
The record tier and the dirty flag read this delta instead of scanning the store. The
dirty flag fires on a non-empty `(key, value)` delta, or when a fingerprint of the
formula set moves, because formula edits move no cell.
Under ADR 007 and its amendment, `CellKey` pairs are `(CategoryId, ItemHandle)`. So a
rename is a display write and emits no delta. The formula-set fingerprint must use the id
form too, or a rename reads as a formula change.
### 4. Aggregates
A group is the cell key minus the None categories. One aggregation circuit serves every
shape key with the same None list, which removes the per-session model clone.
| aggregate | state per group | per delta element |
|---|---|---|
| `Count` | count of values, count of numbers | O(1) |
| `Sum`, `Avg` | an exact accumulator (clause 8), and the count | O(1) |
| `Min`, `Max`, text unanimity | a multiset of values | O(log d) |
A group is read with eager's rules, including NaN and signed zero. `LAST`/`FIRST` and
semi-additive measures (improvise-felo.7, felo.8) reuse the multiset shape, keyed by the
ordering category's order. Those tickets choose the key. It must keep a reorder or an
insert incremental: only the entries under the moved items change, and nothing is
rebuilt.
### 5. Formulas
**Stems** are `distinct` over the delta, keyed by "cell key minus the target and None
categories". A formula counts at a stem where its top-level `WHERE` matches or raw data
exists at the target.
**Inputs.** A formula cell is re-evaluated when one of its inputs moved:
- a `Ref(R)`: the group at the stem with `R`'s category set to `R`;
- an `Agg`: the group at the stem;
- an `Agg` whose inner `WHERE` names `Cat = Item`: the group at the stem with `Cat` set
to `Item`. This is the one join: stems indexed by "stem minus `Cat`".
The affected pairs are evaluated in the formula graph's topological order.
**An aggregate over an expression** (`SUM(Revenue * Price)`) is a pipeline:
1. A per-record stage retracts and reinserts each touched record. It reads the record's
previous cells as the store minus the delta, so it keeps no copy of the values. The
store is the `integrate` operator's trace over the cell delta, read as an input of
the stage, so the stage stays self-contained.
2. A map produces `(context, order key, expr)`. The order key is (context rank, record
key in `CellKey` order), which is eager's concatenation order.
3. The group aggregate, with an error channel that reads the first `Err` in order.
No cross product is built. A formula whose expression reads other records is evaluated
over the **occupied contexts** only, and `Eager` uses the same enumeration. The answer is
unchanged, because a context with no records contributes nothing. One exception: an
inner bare ref that names no item errors at every context, so for it the ref is checked
before enumerating. This replaces `Expansion::members`.
### 6. Circular references
DBSP's least fixed point is not improvise's circularity: `X = X + 1` has none, and the
answer must be `Error("circular")`. The formula graph's strongly connected components
are maintained incrementally as formulas change. Acyclic formulas are evaluated once
each, in topological order. A cyclic component is maintained incrementally too: a change
re-iterates only the cells it reaches. A component that does not converge within the
existing bound of sixteen passes (`MAX_EVAL_DEPTH`) propagates `Error("circular")`, as
today. An optimisation that speeds convergence is allowed only where it is exact: it
gives the same result as the plain iteration, and a test compares the two.
### 7. Operator state follows the request's shape key
A projection request's **shape key** selects operator state, shared by every request
with the same key, from any client. The key is:
- the categories on each axis, by level;
- the None categories;
- the filter keys: each filter category's item, and the hidden items of categories on
None, Page or Filter;
- the page categories, so each page is a maintained slice and paging rebuilds nothing;
- `measure_agg`.
The key uses the canonical parsed id form of the request's query text, so two spellings
of one query share their operator state. Everything else in a request is applied when output is read: hidden and collapsed items
on axis categories, item order, the page selection and prune. Occupancy and column
aggregates are keyed by axis **items**, so a collapse, hide, show, reorder or added item
changes no operator. The visible axis is a remap.
- **A new key** builds its operators as a sliced job under ADR 008's job runner. Deltas
arriving during the build are buffered and added at install. The finished result is the
response to the request that asked for it, not a notice.
- **A formula edit** is a formula delta fed to the circuit.
- **A document swap** rebuilds everything, and a new key's first build is a full pass.
These are the only rebuilds; every other change is a delta.
- **Saving a view definition** touches no operator.
- **Live keys are a cache.** A key is evicted after a period in which no request names
it. The period is a local duration on the host, set by measurement, and never compared
with another machine's clock. Every eviction is announced (clause 9).
**Column totals and hidden rows.** Totals are kept over all rows, because the shape key
is shared and a hide is one client's. Whether hidden items travel in a request is ADR
008's open question. If they do not, a total is read as it is kept. If they do, a total
must leave the hidden rows out, so a read subtracts the request's H hidden cells as a
temporary overlay:
- `Count` subtracts;
- `Min` and `Max` skip keys whose whole count is hidden;
- text counts the distinct keys left;
- `Sum` subtracts exactly through the accumulator.
That is O(H log d) per read for every aggregate, with no view-store state in the model
store.
Aggregate-over-expression totals index their entries by row tuple for the same
subtraction.
### 8. Exact sums
An f64 sum maintained by `+v` and `−v` drifts. Measured by
`crates/improvise-core/examples/retraction_drift.rs`, at 18,000 cells after a million
cycles, it is 1.7e-9 from eager, past `FLOAT_EQ_EPSILON`, so `=` could answer
differently on the two paths. So every `Sum` and `Avg` uses **`ExactSum`**, in both
backends:
- a fixed-point integer covering every finite f64, with counts of NaN and ±∞, which
makes it an abelian group that rounds once, on read;
- measured bit-equal to a fresh sum after a million cycles, and exact against rational
sums on 20,000 cases;
- 568 bytes per group, and about 1.3× the cost of an f64 add per update.
The sign of a zero sum follows IEEE rules for the same addends: it is `-0.0` only if
every addend is `-0.0`, and an exact cancellation gives `+0.0`. It lives in ADR 005's
`numtower` crate. **User-visible change:** displayed sums change in their last bits,
because they become correctly rounded.
### 9. Delivery
The **model store** (the document model and the DBSP graph) and the **view store** (the
cursor, the window and unsaved view definitions) are conceptual boxes, not a network
boundary. Where a piece of state physically lives (a client, the host or a worker) is
an implementation detail. The model store's decisions depend only on the model. A
server-side twin of a client's state, such as its filter, watermark or window, may be
held at the serving end of an edge where that makes the graph genuinely simpler. It is
then part of the edge, not of the model store, and it is disposable: losing it never
gives wrong data. A client that keeps its own watermark recovers with a correction from
it. A client that relied on the twin for its watermark recovers with a full refetch.
**Every edge may cross a network.** An edge, meaning the stream from one operator to the
next, may be in-process or remote. No edge assumes a function call or shared
memory. An operator that sits at a wire boundary is just one implementation of an edge,
an **edge adapter**. So every edge has the same semantics:
- **A wire form.** Every element type that flows on an edge has a serializable form,
defined in `improvise-protocol`. The adapter converts to and from it, and `zset` stays
free of serde. `CellKey` is not `Serialize`, so on the wire a cell is named by its
tuple of ids, its position or an ADR 002 surrogate.
- **The wire form is lossless.** Nothing is lost in transport, numeric precision above
all:
- a value never crosses as a JSON or JS double. It crosses in ADR 005's canonical text:
an inexact finite value as `#i` with its shortest round-trip spelling, which keeps
the sign of zero, infinities as `inf` and `-inf`, and every NaN as the one token
`nan`;
- positions and counts below 2⁵³ are exact as JSON numbers. Wider integers cross as
strings;
- text crosses byte-exact, with no normalisation;
- the delta is the unit that crosses, never `ExactSum`'s limbs. A receiver that needs
a sum gets the host's exact read or rebuilds it from the exact deltas.
The stamp rule detects lost messages, but nothing would detect a lossy encoding. So
every element type has a round-trip property test, where encoding then decoding is
the identity for every value, including `-0.0`, infinities, 28-digit decimals and
`i128`-wide fractions. A NaN round-trips as the canonical NaN.
- **The clock is causal, never a wall clock.** Ordering, staleness and loss detection
use only the stamp sequence (below), a logical clock the host advances once per
committed command. Each process runs its own synchronous sub-circuit, and across an
edge the stamped commands are the clock. No decision compares wall-clock times across
machines, so clock skew cannot affect it. A delta crosses a wire only after its command
commits, so a remote operator never sees a delta that is later rolled back.
- **Ordering and backpressure.** An edge delivers in order. A slow edge coalesces
exactly: deltas are summed, and notices merged. A merged notice carries the range of
publications it covers, `(from, to]`, so the client that merges notices sees no gap
where nothing was lost.
- **An edge is a session, not a connection.** An edge may be one long-lived connection,
or a sequence of requests carrying the same session token or cookie, for example
several HTTP requests. The host decides nothing from the token: it is transport, and
any later authentication is out of scope. What continues across connections is the
client's own watermark and shape keys, which it sends with each request. Every request
can be retried, and a resync from a watermark is idempotent. A client asks for
"notices after my watermark", and the host answers from **one shared, bounded log
of published notices**: a single log, not one per edge or client, holding the notices
themselves (keys and flags, not values). A request older than the log's start is
answered with a correction (below).
- **Failure.** Timeouts are local monotonic durations, never compared across machines.
On a long-lived connection, both ends keep alive, and a host-side timeout is tied to
that open connection. On a sequence of requests, only the
client detects failure: a request that fails or times out. In either case the client
then asks for a correction from its watermark.
- **Correction is Z-set arithmetic.** A correction is `target − held`, where `target`
is the current state for the client's window and `held` is what the client holds.
Applying it with ordinary Z-set addition gives `held + (target − held) = target`, so
the group structure of `ZSet` provides it, with no separate machinery. The client's
position is a **watermark**: the whole stamp `(host instance, document_generation,
publication)` up to which it has applied the stream, so everything at or before it is
applied. A watermark from another host instance or generation is behind every current
stamp. The client holds one per shape key it reads. A watermark is view-store state:
the model store keeps none, so a client sends its own with each request unless a twin
at the serving end holds it.
- where the shared log covers the watermark and the covered range holds only moved
keys, the host sends the current values of those keys, restricted to the client's
window, and the client applies `target − held` over them. The log holds keys and
flags, not values, so the model store keeps no view-store state;
- where the covered range also holds axis changes, relabels or flags, which positions
and labels depend on, the correction carries them too, in order, as a merged notice
does, so the client re-derives positions incrementally;
- only where the log no longer covers the watermark, or the covered range holds an
oversized notice that lists no keys, does the host fall back to `target`
for the whole window, and the client applies `target − held`;
- a full refetch is `target − ∅`.
**Every read sees only the last committed state**, whether it is a correction, a
fetch or `project`, and it is labelled with that state's stamp. A value a failed
command wrote and then rolled back is never read. A correction is labelled with the
publication range it covers, so it is applied only at the watermark it starts from,
which makes it idempotent and resumable. Its elements are `(key, value)` pairs with
weights, in the lossless wire form, so a client and a worker apply it with the same
code.
Aggregation stays on the model side (ADR 008), so the operators that aggregate run on
the model side. Which other edges cross a wire is a deployment choice. The engine's invalidation signal, which ADR 006's replay trigger
relies on, is a correction trigger below, followed from a displayed cell's key back to
the model.
**A request** is a shape key plus a window, given as an absolute cross-product top-left
and dimensions. The window is view-store state, and a request costs O(visible).
- **Positions** are bigints where they exceed 2⁵³, sent as JSON strings.
- **Output** is keyed by tuple identity, with the position attached when sent.
- **An item insert, removal or reorder** is sent as an axis change, and each client
re-derives its positions.
**Pruning support.** Pruning is window-relative: a row or column empty within the
displayed values is hidden. This replaces today's whole-grid `prune_empty`. The model
answers one restricted successor query: "the next k rows after position p with a value
in a column of S", and the same for columns. A per-column sorted index of occupied rows
answers it in O(|S| log n + k|S|), holding O(occupied cells). Whether pruning runs on
the client and which parameters travel are ADR 008's questions, and this design supports
either answer. If the client prunes, it prunes its window, fills forward with that
query, and, when scrolling up or left, finds a new top-left by galloping backwards in
doubling jumps (about 14 to 28 probes at density 10⁻⁴).
**The change notice.** Each command commits on its own (ADR 006), and after each
committed command the host emits one notice, addressed to no client. A frame
whose later command fails keeps its earlier commands, and their notices, published. A
notice carries:
- the stamp `(host instance, document_generation, publication)`. The publication number
increments once per committed command. The host instance changes on restart; it is a
random id, never a wall-clock time;
- axis changes, display relabels as `(id, display)`, and a flag when the formula set,
`measure_agg`, the view definitions or the default view moved. A shared display's
` #2` suffix is never sent: the client computes it from the displays in scope;
- the shape keys evicted since the last notice;
- for each live shape key, the moved output keys (cells, formula cells and totals), as
the circuit computes them. Moved records alone would miss formulas with cross-stem
inputs.
**A client receives only what it needs.** Filtering is part of the DBSP graph, not a
separate mechanism. A client either pulls just the state it needs, or pushes a filter to
the host's side of its connection: a filter operator on its shape's query, written in the
query language, through which its notices pass. Whether the filter also carries the
client's window is an optimisation, since the model store's decisions depend only on the
model. A pushed filter is held by the open connection, or by a twin. Otherwise the
client sends it with each request. While a connection's filter is open, it counts as a request
naming its shape, so its shape is not evicted.
A notice whose key lists exceed a fixed size carries "correct" instead. That is a
fallback, and it costs every reader a correction.
**What a client does.** It intersects the moved keys with its window. A fold row has no
maintained key, so it maps member keys through its own remap first. It then fetches what
moved. It asks for a correction when it sees a publication gap, a new
`document_generation` or host instance, its shape in the eviction list, a notice it
cannot read, or a dropped connection or failed request. A slow client may merge pending notices: the union of
keys, the covered range of publications, and axis changes kept in order. Because the
wire is lossless, outputs keep the circuit's value identity end to end.
**What this replaces.** Today's per-session mirrors and pushed diffs (`Sessions`,
`ClientSession`, `grid_updates`, `Projection`, `Dirty`). The model store keeps no
view-store state, only the shared notice log. `GridExtent` is derived by the client, and in records
mode the record count is fetched. `Projection` also published the ADR 001 basis a write
claims. Whether `project` returns it instead, and what replaces the layout tier, are
ADR 008's questions.
**Edge adapters.** An edge that crosses a wire has one adapter at each end. Toward the
circuit, the client's adapter is a sink for results and notices and a source of requests.
Toward the client it is the mirror image. The client's end keeps the watermark. In
process, an edge is a function call.
### 10. Column width hints
Every column incrementally maintains a width hint, which optimises grid sizing and
reduces jitter. The hint is an integrated aggregate: a histogram from display width to
count, so a retraction is O(log d). Where in the graph the hint is computed is a
scheduling choice, made as an optimisation later.
### 11. Rollback and memory
- **Rollback.** The circuit's state joins ADR 006's per-command snapshot, which is taken
before a command's effects and restored if the command fails. An inverse delta is not
a rewind, because `delay` and `differentiate` keep history. The requirement is that
this snapshot does not cost O(cells) per command against ADR 008's 8 ms budget, and
that writes do not block reads. How is for the implementation epic.
- **Memory.** Resident state beside the store is the multisets, `ExactSum` per summed
group, occupancy per live shape, the successor-query index, the width histograms and
the bounded notice log. The store itself is the integrated base and is never copied. It
is measured at ADR 008's fixture before the store-memory work begins.
## Alternatives rejected
- **Feldera's `dbsp` crate:** 67 required dependencies, including tokio and zstd, and it
fails to build for wasm32.
- **Compiling formulas into operators:** a second evaluator.
- **An ad-hoc dirty-tracking cache:** one invalidation rule per derived fact, each a
second description of evaluation.
- **A topic or subscription mechanism beside the graph:** filtering is an operator in
the graph (clause 9).
## Consequences
**Good:**
- A write's evaluation costs O(delta) plus the formula cells it reaches.
- The store scans and the per-view model clones go.
- ADR 008 gets occupancy, totals and widths from one evaluator.
- Live sources that report changes plug in.
- The model store holds no view-store state.
**Bad:**
- Two backends exist until `Eager` goes.
- Integrated state is O(cells) beside the store.
- `ExactSum` costs 568 bytes per summed group.
- Each command snapshots the circuit's state for rollback.
**Order** (each step is a ticket in the implementation epic, behind the differential
test):
1. `zset`;
2. the `Derived` trait with `Eager`;
3. the cell delta, the record tier and the dirty flag;
4. order-free aggregates and occupancy, measured at ADR 008's fixture;
5. formulas and the aggregate-over-expression pipeline;
6. `Sum`, `Avg` and totals through `ExactSum`;
7. live sources;
8. delivery.
Step 8's required tests:
- a frame `[AddItem, SetCell]` publishes one notice per committed command, numbered
consecutively, and no client asks for a correction;
- a dropped notice is detected;
- a reconnect after a host restart asks for a correction;
- a corrected state equals the target, `held + (target − held) = target`, and
correcting a state that is already correct is the zero Z-set (property tests over
random edit sequences and losses);
- a refused command consumes no publication number, and the commands before it in its
frame stay published;
- a dropped connection or a failed request triggers a correction, and a later request
with the same session resumes from the watermark.
- a zero sum has the sign IEEE gives the same addends: `-0.0` only if every addend is
`-0.0`, and `+0.0` for an exact cancellation (a property test against IEEE).
## Open questions
None.
@@ -0,0 +1,458 @@
//! How far a SUM maintained by insert/retract deltas drifts from the eager fold,
//! in f64, and what an exact accumulator costs. ADR 009 (improvise-felo.4), clause 8.
use std::hint::black_box;
use std::time::Instant;
use improvise_core::format::format_f64;
use improvise_core::model::cell::{CellKey, CellValue};
use improvise_core::model::types::{Aggregation, Model};
const MEASURE: &str = "_Measure";
/// SplitMix64: seeded, so a run is reproducible and a second seed is a new sample.
struct Rng(u64);
impl Rng {
fn next(&mut self) -> u64 {
self.0 = self.0.wrapping_add(0x9E37_79B9_7F4A_7C15);
let mut z = self.0;
z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9);
z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB);
z ^ (z >> 31)
}
fn below(&mut self, n: u64) -> u64 {
self.next() % n
}
}
/// A user-entered decimal as an exact count of 1e-5 units, its magnitude spread
/// over `decades` decades. The f64 the store holds is `typed(units)`.
fn decimal(rng: &mut Rng, decades: u32) -> i64 {
let scale = 10i64.pow(rng.below(u64::from(decades)) as u32);
let cents = rng.below(100_000) as i64;
let sign = if rng.below(4) == 0 { -1 } else { 1 };
sign * cents * scale
}
fn typed(units: i64) -> f64 {
units as f64 / 1e5
}
/// The exact decimal sum, rounded to f64 once: what ADR 005's arithmetic would show.
fn exact(units: &[i64]) -> f64 {
units.iter().map(|&u| i128::from(u)).sum::<i128>() as f64 / 1e5
}
fn key(row: usize) -> CellKey {
CellKey::new(vec![
("Row".to_string(), format!("r{row}")),
(MEASURE.to_string(), "Amount".to_string()),
])
}
/// `cells` rows of one measure; the eager SUM is the one the grid draws.
fn document(rng: &mut Rng, cells: usize, decades: u32) -> (Model, Vec<i64>) {
let mut m = Model::new("drift");
m.add_category("Row").expect("under the limit");
let values: Vec<i64> = (0..cells).map(|_| decimal(rng, decades)).collect();
for (row, v) in values.iter().enumerate() {
m.category_mut("Row")
.expect("just added")
.add_item(format!("r{row}"));
m.set_cell(key(row), CellValue::Number(typed(*v)));
}
(m, values)
}
fn eager(m: &Model) -> f64 {
let total = CellKey::new(vec![(MEASURE.to_string(), "Amount".to_string())]);
m.evaluate_aggregated_f64(&total, Aggregation::OverHidden)
}
/// Neumaier's compensated running sum: the cheap mitigation, measured beside the naive one.
#[derive(Default)]
struct Compensated {
sum: f64,
carry: f64,
}
impl Compensated {
fn add(&mut self, x: f64) {
let t = self.sum + x;
self.carry += match self.sum.abs() >= x.abs() {
true => (self.sum - t) + x,
false => (x - t) + self.sum,
};
self.sum = t;
}
fn value(&self) -> f64 {
self.sum + self.carry
}
}
const LIMB_BITS: u32 = 32;
/// Every finite f64 is `mant * 2^(shift - 1074)` with `shift <= 2045` and `mant < 2^53`.
const LIMBS: usize = 67;
/// The most weight one limb absorbs between carries: each add moves a limb by < 2^32.
const CHUNK: i64 = 1 << 30;
/// An exact sum of f64s: a signed fixed-point integer in units of 2^-1074, wide enough
/// for every finite f64, with non-finite values counted so they can be retracted too.
#[derive(Clone)]
struct ExactSum {
limbs: [i64; LIMBS],
pending: u64,
nan: i64,
pos_inf: i64,
neg_inf: i64,
}
impl ExactSum {
fn new() -> Self {
ExactSum {
limbs: [0; LIMBS],
pending: 0,
nan: 0,
pos_inf: 0,
neg_inf: 0,
}
}
/// Add `x` with any `weight` (a consolidated Z-set weight). Exact: no rounding here.
fn add(&mut self, x: f64, weight: i64) {
match (x.is_nan(), x.is_infinite()) {
(true, _) => self.nan += weight,
(false, true) if x > 0.0 => self.pos_inf += weight,
(false, true) => self.neg_inf += weight,
(false, false) => {
let mut left = weight;
while left != 0 {
let chunk = left.clamp(-CHUNK, CHUNK);
self.add_finite(x, chunk);
left -= chunk;
}
}
}
}
/// Bounded so each limb stays under 2^62: `|weight| <= CHUNK` and `pending`
/// counts magnitude, not calls.
fn add_finite(&mut self, x: f64, weight: i64) {
if self.pending + weight.unsigned_abs() > CHUNK as u64 {
propagate(&mut self.limbs);
self.pending = 0;
}
let bits = x.to_bits();
let biased = ((bits >> 52) & 0x7ff) as u32;
let fraction = bits & ((1 << 52) - 1);
let (mant, shift) = match biased {
0 => (fraction, 0),
e => (fraction | (1 << 52), e - 1),
};
let sign = if x.is_sign_negative() {
-weight
} else {
weight
};
let wide = u128::from(mant) << (shift % LIMB_BITS);
let k = (shift / LIMB_BITS) as usize;
for i in 0..3 {
let part = ((wide >> (LIMB_BITS * i as u32)) & 0xffff_ffff) as i64;
self.limbs[k + i] += sign * part;
}
self.pending += weight.unsigned_abs();
}
/// The exact sum rounded to the nearest f64, ties to even: the one rounding.
fn value(&self) -> f64 {
if self.nan > 0 || (self.pos_inf > 0 && self.neg_inf > 0) {
return f64::NAN;
}
if self.pos_inf > 0 {
return f64::INFINITY;
}
if self.neg_inf > 0 {
return f64::NEG_INFINITY;
}
let mut limbs = self.limbs;
propagate(&mut limbs);
let negative = limbs[LIMBS - 1] < 0;
if negative {
limbs.iter_mut().for_each(|l| *l = -*l);
propagate(&mut limbs);
}
let magnitude = round_to_f64(&digits(&limbs));
match negative {
true => -magnitude,
false => magnitude,
}
}
}
/// Carry each limb into the next, leaving every limb but the top in [0, 2^32).
fn propagate(limbs: &mut [i64; LIMBS]) {
for k in 0..LIMBS - 1 {
let carry = limbs[k] >> LIMB_BITS;
limbs[k] -= carry << LIMB_BITS;
limbs[k + 1] += carry;
}
}
/// A non-negative normalized accumulator as 32-bit digits, the top limb split out.
fn digits(limbs: &[i64; LIMBS]) -> Vec<u64> {
let mut out: Vec<u64> = limbs.iter().map(|&l| l as u64).collect();
let mut top = out.pop().unwrap_or(0);
loop {
out.push(top & 0xffff_ffff);
top >>= LIMB_BITS;
if top == 0 {
break;
}
}
out
}
fn bit(d: &[u64], pos: usize) -> u64 {
(d[pos / 32] >> (pos % 32)) & 1
}
/// Round an exact multiple of 2^-1074 to the nearest f64, ties to even.
fn round_to_f64(d: &[u64]) -> f64 {
let Some(k) = d.iter().rposition(|&x| x != 0) else {
return 0.0;
};
let msb = 32 * k + 63 - d[k].leading_zeros() as usize;
if msb < 53 {
let int = d[0] | (d[1] << 32);
return int as f64 * f64::from_bits(1);
}
let mut mant = (0..53).fold(0u64, |m, i| m | (bit(d, msb - 52 + i) << i));
let round = bit(d, msb - 53) == 1;
let sticky = (0..msb - 53).any(|p| bit(d, p) == 1);
let mut top = msb;
if round && (sticky || mant & 1 == 1) {
mant += 1;
if mant == 1 << 53 {
mant >>= 1;
top += 1;
}
}
let exp = top as i64 - 52 - 1074;
if exp > 971 {
return f64::INFINITY;
}
let scale = match exp >= -1022 {
true => f64::from_bits(((exp + 1023) as u64) << 52),
false => f64::from_bits(1u64 << (exp + 1074)),
};
mant as f64 * scale
}
fn exact_over(values: &[i64]) -> f64 {
let mut acc = ExactSum::new();
values.iter().for_each(|&v| acc.add(typed(v), 1));
acc.value()
}
fn ulps(a: f64, b: f64) -> u64 {
(a.to_bits() as i64).abs_diff(b.to_bits() as i64)
}
struct Sums {
exact: f64,
eager: f64,
naive: f64,
comp: f64,
xs_maintained: f64,
xs_fresh: f64,
}
/// Each running sum against the eager fold (the differential oracle) and against the
/// exact decimal sum. `flip`: the gap reaches `FLOAT_EQ_EPSILON`, so `=` can answer differently.
fn report(label: &str, cycles: usize, s: &Sums) {
let shown = |v: f64| format_f64(v, true, 2);
println!(
"{label:<22} {cycles:>8} {:>10.2e} {:>6} {:>10.2e} {:>6} {:>5} {:>10.2e} {:>10.2e} {:>10.2e} {:>7} {:>10.2e}",
(s.naive - s.eager).abs(),
ulps(s.naive, s.eager),
(s.comp - s.eager).abs(),
ulps(s.comp, s.eager),
(s.naive - s.eager).abs() >= 1e-10,
(s.eager - s.exact).abs(),
(s.naive - s.exact).abs(),
(s.comp - s.exact).abs(),
s.xs_maintained.to_bits() == s.xs_fresh.to_bits(),
(s.xs_maintained - s.exact).abs(),
);
if shown(s.naive) != shown(s.eager) || shown(s.eager) != shown(s.exact) {
println!(
" shown at ,.2: naive {} eager {} exact {} maintained-exact-sum {}",
shown(s.naive),
shown(s.eager),
shown(s.exact),
shown(s.xs_maintained)
);
}
}
/// Each cycle retracts one live value and inserts a fresh one at the same key.
fn random_updates(seed: u64, cells: usize, decades: u32, checkpoints: &[usize]) {
let mut rng = Rng(seed);
let (mut m, mut values) = document(&mut rng, cells, decades);
let mut naive = eager(&m);
let mut comp = Compensated::default();
comp.add(naive);
let mut xs = ExactSum::new();
values.iter().for_each(|&v| xs.add(typed(v), 1));
let mut done = 0;
for &at in checkpoints {
while done < at {
let row = rng.below(cells as u64) as usize;
let new = decimal(&mut rng, decades);
naive = naive - typed(values[row]) + typed(new);
comp.add(-typed(values[row]));
comp.add(typed(new));
xs.add(typed(values[row]), -1);
xs.add(typed(new), 1);
values[row] = new;
m.set_cell(key(row), CellValue::Number(typed(new)));
done += 1;
}
let sums = Sums {
exact: exact(&values),
eager: eager(&m),
naive,
comp: comp.value(),
xs_maintained: xs.value(),
xs_fresh: exact_over(&values),
};
report(&format!("{cells} cells, {decades} dec"), at, &sums);
}
}
/// A large value arrives, small ones accumulate under it, then it is retracted.
fn absorbed_by_a_large_value() {
let mut m = Model::new("absorb");
m.add_category("Row").expect("under the limit");
m.category_mut("Row").expect("just added").add_item("r0");
let mut naive = 0.0;
let mut comp = Compensated::default();
let mut xs = ExactSum::new();
m.set_cell(key(0), CellValue::Number(1e17));
naive += 1e17;
comp.add(1e17);
xs.add(1e17, 1);
for row in 1..=1000 {
m.category_mut("Row")
.expect("just added")
.add_item(format!("r{row}"));
m.set_cell(key(row), CellValue::Number(0.01));
naive += 0.01;
comp.add(0.01);
xs.add(0.01, 1);
}
m.clear_cell(&key(0));
naive -= 1e17;
comp.add(-1e17);
xs.add(1e17, -1);
let mut fresh = ExactSum::new();
(0..1000).for_each(|_| fresh.add(0.01, 1));
let sums = Sums {
exact: 10.0,
eager: eager(&m),
naive,
comp: comp.value(),
xs_maintained: xs.value(),
xs_fresh: fresh.value(),
};
report("1e17, 1000 x 0.01, -1e17", 1001, &sums);
}
/// Nanoseconds per operation over `n` values, best of three runs.
fn per_op(n: usize, mut op: impl FnMut() -> f64) -> f64 {
(0..3)
.map(|_| {
let t = Instant::now();
black_box(op());
t.elapsed().as_secs_f64() * 1e9 / n as f64
})
.fold(f64::INFINITY, f64::min)
}
fn cost() {
let n = 10_000_000;
let mut rng = Rng(7);
let values: Vec<f64> = (0..n).map(|_| typed(decimal(&mut rng, 9))).collect();
let updates = 3 * n;
let plain = per_op(updates, || {
values
.iter()
.fold(0.0, |s, &v| black_box(s + v) - black_box(v) + v)
});
let exact_add = per_op(updates, || {
let mut xs = ExactSum::new();
values.iter().for_each(|&v| {
xs.add(black_box(v), 1);
xs.add(black_box(v), -1);
xs.add(v, 1);
});
xs.limbs[0] as f64
});
let mut xs = ExactSum::new();
values.iter().for_each(|&v| xs.add(v, 1));
let read = per_op(1, || xs.value());
println!();
println!("cost, {n} values x (insert, retract, insert), best of 3, per single update:");
println!(" plain f64 add/sub {plain:>8.2} ns per update");
println!(" exact accumulator {exact_add:>8.2} ns per update");
println!(" exact read (round) {read:>8.0} ns per read");
println!(
" state per group {} bytes, against 8 for an f64",
std::mem::size_of::<ExactSum>()
);
}
/// A consolidated Z-set weight can be large; the accumulator must stay exact at 2^33.
fn large_weight() {
let x = f64::from_bits((1u64 << 32) - 1);
let mut xs = ExactSum::new();
xs.add(x, 1 << 33);
let want = ((1u64 << 32) - 1) as f64 * 2f64.powi(33) * f64::from_bits(1);
println!();
println!(
"large weight: (2^32-1) x 2^-1074 at weight 2^33 reads {:e}, exact {:e}, equal {}",
xs.value(),
want,
xs.value().to_bits() == want.to_bits()
);
}
fn main() {
println!(
"{:<22} {:>8} {:>10} {:>6} {:>10} {:>6} {:>5} {:>10} {:>10} {:>10} {:>7} {:>10}",
"case",
"cycles",
"naive-eag",
"ulp",
"comp-eag",
"ulp",
"flip",
"eag-exact",
"nai-exact",
"cmp-exact",
"xs=fresh",
"xs-exact"
);
let checkpoints = [1_000, 10_000, 100_000, 1_000_000];
for seed in [1, 2] {
random_updates(seed, 1_000, 3, &checkpoints);
random_updates(seed, 18_000, 3, &checkpoints);
random_updates(seed, 18_000, 9, &checkpoints);
}
absorbed_by_a_large_value();
large_weight();
cost();
}
+1
View File
@@ -50,6 +50,7 @@ const SECTIONS = [
["context/adr/005-exact-decimal-arithmetic.md", "adr/005-exact-decimal-arithmetic", "ADR 005 — Exact decimal arithmetic (proposed)"],
["context/adr/007-stable-identities.md", "adr/007-stable-identities", "ADR 007 — Stable identities"],
["context/adr/007-stable-identities-amendment-1.md", "adr/007-stable-identities-amendment-1", "ADR 007, amendment 1"],
["context/adr/009-incremental-evaluation.md", "adr/009-incremental-evaluation", "ADR 009 — Incremental evaluation (proposed)"],
],
},
{