CUDA推力排序或CUB :: DeviceRadixSort

Question

我有一个由float4数组表示的粒子池，其中w分量是粒子当前的生命周期[0，1]。

我需要根据粒子的寿命以降序对这个数组进行排序，这样我就可以为当前有多少粒子“活跃”（寿命大于0）保留一个精确计数器。我需要这个计数器，因为当我需要激活更多的粒子（随机发生）时，它将允许我索引到阵列中的正确位置。

我的微粒数组存储在设备内存中，似乎我应该能够对数组进行排序而不必将数组传输到主机内存。

我在网上找到例子并没有太多的运气，这说明我可以用Thrust或CUB做到这一点。另外，我对使用Thrust犹豫不决，因为我不知道如何防止它退化为合并排序（这比基数排序慢得多），因为我需要根据w组件进行排序。至于CUB，我根本没有找到任何资源来说明我如何做到这一点。

我还希望保存在w组件中的生命周期，因为这使我的生活在我的代码的其他部分更容易。

有没有简单的方法来做到这一点？谢谢你的帮助。

Answer 1

在任一幼崽或推力，我们可以排序的.w“钥匙”只是，那种做一个键值，其中值只是一个线性递增指数：

0, 1, 2, 3, ... 

然后我们可以使用合成重新排列索引序列以便在一个步骤中重新排列原始的float4阵列（有效地按.w排序）。这将允许您保留基数分选速度（无论是幼体还是推力），并且也可能相当有效，因为float4数量只需移动/重新排列一次，而不是在分类操作过程中连续移动。

这里是一个完整的推力示例，在32M元素上，演示了一个“普通”推力排序，使用函数来指定排序.w元素（sort_f4_w），接着是上述方法。在这种情况下，我的特定的设置（Fedora的20，CUDA 7，Quadro5000），第二个方法似乎是约5倍快：

$ cat t686.cu 
#include <iostream> 
#include <vector_types.h> 
#include <stdlib.h> 
#include <thrust/host_vector.h> 
#include <thrust/device_vector.h> 
#include <thrust/sort.h> 
#include <thrust/iterator/transform_iterator.h> 
#include <thrust/iterator/permutation_iterator.h> 
#include <thrust/sequence.h> 
#include <thrust/copy.h> 
#include <thrust/equal.h> 

#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; 
} 

#define DSIZE (32*1048576) 

struct sort_f4_w 
{ 
    __host__ __device__ 
    bool operator()(const float4 &a, const float4 &b) const { 
    return (a.w < b.w);} 
}; 
// functor to extract the .w element from a float4 
struct f4_to_fw : public thrust::unary_function<float4, float> 
{ 
    __host__ __device__ 
    float operator()(const float4 &a) const { 
    return a.w;} 
}; 
// functor to extract the .x element from a float4 
struct f4_to_fx : public thrust::unary_function<float4, float> 
{ 
    __host__ __device__ 
    float operator()(const float4 &a) const { 
    return a.x;} 
}; 


bool validate(thrust::device_vector<float4> &d1, thrust::device_vector<float4> &d2){ 
    return thrust::equal(thrust::make_transform_iterator(d1.begin(), f4_to_fx()), thrust::make_transform_iterator(d1.end(), f4_to_fx()), thrust::make_transform_iterator(d2.begin(), f4_to_fx())); 
} 


int main(){ 
    unsigned long long t1_time, t2_time; 
    float4 *mydata = new float4[DSIZE]; 
    for (int i = 0; i < DSIZE; i++){ 
    mydata[i].x = i; 
    mydata[i].y = i; 
    mydata[i].z = i; 
    mydata[i].w = rand()/(float)RAND_MAX;} 

    thrust::host_vector<float4> h_data(mydata, mydata+DSIZE); 
    // do once as a warm-up run, then report timings on second run 
    for (int i = 0; i < 2; i++){ 
    thrust::device_vector<float4> d_data1 = h_data; 
    thrust::device_vector<float4> d_data2 = h_data; 

    // first time sort using typical thrust approach 
    t1_time = dtime_usec(0); 
    thrust::sort(d_data1.begin(), d_data1.end(), sort_f4_w()); 
    cudaDeviceSynchronize(); 
    t1_time = dtime_usec(t1_time); 
    // now extract keys and create index values, sort, then rearrange 
    t2_time = dtime_usec(0); 
    thrust::device_vector<float> keys(DSIZE); 
    thrust::device_vector<int> vals(DSIZE); 
    thrust::copy(thrust::make_transform_iterator(d_data2.begin(), f4_to_fw()), thrust::make_transform_iterator(d_data2.end(), f4_to_fw()), keys.begin()); 
    thrust::sequence(vals.begin(), vals.end()); 
    thrust::sort_by_key(keys.begin(), keys.end(), vals.begin()); 
    thrust::device_vector<float4> result(DSIZE); 
    thrust::copy(thrust::make_permutation_iterator(d_data2.begin(), vals.begin()), thrust::make_permutation_iterator(d_data2.begin(), vals.end()), result.begin()); 
    cudaDeviceSynchronize(); 
    t2_time = dtime_usec(t2_time); 
    if (!validate(d_data1, result)){ 
     std::cout << "Validation failure " << std::endl; 
     } 
    } 
    std::cout << "thrust t1 time: " << t1_time/(float)USECPSEC << "s, t2 time: " << t2_time/(float)USECPSEC << std::endl; 
} 


$ nvcc -o t686 t686.cu 
$ ./t686 
thrust t1 time: 0.731456s, t2 time: 0.149959 
$