CUDA GPU上的小尺寸二進制搜索

Question

我有一個int64_t類型的大型設備陣列inputValues。此數組中的每個32個元素都按升序排序。我有一個未排序的搜索數組removeValues。

我的意圖是尋找removeValues中inputValues的所有元素，並將它們標記爲-1。什麼是最有效的方法來實現這一點？如果有幫助，我正在使用3.5 cuda設備。

我不是在尋找更高層次的解決方案，即我不想使用推力或幼崽，但是我想用cuda內核編寫它。

我最初的做法是將每個32值加載到線程塊中的共享內存中。每個線程還從removeValues加載單個值，並在共享內存陣列上執行獨立的二進制搜索。如果找到，則根據if條件設置該值。

這種方法不會涉及很多銀行衝突和分支機構分歧嗎？你認爲在實現二分搜索時可以通過使用三元運算符來解決分支分歧嗎？即使解決了這個問題，銀行衝突如何消除？由於排序數組的大小爲32，是否可以使用混洗指令實現二分搜索？這會有幫助嗎？

編輯：我已經添加了一個示例來顯示我打算實現什麼。

假設inputValues是其中每32個元件被排序的矢量：
[2, 4, 6, ... , 64], [95, 97, ... , 157], [1, 3, ... , 63], [...]

此陣列的典型尺寸範圍可以32之間32 * 2〜* 32。值可以從0到INT64_MAX。

的removeValues一個例子是：
[7, 75, 95, 106]

此陣列的典型尺寸範圍可以從1到1024。

操作removeValues之後將是： [-1，75，-1，106]

在inputValues的值保持不變。

Answer 1

我會同意答案（現在刪除）和@harrism評論。由於我對非推力方法付出了一些努力，因此我將介紹我的發現。

我試圖天真地使用__shfl()在warp層次上實現二分搜索，然後在整個數據集上重複該二分搜索，將數據集通過每個32元素組。

這很令人尷尬，但我的代碼比推力要慢20倍左右（實際上如果你用nvprof仔細計時，可能會更糟）。

由於問題中的數據大小非常小，以至於計時不準確，因此我將數據大小稍微大於問題中提出的數據大小。

這裏有2種方式完全樣例：

1. 中有什麼問題大致概括，即使用經洗牌創建一個二進制搜索，最多可以搜索到32元對32元有序陣列。重複這個過程爲32個元素的有序數組，然後通過每個有序數組傳遞整個數據集（希望你現在可以開始看到一些低效率。）

2. 使用推力，基本上與@harrism概述了什麼，即對分組數據集進行排序，然後在其上運行矢量化thrust::binary_search。

這裏的例子：

$ cat t1030.cu 
#include <stdio.h> 
#include <assert.h> 
#include <thrust/host_vector.h> 
#include <thrust/device_vector.h> 
#include <thrust/sort.h> 
#include <thrust/binary_search.h> 

typedef long mytype; 

const int gsize = 32; 
const int nGRP = 512; 
const int dsize = nGRP*gsize;//gsize*nGRP; 

#include <time.h> 
#include <sys/time.h> 
#define USECPSEC 1000000ULL 

unsigned long long dtime_usec(unsigned long long start){ 

    timeval tv; 
    gettimeofday(&tv, 0); 
    return ((tv.tv_sec*USECPSEC)+tv.tv_usec)-start; 
} 

template <typename T> 
__device__ T my_shfl32(T val, unsigned lane){ 
    return __shfl(val, lane); 
} 

template <typename T> 
__device__ T my_shfl64(T val, unsigned lane){ 
    T retval = val; 
    int2 t1 = *(reinterpret_cast<int2 *>(&retval)); 
    t1.x = __shfl(t1.x, lane); 
    t1.y = __shfl(t1.y, lane); 
    retval = *(reinterpret_cast<T *>(&t1)); 
    return retval; 
} 

template <typename T> 
__device__ bool bsearch_shfl(T grp_val, T my_val){ 
    int src_lane = gsize>>1; 
    bool return_val = false; 
    T test_val; 
    int shift = gsize>>2; 
    for (int i = 0; i <= gsize>>3; i++){ 
    if (sizeof(T)==4){ 
     test_val = my_shfl32(grp_val, src_lane);} 
    else if (sizeof(T)==8){ 
     test_val = my_shfl64(grp_val, src_lane);} 
    else assert(0); 
    if (test_val == my_val) return_val = true; 
    src_lane += (((test_val<my_val)*2)-1)*shift; 
    shift>>=1; 
    assert ((src_lane < gsize)&&(src_lane > 0));} 
    if (sizeof(T)==4){ 
    test_val = my_shfl32(grp_val, 0);} 
    else if (sizeof(T)==8){ 
    test_val = my_shfl64(grp_val, 0);} 
    else assert(0); 
    if (test_val == my_val) return_val = true; 
    return return_val; 
} 

template <typename T> 
__global__ void bsearch_grp(const T * __restrict__ search_grps, T *data){ 

    int idx = threadIdx.x+blockDim.x*blockIdx.x; 
    int tid = threadIdx.x; 
    if (idx < gsize*nGRP){ 
    T grp_val = search_grps[idx]; 
    while (tid < dsize){ 
     T my_val = data[tid]; 
     if (bsearch_shfl(grp_val, my_val)) data[tid] = -1; 
     tid += blockDim.x;} 
    } 
} 


int main(){ 

    // data setup 
    assert(gsize == 32); //mandatory (warp size) 
    assert((dsize % 32)==0); //needed to preserve shfl capability 
    thrust::host_vector<mytype> grps(gsize*nGRP); 
    thrust::host_vector<mytype> data(dsize); 
    thrust::host_vector<mytype> result(dsize); 
    for (int i = 0; i < gsize*nGRP; i++) grps[i] = i; 
    for (int i = 0; i < dsize; i++) data[i] = i; 
    // method 1: individual shfl-based binary searches on each group 
    mytype *d_grps, *d_data; 
    cudaMalloc(&d_grps, gsize*nGRP*sizeof(mytype)); 
    cudaMalloc(&d_data, dsize*sizeof(mytype)); 
    cudaMemcpy(d_grps, &(grps[0]), gsize*nGRP*sizeof(mytype), cudaMemcpyHostToDevice); 
    cudaMemcpy(d_data, &(data[0]), dsize*sizeof(mytype), cudaMemcpyHostToDevice); 
    unsigned long long my_time = dtime_usec(0); 
    bsearch_grp<<<nGRP, gsize>>>(d_grps, d_data); 
    cudaDeviceSynchronize(); 
    my_time = dtime_usec(my_time); 
    cudaMemcpy(&(result[0]), d_data, dsize*sizeof(mytype), cudaMemcpyDeviceToHost); 
    for (int i = 0; i < dsize; i++) if (result[i] != -1) {printf("method 1 mismatch at %d, was %d, should be -1\n", i, (int)(result[i])); return 1;} 
    printf("method 1 time: %fs\n", my_time/(float)USECPSEC); 
    // method 2: thrust sort, followed by thrust binary search 
    thrust::device_vector<mytype> t_grps = grps; 
    thrust::device_vector<mytype> t_data = data; 
    thrust::device_vector<bool> t_rslt(t_data.size()); 
    my_time = dtime_usec(0); 
    thrust::sort(t_grps.begin(), t_grps.end()); 
    thrust::binary_search(t_grps.begin(), t_grps.end(), t_data.begin(), t_data.end(), t_rslt.begin()); 
    cudaDeviceSynchronize(); 
    my_time = dtime_usec(my_time); 
    thrust::host_vector<bool> rslt = t_rslt; 
    for (int i = 0; i < dsize; i++) if (rslt[i] != true) {printf("method 2 mismatch at %d, was %d, should be 1\n", i, (int)(rslt[i])); return 1;} 
    printf("method 2 time: %fs\n", my_time/(float)USECPSEC); 

    // method 3: multiple thrust merges, followed by thrust binary search 



    return 0; 
} 

$ nvcc -O3 -arch=sm_35 t1030.cu -o t1030 
$ ./t1030 
method 1 time: 0.009075s 
method 2 time: 0.000516s 
$ 

我是在Linux上運行此，CUDA 7.5，GT640 GPU。顯然，不同GPU上的性能會有所不同，但如果任何GPU顯着縮小差距，我會感到驚訝。

總之，你最好建議使用一個調整好的庫，如推力或幼崽。如果你不喜歡推力的單一性，你可以嘗試幼崽。我不知道Cub是否有二進制搜索，但是對整個分類數據集的單個二分搜索是not a difficult thing to write，這是涉及的時間的較小部分（對於方法2 - 使用nvprof或附加的時間碼可識別）。

由於你的32個元素的分組範圍已經排序，我還考慮了使用多個thrust::merge操作而不是單一排序的想法。我不確定哪個會更快，但由於推力方法已經比32元素洗牌搜索方法快得多，我認爲推力（或幼仔）是明顯的選擇。