|
28 | 28 | #include <raft/linalg/map.cuh> |
29 | 29 | #include <raft/linalg/norm.cuh> |
30 | 30 | #include <raft/linalg/transpose.cuh> |
| 31 | +#include <raft/matrix/gather.cuh> |
31 | 32 | #include <raft/matrix/init.cuh> |
32 | 33 | #include <raft/sparse/convert/coo.cuh> |
33 | 34 | #include <raft/sparse/convert/csr.cuh> |
|
41 | 42 | #include <cuda_fp16.h> |
42 | 43 | #include <rmm/cuda_device.hpp> |
43 | 44 | #include <rmm/device_uvector.hpp> |
| 45 | +#include <thrust/copy.h> |
44 | 46 | #include <thrust/for_each.h> |
| 47 | +#include <thrust/gather.h> |
| 48 | +#include <thrust/iterator/counting_iterator.h> |
45 | 49 |
|
| 50 | +#include <algorithm> |
| 51 | +#include <cmath> |
46 | 52 | #include <cstdint> |
47 | 53 | #include <iostream> |
48 | 54 | #include <optional> |
@@ -581,6 +587,150 @@ void brute_force_search( |
581 | 587 | query_norms ? query_norms->data_handle() : nullptr); |
582 | 588 | } |
583 | 589 |
|
| 590 | +/** |
| 591 | + * @brief The three strategies available for a filtered brute-force search. |
| 592 | + * |
| 593 | + * - `sddmm`: evaluate only the passing (query, row) pairs via a masked matmul. Cost |
| 594 | + * grows with the number of passing rows. |
| 595 | + * - `gather`: compact the passing rows into a dense [n_pass x dim] matrix and search |
| 596 | + * that. Cost grows with the passing fraction, over a fixed setup cost. |
| 597 | + * - `dense`: GEMM against the whole dataset, then mask. Cost is independent of the |
| 598 | + * filter. |
| 599 | + */ |
| 600 | +enum class filtered_search_path { sddmm, gather, dense }; |
| 601 | + |
| 602 | +/** |
| 603 | + * @brief Pick the cheapest search path for a filtered brute-force query. |
| 604 | + * |
| 605 | + * The thresholds are empirical. |
| 606 | + * |
| 607 | + * @param[in] n_dataset rows in the dataset |
| 608 | + * @param[in] dim columns in the dataset |
| 609 | + * @param[in] selectivity fraction of (query, row) pairs the filter passes |
| 610 | + * @param[in] k neighbors requested |
| 611 | + * @param[in] passing rows passing the filter, if the filter passes the same rows for |
| 612 | + * every query. A bitmap passes a different set per query, so it has |
| 613 | + * no such count, and cannot use the gather path. |
| 614 | + */ |
| 615 | +inline filtered_search_path select_filtered_search_path( |
| 616 | + int64_t n_dataset, int64_t dim, double selectivity, int64_t k, std::optional<int64_t> passing) |
| 617 | +{ |
| 618 | + // SDDMM competes against a fixed setup cost, so it stops paying off at an absolute row |
| 619 | + // count rather than at a fraction of the dataset. |
| 620 | + constexpr int64_t kSddmmMaxPassingRows = 2000; |
| 621 | + constexpr double kSddmmMaxSelectivity = 0.03; |
| 622 | + // Gathering trades memory traffic (proportional to dim) against a GEMM (proportional to |
| 623 | + // n_dataset), and cannot beat simply GEMMing everything once most rows pass. |
| 624 | + constexpr double kGatherSelectivityScale = 5.5; |
| 625 | + constexpr double kGatherMaxSelectivity = 0.5; |
| 626 | + // Below this, the gather setup costs more than the entire dense search. |
| 627 | + constexpr double kGatherMinDataset = 200000.0; |
| 628 | + |
| 629 | + const auto n_dataset_d = static_cast<double>(n_dataset); |
| 630 | + |
| 631 | + const bool sddmm_beats_dense = selectivity < kSddmmMaxSelectivity; |
| 632 | + if (!passing) { |
| 633 | + return sddmm_beats_dense ? filtered_search_path::sddmm : filtered_search_path::dense; |
| 634 | + } |
| 635 | + const int64_t n_pass = *passing; |
| 636 | + if (n_pass < kSddmmMaxPassingRows && sddmm_beats_dense) { return filtered_search_path::sddmm; } |
| 637 | + |
| 638 | + // Fraction of the dense search left over after paying the gather setup. |
| 639 | + const double amortized = 1.0 - kGatherMinDataset / n_dataset_d; |
| 640 | + const double gather_max_selectivity = |
| 641 | + amortized * |
| 642 | + std::min(kGatherSelectivityScale / std::sqrt(static_cast<double>(dim)), kGatherMaxSelectivity); |
| 643 | + |
| 644 | + // n_pass > k keeps the compacted dataset big enough to hold k neighbors. |
| 645 | + if (amortized > 0.0 && selectivity < gather_max_selectivity && n_pass > k) { |
| 646 | + return filtered_search_path::gather; |
| 647 | + } |
| 648 | + return filtered_search_path::dense; |
| 649 | +} |
| 650 | + |
| 651 | +/** |
| 652 | + * @brief Filtered search over the rows a bitset selects, by compacting them first. |
| 653 | + * |
| 654 | + * Bitset only: a bitmap selects different rows per query, so there is no single set of |
| 655 | + * rows to compact. |
| 656 | + */ |
| 657 | +template <typename T, typename IdxT, typename BitsT, typename DistanceT = float> |
| 658 | +void brute_force_search_gathered( |
| 659 | + raft::resources const& res, |
| 660 | + const cuvs::neighbors::brute_force::index<T, DistanceT>& idx, |
| 661 | + raft::device_matrix_view<const T, IdxT, raft::row_major> queries, |
| 662 | + const cuvs::core::bitset_view<BitsT, IdxT>& filter_view, |
| 663 | + IdxT n_pass, |
| 664 | + raft::device_matrix_view<IdxT, IdxT, raft::row_major> neighbors, |
| 665 | + raft::device_matrix_view<DistanceT, IdxT, raft::row_major> distances, |
| 666 | + std::optional<raft::device_vector_view<const DistanceT, IdxT>> query_norms = std::nullopt) |
| 667 | +{ |
| 668 | + auto policy = raft::resource::get_thrust_policy(res); |
| 669 | + IdxT n_queries = queries.extent(0); |
| 670 | + IdxT n_dataset = idx.dataset().extent(0); |
| 671 | + IdxT dim = idx.dataset().extent(1); |
| 672 | + IdxT k = neighbors.extent(1); |
| 673 | + |
| 674 | + // 1. enumerate the rows the filter keeps |
| 675 | + auto passing = raft::make_device_vector<IdxT, IdxT>(res, n_pass); |
| 676 | + thrust::copy_if(policy, |
| 677 | + thrust::make_counting_iterator<IdxT>(0), |
| 678 | + thrust::make_counting_iterator<IdxT>(n_dataset), |
| 679 | + passing.data_handle(), |
| 680 | + [filter_view] __device__(IdxT i) { return filter_view.test(i); }); |
| 681 | + |
| 682 | + // 2. compact those rows, and their norms, into a dense dataset |
| 683 | + auto gathered = raft::make_device_matrix<T, IdxT, raft::row_major>(res, n_pass, dim); |
| 684 | + raft::matrix::gather( |
| 685 | + res, |
| 686 | + raft::make_device_matrix_view<const T, IdxT, raft::row_major>( |
| 687 | + idx.dataset().data_handle(), n_dataset, dim), |
| 688 | + raft::make_device_vector_view<const IdxT, IdxT>(passing.data_handle(), n_pass), |
| 689 | + gathered.view()); |
| 690 | + |
| 691 | + auto gathered_norms = |
| 692 | + raft::make_device_vector<DistanceT, IdxT>(res, idx.has_norms() ? n_pass : 0); |
| 693 | + if (idx.has_norms()) { |
| 694 | + thrust::gather(policy, |
| 695 | + passing.data_handle(), |
| 696 | + passing.data_handle() + n_pass, |
| 697 | + idx.norms().data_handle(), |
| 698 | + gathered_norms.data_handle()); |
| 699 | + } |
| 700 | + |
| 701 | + // 3. ordinary unfiltered search over the compacted dataset |
| 702 | + auto compact_neighbors = raft::make_device_matrix<IdxT, IdxT, raft::row_major>(res, n_queries, k); |
| 703 | + std::vector<T*> dataset = {gathered.data_handle()}; |
| 704 | + std::vector<int64_t> sizes = {static_cast<int64_t>(n_pass)}; |
| 705 | + std::vector<DistanceT*> norms; |
| 706 | + if (idx.has_norms()) { norms.push_back(gathered_norms.data_handle()); } |
| 707 | + |
| 708 | + brute_force_knn_impl<int64_t, IdxT, T, DistanceT>( |
| 709 | + res, |
| 710 | + dataset, |
| 711 | + sizes, |
| 712 | + dim, |
| 713 | + const_cast<T*>(queries.data_handle()), |
| 714 | + n_queries, |
| 715 | + compact_neighbors.data_handle(), |
| 716 | + distances.data_handle(), |
| 717 | + k, |
| 718 | + true, |
| 719 | + true, |
| 720 | + nullptr, |
| 721 | + idx.metric(), |
| 722 | + idx.metric_arg(), |
| 723 | + norms.size() ? &norms : nullptr, |
| 724 | + query_norms ? query_norms->data_handle() : nullptr); |
| 725 | + |
| 726 | + // 4. translate compacted row numbers back into the original dataset's numbering |
| 727 | + thrust::gather(policy, |
| 728 | + compact_neighbors.data_handle(), |
| 729 | + compact_neighbors.data_handle() + compact_neighbors.size(), |
| 730 | + passing.data_handle(), |
| 731 | + neighbors.data_handle()); |
| 732 | +} |
| 733 | + |
584 | 734 | template <typename T, typename IdxT, typename BitsT, typename DistanceT = float> |
585 | 735 | void brute_force_search_filtered( |
586 | 736 | raft::resources const& res, |
@@ -619,30 +769,39 @@ void brute_force_search_filtered( |
619 | 769 | const cuvs::core::bitset_view<BitsT, IdxT>>> |
620 | 770 | filter_view; |
621 | 771 |
|
622 | | - IdxT nnz_h = 0; |
623 | | - float sparsity = 0.0f; |
| 772 | + IdxT nnz_h = 0; |
| 773 | + // A bitset passes the same rows for every query, so it has a row count; a bitmap passes a |
| 774 | + // different set per query and has none. |
| 775 | + std::optional<IdxT> n_pass; |
624 | 776 |
|
625 | 777 | const BitsT* filter_data = nullptr; |
626 | 778 |
|
627 | 779 | if (filter_type == cuvs::neighbors::filtering::FilterType::Bitmap) { |
628 | 780 | auto actual_filter = |
629 | 781 | dynamic_cast<const cuvs::neighbors::filtering::bitmap_filter<BitsT, int64_t>*>(filter); |
630 | 782 | filter_view.emplace(actual_filter->view()); |
631 | | - nnz_h = actual_filter->view().count(res); |
632 | | - sparsity = 1.0 - nnz_h / (1.0 * n_queries * n_dataset); |
| 783 | + nnz_h = actual_filter->view().count(res); |
633 | 784 | } else if (filter_type == cuvs::neighbors::filtering::FilterType::Bitset) { |
634 | 785 | auto actual_filter = |
635 | 786 | dynamic_cast<const cuvs::neighbors::filtering::bitset_filter<BitsT, int64_t>*>(filter); |
636 | 787 | filter_view.emplace(actual_filter->view()); |
637 | | - nnz_h = n_queries * actual_filter->view().count(res); |
638 | | - sparsity = 1.0 - nnz_h / (1.0 * n_queries * n_dataset); |
| 788 | + n_pass = actual_filter->view().count(res); |
| 789 | + nnz_h = n_queries * (*n_pass); |
639 | 790 | } else { |
640 | 791 | RAFT_FAIL("Unsupported sample filter type"); |
641 | 792 | } |
642 | 793 |
|
643 | 794 | std::visit([&](const auto& actual_view) { filter_data = actual_view.data(); }, *filter_view); |
644 | 795 |
|
645 | | - if (sparsity < 0.9f) { |
| 796 | + const double selectivity = |
| 797 | + static_cast<double>(nnz_h) / (static_cast<double>(n_queries) * static_cast<double>(n_dataset)); |
| 798 | + const auto path = select_filtered_search_path(n_dataset, dim, selectivity, k, n_pass); |
| 799 | + |
| 800 | + if (path == filtered_search_path::gather) { |
| 801 | + auto bitset_view = std::get<const cuvs::core::bitset_view<BitsT, IdxT>>(*filter_view); |
| 802 | + brute_force_search_gathered<T, IdxT, BitsT, DistanceT>( |
| 803 | + res, idx, queries, bitset_view, *n_pass, neighbors, distances, query_norms); |
| 804 | + } else if (path == filtered_search_path::dense) { |
646 | 805 | raft::resources stream_pool_handle(res); |
647 | 806 | raft::resource::set_cuda_stream(stream_pool_handle, stream); |
648 | 807 | auto idx_norm = idx.has_norms() ? const_cast<DistanceT*>(idx.norms().data_handle()) : nullptr; |
|
0 commit comments