Ordenar 3 matrices por clave en CUDA (usando Thrust quizás)

Question

Tengo 3 matrices del mismo tamaño (más de 300.000 elementos). Una matriz de números flotantes y dos matrices de índices. Entonces, para cada número tengo 2 IDs.

Todas las matrices 3 ya se encuentran en la memoria global de la GPU. Quiero ordenar todos los números con sus IDs en consecuencia.

¿Hay alguna manera en que pueda usar la biblioteca Thrust para realizar esta tarea? ¿Hay alguna manera mejor que la biblioteca Thrust?

Por supuesto, prefiero no copiarlos ay desde la memoria del host un par de veces. Por cierto, son matrices, no vectores.

Gracias por su ayuda con anticipación.

---

solución provisional, pero esto es extremadamente lento. Se tarda casi 4 segundo y mi tamaño de la matriz es con el fin de 300000

thrust::device_ptr<float> keys(afterSum); 
thrust::device_ptr<int> vals0(d_index); 
thrust::device_ptr<int> vals1(blockId); 

thrust::device_vector<int> sortedIndex(numElements); 
thrust::device_vector<int> sortedBlockId(numElements); 

thrust::counting_iterator<int> iter(0); 
thrust::device_vector<int> indices(numElements); 
thrust::copy(iter, iter + indices.size(), indices.begin()); 

thrust::sort_by_key(keys, keys + numElements , indices.begin());  

thrust::gather(indices.begin(), indices.end(), vals0, sortedIndex.begin()); 
thrust::gather(indices.begin(), indices.end(), vals1, sortedBlockId.begin()); 

thrust::host_vector<int> h_sortedIndex=sortedIndex; 
thrust::host_vector<int> h_sortedBlockId=sortedBlockId; 

Answer 1

Por supuesto se puede utilizar de empuje. Primero, debe envolver los punteros de su dispositivo CUDA sin procesar con thrust::device_ptr. Asumiendo que su valores flotantes están en la matriz pkeys y los ID están en las matrices pvals0 y pvals1 y numElements es la longitud de las matrices, algo como esto debería funcionar:

#include <thrust/device_ptr.h> 
#include <thrust/sort.h> 
#include <thrust/gather.h> 
#include <thrust/iterator/counting_iterator.h> 

cudaEvent_t start, stop; 
cudaEventCreate(&start); 
cudaEventCreate(&stop); 

cudaEventRecord(start); 

thrust::device_ptr<float> keys(pkeys); 
thrust::device_ptr<int> vals0(pvals0); 
thrust::device_ptr<int> vals1(pvals1); 

// allocate space for the output 
thrust::device_vector<int> sortedVals0(numElements); 
thrust::device_vector<int> sortedVals1(numElements); 

// initialize indices vector to [0,1,2,..] 
thrust::counting_iterator<int> iter(0); 
thrust::device_vector<int> indices(numElements); 
thrust::copy(iter, iter + indices.size(), indices.begin()); 

// first sort the keys and indices by the keys 
thrust::sort_by_key(keys.begin(), keys.end(), indices.begin()); 

// Now reorder the ID arrays using the sorted indices 
thrust::gather(indices.begin(), indices.end(), vals0.begin(), sortedVals0.begin()); 
thrust::gather(indices.begin(), indices.end(), vals1.begin(), sortedVals1.begin()); 

cudaEventRecord(stop); 
cudaEventSynchronize(stop); 
float milliseconds = 0; 
cudaEventElapsedTime(&milliseconds, start, stop); 
printf("Took %f milliseconds for %d elements\n", milliseconds, numElements); 

Answer 2

me gustaría utilizar zip_iterator para realizar una sort_by_key en ambos vectores de índice al mismo tiempo.

Este sería el siguiente:

typedef typename thrust::tuple<thrust::device_vector<int>::iterator, thrust::device_vector<int>::iterator> IteratorTuple; 
    typedef typename thrust::zip_iterator<IteratorTuple> ZipIterator; 

    // here I suppose your 3 arrays are pointed to by device_ptr as suggested by @harrism 
    thrust::device_vector<float> key(pKey, pKey + numElements); 
    thrust::device_vector<int> val0(pVal0, pVal0 + numElements); 
    thrust::device_vector<int> val1(pVal1, pVal1 + numElements); 

    ZipIterator iterBegin(thrust::make_tuple(val0.begin(), val1.begin())); 
    thrust::sort_by_key(key.begin(), key.end(), iterBegin); 

Answer 3

He comparado los dos enfoques propuestos anteriormente, a saber, que el uso de thrust::zip_iterator y que el uso de thrust::gather. Los he probado en el caso de ordenar dos matrices por clave o tres matrices, según lo solicitado por el afiche. En todos los casos, el enfoque que usa thrust::gather ha demostrado ser más rápido.

EL CASO DE 2 ARRAYS

#include <time.h>  // --- time 
#include <stdlib.h>  // --- srand, rand 

#include <thrust\host_vector.h> 
#include <thrust\device_vector.h> 
#include <thrust\sort.h> 
#include <thrust\iterator\zip_iterator.h> 

#include "TimingGPU.cuh" 

//#define VERBOSE 
//#define COMPACT 

int main() { 

    const int N = 1048576; 
    //const int N = 10; 

    TimingGPU timerGPU; 

    // --- Initialize random seed 
    srand(time(NULL)); 

    thrust::host_vector<int> h_code(N); 
    thrust::host_vector<double> h_x(N); 
    thrust::host_vector<double> h_y(N); 

    for (int k = 0; k < N; k++) {  
     // --- Generate random numbers between 0 and 9 
     h_code[k] = rand() % 10 + 1; 
     h_x[k] = ((double)rand()/(RAND_MAX)); 
     h_y[k] = ((double)rand()/(RAND_MAX)); 
    } 

    thrust::device_vector<int> d_code(h_code); 

    thrust::device_vector<double> d_x(h_x); 
    thrust::device_vector<double> d_y(h_y); 

#ifdef VERBOSE 
    printf("Before\n"); 
    for (int k = 0; k < N; k++) printf("code = %i; x = %f; y = %f\n", h_code[k], h_x[k], h_y[k]); 
#endif 

    timerGPU.StartCounter(); 
#ifdef COMPACT 
    thrust::sort_by_key(d_code.begin(), d_code.end(), thrust::make_zip_iterator(thrust::make_tuple(d_x.begin(), d_y.begin()))); 
#else 

    // --- Initialize indices vector to [0,1,2,..] 
    thrust::counting_iterator<int> iter(0); 
    thrust::device_vector<int> indices(N); 
    thrust::copy(iter, iter + indices.size(), indices.begin()); 

    // --- First, sort the keys and indices by the keys 
    thrust::sort_by_key(d_code.begin(), d_code.end(), indices.begin()); 

    // Now reorder the ID arrays using the sorted indices 
    thrust::gather(indices.begin(), indices.end(), d_x.begin(), d_x.begin()); 
    thrust::gather(indices.begin(), indices.end(), d_y.begin(), d_y.begin()); 
#endif 

    printf("Timing GPU = %f\n", timerGPU.GetCounter()); 

#ifdef VERBOSE 
    h_code = d_code; 
    h_x = d_x; 
    h_y = d_y; 

    printf("After\n"); 
    for (int k = 0; k < N; k++) printf("code = %i; x = %f; y = %f\n", h_code[k], h_x[k], h_y[k]); 
#endif 
} 

EL CASO DE 3 ARRAYS

#include <time.h>  // --- time 
#include <stdlib.h>  // --- srand, rand 

#include <thrust\host_vector.h> 
#include <thrust\device_vector.h> 
#include <thrust\sort.h> 
#include <thrust\iterator\zip_iterator.h> 

#include "TimingGPU.cuh" 

//#define VERBOSE 
//#define COMPACT 

int main() { 

    const int N = 1048576; 
    //const int N = 10; 

    TimingGPU timerGPU; 

    // --- Initialize random seed 
    srand(time(NULL)); 

    thrust::host_vector<int> h_code(N); 
    thrust::host_vector<double> h_x(N); 
    thrust::host_vector<double> h_y(N); 
    thrust::host_vector<double> h_z(N); 

    for (int k = 0; k < N; k++) { 
     // --- Generate random numbers between 0 and 9 
     h_code[k] = rand() % 10 + 1; 
     h_x[k] = ((double)rand()/(RAND_MAX)); 
     h_y[k] = ((double)rand()/(RAND_MAX)); 
     h_z[k] = ((double)rand()/(RAND_MAX)); 
    } 

    thrust::device_vector<int> d_code(h_code); 

    thrust::device_vector<double> d_x(h_x); 
    thrust::device_vector<double> d_y(h_y); 
    thrust::device_vector<double> d_z(h_z); 

#ifdef VERBOSE 
    printf("Before\n"); 
    for (int k = 0; k < N; k++) printf("code = %i; x = %f; y = %f\n", h_code[k], h_x[k], h_y[k]); 
#endif 

    timerGPU.StartCounter(); 
#ifdef COMPACT 
    thrust::sort_by_key(d_code.begin(), d_code.end(), thrust::make_zip_iterator(thrust::make_tuple(d_x.begin(), d_y.begin(), d_z.begin()))); 
#else 

    // --- Initialize indices vector to [0,1,2,..] 
    thrust::counting_iterator<int> iter(0); 
    thrust::device_vector<int> indices(N); 
    thrust::copy(iter, iter + indices.size(), indices.begin()); 

    // --- First, sort the keys and indices by the keys 
    thrust::sort_by_key(d_code.begin(), d_code.end(), indices.begin()); 

    // Now reorder the ID arrays using the sorted indices 
    thrust::gather(indices.begin(), indices.end(), d_x.begin(), d_x.begin()); 
    thrust::gather(indices.begin(), indices.end(), d_y.begin(), d_y.begin()); 
    thrust::gather(indices.begin(), indices.end(), d_z.begin(), d_z.begin()); 
#endif 

    printf("Timing GPU = %f\n", timerGPU.GetCounter()); 

#ifdef VERBOSE 
    h_code = d_code; 
    h_x = d_x; 
    h_y = d_y; 

    printf("After\n"); 
    for (int k = 0; k < N; k++) printf("code = %i; x = %f; y = %f\n", h_code[k], h_x[k], h_y[k]); 
#endif 
} 

de temporización en el caso de 2 arrays para N = 1048576

zip_iterator = 7.34ms 
gather  = 4.27ms 

Momento en el caso de las matrices de 3N = 1048576

zip_iterator = 9.64ms 
gather  = 4.22ms 

Las pruebas realizadas en una GTX 960 tarjeta de NVIDIA.