We have used MPI to divide arithmetic, move arrays, and exchange values with neighbors. Today we will revisit two problems from Chapter 3 of Pacheco and Malensek: the trapezoidal-rule program, which the book improves in stages, and odd-even sorting, which uses repeated neighbor exchanges. We will use shorter versions of the ideas so we can trace what each process does.
The complete programs are here if you want to type along or compare your work:
- mpi_trap_blocks.c — consecutive groups of trapezoids, with broadcast and reduction.
- mpi_odd_even_one_value.c — one number per rank, sorted through neighbor exchanges.
1. How the book improves its trapezoid program
Recall the goal: estimate the area under \(f(x)=x^2\) from 0 to 1. Divide the interval into \(n\) narrow trapezoids, calculate their areas, and add them. The exact area is \(1/3\).
The book does not present one mysterious “final” MPI program. Its mpi_trap1.c through mpi_trap4.c change how input and results move while keeping the numerical problem essentially the same:
| Book version | Main idea | Question it raises |
|---|---|---|
mpi_trap1.c (§3.2) |
Fixed inputs; each nonzero rank sends its partial area to rank 0. | Who owns which trapezoids, and how does rank 0 collect the answer? |
mpi_trap2.c (§3.3) |
Rank 0 reads the inputs and sends them to the other ranks. | How can we avoid entering the same input on every process? |
mpi_trap3.c (§3.4) |
Broadcast the inputs and reduce the partial areas. | Can we express “share with everyone” and “add everyone’s result” directly? |
mpi_trap4.c (§3.5) |
Package several input fields into a derived MPI datatype. | What if we want to send a more complicated record as one unit? |
We will stop at the broadcast-and-reduce idea. The derived datatype is useful later, but it obscures today’s data flow. The collectives make the code shorter and clearer; they do not guarantee that a tiny calculation runs faster.
Predict the work before running it
We will use 12 trapezoids and four ranks. Unlike our earlier cyclic loop, this version gives each rank a consecutive block of three trapezoids. If trapezoids are numbered 0–11:
| Rank | Trapezoid numbers | Part of \([0,1]\) |
|---|---|---|
| 0 | 0, 1, 2 | \([0, 1/4]\) |
| 1 | 3, 4, 5 | \([1/4, 1/2]\) |
| 2 | 6, 7, 8 | \([1/2, 3/4]\) |
| 3 | 9, 10, 11 | \([3/4, 1]\) |
Each rank computes the area for its block. A reduction adds those four partial areas on rank 0. The trapezoids meet at shared endpoints, but no trapezoid is counted twice.
Start the program with the usual MPI setup. On rank 0, set n = 12; broadcast it so every process knows the same n:
int n = 0;
if (rank == 0) n = 12;
MPI_Bcast(&n, 1, MPI_INT, 0, MPI_COMM_WORLD);Set a = 0, b = 1, and h = (b-a)/n on every rank. Since n divides evenly by size, each rank owns local_n = n/size trapezoids. The index of its first trapezoid is rank * local_n:
double local_area = 0.0;
for (int j = 0; j < local_n; j++) {
int i = rank * local_n + j;
double left = a + i * h;
double right = left + h;
local_area += h * (left * left + right * right) / 2.0;
}Pause and predict: What is the first i for rank 2? Why does the loop use j < local_n rather than j < n? Which variable holds only one rank’s answer?
Finish with a reduction and print on rank 0:
double total_area = 0.0;
MPI_Reduce(&local_area, &total_area, 1, MPI_DOUBLE,
MPI_SUM, 0, MPI_COMM_WORLD);
if (rank == 0) printf("Estimated area: %.9f\n", total_area);Compile the complete program and try it:
mpicc -Wall -o mpi_trap_blocks mpi_trap_blocks.c
mpirun -n 4 ./mpi_trap_blocksFor n = 12, expect approximately 0.334490741, a little above \(1/3 \approx 0.333333333\). Change the 12 to 1200, recompile, and run again. Does a finer numerical approximation change the number of ranks or the communication pattern? Would replacing the reduction with rank 0’s explicit receive loop change the mathematics?
Keep MPI_Finalize() in the program. PicoCluster has sometimes stalled during finalization on multi-board runs. To rehearse the example on one board, use /usr/bin/mpirun -hosts pc0 -genv HWLOC_COMPONENTS -gl -n 4 ./mpi_trap_blocks. That is four processes on one board, so do not use it as a four-board speed test.
2. Can MPI sort numbers?
The book’s mpi_odd_even.c sorts a block of numbers on each rank. Each rank sorts its own block, exchanges it with a neighbor, and keeps the appropriate half. The merge steps make that program long. We can see the same communication pattern more clearly with one number per rank.
Suppose ranks 0–3 initially hold:
rank: 0 1 2 3
value: 8 3 6 1In an even phase, pairs (0,1) and (2,3) exchange their values. The lower rank keeps the smaller value; the higher rank keeps the larger. In an odd phase, pair (1,2) does the same. Repeat even, odd, even, odd: after four phases, the values are sorted across the ranks.
| After phase | Pair(s) that compare | Values on ranks 0–3 |
|---|---|---|
| Start | — | 8 3 6 1 |
| 0 (even) | (0,1), (2,3) | 3 8 1 6 |
| 1 (odd) | (1,2) | 3 1 8 6 |
| 2 (even) | (0,1), (2,3) | 1 3 6 8 |
| 3 (odd) | (1,2) | 1 3 6 8 |
Trace it yourself before copying the last column. Why do ranks 0 and 3 sit out the odd phases? Why is one even phase insufficient?
The heart of mpi_odd_even_one_value.c is:
int partner;
if (phase % 2 == 0) {
partner = (rank % 2 == 0) ? rank + 1 : rank - 1;
} else {
partner = (rank % 2 == 0) ? rank - 1 : rank + 1;
}
if (partner >= 0 && partner < size) {
int neighbor_value;
MPI_Sendrecv(&value, 1, MPI_INT, partner, 0,
&neighbor_value, 1, MPI_INT, partner, 0,
MPI_COMM_WORLD, MPI_STATUS_IGNORE);
if (rank < partner && value > neighbor_value) value = neighbor_value;
if (rank > partner && value < neighbor_value) value = neighbor_value;
}MPI_Sendrecv lets each partner exchange the old value before either chooses which one to keep. It pairs the send and receive safely; it does not mean the two transfers happen at precisely the same instant. A rank with no neighbor skips the exchange for that phase. The program uses MPI_Scatter to hand out the four starting numbers and MPI_Gather to show the result on rank 0.
mpicc -Wall -o mpi_odd_even_one_value mpi_odd_even_one_value.c
mpirun -n 4 ./mpi_odd_even_one_valueExpected output:
Before: 8 3 6 1
After: 1 3 6 8Discuss: Which parts of this program depend on neighbors, and which involve the whole group? What would have to change if each rank held a block of values instead of just one? (The book’s answer is: sort each local block, then merge the two blocks during each neighbor exchange so the lower rank keeps the smaller half.)