add some doc about optimisation for query

DevDocMD/architecture/siblings.md

@@ -756,3 +756,129 @@ scan the full unsampled index — never threaded `--subsample`/`--entropy`,
 out of scope here since the reported problem was specifically about the
 `--sankoff`/`--tnt` pipeline's redundant/inconsistent scans, not these
 two standalone flags.
+
+## `query` never benefits from sparse row-major access (found 2026-08-19, not implemented)
+
+Benchmarked `obikmer query` against `global_index_presence` (dense-packed)
+vs. `global_index_presence_sparse` (`pack --sparse`), 100k simulated reads
+× 2 specimens (`benchmark/`, see
+[benchmark_query_testing.md](../implementation/benchmark_query_testing.md)).
+Correctness: 0 mismatches — sparse and dense return bit-identical query
+results. Performance: sparse consistently *slower* than dense (~30-50%,
+reproducible across two runs with warm disk cache), the opposite of
+`pack --sparse`'s stated intent ("faster for single-row access... like
+query").
+
+**Root cause, read from source, not measured in isolation:**
+`KmerPartition::query_partition_with` (`obikpartitionner/src/query_layer.rs:155-220`)
+is architecturally column-major: stage 2 walks `for g in 0..n_cols { for
+slot in hit_slots { layer.col_value(g, slot) } }`, documented (correctly)
+as the right locality strategy for the packed/columnar formats, where
+`col_value` → `PersistentBitMatrix::get` is a genuine O(1) mmap'd column
+read (`persistent.rs:110-113`).
+
+For `Self::Sparse`, that same `get(c, slot)` (`persistent.rs:114-118`)
+allocates a full `n_cols`-wide buffer and calls `fill_row` — materializing
+the *entire row* — just to return one cell. Called from inside the
+column-major double loop, this reconstructs the same row once per genome
+column touched: O(hits × n_cols) full-row rebuilds instead of O(hits).
+`PersistentSparseBitMatrix`'s own native row-major decode
+(`for_each_genome_in_row`, `sparse.rs:164-177`, used correctly by its own
+`row`/`fill_row`/`fill_sub_matrix`) is never reached from the query path
+at all.
+
+**`fill_sub_matrix` (the existing `BinaryMatrix` trait primitive,
+`traits.rs:13-37`) is not the right replacement for `query` either**, even
+once its own dispatch bug is fixed (see next section) — its output shape
+is inherently column-dense: `out[col]` gets an entry for every column,
+including columns with zero hits among the requested slots. On real
+sparse data (a hit typically touching a handful of genomes out of dozens)
+that's still O(n_cols) output regardless of true sparsity. What `query`
+actually wants is the sparse triple stream `(slot, col, value)` it already
+consumes as `QueryHit::Value` — not a materialized sub-matrix.
+
+**Proposed primitive** (design only, not implemented — explicit ask: keep
+count matrices *not excluded*, even though effort right now is
+presence/absence only):
+
+Not a closure-driven `for_each` — a real `Iterator`, one concrete struct
+per matrix format, so the traversal state (current position in the sorted
+slot list, current column, permutation, sparse-row decode cursor…) lives
+in named struct fields instead of being threaded implicitly through
+recursion or a captured closure. RPITIT (stable since 1.75, and this
+workspace is edition 2024) means the trait method can return it without
+naming or boxing the concrete type:
+
+```rust
+/// Yields every nonzero cell among `slots`, in implementation-defined order.
+fn nonzero_iter<'a>(&'a self, slots: &'a [usize]) -> impl Iterator<Item = (usize, usize, u32)> + 'a;
+// item: (idx into `slots`, col, value)
+```
+
+This is the one foundational primitive per format — both `fill_sub_matrix`
+and a `for_each`-style callback become trivial, free consumers of it
+(`.for_each(f)` is already `Iterator::for_each` from std; `fill_sub_matrix`
+becomes "drain the iterator, scatter into `out[][]`"), instead of two
+independently-maintained traversals that can silently diverge (see the bug
+below — this is exactly how it happened).
+
+- **On `PersistentSparseBitMatrix`**: the struct is nearly free to write —
+  it wraps the existing (currently private) `for_each_genome_in_row`
+  per-row decode, advancing to the next `slots` entry on exhaustion. O(Σ
+  row nnz), zero `n_cols`-wide allocation.
+- **On `PersistentBitMatrix::{Packed,Columnar}`**: revised — cheaper than
+  first thought, by reusing the same split already used for
+  `fill_matrix`'s own implementation instead of hand-writing a resumable
+  state machine at the matrix level. The base-vector layer
+  (`BitSliceView`, `views.rs`) already separates the two concerns:
+  `fill_batch_sorted` (`views.rs:55-60`, sorted-slot batch lookup) and a
+  genuine per-bit `Iterator` (`BitSliceIter`, `views.rs:94+`) sit side by
+  side there, one level below the matrix. Adding a
+  "positions among `sorted_slots` where the bit is set" iterator at that
+  same vector level is a `filter` over the existing `get()` — no new
+  state machine, since `std::iter::Filter` already *is* one, generated by
+  the compiler:
+  ```rust
+  fn nonzero_among_sorted<'s>(&'s self, sorted_slots: &'s [usize]) -> impl Iterator<Item = usize> + 's {
+      sorted_slots.iter().copied().filter(move |&slot| self.get(slot))
+  }
+  ```
+  The matrix-level `nonzero_iter` then composes these per column with
+  `flat_map` over `0..n_cols` (each column's hits, tagged with `c`,
+  slot mapped back through the sort permutation `fill_batch`/
+  `fill_batch_sorted` already carry) — again a combinator chain, not a
+  hand-rolled struct. Same algorithm, same mmap/sort locality as today's
+  `fill_sub_matrix`; just assembled from `std` iterator adaptors instead
+  of a loop body writing into a buffer, mirroring the vector/matrix split
+  the codebase already uses for `fill_batch_sorted` rather than
+  introducing a new shape.
+- **`Implicit`**: trivial (`slots.iter().map(|&i| (i, 0, 1))`, one column,
+  always present).
+- **On `PersistentCompactIntMatrix` (counts)**: same treatment as
+  `Packed`/`Columnar` — no sparse count format exists yet ("Explicitly
+  deferred" per `traits.rs:9-12`), so no native low-effort case the way
+  `Sparse` has one, but not excluded either: the iterator's `Item` is
+  already `(usize, usize, u32)`, not `bool`, specifically so presence
+  (`0`/`1`) and counts (arbitrary `u32`) share one primitive instead of a
+  bool/u32 split forcing counts out of the design. Ready for a native
+  sparse-count struct later without a signature change.
+
+Would let `query_partition_with`'s stage 2 collapse to one
+`for (i, g, v) in layer.matrix().nonzero_iter(&hit_slots) { on_event(...) }`
+per layer, format-agnostic, each backend's struct deciding the actual
+traversal.
+
+**This also closes the existing dispatch bug for free, by construction**:
+`PersistentBitMatrix::fill_sub_matrix` (`persistent.rs:190-215`, the enum
+wrapper backing `BinaryMatrix`'s trait impl) today does *not* delegate to
+`PersistentSparseBitMatrix::fill_sub_matrix` for `Self::Sparse` — it
+reimplements the same naive per-(column, slot) `fill_row_bool` loop
+instead, bypassing the efficient native method one file over
+(`sparse.rs:249-258`). `obikphylo::siblings::cache::Mat`
+(`cache.rs:138-145`) independently built its own parallel enum wrapper
+that dispatches correctly — evidence this was worked around rather than
+fixed at the source: two hand-written traversals for the same format,
+free to drift apart, and they did. If `fill_sub_matrix` itself is
+rewritten as "drain `nonzero_iter`, scatter into `out[][]`", there is only
+one traversal per format left to get right — the bug class doesn't just
+get fixed once, it stops being possible to reintroduce.