common_ancestors + children_of: the dense-output batch pair (issue #156 phase 3)

[espg]

Refs #156 — phase 3 of many (not Closes). Phase 1 (mocs_to_orders + the mortie/moc.py domain split) merged as 42db670; this is the dense-output pair.

What this does

Adds the two dense-output batches from the issue #156 audit table. They are one phase because they share a property: their results are fixed-shape arrays, so neither needs an arrow offsets pair on the way out, and both are the cheapest high-value entries in the sweep.

common_ancestors(values, offsets) -> uint64[N]

Segmented reduce of common_ancestor / moc_min. Ragged in (the arrow list layout polygons_to_morton_mocs and mocs_to_orders already use), dense out — one word per group, because the reduction is many→one per group.

The consumer is a per-worker inner loop, not a one-off: zagg's t-digest merge runs for j in np.flatnonzero(~single): at src/zagg/stats/tdigest.py:198 — once per multi-member centroid (δ=512), per cell, per build and per fold, on every Lambda worker, at 65,536 cells/shard in the shipped ATL03 config. That module's own docstring at tdigest.py:177-179 already names it "the same O(n) Python-loop shape issue #279 removed".

One contract decision worth flagging: an empty group is an error, not an empty slot. mocs_to_orders lets an empty item densify to an empty slice, but a many→one reduction over no words has no answer — the scalar refuses empty input, so the batch refuses that group by index.

children_of(words, order, max_cells=None) -> uint64[N, 4**d]

Batch of generate_morton_children, whose wrapper coerces to a scalar parent (tools.py:1167). Every parent must sit at one order p <= order, so each yields exactly 4**d children for d = order - p and the result is a dense (n, 4**d) block.

That shared-order requirement is not new — it is what the consumers already assume, because the loop this replaces ends in np.stack, which raises on rows of unequal width. moczarr wrote the gap down in its own source: src/moczarr/dggs.py:310 is np.stack([generate_morton_children(int(w), level) for w in words]), with the comment at dggs.py:302-305 saying "there is still no vectorized many-parent children kernel". zagg calls the scalar the same way per sub-chunk on every worker (src/zagg/grids/healpix.py:199) and per shard in the shardmap reprojection (src/zagg/catalog/shardmap.py:708).

d == 0 returns the parents verbatim, matching the scalar exactly — that is what preserves a Kind::Point word, which a to_nested/from_nested round trip would re-pack as the order-29 area cell. Pinned by a test.

The result block is allocated fallibly (try_reserve_exact, not vec![0u64; total]), so an order whose (n, 4**d) block the allocator refuses is a catchable ValueError naming the byte count rather than a handle_alloc_error → abort() that kills the interpreter. An opt-in max_cells=None budget (espg's ruling on question (2)) refuses an over-budget n * 4**d result before anything is allocated, which is the only thing that fences an allocation an overcommitting OS accepts and then cannot back. See "Allocation posture" below.

House pattern

Per #153/#154/#160, and unchanged here:

• Plural naming marks many→many; children_of keeps espg's named form from the audit.
• Strict ragged contract on the one ragged input: offsets[0] == 0, offsets[-1] == len(values), exact coverage, both endpoints checked and each naming which one failed.
• Shared scalar params (order), not per-item arrays.
• Lowest-index fail-fast: serial validation pre-pass + per-item capture + index-order scan. No try_reduce short-circuiting (rejected in Batch polygon coverage: polygons_to_morton_mocs (issue #153) #154 for nondeterminism).
• rayon under py.allow_threads, chunked at CHUNK = 2048.

Rust module placement is the house pattern too, not a deviation. src_rust/src/decimal_morton/batch.rs sits exactly where the tree already puts batch submodules — coverage/batch.rs, moc/batch.rs, and PR #158's wkb/batch.rs are all <domain>/batch.rs. That is the axis espg ruled in the plan revision item 1 ("Rust side follows the same axis"), and both scalar kernels this phase wraps (common_ancestor, and the to_nested/from_nested pair the child generator is built from) live in decimal_morton.rs. An earlier revision of this description called it a deviation; it is not.

One place this phase genuinely does differ from its siblings:

• children_of validates entirely up front; common_ancestors does not. children_of must know the row width before it can allocate, so decodability / target order / shared-order all run in the pre-pass. common_ancestors validates only the layout up front and leaves per-group failures to the parallel pass. The reason is cost, not error ordering: the mixed-base-cell check is the reduction (it decodes every word in every group), so a pre-pass that ran it would run the whole kernel twice. Hoisting per se would not break the lowest-index rule — a single index-order pre-pass doing all the checks preserves it, which is exactly what validate_words and moc::batch's validator do. Only a partial hoist would invert (emptiness alone, as its own earlier pass, would report an empty group 7 ahead of a mixed-base-cell group 2). An earlier revision of this description gave the partial-hoist inversion as the reason for the whole design; that was wrong, and the module header now says cost.

The one place the lowest-index promise is qualified

For common_ancestors the rule holds within each pass, not across the two. Layout is a whole-batch serial pre-pass, so a bad offset at a high index outranks a domain failure at a low one:

v, o = ...            # group 0 spans two base cells (domain error at index 0)
o[8] = 99             # group 7's offset is out of bounds (layout error at index 7)
mortie.common_ancestors(v, o)
# ValueError: group 7: offset 99 exceeds value array length 20

That ordering is deliberate and now documented and pinned by a test (test_ancestors_layout_errors_outrank_domain_errors, plus a Rust twin): a layout failure says the offsets array itself is wrong, and the group indices a domain error would be reported by are read out of that same array, so "your offsets are broken" is the failure a caller has to fix first. moc::batch has no such split only because its domain check (a per-MOC cell budget, computable without reducing) is cheap enough to interleave into the index-order layout loop; this op's is not. Among the domain classes themselves the lowest index always wins — verified with offenders straddling the CHUNK seam in both directions, 200 repeats per configuration, one distinct answer every time.

On reuse: there were no pub(crate) helpers in coverage::batch on main to reuse (all its helpers are private, and PR #158 is unmerged). validate_ragged here is the offsets-only spelling of the same contract. It is deliberately not folded into a shared helper with coverage::batch / moc::batch, because each of those interleaves a domain check (the 3-vertex ring minimum; the per-MOC max_cells budget) inside the index-order loop, which is exactly what makes their errors lowest-index across both kinds of failure — a shared helper would have to take that as a callback, and this op has no such per-item pre-check. Raised as a question below.

Allocation posture — children_of must not abort the process

vec![0u64; total] goes through Rust's infallible allocator, whose failure path is handle_alloc_error → abort(). Measured on the pre-fix branch (macOS/arm64, 64 GiB, 10 cores), each row a fresh subprocess:

call	scalar np.stack([generate_morton_children(...)])	batch, before	batch, after
1 order-0 parent → order 29 (2.00 EiB)	MemoryError, exit 0	memory allocation of 2305843009213693952 bytes failed, SIGABRT exit 134	ValueError: 1 parents x 288230376151711744 children each needs 2305843009213693952 bytes; allocation failed, exit 0
12 order-0 parents → order 25 (96 PiB)	MemoryError, exit 0	SIGABRT exit 134	ValueError: ... needs 108086391056891904 bytes; allocation failed, exit 0
64 order-0 parents → order 29 (≥2⁶⁴ elements)	—	ValueError: 64 parents x 288230376151711744 children each overflows	unchanged
100k order-6 parents → order 11 (819 MB)	works	works	works

The existing checked_mul guard only trips at ≥2⁶⁴ elements — ≥64 parents at d=29 — so every reachable blowup sat below it and aborted. The fix is try_reserve_exact on the result block; no new constant, no policy.

What it does not fence. An allocation an overcommitting OS accepts and then cannot back still ends in a kill: 1M order-6 parents → order 14 (d=8, 488 GiB) is try_reserve-clean and dies at exit 137 (SIGKILL) while the block is being written. That is not a batch-specific gap — numpy accepts the same request (np.empty((1_000_000, 65_536), np.uint64) allocates 488 GiB on this box without raising, while (12, 1125899906842624) raises MemoryError), so the scalar loop dies at the same shape. After the fix the two agree everywhere measured: both raise catchably where the allocator refuses, both are killed where the OS overcommits. Refusing the second class needs a policy ceiling — that is what the opt-in max_cells added in 73b0681 is, per espg's ruling on question (2) below.

Cost of the fix. try_reserve_exact + resize loses vec![0; n]'s alloc_zeroed, so the block is memset once before the parallel pass overwrites it. Median of 5 runs, same box, before → after: d=4 14.4x → 12.9x, d=8 8.6x → 7.4x (≈10–15% off the largest results; no measurable change at d ≤ 2, where the result is small). The zero-fill is recoverable by writing into spare_capacity_mut() as MaybeUninit and a single set_len — one unsafe line in a module that has none today, which is why it is not in this diff. Flagged as question (6).

Memory posture — measured, not asserted

#162 exists because merged batches claimed "result + one chunk" while omitting the mandatory input copy (to_vec() is required for a numpy-borrowed slice to cross allow_threads). Both docstrings here state the real model, with numbers.

op	input	result	measured peak	vs input + result	vs result alone
common_ancestors, 5M groups × 3 order-9 words	152.6 MiB	38.1 MiB	191.9 MiB	1.01x	5.0x
children_of, 1M order-6 parents → order 9	7.6 MiB	488.3 MiB	497.1 MiB	1.00x	1.02x

Dense output makes the model short: input copy + result + one chunk of Results, where the chunk term is a flat 2048 × 32 B = 64 KiB regardless of batch size. children_of additionally allocates its block once at the exact final size, so there is no growth-realloc transient and no ragged assembly copy.

common_ancestors has a fourth term, and it scales with the largest group

An earlier revision of this description said "for common_ancestors the input copy is the peak, always". That is false for skewed groups. The scalar kernel builds its own scratch — let mut others: Vec<(u8, u64)> = Vec::with_capacity(rest.len()) (decimal_morton.rs:429), 16 B per non-first word — and holds it for the whole reduction, one per group in flight. So the term is min(threads, n_groups) * 16 B * max_group_size: it grows with the largest single group, which is the one dimension the quoted case (groups of 3) holds constant.

Measured on this branch, every input loaded from disk, same instrument and order across rows:

shape	input	documented model	measured Δpeak	ratio	excess	scratch the formula predicts
5M groups × 3 words, order 9	152.6 MiB	190.8 MiB	191.9 MiB	1.01x	~0	< 1 KiB
40 groups × 1M words, order 9	305.2 MiB	305.2 MiB	458.6 MiB	1.50x	153.4 MiB	152.6 MiB
1 group × 20M words, order 9	152.6 MiB	152.7 MiB	460.0 MiB	3.01x	307.3 MiB	305.2 MiB

The excess is the missing term to within 2 MiB in both skewed rows. The docstrings in decimal_morton/batch.rs and mortie/moc.py now state the four-term model and both regimes: size a worker off input + result for many small groups, and off the largest group when groups are large. For the t-digest consumer (groups of 2–3) the term is noise, which is why the original measurement was accurate and the generalisation was not.

A measurement caveat worth recording, because it bit this PR: sampled RSS under-counts the input copy when the input is constructed in-process, since to_vec is then served out of already-resident freed heap. The same 5M-group case reads as 77.3 MiB that way versus 191.9 MiB when the input is loaded from disk (independently reproduced to 0.1 MiB). The numbers above are the load-from-disk ones.

Benchmarks — honest read

benchmarks/measure_common_ancestors.py and benchmarks/measure_children_of.py, 10 cores, macOS/arm64, ≥100k items except where noted. The children_of scalar side is timed as the consumer actually writes it, np.stack included, since a caller cannot get a dense block out of the scalar without it.

These are medians of 5 runs, with the observed range. An earlier revision of this description published single-run figures (87.8x / 61.7x / 15.4x / 14.6x / 18.7x); re-measuring showed the batch side is a few hundredths of a second at small d, so one shot lands anywhere in a ±20% band and those numbers sat at the top of it. The medians below are lower and reproducible; the ranges are printed so a re-run that disagrees by 10% is not a surprise. Both scripts now carry the table in their module docstring.

case	result	median speedup	range over 5 runs
common_ancestors, 100k groups of 3, order 9	—	17.3x	15.0–18.3
common_ancestors, 500k groups of 2, order 9	—	19.6x	18.1–19.9
children_of, 100k parents, 6 → 7 (d=1)	3.1 MiB	74.2x	68.5–84.8
children_of, 100k parents, 6 → 8 (d=2)	12.2 MiB	50.3x	48.1–55.9
children_of, 100k parents, 6 → 9 (d=3)	48.8 MiB	28.2x	26.3–28.8
children_of, 100k parents, 6 → 10 (d=4)	195.3 MiB	12.9x	12.3–13.1
children_of, 2k parents, 6 → 14 (d=8, the shipped zagg shape)	1000.0 MiB	7.4x	7.1–7.5

These land high, unlike phase 1's mocs_to_orders (2.96–5.80x), and the reason is the same reason phase 1 landed low: what dominates. Phase 1 was bandwidth-bound on writing a flat result far larger than its input, so removing the boundary cost bought little. Here the per-item work is tiny (a shared-prefix walk over 2-3 words; a contiguous nested run) and the result is small or written once at exact size, so the Python boundary is the wall — which is what a batch removes.

The children_of decay does not flatten inside the published range, and d=8 is below its low end. 74x at d=1 falls to 12.9x at d=4 as the result turns it back into a bandwidth measurement, and keeps falling to 7.4x at d=8 — the 65,536 cells/shard shape zagg ships, and the only row a zagg or moczarr reader will actually look for. An earlier revision predicted d=8 would "sit nearer the low end" of a d≤4 table; it sits below it. The row is now measured, not predicted, and the benchmark docstring says so. (Independently measured at 8.71x pre-fix on real order-6 shard keys from the in-tree Antarctic drainage-basin fixture coarsened with clip2order(6, ...), 644 distinct keys — the random-parent corpus and the real-key corpus agree.)

Testing

• mortie/tests/test_dense_batch.py — 51 tests, all passing (42 in phase 3a/3b, 4 folding review, 5 folding the max_cells ruling).
• 21 Rust unit tests in decimal_morton::batch; cargo test --lib is 297 passed / 0 failed / 1 ignored.
• Full suite: 1119 passed, 12 skipped.

Coverage of the acceptance list:

• Per-item byte parity with the scalar twin, over randomized inputs and the in-tree Antarctic basin fixtures (Ant_Grounded_DrainageSystem_Polygons.txt) — for common_ancestors grouped by split_base_cells, for children_of at real order-9 shard keys refined the way zagg's two-step does it. Plus test_children_stacked_scalar_loop_is_reproduced_exactly, which is moczarr's np.stack([...]) expression verbatim.
• Edges: order-0 parents and base-cell groups; order-29 / max-order targets; single item; empty batch ((0,) and (0, 1)); a group of one word (returns that word); children_of with d = 0 (identity, and the point-word case); a word finer than order (refused, index named); mixed parent orders (refused, index named); undecodable words (refused, index named).
• Allocation: an unservable result is a ValueError a plain except ValueError catches (the PR Batch MOC densify: mocs_to_orders + the mortie/moc.py extraction (issue #156) #160 lesson, where a PanicException escaped even except Exception) — and, because pytest.raises cannot observe an abort(), a subprocess test asserts the interpreter survives it (exit 0). Plus a 819 MB result that must still go through, so the fallible allocation is not a de facto budget.
• Error ordering: the layout-before-domain asymmetry is pinned in both Python and Rust, so it cannot drift silently.
• Determinism: test_ancestors_lowest_index_holds_across_chunk_boundaries and its children_of twin are parametrized over offenders at indices 0, 1, 7, CHUNK-1, CHUNK, CHUNK+1, 2*CHUNK+5, 3*CHUNK-1 and repeat each 5x; the ..._wins_over_later_offenders pair plants five simultaneous offenders and heals them lowest-first, asserting the named index each time (3x per step). Plus 10-run byte-equality checks on both ops.
• GIL genuinely released: a pure-Python counter thread ticks during a large call, same instrument as test_threading.test_gil_released_during_rust_compute.

Gates

cargo fmt --check clean · cargo clippy --lib --benches — no new warnings (the remaining 7 are the pre-existing coverage/tests.rs, geo2mort.rs, prefix_trie.rs set) · cargo test --lib 297 passed / 1 ignored · pytest 1119 passed / 12 skipped · flake8 mortie --select=E9,F63,F7,F82 clean · flake8 --max-line-length=88 clean on every touched file (the three tools.py hits are pre-existing, at lines 958/975/986 — the on_antimeridian one is #151 item 4) · numpydoc lint clean on moc.py, tools.py, __init__.py · doctests pass on both touched modules.

Phases

• Phase 3a — common_ancestors: Rust kernel + binding + mortie/moc.py wrapper, beside its scalar twin.
• Phase 3b — children_of: Rust kernel + binding + mortie/tools.py wrapper, beside its scalar twin.
• Tests, benchmarks, USAGE.md, docs/api/moc.md, docs/api/tools.md, __init__ exports.
• Adversarial self-review folded (4 findings: fallible allocation, the scratch term, the hoisting justification, the bench table).
• espg's rulings on questions (2) and (3) folded: opt-in max_cells, and (0, 1) kept with the reasoning recorded.

Landed as one commit because the two share a Rust module and a test file; splitting them would have been file-splitting, not phasing.

Questions for review

1. tools.py is now 1,451 lines, against this repo's ~1,000-line aim (CLAUDE.md §4). It was already over at 1,354 before this PR. The instruction for this phase was explicit: put children_of beside its scalar twin and note that Split tools.py and geometry.py into domain modules mirroring the Rust tree #159 will move the hierarchy/order functions to orders.py later, and do not start the Split tools.py and geometry.py into domain modules mirroring the Rust tree #159 split here (it is blocked). I followed that, but the rule says stop and raise when a file would cross, so raising it.

   The remedy has a deadline that is earlier than "after the sweep." children_of genuinely belongs beside generate_morton_children, so 1,451 is the right place to land this phase — but every remaining phase lands in the same file too: phase 4's words_to_decimals / hive_paths, phase 5's mort2polygons / mort2bboxes. That puts tools.py near ~1,700 by the end of the sweep. So the ask is to unblock Split tools.py and geometry.py into domain modules mirroring the Rust tree #159 before phase 4, not after phase 5, so the split absorbs the batch additions instead of chasing them. mortie/moc.py is fine at 644; src_rust/src/decimal_morton/batch.rs is 702.

   Second, smaller point: mortie's CLAUDE.md has no 1,200-line raise — that is zagg's standing ruling (its issues #351/#358), not this repo's, and it should not be inherited silently. If you want the same allowance here it belongs in mortie's own CLAUDE.md as its own decision.

2. Should children_of carry a max_cells budget? RULED — option (b), implemented in 73b0681. espg ruled the opt-in form: "it costs adopters nothing and gives a Lambda worker a way to fail catchably instead of being killed."

   children_of(words, order, max_cells=None). Unset — the default — behaviour is byte-identical to before: same SHA-256 over children_of output across six shapes (d=0, 3, 5, 8 spans, order-0 and order-29 edges, plus the empty batch) on the pre- and post-change builds, and the try_reserve path and its documented overcommit residual are untouched. Set, the result's n * 4**d cell count is checked before anything is allocated and refused with a ValueError in moc_to_order's wording, so the two read alike:

   ValueError: children_of would generate 65536000000 cells (1000000 parents x 65536
   children each) at order 14, exceeding max_cells=134217728. Pass a larger max_cells,
   or max_cells=None to proceed (risking OOM), or refine to a coarser order.
   

   That is the case that motivated the ruling and it is the acceptance test: zagg's shipped d=8 shape at fleet scale (1M order-6 parents → order 14, 488 GiB) is the one try_reserve cannot refuse, because the OS accepts the reservation and kills the process mid-write. With a worker's own ceiling set it comes back catchable.

   Why the default is the opposite of moc_to_order's, stated in the docstring so the API does not read as inconsistent: moc_to_order defaults its budget on because a densify explodes from a tiny input (Σ 4**(order - depth), pre-emptive flat-cover size guard #80) and the caller cannot cheaply predict the output; children_of's output is exactly n * 4**d, computable from the arguments before the call, so a default guard would refuse calls the caller already knows are fine. Same parameter name, opposite default, for a stated reason. The None polarity inverts with it — here None means no budget, not disable a default — and that is written down too.

   Follows Consolidated follow-ups from the June/July sweeps (#73, #88, #97) #108's ruled semantics where they apply: ValueError (a caller-controlled argument condition, not a failed allocation), caller-settable, negative rejected, and clamped at u64::MAX the way mocs_to_orders clamps.

   One ordering decision, made deliberately and pinned: the budget is compared in u128 and checked ahead of the checked_mul element-count guard, so an explicit budget answers even for a request too large to represent as a usize. 64 order-0 parents at order 29 gives ... children each overflows with no budget and ... exceeding max_cells=1048576 with one; likewise the budget outranks the allocation failed refusal. Rationale: an explicit argument condition is more actionable than an internal representability diagnostic the caller cannot act on (Consolidated follow-ups from the June/July sweeps (#73, #88, #97) #108's framing). Pinned in children_max_cells_outranks_the_overflow_guard (Rust) and test_children_max_cells_outranks_the_overflow_guard (Python).

   The residual stands and is why (b) beat (a): an unset budget still cannot close the overcommit class, and the scalar loop has the identical gap (np.empty((1_000_000, 65_536), np.uint64) allocates 488 GiB on this box without raising). A caller-set ceiling is the only thing that closes it, which is exactly what (b) provides without changing anyone's behaviour by default.

3. children_of on an empty batch returns shape (0, 1). RULED — keep (0, 1); no parent_order= parameter. No behaviour change; the reasoning is now recorded in the docstring (722ec03) so it is not re-opened. With no words there is no source order to derive 4**(order - parent_order) from, so a derived width would require callers to pass parent_order — repeating, in the empty case alone, what the data carries in every other case. Consumers that need a typed empty already special-case it: moczarr's dggs.py builds (0, 4**(level - order)) on its own empty branch precisely because it knows its source order at that point. So (0, 1) is the honest shape for "no rows, width unknown", and the width belongs in the consumer-side special case.

4. The three ragged validators are now structurally parallel and not shared (coverage::batch, moc::batch, decimal_morton::batch). I deliberately did not refactor: the other two interleave a domain check inside the index-order loop, so unifying them needs a callback, and touching coverage/batch.rs right now would collide with PR Rust WKB reader: backend-free geometry ingest, plus from_wkbs (issue #157) #158. Worth a follow-up issue after Rust WKB reader: backend-free geometry ingest, plus from_wkbs (issue #157) #158 merges, or leave the three as they are?

5. Should common_ancestors grow a split_base_cells-style escape? Today a group spanning more than one base cell refuses the whole call by index, matching the scalar. That is the right default, but the real Antarctic fixtures do straddle base cells — the parity test has to run split_base_cells first to build groups the reduction accepts. If a consumer wants "reduce each group, and partition the ones that straddle", that is a different op and I would rather it be filed than smuggled in here.

6. Should children_of recover the zero-fill the fallible allocation costs? try_reserve_exact + resize memsets the block before the parallel pass overwrites every word of it — measured at ~10–15% on the large-d shapes (d=4 14.4x → 12.9x, d=8 8.6x → 7.4x; no measurable change at d ≤ 2). Writing into spare_capacity_mut() as [MaybeUninit<u64>] and calling set_len after the last chunk succeeds removes it entirely and is provably fully-initialised (words.len() rows × width = total, and an early return drops the vec at len 0). It costs one unsafe line in a module that has none today — unsafe in this crate is currently confined to arrow_ffi.rs and the C ABI export. I left it out rather than expand the scope of a correctness fold; say the word and it is a two-line follow-up.

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 95.31%. Comparing base (42db670) to head (722ec03).

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #164      +/-   ##
==========================================
+ Coverage   95.27%   95.31%   +0.04%     
==========================================
  Files          12       12              
  Lines        1587     1601      +14     
==========================================
+ Hits         1512     1526      +14     
  Misses         75       75              

Flag	Coverage Δ
unittests	95.31% <100.00%> (+0.04%)	⬆️

Flags with carried forward coverage won't be shown. Click here to find out more.

Files with missing lines	Coverage Δ
mortie/__init__.py	89.28% <ø> (ø)
mortie/moc.py	100.00% <100.00%> (ø)
mortie/tools.py	98.10% <100.00%> (+0.05%)	⬆️

---

Continue to review full report in Codecov by Harness.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 42db670...722ec03. Read the comment docs.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[espg]

🤖 from Claude (review)

[high] children_of aborts the whole process on a large-d allocation, where the scalar loop it replaces raises a catchable MemoryError — and the checked_mul guard never fires below that point.

Measured on this branch (rebuilt with maturin develop --release, macOS/arm64, 64 GiB, 10 cores), each row a fresh subprocess:

call	scalar np.stack([generate_morton_children(...)])	batch children_of
1 order-0 parent → order 29 (2.00 EiB)	MemoryError: Unable to allocate 2.00 EiB ..., exit 0	memory allocation of 2305843009213693952 bytes failed → SIGABRT, exit 134
12 order-0 parents → order 25 (96 PiB)	MemoryError: Unable to allocate 8.00 PiB ..., exit 0	SIGABRT, exit 134
1M order-6 parents → order 14 (d=8, the shipped zagg shape, 488 GiB)	—	SIGKILL, exit 137

vec![0u64; total] on this line goes through Rust's infallible allocator, so an allocation failure is handle_alloc_error → abort(). There is no Python exception, no traceback, and no except that can see it — the interpreter is gone. On a Lambda worker that is a silent worker death rather than an error a caller can catch and retry at a coarser order.

The checked_mul above is decorative in practice: it only trips at n * 4**d >= 2**64 elements, i.e. ≥64 parents at d=29. I confirmed that case does return cleanly —

ValueError: 64 parents x 288230376151711744 children each overflows

— but every reachable blowup (2 EiB, 96 PiB, 488 GiB) sits below the guard and aborts instead.

This also inverts the premise of PR question (2), option (a): "adding one would refuse calls the Python loop it replaces completes". The Python loop does not complete them — np.arange(span, dtype=np.uint64) in generate_morton_children raises a catchable MemoryError at exactly the sizes where the batch aborts. A guard here would therefore move the batch closer to the scalar's observable behaviour, not away from it. That is the same reasoning #108 used to make moc_to_order's ceiling a catchable ValueError, and the same reasoning PR #160 used to fence a PanicException the batch could not otherwise surface.

The budget question (a/b/c) is espg's to rule, but the abort is independent of it: even with no max_cells, replacing vec![0u64; total] with a try_reserve-based allocation (or a total * 8 ceiling check returning Err(String)) turns this into a ValueError and needs no new constant.

Repro:

import numpy as np, mortie
from mortie import _rustie
w = np.asarray(_rustie.rust_nested2mort(
    np.arange(12, dtype=np.uint64), np.zeros(12, np.uint8)), np.uint64)
mortie.children_of(w, 25)     # -> SIGABRT, exit 134, no Python exception

[espg]

🤖 from Claude

Fixed in 48e230d. vec![0u64; total] is now try_reserve_exact + resize, so the failure comes back as Err(String) and the binding turns it into ValueError. Verified on this branch, each row a fresh subprocess, rebuilt with maturin develop --release:

call	before	after
1 order-0 parent → order 29 (2.00 EiB)	memory allocation of 2305843009213693952 bytes failed, exit 134	ValueError: 1 parents x 288230376151711744 children each needs 2305843009213693952 bytes; allocation failed, exit 0
12 order-0 parents → order 25 (96 PiB)	exit 134	ValueError: ... needs 108086391056891904 bytes; allocation failed, exit 0
64 order-0 parents → order 29	ValueError: ... children each overflows	unchanged
100k order-6 parents → order 11 (819 MB)	works	works — the guard is not a de facto budget

Three tests, all in 48e230d:

• test_children_oversized_result_raises_instead_of_aborting — asserts a plain except ValueError sees it (the PR Batch MOC densify: mocs_to_orders + the mortie/moc.py extraction (issue #156) #160 pattern you cite), and that the checked_mul overflow message is still separately reachable.
• test_children_oversized_result_leaves_the_process_alive — the one that is actually the point: pytest.raises cannot observe an abort(), so this runs the repro in a subprocess and asserts exit 0 plus ALIVE on stdout.
• children_oversized_result_errors_instead_of_aborting on the Rust side, pinning the byte count in the message.

Your point about PR question (2)'s premise is taken, and the question is rewritten. It now leads with the correction — the loop does not complete those calls, np.arange raises catchably — and I moved my own lean to (b) max_cells=None on the strength of it.

One thing your table's third row does not fix, and I want it on the record rather than implied away. The 1M order-6 → order 14 case (d=8, 488 GiB) is still exit 137. try_reserve returns Ok there — macOS accepts the reservation and the kernel kills the process later, while the block is being written. That is not a batch-vs-scalar gap: I measured numpy doing the same thing, which is the honest comparison your — left open:

np.empty((1_000_000, 65_536), np.uint64)      -> ALLOCATED ok, nbytes 524288000000
np.empty((12, 1125899906842624), np.uint64)   -> MemoryError: Unable to allocate 96.0 PiB ...

So after this fix the batch and the scalar agree at every size measured: both raise catchably where the allocator refuses, both are killed where the OS overcommits and then cannot back it. Closing the second class needs a policy ceiling, which is exactly what question (2) now asks you to rule on — and I've said in the PR body that the abort fix was independent of and prior to that ruling, as you framed it.

Cost, disclosed rather than buried. try_reserve_exact + resize gives up vec![0; n]'s alloc_zeroed, so the block is memset once before the parallel pass overwrites every word of it. Median of 5 runs, same box, before → after: d=4 14.4x → 12.9x, d=8 8.6x → 7.4x; no measurable change at d ≤ 2. It is fully recoverable by writing into spare_capacity_mut() as [MaybeUninit<u64>] with a single set_len after the last chunk succeeds, but that is one unsafe line in a module that has none — unsafe in this crate is confined to arrow_ffi.rs and the C ABI export today — so I did not expand a correctness fold into it. Raised as new PR question (6) for espg.

[espg]

🤖 from Claude (review)

[medium] The common_ancestors memory model omits the scalar kernel's own per-group scratch, so "the input copy is therefore the peak, always" is false for skewed groups — measured 3.00x the stated model.

common_ancestor (decimal_morton.rs:429) builds let mut others: Vec<(u8, u64)> = Vec::with_capacity(rest.len()) — 16 B per non-first word in the group — and holds it for the whole reduction, once per group in flight. The model documented here (input copy + result + one 64 KiB chunk) has no term for it.

Independently measured on this branch, every input loaded from disk so the input copy is not masked by freed heap. The caveat four lines below is real and I reproduced it to 0.1 MiB — the 5M x 3 case reads 77.2 MiB built in-process versus 192.0 MiB loaded from disk, against the 77.3 / 191.9 stated here — so the load-from-disk numbers below are the comparable ones:

shape	input	documented model	measured Δpeak	ratio	excess
5M groups x 3 words, order 9 (the case quoted here)	152.6 MiB	190.8 MiB	191.5 MiB	1.00x	0.7 MiB
40 groups x 1M words, order 9	305.2 MiB	305.2 MiB	458.4 MiB	1.50x	153.2 MiB
1 group x 20M words, order 9	152.6 MiB	152.7 MiB	458.6 MiB	3.00x	305.9 MiB

The excess is exactly the missing term, to within a MiB:

• single 20M-word group: 16 B x 19,999,999 = 305.2 MiB (measured excess 305.9)
• 40 groups of 1M: 10 cores x 16 B x 1,000,000 = 152.6 MiB (measured excess 153.2)

So the bound is min(threads, n_groups) * 16 B * max_group_size, not a flat 64 KiB — it scales with the largest group, which is the one dimension the documented case (groups of 3) holds constant.

The quoted measurement is accurate and for the t-digest consumer this term is noise. The problem is the generalisation: "always" is what a caller sizing a worker will read, and #162 exists because a memory claim was carried past the shape it was measured on. Suggest scoping the sentence to small groups and stating the scratch term alongside the chunk term; the same sentence is mirrored in mortie/moc.py's docstring ("the input copy is the peak").

[espg]

🤖 from Claude

Fixed in 6555542. "The input copy is therefore the peak, always" is gone; both docstrings (decimal_morton/batch.rs module header and mortie/moc.py) now state a four-term model with the scratch named and bounded:

peak is input copy + result + one 64 KiB chunk + the reduction's own scratch, where the scratch is common_ancestor's others buffer — 16 B per non-first word in the group it is reducing, held for that whole reduction, one per group in flight — so the term is min(threads, groups) * 16 B * max_group_size. It scales with the largest single group, not with n.

That last sentence is the part I think matters most and it is your finding, not mine: the term is invisible in the consumer shape and dominant in a skewed one precisely because it keys off the one dimension the quoted case holds constant.

Re-measured independently on this branch, every input loaded from disk, same instrument across rows:

shape	input	documented model	measured Δpeak	ratio	excess	formula predicts
5M groups × 3 words, order 9	152.6 MiB	190.8 MiB	191.9 MiB	1.01x	~0	< 1 KiB
40 groups × 1M words, order 9	305.2 MiB	305.2 MiB	458.6 MiB	1.50x	153.4 MiB	152.6 MiB
1 group × 20M words, order 9	152.6 MiB	152.7 MiB	460.0 MiB	3.01x	307.3 MiB	305.2 MiB

Within 2 MiB of your 458.4 / 458.6 on both skewed rows, and the excess matches the formula to within 2 MiB in both — so the term is identified, not just observed. The PR body carries the same table.

Both docstrings now close on the two regimes rather than one rule: size a worker off input + result for many small groups, and off the largest group when groups are large.

I did not change the kernel to eliminate the buffer (it is a real second pass over the group's decoded words, so removing it is a rewrite of common_ancestor, not a batch-layer change) — that would be a scalar-side optimisation and is out of this diff's scope. Say if you want it filed.

[espg]

🤖 from Claude (review)

[low] The lowest-index guarantee holds across domain classes (verified) but not across the layout/domain split — and the reason given here for the asymmetry only holds for a partial hoist.

Three things, all about this paragraph.

(1) The domain-class behaviour is correct and I could not break it. With offenders straddling the CHUNK seam and 200 repeats per configuration, the lowest index always wins across differing failure classes:

empty group at	mixed-base-cells at	reported
4000	3	group 3: inputs span multiple base cells ...
3	4000	group 3: empty input has no common ancestor
2049	2047	group 2047: ... multiple base cells
2047	2049	group 2047: empty input ...
0	5000	group 0: empty input ...
5000	0	group 0: ... multiple base cells

200/200 runs of the first two rows each produced a single distinct answer. So the shipped behaviour is what this paragraph promises.

(2) Layout errors are in a different pass, so they outrank domain errors regardless of index — the same inversion this paragraph says it is avoiding, on the other axis:

v, o = ...            # group 0 spans two base cells (domain error at index 0)
o[8] = 99             # group 7's offset is out of bounds (layout error at index 7)
mortie.common_ancestors(v, o)
# ValueError: group 7: offset 99 exceeds value array length 20

moc::batch::validate_batch does not have this gap — its domain check (the per-MOC max_cells budget) is interleaved into the same index-order loop as the offsets checks, which is exactly what line 53 credits it with. Here validate_ragged has no domain check to interleave, so the guarantee really is weaker; but the module header ("naming the lowest-index offending item"), the common_ancestors doc comment, and mortie/moc.py's Raises section all list layout and domain failures under one undifferentiated lowest-index promise.

(3) The justification as written is not the load-bearing one. Hoisting all the per-group checks into a single index-order pre-pass would preserve lowest-index, not break it — that is precisely what validate_words does with its three checks 150 lines below, and what moc::batch does. Only a partial hoist (emptiness alone, as its own earlier pass) produces the "empty group 7 ahead of mixed-base-cell group 2" inversion described here. The real reason not to hoist is cost: the mixed-base-cell check is the reduction, so a pre-pass would run the whole kernel twice. That is a stronger argument than the one written, and worth correcting because later phases will copy this paragraph.

[espg]

🤖 from Claude

All three points taken; fixed in 826e793.

(3) — the justification. You are right that the written reason is false as stated, and right about why: a full index-order pre-pass preserves lowest-index, and validate_words 150 lines below is the counterexample in this same file. The module header now says the load-bearing reason:

The reason is cost, not ordering: the mixed-base-cell check is the reduction — it decodes every word in every group — so a pre-pass that ran it would run the whole kernel twice. Hoisting per se would not break the lowest-index rule: a single index-order pre-pass doing all the checks preserves it, which is exactly what validate_words and crate::moc::batch's validator do. Only a partial hoist would invert — emptiness alone, as its own earlier pass, would report an empty group 7 ahead of a mixed-base-cell group 2.

Your "later phases will copy this paragraph" is the reason I rewrote rather than deleted it. The partial-hoist inversion is kept as the parenthetical it always should have been.

(2) — the layout/domain inversion. Reproduced exactly:

domain only  -> group 0: inputs span multiple base cells and have no common ancestor; ...
layout at 7  -> group 7: offset 99 exceeds value array length 20

Now documented as deliberate rather than left implicit, with the reason stated: a layout failure says the offsets array itself is wrong, and the group indices a domain error would be reported by are read out of that same array — so "your offsets are broken" is the failure a caller has to fix first. And I took your moc::batch contrast into the text: it has no such split only because its domain check (a per-MOC cell budget, computable without reducing) is cheap enough to interleave; this op's is not, which is the same cost argument as (3).

Documented and pinned, so it cannot drift silently now that it is a stated contract:

• ancestors_layout_errors_outrank_domain_errors (Rust) — group 0 domain error on its own, then a group-2 layout error added, which outranks it.
• test_ancestors_layout_errors_outrank_domain_errors (Python) — your exact repro shape, o[8] = 99 against a mixed-base-cell group 0.

The undifferentiated promise you flagged in three places is split in all three: the module header, common_ancestors's # Errors section, and mortie/moc.py's Raises now all say the lowest-index rule holds within a pass — layout first for the whole batch, then domain in index order — instead of implying it holds across both.

(1) — nothing to do, and I did not touch it. Your 200-repeat sweep is the evidence that the shipped behaviour is correct; I only corrected the comment in ancestors_errors_name_lowest_index_group, which claimed the mixed-base-cell-beats-empty case was "the whole reason emptiness is not hoisted into the pre-pass". Both are domain failures found in the same pass, so the lower index wins regardless of which check is cheaper — which is your point, restated where a future reader will hit it.

[espg]

🤖 from Claude (review)

[low] The d=8 shape zagg actually ships measures 8.7x — below the table's 15.4x low end, not "nearer the low end" of it.

Measured on this branch with real order-6 shard keys (the in-tree Antarctic drainage-basin fixture coarsened with clip2order(6, ...), 644 distinct keys), refined to order 14 — the 65,536 cells/shard shape the PR body cites:

n=2000 parents, 6 -> 14 (d=8, 65,536 children/parent), result 1000 MiB
scalar np.stack loop : 2.79 s
batch children_of    : 0.32 s
speedup              : 8.71x

For calibration, my reruns of the published rows on the same box (10 cores) also come in under the table — common_ancestors 12.8x (100k x 3) and 13.4x (500k x 2); children_of 63.8x (d=1) and 13.5x (d=4) — so part of the gap is machine load. But the decay is monotone and continues past d=4 rather than flattening inside the published range, which is the substantive point: d=8 is the only row a zagg or moczarr reader will look for, and it is the one row that is predicted rather than measured. Either add the measured d=8 row to the table, or reword to "below the low end".

[espg]

🤖 from Claude

Taken on both counts; fixed in 41e4508. The optimistic set is gone — I did not keep it and explain the gap, because re-measuring showed there is nothing to explain except run-to-run variance, and the published numbers sat at the top of it.

Republished as medians of 5 runs with the observed range, same box (10 cores, macOS/arm64), post-fold build. Both benchmark scripts now carry this table in their module docstring, so it is reproducible from the file rather than only from the PR:

case	result	median	range over 5 runs	was published
common_ancestors, 100k groups of 3, order 9	—	17.3x	15.0–18.3	14.6x
common_ancestors, 500k groups of 2, order 9	—	19.6x	18.1–19.9	18.7x
children_of, 100k, 6 → 7 (d=1)	3.1 MiB	74.2x	68.5–84.8	87.8x
children_of, 100k, 6 → 8 (d=2)	12.2 MiB	50.3x	48.1–55.9	61.7x
children_of, 100k, 6 → 9 (d=3)	48.8 MiB	28.2x	26.3–28.8	— (new)
children_of, 100k, 6 → 10 (d=4)	195.3 MiB	12.9x	12.3–13.1	15.4x
children_of, 2k, 6 → 14 (d=8)	1000.0 MiB	7.4x	7.1–7.5	predicted only

The range column is the point: at small d the batch side is a few hundredths of a second, so one shot lands anywhere in a ±20% band — d=1 alone spanned 68.5–84.8 across five runs. That is where 87.8x came from, and why your reruns and mine both came in under it. The docstrings now say to quote the median and the range.

Your substantive point — the decay does not flatten inside the published range — is now measured, not reworded. I added the d=3 row so the curve is visible (74 → 50 → 28 → 12.9 → 7.4) and the d=8 row is measured rather than predicted. The "expect d=8 to sit nearer the low end" sentence is replaced by a statement that it sits below it, and measure_children_of.py's docstring now warns in its own words: a figure quoted at d ≤ 4 does not bound d = 8; measure the d you intend to run.

Your 8.71x on the real order-6 shard keys and my 7.4x on the random-parent corpus are the same measurement either side of the fold: 8.71 vs my 8.64 pre-fix median on the same shape, and the drop to 7.4 is the memset the try_reserve fix in 48e230d costs (disclosed on that thread and in the PR body). So the two corpora agree, which is worth having established — the random-parent bench is a fair stand-in for the shipped shape.

[codspeed-hq]

Merging this PR will degrade performance by 14.16%

❌ 1 regressed benchmark
✅ 70 untouched benchmarks
⏩ 1 skipped benchmark¹

Warning

Please fix the performance issues or acknowledge them on CodSpeed.

Performance Changes

	Benchmark	BASE	HEAD	Efficiency
❌	coverage_square[4]	147.2 µs	171.5 µs	-14.16%

Tip

Investigate this regression by commenting @codspeedbot fix this regression on this PR, or directly use the CodSpeed MCP with your agent.

---

_{Comparing claude/156-common-ancestors (722ec03) with main (42db670)}

Footnotes

1. 1 benchmark was skipped, so the baseline result was used instead. If it was deleted from the codebase, click here and archive it to remove it from the performance reports. ↩

[espg]

🤖 from Claude

espg ruling (in-session, 2026-08-08): question (2) resolves to option (b) — add max_cells=None to children_of as an opt-in budget with no default, so callers can ask for a catchable refusal without changing anyone's behaviour. Rationale on the record: it costs adopters nothing and gives a Lambda worker a way to fail catchably instead of being SIGKILLed by an overcommitting OS.

Implementing now; waiting cleared while that lands, and will go back on with the PR when it does.

Two things to carry into the implementation, recorded so they are not rediscovered:

• The parameter's default is deliberately opposite to moc_to_order's, and the docstring must say why. moc_to_order(max_cells=_FLAT_COVER_WARN_THRESHOLD) defaults on because a densify can explode from a small input — the caller cannot cheaply predict the output size. children_of's output size is exactly n * 4**d cells, computable from the arguments before the call, so a default guard would refuse calls the caller already knows are fine. Same name, opposite default, for a stated reason.
• This does not close the overcommit class by itself — it is the only thing that can close it, which is why (b) was ruled, but the residual recorded in the abort fold stands: try_reserve returns Ok under macOS overcommit and the kernel kills during the write. With max_cells unset the behaviour is unchanged, including that residual, and the scalar has the same gap (np.empty((1_000_000, 65_536), np.uint64) also allocates 488 GiB without raising).

[espg]

🤖 from Claude

espg ruling (in-session, 2026-08-08): question (3) resolves to the stated lean — keep (0, 1), no parent_order= parameter.

Reasoning on the record: deriving the width from the data beats making callers repeat what the data already says, and the empty case is already special-cased consumer-side — moczarr's dggs.py builds (0, 4**(level - order)) on its own empty branch precisely because it knows its source order, which the batch does not when there are no words to read it from.

No code change. A one-line note is being added to the docstring recording why the width is 1 rather than derived, so the question is not re-litigated by the next reader.

That closes (3). Still open on this PR: (1) the tools.py overage at 1,451 (with the reviewer's and the fold's independent recommendation to unblock #159 before phase 4, since phases 4 and 5 also land there and would push it toward ~1,700), and (6) whether to take the unsafe spare_capacity_mut() + set_len recovery for the memset cost the abort fix introduced (d=4 14.4x → 12.9x, d=8 8.6x → 7.4x). (2) is ruled and being implemented as opt-in max_cells=None.

[espg]

🤖 from Claude

espg ruling (in-session, 2026-08-08): question (6) resolves NO — do not take the unsafe recovery. Recording the reasoning so it is not re-opened.

The memset is the price of making the allocation fallible: vec![0u64; total] uses alloc_zeroed (free zero pages from the OS on large blocks) but aborts on failure; try_reserve_exact + resize is fallible but must write the zeros. There is no safe way to have both — fallible allocation yields uninitialised capacity, and safely reaching len == n means writing n values. Fallibility was the right trade.

Why not buy the 14% back with set_len:

1. It changes the crate's unsafe story from a defensible one-liner to a judgement call. unsafe here is currently confined to arrow_ffi.rs and the C ABI export — places FFI forces it. "Unsafe only where the boundary requires it" is auditable; "unsafe wherever it is faster" is not, and the boundary does not move back.
2. The proof obligation is not one-time, and this PR already demonstrates why. set_len after a parallel fill obliges every future change to re-establish full initialisation. The max_cells opt-in ruled on this same PR adds exactly such an early-return path — had the unsafe landed first, that ruling would have silently invalidated the proof, with nothing in the type system to catch it.
3. The trade is small where it lands. 7.4x vs 8.6x against the scalar loop are both large wins, and at d=8 the dominant cost is materialising 65,536 children per parent at all. If that shape's throughput ever matters, the lever is not set_len — it is not requesting 488 GiB in one call, which max_cells now lets a caller express.

Question (6) closed. The measured cost stays documented in the benchmark table so the trade is visible rather than lost.

[espg]

🤖 from Claude

Both rulings are landed and pushed — 73b0681 (opt-in max_cells) and 722ec03 (the empty-width rationale). CI: 16 SUCCESS, 0 failing; only the CodSpeed "Rust benchmarks" job is outstanding, which runs ~27 min on every branch including main and is a perf-tracking job, not a correctness gate. Reapplying waiting.

This PR is now genuinely ready. Every phase, review finding, and ruling is folded. Two things from the implementation worth a reviewer's eye:

The checked_mul ordering needed a real decision. max_cells now fires first, compared in u128, so it stays reachable at sizes where the element count cannot be represented as a usize at all. The rationale, recorded in the code: a caller-controlled argument condition (#108's framing) is more actionable than an internal representability diagnostic the caller cannot act on. Consequence: the checked_mul overflow guard is now only reachable with max_cells=None. Both sides are pinned.

The acceptance case is the one that motivated the ruling. 1M order-6 parents → order 14 — the shipped d=8 shape, 488 GiB, and the exact case try_reserve cannot refuse under macOS overcommit — now returns a catchable ValueError naming 65536000000 cells when a worker sets a ceiling. The same shape under the budget still returns (1000, 65536), so it is a ceiling rather than a refusal of the shape.

Unset is byte-identical, verified rather than asserted: the same SHA-256 over children_of output across six shapes (d=0 through d=5, plus both order edges) and the empty batch, computed on builds of 41e4508 and 722ec03.

Standing for espg: (1) tools.py at 1,451, with the reviewer's and fold's independent recommendation to unblock #159 before phase 4 rather than after phase 5; (4) whether to file the shared-ragged-validator follow-up (my read: leave the three as they are — the duplication is a contract check, and a callback abstraction would have one degenerate user); (5) and (6) are closed.