En el post anterior quedamos en el fork-join, memoria compartida vs distribuida, y un primer “hola mundo” con OpenMP. Ahora toca la parte que de verdad vas a usar todos los días en este mundo de programación paralela, repartir un bucle entre threads, y no volarte el pie en el intento, porque créeme si pasa 😭
La directiva for 🔁
#pragma omp for [cláusulas]
for-loopEsto le dice a OpenMP que reparta las iteraciones del bucle siguiente entre los threads que ya existen. Ojo con ese detalle, for no crea threads nuevos, tiene que estar adentro de una región paralela para que haya algo entre lo cual repartir. No hay barrera implícita al entrar, pero sí una al salir.
#pragma omp parallel
{
#pragma omp for
for (i = 0; i < 1000; i++)
A[i] = B[i];
}Con dos threads, ese bucle de 1000 iteraciones queda repartido automáticamente en algo como i=0..500 para el Thread 1 e i=500..1000 para el Thread 2. Las variables A y B son compartidas, pero i es local a cada thread (recordá, el índice de un for asociado siempre es private), y toma un valor inicial y final distinto según el thread.
for contra calcularlo vos mismo 🧮
Estos dos códigos son equivalentes.
/* usando la directiva for */
#pragma omp parallel
#pragma omp for
for (i = 0; i < n; i++)
z[i] = a*x[i] + y;/* repartiendo a mano, solo con parallel */
#pragma omp parallel private(id, num, istart, iend)
{
id = omp_get_thread_num();
num = omp_get_num_threads();
istart = id * n / num;
iend = min(n, (id + 1) * n / num);
for (i = istart; i < iend; i++)
z[i] = a*x[i] + y;
}La diferencia es puro trabajo manual sin for, somos nosotros quienes vamos a calcular qué rango de i le toca a cada thread usando su id y el total de threads. Con for, OpenMP hace ese cálculo por vos. En la práctica siempre vas a preferir for, salvo que necesites un control muy fino del reparto que la cláusula schedule no te dé (más sobre eso más abajo).
Ejercicio, antes de seguir. ¿Qué diferencia hay entre estos dos programas?
/* programa A */
#pragma omp parallel private(i)
{
for (i = 0; i < 10; i++)
printf("Hello world %d\n", i);
}/* programa B */
#pragma omp parallel private(i)
{
#pragma omp for
for (i = 0; i < 10; i++)
printf("Hello world %d\n", i);
}Antes de bajar y mirar la respuesta, porfa piénsalo un poco, no leas lo que escribiré abajo, o bueno, también es válido a veces una ayuda, pero sí sería bastante bueno que lo trataras de resolver tú mismo.
La respuesta. En el programa A, cada thread ejecuta el bucle completo, las 10 iteraciones enteras, porque no hay ninguna directiva de reparto de trabajo adentro de la región paralela. Con 4 threads, vas a ver el mensaje “Hello world 0” hasta “Hello world 9” impreso 4 veces, una tanda por thread. En el programa B, el #pragma omp for reparte esas 10 iteraciones entre los threads disponibles, así que cada línea se imprime una sola vez en total, repartida entre threads. Es la diferencia entre “cada thread hace todo el trabajo” y “el trabajo se divide entre threads”, y es exactamente el tipo de error que se comete cuando te olvidás el for adentro de un parallel.
No todo bucle se puede paralelizar ⚠️
Para que for funcione, el bucle tiene que cumplir tres restricciones.
- Ser un bloque estructurado. Nada de
breaknigotopara salir antes de tiempo. - El número de iteraciones tiene que poder calcularse de antemano, antes de entrar al bucle.
- No puede haber dependencias entre iteraciones.
Las dos primeras son mecánicas, pero la tercera es la que de verdad importa entender. La regla para detectar dependencias es simple de enunciar. Si ejecutás el bucle en orden inverso y el resultado cambia, el bucle no es paralelo. Si da lo mismo en cualquier orden, sí lo es.
flowchart LR
subgraph OK["✅ A[i] = A[i] + B[i] → paralelo"]
direction LR
a0[A0] -.-> a0
a1[A1] -.-> a1
a2[A2] -.-> a2
end
classDef ok fill:#dcefdd,stroke:#1f5c2e,color:#1f5c2e,font-weight:bold;
class a0,a1,a2 ok
style OK fill:#f2faf3,stroke:#4a8f57flowchart LR
subgraph FWD["🚫 A[i] = A[i+1] + B[i] → no paralelo, dependencia hacia adelante"]
direction RL
b1[A1] --> b0[A0]
b2[A2] --> b1
b3[A3] --> b2
end
classDef bad fill:#fbd8d3,stroke:#8a241b,color:#8a241b,font-weight:bold;
class b0,b1,b2,b3 bad
style FWD fill:#fdf3f2,stroke:#c0392bflowchart LR
subgraph BACK["🚫 A[i] = A[i-1] + B[i] → no paralelo, dependencia hacia atrás"]
direction LR
c0[A0] --> c1[A1]
c1 --> c2[A2]
c2 --> c3[A3]
end
classDef bad fill:#fbd8d3,stroke:#8a241b,color:#8a241b,font-weight:bold;
class c0,c1,c2,c3 bad
style BACK fill:#fdf3f2,stroke:#c0392bCon A[i] = A[i] + B[i], cada iteración lee y escribe únicamente su propio elemento. No le importa nada de lo que hagan las demás, así que las podés ejecutar en cualquier orden. Con A[i] = A[i+1] + B[i], la iteración i necesita el valor original de A[i+1], que la iteración i+1 está a punto de sobrescribir. Si i+1 corre primero, i lee el valor equivocado. Y con A[i] = A[i-1] + B[i], cada iteración necesita el resultado que acaba de escribir la iteración anterior, así que directamente tienen que correr en orden, una detrás de la otra. Ninguno de estos dos últimos casos se puede repartir con for tal cual está.
Segundo ejercicio. Yo sé que me estoy poniendo como un profesor, con bastantes ejercicios por lección, pero de verdad tenemos que ejercitar un poco la mente con este tema para poder entenderlo ya que puede parecer un poco magia a veces, ahora sí, el ejercicio, analizá si estos tres fragmentos son paralelos, mirando qué elementos de A se leen y se escriben en cada iteración.
/* fragmento 1 */
for (i = 0; i < n; i++)
A[i] = A[i+n] + B[i];/* fragmento 2 */
for (i = 0; i < n; i++)
A[i] = A[i+m] + B[i];/* fragmento 3 */
for (i = 0; i < n; i++) {
C[i] = A[3*i+1] + 1;
A[2*i+7] = B[i] - 3;
}El fragmento 1 sí es paralelo. Como i va de 0 a n-1, el índice i+n siempre cae fuera del rango que el propio bucle escribe (A[i+n] nunca coincide con ningún A[i] que otra iteración esté modificando), así que no hay conflicto real. El fragmento 2 depende de m. Si m es mayor o igual que n, es el mismo caso que el fragmento 1 y es paralelo. Si m es más chico, A[i+m] puede solaparse con algún A[j] que otra iteración escribe, y ahí ya no podés garantizar nada sin conocer el valor exacto. El fragmento 3 es el más interesante. Escribe en A[2*i+7] y lee de A[3*i+1], dos índices con progresiones distintas. Para que haya dependencia real, tendría que existir algún par de iteraciones i, j donde 2*i+7 == 3*j+1, en general, para valores arbitrarios de n, esa coincidencia sí puede ocurrir, así que el fragmento no es seguro de paralelizar sin analizar el rango concreto de i. Este último es justo el tipo de caso donde “parece paralelo a simple vista” pero hay que hacer la cuenta antes de confiar.
Yendo más profundo, repartir trabajo sin pisarse 🕵️
Ya sabemos cómo identificar qué se puede paralelizar. Ahora veamos las herramientas para hacerlo bien cuando el reparto no es tan directo como sumar dos vectores.
Veamos el siguiente ejemplo práctico, tal vez lo tengas que leer varias veces y hacer por tu cuenta, yo no lo entendí a la primera, si te soy honesto.
reduction, sumar 100.000 números en una línea ➕
#include <stdio.h>
#include <omp.h>
#define N 100000
int main(void) {
double A[N], suma = 0.0;
for (int i = 0; i < N; i++) A[i] = 1.0;
#pragma omp parallel for reduction(+:suma)
for (int i = 0; i < N; i++)
suma += A[i];
printf("Suma total = %.0f\n", suma); /* 100000 */
return 0;
}Por debajo, reduction(+:suma) hace tres cosas. Le da a cada thread su propia copia privada de suma, inicializada en el elemento neutro de la operación (0 para +, 1 para *). Cada thread acumula en su copia solo los elementos que le tocaron, sin pisar la copia de nadie más (por eso no hace falta ni critical ni atomic adentro del bucle). Y al final de la región combina todas las copias en la variable original, con un esquema interno de reducción en árbol, en cada paso se suman parejas de valores y el número de valores vivos se reduce a la mitad, así que con P threads hacen falta del orden de log₂(P) pasos en vez de P-1 sumas una por una.
flowchart TD
S0[suma₀] --> R1
S1[suma₁] --> R1[+]
S2[suma₂] --> R2
S3[suma₃] --> R2[+]
R1 --> R3[+]
R2 --> R3
R3 --> T[total]
classDef leaf fill:#dcefdd,stroke:#1f5c2e,color:#1f5c2e,font-weight:bold;
classDef op fill:#4a8f57,stroke:#1f5c2e,color:#ffffff,font-weight:bold;
classDef final fill:#c0392b,stroke:#8a241b,color:#ffffff,font-weight:bold;
class S0,S1,S2,S3 leaf
class R1,R2,R3 op
class T finalreduction también admite otros operadores, * (producto), max/min, &&/||, &/|/^, cada uno con su elemento neutro correspondiente. Un matiz para tener en cuenta. La suma en coma flotante no es estrictamente asociativa por el redondeo, así que sumar los mismos números en distinto orden (como hace reduction al repartirlos entre threads, comparado con sumarlos uno a uno en secuencial) puede dar una diferencia mínima en los últimos decimales. No es un bug, es el precio de paralelizar algo que en matemática pura es asociativo pero en aritmética real no lo es del todo.
La suma de 100.000 números, paso a paso 🌳
reduction es cómodo, pero se siente casi mágico si no viste nunca qué hace por debajo. Así que vamos a desplegar el mismo problema a mano, elemento por elemento, hasta que quede clarísimo por qué el código de arriba tiene exactamente esa forma.
El problema. 10 threads suman esos mismos 100.000 números guardados en el vector A. Cada thread Pn (Pn = 0…9) se encarga de su porción de 10.000 elementos consecutivos y guarda el resultado parcial en sum[Pn]. Hasta acá no hay nada que coordinar, cada uno suma su trozo de forma completamente independiente.
sum[Pn] = 0;
for (i = 10000*Pn; i < 10000*(Pn+1); i++)
sum[Pn] = sum[Pn] + A[i]; /* suma de la porción local, sin coordinación */Lo interesante empieza ahora. Hay que combinar esos 10 resultados parciales en un único total. La forma ingenua sería que el thread 0 sumara sum[1]+sum[2]+…+sum[9] uno a uno, 9 sumas secuenciales, sin aprovechar que hay 10 threads libres para ayudar. El truco más inteligente es una reducción en árbol. En cada paso, la mitad de los valores todavía vivos se suma con su pareja, y el número de valores vivos se reduce aproximadamente a la mitad. Con P valores hacen falta del orden de log₂(P) pasos en vez de P-1 sumas una detrás de otra.
Con 10 threads la diferencia (3 pasos contra 9) no impresiona demasiado, pero se nota de verdad al escalar. Con P = 1.024 threads, la suma ingenua necesitaría 1.023 pasos, uno detrás del otro, mientras que el árbol lo hace en 10 pasos (log₂(1024) = 10), porque en cada paso todos los threads todavía vivos trabajan a la vez en vez de esperar turno. Esa diferencia, lineal contra logarítmica, es la razón de fondo por la que ninguna biblioteca paralela suma los resultados parciales uno a uno.
flowchart TB
subgraph I1["Inicio: 10 sumas parciales"]
direction LR
s0[S0] & s1[S1] & s2[S2] & s3[S3] & s4[S4] & s5[S5] & s6[S6] & s7[S7] & s8[S8] & s9[S9]
end
subgraph I2["Iteración 1, half: 10 → 5"]
direction LR
t0[S0+S5] & t1[S1+S6] & t2[S2+S7] & t3[S3+S8] & t4[S4+S9]
end
subgraph I3["Iteración 2, half impar: S4 se absorbe, luego 5 → 2"]
direction LR
u0["S0+S4+S2"] & u1["S1+S3"]
end
subgraph I4["Iteración 3, half: 2 → 1"]
v0["total en sum[0]"]
end
I1 --> I2 --> I3 --> I4
classDef alive fill:#dcefdd,stroke:#1f5c2e,color:#1f5c2e,font-weight:bold;
classDef dead fill:#f2f2f2,stroke:#9a9a9a,color:#9a9a9a;
classDef total fill:#c0392b,stroke:#8a241b,color:#ffffff,font-weight:bold;
class s0,s1,s2,s3,s4,t0,t1,t2,t3,t4,u0,u1 alive
class s5,s6,s7,s8,s9 dead
class v0 totalEl código completo, con la variable half llevando la cuenta de cuántos valores siguen vivos en cada momento.
half = 10; /* número de threads */
repeat
synch(); /* barrera: espera a que termine la suma parcial */
if (half % 2 != 0 && Pn == 0)
sum[0] = sum[0] + sum[half-1];
half = half / 2;
if (Pn < half) sum[Pn] = sum[Pn] + sum[Pn+half];
until (half == 1); /* el resultado final queda en sum[0] */synch() es una barrera. Nadie pasa a la siguiente iteración hasta que todos los threads llegaron a ese punto. Y no está de adorno. Sin ella, nada impide que un thread rápido entre a la iteración 2 mientras otro más lento todavía no terminó de escribir su resultado de la iteración 1. Un entrelazamiento concreto que rompe el resultado, sin synch(), en la primera iteración (half pasando de 10 a 5, sum0 += sum5).
| Instante | P0 (calcula sum0 += sum5) | P5 (todavía terminando su suma local) |
|---|---|---|
| t1 | lee sum5 → valor viejo, sin actualizar |
sigue en el bucle sumando su porción de A |
| t2 | escribe sum0 = sum0 + (valor viejo de sum5) |
termina de sumar y recién ahí escribe el sum5 correcto |
El resultado. sum0 queda con un total que no incluye todos los elementos que P5 debía sumar. Un error silencioso, sin ningún mensaje ni caída del programa, y encima no reproducible siempre igual (a veces P5 sí llega a tiempo y el resultado sale bien, dependiendo de qué tan cargada esté la máquina). Eso es justo lo que la barrera evita, obliga a los 10 threads a llegar todos al mismo punto antes de que ninguno lea el resultado de otro.
La traza completa, iteración a iteración.
| Iteración | half antes → después | Operación tras la barrera |
|---|---|---|
| 1 | 10 → 5 | sum0+=sum5; sum1+=sum6; sum2+=sum7; sum3+=sum8; sum4+=sum9 |
| 2 | 5 → 2 | half es impar, primero sum0+=sum4 (corrección), luego sum0+=sum2; sum1+=sum3 |
| 3 | 2 → 1 | sum0+=sum1 → fin, total en sum[0] |
Y ahí está, eso, exactamente eso, es lo que la cláusula reduction(+:suma) te ahorra escribir a mano.
La condición de carrera más típica 🏁
Por defecto, una variable declarada afuera de la región paralela es shared. Escribir sobre ella desde varios threads sin protección es una condición de carrera. El resultado depende de en qué orden se entrelacen las ejecuciones, y puede cambiar de una corrida a otra.
int maximo = a[0];
#pragma omp parallel for
for (int i = 1; i < N; i++) {
if (a[i] > maximo) /* ¡incorrecto! lectura+escritura sin proteger */
maximo = a[i]; /* sobre una variable compartida */
}Un entrelazamiento concreto que da un resultado equivocado, con maximo en 5 y dos threads procesando a[i]=8 y a[j]=6 casi al mismo tiempo.
| Instante | Thread A (procesa a[i]=8) | Thread B (procesa a[j]=6) |
|---|---|---|
| t1 | lee maximo → 5 |
lee maximo → 5 |
| t2 | compara 8 > 5 → cierto | compara 6 > 5 → cierto |
| t3 | escribe maximo = 6 |
|
| t4 | escribe maximo = 8 |
Con este orden el resultado sale bien de casualidad (maximo = 8), pero si las escrituras de t3 y t4 se intercambian, el Thread A escribe 8 primero y el Thread B lo pisa con 6 después. El valor 8 se pierde, aunque ambos threads “vieron” en algún momento que 8 era mayor. Esto se llama lost update. El problema no es que un thread lea basura, es que la secuencia leer-comparar-escribir de uno se mete en el medio de la de otro sin que ninguno se entere del progreso del otro.
El arreglo correcto no es declarar maximo como private (perdería el valor al salir del bucle), sino usar reduction también acá, con el operador max.
maximo = a[0];
#pragma omp parallel for reduction(max:maximo)
for (int i = 1; i < N; i++)
if (a[i] > maximo) maximo = a[i]; /* ahora sí, correcto */private(x) sigue siendo necesaria para variables auxiliares que cada iteración recalcula desde cero (un índice temporal, un acumulador local que no se lee de una iteración a otra) y que precisamente por eso no deberían compartirse entre threads.
critical y atomic, proteger una actualización compartida 🔒
atomic protege una sola operación de memoria (lectura, modificación y escritura de una variable) y suele traducirse a una instrucción de hardware atómica. Es la opción más rápida cuando aplica. critical protege un bloque de código arbitrario, con más overhead pero más flexibilidad.
long contador = 0;
#pragma omp parallel for
for (int i = 0; i < 1000000; i++) {
#pragma omp atomic
contador++; /* una sola operación → atomic alcanza */
}
printf("contador = %ld (esperado: 1000000)\n", contador);#pragma omp parallel for
for (int i = 0; i < n; i++) {
double local = f(i); /* cómputo caro, afuera de la sección crítica */
#pragma omp critical
{ /* if + dos escrituras: no entra en un atomic */
if (local > mejor_valor) {
mejor_valor = local;
mejor_indice = i;
}
}
}atomic solo admite formas muy puntuales de actualizar una única variable escalar, x++, x--, x op= expr, o x = x op expr. Por eso el segundo ejemplo no puede usar atomic. Hay un if y dos variables distintas a la vez, necesita critical. Y una advertencia de rendimiento que vale la pena recordar. Meter un critical adentro de un bucle muy iterado serializa justo esa parte. Si la mayor parte del tiempo del bucle se va en la sección crítica, el paralelismo no te sirve de mucho, por rápido que sea el resto del cuerpo.
Con esto ya tenemos lo esencial de OpenMP, cómo repartir un bucle, cómo saber si se puede repartir, y cómo protegerte de pisarte con vos mismo cuando varios threads tocan el mismo dato.
La canción del post
Ghost Town de Kanye West. La verdad no tiene nada de especial, pero la escuché muchísimo últimamente y es una canción bastante buena, de hecho la escuché varias veces mientras redactaba este post jajaja.