feat(svar2)!: plan conversion concurrency under core and memory budgets (+ scale-bench harness) - #141
feat(svar2)!: plan conversion concurrency under core and memory budgets (+ scale-bench harness)#141d-laub wants to merge 53 commits into
Conversation
Adds a reproducible sweep over the sub-contig sharded VCF reader budget (reader_workers x per-shard HTSlib threads x OVERSHARD_FACTOR) to evaluate PR #140, plus BENCH-ONLY env hooks in orchestrator.rs so one build covers the whole space (GENORAY_READER_WORKERS / GENORAY_SHARD_HTSLIB / GENORAY_OVERSHARD). Unset, they leave the planner untouched. Every configuration is gated by a whole-store byte-identity oracle; all 60+ swept configurations matched, confirming sharding does not perturb output at any worker count. Findings are written up in scripts/bench_sharded_vcf/README.md. Co-Authored-By: Claude Opus 5 <[email protected]>
Defines a harness that predicts biobank-scale (500k sample) conversion behaviour from affordable measurements, and discriminates three hypotheses for what a reader-budget autotuner should key on: a static cap, a formula in cohort size, or a byte-bounded global reader pool. Co-Authored-By: Claude Opus 5 <[email protected]>
Eight tasks across five waves; waves 1, 2 and 5 hold pairs of tasks on disjoint files that run concurrently. Co-Authored-By: Claude Opus 5 <[email protected]>
Implements the pure-math core of the scale-bench harness: fits V/cost/RAM scaling laws from measured records and turns the design spec's falsifiable criteria into an H1/H2/H3/none hypothesis verdict. Fixed one arithmetic inconsistency in the source brief's worked test example: the planted exponent 1.6 in test_v_law_reports_low_r2_on_nonlinear_data produces r2=0.983 against the brief's own reference fit_v_law, which fails the brief's own `r2 < 0.98` assertion (verified independently outside the test suite). Bumped the planted exponent to 2.0 (r2~=0.962, clear margin) without touching V_LAW_MIN_R2 or any other threshold/production code.
check_oracle pooled digests across the whole results file, so a legitimate multi-corpus sweep (different corpora hold different variant data, so their store digests differ by design) hit a guaranteed false-positive RuntimeError the moment a second corpus's results landed. Group the check by corpus (via a point_id -> corpus map rebuilt from the plan) and run it after every point instead of only at the end, so a real same-corpus divergence surfaces immediately -- append_ndjson already fsyncs each record, so nothing is lost by failing fast.
…ADME wording
The Slurm driver never entered the repo root, so `python -m scripts.bench_svar2.*`
would fail with ModuleNotFoundError unless submitted from the worktree itself;
add cd "$WT" so the script is location-independent. The contig sweep plan's
concurrency loop iterated the literal tuple (1, min(c, 4)), which collapses to
(1, 1) at c=1 and runs that configuration twice (sweep.py snapshots pending
points once, so both duplicates execute, not just the second on restart);
iterate sorted({1, min(c, 4)}) instead so c=1 yields one point, with a comment
on why a high-split counterfactual is inherently absent at a single contig.
Also reword the README to stop conflating seed-deterministic corpora with
runtime-measured (non-reproducible) results.
Adds the ~2-minute regression tier: a tiny 200-sample/20k-variant corpus
run at reader_workers in {1,3,7}, checked against committed wall_s/maxrss_mb
baselines with a one-sided 25% tolerance band (faster/leaner never fails,
a missing baseline or failed run is always reported).
Co-Authored-By: Claude Opus 5 <[email protected]>
…icated allocation Fix round 1 for the fast regression tier. Finding 1 (critical): baselines were keyed by SweepPoint.point_id, which hashes the corpus path under $CLAUDE_JOB_DIR -- a per-session ephemeral directory. Any other session, same box and cores, reported "no baseline recorded" for every point and exited 1. The committed file is now keyed by reader_workers and re-keyed onto the session's point_ids. Finding 2 (important): wall_s is no longer a hard gate. A 25% band over a contention-inflated baseline absorbs a real 30-60% slowdown and still passes. Measurement shows a dedicated allocation does not fix this either -- an immediate re-check on an idle node swung -24%/-9%/+13% -- so the comment no longer suggests promoting it there. maxrss_mb remains the hard gate. Finding 3 (important): the corpus shrank 20_000 -> 2_000 variants and the runtime claim is now measured, not asserted: 65 s cold including corpus generation, 58 s warm, on a dedicated 8-CPU allocation. Also adds regression_record.sbatch so baselines are recorded reproducibly, restores the warm-up rep (recording without it baked a 4x first-touch transient into the baselines), and guards against a stale cached corpus -- the gate is one-sided, so baselines from a larger corpus checked against a smaller one passed vacuously.
…widths Fix round 2 for the fast regression tier, from the scoped re-review. The theme is the same in every case: the hard gate is one-sided, so anything that lets a wrong baseline through does not fail loudly -- it silently stops gating. N1 (important): --record wrote baselines even when a point failed. run_point returns ok=False records that still carry a maxrss_mb from the crashed child, and recording runs unattended under sbatch, so a crash-inflated number would be committed with nobody watching and could never be exceeded again. Recording now exits 1 without touching the baseline file if any point failed. N2 (important): threads (allocation width) reaches the conversion as `-@ N` and sizes a rayon pool, so it moves the maxrss_mb that gates -- but the re-keying in round 1 removed the accidental protection point_id hashing used to provide. The file now records the width it was taken at and a check at a different width is refused with an actionable message instead of silently compared. N3/N4/N5 (minor): corpus_is_current now also invalidates on a GENERATOR_VERSION bump, requires the corpus the manifest describes to still exist (the default workdir is scratch that gets cleaned), and compares against generate()'s floored variant count so a future non-divisible multi-contig CORPUS does not regenerate on every invocation forever. N6 (minor): the baseline pairing reads reader_workers off each point instead of zipping against a parallel WORKERS sequence -- a second source of truth whose mispairings would be silent for the same one-sided reason. N7 (minor): the HARD_METRICS comment called the recording node "idle"; the job log shows loadavg ~22-31 with isolated cores. It now cites both recordings, which disagree by -27% on identical code -- a stronger argument for not gating wall_s than the one it replaces. Baselines re-recorded by job 13332630 on a dedicated 8-CPU allocation; measured runtime in the docstring updated to that job (50 s cold, 50 s warm).
The whole-branch review found the harness ran cleanly but could not
discriminate the three hypotheses it exists to decide between. Engineering
quality was not the problem; three gaps between the code and the design spec
were.
C1 (critical): the reorder-backlog gauge sampled AFTER inserting the just-read
chunk, so pending_hw was never below 1 in any sharded run. decide()'s H3 gate
is `pending_hw / workers >= 0.5`, and every scale point runs at workers=1, so
H3 fired on the first row of every sweep. Planting synthetic H1 and H2 data
returned H3 both times. Since H3 is also the spec's stated prior, the overnight
job would have produced an authoritative-looking confirmation of the prior with
no evidence behind it. The gauge now samples BEFORE the insert, via a
PendingBacklog wrapper whose single insert_observing method makes
insert-without-observe unrepresentable; a chunk released the instant it arrives
contributes 0. Remove sites deliberately do not observe -- fetch_max ignores
lower values, so it would be a provable no-op.
C2 (critical): the V-linearity ladder had no corpora and no plan points. The
spec calls V-linearity the only reason 500,000 samples is reachable and gates
on R^2 >= 0.98, but the S-ladder holds the cell budget fixed (V = 1.4e9 / S),
so phase1_s is roughly constant while V varies 2000x -- fitting it yields
R^2 ~ 0 and the harness correctly refuses to project, producing no headline
result. Adds a vlinear plan (V in {25k,50k,100k,200k} at fixed S=250 and fixed
chunk_size, only the corpus varying) plus its four corpora.
I1 (important): model.py had no entry point and was imported by nothing but its
test, so the sweep ended at three NDJSON files. Adds
`python -m scripts.bench_svar2.model` to read results, fit the laws, print the
verdict and extrapolation, and gate the hold-out on the spec's 25% threshold.
Also: extrapolate now projects the pending term fit_ram_law was fitted with
(dropping it understates peak RSS by the reorder factor, which is exactly what
H3 is about); extrapolation_factor was variants/n_points, a variant count over
a point count; predicted_wall_s scales the per-variant slope by the fitted
read-cost exponent instead of applying a small-cohort slope at 500,000 samples;
degenerate _linfit fits report an unbounded rather than zero-width CI; the
scale plan gained production-default chunk_size=25_000 points, without which
the RAM law is fitted at essentially one chunk_bytes and the OOM deliverable
is unobtainable; probe.py attributes oom_at_rss_mb only on an actual allocation
failure rather than on any nonzero exit, and pins MALLOC_ARENA_MAX so the
RLIMIT_AS ceiling tracks RSS; the regression tier gates on the
worker-attributable RSS delta rather than absolute RSS, which the fixed
footprint otherwise swamps.
chunk_size in the regression tier stays at the production 25_000: shrinking it
to 128 was proposed on the theory that a large chunk over a small corpus never
fills, but read_next_chunk allocates the full grid up front and truncates only
after EOF, so shrinking it halved the signal the delta gate reads (measured
21.4 -> 7.8 MB at w=7). Baselines re-recorded on a dedicated allocation.
Verified: 97 bench tests pass (was 66); ruff clean; cargo check and cargo test
--no-default-features --features conversion both pass (Slurm job 13332809).
…isoned
Second whole-branch review round. Both criticals are cases where a plausible
mechanism was wrong and only measurement caught it.
The H3 gate still could not discriminate. The previous round removed the +1
gauge floor that made H3 fire at w=1, but not the architectural reason it fires
at w>1: ReorderBuffer::push releases a chunk on arrival only when its ordinal is
the head, so the w-1 units ahead of the head keep everything they produce
buffered until the head unit finishes. Roughly (w-1)*chunks_per_unit chunks are
resident even with perfectly balanced readers and zero skew -- a real 12-unit,
w=3, overshard=4 probe log in this repo sustains pending=5. Every planned sweep
row (scale w in {2,3,5,7,11}, contig w in {12,6,3}) would have tripped
`pending_hw >= w/2` and returned H3 before H1 or H2 were evaluated. The planted
H1/H2 tests passed only because they hand-planted pending=0, which the harness
cannot produce at w>1. H3(a) now gates on the backlog's BYTE share of measured
peak RSS (kappa*pending*chunk_bytes / maxrss_mb >= 0.25), which tests the
hypothesis's actual claim -- bytes, not worker count, set peak RSS -- rather
than a count that grows structurally. It stays reachable: a 1 GB chunk fires it
at w=1, where "structural" is not an argument. knee_spread and beta_diff_ci95
now compute above both H3 returns, so an H3 verdict no longer ships evidence
with the H1/H2 terms missing.
The production-default points poisoned the RAM law. BitGrid3::zeros is
vec![0u64; n] -> alloc_zeroed -> calloc, so untouched pages never become
resident: a 3 GB zeroed allocation adds 0 MB to ru_maxrss, measured on this
node. The earlier justification for those points ("a small corpus still pays the
large allocation") is true of address space and false of maxrss_mb, which is the
metric recorded. Since the sweep fixes the cell budget, large-S corpora are tiny
in V, so a 25,000-variant chunk never fills -- nominal chunk_bytes of 1,562 and
3,125 MB against ~350 MB actually touched, while every other row in the sweep
sits under ~120. Two enormous-leverage points with a local slope of ~0.3 dragged
OLS kappa about 10x below its true ~3, which would have projected ~1.3 GB
instead of ~9.8 GB at S=500,000 and reported the hardcoded chunk_size=25_000 as
SAFE at biobank scale -- the reverse of the spec's headline arithmetic. The
regressor is now bounded by the corpus: a chunk cannot hold more variants than
exist. extrapolate is deliberately unchanged, since at V=1e9 the chunk does fill.
Also: the cohort correction was built on cpu_shard_pct, a UTILIZATION percentage
that pegs near 100% at every S on a reader-bound pipeline, so beta_read ~ 0 and
the correction collapsed to a no-op -- predicted_wall_s at 500,000 samples was
the 250-sample wall time. It now fits an absolute per-variant cost against S, and
extrapolation is skipped with a reason rather than silently assuming cohort size
does not matter. The regression tier's delta band gained an 8 MB floor: two
dedicated-allocation recordings of identical code disagree by up to 5.98 MB, so
a pure 25% band would have failed unchanged code by 60%. rss_ceiling_mb (and
therefore MALLOC_ARENA_MAX=1, measured at 73% slower in a multithreaded
conversion regime) is now set only on the points that exist to probe the
ceiling, leaving every law-fitting point on the production allocator.
Retracts the disproved allocated-up-front mechanism from the regression tier's
comment; the chunk_size=25_000 choice stands on its measurement (worker deltas
4.2/21.4 MB vs 3.6/7.8 MB at 128) and is no longer explained by it.
Verified: 105 bench tests pass (was 97); ruff clean; records.py untouched.
4b69b28 to
789a5fb
Compare
Third and final review round: the six Minor findings plus the deferred Verdict round-trip item. Two of these are data-collection bugs and had to land before the cluster sweep started, not after. probe.py, cpu column alignment. `cpu_shard_pct` and `cpu_exec_pct` were appended to independent lists, but `model._median_costs` ZIPS them --- so a single `n/a` in one column (the single-reader fallback path emits one) did not merely drop a sample, it shifted every later `cpu_exec` sample against a different tick's `cpu_shard` value. The corruption is silent because both tuples still look well-formed, and it lands directly on `c_read / c_exec`, which is the knee. Both fields come from the same `pipeline sampler` line and nothing downstream reads either series alone, so a tick missing one is now dropped from both. probe.py, OOM attribution. `_is_oom_failure`'s own docstring lists "a preemption signal" among the cases it exists to exclude, while its SIGKILL branch returned True unconditionally --- readmitting exactly that. Slurm ends a preempted or time-limited job with SIGKILL. Worse, under this harness's own configuration the bare signal is near-certainly NOT an OOM: `rss_ceiling_mb` is installed as RLIMIT_AS (60 GB) while the sweep's cgroup allows 120 GB, so genuine exhaustion trips RLIMIT_AS first and dies via SIGABRT carrying the allocator message the regex already catches; the cgroup OOM killer cannot fire until twice the ceiling. SIGKILL now requires the process to have reached at least half the ceiling, which keeps the branch meaningful if that relationship is ever reconfigured without minting "OOMs at scale" data out of cluster scheduling. model.py, hold-out scoring. The V-law is fitted `phase1_s ~ a + b*V`, so the projection is a phase-1 time, but the hold-out check scored it against `wall_s` --- which also carries the reader-independent rayon merge tail and process startup and is therefore always larger. Every hold-out error was inflated one-sidedly, feeding a 25% gate whose documented meaning is "this invalidates the model, not just this point". The key is renamed `predicted_wall_s` -> `predicted_phase1_s`: the name invited exactly one wrong comparison and got it, and there is no correct way to compare that quantity against a wall time. When `phase1_s` is 0 the time half of the gate is now skipped and said to be skipped rather than falling back to `wall_s`. The V-ladder and hold-out corpora are both single-contig, so both sides are one uncontended span. model.py, H1 support floor. `decide` had no minimum-points guard, so a single surviving cohort size gave `knee_spread == 0` --- spread over one point is 0 by definition --- and returned a confident "a static cap suffices, no autotuner needed". That is the shape a partly-failed sweep takes: rows that yield no usable cpu ticks get dropped, `knees` shrinks, and the verdict grew MORE confident as evidence disappeared. H1 now requires 3 distinct cohort sizes (the scale plan supplies 7), `knee_points` rides along on every verdict so a spread of 0 can be told apart from a spread of 0 with support, and the under-powered case says the flatness is unevaluable instead of silently reading as a finding. regression.py, uniform-shift blind spot. The delta gate subtracts the `reader_workers=1` point off every other point, so it cannot see a regression that moves all worker counts equally --- a bigger shared buffer, an extra copy of the dense grid, a leak scaling with variants rather than readers. With absolute `maxrss_mb` reporting-only since Finding I8, such a change passed the hard gate at any magnitude. Adds an absolute backstop at a deliberately loose 50% band: a doubling detector, not a drift detector, so it does not reintroduce the I8 sensitivity problem. The new test confirms a uniform +400 MB shift leaves the delta gate silent and is caught only by the backstop. monitor.rs: `PendingGauge`'s doc comment now actually documents the before-insert counting convention that `shard_exec.rs` cross-references it for, including why an off-by-one there reads downstream as a real second peak-RSS term. sweep_scale.sbatch rebuilds the extension before measuring. `pixi run` does not, so the sweep could otherwise spend 24 hours measuring whichever .so happened to be installed --- including one predating the PendingGauge counting change, where `pending_highwater` means something different. Also adds a characterization test pinning that `Verdict.evidence` does NOT survive a JSON round-trip (tuples degrade to lists, int keys to strings) because `_tuple_fields` only inspects top-level field annotations and `evidence` is `dict[str, Any]`. Verdict is print-only today, which is the only reason that is harmless; records.py is frozen, so this documents the limit for whoever adds persistence. Gitignores the sbatch logs that land in the repo root. Verified on a dedicated allocation: 117 bench tests pass (was 105); ruff check and format clean; records.py untouched. Preflight job 13335579 confirms the fixes on real trace data --- pending_highwater=2 at w=3 with 6.66 MB of backlog bytes, and cpu_shard/cpu_exec parsed index-aligned. That w=3 run is also a live instance of the structural floor the previous round identified: the old count-based H3 gate would have fired on it, and the byte-share gate correctly does not.
The 24h scale sweep deadlocked 3h15m in with a 0-byte hold-out corpus: every
process at 0% CPU, bgzip with nothing to compress, and pool workers respawning
at staggered times. Three defects compounded.
Block sizing ignored cohort size. `BLOCK_VARIANTS = 2_000` bounds a task by
VARIANTS alone, so a block holds `2000 * n_samples` cells -- 8e6 at the contig
axis's 4,000 samples, but 2e8 at the hold-out's 100,000, with 14 such blocks
formatting concurrently under `--procs 16`. The large-S scale corpora survived
only by accident: they have 2-3 blocks total, so at most 3 were ever live.
`.astype(str)` on an int64 array returns dtype `<U21` -- numpy sizes the result
for the widest possible int64, 84 bytes per element. The four dp/gq/ad/(ad*2)
copies alone are ~67 GB over a 2e8-cell block, which is why only the hold-out
died: it is the sole corpus using `--format-fields`. Every value here is
bounded under 1000, so `.astype("U3")` is lossless and the emitted text is
byte-identical.
`mp.Pool` turned the resulting OOM kills into a silent hang. It quietly
repopulates dead workers while `imap` waits forever on results that will never
arrive -- no error, no output, no progress. That is what converted a crash
into three lost hours. `ProcessPoolExecutor` raises `BrokenProcessPool`, so
`set -e` now kills the job in seconds with a diagnosis.
`plan_blocks` bounds peak pool memory, with the two knobs deliberately
asymmetric. `procs` is capped freely: concurrency alone cannot affect output
because `.map` preserves order and `_format_block` seeds per block. Block size
is reduced only on the FORMAT path and only as a function of corpus shape,
because it sets both the position striping and the per-block seed. The GT-only
path keeps cutting at `BLOCK_VARIANTS` so the 11 corpora already generated for
this sweep, and the regression baselines recorded against one of them, stay
byte-reproducible; changing that is a GENERATOR_VERSION bump, not a bug fix.
Two follow-on bugs, both caught only by end-to-end verification:
`_block_positions` striped at the BLOCK_VARIANTS constant rather than the
actual block size. With smaller blocks each stripe sat mostly empty and the
last block began at `n_blocks * BLOCK_VARIANTS * stride`, ~20x past the
declared contig length. Positions still came out sorted, so nothing looked
wrong until tabix rejected the finished file -- after 21 minutes and 12 GB had
been written. The stripe now uses the block size actually cut, restoring the
`per_contig * stride` bound the surrounding comment already claimed.
The first cut of the memory fix derived the per-block budget as a `procs`
share of the pool budget, which made `--procs` change block partitioning and
therefore output bytes. That silently broke the determinism `_format_block`'s
per-block seeding exists to guarantee, and every unit-level check still passed
-- only a byte comparison of the same corpus at two worker counts caught it.
The budget is now a per-block constant.
Also: the sbatch cached corpora on `[ -f "$OUT" ]`, which a 0-byte .vcf.gz
satisfies, so a resubmit would have SKIPPED regeneration and failed later. It
now caches on the manifest, which `generate` writes last, after tabix and the
record-count check.
Verified: 17 corpus tests pass (was 12). New coverage asserts blocks actually
shrink before testing them (the first version of that test was vacuous and
passed), that block size is invariant across procs in {1,2,4,8,16,48}, that a
sub-2000-block corpus is tabix-indexable with POS inside the contig and
strictly increasing, that the same request at procs=1 and procs=16 is
byte-identical, that every GT-only shape on the sample axis still cuts at
2000, and that an unfittable block raises instead of hanging.
The V-linearity ladder was falsified by data the sweep had already collected, at the same cohort size, so no cohort scaling was even involved in the refutation. Per-variant phase-1 cost fell monotonically across every rung -- 3.04e-4, 2.94e-4, 1.81e-4, 1.61e-4 s/variant -- and reached 1.12e-5 at the scale sweep's own S=250 point (V=5,600,000), still falling. The ladder never left the fixed-cost-dominated regime, so fitting a LINE through it and stretching that to V=1e9 predicted 740.5s at V=5,600,000 where 62.9s was measured: 11.8x high. Downstream, that one slope was the difference between a usable model and a useless one. The driver reported R^2=0.9738 against its own 0.98 gate, a 2820% hold-out error against a 25% gate, and a predicted phase-1 time of 4.0e9 seconds -- 128 years -- for the target regime. The mechanism was per-CHUNK cost, not per-variant work. `VLINEAR_CHUNK_SIZE` was derived from `min(VLINEAR_VARIANTS)`, which pinned the whole ladder at chunk_size=781 while `extrapolate` targets PROD_CHUNK_SIZE=25_000, 32x larger. At S=250 that is decisive: the V=200_000 rung is 257 chunks/32.1s, and the V=5,600,000 point at chunk_size=25_000 is 224 chunks/62.9s -- near-identical chunk counts, 28x the data, 2x the time. The fitted "per-variant slope" was a per-chunk cost divided by 781, then applied at a 32x larger chunk size. Both halves of the confound are fixed. The chunk size is pinned to PROD_CHUNK_SIZE, stated directly rather than derived, because a V-law is only usable at the chunk size it was fitted under. The V range moves into the many-chunk regime: 800_000 (the smallest V clearing the >=MIN_CHUNKS floor at that chunk size) through 5,600,000, every rung >=32 chunks. A raise-on-violation guard replaces the arithmetic coincidence that `_chunk_size_for(800_000)` happens to return 25_000, so a future V change cannot silently reintroduce the confound. Two things improve for free. The top rung cuts the V-law's extrapolation stretch to V=1e9 from ~5000x to ~179x -- the step the driver itself flagged as least-supported -- and it matches the s250 corpus shape exactly, so it cross-checks against an independently measured point. Cost is unchanged in the axis that matters: S stays pinned at 250, so the top rung is 1.4e9 cells, the same as the largest scale corpus, and peak RSS on this ladder was ~450MB. No stale-data hazard: `load_sweep` keeps only records resolvable to a plan point, so the four records measured at chunk_size=781 are dropped by point_id and named in `excluded` rather than pooled with the new ones. Verified against rebuilt plans -- the driver reports 4 exclusions, skips the V-law, and skips the extrapolation instead of emitting a number it cannot defend. Also: sweep_scale.sbatch derives `--threads` from `os.sched_getaffinity(0)` instead of hardcoding 48. `-@` is the conversion's TOTAL budget and reader workers are carved out of it, so lowering --cpus-per-task while leaving --threads at 48 would have oversubscribed the allocation and measured the reader knee against a pool the job does not own. --mem drops 120G -> 64G on measured evidence (peak MaxRSS 4.58GB across 20 scale points; the OOM probes bound ADDRESS SPACE via RLIMIT_AS=60GB, which trips before 64GB resident). --time 24h -> 72h since the partition allows 14 days and the s500000 ladder alone projects to ~18h. Plan order is now cheapest-and-most-load-bearing first. Running `scale` first nearly cost the study its answer: vlinear is ~10 minutes and fits the V-law, the V-law is the only reason 500,000 samples is reachable, and a timeout inside scale would have left the model unable to produce ANY verdict after ~24h of measurement. Safe to reorder because `pending_points` resumes.
The sharded VCF reader built one `VcfRecordSource` per shard via `::new`, which resolves every requested sample name through `HeaderView::sample_id` -- and that is `samples().iter().position(..)`, a linear scan. Resolving S samples against an S-sample header is therefore O(S^2), and the sharded path paid it once per shard while shard count is `reader_workers * OVERSHARD_FACTOR`. Total cost: O(reader_workers * S^2). At biobank scale that was the entire conversion. Regressing measured CPU seconds on shard count across a 500,000-sample corpus gives ~1,150s per shard with an intercept of -973 -- statistically zero -- meaning essentially none of phase-1 CPU was record decoding. It also inverted the whole point of the sharded reader: adding readers made conversion monotonically SLOWER, 4,682s at w=1 to 7,810s at w=11 (1.67x), while burning a full extra core per worker (cpu_s 4,591 -> 60,239). The regime flip with cohort size is the tell. The same ladder on small cohorts behaves as designed -- at S=250 readers give a 3.1x speedup with FLAT cpu_s (63s at both w=1 and w=11), and at S=16,000 a 2.7x speedup -- because S^2 is negligible there. Only once S^2 dominates does adding readers become pure cost. Two independent fixes, both needed: 1. `resolve_against_header` builds a name->index HashMap once, making one resolution O(k + S) instead of O(k * S). `new` and `resolve_sample_indices` now share it, so they cannot drift. Duplicate names keep the FIRST index (`or_insert`), matching `position`'s semantics, and the not-found message is byte-identical, so the existing error-parity test between the two still holds. 2. The orchestrator hoists resolution out of the per-shard closure and passes the resolved vector to `with_sample_indices`, dropping the `* reader_workers` factor along with one file open and header parse per shard. This is exactly the reuse `with_sample_indices`'s docstring was written for; the sharded path simply never adopted it. The hoist is only sound if the mapping is shard-invariant, so `sharded_permuted_sample_subset_matches_unsharded_bytes` pins it: a REORDERED, non-contiguous subset (header S0..S3, requesting S2, S0, S3) must produce byte-identical output sharded and unsharded. The pre-existing sharded test passes the full cohort in header order, where the mapping is the identity and a resolution bug is invisible. Verified sensitive by mutation: reversing the hoisted vector fails the new test, and it passes once reverted. 413 tests pass (412 + the new guard), clippy clean under -D warnings, and `cargo check --no-default-features` still builds the query-only core. No public API change, so no `skills/genoray-api/SKILL.md` update: `resolve_against_header` is private, and both public constructors keep their signatures and their observable behavior. Nothing reachable from `import genoray` moves.
The harness found a real O(w·S²) defect in the sharded reader (f17dcac)The completed scale sweep (27/27 points, 1d21h) showed reader workers making conversion monotonically slower at biobank scale — the opposite of PR #140's intent:
Cost per shard is flat across an 11× range, and regressing Root cause: the sharded path built one This is why the ladder flips sign with cohort size: at S=250 readers give a 3.1× speedup with flat Fix: (1) Guarded by |
Peak RSS carries a large term proportional to cohort size that has nothing to do with chunk bytes, and the law had nowhere to put it. Holding chunk_bytes pinned at 10.9MB and varying only the cohort, measured peak RSS runs 789MB at S=4,000 to 5,061MB at S=500,000 -- a 6.4x spread against a regressor that did not move. Above a ~780MB floor the growth is linear in S (2,161 -> 4,281MB for a 2x cohort), i.e. ~11KB/sample of per-sample accumulation buffer that exists whether or not a chunk is in flight. Modelling that as a bare constant intercept left the fit at R^2=0.057 across a real 39-point sweep -- the regressor explained essentially nothing -- and pushed kappa to 2.94 against per-worker slopes implying far less. Adding the term takes the SAME data to R^2=0.913 with kappa 1.11. That is not cosmetic, because the H3 verdict is computed from kappa. H3(a) asks what share of peak RSS the reorder backlog accounts for, and with kappa nearly 3x too high the backlog looked 3x more important than it is: the share was 45% against a 25% gate. Corrected, it is 17% -- BELOW the gate. So the criterion this harness had reported H3 on for three consecutive runs was an artifact of the mis-specified law. H3 still stands, but now on the independent multi-contig evidence (same total readers split differently across contigs differ by 199% wall time), which is a different and better supported claim. The biobank RSS projection drops 6.756e6 -> 2.552e6 MB. Both regressors are solved jointly rather than by fitting one and regressing the other on the residual: real sweeps correlate them (bigger cohorts get smaller chunk_size), so sequential fitting hands the shared variance to whichever goes first. A single-cohort sweep makes the `samples` column a constant multiple of the intercept column, leaving the two unidentifiable -- least squares then returns a minimum-norm split that is pure arithmetic artifact, and `extrapolate` multiplies that coefficient by 500,000. `fit_ram_law` now detects the degenerate case and drops the regressor so `base_mb` owns the constant, which is all a one-cohort sweep can support. Rows become a `RamRow` NamedTuple. The payload is four ints and a float, and swapping `chunk_bytes` for `samples` positionally still fits, just wrongly, with nothing to catch it; every consumer now reads by name. Tests: `test_ram_law_recovers_planted_cohort_term` plants both regressors CORRELATED the way a real sweep correlates them, so it also pins the joint solve. `test_ram_law_without_cohort_term_cannot_fit_cohort_scaling` plants RSS varying only with cohort at fixed chunk bytes -- the exact shape that exposed the defect -- so dropping the regressor later fails loudly instead of returning a confident, wrong kappa. The nominal-chunk regression test keeps its intent but restates its assertion: with a cohort regressor present nominal chunk_bytes is near-collinear with `samples`, so the estimate blows up rather than collapsing to ~0.23. It is not 3, which is the claim. 129 bench tests pass, ruff clean. Known-remaining and NOT addressed here: the hold-out still fails its 25% gate (phase1 63%, rss 65%) because it is the only corpus with FORMAT fields -- every law-fitting corpus is F=0, so the model is fitted on F=0 and validated on F=3. That is a harness design flaw, tracked separately.
The hold-out gate reported "MODEL FAILURE: this invalidates the model" at 63% phase-1 error, and the cause was not the extrapolation. The hold-out is the ONLY corpus in the harness carrying FORMAT fields (DP, GQ, AD); every law-fitting corpus is F=0. No cost law carries an F term at all -- `extrapolate` threads `format_fields` into chunk_bytes for the RSS side only -- so the model was fitted on F=0 and validated on F=3, and the error it was condemned for is dominated by an unmodelled FORMAT-decode cost. Decoding four per-sample values instead of one is the right order for the 2.74x observed. That is a harness design flaw, not a model failure, and the difference matters: the gate's whole job is to say whether the S,V extrapolation to biobank scale can be trusted, and it was answering a question about F. Two changes, and deliberately not a third. `decide`'s caller now computes which F values the laws were actually fitted from and reports a hold-out outside that set as OUT-OF-DOMAIN rather than as a model failure. It still prints the errors -- silence would be worse than a wrong verdict -- but does not condemn the extrapolation on evidence that cannot bear it, and it names the fix (fit an F law, or hold out in-domain). `HOLDOUT_F0` adds a second hold-out at the same S and V with no FORMAT fields, so the gate has something IN the fitted domain to check and its error is attributable to S,V extrapolation and nothing else. The F=3 hold-out is KEPT and still runs. Deleting it would make the gate pass by removing the dimension the model cannot handle, which is how a benchmark becomes dishonest. Real biobank VCFs carry FORMAT fields, so the standing record of how far off the model is on data it does not cover is worth more than a green gate. Fitting an actual F law needs corpora at several F across several S -- a new sweep dimension, tracked separately, not smuggled in here. Tests: `test_main_end_to_end`'s fixture had the same confound the real sweep did -- it drove the gate through an F=3 hold-out while fitting on F=0, so it asserted MODEL FAILURE while exercising the path where that must not fire. It is now F=0 and additionally asserts OUT-OF-DOMAIN is absent, pinning the in-domain half. `test_holdout_outside_the_fitted_format_field_domain_is_not_a _model_failure` pins the other half with a deliberately absurd record, so the gate WOULD fire if the point were in-domain -- otherwise the assertion is vacuous. 130 bench tests pass, ruff clean.
… pin it The cohort law was fitted from the scale ladder, which holds S*V = CELLS_BUDGET at every rung. That makes its regressand `log(phase1/V)` identically `log(phase1) + log(S) - log(cells)`, so the fitted slope is `1 + dlog(phase1)/dlog(S)` -- and a constant-cells ladder is BUILT so every rung does the same total work, leaving phase1 flat (36-44s across a 2000x cohort range) and the slope collapsed to 1. It reported beta=1.0020 with a 95% CI of [0.9689, 1.0352]: tight enough to read as a solid measurement, from a design that could not have returned anything else. `extrapolate` then raised that number to a 2000x cohort ratio. Fit beta from two V-ladders instead. Within a ladder S is fixed and V varies, so each fitted slope IS a per-variant cost at that cohort size, and beta is the log-ratio of the two. Measured that way beta = 0.9860, CI [0.9808, 0.9912] -- OUTSIDE the forced fit's own CI, so that fit was not merely unidentified but wrong and falsely confident. - `fit_cohort_beta_from_ladders` propagates the two slope standard errors through the log-ratio rather than reporting the zero residual two points always leave, which would collapse the CI to a point and read as certainty. - `cohort_beta_is_design_forced` detects the degenerate ladder and names the remedy, so a run without the second ladder says so instead of reporting a confident-looking number. - `VLINEAR2_SAMPLES` is deliberately NOT the hold-out's cohort size: a ladder at the hold-out's S would make the hold-out an interpolation inside the fitted data and the gate would go quiet for the wrong reason. The two ladders (S=250 and S=250,000) bracket the hold-out's S=100,000, with tests pinning that invariant. Also record the machine. `point_id` hashes every field of `SweepPoint` but deliberately not the node, so a resumed sweep skips work already paid for -- which also means two records can share a point_id and disagree wildly. Measured: point 2ac9bbbfbe0dc691 ran 151.9s on carter-cn-03 and 73.2s on carter-cn-04, a 2.08x spread against a 25% gate, while same-node controls reproduced within 1.9% and a repeated control was identical. That one cross-node record was the entire reported "40% MODEL FAILURE". `model.py` now reports the numbers and declines to charge a node gap to the model. Together these take the in-domain F=0 hold-out from a 40% MODEL FAILURE to 13% phase-1 / 3% RSS, both under the gate, and move the S=500,000 V=1e9 projection from 4539h to 4020h. beta is not a true constant -- per-cell cost is mildly U-shaped in S, so the S=250/S=100,000 pair gives 0.9592 and the S=250/S=250,000 pair 0.9860, with non-overlapping CIs. A single power law is an approximation, fitted across the range that brackets the hold-out. Recorded in the README. Co-Authored-By: Claude Opus 5 <[email protected]>
The scale-bench harness returned H3 -- worker count is the wrong invariant, in-flight bytes is the right one -- and this is the design that acts on it. `plan_thread_budget` charges every contig MIN_THREADS_PER_CHROM = 6 cores: four pipeline threads plus two HTSlib decode threads. Both halves are wrong on the sharded VCF path. The HTSlib pool does not exist there (SHARDED_VCF_HTSLIB_THREADS_PER_READER is 0 -- shard readers decompress inline), and of the four pipeline threads only the executor is CPU-bound. Measured, the 22-contig cc=1 w=12 row runs 16 threads on 2.02 cores. So the planner reserves ~6 cores for something that consumes 0.4-1.0, and plan_thread_budget(48, 22) returns cc=7 where the machine could run all 22. What that costs, holding total reader threads fixed at 12 and varying only the split: 22 contigs go 35.11s at cc=1,w=12 versus 11.74s at cc=4,w=3 -- 2.99x for identical work (both sum to phase1 = 33.0s of per-contig spans). The best of those still leaves 43 of 48 cores idle. More readers cannot substitute: `run_compute_engine` is a serial recv() loop, so a large contig runs its executor at 99-101% with the dense channel full at its cap, and many small contigs leave every executor at 36-43% with the channel empty. Neither regime improves with w; both improve with cc. The design plans cc under explicit core and memory constraints (the fitted RAM law becomes the byte budget), orders contigs by cost estimates read from index metadata already on disk, and derives w from a probed t_read/t_exec ratio rather than a hardcoded knee -- node speed here varies 2.08x, which a fitted knee cannot see. Two things called out rather than left implicit: max_mem=None means a detected budget, not unbounded, which is a default behavior change; and the budget must come from the cgroup limit, not /proc/meminfo, since under Slurm the latter reports the node rather than the job -- on exactly the allocations where the planner matters most. Sub-contig executors are deliberately out of scope. At cc = n_contigs the makespan floor is the largest contig's serial-executor span (~1.8x optimal for a whole-genome VCF) and removing it needs output-assembly changes that belong in their own spec. Co-Authored-By: Claude Opus 5 <[email protected]>
Seven tasks against the approved spec. Tasks 1-4 (cgroup memory detection, the constrained planner, index-derived contig costs, the rate probe) are mutually independent and meant to be dispatched in parallel; Task 5 wires them, and Tasks 6-7 (the digest gate, docs plus a bench axis) follow. Two things the plan resolves rather than leaves to the implementer. The spec's three-tier cost fallback becomes two: the middle tier, the linear index's compressed byte extent, needs more unsafe CSI walking than the rest of the module for what is only a sort key. And tier 1 (hts_idx_get_stat) is documented for BAM with no guarantee CSI-over-VCF populates it, so the task is written to let its test decide -- with an explicit instruction to delete tier 1 rather than weaken the assertion if the counts come back equal. The digest-invariance test is the gate the change rides on: cc, w, and contig order all move, and each can perturb chunk ordinals, per-chunk ledgers, or long-allele bank offsets. Its fixture uses eight contigs with deliberately unequal record counts, since with equal contigs longest-first ordering is a no-op and the test would prove nothing. Co-Authored-By: Claude Opus 5 <[email protected]>
Under Slurm /proc/meminfo reports the node while the job is capped by its cgroup, so reading the node hands a planner a budget it does not have -- on exactly the allocations where planning matters. Prefer the cgroup limit, treat cgroup v2's "max" and v1's PAGE_COUNTER_MAX sentinel as absent, and apply a 0.8 fraction so the first prediction error is not an OOM kill.
plan_thread_budget charges every contig MIN_THREADS_PER_CHROM = 6 cores: four pipeline threads plus two HTSlib decode threads. On the sharded path the HTSlib pool does not exist (readers decompress inline) and only the executor is CPU-bound -- a measured 22-contig run put 16 threads on 2.02 cores. So it reserved ~6 cores for something consuming 0.4-1.0, and returned cc=7 on a 48-core box with 22 contigs. plan_sharded charges 1 + reader_workers and bounds concurrency by the fitted RAM law as well, which is what the scale bench's H3 verdict asked for: the invariant is in-flight bytes, not worker count. A budget too small for one contig is an error rather than a cc=0 plan that silently writes nothing.
Contigs are wildly unequal and were dispatched in whatever order the caller supplied, so the longest could start last. These estimates order them longest-first; rayon's work stealing does the rest. Only ratios matter -- the values are a sort key and nothing else -- so a coarse fallback tier is acceptable and the absolute unit is allowed to differ between tiers. An unestimated contig sorts first: guessing high costs a slightly worse order, guessing low risks the failure the ordering exists to prevent. Ties break by name so dispatch order stays deterministic.
w readers keep one executor fed when w/t_read >= 1/t_exec, so w is ceil(t_read/t_exec) -- rounding up, since rounding down starves the serial stage this exists to keep busy. Clamped at W_MAX=16, well clear of the knee of 3-7 the scale bench observed, so reaching the clamp flags a bad probe rather than silently over-provisioning.
Times one shard worker's chunk production against dense2sparse_vk's chunk consumption on a bounded prefix, so the reader count comes from a measurement on this input and this machine rather than a knee fitted elsewhere. One reader and no shards, or t_read is not one worker's rate and the ratio it feeds means nothing. Writes no user-visible bytes.
Concurrency now comes from plan_sharded under core and memory constraints rather than plan_thread_budget's core arithmetic, contigs dispatch longest-first so the longest cannot start last, and reader_workers can be probed from the largest contig with --tune. max_mem=None means a DETECTED budget, not unbounded -- a default behavior change. Unbounded preserves exactly the biobank-scale OOM exposure the byte-budgeted planner exists to remove. chroms keeps its original order everywhere except dispatch: finalize_fields and write_meta consume it as store layout, so reordering it there would move output bytes. Co-Authored-By: Claude Opus 5 <[email protected]>
detect_memory_budget() raised RuntimeError with no cgroup limit and no readable /proc/meminfo (every macOS run), turning from_vcf's default max_mem=None into a hard failure before any byte was read. Catch it, warn, and fall back to core-bound-only planning instead -- detection is an optimization, not a requirement, same principle as the tune probe's own failure path. Cross-reference from_vcf's and from_vcf_list's max_mem docstrings: they are different quantities (whole-process planning budget vs. one dense chunk) on the same class, and a value tuned for one silently means something else on the other. Also: filter src/logging.rs's two global-subscriber tests to their own tracing target, since CURRENT_SINK is process-global and was capturing tracing events emitted by contig_cost's and tune's own tests running concurrently on other threads (10/10 clean default-parallelism runs after the fix, vs. an observed 36-events-captured-instead-of-1 flake before); hoist the longest-first cost estimate above the tune branch so it is computed once instead of twice; and note at the dispatch site that `results` comes back in dispatch order, not chroms order, and why that is still safe. Co-Authored-By: Claude Opus 5 <[email protected]>
concurrent_chroms, reader_workers, and contig dispatch order all move under the tuned planner, and each can perturb chunk ordinals, per-chunk ledgers, or long-allele bank offsets. Pins the store digest across four schedules spanning the corners the planner can now reach, and across tune on/off.
The contig axis only ever compared cc=1 against cc=4. The new planner reaches cc=15 on a 48-core box, so the sweep needs points there to say whether its choice is the good one. SKILL.md states that max_mem=None is a detected budget rather than unbounded, since a reader will otherwise assume the opposite. The cc=4 corner is dropped from the new axis: at this corpus/reader_workers/ chunk_size it is byte-identical to the existing contig axis's c=22, concurrent=4 point (same point_id), so keeping it would silently duplicate a measurement under a second name.
The schedule-invariance fixture used chunk_size=64 with a max of 32 records per contig, so every contig fit in one chunk and the long-allele bank (spill threshold 13 bytes) was never written -- two of the three risks the gate claims to cover (chunk ordinals across chunks, bank offsets) went untested. Lower chunk_size to 8 so the largest contig spans four chunks, and give four of the eight contigs one over-threshold ALT each, planted mid-contig so bank offsets can interleave across schedules. This also gives the reader-rate probe (PROBE_CHUNKS=2) its intended two-chunk average instead of degenerating to one. The digest still holds across all four schedules with the strengthened fixture. Also: assert the max_mem-rejection test leaves no store on disk (it only checked for the exception before), and add scripts/bench_svar2/probe.py's reciprocal pointer to tests/_oracle.py::store_digest so the duplication obligation is visible from both sides. Switch schedule env hooks from os.environ.update/pop to monkeypatch.setenv so a pre-existing value would be restored rather than deleted.
build_plans.py already emitted a concurrency.json plan, but nothing ran it (sweep_scale.sbatch's dispatch loop didn't name it) and nothing loaded it (model.py's _SWEEP_NAMES didn't name it) -- the axis this task exists to add produced no data. Wired concurrency into both: added to _SWEEP_NAMES (with a docstring update and a main() fix, since the old positional generator-unpack over _SWEEP_NAMES silently assumed exactly five names) and to the sbatch dispatch loop. No new corpus generation needed -- every concurrency point reuses the s4000_c22 corpus the contig axis already generates. Also closes the identical omission for vlinear2 (added in 407a56c, never wired up): added to the sbatch dispatch loop, and -- the more serious half of the gap -- added the vlinear2 corpus-generation block that never existed, without which dispatching it would hit run_sweep's uncaught FileNotFoundError on a missing manifest and abort the whole sweep under set -e. vlinear2 was already present in model.py's _SWEEP_NAMES. Routed the concurrency loop through the existing _point() helper instead of constructing SweepPoint directly (verified byte-identical point_ids). Added a covering test asserting build_plans()'s axes match _SWEEP_NAMES exactly, so a future axis added to one without the other fails a fast unit test instead of shipping silently inert. holdout_f0's corpus is also never generated in sweep_scale.sbatch (a separate, pre-existing gap discovered in passing) -- left unfixed and flagged in the task report for follow-up.
…timates from over-tightening the memory budget `tune=True`'s probe can pick a reader-worker count up to 6.2x the default's per-contig memory cost, so a probe that SUCCEEDED could turn a conversion that plans fine untuned into InsufficientMemory -- contradicting `tune`'s documented "pure optimization" contract. `plan_sharded_tuned` retries once at the default worker count (and warns) before letting the error surface. Separately, the planner fed the NOMINAL chunk_size * per_variant_bytes into the RAM law, which was fitted against RESIDENT chunk bytes (min(chunk_size, variants) -- BitGrid3::zeros is a calloc, so an oversized chunk_size costs address space, not RSS). `estimate_contig_costs` now tags its return value (`ContigCosts`) with which tier produced it, so the planner can safely narrow by the largest contig's true record count only when those counts are exact (the index tier), never when they're header-length estimates in a different unit.
… self-explaining at runtime `max_mem` means a whole-process planning budget on `from_vcf` but a single in-flight dense-chunk cap on `from_vcf_list` -- deliberately not renamed (public API, plan-specified), but an order-of-magnitude mix-up between the two was easy to make silently. `PlanError::InsufficientMemory`'s message now points at `from_vcf_list` when it fires (a `from_vcf_list`-sized value passed to `from_vcf` lands under its ~1.2 GB practical floor and hits this error anyway). `from_vcf_list` now warns when `max_mem` is far above a single dense chunk's target size, catching the reverse direction, which otherwise silently derives an enormous `chunk_size` instead of raising.
…the root cgroup detect_memory_budget read /sys/fs/cgroup/memory.max / .../memory.limit_in_bytes directly -- the ROOT cgroup's own limit, which is unlimited on any host without a cgroup namespace. Under Slurm the real job limit lives several levels below the root (/slurm/uid_<uid>/job_<id>/...), so this silently planned against the whole node instead of the job's actual allocation: on this cluster, a real 64 GiB job cgroup read as an unbounded ~866 GB budget via the /proc/meminfo fallback, a 15.7x over-estimate with no warning. Resolves the process's actual cgroup path from /proc/self/cgroup first (v2: walks every ancestor up to the root and takes the minimum, since an ancestor can be more restrictive than the leaf; v1: reads the leaf's memory.limit_in_bytes, rejecting the uncapped sentinel), falling back to the old fixed root paths -- then /proc/meminfo -- only when self-discovery finds nothing. Verified against this exact job's real cgroup v1 limit (68,719,476,736 bytes): detect_memory_budget() now returns ~54.98 GB (0.8 x 64 GiB) instead of the prior ~865.7 GB.
…-invariance gate The digest-invariance gate's fixture plants a long ALT on chr2/chr4/chr6/chr8 specifically to exercise the long-allele bank, but nothing asserted the bank was actually non-empty -- digest-invariance can't distinguish "correctly empty" from "incorrectly empty", so a future edit shortening _LONG_ALT below MAX_INLINE_ALT_LEN would silently empty the bank and this gate, the branch's only end-to-end safety net, would still pass green. Asserts chr8/indel/long_alleles.bin (confirmed against src/layout.rs) exists and is non-empty on one representative store, which the digest check above already proved is byte-identical across every schedule.
…al in the concurrency axis The concurrency axis hardcoded chunk_size=10_937 even though the identical derived value (cs = size_corpus(4_000, CELLS_BUDGET)) was already in scope three lines above for the contig axis. If CELLS_BUDGET or MIN_CHUNKS moves, the contig axis follows and the concurrency axis silently doesn't -- which also invalidates that block's comment justifying the omitted cc=4 corner, since that argument depends on byte-identical point_ids between the two axes.
…rong warning message Re-review found the memory-budget test rewrite had lost its headline claim: after the self-discovery fix, no remaining test had both a resolvable cgroup limit AND a present /proc/meminfo, so a regression that reordered meminfo ahead of the cgroup -- i.e. reintroduced the Critical this file exists to guard -- would have passed every test in the file. Restores that precedence assertion for both the v2 and v1 tiers (v1 especially, since that's the tier that broke on this cluster), restores direct coverage of "v2 max falls through to v1", and adds the missing multi-controller v1 line fixture (9:cpu,memory:/path). Also fixes from_vcf_list's new large-max_mem warning, which unconditionally claimed "This will silently derive a very large chunk_size" even when an explicit chunk_size makes max_mem_bytes dead code instead -- the message now says which of the two actually happens. Documents both of Finding 2's fixes for `tune` (a successful-but-oversized probe also can't fail a conversion, not just a failed one) and the new from_vcf_list max_mem sanity warning in skills/genoray-api/SKILL.md, per this repo's public-API-doc-in-the-same-PR rule.
`SparseVar2.from_vcf_list`'s `max_mem` used to cap the bytes of one
in-flight dense chunk directly. `SparseVar2.from_vcf`'s `max_mem` (same
name) means a whole-process planning budget. Having one name mean two
orders-of-magnitude-different things on sibling methods was confusing and
error-prone -- the previous commit only papered over it with a sanity
warning. This makes `max_mem` consistently mean a whole-process budget
across both methods, per owner decision.
`from_vcf_list` has no fitted concurrency planner the way the sharded
`from_vcf` reader does (its contigs convert strictly sequentially, see
`orchestrator::run_vcf_list`), so it derives a per-chunk byte target from
the budget instead of planning concurrency:
chunk_target = min(_DENSE_CHUNK_TARGET_BYTES,
max_mem // (concurrent_jobs * in_flight_chunks_per_job))
`concurrent_jobs=1` (contigs run one at a time, confirmed in
`orchestrator.rs`). `in_flight_chunks_per_job=8` is a fixed constant read
directly off the Rust pipeline's architecture (a bounded reader/executor
channel of capacity 6, plus one chunk each thread may hold outside the
channel) -- not a fitted memory law, since none exists for this pipeline
(the sharded path's RAM_BASE_MB/RAM_PER_SAMPLE_MB/RAM_KAPPA were fitted
elsewhere and must not be reused here). `baseline=0` for the same reason.
`_DENSE_CHUNK_TARGET_BYTES` (~256 MiB) stays a ceiling so a large budget
can't grow chunks past today's size -- only a tight budget shrinks them.
`max_mem=None` now detects the budget the same way `from_vcf` does (80%
of the cgroup limit, or /proc/meminfo outside a cgroup), degrading to the
historical fixed dense-chunk target with a warning if detection fails.
Removes the now-obsolete `_MAX_MEM_SANITY_MULTIPLE` mixup warning (and its
test) and every doc/comment describing the two `max_mem`s as different
quantities, across from_vcf/from_vcf_list docstrings, src/budget.rs's
InsufficientMemory message, and skills/genoray-api/SKILL.md.
BREAKING CHANGE: `SparseVar2.from_vcf_list(max_mem=...)` no longer caps a
single dense chunk; it is now a whole-process byte budget, matching
`SparseVar2.from_vcf`. A call like `max_mem="512MiB"` that used to mean
"let one dense chunk use up to 512 MiB" now means "the whole process
should use at most 512 MiB", which derives a much smaller chunk_size.
Callers relying on the old per-chunk-cap meaning must re-tune their
`max_mem` value.
Co-Authored-By: Claude Opus 5 <[email protected]>
Follow-up to e5f2b34, fixing three Important findings from code review: 1. `_VCF_LIST_CONCURRENT_JOBS`/`_VCF_LIST_IN_FLIGHT_CHUNKS_PER_JOB` were Python literals duplicating Rust architecture (contig sequencing, reader/executor channel capacity) with only a one-way doc pointer. Export both as named Rust constants (`orchestrator::VCF_LIST_CONCURRENT_CHROMS`, `orchestrator::VCF_LIST_DENSE_CHANNEL_CAP`) through the `_core` PyO3 module so Python reads them live instead of hardcoding copies that could silently drift (in the unsafe, under-estimating direction) if the orchestrator's concurrency or channel tuning ever changes. 2. Fix a reintroduced bug: the max_mem detect/parse block ran unconditionally ahead of the `chunk_size is None` guard, so a caller passing an explicit `chunk_size` on a host where memory-budget detection fails got a spurious "could not detect a memory budget ... deriving chunk_size ..." warning for a derivation that wasn't happening. Move the block inside the guard so detection is only attempted when chunk_size is actually being derived. Adds a regression test. 3. State plainly, in the paragraph that makes the whole-process promise (not just the chunk_size paragraph), that max_mem only sizes the dense-chunk term -- it does not cap from_vcf_list's per-input-file reader overhead, which scales with cohort size independent of max_mem. Also tightens a budget.rs test to assert the actual remedy phrase instead of two independently-satisfiable substrings, and clarifies that the in-flight dense-chunk accounting excludes the separate (tiny) SparseChunk queue. Co-Authored-By: Claude Opus 5 <[email protected]>
…line comment The comment explaining why EventSink's per-chrom progress buffer is safe under concurrent dispatch still attributed `concurrent_chroms` to `plan_thread_budget`. Since the tuned load-balancing change it comes from `plan_sharded`, which plans against a different per-contig cost model. Comment-only; no behavior change. Co-Authored-By: Claude Opus 5 <[email protected]>
Two stacked efforts on one branch: the conversion scale-bench harness (earlier commits), and the tuned load balancing it was built to justify (this plan).
Why
plan_thread_budgetcharged every contig 6 cores — four mostly-blocked pipeline threads plus an HTSlib decode pool the sharded path never allocates. A measured 22-contig run put 16 threads on 2.02 cores, so it returnedcc=7on a 48-core box with 22 contigs, and nothing bounded memory at all.What changed
budget::plan_sharded1 + reader_workers) and the fitted RAM lawcontig_costtunet_read/t_execto size readers from measurement, not a fitted knee_utils.detect_memory_budget/proc/meminforeports the node, the cgroup reports the job)tests/test_svar2_schedule_invariance.py(cc, w, order)SparseVar2.from_vcf_list(max_mem=)changes meaning: it was a cap on one in-flight dense chunk; it is now a whole-process budget, matching the newfrom_vcf(max_mem=). An existingfrom_vcf_list(max_mem="512MiB")caller now says "the whole process may use 512 MiB" and gets a correspondingly smallerchunk_size.max_mem=Nonemeans a detected budget (cgroup limit × 0.8), not unbounded. Detection failure warns and degrades to core-bound planning rather than raising.Measured result
Run on carter-cn-04,
--cpus-per-task=48 --mem=64G, at commit7fbdc49. Corpus: S=4,000, 22 contigs, 349,998 variants,chunk_size=10,937. Best-of-3 after a page-cache warm-up rep (probe.run_point).The planner, given no overrides, chooses
concurrent_chroms=11, reader_workers=3and converts in 8.05 s using 10.87 of 48 cores (3,044 MB peak RSS).cc=1, w=12cc=4, w=3That is 1.46× faster than the 11.74 s target and 4.36× faster than
cc=1.The
concurrencysweep axis, measured on the same node and corpus at fixedw=3:Three things worth stating plainly, because none of them match the plan's prediction:
cc=15; the planner pickscc=11—(48−1)/(1+3). The prediction was wrong in the lucky direction:cc=15measures 9.21 s, so following it would have been slower than what the planner actually does.cc. It bottoms out at or belowcc=8(7.59 s) and degrades after. Beyondcc≈8the extra contigs buy no wall time — measured core usage saturates around 10 regardless — whilecpu_sclimbs from 76 to 98 and peak RSS more than doubles. The planner'scc=11is on the flat shoulder, ~6% off the best sampled point, not on the optimum. It clears the bar it was built to clear; it is not a tuner.cc=1row is 34.10 s against the baseline'scc=1row at 35.11 s — 2.9% apart. Node speed on this cluster varies by 2.08×, so without that control a cross-run comparison would be assumption, not measurement.Schedule invariance holds on real data: all four fixed-
ccpoints and the planner's own run produce the byte-identical store digest8dbfe4e829fbf0e8— the invarianttests/test_svar2_schedule_invariance.pygates in CI, here at 22 contigs and 350k variants.The memory-budget fix is confirmed live: inside the 64 GiB Slurm cgroup,
detect_memory_budget()returns 54,975,581,388 bytes = exactly 0.8 × 64 GiB. Before the fix it read the root cgroup and returned 865.7 GB.Findings worth flagging
Three real bugs surfaced during review, none of them in the original plan:
contig_costpassed a BCF-headerridtohts_idx_get_stat, but for bgzipped VCF the.csiis tabix-numbered over covered contigs only, and that function does no bounds checking. A test written before the fix showed a zero-record header contig silently misattributing another contig's count and reading past the array. Now resolved viatbx_name2idwith every id bounded byhts_idx_nseq.detect_memory_budgetread the root cgroup, returning 865.7 GB against a real 64 GiB job limit — a 15.7× silent over-estimate, on exactly the platform its docstring names. Now resolves the process's own cgroup (v2 walks ancestors and takes the minimum).Verification
cargo check --no-default-features(query-core gate gvl links against) — cleancargo test --no-default-features --features conversion— 439 passed, 0 failedpixi run pytest tests/ -q -m "not network"— 963 passed, 5 skipped, 16 xfailed, 0 failedKnown follow-ups
Filed rather than folded in:
sweep_scale.sbatchnever generates theholdout_f0corpus it dispatches a plan point for (sameFileNotFoundError-under-set -emechanism fixed here forvlinear2).concurrent_chroms=11, htslib_threads=2, monolithic_reader_active=66against a 48-core budget.run_conversion_pipelinenow takes 20 parameters, including two structurally indistinguishable field-manifest vectors.🤖 Generated with Claude Code