[C++] Geen speedup na parellisatie met OpenMP

[Verwijderd]

Ik ben bezig aan een algoritme dat ik op eight-core machines wil draaien om het allemaal wat sneller te laten gaan. Ik ben daarom mijn code aan het paralleliseren met OpenMP (ik gebruik GCC 4.2.2). Tot nu toe is m'n geparalleliseerde code echter eerder trager dan sneller. Een voorbeeld:

C++:

    int i, j, k;
    double val, totalSum, *rowSums;
    rowSums = (double*) calloc(n, sizeof(double));
    #pragma omp barrier

    // Compute the Gaussian kernel and corresponding row sums
    #pragma omp parallel for schedule(static) private(j,k,val) reduction(+:totalSum)
    for(i = 0; i < n; i++) {
        for(j = i + 1; j < n; j++) {
        
            // Compute the value
            val = 0.0;
            for(k = 0; k < d; k++) {
                val += pow(*(data + i * d + k) - *(data + j * d + k), 2.0);
            }
            val = exp(-val / (2 * pow(sigma, 2.0)));                
            
            // Store the value and update the row sums
            *(P + i * n + j) = val;
            *(P + j * n + i) = val;
            *(rowSums + i) += val;
            *(rowSums + j) += val;
            totalSum += val;
        }
    }
    #pragma omp barrier

Als ik de #pragma's uitvink is de code ongeveer 20% sneller op een quad-core machine (i.p.v. langzamer) voor n=5000. Ik heb dynamic scheduling geprobeerd met verschillende chunk sizes, maar dat levert weinig op in dit geval. Treden er lock-ups op in het sommeren van val of het schrijven naar P of rowSums ofzo? Hoe kan ik daar iets aan doen?


De double-arrays P en data, en double sigma komen trouwens van buiten deze code; maar het idee is duidelijk denk ik...

[.oisyn]

n=5000 is ook wel erg veel. Wat is d in dit geval? De code lijkt namelijk niet heel erg veel te doen, en dan krijg je al snel dat threading overhead de overhand krijgt. Een CPU kun je het beste lekker aan het werk zetten, ipv meer bezig laten zijn met het schedulen van korte threads.

[Verwijderd]

Dat n groot is, is toch juist mooi! Iedere core krijgt een stuk of 1250 van de iteraties voor z'n kiezen van de static scheduler (die maar een overhead van een paar milliseconden heeft?).

Normaal gesproken is d heel klein, 2 of 3 ofzo. De meeste rekentijd zit denk ik in het 0.5 * 5000² uitvoeren van de exp() functie, dat duurt ongeveer een seconde ofzo. Grotere stukken code paralleliseren is voor mij geen oplossing...

[leuk_he]

Dit is nieuw voor mij, dus ignore me als ik domme dingen zeg.

je maakt een static schedule aan. Maar niet alle chunks lijken me even groot. de for J loop i bij i=1 5000 elementen groot, en bij i =4999 1 element groot. Dynaimc schedule? Ik neem aan dat alleen de for (i=...) loop parallel gemaakt wordt.

[Verwijderd]

Klopt als een bus; ik heb dus zeker idle CPU ticks. De dynamic scheduler heeft echter een stuk meer overhead (ik heb verschillende chunk sizes geprobeerd), dus levert geen snelheidswinst op.
Dat ik wat idle ticks heb verklaart m.i. echter nog steeds niet waarom de code trager is op 4 cores dan op 1.

[.oisyn]

Verwijderd schreef op vrijdag 25 januari 2008 @ 17:58:
Dat n groot is, is toch juist mooi! Iedere core krijgt een stuk of 1250 van de iteraties voor z'n kiezen van de static scheduler (die maar een overhead van een paar milliseconden heeft?).

Ah ja natuurlijk, ik zat even te denken in de richting van 5000 threads ofzo. Maar ze worden natuurlijk gewoon evenredig verdeeld over de cores.

[leuk_he]

Verwijderd schreef op vrijdag 25 januari 2008 @ 17:11:
Als ik de #pragma's uitvink is de code ongeveer 20% sneller op een quad-core machine (i.p.v. langzamer) voor n=5000. Ik heb dynamic scheduling geprobeerd met verschillende chunk sizes, maar dat levert weinig op in dit geval. Treden er lock-ups op in het sommeren van val of het schrijven naar P of rowSums ofzo? Hoe kan ik daar iets aan doen?

Dan haal je P en Rowsum er toch even uit om te benchmarken of de locking die daarop optreedt in de weg zit? (ook al is het resultaat dan niet meer zinnig)

Duidelijk lijkt me wel dat als je over multiple dies(neem aan dat acht cores 2 quadcores is) praat alle access van P en rowsum via main memory gaat, en dus de access daarop relatief duur is omdat die hele tijd de cache geflushed moet worden tussen beide fysieke CPU's

Ik zie allen niet zo snel een oplossing daarvoor.

[Verwijderd]

Als ik de P en rowSums regels eruit haal en in de loop die val sommeert de data even vervang door een vast getalletje (zodat ik zeker nergens lockups heb), zijn de single- en multicore versie precies even snel. Het zal dus wel een combinatie van lockups en cache zijn die de problemen veroorzaakt.

De P-matrix hoeft echter niet gelockt te worden, omdat ik nooit op dezelfde plek iets opsla. Kan ik niet op één of andere manier aangeven dat P niet gelockt moet worden? En is het niet mogelijk om op de rowSums ook een soort + reductie te doen, net als bij de totalSum?

In de standaardvoorbeeldjes van OpenMP doen ze net of dit soort dingen vanzelf goed zouden gaan

C++:

#define N 100000
int main(int argc, char *argv[])
{
  int i, a[N];
  #pragma omp parallel for
  for (i=0;i<N;i++) 
     a[i]= 2*i;
  return 0;
}

[Nick The Heazk]

Je performantie kan natuurlijk afnemen door de interactie overhead die je krijgt. Zo update je bijvoorbeeld steeds

C++:

*(P + i * n + j) = val;
*(P + j * n + i) = val;
*(rowSums + i) += val;
*(rowSums + j) += val; 

bij elke kolom. Mijn inziens update je hier gedeelde variabelen, waardoor je een immense interactie overhead krijgt. Een immense overhead krijg je allesinds op een disjunct memory space architectuur. Hetgeen je waarschijnlijk niet hebt?. Zo ga je trouwens ook last krijgen van false sharing; je cache van de andere processoren/cores zal steeds invalidated worden. (Daar opeenvolgende elementen van deze rij ook steeds in de cache worden geladen).

die barrier aan het einde van je lus is overigens niet nodig. Alle threads zullen standaard synchroniseren aan het einde van de omp lus.

Edit: Het feit dat je efficiëntie afneemt als je meer processoren inschakelt, met een constante probleem-grootte is overigens niet zo verwonderlijk. Dat fenomeen staat ook wel bekend als Amdahl's Law. (Al heb jij in dit geval wel te maken met een enorme afname in efficiëntie. Dat is waarschijnlijk te wijten aan de false sharing die optreedt.)

[ Voor 63% gewijzigd door Nick The Heazk op 25-01-2008 19:37 ]

[Soultaker]

Ik denk ook dat het misgaat door het sommeren. Sowieso zou ik de code altijd opsplitsen in verschillende loops zodat je binnen de loop maar één ding doet; dat maakt het een stuk simpeler en zorgt dat je bijvoorbeeld meer relevante data in je cache kunt houden.

Dan krijg je dus iets als:

C++:

    #pragma omp parallel for private(j,k,val)
    for(i = 0; i < n; i++) {
        for(j = i + 1; j < n; j++) {
            val = 0.0;
            for(k = 0; k < d; k++) {
                val += pow(*(data + i * d + k) - *(data + j * d + k), 2.0);
            }
            val = exp(-val / (2 * pow(sigma, 2.0)));
            *(P + i * n + j) = *(P + j * n + i) = val;
        }
    }

    totalSum = 0;
    #pragma omp parallel for
    for(i = 0; i < n; i++) {
        /* Als de som eerst in een lokale variable berekend, lijkt me dat je
           synchronisatie in de binnenste loop sowieso vermijdt. */
        double sum = 0;
        for(j = 0; j < n; ++j)
           sum += *(P + i*n + j);
        /* Hier vindt waarschijnlijk wel synchronisatie plaats */
        rowSums[i] = sum;
        totalSum += sum;
    }

Verder is het altijd zinnig bij dit soort code om te kijken naar de assembly uitvoer die gegenereert wordt. Die OpenMP directives worden ook maar gewoon naar functional calls gemapt.

[leuk_he]

Eerst een klein stukje disassemby: (vs2008 )

code:

double var;
double P[100];
int i;
...
        P[i]+=var;
0040103C  fadd        qword ptr [esp+ecx*8] 
0040103F  mov         dword ptr [esp+4],ecx 
00401043  fstp        qword ptr [esp+ecx*8]

Ofwel, het updaten van die double is geen atomaire actie, dus hier wordt voor parrellelisatie locking toegevoegd.

Soultaker schreef op vrijdag 25 januari 2008 @ 19:40:
Ik denk ook dat het misgaat door het sommeren. Sowieso zou ik de code altijd opsplitsen in verschillende loops zodat je binnen de loop maar één ding doet; dat maakt het een stuk simpeler en zorgt dat je bijvoorbeeld meer relevante data in je cache kunt houden.

Dan krijg je dus iets als:

C++:

    #pragma omp parallel for private(j,k,val)
            *(P + i * n + j) = *(P + j * n + i) = val;
        }

Ik zou
*(P + i * n + j) houden, en de 2e P in de 2e loop toewijzen. Je doet 2 toewijzingen(2x cache flush, wellicht 2x een lock , dus 1 dure interactie kun je naar de 2e loop verschuiven.
....

#pragma omp parallel for
totalSum = 0;

Wellicht kun je beter geen parallelisatie doen in de 2e loop.

[ Voor 0% gewijzigd door leuk_he op 25-01-2008 23:01 . Reden: parallelisatie parallelisere geparalleliseerde ... wat een rot woord ]

[Soultaker]

leuk_he schreef op vrijdag 25 januari 2008 @ 23:00:
Wellicht kun je beter geen parallelisatie doen in de 2e loop.

Mits n groot genoeg is, waarom niet?

[H!GHGuY]

kun je i,j en k ook niet lokaler declareren?

je declareert i, j en k waarna je een barrier plaatst.
kan OpenMP dan niet besluiten dat deze global zijn en dus geshared moeten worden met locking?

(ik ken eigenlijk geen OpenMP)

[Confusion]

Heb je al gegoogled op parallel algorithm gaussian kernel? Ik vermoed dat er in dit topic een wiel opnieuw uitgevonden wordt.

[Verwijderd]

Dat gaat vooral over convoluties met Gaussian kernels, dat is toch wat anders. Ik vind één artikel over parallelisatie van Gaussian kernels met MPI, waar ze ook geen speed-up haalden met meer nodes (maar goed, dat was geen shared memory architecture).

Die integers j en k kan ik niet binnen de #pragma pas declareren; als je dat doet snapt hij niet waar die private(j,k) in de #pragma clause over gaat.

Ik heb nu alle operaties zoveel mogelijk uit elkaar gehaald, en ook een reduction toegevoegd voor bijvoorbeeld die val. Als ik dan de schedule van de duurste loop (die met de pow() en de exp()) op dynamic zet met een redelijke chunksize (heb nu 100), heb ik op een dualcore machine iig geen slowdown meer. Als er nou nog een manier was om te zeggen dat hij P niet hoeft te locken in onderstaande code...

C++:

// Compute the Gaussian kernel
    #pragma omp parallel for schedule(dynamic, 50) private(j,k) reduction(+:val)
    for(int i = 0; i < n; i++) {
        for(j = i + 1; j < n; j++) {
        
            // Compute the value
            val = 0.0;
            for(k = 0; k < d; k++) {
                val += pow(data[i * d + k] - data[j * d + k], 2.0);
            }
            val = exp(-val / (2 * pow(sigma, 2.0)));                
            
            // Store the value
            P[i * n + j] = val;
        }
    }

[H!GHGuY]

vervangen van P[i*n+j] door P[i][j]?

Is voor de compiler wsch een stuk duidelijker.

[Reptile209]

Ik heb de ballen verstand van parallel programming, dus handle this post with care .
Kan je niet op een of andere mantier de matrix voor de output ook lokaal maken voor iedere processor? Na je parallele rekenwerk voeg je de matriches dan samen in één grote copy-operatie die vast minder duur is dan al het rekenwerk.

[Nick The Heazk]

Verwijderd schreef op zondag 27 januari 2008 @ 01:10:
. Als er nou nog een manier was om te zeggen dat hij P niet hoeft te locken in onderstaande code...

Dat kun je bereiken door je berekende waarde niet meteen in de globale rij op te slaan. Declareer een lokale rij voor elk van de processors. Deze kan je dan helemaal aan het einde gaan kopiëren in de globale structuur.

Dat zorgt voor een veel hogere lokaliteit van data-toegangen. Bovendien verminder je ook het aantal keren dat je een lock aanvraagt; het is niet onwaarschijnlijk dat de compiler voor het kopiëren van de rij met resultaten in de globale structuur slechts een keer een lock hoeft aan te vragen. Dit in tegenstelling natuurlijk met die binnenste lus, alwaar dat niet mogelijk is; daar zal je steeds een lock aanvragen om die waarde op te slaan.

Zoals hier ook al boven mij wordt aangegeven

[ Voor 31% gewijzigd door Nick The Heazk op 27-01-2008 17:43 ]

[Verwijderd]

Misschien een domme vraag, maar hoe declareer ik één lokale array voor iedere thread? Als ik in iedere i-iteratie een nieuwe malloc() moet doen wordt het natuurlijk niet bepaald sneller:

C++:

    // Compute the Gaussian kernel
    #pragma omp parallel for schedule(dynamic, 50) private(j,k) reduction(+:val)
    for(int i = 0; i < n; i++) {
    
        // Initialize local row
        double* rowP = (double*) malloc(n * sizeof(double));
    
        for(j = i + 1; j < n; j++) {
        
            // Compute the value
            val = 0.0;
            for(k = 0; k < d; k++) {
                val += pow(data[i * d + k] - data[j * d + k], 2.0);
            }
            val = exp(-val / (2 * pow(sigma, 2.0)));                
            
            // Store the value in row
            rowP[j] = val;
        }
        
        // Copy to global row
        memcpy(P + i * n, rowP, n * sizeof(double));
        free(rowP);
    }

[Verwijderd]

Die integers j en k kan ik niet binnen de #pragma pas declareren; als je dat doet snapt hij niet waar die private(j,k) in de #pragma clause over gaat.

Als een variable gedefineerd word in de scope van een omp directive, dan is hij automatisch private, dus je kan dat hele private weglaten.

Verwijderd schreef op zondag 27 januari 2008 @ 21:33:
Misschien een domme vraag, maar hoe declareer ik één lokale array voor iedere thread? Als ik in iedere i-iteratie een nieuwe malloc() moet doen wordt het natuurlijk niet bepaald sneller:

C++:

#pragma omp parallel
{
  double* rowP = (double*) malloc((n / omp_get_num_threads() + 1) * sizeof(double));

  // geen "parallel" in dit directive, de forking is al gedaan door het 
  // omp parallel directive hier boven. Dit for directive geeft aan 
  // dat elk van de bestaande threads, slechts een deel van de loop uit moet voeren.
  #pragma omp for
  for(int i = 0; i < n; i++)
  {
  }

  free(rowP);
}

[Verwijderd]

Ah op die fiets! Ik had het zelf kunnen bedenken
Het draait nu als een zonnetje op m'n dualcore (70% speedup), dus ik ben benieuwd hoe het morgen op de eight-core loopt. Voor de geïnterresseerden de uiteindelijke code:

C++:

    // Compute the Gaussian kernel
    #pragma omp parallel
    {               
        // Initialize memory to store row of P
        double* rowP = (double*) malloc(n * sizeof(double));
        double val;
    
        #pragma omp for schedule(dynamic, 50)
        for(int i = 0; i < n; i++) {        
            for(int j = i + 1; j < n; j++) {
            
                // Compute the value
                val = 0.0;
                for(int k = 0; k < d; k++) {
                    val += pow(data[i * d + k] - data[j * d + k], 2.0);
                }
                val = exp(-val / (2 * pow(sigma, 2.0)));                
                
                // Store the value in row
                rowP[j] = val;
            }
            
            // Copy to global row
            memcpy(P + i * n + i + 1, rowP + i + 1, (n - i - 1) * sizeof(double));
        }
        
        // Free the memory
        free(rowP);
    }

edit: Als ik de data-array voor de parallele for-loop kopieer in een private array wordt het zelfs nog iets sneller, omdat je dan ook van de lockups in het ophalen van de data af bent...

[ Voor 9% gewijzigd door Verwijderd op 28-01-2008 05:01 ]

[pkuppens]

Niet relevante opmerking, maar wat gebeurt er als je deze constante berekening
uit de dubbele for loop haalt, en door een sigma_squared_doubled vervangt?

C++:

sigma_squared_doubled = (2 * pow(sigma, 2.0))
for (blahi)
  for(blahj)
    dingetje = f(iets/sigma_squared_doubled)

[Verwijderd]

Je zou zeggen dat het dan sneller zou moeten worden, maar in de praktijk levert het niks op. Misschien dat de compiler dat oplost? En de exp() is natuurlijk een veel duurdere bewerking...

[mineral]

Probeer de deling er ook eens uit te halen. Dus f(iets*inverse_sigma_squared_doubled). Delingen zijn vele malen duurder dan vermenigvuldigingen.

Verder is mijn ervaring met pow(x, 2.0) dat dit niet altijd goed geoptimaliseerd wordt. Het wil wel eens helpen om de vermenigvuldiging expliciet uit te schrijven.

[Verwijderd]

Goed idee ja

De pow() is inderdaad langzamer dan de vermenigvuldiging gewoon uitschrijven, ik denk door de extra function-call. Maar ik bereken de macht nu nog maar één keer, dus dat zal niet veel uitmaken...