# Help how to implement this algorithm in CUDA

Hi guys,

I am doing a comparison between CUDA, MPI, OpenMP and Pthreads.

I implemented the algorithm bitonic MergeSort in Pthread, MPI and OpenMP.

I try to implement it in CUDA, but I no have sucess, because device and global methods in CUDA have restrictions ( per example memcpy function in device and the limit of 16KB per block of threads)

Can you help me?

Bellow the algorithm in Pthreads and more bellow in CUDA

[codebox]#include <stdio.h>

#include <stdlib.h>

#include <sys/time.h> /* tempo de execuÃ§Ã£o */

#include <string.h> /* memcpy */

/*#define LARGURA 524288 //2^19

#define ELEMENTOS 131072

#define LARGURA 262144 //2^18

#define ELEMENTOS 65536

#define LARGURA 131072 //2^17

#define ELEMENTOS 32768

#define LARGURA 65536 //2^16

#define ELEMENTOS 16384

#define LARGURA 32768 //2^15

#define ELEMENTOS 8192

#define LARGURA 16384 //2^14

#define ELEMENTOS 4096

#define LARGURA 8192 //2^13

#define ELEMENTOS 2048

#define LARGURA 4096 //2^12

#define ELEMENTOS 1024

#define LARGURA 2048 //2^11

#define ELEMENTOS 512

#define LARGURA 1024 //2^10

#define ELEMENTOS 256

*/

// Tamanho para teste simples

#define LARGURA 8

#define ELEMENTOS 2

// Debug

#define DEBUG 0

int vetor[LARGURA];

void quickSort_it (int v, int p, int r)

{

``````int j, *pilhap, *pilhar, t;

pilhap = malloc ((r-p+1) * sizeof (int));

pilhar = malloc ((r-p+1) * sizeof (int));

pilhap[0] = p;

pilhar[0] = r;

t = 0;

while (t >= 0)

{

p = pilhap[t];

r = pilhar[t];

t--;

if (p < r)

{

j = separa (v, p, r);

t++;

pilhap[t] = p;

pilhar[t] = j-1;

t++;

pilhap[t] = j+1;

pilhar[t] = r;

}

}
``````

}

int separa (int v, int p, int r)

{

``````int c = v[p], i = p+1, j = r, t;

while ( i <= j)

{

if (v[i] <= c) ++i;

else if (c < v[j]) --j;

else

{

t = v[i], v[i] = v[j], v[j] = t;

++i; --j;

}

}

v[p] = v[j], v[j] = c;

return j;
``````

}

void carregaVetor()

{

``````FILE *arquivo;

int i=0, num;

arquivo = fopen("randomicos.txt", "r");

for(i=0; i < LARGURA; i++)

{

fscanf (arquivo, "%d", &num);

vetor[i] = num;

}

fclose(arquivo);
``````

}

void imprimeVetor(){

``````int i;

for(i=0; i < LARGURA; i++){

printf("%d\n", vetor[i]);

}
``````

}

void *bitonicMergesort(void *args)

{

``````int i, j, k, l, numprocs, myid, aux;

numprocs = (LARGURA / ELEMENTOS);

myid = (int)args;

int vetorLocal[ELEMENTOS * 2];

quickSort_it(vetor, (myid * ELEMENTOS), (myid * ELEMENTOS) + ELEMENTOS -1);

for (k = 1; k < numprocs / 2; k = k * 2)

{//Controle da fase

if ((myid % (k * 4)) < (k * 2))

{//ORDENA CRESCENTE

if (myid % (k * 2) < k)

{ //Processa

memcpy(vetorLocal,&vetor[myid * ELEMENTOS],sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS],&vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0,(ELEMENTOS *2)-1);

memcpy(&vetor[myid * ELEMENTOS],vetorLocal ,sizeof(int) * ELEMENTOS);

memcpy(&vetor[(myid + k) * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

}

else

{ //ORDENA DECRESCENTE

if (myid % (k * 2) < k)

{ //Processa

memcpy(vetorLocal,&vetor[myid * ELEMENTOS],sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS],&vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0,(ELEMENTOS *2)-1);

memcpy(&vetor[(myid + k) * ELEMENTOS],vetorLocal ,sizeof(int) * ELEMENTOS);

memcpy(&vetor[myid * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

}

if (k > 1)

{ //Acerta os Valores

for (l = k / 2; l >= 1;  l = l / 2)

{

if (myid % (k * 4) < (k * 2))

{ //Crescente

if (myid % (l * 2) < l)

{ //Processa

memcpy(vetorLocal,&vetor[myid * ELEMENTOS],sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS], &vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0, (ELEMENTOS * 2) -1);

memcpy(&vetor[myid * ELEMENTOS],vetorLocal ,sizeof(int) * ELEMENTOS);

memcpy(&vetor[(myid + k) * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

}

else{ //Decrescente

if (myid % (l * 2) < l)

{ //Processa

memcpy(vetorLocal,&vetor[myid * ELEMENTOS],sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS], &vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0, (ELEMENTOS * 2) -1);

memcpy(&vetor[(myid + k) * ELEMENTOS],vetorLocal ,sizeof(int) * ELEMENTOS);

memcpy(&vetor[myid * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

}

}

}

}

//FASE 2 - Bitonico

for (k = numprocs/2; k >= 1; k = k / 2)

{ //Controle da fase

if (myid % (k * 2) < k)

{ //Processa

memcpy(vetorLocal,&vetor[myid * ELEMENTOS],sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS],&vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0, (ELEMENTOS * 2) -1);

memcpy(&vetor[myid * ELEMENTOS],vetorLocal ,sizeof(int) * ELEMENTOS);

memcpy(&vetor[(myid + k) * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

}
``````

}

void escreveTempos(double tempo[5]){

``````FILE *arquivo;

int i;

for(i=0; i < 5; i++)

fprintf(arquivo, "-> %d itens: %lf s\n",LARGURA, tempo[i]);

fclose(arquivo);
``````

}

int main ()

{

``````double diferenca[5];

struct timeval tvAntes, tvDepois;

int i,j;

for (j=0; j < 5; j++)

{

carregaVetor();

if(gettimeofday(&tvAntes,NULL)<0)

perror("gettimeofday");

for (i = 0; i < (LARGURA / ELEMENTOS); i++)

for (i = 0; i < (LARGURA / ELEMENTOS); i++)

if(gettimeofday(&tvDepois,NULL)<0)

perror("gettimeofday");

diferenca[j] = (((double)tvDepois.tv_sec) - ((double)tvAntes.tv_sec)) + ((((double)tvDepois.tv_usec) - ((double)tvAntes.tv_usec)) / 1000000);

}

escreveTempos(diferenca);

#if DEBUG == 1

imprimeVetor();

#endif

return 0;
``````

}

[/codebox]

My CUDA implementation

[codebox]

#ifndef BITONIC_KERNEL_CU

#define BITONIC_KERNEL_CU

// Tamanho dos vetores utilizados nos testes

/*#define LARGURA 524288 //2^19

#define ELEMENTOS 131072

#define LARGURA 262144 //2^18

#define ELEMENTOS 65536

#define LARGURA 131072 //2^17

#define ELEMENTOS 32768

#define LARGURA 65536 //2^16

#define ELEMENTOS 16384

#define LARGURA 32768 //2^15

#define ELEMENTOS 8192

#define LARGURA 16384 //2^14

#define ELEMENTOS 4096

#define LARGURA 8192 //2^13

#define ELEMENTOS 2048

#define LARGURA 4096 //2^12

#define ELEMENTOS 1024

#define LARGURA 2048 //2^11

#define ELEMENTOS 512 */

//#define LARGURA 512//2^10

//#define ELEMENTOS 128

//#define LARGURA 256//2^10

//#define ELEMENTOS 64

//#define LARGURA 8 //2^10

//#define ELEMENTOS 2

#define LARGURA 128 //2^10

#define ELEMENTOS 32

//#define LARGURA 8 //2^10

//#define ELEMENTOS 2

/**

*/

device static int separa (int v, int p, int r)

{

``````int c = v[p], i = p+1, j = r, t;

while (i <= j)

{

if (v[i] <= c) ++i;

else if (c < v[j]) --j;

else {

t = v[i], v[i] = v[j], v[j] = t;

++i; --j;

}

}

v[p] = v[j], v[j] = c;

return j;
``````

}

/**

• Algoritmo quicksort iterativo

*/

device static void quickSort_it(int v, int p, int r)

{

``````int j, t;

extern __shared__ int pilhap[];

extern __shared__ int pilhar[];

pilhap[0] = p;

pilhar[0] = r;

t = 0;

while (t >= 0) {

p = pilhap[t];

r = pilhar[t];

t--;

if (p < r) {

j = separa(v, p, r);

t++;

pilhap[t] = p;

pilhar[t] = j - 1;

t++;

pilhap[t] = j + 1;

pilhar[t] = r;

}

}

pilhap[0] = p;

pilhar[0] = r;

t = 0;

while (t >= 0) {

p = pilhap[t];

r = pilhar[t];

t--;

while (p < r) {

j = separa(v, p, r);

t++;

pilhap[t] = p;

pilhar[t] = j - 1;

p = j + 1;

}

}
``````

}

/**

• Bitonic MergeSort

*/

global static void bitonicMergeSort(int * vetor)

{

``````const unsigned int myid = threadIdx.x;

unsigned int k, l, numprocs;

unsigned int vetorLocal[ELEMENTOS * 2];

extern __shared__ int vetorLocal[];

numprocs = (LARGURA / ELEMENTOS);

quickSort_it(vetor, (myid * ELEMENTOS), (myid * ELEMENTOS) + ELEMENTOS - 1);

for (k = 1; k < numprocs / 2; k = k * 2)

{

//Controle da fase

if ((myid % (k * 4)) < (k * 2))

{

//ORDENA CRESCENTE

if (myid % (k * 2) < k)

{

//Processa

memcpy(vetorLocal, &vetor[myid * ELEMENTOS], sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS], &vetor[(myid + k) * ELEMENTOS],sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0, (ELEMENTOS * 2) - 1);

memcpy(&vetor[myid * ELEMENTOS], vetorLocal, sizeof(int) * ELEMENTOS);

memcpy(&vetor[(myid + k) * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

} else

{

//ORDENA DECRESCENTE

if (myid % (k * 2) < k)

{

//Processa

memcpy(vetorLocal, &vetor[myid * ELEMENTOS], sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS], &vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0, (ELEMENTOS * 2) - 1);

memcpy(&vetor[(myid + k) * ELEMENTOS], vetorLocal, sizeof(int) * ELEMENTOS);

memcpy(&vetor[myid * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

}

if (k > 1)

{

//Acerta os Valores

for (l = k / 2; l >= 1; l = l / 2)

{

if (myid % (k * 4) < (k * 2))

{

//Crescente

if (myid % (l * 2) < l)

{

//Processa

memcpy(vetorLocal, &vetor[myid * ELEMENTOS], sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS], &vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0, (ELEMENTOS * 2) - 1);

memcpy(&vetor[myid * ELEMENTOS], vetorLocal, sizeof(int) * ELEMENTOS);

memcpy(&vetor[(myid + k) * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

} else {

//Decrescente

if (myid % (l * 2) < l) {

//Processa

memcpy(vetorLocal, &vetor[myid * ELEMENTOS], sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS], &vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetorLocal, 0, (ELEMENTOS * 2) - 1);

memcpy(&vetor[(myid + k) * ELEMENTOS], vetorLocal, sizeof(int) * ELEMENTOS);

memcpy(&vetor[myid * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}

}

}

}

}

//FASE 2 - Bitonico

for (k = numprocs / 2; k >= 1; k = k / 2) {

//Controle da fase

if (myid % (k * 2) < k)

{

//Processa

memcpy(vetorLocal, &vetor[myid * ELEMENTOS], sizeof(int) * ELEMENTOS);

memcpy(&vetorLocal[ELEMENTOS], &vetor[(myid + k) * ELEMENTOS], sizeof(int) * ELEMENTOS);

quickSort_it(vetor, 0, (ELEMENTOS * 2) - 1);

memcpy(&vetor[myid * ELEMENTOS], vetorLocal, sizeof(int) * ELEMENTOS);

memcpy(&vetor[(myid + k) * ELEMENTOS], &vetorLocal[ELEMENTOS], sizeof(int) * ELEMENTOS);

}