En el artículo anterior se mostraba el cálculo de números primos con la librería GMP y en el anterior a éste, el cálculo con SQL y Perl.
Ahora vamos un paso más allá, distribuyendo el cálculo en procesos paralelos. Para ello, hay múltiples soluciones, como por ejemplo:
- Threading Building Blocks de Intel.
- Un Tesla de NVidia, que se programa con una variante de C llamada CUDA.
- Threads de Java sobre HotSpots.
- El lenguaje X10, basado en Java.
- etc.
Sin embargo, la opción más sencilla probablemente sea OpenMP, un estándar abierto que ya se encuentra integrado en el compilador GCC.
El siguiente código, que se compila con "
gcc -fopenmp -o hola hola.c", es un ejemplo muy sencillo del uso de OpenMP:hola.c
#include <omp.h>
#include <stdio.h>
// @author José Cotrino (http://www.cotrino.com)
int main() {
#pragma omp parallel
printf("- Hola caracola desde el hilo %d, hay %d hilos\n", omp_get_thread_num(), omp_get_num_threads());
const int N = 100000;
const int MAX_THREADS = 20;
int i, a[N], p[MAX_THREADS];
for (i=0;i<MAX_THREADS;i++) p[i]= 0;
printf("\n- Calculando un bucle 'for' de %d iteraciones en paralelo: \n + cada iteración puede ser asignada teóricamente a un procesador\n + el array sobre el que se calcula está en memoria compartida\n",N);
#pragma omp parallel for
for (i=0;i<N;i++) {
a[i]= 2*i;
p[omp_get_thread_num()] ++;
}
printf("\n- Yatá!\n");
for (i=0;i<MAX_THREADS;i++) {
if (p[i] != 0)
printf(" + Hilo %d se ha ejecutado %d veces\n", i, p[i]);
}
}
En un Intel Core 2 Duo obtenemos la siguiente salida:
- Hola caracola desde el hilo 1, hay 2 hilos
- Hola caracola desde el hilo 0, hay 2 hilos
- Calculando un bucle 'for' de 100000 iteraciones en paralelo:
+ cada iteración puede ser asignada teóricamente a un procesador
+ el array sobre el que se calcula está en memoria compartida
- Yatá!
+ Hilo 0 se ha ejecutado 50000 veces
+ Hilo 1 se ha ejecutado 50000 veces
Volviendo al código que teníamos para calcular primos con GMP, lo que haremos es simplemente dividir todo el rango en varios trozos, de forma que cada procesador calcule los números primos dentro de un subrango. Por ejemplo, si tenemos 2 núcleos y queremos calcular todos los primos entre 0 y 1000, en un procesador calculamos entre 0 y 500, y en el otro entre 501 y 1000. Esto es posible porque el test de primalidad que utiliza GMP no necesita conocer los primos anteriores.
El código queda así:
paraprimes.cpp
#include <stdio.h>
#include <time.h>
#include <gmp.h>
#include <omp.h>
// @author José Cotrino (http://www.cotrino.com)
void calculatePrimes ( mpz_t N_max ) {
time_t start,end;
double elapsed;
mpz_t N_primes, N_range;
mpz_init(N_primes);
mpz_init(N_range);
mpz_set_ui(N_primes, 0);
mpz_set(N_range, N_max);
time(&start);
// calculate N/threads
int N_cores = omp_get_num_procs();
mpz_cdiv_q_ui(N_range, N_range, N_cores);
// each thread calculates the primes between N/t*i and N/t*(i+1), where:
// - N is the total range where to search for primes
// - t is the amount of available processors
// - i is the number of the current thread/processor (0, 1, 2, 3...)
#pragma omp parallel shared(N_range, N_primes)
{
int thread = omp_get_thread_num();
mpz_t prime, N_max_thread, N_primes_thread;
mpz_init(prime);
mpz_init(N_max_thread);
mpz_init(N_primes_thread);
mpz_set(prime, N_range);
mpz_mul_ui(prime, prime, thread);
mpz_add_ui(prime, prime, 1);
mpz_set(N_max_thread, N_range);
mpz_mul_ui(N_max_thread, N_max_thread, thread+1);
mpz_set_ui(N_primes_thread, 0);
//printf("Thread %d searching in range %s, %s\n", thread, mpz_get_str(NULL,10,prime), mpz_get_str(NULL,10,N_max_thread));
while( mpz_cmp(prime, N_max_thread) < 0 ) {
mpz_nextprime(prime, prime);
//printf("%s, ", mpz_get_str(NULL,10,prime));
mpz_add_ui(N_primes_thread, N_primes_thread, 1);
}
// skip the last prime, which is over the range
mpz_sub_ui(N_primes_thread, N_primes_thread, 1);
// add the primes found by this thread to the whole amount
mpz_add(N_primes, N_primes, N_primes_thread);
}
time(&end);
elapsed = difftime(end, start);
printf("%s primes found below %s in %.0f seconds\n",
mpz_get_str(NULL,10,N_primes),
mpz_get_str(NULL,10,N_max),
elapsed);
//printf("%s %f\n",mpz_get_str(NULL,10,N_max),elapsed);
}
int main() {
int i;
mpz_t N_max;
mpz_init(N_max);
mpz_set_ui(N_max, 100);
printf("%d cores available\n", omp_get_num_procs());
for(i=0; i<8; i++) {
calculatePrimes(N_max);
mpz_mul_ui(N_max, N_max, 10);
}
}
Y se compila con:
g++ -lgmp -fopenmp paraprimes.cpp -o paraprimesLos resultados en un Intel Core 2 Duo T5750 son impresionantes:
| Máximo número | Primos | Duración |
|---|---|---|
| 100 | 25 | 0s |
| 1.000 | 168 | 0s |
| 10.000 | 1.229 | 0s |
| 100.000 | 9.592 | 0s |
| 1.000.000 | 78.498 | 0s |
| 10.000.000 | 664.579 | 7s |
| 100.000.000 | 5.761.455 | 66s (~1.1min) |
| 1.000.000.000 | 50.847.534 | 703s (~11.7min) |
Comparamos de nuevo con las implementaciones anteriores para calcular los primos por debajo de 100 millones:
- El algoritmo en C utilizando GMP sin paralelizar necesitaba 126 segundos. Con dos procesadores necesitamos 66 segundos, es decir, es 1.91 veces más rápido.
- El algoritmo en Perl necesitaba 5.929 segundos, que comparados con los 66 segundos actuales nos da que este algoritmo en paralelo es 89.83 veces más rápido.
Teóricamente, con este algoritmo y si la ley de Amdahl lo permite, con un quad-core se podrían calcular todos los primos por debajo de 100 millones en apenas 35 segundos, con un octo-core en 17 segundos, y con un Tesla en ná y menos. Me mantengo a la espera de que alguien me regale uno para poder hacer la prueba.
Entradas
