Test Millera-Rabina - wykorzystanie cuda.

Tak ostatnio zabrałem się za wykorzystanie rdzeni cuda w obliczeniach. Poczytałem na ten temat i postanowiłem przetestować, to co wiem na przykładzie algorytmu sprawdzającego czy liczba jest liczbą pierwszą. Teoretycznie wszystko działa, ale wzrost wydajności praktycznie żaden. Oryginalny algorytm został zaczerpnięty stąd: http://edu.i-lo.tarnow.pl/inf/alg/001_search/0019.php

Moje przełożenie algorytmu, aby wykorzystał GPU:

#include "cuda.h"
#include "../Cuda_Examples/common/book.h"
#include <time.h>


#define 	N	4000
#define 	B	4


__global__ void setAll(bool *b){
	b[threadIdx.x + blockIdx.x * blockDim.x] = true;
}

__host__ ulong Random(ulong a, ulong b){
	ulong w;
	for(int i = 1; i<=8; i++){
		w <<= 8;
		w |= rand() % 256;
	}
	return a + (w % (b-a));
}


__device__ ulong MulModulo(ulong a, ulong b, ulong n){
	ulong w = 0;
	for(ulong m = 1; m; m<<=1){
		if(b & m) w = (w + a) % n;
		a = (a << 1) % n;
	}
	return w;
}


__device__ ulong PowModulo(ulong a, ulong e, ulong n){
	ulong p, w;
	p = a; w = 1;
	for(ulong m = 1; m; m<<=1){
		if(e & m) w = MulModulo(w, p ,n);
		p = MulModulo(p,p,n);
	}
	return w;
}


__global__ void check(ulong v, ulong *r, bool *b){
	ulong d, x;
	int s = 0;
	int index = threadIdx.x + blockIdx.x * blockDim.x;

	for(d = v - 1; d % 2 == 0; s++) d /= 2;

	x = PowModulo(r[index], d, v);
	if(!((x == 1)||(x == v - 1))){
		for(int j = 1; (j < s) && (x != v - 1); j++){
			x = MulModulo(x, x, v);
			if(x == 1){
				b[index] = false; break;
			}
		}
		if(x != v - 1){
			b[index] = false;
		}
	}
}


int main(void){
	cudaEvent_t start, stop;
	cudaEventCreate(&start);
	cudaEventCreate(&stop);

	ulong value = 0;

	printf("Wprowadź liczbę: ");
	scanf("%lu", &value);

	cudaEventRecord(start, 0);

	srand((unsigned)time(NULL));
	ulong *random = (ulong*)malloc(sizeof(ulong)*N);
	bool *bool_f = (bool*)malloc(sizeof(bool)*N);


	bool *dev_b;
	ulong *dev_r;

	for(int i = 0; i<N; i++){
		random[i] = Random(2, value - 2);
	}

	HANDLE_ERROR(cudaMalloc((void**)&dev_r, N * sizeof(ulong)));
	HANDLE_ERROR(cudaMalloc((void**)&dev_b, N * sizeof(bool)));
	HANDLE_ERROR(cudaMemcpy(dev_r, random, N * sizeof(ulong), cudaMemcpyHostToDevice));
	HANDLE_ERROR(cudaMemcpy(dev_b, bool_f, N * sizeof(bool), cudaMemcpyHostToDevice));

	setAll<<<B, N/B>>>(dev_b);

	check<<<B, N/B>>>(value, dev_r, dev_b);

	HANDLE_ERROR(cudaMemcpy(bool_f, dev_b, N * sizeof(bool), cudaMemcpyDeviceToHost));

	for(int i = 0; i < N; i++){
		if((bool_f[i]==true) && (bool_f[0]==true))
			bool_f[0] = true;
		else
			bool_f[0] = false;
	}
	if(bool_f[0]==true)
		printf("TAK \n");
	else
		printf("NIE \n");


	cudaFree(dev_b);
	cudaFree(dev_r);
	free(random);
	free(bool_f);

	cudaEventRecord(stop, 0);
	cudaEventSynchronize(stop);

	float elapsedTime;

	HANDLE_ERROR(cudaEventElapsedTime(&elapsedTime, start, stop));

	printf("Średni czas: %3.1f ms \n", elapsedTime);

	HANDLE_ERROR(cudaEventDestroy(start));
	HANDLE_ERROR(cudaEventDestroy(stop));

}

Przerobiłem również algorytm CPU, aby wykorzystywał dwa rdzenie CPU w sposób podobny jak program na GPU.

Moje wyniki [9223372036854775783, 4000]:

• CPU [1 thread] - 658 ms
• CPU [2 threads] - 284 ms
• GPU [4 x 1000 threads] - 200 ms

Wiem, że algorytm dla GPU można jeszcze zoptymalizować, ale samo przenoszenie danych pomiędzy CPU a GPU nie wpływa znacząco na wynik. Wiem też, że karta jednocześnie wykona 32 x 2 wątki, co powinno dać teoretycznie 64 krotny wzrost wydajności, a tymczasem jest nieco ponad 3 krotny. Oczywiście podczas uruchomienia programu sterownik nvidii pokazuje użycie karty w 99%, a wykorzystanie szyny danych 20%, pamięci: około 450 MB.

Platforma:

• CPU: Intel C2D E8400 @3000 MHz
• RAM: 3 GB DDR2 800 MHz
• GPU: NVIDIA GTX 650 1 GB DDR5
• SYSTEM: Ubuntu 14.04 LTS x64

Sugerujesz, że karta graficzna ma 64 rdzenie o takiej samej wydajności co rdzenie procesorów Intela ;)?