Soma cumulativa paralela (prefixo) no OpenMP: comunicando valores entre threads
Suponha que eu tenha uma funçãof(i)
que depende de um índicei
(entre outros valores que não podem ser pré-computados). Eu quero preencher uma matriza
de modo aa[n] = sum(f(i)) from i=0 to n-1
.
Editar: Após um comentário de Hristo Iliev, percebi que o que estou fazendo ésoma cumulativa / prefixo.
Isso pode ser escrito em código como
float sum = 0;
for(int i=0; i<N; i++) {
sum += f(i);
a[i] = sum;
}
Agora eu quero usar o OpenMP para fazer isso em paralelo. Uma maneira que eu poderia fazer isso com o OpenMP é escrever os valores paraf(i)
em paralelo e, em seguida, cuidar da dependência em série. E sef(i)
é uma função lenta, então isso pode funcionar bem, já que o loop não paralelo é simples.
#pragma omp parallel for
for(int i=0; i<N; i++) {
a[i] = f(i);
}
for(int i=1; i<N; i++) {
a[i] += a[i-1];
}
Mas é possível fazer isso sem o loop não paralelo com o OpenMP. A solução, no entanto, que eu inventei é complicada e talvez hackeada. Então, minha pergunta é se existe uma maneira mais simples e menos complicada de fazer isso com o OpenMP?
O código abaixo basicamente executa o primeiro código que listei para cada thread. O resultado é que os valores dea
em um determinado segmento estão corretos até uma constante. Eu salvo a soma de cada thread para um arraysuma
comnthreads+1
elementos. Isso me permite comunicar entre threads e determinar o deslocamento constante para cada thread. Então eu corrijo os valores dea[i]
com o deslocamento.
float *suma;
#pragma omp parallel
{
const int ithread = omp_get_thread_num();
const int nthreads = omp_get_num_threads();
const int start = ithread*N/nthreads;
const int finish = (ithread+1)*N/nthreads;
#pragma omp single
{
suma = new float[nthreads+1];
suma[0] = 0;
}
float sum = 0;
for (int i=start; i<finish; i++) {
sum += f(i);
a[i] = sum;
}
suma[ithread+1] = sum;
#pragma omp barrier
float offset = 0;
for(int i=0; i<(ithread+1); i++) {
offset += suma[i];
}
for(int i=start; i<finish; i++) {
a[i] += offset;
}
}
delete[] suma;
Um teste simples é apenas para definirf(i) = i
. Então a solução éa[i] = i*(i+1)/2
(e no infinito é-1/12).