From 2a66f9fcc6a54fd0d4dcc6bc2991d83f55890f27 Mon Sep 17 00:00:00 2001
From: mikessh <mikhail.shugay@gmail.com>
Date: Sun, 18 Oct 2026 23:01:32 +0300
Subject: [PATCH 1/5] Stream native search hits in bounded complete-query
 buffers

Co-Authored-By: Codex <noreply@openai.com>
---
 docs/api.rst                        |  13 +++
 include/seqtree/seqtree.hpp         |   8 ++
 src/_bindings.cpp                   |  33 +++++++
 src/index.cpp                       |  59 ++++++++++++
 tests/python/test_search_consume.py | 134 ++++++++++++++++++++++++++++
 5 files changed, 247 insertions(+)
 create mode 100644 tests/python/test_search_consume.py

diff --git a/docs/api.rst b/docs/api.rst
index f4a940c..af966cf 100644
--- a/docs/api.rst
+++ b/docs/api.rst
@@ -47,6 +47,19 @@ and exact constrained dynamic programming; other profiles use the generic search
 worker. Existing eligible count paths are retained when the bound is redundant.
 Explicit ``seqtrie`` remains unsupported.
 
+For full hits with bounded output, use
+``Index.search_batch_consume(queries, params, consume, threads=0, chunk_size=32, max_hit_bytes=64*1024**2)``.
+It searches the complete population with one native pool. Serialized callbacks,
+which may run on native worker threads,
+receive ``(start, stop, owned_buffer)`` with read-only uint32 columns
+``query_index, ref_id, score, n_subs, n_ins, n_dels``. Every query is complete;
+block order may vary across threads, while hit order within each query is unchanged.
+Each output buffer obeys the explicit byte ceiling. An oversized single query
+raises rather than truncating its hits. Workers also retain the current block's
+native hit scratch. Consume and release buffers to bound total memory; buffers
+explicitly retained by a callback remain owned. Parameters and matrices are
+snapshotted, callbacks are serialized, and workers join before errors propagate.
+
 Floating alignment target reductions
 ------------------------------------
 
diff --git a/include/seqtree/seqtree.hpp b/include/seqtree/seqtree.hpp
index 663a228..a1c6073 100644
--- a/include/seqtree/seqtree.hpp
+++ b/include/seqtree/seqtree.hpp
@@ -2,6 +2,7 @@
 #include "seqtree/types.hpp"
 
 #include <cstdint>
+#include <functional>
 #include <memory>
 #include <optional>
 #include <string>
@@ -13,6 +14,7 @@ namespace seqtree {
 
 struct Trie;     // internal, src/trie.hpp
 struct Scratch;  // internal, src/engines.hpp
+struct CountBuffer;
 
 // Maps characters to small integer codes (<=32) and back. Case-insensitive on
 // input. Invalid characters encode to kInvalid.
@@ -134,6 +136,12 @@ class Index {
     std::vector<std::vector<Hit>> search_batch(const std::vector<std::string>& queries,
                                                const SearchParams&, int threads = 0) const;
 
+    // One pool, serialized owned primitive blocks; each query is emitted completely.
+    // Columns: query index, ref_id, score, n_subs, n_ins, n_dels (all nonnegative).
+    void search_batch_consume(const std::vector<std::string>& queries, const SearchParams&,
+        const std::function<void(size_t, size_t, CountBuffer&&)>& consume,
+        int threads = 0, size_t chunk_size = 32, size_t max_hit_bytes = 64*1024*1024) const;
+
     // Exact edit-count bins 0..max_total_edits; one reusable hit buffer per worker.
     std::vector<std::vector<uint64_t>> edit_histogram_batch(
         const std::vector<std::string>& queries, const SearchParams&,
diff --git a/src/_bindings.cpp b/src/_bindings.cpp
index a30c163..c76fccb 100644
--- a/src/_bindings.cpp
+++ b/src/_bindings.cpp
@@ -12,6 +12,7 @@
 #include <nanobind/stl/vector.h>
 
 #include <algorithm>
+#include <exception>
 #include <optional>
 #include <string>
 #include <vector>
@@ -186,6 +187,29 @@ nb::list py_search_batch(const Index& idx, const std::vector<std::string>& queri
     return out;
 }
 
+void py_search_batch_consume(const Index& idx, const std::vector<std::string>& queries,
+        const PyParams& pp, nb::object consume, int threads, size_t chunk_size, size_t max_hit_bytes) {
+    if (!PyCallable_Check(consume.ptr())) throw nb::value_error("consumer must be callable");
+    // A callback may mutate the original Python parameters; native workers read
+    // only these owned matrix/positional snapshots until the whole pool joins.
+    const PyParams owned = pp;
+    auto mat = make_matrix(owned, idx.alphabet());
+    SearchParams cp = to_cpp(owned, mat ? &*mat : nullptr);
+    std::exception_ptr error;
+    {
+        nb::gil_scoped_release release;
+        try {
+            idx.search_batch_consume(queries, cp,
+                [&](size_t start, size_t stop, CountBuffer&& buffer) {
+                    nb::gil_scoped_acquire acquire;
+                    consume(start, stop, nb::cast(std::move(buffer)));
+                }, threads, chunk_size, max_hit_bytes);
+        } catch (...) { error = std::current_exception(); }
+    }
+    // Restore the caller's GIL before Python exceptions are restored/destroyed.
+    if (error) std::rethrow_exception(error);
+}
+
 std::vector<std::vector<uint64_t>> py_edit_histogram_batch(
         const Index& idx, const std::vector<std::string>& queries,
         const PyParams& pp, int threads, bool exclude_exact, int max_indels) {
@@ -718,6 +742,15 @@ NB_MODULE(_core, m) {
              nb::arg("threads") = 0,
              "Search many queries in parallel (releases the GIL). ``threads=0`` uses all "
              "cores. Returns one hit list per query, in input order.")
+        .def("search_batch_consume", &py_search_batch_consume,
+             nb::arg("queries"), nb::arg("params"), nb::arg("consume"), nb::arg("threads") = 0,
+             nb::arg("chunk_size") = 32, nb::arg("max_hit_bytes") = 64*1024*1024,
+             "One population/pool, serialized consume(start, stop, owned_buffer) callbacks. "
+             "Read-only uint32 rows: query_index, ref_id, score, n_subs, n_ins, n_dels. "
+             "Every query is complete; block order may vary across threads. Output bytes per "
+             "block never exceed max_hit_bytes; an oversized query raises, never truncates. "
+             "Native worker scratch also holds the current block of hits. Consumer-retained "
+             "buffers remain owned; callers control their retention. Parameters are snapshotted.")
         .def("edit_histogram_batch", &py_edit_histogram_batch,
              nb::arg("queries"), nb::arg("params"), nb::arg("threads") = 0,
              nb::arg("exclude_exact") = false, nb::arg("max_indels") = -1,
diff --git a/src/index.cpp b/src/index.cpp
index 559e033..3231bb5 100644
--- a/src/index.cpp
+++ b/src/index.cpp
@@ -5,9 +5,12 @@
 #include "count_unit.hpp"
 
 #include <algorithm>
+#include <atomic>
 #include <cstring>
 #include <fstream>
+#include <iterator>
 #include <limits>
+#include <mutex>
 #include <stdexcept>
 #include <vector>
 
@@ -188,6 +191,62 @@ std::vector<std::vector<Hit>> Index::search_batch(const std::vector<std::string>
     return results;
 }
 
+void Index::search_batch_consume(const std::vector<std::string>& queries,
+        const SearchParams& p,
+        const std::function<void(size_t, size_t, CountBuffer&&)>& consume,
+        int threads, size_t chunk_size, size_t max_hit_bytes) const {
+    constexpr size_t columns = 6;
+    const size_t max_rows = max_hit_bytes / (columns * sizeof(uint32_t));
+    if (!chunk_size || !max_rows || !consume)
+        throw std::invalid_argument("positive chunk_size/max_hit_bytes and a consumer required");
+    if (queries.size() > uint64_t(std::numeric_limits<uint32_t>::max()) + 1)
+        throw std::overflow_error("query indices exceed the uint32 output buffer");
+    // Bound both the primitive output and the per-worker Hit vectors. A single
+    // query may still exceed the explicit output cap; detect that without truncation.
+    chunk_size = std::min(chunk_size, std::max(size_t(1), max_rows / std::max(size_t(1), size_t(size()))));
+    const size_t blocks = queries.size() / chunk_size + (queries.size() % chunk_size != 0);
+    std::atomic<bool> cancelled{false};
+    std::mutex callback_mutex;
+    struct Local { Searcher searcher; std::vector<std::vector<Hit>> hits; };
+    parallel_for(blocks, threads,
+        [&] { return Local{Searcher(*this), std::vector<std::vector<Hit>>(chunk_size)}; },
+        [&](size_t block, Local& worker) {
+            if (cancelled.load()) return;
+            try {
+                const size_t start = block * chunk_size;
+                const size_t stop = start + std::min(chunk_size, queries.size() - start);
+                size_t rows = 0;
+                for (size_t i = start; i < stop; ++i) {
+                    if (cancelled.load()) return;
+                    auto& hits = worker.hits[i-start];
+                    worker.searcher.search_into(queries[i], p, hits);
+                    if (hits.size() > max_rows - rows)
+                        throw std::overflow_error("complete query hits exceed max_hit_bytes; increase the explicit buffer limit");
+                    rows += hits.size();
+                }
+                CountBuffer buffer;
+                buffer.shape = {rows, columns};
+                buffer.data.resize(rows * columns);
+                size_t offset = 0;
+                for (size_t i = start; i < stop; ++i)
+                    for (const auto& hit : worker.hits[i-start]) {
+                        const uint32_t values[] = {uint32_t(i), hit.ref_id, uint32_t(hit.score),
+                            hit.n_subs, hit.n_ins, hit.n_dels};
+                        std::copy(std::begin(values), std::end(values), buffer.data.begin() + offset);
+                        offset += columns;
+                    }
+                std::lock_guard<std::mutex> lock(callback_mutex);
+                if (!cancelled.load()) {
+                    try { consume(start, stop, std::move(buffer)); }
+                    catch (...) { cancelled.store(true); throw; }
+                }
+            } catch (...) {
+                cancelled.store(true);
+                throw;
+            }
+        });
+}
+
 std::vector<std::vector<uint64_t>> Index::edit_histogram_batch(
         const std::vector<std::string>& queries, const SearchParams& p,
         int threads, bool exclude_exact, int max_indels) const {
diff --git a/tests/python/test_search_consume.py b/tests/python/test_search_consume.py
new file mode 100644
index 0000000..c3800dc
--- /dev/null
+++ b/tests/python/test_search_consume.py
@@ -0,0 +1,134 @@
+"""Bounded primitive search output retains the ordinary batch's exact hits."""
+import pytest
+from seqtree import Index, SearchParams
+
+
+@pytest.mark.parametrize('threads', [1, 2])
+@pytest.mark.parametrize('weighted', [False, True])
+def test_consumer_equals_full_batch_and_keeps_owned_complete_buffers(threads, weighted):
+    refs = ['AAAA', 'AAAC', 'AACC', 'CCCC', 'AAA', 'AAAAA', 'AAAA']
+    queries = ['AAAA', 'CCCC', 'GGGG', 'AAAC'] * 9
+    index = Index.build(refs, alphabet='nt')
+    params = SearchParams(max_subs=2, max_ins=1, max_dels=1, max_total_edits=2)
+    if weighted:
+        from seqtree import SubstitutionMatrix
+        params.matrix = SubstitutionMatrix.from_similarity([[2 if i == j else -1 for j in range(4)] for i in range(4)])
+        params.gap_open = params.gap_extend = 3
+    expected = [[(h.ref_id, h.score, h.n_subs, h.n_ins, h.n_dels) for h in hits]
+        for hits in index.search_batch(queries, params, threads)]
+    blocks = []
+    def consume(start, stop, buffer):
+        view = memoryview(buffer)
+        assert view.readonly and view.shape[1] == 6
+        assert view.nbytes <= 24 * len(refs) * 2
+        assert 0 < stop-start <= 2
+        blocks.append((start, stop, buffer))  # owns data after callback return
+    index.search_batch_consume(queries, params, consume, threads=threads,
+        chunk_size=2, max_hit_bytes=24*len(refs)*2)
+    observed = [[] for _ in queries]
+    covered = []
+    for start, stop, buffer in sorted(blocks, key=lambda b: b[0]):
+        covered.extend(range(start, stop))
+        for row in memoryview(buffer).tolist():
+            assert start <= row[0] < stop
+            observed[row[0]].append(tuple(row[1:]))
+    assert covered == list(range(len(queries)))
+    assert observed == expected
+
+
+@pytest.mark.parametrize('base_exception', [False, True])
+def test_consumer_overflow_and_callback_errors_join_workers_and_preserve_index(base_exception):
+    index = Index.build(['AAAA'] * 8, alphabet='nt')
+    params = SearchParams(max_subs=1)
+    with pytest.raises(OverflowError, match='max_hit_bytes'):
+        index.search_batch_consume(['AAAA']*20, params, lambda *_: None,
+            threads=2, max_hit_bytes=24)
+    class CallbackError(BaseException if base_exception else Exception):
+        pass
+    calls = []
+    def fail(*args):
+        calls.append(1)
+        raise CallbackError('consumer stopped')
+    with pytest.raises(CallbackError, match='consumer stopped'):
+        index.search_batch_consume(['AAAA']*20, params, fail, threads=2, chunk_size=1)
+    assert calls == [1]
+    assert len(index.search('AAAA', params)) == 8
+
+
+def test_consumer_snapshots_params_and_bounds_dense_output():
+    import threading
+    from seqtree import PositionalMatrix, SubstitutionMatrix
+    index = Index.build(['AAAA']*100, alphabet='nt')
+    params = SearchParams()
+    params.pos_matrix = PositionalMatrix.from_weights(SubstitutionMatrix.unit(4), [1, 1, 1, 0])
+    calls = []
+    workers = set()
+    def consume(start, stop, buffer):
+        assert stop-start == 1 and memoryview(buffer).nbytes == 2400
+        params.max_subs = 0  # mutable Python parameters cannot change this batch
+        params.pos_matrix = PositionalMatrix.from_weights(SubstitutionMatrix.unit(4), [1, 1, 1, 1])
+        assert memoryview(buffer).shape == (100, 6)
+        calls.append(start)
+        workers.add(threading.get_native_id())
+    index.search_batch_consume(['AAAC']*1000, params, consume, threads=2,
+        chunk_size=32, max_hit_bytes=2400)
+    assert sorted(calls) == list(range(1000))
+    assert len(workers) <= 2  # one native pool, not one new pool per output block
+    assert index.search('AAAC', params) == []
+
+
+def test_consumer_emits_empty_queries_and_native_errors_leave_index_usable():
+    index = Index.build(['AAAA'], alphabet='nt')
+    blocks = []
+    index.search_batch_consume(['CCCC'], SearchParams(),
+        lambda start, stop, buffer: blocks.append((start, stop, memoryview(buffer).shape)))
+    assert blocks == [(0, 1, (0, 6))]
+    with pytest.raises(ValueError, match='invalid symbol'):
+        index.search_batch_consume(['AAAA', '!!!!'], SearchParams(), lambda *_: None, threads=2)
+    assert len(index.search('AAAA', SearchParams())) == 1
+
+
+def test_callbacks_remain_serial_when_python_releases_the_gil():
+    import time
+    active = [0]
+    calls = []
+    def consume(start, stop, buffer):
+        active[0] += 1
+        assert active[0] == 1
+        time.sleep(.001)
+        calls.append(start)
+        active[0] -= 1
+    Index.build(['AAAA'], alphabet='nt').search_batch_consume(
+        ['AAAA']*32, SearchParams(), consume, threads=2, chunk_size=1)
+    assert sorted(calls) == list(range(32)) and active == [0]
+
+
+def test_consumer_peak_memory_does_not_accumulate_dense_cohort_hits():
+    import json
+    import os
+    import subprocess
+    import sys
+    if sys.platform == 'win32':
+        pytest.skip('resource.ru_maxrss is unavailable on Windows; buffer/pool bounds are tested above')
+    code = '''import json,resource,sys
+from seqtree import Index,SearchParams
+index=Index.build(['AAAA']*1000,alphabet='nt')
+count=[0]
+def consume(start,stop,buffer):
+    count[0]+=memoryview(buffer).shape[0]
+index.search_batch_consume(['AAAA']*int(sys.argv[1]),SearchParams(),consume,
+    threads=2,chunk_size=32,max_hit_bytes=96000)
+peak=resource.getrusage(resource.RUSAGE_SELF).ru_maxrss
+print(json.dumps({'hits':count[0],'peak':peak*(1 if sys.platform=='darwin' else 1024)}))
+'''
+    small, large = [json.loads(subprocess.check_output([sys.executable, '-c', code, str(n)],
+        text=True, env=os.environ, timeout=30)) for n in (200, 20000)]
+    assert small['hits'] == 200000 and large['hits'] == 20000000
+    assert large['peak'] <= small['peak'] + 32*1024**2
+
+
+@pytest.mark.parametrize('kwargs', [{'chunk_size':0}, {'max_hit_bytes':0}])
+def test_consumer_rejects_invalid_output_limits(kwargs):
+    with pytest.raises(ValueError):
+        Index.build(['AAAA'], alphabet='nt').search_batch_consume(
+            ['AAAA'], SearchParams(), lambda *_: None, **kwargs)

From ef2e5f6141ed28967ad2ed4d0e49168452c56e08 Mon Sep 17 00:00:00 2001
From: mikessh <mikhail.shugay@gmail.com>
Date: Sun, 18 Oct 2026 23:35:37 +0300
Subject: [PATCH 2/5] Coordinate bounded two-index search in one native pool

Co-Authored-By: Codex <noreply@openai.com>
---
 docs/api.rst                        |  9 ++++
 include/seqtree/seqtree.hpp         |  7 +++
 src/_bindings.cpp                   | 32 +++++++++++
 src/index.cpp                       | 67 +++++++++++++++++++++++
 tests/python/test_paired_consume.py | 83 +++++++++++++++++++++++++++++
 5 files changed, 198 insertions(+)
 create mode 100644 tests/python/test_paired_consume.py

diff --git a/docs/api.rst b/docs/api.rst
index af966cf..ccdcded 100644
--- a/docs/api.rst
+++ b/docs/api.rst
@@ -60,6 +60,15 @@ native hit scratch. Consume and release buffers to bound total memory; buffers
 explicitly retained by a callback remain owned. Parameters and matrices are
 snapshotted, callbacks are serialized, and workers join before errors propagate.
 
+For coordinated searches over two same-alphabet indexes, use
+``Index.search_paired_batch_consume(other, queries, other_queries, params, consume, threads=0, chunk_size=32, max_hit_bytes=64*1024**2)``.
+Both equal-length query populations run in one pool with two reusable Searchers
+per worker. Serialized callbacks receive ``(start, stop, owned_buffer,
+other_owned_buffer)`` for the same complete query block. The byte ceiling applies
+to the sum of both output buffers; ownership, snapshots and error cleanup follow
+the single-index consumer contract. This API provides no pairing or score fusion;
+the caller owns those interpretations.
+
 Floating alignment target reductions
 ------------------------------------
 
diff --git a/include/seqtree/seqtree.hpp b/include/seqtree/seqtree.hpp
index a1c6073..7a3c2ca 100644
--- a/include/seqtree/seqtree.hpp
+++ b/include/seqtree/seqtree.hpp
@@ -142,6 +142,13 @@ class Index {
         const std::function<void(size_t, size_t, CountBuffer&&)>& consume,
         int threads = 0, size_t chunk_size = 32, size_t max_hit_bytes = 64*1024*1024) const;
 
+    // Coordinate two indexes in one pool; both owned buffers cover the same queries.
+    void search_paired_batch_consume(const Index& other,
+        const std::vector<std::string>& queries, const std::vector<std::string>& other_queries,
+        const SearchParams&,
+        const std::function<void(size_t, size_t, CountBuffer&&, CountBuffer&&)>& consume,
+        int threads = 0, size_t chunk_size = 32, size_t max_hit_bytes = 64*1024*1024) const;
+
     // Exact edit-count bins 0..max_total_edits; one reusable hit buffer per worker.
     std::vector<std::vector<uint64_t>> edit_histogram_batch(
         const std::vector<std::string>& queries, const SearchParams&,
diff --git a/src/_bindings.cpp b/src/_bindings.cpp
index c76fccb..d9070ff 100644
--- a/src/_bindings.cpp
+++ b/src/_bindings.cpp
@@ -210,6 +210,30 @@ void py_search_batch_consume(const Index& idx, const std::vector<std::string>& q
     if (error) std::rethrow_exception(error);
 }
 
+void py_search_paired_batch_consume(const Index& idx, const Index& other,
+        const std::vector<std::string>& queries, const std::vector<std::string>& other_queries,
+        const PyParams& pp, nb::object consume, int threads, size_t chunk_size, size_t max_hit_bytes) {
+    if (!PyCallable_Check(consume.ptr())) throw nb::value_error("consumer must be callable");
+    // A callback may mutate the original Python parameters; native workers read
+    // only these owned matrix/positional snapshots until the whole pool joins.
+    const PyParams owned = pp;
+    auto mat = make_matrix(owned, idx.alphabet());
+    SearchParams cp = to_cpp(owned, mat ? &*mat : nullptr);
+    std::exception_ptr error;
+    {
+        nb::gil_scoped_release release;
+        try {
+            idx.search_paired_batch_consume(other, queries, other_queries, cp,
+                [&](size_t start, size_t stop, CountBuffer&& buffer, CountBuffer&& other_buffer) {
+                    nb::gil_scoped_acquire acquire;
+                    consume(start, stop, nb::cast(std::move(buffer)), nb::cast(std::move(other_buffer)));
+                }, threads, chunk_size, max_hit_bytes);
+        } catch (...) { error = std::current_exception(); }
+    }
+    // Restore the caller's GIL before Python exceptions are restored/destroyed.
+    if (error) std::rethrow_exception(error);
+}
+
 std::vector<std::vector<uint64_t>> py_edit_histogram_batch(
         const Index& idx, const std::vector<std::string>& queries,
         const PyParams& pp, int threads, bool exclude_exact, int max_indels) {
@@ -751,6 +775,14 @@ NB_MODULE(_core, m) {
              "block never exceed max_hit_bytes; an oversized query raises, never truncates. "
              "Native worker scratch also holds the current block of hits. Consumer-retained "
              "buffers remain owned; callers control their retention. Parameters are snapshotted.")
+        .def("search_paired_batch_consume", &py_search_paired_batch_consume,
+             nb::arg("other"), nb::arg("queries"), nb::arg("other_queries"),
+             nb::arg("params"), nb::arg("consume"), nb::arg("threads") = 0,
+             nb::arg("chunk_size") = 32, nb::arg("max_hit_bytes") = 64*1024*1024,
+             "One pool with two Searchers; consume(start, stop, owned_buffer, other_owned_buffer). "
+             "Both six-column buffers cover the same complete query block. Combined output "
+             "bytes are bounded by max_hit_bytes; block order may vary. Same-alphabet indexes "
+             "and equal-length query populations required. Parameters are snapshotted.")
         .def("edit_histogram_batch", &py_edit_histogram_batch,
              nb::arg("queries"), nb::arg("params"), nb::arg("threads") = 0,
              nb::arg("exclude_exact") = false, nb::arg("max_indels") = -1,
diff --git a/src/index.cpp b/src/index.cpp
index 3231bb5..e906258 100644
--- a/src/index.cpp
+++ b/src/index.cpp
@@ -247,6 +247,73 @@ void Index::search_batch_consume(const std::vector<std::string>& queries,
         });
 }
 
+void Index::search_paired_batch_consume(const Index& other,
+        const std::vector<std::string>& queries, const std::vector<std::string>& other_queries,
+        const SearchParams& p,
+        const std::function<void(size_t, size_t, CountBuffer&&, CountBuffer&&)>& consume,
+        int threads, size_t chunk_size, size_t max_hit_bytes) const {
+    if (queries.size() != other_queries.size())
+        throw std::invalid_argument("paired query populations must have equal lengths");
+    if (alphabet() != other.alphabet())
+        throw std::invalid_argument("paired indexes must use the same alphabet");
+    constexpr size_t columns = 6;
+    const size_t max_rows = max_hit_bytes / (columns * sizeof(uint32_t));
+    if (!chunk_size || !max_rows || !consume)
+        throw std::invalid_argument("positive chunk_size/max_hit_bytes and a consumer required");
+    if (queries.size() > uint64_t(std::numeric_limits<uint32_t>::max()) + 1)
+        throw std::overflow_error("query indices exceed the uint32 output buffer");
+    // Bound both the primitive output and the per-worker Hit vectors. A single
+    // query may still exceed the explicit output cap; detect that without truncation.
+    chunk_size = std::min(chunk_size, std::max(size_t(1), max_rows / std::max(size_t(1), size_t(size()) + size_t(other.size()))));
+    const size_t blocks = queries.size() / chunk_size + (queries.size() % chunk_size != 0);
+    std::atomic<bool> cancelled{false};
+    std::mutex callback_mutex;
+    struct Local { Searcher searcher, other_searcher; std::vector<std::vector<Hit>> hits, other_hits; };
+    parallel_for(blocks, threads,
+        [&] { return Local{Searcher(*this), Searcher(other), std::vector<std::vector<Hit>>(chunk_size), std::vector<std::vector<Hit>>(chunk_size)}; },
+        [&](size_t block, Local& worker) {
+            if (cancelled.load()) return;
+            try {
+                const size_t start = block * chunk_size;
+                const size_t stop = start + std::min(chunk_size, queries.size() - start);
+                size_t rows = 0;
+                for (size_t i = start; i < stop; ++i) {
+                    if (cancelled.load()) return;
+                    auto& hits = worker.hits[i-start];
+                    worker.searcher.search_into(queries[i], p, hits);
+                    auto& other_hits = worker.other_hits[i-start];
+                    worker.other_searcher.search_into(other_queries[i], p, other_hits);
+                    if (hits.size() > max_rows - rows || other_hits.size() > max_rows - rows - hits.size())
+                        throw std::overflow_error("complete query hits exceed max_hit_bytes; increase the explicit buffer limit");
+                    rows += hits.size() + other_hits.size();
+                }
+                auto pack = [&](const auto& groups) {
+                    CountBuffer buffer;
+                    size_t count = 0;
+                    for (size_t i = start; i < stop; ++i) count += groups[i-start].size();
+                    buffer.shape = {count, columns};
+                    buffer.data.reserve(count * columns);
+                    for (size_t i = start; i < stop; ++i)
+                        for (const auto& hit : groups[i-start]) {
+                            const uint32_t values[] = {uint32_t(i), hit.ref_id, uint32_t(hit.score),
+                                hit.n_subs, hit.n_ins, hit.n_dels};
+                            buffer.data.insert(buffer.data.end(), std::begin(values), std::end(values));
+                        }
+                    return buffer;
+                };
+                auto buffer = pack(worker.hits), other_buffer = pack(worker.other_hits);
+                std::lock_guard<std::mutex> lock(callback_mutex);
+                if (!cancelled.load()) {
+                    try { consume(start, stop, std::move(buffer), std::move(other_buffer)); }
+                    catch (...) { cancelled.store(true); throw; }
+                }
+            } catch (...) {
+                cancelled.store(true);
+                throw;
+            }
+        });
+}
+
 std::vector<std::vector<uint64_t>> Index::edit_histogram_batch(
         const std::vector<std::string>& queries, const SearchParams& p,
         int threads, bool exclude_exact, int max_indels) const {
diff --git a/tests/python/test_paired_consume.py b/tests/python/test_paired_consume.py
new file mode 100644
index 0000000..dba47f7
--- /dev/null
+++ b/tests/python/test_paired_consume.py
@@ -0,0 +1,83 @@
+"""Coordinated full-population search preserves each chain's ordinary hits."""
+import pytest
+from seqtree import Index, SearchParams, PositionalMatrix, SubstitutionMatrix
+
+@pytest.mark.parametrize('threads',[1,2])
+@pytest.mark.parametrize('weighted',[False,True])
+def test_paired_buffers_equal_legacy(threads,weighted):
+    a=Index.build(['AAAA','AAAC','AAAA','CCCC'],alphabet='nt')
+    b=Index.build(['CCCC','CCCA','AAAA'],alphabet='nt')
+    qa=['AAAA','GGGG','AAAC']*25;qb=['CCCC','GGGG','AAAA']*25
+    p=SearchParams(max_subs=2,max_ins=1,max_dels=1,max_total_edits=2)
+    if weighted:
+        p.matrix=SubstitutionMatrix.unit(4);p.gap_open=p.gap_extend=2
+    expected=[[[tuple(h) for h in row] for row in index.search_batch(q,p,threads)] for index,q in [(a,qa),(b,qb)]]
+    blocks=[]
+    a.search_paired_batch_consume(b,qa,qb,p,lambda *args:blocks.append(args),threads=threads,chunk_size=2,max_hit_bytes=24*7*2)
+    actual=[[[] for _ in qa] for _ in range(2)];covered=[]
+    for start,stop,*buffers in sorted(blocks):
+        covered.extend(range(start,stop))
+        assert 0<stop-start<=2
+        assert sum(memoryview(buf).nbytes for buf in buffers)<=24*7*2
+        for chain,buf in enumerate(buffers):
+            assert memoryview(buf).readonly
+            for qi,*hit in memoryview(buf).tolist():actual[chain][qi].append(tuple(hit))
+    assert covered==list(range(len(qa))) and actual==expected
+
+@pytest.mark.parametrize('base',[False,True])
+def test_paired_failure_cleanup_and_snapshot(base):
+    a=Index.build(['AAAA']*10,alphabet='nt');b=Index.build(['AAAC']*10,alphabet='nt')
+    p=SearchParams(max_subs=1)
+    p.pos_matrix=PositionalMatrix.from_weights(SubstitutionMatrix.unit(4),[1,1,1,0])
+    calls=[]
+    def callback(start,stop,x,y):
+        p.max_subs=0;p.pos_matrix=PositionalMatrix.from_weights(SubstitutionMatrix.unit(4),[1,1,1,1])
+        assert memoryview(x).shape[0]==10 and memoryview(y).shape[0]==10
+        calls.append(start)
+    a.search_paired_batch_consume(b,['AAAA']*20,['AAAA']*20,p,callback,threads=2,chunk_size=1)
+    assert sorted(calls)==list(range(20))
+    class Stop(BaseException if base else Exception):pass
+    calls.clear()
+    def fail(*_):calls.append(1);raise Stop('stop')
+    with pytest.raises(Stop,match='stop'):
+        a.search_paired_batch_consume(b,['AAAA']*20,['AAAA']*20,p,fail,threads=2,chunk_size=1)
+    assert calls==[1] and len(a.search('AAAA',p))==10
+    with pytest.raises(OverflowError,match='max_hit_bytes'):
+        a.search_paired_batch_consume(b,['AAAA'],['AAAA'],p,lambda *_:None,max_hit_bytes=24)
+    with pytest.raises(ValueError,match='equal lengths'):
+        a.search_paired_batch_consume(b,['AAAA'],[],p,lambda *_:None)
+    with pytest.raises(ValueError,match='same alphabet'):
+        a.search_paired_batch_consume(Index.build(['AAAA']),[],[],p,lambda *_:None)
+
+
+def test_paired_callbacks_serial_and_invalid_query_cleanup():
+    import threading,time
+    a=Index.build(['AAAA']*3,alphabet='nt');b=Index.build(['CCCC']*4,alphabet='nt')
+    p=SearchParams(max_subs=4);active=0;workers=set()
+    def consume(start,stop,x,y):
+        nonlocal active
+        active+=1;assert active==1
+        workers.add(threading.get_native_id())
+        time.sleep(.001)
+        assert active==1;active-=1
+    a.search_paired_batch_consume(b,['AAAA']*40,['CCCC']*40,p,consume,threads=2,chunk_size=1)
+    assert len(workers)<=2
+    with pytest.raises(ValueError):
+        a.search_paired_batch_consume(b,['AAAA']*20,['***']*20,p,consume,threads=2,chunk_size=1)
+    assert len(b.search('CCCC',p))==4
+
+
+def test_paired_peak_memory_does_not_retain_first_chain_population():
+    import json,os,subprocess,sys
+    if sys.platform=='win32':
+        pytest.skip('resource.ru_maxrss unavailable; primitive ownership and caps tested above')
+    code='''import json,resource,sys
+from seqtree import Index,SearchParams
+a=Index.build(['AAAA']*1000,alphabet='nt');b=Index.build(['CCCC']*1000,alphabet='nt');count=[0]
+def consume(start,stop,x,y):count[0]+=memoryview(x).shape[0]+memoryview(y).shape[0]
+a.search_paired_batch_consume(b,['AAAA']*int(sys.argv[1]),['CCCC']*int(sys.argv[1]),SearchParams(),consume,threads=2,chunk_size=32,max_hit_bytes=96000)
+print(json.dumps({'hits':count[0],'peak':resource.getrusage(resource.RUSAGE_SELF).ru_maxrss*(1 if sys.platform=='darwin' else 1024)}))'''
+    small,large=[json.loads(subprocess.check_output([sys.executable,'-c',code,str(n)],
+        env=os.environ,text=True,timeout=30)) for n in (200,20000)]
+    assert small['hits']==400000 and large['hits']==40000000
+    assert large['peak']<=small['peak']+32*1024**2

From 59e2e28998e2c503d54d04709cb93e42b0aac742 Mon Sep 17 00:00:00 2001
From: mikessh <mikhail.shugay@gmail.com>
Date: Mon, 19 Oct 2026 00:00:14 +0300
Subject: [PATCH 3/5] Allow explicitly empty axes in coordinated search

Co-Authored-By: Codex <noreply@openai.com>
---
 docs/api.rst                        |  6 ++++--
 include/seqtree/seqtree.hpp         |  2 +-
 src/_bindings.cpp                   |  9 +++++----
 src/index.cpp                       |  8 +++++---
 tests/python/test_paired_consume.py | 15 +++++++++++++++
 5 files changed, 30 insertions(+), 10 deletions(-)

diff --git a/docs/api.rst b/docs/api.rst
index ccdcded..18bcd95 100644
--- a/docs/api.rst
+++ b/docs/api.rst
@@ -61,13 +61,15 @@ explicitly retained by a callback remain owned. Parameters and matrices are
 snapshotted, callbacks are serialized, and workers join before errors propagate.
 
 For coordinated searches over two same-alphabet indexes, use
-``Index.search_paired_batch_consume(other, queries, other_queries, params, consume, threads=0, chunk_size=32, max_hit_bytes=64*1024**2)``.
+``Index.search_paired_batch_consume(other, queries, other_queries, params, consume, threads=0, chunk_size=32, max_hit_bytes=64*1024**2, skip_empty=False)``.
 Both equal-length query populations run in one pool with two reusable Searchers
 per worker. Serialized callbacks receive ``(start, stop, owned_buffer,
 other_owned_buffer)`` for the same complete query block. The byte ceiling applies
 to the sum of both output buffers; ownership, snapshots and error cleanup follow
 the single-index consumer contract. This API provides no pairing or score fusion;
-the caller owns those interpretations.
+the caller owns those interpretations. ``skip_empty=True`` emits no hits for an
+empty query axis while retaining searches on its partner; the default preserves
+ordinary empty-string search semantics.
 
 Floating alignment target reductions
 ------------------------------------
diff --git a/include/seqtree/seqtree.hpp b/include/seqtree/seqtree.hpp
index 7a3c2ca..56598dc 100644
--- a/include/seqtree/seqtree.hpp
+++ b/include/seqtree/seqtree.hpp
@@ -147,7 +147,7 @@ class Index {
         const std::vector<std::string>& queries, const std::vector<std::string>& other_queries,
         const SearchParams&,
         const std::function<void(size_t, size_t, CountBuffer&&, CountBuffer&&)>& consume,
-        int threads = 0, size_t chunk_size = 32, size_t max_hit_bytes = 64*1024*1024) const;
+        int threads = 0, size_t chunk_size = 32, size_t max_hit_bytes = 64*1024*1024, bool skip_empty = false) const;
 
     // Exact edit-count bins 0..max_total_edits; one reusable hit buffer per worker.
     std::vector<std::vector<uint64_t>> edit_histogram_batch(
diff --git a/src/_bindings.cpp b/src/_bindings.cpp
index d9070ff..8d23730 100644
--- a/src/_bindings.cpp
+++ b/src/_bindings.cpp
@@ -212,7 +212,7 @@ void py_search_batch_consume(const Index& idx, const std::vector<std::string>& q
 
 void py_search_paired_batch_consume(const Index& idx, const Index& other,
         const std::vector<std::string>& queries, const std::vector<std::string>& other_queries,
-        const PyParams& pp, nb::object consume, int threads, size_t chunk_size, size_t max_hit_bytes) {
+        const PyParams& pp, nb::object consume, int threads, size_t chunk_size, size_t max_hit_bytes, bool skip_empty) {
     if (!PyCallable_Check(consume.ptr())) throw nb::value_error("consumer must be callable");
     // A callback may mutate the original Python parameters; native workers read
     // only these owned matrix/positional snapshots until the whole pool joins.
@@ -227,7 +227,7 @@ void py_search_paired_batch_consume(const Index& idx, const Index& other,
                 [&](size_t start, size_t stop, CountBuffer&& buffer, CountBuffer&& other_buffer) {
                     nb::gil_scoped_acquire acquire;
                     consume(start, stop, nb::cast(std::move(buffer)), nb::cast(std::move(other_buffer)));
-                }, threads, chunk_size, max_hit_bytes);
+                }, threads, chunk_size, max_hit_bytes, skip_empty);
         } catch (...) { error = std::current_exception(); }
     }
     // Restore the caller's GIL before Python exceptions are restored/destroyed.
@@ -778,11 +778,12 @@ NB_MODULE(_core, m) {
         .def("search_paired_batch_consume", &py_search_paired_batch_consume,
              nb::arg("other"), nb::arg("queries"), nb::arg("other_queries"),
              nb::arg("params"), nb::arg("consume"), nb::arg("threads") = 0,
-             nb::arg("chunk_size") = 32, nb::arg("max_hit_bytes") = 64*1024*1024,
+             nb::arg("chunk_size") = 32, nb::arg("max_hit_bytes") = 64*1024*1024, nb::arg("skip_empty") = false,
              "One pool with two Searchers; consume(start, stop, owned_buffer, other_owned_buffer). "
              "Both six-column buffers cover the same complete query block. Combined output "
              "bytes are bounded by max_hit_bytes; block order may vary. Same-alphabet indexes "
-             "and equal-length query populations required. Parameters are snapshotted.")
+             "and equal-length query populations required. Parameters are snapshotted. "
+             "skip_empty=True emits no hits for empty query axes; other axes remain searched.")
         .def("edit_histogram_batch", &py_edit_histogram_batch,
              nb::arg("queries"), nb::arg("params"), nb::arg("threads") = 0,
              nb::arg("exclude_exact") = false, nb::arg("max_indels") = -1,
diff --git a/src/index.cpp b/src/index.cpp
index e906258..7d2f2bc 100644
--- a/src/index.cpp
+++ b/src/index.cpp
@@ -251,7 +251,7 @@ void Index::search_paired_batch_consume(const Index& other,
         const std::vector<std::string>& queries, const std::vector<std::string>& other_queries,
         const SearchParams& p,
         const std::function<void(size_t, size_t, CountBuffer&&, CountBuffer&&)>& consume,
-        int threads, size_t chunk_size, size_t max_hit_bytes) const {
+        int threads, size_t chunk_size, size_t max_hit_bytes, bool skip_empty) const {
     if (queries.size() != other_queries.size())
         throw std::invalid_argument("paired query populations must have equal lengths");
     if (alphabet() != other.alphabet())
@@ -280,9 +280,11 @@ void Index::search_paired_batch_consume(const Index& other,
                 for (size_t i = start; i < stop; ++i) {
                     if (cancelled.load()) return;
                     auto& hits = worker.hits[i-start];
-                    worker.searcher.search_into(queries[i], p, hits);
+                    if (skip_empty && queries[i].empty()) hits.clear();
+                    else worker.searcher.search_into(queries[i], p, hits);
                     auto& other_hits = worker.other_hits[i-start];
-                    worker.other_searcher.search_into(other_queries[i], p, other_hits);
+                    if (skip_empty && other_queries[i].empty()) other_hits.clear();
+                    else worker.other_searcher.search_into(other_queries[i], p, other_hits);
                     if (hits.size() > max_rows - rows || other_hits.size() > max_rows - rows - hits.size())
                         throw std::overflow_error("complete query hits exceed max_hit_bytes; increase the explicit buffer limit");
                     rows += hits.size() + other_hits.size();
diff --git a/tests/python/test_paired_consume.py b/tests/python/test_paired_consume.py
index dba47f7..3e504c9 100644
--- a/tests/python/test_paired_consume.py
+++ b/tests/python/test_paired_consume.py
@@ -81,3 +81,18 @@ def consume(start,stop,x,y):count[0]+=memoryview(x).shape[0]+memoryview(y).shape
         env=os.environ,text=True,timeout=30)) for n in (200,20000)]
     assert small['hits']==400000 and large['hits']==40000000
     assert large['peak']<=small['peak']+32*1024**2
+
+
+def test_skip_empty_preserves_other_axis_and_original_indices():
+    a=Index.build(['A'],alphabet='nt');b=Index.build(['C'],alphabet='nt')
+    p=SearchParams(max_ins=1,max_total_edits=1)
+    qa=['A','','A',''];qb=['','C','C',''];blocks=[]
+    a.search_paired_batch_consume(b,qa,qb,p,lambda *x:blocks.append(x),threads=2,chunk_size=1,skip_empty=True)
+    observed=[[],[]]
+    for start,stop,x,y in blocks:
+        for axis,buffer in enumerate((x,y)):
+            observed[axis].extend(row[0] for row in memoryview(buffer).tolist())
+    assert sorted(observed[0])==[0,2] and sorted(observed[1])==[1,2]
+    legacy=[]
+    a.search_paired_batch_consume(b,[''],[''],p,lambda *x:legacy.append(x))
+    assert memoryview(legacy[0][2]).shape[0]==memoryview(legacy[0][3]).shape[0]==1

From a999806a735e4570f3d190a85dac20b9abd7a856 Mon Sep 17 00:00:00 2001
From: mikessh <mikhail.shugay@gmail.com>
Date: Mon, 19 Oct 2026 12:29:12 +0300
Subject: [PATCH 4/5] Document bounded hit consumers and coordinated empty axes

Co-Authored-By: Codex <noreply@openai.com>
---
 CHANGELOG.md | 9 +++++++++
 1 file changed, 9 insertions(+)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index 0e93a5e..d35c4ef 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -3,6 +3,15 @@
 All notable changes to `seqtree`. Dates are release dates. From 1.0.0 the project follows semantic
 versioning: breaking changes need a **major** bump.
 
+## [Unreleased]
+
+### Added
+
+- Bounded complete-query hit consumers for single and coordinated two-index
+  searches, with owned numeric buffers, explicit byte budgets and worker cleanup.
+- Optional empty query axes in coordinated search retain original query offsets
+  and supported marginal evidence without creating a query for the missing axis.
+
 ## [1.2.0] — 2026-10-17
 
 ### Fixed

From 43f6d941e4e90517107ef0015acf9a9cad67f464 Mon Sep 17 00:00:00 2001
From: mikessh <mikhail.shugay@gmail.com>
Date: Mon, 19 Oct 2026 12:42:56 +0300
Subject: [PATCH 5/5] Prepare seqtree 1.3.0 bounded consumer release

Co-Authored-By: Codex <noreply@openai.com>
---
 CHANGELOG.md   | 2 +-
 pyproject.toml | 2 +-
 2 files changed, 2 insertions(+), 2 deletions(-)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index d35c4ef..e6b740b 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -3,7 +3,7 @@
 All notable changes to `seqtree`. Dates are release dates. From 1.0.0 the project follows semantic
 versioning: breaking changes need a **major** bump.
 
-## [Unreleased]
+## [1.3.0] — 2026-10-19
 
 ### Added
 
diff --git a/pyproject.toml b/pyproject.toml
index b13195d..9a00cc6 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -4,7 +4,7 @@ build-backend = "scikit_build_core.build"
 
 [project]
 name = "seqtree"
-version = "1.2.0"
+version = "1.3.0"
 description = "Fast fuzzy search over biological sequences (C++ core, Python bindings)"
 readme = "README.md"
 requires-python = ">=3.10"