3/4 perf: external memory, 24.2 s -> 1.56 s on chr21 FASTA - #5
3/4 perf: external memory, 24.2 s -> 1.56 s on chr21 FASTA#5BenjaminDEMAILLE wants to merge 22 commits into
Conversation
Two prerequisites for the performance work that follows, neither of which changes behaviour. The crate had no test covering the LCP array. That is the riskiest possible gap for this algorithm: the public entry points discard the array, but the *next* merge level consumes it in the three-case decision, so a single wrong LCP entry silently reorders suffixes one level up and the SA comes out subtly wrong. Add four tests that check `lcp[0] == 0` and `lcp[i] == lcp(text[sa[i-1]..], text[sa[i]..])` against a naive oracle, over fixtures, random texts across four alphabet sizes, long runs and periodic text, and finite `max_context`. `bench/README.md` claims the published numbers were taken with fat LTO and one codegen unit, supplied by a parent workspace. That workspace is not in this repo, so every build made from it since the crate went standalone has used `lto = false, codegen-units = 16`. Pin the profile here. In a library crate `[profile.release]` applies only when this crate is the workspace root, so it affects this repo's own tests, examples and benches and is invisible to downstream consumers. Measured on Apple M4 Max, 12 threads, chr21 FASTA (47.5 MB): 27.8 s -> 24.0 s wall. Neutral on N-free DNA. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The CaPS-SA merge kernel stays the general path, but it is not the right
algorithm for the most common request: the standard lexicographic suffix
array of a byte text, unsegmented, with no context bound. Add
`src/radix.rs` and route that case through it.
Profiling separated two distinct costs, and the merge kernel pays both:
* Step count. `n log n` merge steps, most resolved by a symbol
comparison at a random text address. 80 MB of N-free DNA is ~2.1e9
steps at ~13 ns each.
* Scan length. Every leaf merge starts at `m = 0`, so two suffixes
sharing a long prefix cost a scan proportional to that prefix. Genome
FASTA carries megabyte-scale runs of `N` (period-61 once line wrapping
is included), where one comparison scans millions of bytes. That drives
the cost per merge step from 13 ns to 222 ns, a 16x penalty which is
entirely scan time. This is what made real chr21 20x slower than
N-free DNA of comparable size.
The new path removes both. It sorts by a packed fixed-depth key, then
resolves the remainder by prefix doubling on ranks. The packing picks the
narrowest field width in {1,2,4,8} bits that holds the alphabet, so DNA
over {0,1,2,3} resolves 32 symbols per key rather than the 8 a raw byte
key gives. After the seed no comparison reads the text again, so a
megabyte run of `N` costs exactly what random DNA costs.
The seed sorts by `(key, min(n - p, k))`. The second component is
required, not cosmetic: zero-padding makes a short suffix share a key
with any suffix continuing in zeros, and `0` is a real symbol in every
DNA encoding. Ordering by visible length puts the proper prefix first,
which is the crate's shorter-is-smaller convention. Without it `[0, 0]`
leaves two positions permanently tied and doubling cannot terminate.
Guards are soundness conditions, not heuristics, and all three default
to declining:
* `max_context` must be unbounded; a finite bound makes the merge's
comparator fall through to `boundary_order`, which compares lengths,
so it is not lexicographic.
* `LimitProvider::plain_lex_len` must report the full text. New method,
defaulting to `None`, overridden only by `PlainText`. An impl that
delegates `lim_at` to `PlainText` but overrides `boundary_order` for a
different convention (STAR's spacer-as-largest) inherits `None` and
stays on the merge kernel without changing a line.
* `S` must be exactly `u8`. Packing wider symbols into an order-
preserving key is endianness-dependent: for `u16` on a little-endian
host `0x0100 > 0x0001` as values but their byte views compare the
other way. The rest of the crate avoids this only because
`LcpDispatch` resolves equality over bytes and recovers ordering
through `S: Ord`.
Tests: exhaustive over every binary text to length 10 and every ternary
text to length 6, random texts across seven alphabet widths, texts where
a real `0` collides with padding, long runs, periodic text, and the
wrapped-FASTA `N`-block shape.
Measured on Apple M4 Max, 12 threads, against the previous kernel, with
byte-identical suffix arrays on both real inputs:
chr21 fwd+revcomp, N-free, 80 MB 6.08 s -> 1.14 s CPU 28.1 s -> 5.0 s
chr21 FASTA, 47.5 MB 27.8 s -> 1.16 s CPU 283.5 s -> 5.2 s
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Checking a new construction path against the old one only shows the two agree. Checking adjacent suffixes directly is O(n · lcp), which is unusable on exactly the repetitive inputs that need checking most: on a chr21 FASTA a single adjacent pair can share megabytes. Use the fixpoint characterisation instead. With `rank` the inverse of `sa` and `f(p) = (text[p], rank[p + 1])`, taking `rank[n]` as less than every real rank, a permutation of `0..n` is the suffix array of `text` if and only if `f` is strictly increasing along it. That is one pass to invert plus one pass to compare, independent of any construction algorithm and independent of LCP length. Exposed as `caps_sa::verify_sa` and wired to a `--verify` flag on the bench CLI, off by default so it never contaminates a timing run. Full-scale results on Apple M4 Max, 12 threads: chr21 fwd+revcomp, N-free, 80,177,238 entries verify OK in 1.14 s chr21 FASTA, 47,488,540 entries verify OK in 0.57 s Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Phase timings (`CAPS_SA_PROFILE=1`) showed the two rank-scatter passes
were fully sequential and had become the largest single cost after the
seed sort: 0.34 s of a 1.17 s build on 80 MB of DNA.
Both scatter through a permutation -- the target index is `sa[i]`, not
`i` -- so the writes are not expressible as disjoint sub-slices and
`split_at_mut` does not apply. They are nonetheless disjoint: `sa` is a
permutation and the ranges being processed partition its index space, so
every slot is written exactly once. Introduce a small `Scatter` wrapper
that encodes precisely that contract in its `unsafe fn set`, and drive
both passes with rayon.
Grouping now has each index decide for itself whether it starts a group;
the index that does owns the group, walks it to find the end, and writes
its members' ranks. Exactly one owner per group, and `collect` on an
indexed parallel iterator preserves order, so the group list still comes
out sorted, which is what `split_disjoint` relies on.
Also materialise the successor ranks before sorting each group.
`sort_unstable_by_key` re-evaluates its key function O(len log len)
times and every evaluation was a random probe into `rank`; paying once
per element makes the sort's memory traffic sequential.
Apple M4 Max, 12 threads, suffix arrays byte-identical to the previous
kernel and independently `--verify`-checked:
chr21 fwd+revcomp, N-free, 80 MB grouping 0.341 s -> 0.075 s
total 1.17 s -> 0.89 s
chr21 FASTA, 47.5 MB grouping 0.172 s -> 0.037 s
total 1.10 s -> 1.04 s
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The repository had no Rust CI at all: the only workflow deployed the docs site, while 73 tests sat in the tree with nothing running them. That is not a safe baseline for changing the sorting kernel. Covers both architectures that matter here, since the LCP kernel and the pooled external-memory bucket path are the parts that diverge per platform: macOS is aarch64/NEON, Ubuntu is x86_64/AVX2. Runs the tests in debug as well as release, because debug is what exercises the `debug_assert`s guarding the unchecked scatter in `radix.rs` and the buffer-length invariants in the merge kernel. Adds fmt, clippy with warnings denied, a check against the declared 1.89 MSRV, and rustdoc with broken links denied. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
`bench/run.sh` compares against upstream C++ and needs both binaries
prebuilt. Add a self-contained script for the case that prompted this
work: it fetches hg38 chr21 and prepares *both* inputs, which is the
distinction that made the original slowdown report hard to interpret.
chr21.0123 forward ++ revcomp, one byte per base, codes 0..=3,
ambiguous bases dropped. ~80 MB, alphabet 4, no long runs.
The input libsais is normally benchmarked on.
chr21.fa the raw FASTA, still carrying its ~6.6 Mb of `N`. Wrapped
at 60 columns, so the `N` blocks are a period-61 repeat
rather than a plain run.
Benchmarking one implementation on the first and another on the second
compares two different problems. The second is the realistic input and
the one that used to be pathological.
Builds with `-C target-cpu=native` and runs each case through
`--verify`, so the harness reports correctness alongside timing.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Adds to README: the in-memory fast path with measured numbers, a "Choosing a path" table giving the three soundness conditions and why each one exists, and `verify_sa` usage. Adds to bench/README: the chr21 section, including the per-merge-step arithmetic that separates the two inputs (13 ns/step on N-free DNA against 222 ns/step on the FASTA, same kernel and same machine), and the phase breakdown of the fast path. Corrects two claims that were misleading: The build paragraph attributed `lto = "fat"` and `codegen-units = 1` to a parent workspace. No such workspace is in the repository, so from the commit that made the crate standalone until the profile was added, every build made from this repo used `lto = false, codegen-units = 16`. Numbers taken in that window are not comparable with numbers taken now. The 97.54% `lcp_u8_avx2` profile was read as "LCP scanning is expensive", which motivated widening the scan through AVX2, AVX-512 and the hybrid. The LCP kernel is also where the two random text loads happen, so on short-LCP input those samples are load stalls and a wider vector cannot help. The existing AVX-512 ablation already showed this: the 64-byte-only variant was 16% slower on rand100m and only the long-LCP human slice gained. Both readings are now stated. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Public documentation linked to `FilteredSource`, which is private, so `cargo doc` fails under `RUSTDOCFLAGS=-D warnings`. Pre-existing, but it blocks the rustdoc job added in the previous commit. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The seed was a `par_sort_unstable` over a materialised `Vec<(u64, u32, I)>`, which was both the largest remaining cost and the peak-memory driver. Three things change. The source is never materialised: keys are recomputed from `text` in the histogram pass and again in the scatter, which trades a random read of an n-element key array for a sequential read of the text. Peak memory drops to the two destination buffers, 12 bytes per position at `I = u32` against the 16 a `(u64, u32, I)` record costs. And the top-level partition becomes a counting pass, which parallelises evenly, where a parallel comparison sort's first partitioning steps are close to serial. 11 bits (2048 buckets) keeps the write-combining state near 512 KB and inside a core's private cache. 16 bits would need 16 MB of open write lines and thrash the TLB instead. The visible-length tie-break no longer needs storing. Only the last `k - 1` positions can have a visible length below `k`, so it is a function of the position alone and is applied in the per-bucket sort and in the group scan. Apple M4 Max, 12 threads, chr21 fwd+revcomp 80 MB, output unchanged and `--verify` clean: peak RSS 2.83 GB -> 2.21 GB (-22%) seed 0.341 s -> 0.289 s total 0.892 s -> 0.844 s Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Closes two of the three gaps listed as limitations when the fast path landed. **Subsets.** Doubling cannot be restricted to a subset directly, because a round compares `rank[p + d]` and that successor is generally not in the subset, so ranks must be defined for every position in the text. Building the whole array and filtering it sidesteps that, and the filter is one O(n) pass since the full array is already ordered. Worth it when the subset is a real fraction of the text, which is what this API exists for: STAR-style indexing keeps every ACGT position and drops only spacers. Below one eighth of the text it declines, because O(n) to build and discard would dwarf the O(m log m) the merge kernel needs. That ratio is a performance heuristic; the guards it sits behind remain correctness conditions. It also declines on duplicate or out-of-range positions. The output is a permutation of the input *multiset*, which a membership filter cannot reproduce. **In-memory sample sort.** `build_in_memory_sample_sort` exists to sort in RAM, so where the doubling path applies it is strictly better: same output, no bucket machinery, and none of the scan cost on repeat-heavy text. `build_ext_mem` deliberately does *not* get this. Its purpose is to bound peak memory, and routing it through an in-memory algorithm would defeat exactly that. It stays on the merge kernel, and the remaining limitation is now stated as a deliberate choice rather than an omission. Apple M4 Max, 12 threads, chr21 fwd+revcomp 80 MB. `--in-mem-ss` output verified identical to the in-memory path's: --in-mem-ss 3.41 s -> 1.18 s wall, 34.5 s -> 7.2 s CPU Tests: duplicate positions keep their multiplicity, a tiny subset of a large text still matches brute force through the merge kernel, and out-of-range positions still panic rather than silently returning a wrong answer. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Records the counting-sort seed (peak RSS 2.83 GB -> 2.21 GB, total 0.89 s -> 0.84 s on the 80 MB input), adds the `--in-mem-ss` row, and replaces the "everything else uses the merge kernel" line with what is now actually true: subsets and in-memory sample sort are covered, and `build_ext_mem` stays on the merge kernel deliberately, because bounding peak memory is the whole point of that path. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
`build_ext_mem` was left on the merge kernel because bounding peak
memory is its purpose and prefix doubling needs a rank for every
position in the text. That left it paying the full scan cost on
repeat-heavy input. Profiling put 94% of an ext-mem run on chr21 FASTA
in phase 1:
chr21.0123 (N-free, 80 MB) phase1 1.15 s of 3.49 s total
chr21.fa (47.5 MB) phase1 22.85 s of 24.19 s total
Normalised that is 0.014 s/MB against 0.48 s/MB, the same 34x scan
penalty the in-memory path had, and for the same reason: two suffixes
inside a long repeat agree for as far as it continues, so one comparison
scans megabytes.
Fix it in the comparator rather than the algorithm, which keeps the
memory bound intact. If `text[s..e)` has period `q` and two suffixes
start at `a < b` inside it with `(b - a) % q == 0`, they agree until the
later one reaches `e`, so `lcp(a, b) >= e - b` is known in O(1) from the
run's bounds with no scanning. When the phase does not match, the two
must differ within `q` symbols and the ordinary scan is already short.
The scan is additionally bounded so it stops at a run's start rather
than traversing it.
Detecting only single-symbol runs would have missed the case that
actually occurs: in wrapped FASTA an `N` block is 60 `N`s then a
newline, which is period 61, not period 1. Periods up to 64 are
considered.
Detection is two-stage so texts without runs pay almost nothing. A
sampling pass looks for any periodic window and collects the periods
that occur; the full scan runs only for those. On N-free DNA the sample
finds nothing, the table is empty, and every query short-circuits on a
slice-empty check. The table itself is a few dozen entries, so the
memory bound is untouched.
`Cmp` bundles the SIMD dispatch with the run table and replaces the bare
`LcpDispatch` threaded through the merge kernel, so phase 1, the phase-3
pivot searches and the phase-4 cascade all benefit. It stays `Copy` and
still travels through the recursion in registers.
Apple M4 Max, 12 threads. Ext-mem output verified identical to the
in-memory suffix array on both inputs:
chr21.fa 24.19 s -> 2.48 s wall, 268 s -> 23.3 s CPU
phase 1 22.85 s -> 0.95 s
peak RSS 147 MB -> 151 MB
chr21.0123 3.49 s -> 3.55 s wall (unchanged; no runs to find)
Tests: the run-aware LCP is checked against a byte-at-a-time oracle over
sampled position pairs on homopolymers, wrapped-FASTA blocks and
multi-period texts, plus detection shape (sorted, disjoint, period claim
actually holds) and `max_context` behaviour inside a run.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Adds the ext-mem phase breakdown that located the cost (94% of the chr21 FASTA run in phase 1), the before/after table, and the detail worth keeping: a homopolymer detector would not have worked, because in 60-column wrapped FASTA the longest single-byte run is 60. Each line of `N`s ends in a newline, so the real structure is a period-61 repeat spanning 6.6 Mb. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Phase 4 had become the largest ext-mem cost, and its merge CPU (19.2 s
on 80 MB of DNA) exceeded the entire CPU of the in-memory path. The
merge is latency-bound: the tied branch dereferences the text at two
random addresses, and the address for step i+1 is not known until step i
retires, so there is no memory-level parallelism and the hardware
prefetcher cannot see the pattern.
The candidate positions themselves live in the two index arrays, which
are sequential and already in cache, so the addresses several steps
ahead are known even though the dependent loads are not. Issue them as
prefetches, offset by the current boundary LCP `m`, which estimates
where the next scans start.
This is not the prefetch recorded as a negative result in `lcp.rs`. That
one sat inside the strided scan loop, which the hardware prefetcher
already covers. This one targets the random access, which it cannot.
Apple M4 Max, 12 threads, ext-mem, output verified identical to the
in-memory suffix array on both inputs:
chr21.0123, 80 MB phase4 merge CPU 19.20 s -> 12.42 s
phase4 wall 2.17 s -> 1.46 s
total 3.55 s -> 2.90 s CPU 33.8 -> 27.0 s
chr21.fa, 47.5 MB phase4 merge CPU 7.99 s -> 6.76 s
total 2.56 s -> 1.99 s CPU 23.3 -> 18.0 s
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Sweeping `p` directly showed the previous default sitting on the wrong
side of a flat region. Total work is `n log n` either way, since a
smaller `p` moves levels out of phase 4's per-partition cascade and into
phase 1's merge sort, but the constants differ: phase 3 shrinks
quadratically in `p`, and phase 4's cascade does a full pass over its
partition per level.
Peak RSS is set by phase 4 holding `4 x threads` partitions of `n / p`
records at once, so it only starts growing once `p` is small enough for
that product to rival the text. Measured on chr21 forward ++ revcomp
(80 MB), 12 threads:
p total peak RSS
48 2.58 s 987 MB
96 2.49 s 538 MB
306 2.85 s 282 MB
612 2.80 s 202 MB <- 131072
1224 3.05 s 205 MB <- 65536 (previous default)
128Ki is the largest step that costs nothing in memory: same peak RSS,
~8% less wall. Going further trades real memory for speed, which is the
opposite of what this path is for, so it stays available through
`ExtMemOpts::subproblem_count` rather than becoming the default.
At genome scale this changes nothing. `PHASE1_MAX_PARTITIONS` already
binds for any n above ~1 GB, so GRCh38 still gets p = 8192.
Output verified identical to the in-memory suffix array on both inputs.
chr21.0123, 80 MB 3.01 s -> 2.77 s
chr21.fa, 47.5 MB 2.06 s -> 1.86 s
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Phase 1 had become the largest external-memory phase (1.20 s of 2.54 s on 80 MB of DNA) and was already parallelising at ~97%, so the way forward was less work rather than more threads. It was sorting every subarray from singletons, which is the case the merge kernel handles worst: each leaf merge starts at `m = 0` and orders two suffixes by scanning the text at two random addresses. Sorting by the packed key first resolves the leading `k` symbols with no text access at all (32 symbols for DNA at 2 bits each), and yields the LCP between adjacent runs for free from `(key_a ^ key_b).leading_zeros()`. Only suffixes agreeing through all `k` symbols reach the merge kernel, on the short slice they occupy. This is the bounded-memory counterpart to the prefix doubling in `build_in_memory`. Doubling itself is not available here: it needs a rank for every position in the text, which is exactly the memory this path refuses to spend. A fixed-depth key needs none. The cross-run LCP is capped by both suffixes' lengths. Padding can agree with a real `0` symbol past the end of the shorter suffix, so the raw `leading_zeros` count can overstate it, and a wrong LCP would silently corrupt the order at the next merge level. Gated on the same conditions as the other fast paths (`u8` symbols, plain lexicographic comparator, unbounded `max_context`) and falls back to `merge_sort` otherwise. The alphabet scan that picks the field width runs once per build, not once per subarray. Apple M4 Max, 12 threads, output verified identical to the in-memory suffix array on both inputs: chr21.0123, 80 MB phase1 1.201 s -> 0.316 s total 2.54 s -> 2.11 s chr21.fa, 47.5 MB phase1 0.898 s -> 0.415 s total 1.86 s -> 1.65 s Peak RSS 205 MB -> 233 MB: the key vector is 16 bytes per record over one subarray per worker. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
`CascadeWorkspace::merge_one_level` walked its run pairs sequentially, even though the pairs at a level are independent and write to disjoint destination ranges. The only thing forcing the order was the running `src_off` / `dst_off`, and both are prefix sums, so they can be computed up front and each pair handed its own sub-slices. This is what capped phase 4's parallel efficiency. Partition-level parallelism (`4 x threads` at once) hides it while there are many partitions in flight, but each partition's cascade ends in a single 2-way merge over the whole partition, and those tails serialise. Apple M4 Max, 12 threads, output verified identical to the in-memory suffix array on both inputs: chr21.0123, 80 MB 2.11 s -> 2.02 s chr21.fa, 47.5 MB 1.65 s -> 1.55 s A grain-size threshold was tried, on the theory that the small early levels would not pay for their rayon tasks. It measured worse on both inputs (2.22 s and 1.63 s), so the pairs are merged in parallel at every level. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Emitting is Theta(n) and single-threaded by construction: the caller's closure is `FnMut` and its ordering is the point. But it also did not overlap anything, so once per chunk every worker sat idle while the main thread drained the merged results. A scoped producer now merges chunk c+1, itself rayon-parallel, while the main thread emits chunk c. The channel bound of one keeps at most two chunks resident, so the transient cost is one extra chunk of merged positions rather than the unbounded queue an unsynchronised producer would build. This is not the existing `ordered_phase4_emit` path, which coordinates at *partition* granularity through an mpsc channel and a `BTreeMap` and measured slower than plain collect-then-emit. Here the producer hands over whole chunks that are already in order, so the consumer only drains them and no reordering structure is needed. That path is left untouched behind its opt-in flag. Apple M4 Max, 12 threads, output verified identical to the in-memory suffix array on both inputs: chr21.0123, 80 MB 2.02 s -> 1.82 s peak RSS 233 -> 240 MB chr21.fa, 47.5 MB 1.55 s -> 1.45 s peak RSS 182 -> 194 MB Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
A partition arrives as `p` sorted sub-subarrays, and the cascade merged
them pairwise in `log2(p)` levels, each a full LCP-enhanced pass. That
was the largest single cost left in the ext-mem build: 15.4 CPU-seconds
of a 15.1-second run on 80 MB of DNA.
When a packed key applies, discarding that sortedness and re-sorting the
partition outright is much cheaper. One key sort resolves the leading
`k` symbols with no text access, and only suffixes agreeing through all
of them reach the merge kernel, so `log2(p)` passes collapse into one.
Gated on the run table being empty, which is the point worth recording.
A long periodic run is exactly a stretch where a fixed-depth key
resolves nothing, since every suffix inside it shares the whole key, so
the re-sort would hand the merge kernel one enormous tied group and
throw away ordering phase 1 had already established. Measured
unconditionally it was a large regression on the `N`-heavy input:
cascade key re-sort
chr21.0123 (no runs) 1.95 s 1.65 s
chr21.fa (6.6 Mb N) 1.45 s 2.23 s
So the run table, which already exists to make scans cheap, doubles as
the predicate for whether a fixed-depth key can be expected to pay off.
Apple M4 Max, 12 threads, output verified identical to the in-memory
suffix array on both inputs:
chr21.0123, 80 MB 1.85 s -> 1.65 s phase4 merge CPU 15.4 -> 7.8 s
chr21.fa, 47.5 MB 1.45 s -> 1.47 s (unchanged; takes the cascade)
Peak RSS on the run-free input rises 240 MB -> 285 MB: the key vector
and the sort buffers are live for each of the `4 x threads` partitions
in flight. Still an order of magnitude under the in-memory path's
2.2 GB, which is the comparison that matters for choosing this path.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
|
Part of the stack tracked in #7. |
|
Status against 0.7.0: this does not rebase as it stands, and most of it should What survives that is the phase-1 packed-key seed, re-derived against 0.7.0 in Leaving this open but stale rather than force-pushing something unreviewable. |
build_ext_memstays on the merge kernel by design — prefix doubling needs a rank for every text position, which would defeat the memory bound the path exists to provide. So it was still paying the full scan cost. Profiling put 94% of a chr21 FASTA run in phase 1:Six changes, each measured separately:
The interesting ones
Run skipping. If
text[s..e)has periodqand two suffixes start ata < binside it with(b - a) % q == 0, they agree until the later reachese, solcp(a, b) >= e - bis known inO(1)with nothing scanned. A homopolymer detector would not have worked: in 60-column wrapped FASTA the longest single-byte run is 60, because eachNline ends in a newline. The real structure is a period-61 repeat spanning 6.6 Mb. Detection is two-stage, so run-free texts pay one sampling sweep and get an empty table.Prefetch. The merge is latency-bound: the tied branch dereferences the text at two random addresses and the next address is unknown until the current step retires. The candidate positions live in the index arrays, which are sequential and cached, so addresses several steps ahead are known. This is not the prefetch recorded as a negative result in
lcp.rs— that one sat inside the strided scan, which hardware already covers.Subarray target. Sweeping
pshowed the old default on the wrong side of a flat region. 128Ki is the largest step that costs nothing in memory; going further trades real memory for speed, so it stays caller-controlled.Partition re-sort, gated on the run table. Unconditionally it was a large regression on the
N-heavy input (1.45 s → 2.23 s): a long run is exactly where a fixed-depth key resolves nothing, so the re-sort hands the merge kernel one enormous tied group and discards ordering phase 1 had established. The run table, which exists to make scans cheap, doubles as the predicate for whether a fixed-depth key will pay off.Negative result, recorded
A grain-size threshold on the parallel cascade levels measured worse on both inputs (2.22 s and 1.63 s), so every level goes parallel.
Peak RSS 214 MB → 285 MB, still an order of magnitude under the in-memory path. Ext-mem output verified identical to the in-memory suffix array on both inputs at every step.
🤖 Generated with Claude Code