diff --git a/CONTINUITY.md b/CONTINUITY.md index 52b3354..c67f0b7 100644 --- a/CONTINUITY.md +++ b/CONTINUITY.md @@ -65,7 +65,10 @@ These algorithms do not all share the same dependency profile or feature set. - Top-k search is the baseline capability. - Range search is modeled through `SearchConfig` and is not universally implemented. -- Streaming insert is not a generic capability; default behavior throws. +- Streaming insert is implemented by `BruteForceSearch`, `LbBruteforce`, and `Coconut`; + the base implementation still throws for algorithms that do not support it. +- Bruteforce streaming grows the owned in-memory database incrementally. LbBruteforce also + computes SAX summaries incrementally, using the breakpoint set fixed by the initial build. - `setNormalized(bool)` is a declaration about the input data, not a preprocessing step. - Data can come from in-memory arrays or file-backed sources via `DataSource`. diff --git a/README.md b/README.md index 902b3dc..124d7cb 100644 --- a/README.md +++ b/README.md @@ -46,8 +46,8 @@ The following table summarizes the key features of each algorithm: | Algorithm | Description | |-----------|-------------| -| **Bruteforce** | Naive parallel similarity search implementation | -| **Lower Bound Bruteforce** | Optimized bruteforce with lower bounding for the distance calculations | +| **Bruteforce** | Naive parallel similarity search implementation with incremental streaming inserts | +| **Lower Bound Bruteforce** | Optimized bruteforce with lower bounding and incremental streaming inserts | | **[MESSI](https://helios2.mi.parisdescartes.fr/~themisp/messi/)** | In-memory parallel similarity search | | **[PARIS](https://helios2.mi.parisdescartes.fr/~themisp/paris/)** | Disk-based parallel similarity search | | **[SING](https://helios2.mi.parisdescartes.fr/~themisp/sing/)** | GPU-accelerated in-memory parallel similarity search | @@ -58,6 +58,23 @@ The following table summarizes the key features of each algorithm: | **[FreSH](http://publications.ics.forth.gr/tech-reports/2023/2023.TR489_FreSh_A_LockFree_Data_Series_Index.pdf)** | In-memory lock-free parallel similarity search using an iSAX index (SRDS 2023) | | **[COCONUT](http://www.vldb.org/pvldb/vol11/p677-kondylakis.pdf)** | Sortable-SAX index built bottom-up; supports both static datasets and **streaming** (incremental) inserts (PVLDB 2018) | +### Incremental streaming inserts + +`BruteForceSearch`, `LbBruteforce`, and `Coconut` implement the common streaming API. Build +the initial index once, then append one series or a contiguous batch without rebuilding: + +```cpp +daisy::BruteForceSearch search(daisy::DistanceType::L2_SQUARED); +search.buildIndex(initial_data, initial_size, dim); +search.insert(one_series); +search.insertBatch(batch_data, batch_size); +``` + +Inserted series receive consecutive IDs beginning at the size of the initial database and are +immediately visible to top-k and range searches. `LbBruteforce` computes a SAX summary for each +insert using the breakpoints established during the initial build. Inserts can reallocate the +owned database, so callers should not retain a pointer returned by `getDatabase()` across them. + ## Quickstart @@ -206,5 +223,3 @@ For questions and suggestions through mail, you can contact us at [manos.chatzak - - diff --git a/demos/CMakeLists.txt b/demos/CMakeLists.txt index f4cc45f..6778c9d 100644 --- a/demos/CMakeLists.txt +++ b/demos/CMakeLists.txt @@ -37,6 +37,13 @@ if(BUILD_DEMO) message(STATUS "Include directories added for demo_bruteforce_L2Square.") endif() + add_executable(demo_Bruteforce_Streaming demo_Bruteforce_Streaming.cpp) + target_link_libraries(demo_Bruteforce_Streaming PRIVATE dino_lib commons_lib) + target_include_directories(demo_Bruteforce_Streaming PRIVATE + ${CMAKE_CURRENT_SOURCE_DIR}/../lib + ${CMAKE_CURRENT_SOURCE_DIR}/../commons + ) + # ////// COCONUT (static + streaming) ////// if(BUILD_COCONUT) if(DEBUG_MSG) @@ -124,6 +131,13 @@ if(BUILD_DEMO) message(STATUS "Include directories added for demo_LbBruteforce_L2Square.") endif() + add_executable(demo_LbBruteforce_Streaming demo_LbBruteforce_Streaming.cpp) + target_link_libraries(demo_LbBruteforce_Streaming PRIVATE dino_lib commons_lib) + target_include_directories(demo_LbBruteforce_Streaming PRIVATE + ${CMAKE_CURRENT_SOURCE_DIR}/../lib + ${CMAKE_CURRENT_SOURCE_DIR}/../commons + ) + # ////// LBBRUTEFORCE DTW ////// if(DEBUG_MSG) message(STATUS "---") diff --git a/demos/demo_Bruteforce_Streaming.cpp b/demos/demo_Bruteforce_Streaming.cpp new file mode 100644 index 0000000..e1e0da4 --- /dev/null +++ b/demos/demo_Bruteforce_Streaming.cpp @@ -0,0 +1,40 @@ +// Bruteforce streaming: build once, then append new series without rebuilding. + +#include "../commons/dataloaders.hpp" +#include "../lib/daisy.hpp" + +#include + +int main() +{ + const daisy::idx_t dim = 96; + const daisy::idx_t initial = 50000; + const daisy::idx_t batch = 25000; + const daisy::idx_t n_query = 5; + const daisy::idx_t k = 5; + + float *stream = loadRandomData(initial + batch, dim, 100, true); + float *query = loadRandomData(n_query, dim, 50, true); + + daisy::BruteForceSearch search(daisy::DistanceType::L2_SQUARED); + search.buildIndex(stream, initial, dim); + search.insert(stream + initial * dim); + search.insertBatch(stream + (initial + 1) * dim, batch - 1); + + daisy::idx_t *indices = new daisy::idx_t[n_query * k]; + float *distances = new float[n_query * k]; + search.searchIndex(query, n_query, k, indices, distances); + + std::printf("Bruteforce now contains %llu series. Query 0 kNN: ", + search.getNDatabase()); + for (daisy::idx_t j = 0; j < k; ++j) + std::printf("%llu(%.3f) ", indices[j], distances[j]); + std::printf("\n"); + + delete[] stream; + delete[] query; + delete[] indices; + delete[] distances; + return 0; +} + diff --git a/demos/demo_Bruteforce_Streaming.py b/demos/demo_Bruteforce_Streaming.py new file mode 100644 index 0000000..2feff1e --- /dev/null +++ b/demos/demo_Bruteforce_Streaming.py @@ -0,0 +1,40 @@ +import numpy as np + +from daisy import BruteForceSearch, DistanceType, LbBruteforce + + +def run_streaming(index, name, initial, single, batch, query): + index.buildIndex(initial) + index.insert(single) + index.insertBatch(batch) + indices, distances = index.searchIndex(query, 3) + print(f"{name} query 0 indices:", indices[0]) + print(f"{name} query 0 distances:", distances[0]) + + +def main(): + rng = np.random.default_rng(100) + stream = rng.normal(size=(1000, 32)).astype(np.float32) + query = rng.normal(size=(5, 32)).astype(np.float32) + + run_streaming( + BruteForceSearch(DistanceType.L2_SQUARED), + "BruteForceSearch", + stream[:500], + stream[500], + stream[501:], + query, + ) + run_streaming( + LbBruteforce(DistanceType.L2_SQUARED), + "LbBruteforce", + stream[:500], + stream[500], + stream[501:], + query, + ) + + +if __name__ == "__main__": + main() + diff --git a/demos/demo_LbBruteforce_Streaming.cpp b/demos/demo_LbBruteforce_Streaming.cpp new file mode 100644 index 0000000..881187b --- /dev/null +++ b/demos/demo_LbBruteforce_Streaming.cpp @@ -0,0 +1,40 @@ +// Lower-bound Bruteforce streaming: raw series and their SAX summaries are appended together. + +#include "../commons/dataloaders.hpp" +#include "../lib/daisy.hpp" + +#include + +int main() +{ + const daisy::idx_t dim = 96; + const daisy::idx_t initial = 50000; + const daisy::idx_t batch = 25000; + const daisy::idx_t n_query = 5; + const daisy::idx_t k = 5; + + float *stream = loadRandomData(initial + batch, dim, 100, true); + float *query = loadRandomData(n_query, dim, 50, true); + + daisy::LbBruteforce search(daisy::DistanceType::L2_SQUARED); + search.buildIndex(stream, initial, dim); + search.insert(stream + initial * dim); + search.insertBatch(stream + (initial + 1) * dim, batch - 1); + + daisy::idx_t *indices = new daisy::idx_t[n_query * k]; + float *distances = new float[n_query * k]; + search.searchIndex(query, n_query, k, indices, distances); + + std::printf("LbBruteforce now contains %llu series. Query 0 kNN: ", + search.getNDatabase()); + for (daisy::idx_t j = 0; j < k; ++j) + std::printf("%llu(%.3f) ", indices[j], distances[j]); + std::printf("\n"); + + delete[] stream; + delete[] query; + delete[] indices; + delete[] distances; + return 0; +} + diff --git a/docs/demos-guide.md b/docs/demos-guide.md index 32b3315..f825f2a 100644 --- a/docs/demos-guide.md +++ b/docs/demos-guide.md @@ -4,9 +4,9 @@ The demos module provides practical examples of how to use the DaiSy library's a Each demo illustrates a specific algorithm or specific distance metric. This module includes both C++ and Python implementations for various algorithms and use cases. Most demos follow the same batch pattern: `buildIndex(...)` once, then `searchIndex(...)`. -The **Coconut** algorithm additionally supports **streaming**: `demo_Coconut_L2Square` shows -the static (batch) build, while `demo_Coconut_Streaming` builds on an initial batch and then -`insert`/`insertBatch`es new series into the live index, querying after each step. +**Bruteforce**, **LbBruteforce**, and **Coconut** additionally support streaming through +`insert(...)` and `insertBatch(...)`. See `demo_Bruteforce_Streaming`, +`demo_LbBruteforce_Streaming`, and `demo_Coconut_Streaming` for live-index examples. ## Demo Program Structure diff --git a/docs/how-to-contribute.md b/docs/how-to-contribute.md index 0182b02..adc1939 100644 --- a/docs/how-to-contribute.md +++ b/docs/how-to-contribute.md @@ -12,7 +12,7 @@ Here is a (non-exhaustive) mockup of our future and ongoing goals: - Extension of DaiSy for subsequence similarity search - Extension of DaiSy for more algorithms (e.g., SFA, Hercules, Dumpy, etc.) -- Streaming / updatable indexing for more algorithms (currently supported by Coconut) +- Streaming / updatable indexing for more algorithms (currently supported by Bruteforce, LbBruteforce, and Coconut) - Implementation of a DaiSy autotuner to automatically optimize indexing and search parameters - Extension of DaiSy to support learned optimization approaches, e.g., LeaFi and ProS diff --git a/lib/algos/Bruteforce.cpp b/lib/algos/Bruteforce.cpp index 8697c13..8b2677b 100644 --- a/lib/algos/Bruteforce.cpp +++ b/lib/algos/Bruteforce.cpp @@ -1,5 +1,11 @@ #include "Bruteforce.hpp" +#include +#include +#include +#include +#include + namespace daisy { @@ -31,33 +37,124 @@ namespace daisy void BruteForceSearch::buildIndex(DataSource *data_source) { - this->dim = data_source->getDim(); - this->n_database = data_source->getTotalRecords(); + if (data_source == nullptr) + throw std::invalid_argument("BruteForceSearch::buildIndex received a null data source"); - if (this->n_database == 0) - { + const idx_t new_dim = data_source->getDim(); + idx_t new_n_database = data_source->getTotalRecords(); + if (new_dim == 0) + throw std::invalid_argument("BruteForceSearch::buildIndex requires a positive dimension"); + + if (new_n_database == 0) + { data_source->reset(); idx_t count = 0; - float *dummy = new float[this->dim]; - while (data_source->nextRecord(dummy)) + std::vector dummy(new_dim); + while (data_source->nextRecord(dummy.data())) { count++; } - delete[] dummy; - this->n_database = count; + new_n_database = count; data_source->reset(); } - this->database = new float[this->n_database * this->dim]; - float *record = new float[this->dim]; + if (new_n_database > std::numeric_limits::max() / new_dim) + throw std::length_error("BruteForceSearch database is too large"); + + const idx_t new_capacity = std::max(new_n_database, 1); + std::unique_ptr new_database( + new float[static_cast(new_capacity) * static_cast(new_dim)]); + std::vector record(new_dim); idx_t idx = 0; - while (data_source->nextRecord(record)) + while (idx < new_n_database && data_source->nextRecord(record.data())) { - std::copy(record, record + this->dim, this->database + idx * this->dim); + std::copy(record.begin(), record.end(), + new_database.get() + static_cast(idx) * new_dim); idx++; } - delete[] record; + + delete[] this->database; + this->database = new_database.release(); + this->dim = new_dim; + this->n_database = idx; + this->database_capacity = new_capacity; + } + + void BruteForceSearch::reserveDatabase(idx_t required_capacity) + { + if (required_capacity <= this->database_capacity) + return; + + idx_t new_capacity = std::max(this->database_capacity, 1); + while (new_capacity < required_capacity) + { + if (new_capacity > std::numeric_limits::max() / 2) + { + new_capacity = required_capacity; + break; + } + new_capacity *= 2; + } + + if (new_capacity > std::numeric_limits::max() / this->dim) + throw std::length_error("BruteForceSearch database is too large"); + + std::unique_ptr grown_database( + new float[static_cast(new_capacity) * static_cast(this->dim)]); + std::copy_n(this->database, + static_cast(this->n_database) * static_cast(this->dim), + grown_database.get()); + + delete[] this->database; + this->database = grown_database.release(); + this->database_capacity = new_capacity; + } + + void BruteForceSearch::insert(const float *series) + { + insertBatch(series, 1); + } + + void BruteForceSearch::insertBatch(const float *data, idx_t n) + { + if (n == 0) + return; + if (this->database == nullptr || this->dim == 0) + throw std::runtime_error("BruteForceSearch::insertBatch requires an initial buildIndex first"); + if (data == nullptr) + throw std::invalid_argument("BruteForceSearch::insertBatch received null data"); + if (n > std::numeric_limits::max() - this->n_database) + throw std::length_error("BruteForceSearch database size overflow"); + if (n > std::numeric_limits::max() / this->dim) + throw std::length_error("BruteForceSearch insert batch is too large"); + + const size_t current_values = + static_cast(this->n_database) * static_cast(this->dim); + const size_t inserted_values = static_cast(n) * static_cast(this->dim); + const uintptr_t database_begin = reinterpret_cast(this->database); + const uintptr_t database_end = database_begin + current_values * sizeof(float); + const uintptr_t data_address = reinterpret_cast(data); + const bool aliases_database = data_address >= database_begin && data_address < database_end; + size_t source_offset = 0; + if (aliases_database) + { + const uintptr_t byte_offset = data_address - database_begin; + if (byte_offset % sizeof(float) != 0) + throw std::invalid_argument("BruteForceSearch::insertBatch received an unaligned database pointer"); + source_offset = static_cast(byte_offset / sizeof(float)); + if (inserted_values > current_values - source_offset) + throw std::invalid_argument("BruteForceSearch::insertBatch source exceeds the live database"); + } + + const idx_t required_capacity = this->n_database + n; + reserveDatabase(required_capacity); + const float *source = aliases_database ? this->database + source_offset : data; + + std::copy_n(source, + inserted_values, + this->database + static_cast(this->n_database) * this->dim); + this->n_database = required_capacity; } void BruteForceSearch::searchIndexL2Squared(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D) @@ -194,6 +291,7 @@ namespace daisy I.resize(n_query); D.resize(n_query); + const float abandon_bound = std::nextafter(config.r, FLT_MAX); #pragma omp parallel num_threads(num_threads) { @@ -209,7 +307,7 @@ namespace daisy float dist = distance_computer->compute_dist(const_cast(q_vec), const_cast(db_vec), dim, - config.r); + abandon_bound); if (dist <= config.r) hits.emplace_back(dist, dbi); } @@ -230,5 +328,6 @@ namespace daisy BruteForceSearch::~BruteForceSearch() { delete[] database; + database = nullptr; } -} \ No newline at end of file +} diff --git a/lib/algos/Bruteforce.hpp b/lib/algos/Bruteforce.hpp index 645b8d1..a7f0d6c 100644 --- a/lib/algos/Bruteforce.hpp +++ b/lib/algos/Bruteforce.hpp @@ -13,6 +13,9 @@ namespace daisy class BruteForceSearch : public SimilaritySearchAlgorithm { private: + idx_t database_capacity = 0; + + void reserveDatabase(idx_t required_capacity); void searchIndexL2Squared(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D); void searchIndexDTW(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D); @@ -30,15 +33,19 @@ namespace daisy throw std::runtime_error("BruteForceSearch requires in-memory data. Use buildIndex(database, n_database, dim) instead."); } + // Append owned copies of new series to the live in-memory database. + void insert(const float *series) override; + void insertBatch(const float *data, idx_t n) override; + void searchIndex(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D) override; void searchIndex(const float *query, idx_t n_query, const SearchConfig &config, std::vector> &I, std::vector> &D) override; - ~BruteForceSearch(); + ~BruteForceSearch() override; }; } -#endif \ No newline at end of file +#endif diff --git a/lib/algos/LbBruteforce.cpp b/lib/algos/LbBruteforce.cpp index 218448f..80359f1 100644 --- a/lib/algos/LbBruteforce.cpp +++ b/lib/algos/LbBruteforce.cpp @@ -1,6 +1,12 @@ #include "LbBruteforce.hpp" #include "../isax/iSAXIndex.hpp" +#include +#include +#include +#include +#include + namespace daisy { @@ -37,9 +43,15 @@ namespace daisy void LbBruteforce::buildIndex(DataSource *data_source) { + if (data_source == nullptr) + throw std::invalid_argument("LbBruteforce::buildIndex received a null data source"); + this->dim = data_source->getDim(); this->n_database = data_source->getTotalRecords(); + if (this->dim == 0) + throw std::invalid_argument("LbBruteforce::buildIndex requires a positive dimension"); + if (this->n_database == 0) { @@ -55,7 +67,11 @@ namespace daisy data_source->reset(); } - this->database = new float[this->n_database * this->dim]; + this->database_capacity = std::max(this->n_database, 1); + if (this->database_capacity > std::numeric_limits::max() / this->dim) + throw std::length_error("LbBruteforce database is too large"); + + this->database = new float[static_cast(this->database_capacity) * this->dim]; float *record = new float[this->dim]; idx_t idx = 0; while (data_source->nextRecord(record)) @@ -87,7 +103,10 @@ namespace daisy this->index = isax_index_init_inmemory(this->index_settings); - this->db_sax_representations = (sax_type **)malloc(n_database * sizeof(sax_type *)); + this->db_sax_representations = + (sax_type **)malloc(static_cast(this->database_capacity) * sizeof(sax_type *)); + if (this->db_sax_representations == nullptr) + throw std::bad_alloc(); #pragma omp parallel for num_threads(num_threads) for (idx_t dbi = 0; dbi < n_database; dbi++) @@ -108,6 +127,130 @@ namespace daisy } } + void LbBruteforce::reserveDatabase(idx_t required_capacity) + { + if (required_capacity <= this->database_capacity) + return; + + idx_t new_capacity = std::max(this->database_capacity, 1); + while (new_capacity < required_capacity) + { + if (new_capacity > std::numeric_limits::max() / 2) + { + new_capacity = required_capacity; + break; + } + new_capacity *= 2; + } + + if (new_capacity > std::numeric_limits::max() / this->dim || + new_capacity > std::numeric_limits::max() / sizeof(sax_type *)) + throw std::length_error("LbBruteforce database is too large"); + + std::unique_ptr grown_database( + new float[static_cast(new_capacity) * static_cast(this->dim)]); + std::copy_n(this->database, + static_cast(this->n_database) * static_cast(this->dim), + grown_database.get()); + + sax_type **grown_sax = + (sax_type **)malloc(static_cast(new_capacity) * sizeof(sax_type *)); + if (grown_sax == nullptr) + throw std::bad_alloc(); + std::copy_n(this->db_sax_representations, + static_cast(this->n_database), grown_sax); + + delete[] this->database; + free(this->db_sax_representations); + this->database = grown_database.release(); + this->db_sax_representations = grown_sax; + this->database_capacity = new_capacity; + } + + void LbBruteforce::insert(const float *series) + { + insertBatch(series, 1); + } + + void LbBruteforce::insertBatch(const float *data, idx_t n) + { + if (n == 0) + return; + if (this->database == nullptr || this->index == nullptr || + this->index_settings == nullptr || this->dim == 0) + throw std::runtime_error("LbBruteforce::insertBatch requires an initial buildIndex first"); + if (data == nullptr) + throw std::invalid_argument("LbBruteforce::insertBatch received null data"); + if (n > std::numeric_limits::max() - this->n_database) + throw std::length_error("LbBruteforce database size overflow"); + if (n > std::numeric_limits::max() / this->dim) + throw std::length_error("LbBruteforce insert batch is too large"); + + const idx_t required_capacity = this->n_database + n; + const size_t current_values = + static_cast(this->n_database) * static_cast(this->dim); + const size_t inserted_values = static_cast(n) * static_cast(this->dim); + const uintptr_t database_begin = reinterpret_cast(this->database); + const uintptr_t database_end = database_begin + current_values * sizeof(float); + const uintptr_t data_address = reinterpret_cast(data); + const bool aliases_database = data_address >= database_begin && data_address < database_end; + size_t source_offset = 0; + if (aliases_database) + { + const uintptr_t byte_offset = data_address - database_begin; + if (byte_offset % sizeof(float) != 0) + throw std::invalid_argument("LbBruteforce::insertBatch received an unaligned database pointer"); + source_offset = static_cast(byte_offset / sizeof(float)); + if (inserted_values > current_values - source_offset) + throw std::invalid_argument("LbBruteforce::insertBatch source exceeds the live database"); + } + + reserveDatabase(required_capacity); + const float *source_data = aliases_database ? this->database + source_offset : data; + + activateBreakpoints(); + + std::vector new_sax_records; + new_sax_records.reserve(static_cast(n)); + try + { + for (idx_t i = 0; i < n; ++i) + { + sax_type *sax = (sax_type *)malloc( + sizeof(sax_type) * this->index->settings->paa_segments); + if (sax == nullptr) + throw std::bad_alloc(); + + const float *series = source_data + static_cast(i) * this->dim; + if (!this->distance_computer->compute_sax_from_ts( + series, + sax, + this->index->settings->ts_values_per_paa_segment, + this->index->settings->paa_segments, + this->index->settings->sax_alphabet_cardinality, + this->index->settings->sax_bit_cardinality)) + { + free(sax); + throw std::runtime_error("LbBruteforce::insertBatch failed to compute SAX representation"); + } + new_sax_records.push_back(sax); + } + + std::copy_n(source_data, + inserted_values, + this->database + static_cast(this->n_database) * this->dim); + for (idx_t i = 0; i < n; ++i) + this->db_sax_representations[this->n_database + i] = new_sax_records[i]; + this->n_database = required_capacity; + } + catch (...) + { + for (sax_type *sax : new_sax_records) + free(sax); + throw; + } + } + void LbBruteforce::searchIndex(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D) { if (this->distance_type == DistanceType::L2_SQUARED) @@ -127,6 +270,10 @@ namespace daisy void LbBruteforce::searchIndexL2Squared(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D) { + if (!validateSearchParams(k, n_query)) + return; + activateBreakpoints(); + #pragma omp parallel num_threads(num_threads) { #pragma omp for @@ -196,6 +343,10 @@ namespace daisy void LbBruteforce::searchIndexDTW(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D) { + if (!validateSearchParams(k, n_query)) + return; + activateBreakpoints(); + #pragma omp parallel num_threads(num_threads) { #pragma omp for @@ -297,6 +448,12 @@ namespace daisy return; } + if (database == nullptr || index == nullptr) + throw std::runtime_error("LbBruteforce index must be built before searching"); + if (n_query == 0) + throw std::invalid_argument("n_query must be greater than 0"); + activateBreakpoints(); + float r = config.r; I.assign(n_query, {}); D.assign(n_query, {}); @@ -364,17 +521,31 @@ namespace daisy free(db_sax_representations[dbi]); } free(db_sax_representations); + db_sax_representations = nullptr; if (this->index) { + if (this->index->fbl) + { + destroy_fbl(this->index->fbl); + this->index->fbl = nullptr; + } + free(this->index->answer); + this->index->answer = nullptr; free(this->index); this->index = nullptr; } if (this->index_settings) { + if (daisy_active_breakpoints == this->index_settings->breakpoints) + set_active_breakpoints(nullptr, nullptr); + free(this->index_settings->max_sax_cardinalities); + free(this->index_settings->bit_masks); + free(this->index_settings->breakpoints_owned); + free(this->index_settings->breakpoints_max_owned); free(this->index_settings); this->index_settings = nullptr; } } -} \ No newline at end of file +} diff --git a/lib/algos/LbBruteforce.hpp b/lib/algos/LbBruteforce.hpp index 56e7e7f..f7600fc 100644 --- a/lib/algos/LbBruteforce.hpp +++ b/lib/algos/LbBruteforce.hpp @@ -17,6 +17,11 @@ namespace daisy class LbBruteforce : public SimilaritySearchAlgorithm { + private: + idx_t database_capacity = 0; + + void reserveDatabase(idx_t required_capacity); + public: LbBruteforce(DistanceType distance_type); void setNumThreads(int num_threads); @@ -51,6 +56,10 @@ namespace daisy throw std::runtime_error("LbBruteforce requires in-memory data. Use buildIndex(database, n_database, dim) instead."); } + // Append owned copies and compute their SAX summaries with the build-time breakpoints. + void insert(const float *series) override; + void insertBatch(const float *data, idx_t n) override; + void searchIndex(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D) override; void searchIndex(const float *query, idx_t n_query, const SearchConfig &config, std::vector> &I, @@ -58,9 +67,9 @@ namespace daisy void searchIndexL2Squared(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D); void searchIndexDTW(const float *query, const idx_t n_query, const idx_t k, idx_t *I, float *D); - ~LbBruteforce(); + ~LbBruteforce() override; }; } -#endif \ No newline at end of file +#endif diff --git a/lib/algos/SimilaritySearchAlgorithm.hpp b/lib/algos/SimilaritySearchAlgorithm.hpp index a41d051..defe691 100644 --- a/lib/algos/SimilaritySearchAlgorithm.hpp +++ b/lib/algos/SimilaritySearchAlgorithm.hpp @@ -104,7 +104,7 @@ namespace daisy } // Streaming API: incrementally add series to a live index. Only some algorithms - // support it (currently Coconut); the default throws. + // support it (BruteForceSearch, LbBruteforce, and Coconut); the default throws. virtual void insert(const float *series) { (void)series; diff --git a/lib/isax/SAX.cpp b/lib/isax/SAX.cpp index c5286f0..2116c06 100644 --- a/lib/isax/SAX.cpp +++ b/lib/isax/SAX.cpp @@ -637,6 +637,14 @@ namespace daisy __m256i region_lowerv_1_offset = _mm256_add_epi32(region_lowerv_1, vloweroffset); __m256i region_upperv_0_offset = _mm256_add_epi32(region_upperv_0, vupperoffset); __m256i region_upperv_1_offset = _mm256_add_epi32(region_upperv_1, vupperoffset); + // The upper edge of the last SAX region is represented by max_val rather than a + // table entry. Clamp that lane before the unconditional AVX2 gather; the gathered + // value is discarded by upper_juge_maxv_* below, but reading index tri_size would + // otherwise be one float past the breakpoint table. + const __m256i last_breakpoint_offset = + _mm256_set1_epi32(offset + max_cardinality - 2); + region_upperv_0_offset = _mm256_min_epi32(region_upperv_0_offset, last_breakpoint_offset); + region_upperv_1_offset = _mm256_min_epi32(region_upperv_1_offset, last_breakpoint_offset); _mm256_storeu_si256((__m256i_u *)&(region_lower[0]), region_lowerv_0); _mm256_storeu_si256((__m256i_u *)&(region_lower[8]), region_lowerv_1); _mm256_storeu_si256((__m256i_u *)&(region_upper[0]), region_upperv_0); @@ -1138,4 +1146,4 @@ namespace daisy distance = ratio_sqrt * distance; return distance; } -} \ No newline at end of file +} diff --git a/pybinds/setup.cpp b/pybinds/setup.cpp index 480795f..c781f0c 100644 --- a/pybinds/setup.cpp +++ b/pybinds/setup.cpp @@ -82,6 +82,25 @@ PYBIND11_MODULE(_core, m) daisy::InMemoryDataSource data_source(static_cast(buf.ptr), n, d); self.buildIndex(&data_source); }, "Build the index from a 2D float32 numpy array") + // Streaming: append one vector or a batch to the live in-memory database. + .def("insert", [](daisy::BruteForceSearch &self, pybind11::array_t series) + { + pybind11::buffer_info buf = series.request(); + if (buf.ndim != 1) + throw std::runtime_error("insert expects a 1D float32 array"); + if (self.getDim() != 0 && static_cast(buf.shape[0]) != self.getDim()) + throw std::runtime_error("insert series dimension does not match the index dimension"); + self.insert(static_cast(buf.ptr)); }, "Incrementally insert one series into the live index") + + .def("insertBatch", [](daisy::BruteForceSearch &self, pybind11::array_t batch) + { + pybind11::buffer_info buf = batch.request(); + if (buf.ndim != 2) + throw std::runtime_error("insertBatch expects a 2D float32 array"); + if (self.getDim() != 0 && static_cast(buf.shape[1]) != self.getDim()) + throw std::runtime_error("insertBatch series dimension does not match the index dimension"); + self.insertBatch(static_cast(buf.ptr), buf.shape[0]); }, "Incrementally insert a batch of series into the live index") + // Bind method to perform similarity search .def("searchIndex", [](daisy::BruteForceSearch &self, pybind11::array_t query, daisy::idx_t k) { @@ -165,6 +184,25 @@ PYBIND11_MODULE(_core, m) daisy::InMemoryDataSource data_source(static_cast(buf.ptr), n, d); self.buildIndex(&data_source); }, "Build the index from a 2D float32 numpy array") + // Streaming: append raw series and their SAX summaries to the live index. + .def("insert", [](daisy::LbBruteforce &self, pybind11::array_t series) + { + pybind11::buffer_info buf = series.request(); + if (buf.ndim != 1) + throw std::runtime_error("insert expects a 1D float32 array"); + if (self.getDim() != 0 && static_cast(buf.shape[0]) != self.getDim()) + throw std::runtime_error("insert series dimension does not match the index dimension"); + self.insert(static_cast(buf.ptr)); }, "Incrementally insert one series and its SAX summary into the live index") + + .def("insertBatch", [](daisy::LbBruteforce &self, pybind11::array_t batch) + { + pybind11::buffer_info buf = batch.request(); + if (buf.ndim != 2) + throw std::runtime_error("insertBatch expects a 2D float32 array"); + if (self.getDim() != 0 && static_cast(buf.shape[1]) != self.getDim()) + throw std::runtime_error("insertBatch series dimension does not match the index dimension"); + self.insertBatch(static_cast(buf.ptr), buf.shape[0]); }, "Incrementally insert a batch of series and their SAX summaries into the live index") + // Search the index using a query array and return (indices, distances) .def("searchIndex", [](daisy::LbBruteforce &self, pybind11::array_t query, daisy::idx_t k) { diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 9136ae7..f3c6f08 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -82,6 +82,14 @@ if(DEBUG_MSG) message(STATUS "Tests discovered for test_bruteforce_L2Square.") endif() +# ////// BRUTEFORCE Streaming ////// +add_executable(test_Bruteforce_Streaming test_Bruteforce_Streaming.cpp) +target_link_libraries(test_Bruteforce_Streaming PRIVATE GTest::gtest_main dino_lib commons_lib stdc++fs) +target_include_directories(test_Bruteforce_Streaming PRIVATE + ${CMAKE_CURRENT_SOURCE_DIR}/../lib + ${CMAKE_CURRENT_SOURCE_DIR}/../commons) +gtest_discover_tests(test_Bruteforce_Streaming WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}) + # # ////// BRUTEFORCE DTW ////// # if(DEBUG_MSG) # message(STATUS "---") @@ -186,6 +194,14 @@ if(DEBUG_MSG) message(STATUS "Tests discovered for test_LbBruteforce_L2Square.") endif() +# ////// LBBRUTEFORCE Streaming ////// +add_executable(test_LbBruteforce_Streaming test_LbBruteforce_Streaming.cpp) +target_link_libraries(test_LbBruteforce_Streaming PRIVATE GTest::gtest_main dino_lib commons_lib stdc++fs) +target_include_directories(test_LbBruteforce_Streaming PRIVATE + ${CMAKE_CURRENT_SOURCE_DIR}/../lib + ${CMAKE_CURRENT_SOURCE_DIR}/../commons) +gtest_discover_tests(test_LbBruteforce_Streaming WORKING_DIRECTORY ${CMAKE_SOURCE_DIR}) + # # ////// LBBRUTEFORCE DTW ////// # if(DEBUG_MSG) # message(STATUS "---") @@ -748,6 +764,15 @@ if(BUILD_PYTHON) set_tests_properties(test_range_search_py PROPERTIES ENVIRONMENT "PYTHONPATH=${CMAKE_SOURCE_DIR}" ) + + add_test( + NAME test_streaming_py + COMMAND ${Python3_EXECUTABLE} ${CMAKE_CURRENT_SOURCE_DIR}/test_streaming.py + WORKING_DIRECTORY ${CMAKE_SOURCE_DIR} + ) + set_tests_properties(test_streaming_py PROPERTIES + ENVIRONMENT "PYTHONPATH=${CMAKE_SOURCE_DIR}" + ) endif() endif() diff --git a/tests/test_Bruteforce_Streaming.cpp b/tests/test_Bruteforce_Streaming.cpp new file mode 100644 index 0000000..3fda828 --- /dev/null +++ b/tests/test_Bruteforce_Streaming.cpp @@ -0,0 +1,110 @@ +#include + +#include "../lib/algos/Bruteforce.hpp" + +#include +#include + +namespace +{ + constexpr int DIM = 32; + + std::vector makeSeries(int n, int first_phase = 0) + { + std::vector data(static_cast(n) * DIM); + for (int i = 0; i < n; ++i) + { + float *series = data.data() + static_cast(i) * DIM; + const float phase = static_cast(first_phase + i) * 0.37f; + double mean = 0.0; + for (int j = 0; j < DIM; ++j) + { + series[j] = std::sin(0.21f * j + phase) + + 0.35f * std::cos(0.47f * j - phase); + mean += series[j]; + } + mean /= DIM; + + double variance = 0.0; + for (int j = 0; j < DIM; ++j) + variance += (series[j] - mean) * (series[j] - mean); + const double stddev = std::sqrt(variance / DIM); + for (int j = 0; j < DIM; ++j) + series[j] = static_cast((series[j] - mean) / stddev); + } + return data; + } +} + +TEST(BruteforceStreamingTest, RequiresBuildAndRejectsNullInput) +{ + daisy::BruteForceSearch search(daisy::DistanceType::L2_SQUARED); + EXPECT_THROW(search.insert(nullptr), std::runtime_error); + + auto initial = makeSeries(2); + search.buildIndex(initial.data(), 2, DIM); + EXPECT_THROW(search.insert(nullptr), std::invalid_argument); +} + +TEST(BruteforceStreamingTest, SingleAndBatchInsertsAreImmediatelySearchable) +{ + auto all = makeSeries(8); + daisy::BruteForceSearch search(daisy::DistanceType::L2_SQUARED); + search.buildIndex(all.data(), 3, DIM); + + daisy::SimilaritySearchAlgorithm *streaming = &search; + streaming->insert(all.data() + 3 * DIM); + streaming->insertBatch(all.data() + 4 * DIM, 4); + + ASSERT_EQ(search.getNDatabase(), 8u); + for (int j = 0; j < DIM; ++j) + EXPECT_FLOAT_EQ(search.getDatabase()[7 * DIM + j], all[7 * DIM + j]); + + daisy::idx_t index = 0; + float distance = -1.0f; + search.searchIndex(all.data() + 7 * DIM, 1, 1, &index, &distance); + EXPECT_EQ(index, 7u); + EXPECT_FLOAT_EQ(distance, 0.0f); + + daisy::SearchConfig range; + range.type = daisy::QueryType::RANGE; + range.r = 0.0f; + std::vector> indices; + std::vector> distances; + search.searchIndex(all.data() + 3 * DIM, 1, range, indices, distances); + ASSERT_EQ(indices.size(), 1u); + ASSERT_EQ(indices[0].size(), 1u); + EXPECT_EQ(indices[0][0], 3u); + EXPECT_FLOAT_EQ(distances[0][0], 0.0f); +} + +TEST(BruteforceStreamingTest, DtwSearchIncludesInsertedSeries) +{ + auto all = makeSeries(6, 20); + daisy::BruteForceSearch search(daisy::DistanceType::DTW); + search.buildIndex(all.data(), 3, DIM); + search.insert(all.data() + 3 * DIM); + search.insertBatch(all.data() + 4 * DIM, 2); + + daisy::idx_t index = 0; + float distance = -1.0f; + search.searchIndex(all.data() + 5 * DIM, 1, 1, &index, &distance); + EXPECT_EQ(index, 5u); + EXPECT_NEAR(distance, 0.0f, 1e-6f); +} + +TEST(BruteforceStreamingTest, CanInsertFromItsOwnDatabaseAcrossReallocation) +{ + auto initial = makeSeries(4, 150); + daisy::BruteForceSearch search(daisy::DistanceType::L2_SQUARED); + search.buildIndex(initial.data(), 4, DIM); + + search.insert(search.getDatabase() + 2 * DIM); + ASSERT_EQ(search.getNDatabase(), 5u); + + daisy::idx_t index = 0; + float distance = -1.0f; + search.searchIndex(search.getDatabase() + 4 * DIM, 1, 1, &index, &distance); + EXPECT_EQ(index, 2u); + EXPECT_FLOAT_EQ(distance, 0.0f); +} diff --git a/tests/test_LbBruteforce_Streaming.cpp b/tests/test_LbBruteforce_Streaming.cpp new file mode 100644 index 0000000..9c5115d --- /dev/null +++ b/tests/test_LbBruteforce_Streaming.cpp @@ -0,0 +1,164 @@ +#include + +#include "../lib/algos/Bruteforce.hpp" +#include "../lib/algos/LbBruteforce.hpp" + +#include +#include + +namespace +{ + constexpr int DIM = 32; + + std::vector makeSeries(int n, int first_phase = 0) + { + std::vector data(static_cast(n) * DIM); + for (int i = 0; i < n; ++i) + { + float *series = data.data() + static_cast(i) * DIM; + const float phase = static_cast(first_phase + i) * 0.31f; + double mean = 0.0; + for (int j = 0; j < DIM; ++j) + { + series[j] = std::sin(0.19f * j + phase) + + 0.4f * std::cos(0.43f * j - 0.5f * phase); + mean += series[j]; + } + mean /= DIM; + + double variance = 0.0; + for (int j = 0; j < DIM; ++j) + variance += (series[j] - mean) * (series[j] - mean); + const double stddev = std::sqrt(variance / DIM); + for (int j = 0; j < DIM; ++j) + series[j] = static_cast((series[j] - mean) / stddev); + } + return data; + } + + void expectMatchesBruteforce(daisy::LbBruteforce &search, + const std::vector &data, + int n_database, + const std::vector &queries, + int n_query, + int k) + { + daisy::BruteForceSearch ground_truth(daisy::DistanceType::L2_SQUARED); + ground_truth.buildIndex(const_cast(data.data()), n_database, DIM); + + std::vector expected_indices(static_cast(n_query) * k); + std::vector expected_distances(static_cast(n_query) * k); + std::vector actual_indices(static_cast(n_query) * k); + std::vector actual_distances(static_cast(n_query) * k); + + ground_truth.searchIndex(queries.data(), n_query, k, + expected_indices.data(), expected_distances.data()); + search.searchIndex(queries.data(), n_query, k, + actual_indices.data(), actual_distances.data()); + + for (size_t i = 0; i < actual_indices.size(); ++i) + { + EXPECT_EQ(actual_indices[i], expected_indices[i]); + EXPECT_NEAR(actual_distances[i], expected_distances[i], 1e-4f); + } + } +} + +TEST(LbBruteforceStreamingTest, RequiresBuildAndRejectsNullInput) +{ + daisy::LbBruteforce search(daisy::DistanceType::L2_SQUARED); + EXPECT_THROW(search.insert(nullptr), std::runtime_error); + + auto initial = makeSeries(2); + search.buildIndex(initial.data(), 2, DIM); + EXPECT_THROW(search.insert(nullptr), std::invalid_argument); +} + +TEST(LbBruteforceStreamingTest, SingleAndBatchInsertsUpdateDataAndSaxSummaries) +{ + auto all = makeSeries(15); + auto queries = makeSeries(4, 40); + + daisy::LbBruteforce search(daisy::DistanceType::L2_SQUARED); + search.setNumThreads(1); + search.buildIndex(all.data(), 8, DIM); + expectMatchesBruteforce(search, all, 8, queries, 4, 4); + + daisy::SimilaritySearchAlgorithm *streaming = &search; + streaming->insert(all.data() + 8 * DIM); + expectMatchesBruteforce(search, all, 9, queries, 4, 4); + + streaming->insertBatch(all.data() + 9 * DIM, 6); + ASSERT_EQ(search.getNDatabase(), 15u); + expectMatchesBruteforce(search, all, 15, queries, 4, 4); + + daisy::idx_t index = 0; + float distance = -1.0f; + search.searchIndex(all.data() + 14 * DIM, 1, 1, &index, &distance); + EXPECT_EQ(index, 14u); + EXPECT_FLOAT_EQ(distance, 0.0f); + + daisy::SearchConfig range; + range.type = daisy::QueryType::RANGE; + range.r = 0.0f; + std::vector> indices; + std::vector> distances; + search.searchIndex(all.data() + 8 * DIM, 1, range, indices, distances); + ASSERT_EQ(indices.size(), 1u); + ASSERT_EQ(indices[0].size(), 1u); + EXPECT_EQ(indices[0][0], 8u); + EXPECT_FLOAT_EQ(distances[0][0], 0.0f); +} + +TEST(LbBruteforceStreamingTest, DtwSearchIncludesInsertedSaxRecords) +{ + auto all = makeSeries(7, 60); + daisy::LbBruteforce search(daisy::DistanceType::DTW); + search.setNumThreads(1); + search.buildIndex(all.data(), 4, DIM); + search.insert(all.data() + 4 * DIM); + search.insertBatch(all.data() + 5 * DIM, 2); + + daisy::idx_t index = 0; + float distance = -1.0f; + search.searchIndex(all.data() + 6 * DIM, 1, 1, &index, &distance); + EXPECT_EQ(index, 6u); + EXPECT_NEAR(distance, 0.0f, 1e-6f); +} + +TEST(LbBruteforceStreamingTest, EquidepthInsertsReuseInitialBreakpoints) +{ + auto all = makeSeries(14, 90); + auto queries = makeSeries(3, 120); + + // Scaling and offset make the input intentionally non-z-normalized. + for (float &value : all) + value = 25.0f + 7.0f * value; + for (float &value : queries) + value = 25.0f + 7.0f * value; + + daisy::LbBruteforce search(daisy::DistanceType::L2_SQUARED); + search.setNormalized(false); + search.setNumThreads(1); + search.buildIndex(all.data(), 8, DIM); + search.insertBatch(all.data() + 8 * DIM, 6); + + expectMatchesBruteforce(search, all, 14, queries, 3, 4); +} + +TEST(LbBruteforceStreamingTest, CanInsertFromItsOwnDatabaseAcrossReallocation) +{ + auto initial = makeSeries(4, 180); + daisy::LbBruteforce search(daisy::DistanceType::L2_SQUARED); + search.setNumThreads(1); + search.buildIndex(initial.data(), 4, DIM); + + search.insert(search.getDatabase() + DIM); + ASSERT_EQ(search.getNDatabase(), 5u); + + daisy::idx_t index = 0; + float distance = -1.0f; + search.searchIndex(search.getDatabase() + 4 * DIM, 1, 1, &index, &distance); + EXPECT_EQ(index, 1u); + EXPECT_FLOAT_EQ(distance, 0.0f); +} diff --git a/tests/test_streaming.py b/tests/test_streaming.py new file mode 100644 index 0000000..ae0ea65 --- /dev/null +++ b/tests/test_streaming.py @@ -0,0 +1,36 @@ +import unittest + +import numpy as np + +from daisy import BruteForceSearch, DistanceType, LbBruteforce + + +class StreamingBindingsTest(unittest.TestCase): + def test_bruteforce_and_lb_bruteforce(self): + rng = np.random.default_rng(123) + data = rng.normal(size=(8, 32)).astype(np.float32) + + for algorithm in (BruteForceSearch, LbBruteforce): + with self.subTest(algorithm=algorithm.__name__): + index = algorithm(DistanceType.L2_SQUARED) + index.buildIndex(data[:3]) + index.insert(data[3]) + index.insertBatch(data[4:]) + + indices, distances = index.searchIndex(data[7:8], 1) + self.assertEqual(int(indices[0, 0]), 7) + self.assertAlmostEqual(float(distances[0, 0]), 0.0, places=6) + + def test_dimension_validation(self): + index = BruteForceSearch(DistanceType.L2_SQUARED) + index.buildIndex(np.zeros((2, 32), dtype=np.float32)) + + with self.assertRaises(RuntimeError): + index.insert(np.zeros(31, dtype=np.float32)) + with self.assertRaises(RuntimeError): + index.insertBatch(np.zeros((2, 31), dtype=np.float32)) + + +if __name__ == "__main__": + unittest.main() +