Autor Wątek: CUDA i działanie kernela  (Przeczytany 1422 razy)

Offline bartekm

  • Użytkownik
    • Wordpress Blog

# Kwiecień 21, 2010, 15:11:12
Witam, mam taki problem, mam sobie taki kawałek kodu, jest to kernel CUDA, czyli funkcja wywoływana na GPU. Ten kod powinien obliczać przyspieszenie działające pomiędzy ciałem "i-tym" i "j-tym" ciałem. Jednak nie oblicza mi poprawnie tzn, zwraca a = 0. I moje pytanie co robię nie tak.


__global__ void ObliczPrzysp(float3 *pol, float3 *a, float *mass, int M, float GG, float pp)
{
int i = blockIdx.x * blockDim.x + threadIdx.x;

int j = blockIdx.y * blockDim.y + threadIdx.y;

float4 delta;
float Q;

if(i<M && j<M)
{
if(i!=j)
{
delta.x = pol[i].x - pol[j].x;
delta.y = pol[i].y - pol[j].y;
delta.z = pol[i].z - pol[j].z;
delta.w = sqrtf(delta.x*delta.x + delta.y*delta.y + delta.z*delta.z + pp);
Q = GG*mass[j]/(delta.w*delta.w*delta.w);
a[i].x += Q*delta.x;
a[i].y += Q*delta.y;
a[i].z += Q*delta.z;
}
}
}


Ten sam kod napisany pod jedną iterację, tzn:

__global__ void ObliczPrzysp(float3 *pol, float3 *a, float *mass, int M, float GG, float pp)
{
int i = blockIdx.x * blockDim.x + threadIdx.x;

float4 delta;
float Q;

for(j=0;j<M;j++)
{
if(i!=j)
{
delta.x = pol[i].x - pol[j].x;
delta.y = pol[i].y - pol[j].y;
delta.z = pol[i].z - pol[j].z;
delta.w = sqrtf(delta.x*delta.x + delta.y*delta.y + delta.z*delta.z + pp);
Q = GG*mass[j]/(delta.w*delta.w*delta.w);
a[i].x += Q*delta.x;
a[i].y += Q*delta.y;
a[i].z += Q*delta.z;
}
}
}


działa poprawnie, tzn. oblicza mi "a" i wszystko przebiega dobrze, ale wydajność nie jest taka jaką bym chciał, ze względu na umiejscowioną pętlę "for" wewnątrz kernela. Dla dużych M, wątek musi wykonać naprawdę dużo obliczeń, co znacznie spowalnia aplikację i zużywa grafikę tak, że nie da się normalnie pracować. Czy da się to jakoś zooptymalizować? Jak powinien wyglądać kod, żeby to działało pod "dwuwymiarowe bloki" tzn gdy:
int i = blockIdx.x * blockDim.x + threadIdx.x;

int j = blockIdx.y * blockDim.y + threadIdx.y;

Offline Mr. Spam

  • Miłośnik przetworów mięsnych

Offline _OskaR

  • Użytkownik

# Kwiecień 21, 2010, 16:18:29
Specjalistą w tej kwestii nie jestem, więc napiszę tylko to, co przychodzi mi do głowy:
1. blockIdx.y itd. - próbowałeś w ogóle wyświetlić iteratory? Tworzysz w ogóle właściwą liczbę wątków?
2. Możesz wygenerować tablicę ze wszystkimi kombinacjami iteratorów. Pomijając powtórzenia - dla 100 elementów masz 10k kombinacji. Rozbijasz na 10k wątków itd. Pojawia się jednak problem zależności (w wielu wątkach zapisujesz do tej samej zmiennej) - nie wiem, jak CUDA sobie z tym radzi. Zawsze można zrobić 10k wyników, a potem osobną funkcję, 100 wątków i jedziesz z a.x += Q*delta.
3. Robisz owe 10k wątków i "odkodowujesz" i,j - blockIdx.x * blockDim.x + threadIdx.x; traktujesz jako 100*i+j - odwracasz działanie, korzystając m.in. z dzielenia modulo i masz.
4. Może pierwszy wariant działa niepoprawnie, tj.  masz np. 100 razy mniej obliczeń i dlatego uważasz, że drugi wariant muli.

Offline lux

  • Użytkownik
    • homepage

# Kwiecień 21, 2010, 23:51:19
Problemem jest ten kawałek kodu:

a[i].x += Q*delta.x;
a[i].y += Q*delta.y;
a[i].z += Q*delta.z;

Zapisywać do a(i) może wiele wątków jednocześnie i nie masz żadnej gwarancji, że otrzymasz tam poprawną sumę. Najprostszym rozwiązaniem jest zapisywanie sobie sum częsciowych w shared memory (po jednej komórce dla każdego wątków w bloku) i na koniec posumowanie za pomocą operacji atomowych jednym wątkiem tego co jest w shared memory.

Możesz też spróbować zamienić i, j rolami. Najlepiej by było gdyby w obrębie bloku każdy wątek j zapisywał (dodawał) do innego i (wtedy nie ma tego problemu z dodawaniem do a).

Offline bartekm

  • Użytkownik
    • Wordpress Blog

# Kwiecień 22, 2010, 13:28:46
Najlepiej by było gdyby w obrębie bloku każdy wątek j zapisywał (dodawał) do innego i (wtedy nie ma tego problemu z dodawaniem do a).

Właśnie to chciałbym uzyskać, tylko nie wiem jak. Jeśli:
int i = blockIdx.x * blockDim.x + threadIdx.x;
// Jakis kod

jest to samo co:
for(i=0;i<N;i++)
{
   // Jakis kod
}
to myślałem, że dając:
int i = blockIdx.x * blockDim.x + threadIdx.x;

int j = blockIdx.y * blockDim.y + threadIdx.y;

// Jakis kod

będzie to samo co:
for(i=0;i<N;i++)
{
   for(j=0;j<N;j++)
   {
       // Jakis kod
   }
}
tak przynajmniej wynika z przykładów w dokumentacji, ale widzę, że nie ma tak lekko. Będę próbował coś zdziałać w tym kierunku, o którym pisałeś, choć nie znam na razie zasady działania shared memory.

@_OskaR
Liczba wątków jest chyba ok, jest określona przez wyrażenie
dim3 dimBlock(Liczba.watkow.x, Liczba.watkow.y);
dim3 dimGrid(N/Liczba.watkow.x, N/Liczba.watkow.y);

A potem kernel jest odpalany przez:

Oblicz <<< dimGrid, dimBlock >>> ( // Jakies zmienne );

Z pierwszym wariantem jest na pewno coś nie tak, o ile potrafię sobie poradzić z jednowymiarowym kernelem to dwuwymiarowe wogóle mi nie wychodzą. Zrobienie tablicy z wszystkimi możliwościami raczej nie wchodzi w rachubę ze względu na duży rozmiar takiej tablicy i ograniczonej pamięci karty graficznej. Dla 100 elementów może by przeszło, ale dla np. 100 000 już nie.

Offline _OskaR

  • Użytkownik

# Kwiecień 22, 2010, 13:51:27
Z pierwszym wariantem jest na pewno coś nie tak, o ile potrafię sobie poradzić z jednowymiarowym kernelem to dwuwymiarowe wogóle mi nie wychodzą. Zrobienie tablicy z wszystkimi możliwościami raczej nie wchodzi w rachubę ze względu na duży rozmiar takiej tablicy i ograniczonej pamięci karty graficznej. Dla 100 elementów może by przeszło, ale dla np. 100 000 już nie.
Dlatego napisałem jeszcze punkt 3. - bez tablicy, bez zużywania dużej ilości pamięci.
Myślę, że Twój problem da się rozwiązać, tylko trzeba pokombinować. Wstępna wersja algorytmu - robisz punkt 3., czyli 10k wątków dla 100 elementów (100*100) - lub odpowiednio więcej/mniej. W każdym wątku zapisujesz wynik do tablicy (10k elementów, więc nie będzie zależności). To koniec pierwszej funkcji. Masz więc wyniki częściowe. Teraz sumowanie (druga funkcja) - z tego, co zrozumiałem, dla 10k elementów częściowych, wyjdzie 100 wyników końcowych - no to dajesz tablicę 100 elementów, 100 wątków i w każdym pętla, gdzie sumujesz elementy z tablicy 10k - pierwszy wątek - suma elem. od 0 do 99, drugi od 100 do 199 itd.
Wspominałeś o pamięci - można zmniejszyć jej zużycie, realizując algorytm w kilku etapach - zamiast 10k wątków i tablicy 10k, dajesz np. 1k wątków i tablicę 1k. Wykonujesz 2 funkcje i powtarzasz to 10 razy (pamiętając, że trzeba wtedy trochę inaczej odkodować i,j - w każdej kolejnej iteracji "przesuwasz" przedział). Dodatkowo, w drugiej funkcji (robisz 100 wątków - bo tyle będzie wyników końcowych - jak wcześniej), w każdej iteracji, w każdym wątku, po prostu obliczasz sumę 10 elementów - czyli 100 wątków * 10 iteracji * 10 elementów - 10k wyników częściowych (czyli ilość się zgadza).

Offline bartekm

  • Użytkownik
    • Wordpress Blog

# Kwiecień 25, 2010, 15:46:46
Problem rozwiązany. Błędem było zadeklarowanie bloku wątków tzn.
dim3 dimBlock(Liczba.watkow.x, Liczba.watkow.y);
dim3 dimGrid(N/Liczba.watkow.x, N/Liczba.watkow.y);

Do poprawnego działania Liczba.watkow.x*Liczba.watkow.y musi być mniejsza lub równa od maksymalnej ilości wątków per blok (u mnie 512), po zmianie na:
dim3 dimBlock(16, 16);

Wszystko działa ok, bo 16*16<512. Szkoda tylko, że nie wpłynęło to na wydajność, a nawet jest gorsza niż przy jednowymiarowym bloku:(