Fix MESSI top-k returning duplicate series

lib/algos/Messi.cpp

@@ -137,7 +137,7 @@ namespace daisy
                         pos = *node->buffer->full_position_buffer[i] / index->settings->timeseries_size;
                     }
                     pthread_rwlock_wrlock(lock_queue);
-                    pqueue_bsf_insert(pq_bsf, dist, static_cast<long int>(pos), node);
+                    pqueue_bsf_insert_no_dup_position(pq_bsf, dist, static_cast<long int>(pos), node);
                     pthread_rwlock_unlock(lock_queue);
                 }
             }
@@ -156,7 +156,7 @@ namespace daisy
                         pos = *node->buffer->tmp_full_position_buffer[i] / index->settings->timeseries_size;
                     }
                     pthread_rwlock_wrlock(lock_queue);
-                    pqueue_bsf_insert(pq_bsf, dist, static_cast<long int>(pos), node);
+                    pqueue_bsf_insert_no_dup_position(pq_bsf, dist, static_cast<long int>(pos), node);
                     pthread_rwlock_unlock(lock_queue);
                 }
             }
@@ -178,7 +178,7 @@ namespace daisy
                     {
                         pthread_rwlock_wrlock(lock_queue);
 
-                        pqueue_bsf_insert(pq_bsf, dist, *node->buffer->partial_position_buffer[i] / index->settings->timeseries_size, node);
+                        pqueue_bsf_insert_no_dup_position(pq_bsf, dist, *node->buffer->partial_position_buffer[i] / index->settings->timeseries_size, node);
 
                         pthread_rwlock_unlock(lock_queue);
                     }
@@ -229,7 +229,7 @@ namespace daisy
                 if (dist <= pq_bsf->knn[pq_bsf->k - 1])
                 {
                     pthread_rwlock_wrlock(lock_queue);
-                    pqueue_bsf_insert(pq_bsf, dist, *node->buffer->partial_position_buffer[i] / index->settings->timeseries_size, node);
+                    pqueue_bsf_insert_no_dup_position(pq_bsf, dist, *node->buffer->partial_position_buffer[i] / index->settings->timeseries_size, node);
                     pthread_rwlock_unlock(lock_queue);
                 }
             }

tests/CMakeLists.txt

@@ -370,6 +370,29 @@ gtest_discover_tests(test_LbBruteforce_Streaming WORKING_DIRECTORY ${CMAKE_SOURC
      message(STATUS "Tests discovered for test_Messi_L2Square.")
  endif()
 
+# ////// MESSI top-k distinct results //////
+add_executable(
+    test_Messi_TopK_Distinct
+    test_Messi_TopK_Distinct.cpp
+)
+target_link_libraries(
+    test_Messi_TopK_Distinct
+    PRIVATE
+    GTest::gtest_main
+    dino_lib
+    commons_lib
+    stdc++fs
+)
+target_include_directories(test_Messi_TopK_Distinct
+    PRIVATE
+    ${CMAKE_CURRENT_SOURCE_DIR}/../lib
+    ${CMAKE_CURRENT_SOURCE_DIR}/../commons
+)
+gtest_discover_tests(
+    test_Messi_TopK_Distinct
+    WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}
+)
+
 # ////// MESSI streaming //////
 add_executable(
     test_Messi_Streaming

tests/test_Messi_TopK_Distinct.cpp

@@ -0,0 +1,79 @@
+#include <gtest/gtest.h>
+
+#include <algos/Bruteforce.hpp>
+#include <algos/Messi.hpp>
+
+#include <algorithm>
+#include <random>
+#include <set>
+#include <vector>
+
+// Regression test: MESSI's approximate phase scans one leaf to seed the
+// best-so-far queue, then the parallel phase scans every leaf again. Without
+// de-duplication by position the seeded leaf's series were inserted twice, so
+// the k-NN list contained the same series several times and displaced real
+// neighbours.
+namespace
+{
+    constexpr daisy::idx_t kRows = 20000;
+    constexpr daisy::idx_t kDim = 128;
+    constexpr daisy::idx_t kQueries = 50;
+    constexpr daisy::idx_t kK = 10;
+
+    // 50 Gaussian clusters with isotropic noise: dense neighbourhoods, so the
+    // leaf scanned by the approximate phase holds members of the final top-k.
+    std::vector<float> clusteredData(daisy::idx_t rows, unsigned seed)
+    {
+        std::mt19937 centre_rng(42);
+        std::normal_distribution<float> normal;
+        std::vector<float> centres(50 * kDim);
+        for (auto &v : centres)
+            v = 6.0F * normal(centre_rng);
+        std::mt19937 rng(seed);
+        std::uniform_int_distribution<int> pick(0, 49);
+        std::vector<float> data(static_cast<size_t>(rows) * kDim);
+        for (daisy::idx_t r = 0; r < rows; ++r)
+        {
+            const float *centre = centres.data() + pick(rng) * kDim;
+            for (daisy::idx_t d = 0; d < kDim; ++d)
+                data[static_cast<size_t>(r) * kDim + d] = centre[d] + normal(rng);
+        }
+        return data;
+    }
+}
+
+TEST(MessiTopKDistinctTest, ReturnsDistinctIdsAndMatchesBruteForce)
+{
+    auto database = clusteredData(kRows, 1);
+    auto queries = clusteredData(kQueries, 2);
+
+    daisy::BruteForceSearch brute(daisy::DistanceType::L2_SQUARED);
+    brute.buildIndex(database.data(), kRows, kDim);
+    std::vector<daisy::idx_t> expected_ids(kQueries * kK);
+    std::vector<float> expected_distances(kQueries * kK);
+    brute.searchIndex(queries.data(), kQueries, kK, expected_ids.data(),
+                      expected_distances.data());
+
+    daisy::MessiConfig config;
+    config.search_workers = 4;
+    config.index_workers = 2;
+    daisy::Messi messi(daisy::DistanceType::L2_SQUARED, config);
+    messi.setNormalized(false);
+    messi.buildIndex(database.data(), kRows, kDim);
+    std::vector<daisy::idx_t> ids(kQueries * kK);
+    std::vector<float> distances(kQueries * kK);
+    messi.searchIndex(queries.data(), kQueries, kK, ids.data(), distances.data());
+
+    size_t found = 0;
+    for (daisy::idx_t q = 0; q < kQueries; ++q)
+    {
+        std::set<daisy::idx_t> returned(ids.begin() + q * kK, ids.begin() + (q + 1) * kK);
+        EXPECT_EQ(returned.size(), static_cast<size_t>(kK))
+            << "query " << q << " returned duplicate ids";
+        EXPECT_EQ(ids[q * kK], expected_ids[q * kK]) << "nearest neighbour differs, query " << q;
+        for (daisy::idx_t j = 0; j < kK; ++j)
+            found += returned.count(expected_ids[q * kK + j]);
+    }
+    const double recall = static_cast<double>(found) / (kQueries * kK);
+    EXPECT_GE(recall, 0.99) << "MESSI top-" << kK << " recall vs brute force";
+}