Parallel Computing [BCS702]

Module-3
Distributed memory programming with MPI – MPI functions, The trapezoidal rule in MPI,
Dealing with I/O, Collective communication, MPI-derived datatypes, Performance evaluation of
MPI programs, A parallel sorting algorithm.
Introduction
In parallel MIMD systems, there are two main types: distributed-memory and shared-memory
systems. In a distributed-memory model, each core has its own private memory and
communicates via a network (see Fig. 3.1). In contrast, shared-memory systems allow all cores
to access a common global memory (see Fig. 3.2). This chapter begins exploring how to program
distributed-memory systems using message-passing.

Fig 3.1: A distributed memory system

Fig 3.2: A shared memory system

In message-passing programs, each core-memory pair typically runs a separate process, and
communication occurs through send and receive function calls. We'll be using MPI (Message-
Passing Interface), a standardized library—not a new language—that enables such
communication in C and Fortran programs. MPI provides various send and receive functions,
which we’ll explore further. We'll explore collective communication functions in MPI, which
involve multiple processes. Alongside learning MPI functions, we'll cover key topics like
data partitioning, distributed-memory I/O, and parallel program performance.

3.1 Getting Started

Many of us were first introduced to programming through a version of the "hello, world" program
from Kernighan and Ritchie's classic book.

Dept. of CSE- GMIT PROF. SMITHA M S 1

 Parallel Computing [BCS702]

We’ll create a simple MPI-based “hello, world” program where only process 0 prints the output.
The other processes, identified by ranks from 0 to p−1, send messages to process 0, which prints
them. This setup reflects a common pattern in parallel programming.

Program 3.1: MPI program that prints greetings from the processes

3.1.1 Compilation and execution

Compiling and running the program may vary by system, so it's best to consult a local expert.
Typically, we assume using a text editor for writing code and the mpicc command in the terminal
to compile MPI programs.

Dept. of CSE- GMIT PROF. SMITHA M S 2

 Parallel Computing [BCS702]

The mpicc command is usually a wrapper script for the C compiler. It simplifies compilation by
automatically including the necessary header files and linking the required libraries.
Many systems also support program startup with mpiexec:

So to run the program with one process, we’d type

and to run the program with four processes, we’d type

With one process, the program’s output would be

and with four processes, the program’s output would be

The <mpiexec> command starts <number of processes> instances of the <mpi_hello> program
and may assign each to a specific core. Once started, the MPI system manages communication
between the processes.
3.1.2 MPI programs
This is a C program that includes standard header files like <stdio.h> and <string.h>, along with
a main function. Additionally, Line 3 includes <mpi.h>, which provides all necessary MPI
definitions, prototypes, and macros required to compile an MPI program.
All MPI identifiers begin with MPI_, making them easy to distinguish. MPI function names and
types use a capital letter after the underscore (e.g., MPI_Comm_rank), while constants and
macros are written in all capital letters (e.g., MPI_COMM_WORLD).
3.1.3 MPI_Init and MPI_Finalize
In Line 12, The MPI_Init function initializes the MPI environment, performing setup tasks like
allocating message buffers and assigning ranks to processes. It must be the first MPI function
called in the program, and its syntax includes pointers to argc and argv.

The arguments argc_p and argv_p in MPI_Init are typically the pointers to main's arguments, but
if not used, NULL can be passed. MPI_Init returns an integer error code, though it’s often
ignored to keep the code clean. The function MPI_Finalize, usually called at the end, tells MPI
that the program is done and allows MPI to release its resources.

Dept. of CSE- GMIT PROF. SMITHA M S 3

 Parallel Computing [BCS702]

In Line 31, The MPI_Finalize call signals that the program is done using MPI and allows the
system to free any allocated resources.

No MPI functions should be called after MPI_Finalize, making it the final call in a typical MPI
program structure.

It's not required to pass pointers to argc and argv to MPI_Init, nor must MPI_Init and
MPI_Finalize be called from main.
3.1.4 Communicators, MPI_Comm_size, and MPI_Comm_rank
In MPI, a communicator is a group of processes that can exchange messages, and MPI_Init
defines one such communicator called MPI_COMM_WORLD, which includes all user-started
processes. Function calls in Line 13 and 14 often use MPI_COMM_WORLD to obtain
information about the participating processes.

Both MPI_Comm_size and MPI_Comm_rank take a communicator of type MPI_Comm as their

first argument. MPI_Comm_size gives the total number of processes (commonly stored in
comm_sz), and MPI_Comm_rank gives the rank of the calling process (commonly stored in
my_rank) within that communicator.
3.1.5 SPMD programs
In parallel programming with MPI, it's common to write a single program for all processes, even
though each process may perform different tasks. For example, in our case, process 0 receives
and prints messages, while other processes create and send them. This is done using conditional
branching based on the process rank. This style is known as Single Program, Multiple Data
(SPMD), and in our program, the if-else structure in Line 16-29 helps implement this model.
It's important to note that our program is designed to run with any number of processes. While
we've tested it with 1 or 4 processes, it could also work with 1000 or even 100,000, depending

Dept. of CSE- GMIT PROF. SMITHA M S 4

 Parallel Computing [BCS702]

on system resources. Although MPI doesn't enforce this flexibility, it's a common practice to
write scalable programs, since the exact number of available cores can change. For instance, we
might use a 20-core system today and access a 500-core system tomorrow.
3.1.6 Communication
In lines 17–18, all processes except process 0 create a message to send using sprintf, which
formats the message into a string. Lines 19–20 handle the sending of these messages to process
0. Meanwhile, process 0 prints its own message using printf and then enters a loop to receive and
print messages from processes 1 through comm_sz − 1. Lines 25–26 are responsible for receiving
messages from each process one by one.
3.1.7 MPI_Send
The sends performed by processes 1 through comm_sz − 1 involve some complexity, so it's
helpful to examine them in detail. Each send operation is executed using the MPI_Send function,
which follows the syntax shown above.

The first three arguments of MPI_Send—msg_buf_p, msg_size, and msg_type—specify what

data is being sent. The other arguments—dest, tag, and communicator—indicate where the
message should go.
The first argument msg_buf_p points to the memory block holding the message content, which
in our program is the string greeting. The second and third arguments, msg_size and msg_type,
specify the size and type of data to send—in this case, strlen(greeting) + 1 characters of type
MPI_CHAR. Together, these inform MPI of the message’s length and data type.
MPI uses a special type called MPI_Datatype for the msg_type argument because C types like
int and char can’t be directly passed to functions. MPI provides predefined constants for common
types, as shown in Table 3.1.
The message size (msg_size) doesn't have to match the full size of the buffer (like greeting), but
it must not exceed it. For example, if the message is 31 characters and the buffer is 100 characters,
that's acceptable. The key is ensuring the message fits within the allocated storage.
Fourth, dest argument in MPI_Send specifies the rank of the receiving process, while the tag
helps distinguish between similar messages. For example, process 1 can send different types of
data to process 0 using different tags—such as tag 0 for data to print and tag 1 for data to
compute—so the receiver knows how to handle each message.
Table 3.1: Some predefined MPI datatypes.

Dept. of CSE- GMIT PROF. SMITHA M S 5

 Parallel Computing [BCS702]

The Fifth argument in MPI_Send is a communicator, which defines the group of processes that
can communicate with each other. Messages sent within one communicator can’t be received by
processes in a different one. This helps prevent accidental message reception and supports safer
communication in complex MPI programs.
3.1.8 MPI_Recv
The first six arguments to MPI_Recv correspond to the first six arguments of MPI_Send:

The first three arguments of MPI_Recv define where to store the received message: msg_buf_p
points to the memory, buf_size is the number of items it can hold, and buf_type specifies the data
type. The next three—source, tag, and communicator—identify the message’s origin, label, and
communication group. The last argument, status_p, often isn’t used and can be replaced with
MPI_STATUS_IGNORE if the status info isn’t needed.
3.1.9 Message matching
Suppose process q calls MPI_Send with

Also suppose that process r calls MPI_Recv with

Dept. of CSE- GMIT PROF. SMITHA M S 6

 Parallel Computing [BCS702]

Then the message sent by q with the above call to MPI_Send can be received by r with the call
to MPI_Recv if

For a message to be received successfully, the arguments send_buf_p/recv_buf_p,

send_buf_sz/recv_buf_sz, and send_type/recv_type must define compatible buffers. This means
the receive buffer should be large enough and of the correct type to match the sent data. Most of
the time, following this compatibility rule is sufficient.
If recv_type = send_type and recv_buf_sz ≥ send_buf_sz, then the message sent by q can be
successfully received by r.
Sometimes a process may receive messages from multiple processes without knowing the order
in which the messages will arrive.
To avoid this problem MPI provides a special constant MPI_ANY_SOURCE that can be passed
to MPI_Recv. Then, if process 0 executes the following code, it can receive the results in the
order in which the processes finish:

If a process is receiving several messages with different tags and doesn’t know the order they’ll
arrive, it can use the special constant MPI_ANY_TAG in the MPI_Recv function. This allows
the process to accept a message with any tag, regardless of the sending order.
Here are two key points about wildcard arguments in MPI:
1. Only the receiver can use wildcard values like MPI_ANY_SOURCE or
MPI_ANY_TAG; senders must always specify the exact destination and tag. This means
MPI communication is "push"-based—senders initiate the data transfer.
2. Wildcards cannot be used for communicator arguments—both sender and receiver must
use the same communicator.
3.1.10 The status_p argument
If you think about these rules for a minute, you’ll notice that a receiver can receive a message
without knowing
1. the amount of data in the message,
2. the sender of the message, or

Dept. of CSE- GMIT PROF. SMITHA M S 7

 Parallel Computing [BCS702]

3. the tag of the message.

To find out details like the sender’s rank or the tag of a received message, the receiver can use
the last argument of MPI_Recv, which is a pointer to an MPI_Status struct. This struct includes
fields such as MPI_SOURCE, MPI_TAG, and MPI_ERROR, which store information about the
received message.

Then after a call to MPI_Recv, in which &status is passed as the last argument, we can determine
the sender and tags by examining the two members:

The exact amount of data received isn't directly available from the MPI_Status struct. Instead,
you can use the MPI_Get_count function to retrieve it. For example, if your receive buffer type
is recv_type and you passed &status to MPI_Recv, you can call MPI_Get_count to find out how
many items were received.

will return the number of elements received in the count argument. In general, the syntax of
MPI_Get_count is

The count isn't directly stored in the MPI_Status variable because it depends on the type of data
received and would require a calculation (like dividing bytes received by bytes per item). If this
information isn't needed, it's more efficient to skip calculating it.
3.1.11 Semantics of MPI_Send and MPI_Recv
When sending a message, the process creates the message with both data and metadata
(destination rank, sender rank, tag, communicator, and message size).
 There are two ways MPI can handle the message: buffering or blocking.
 In buffering, the message is stored internally by MPI and MPI_Send returns
immediately.
 In blocking, the function waits until the message transmission can begin, possibly
delaying the return of MPI_Send.
 When MPI_Send returns, it doesn't guarantee the message has been sent—only that the
send buffer can be reused.
 If confirmation of message delivery or immediate return is needed, MPI provides other
specialized send functions.
Here is a simplified, pointwise version of the passage:

Dept. of CSE- GMIT PROF. SMITHA M S 8

 Parallel Computing [BCS702]

 The behavior of MPI_Send depends on the MPI implementation.

 Most implementations use a cutoff message size:

o If the message is smaller than the cutoff, it is buffered.

o If the message is larger, MPI_Send will block until it can proceed.

 In contrast, MPI_Recv always blocks until a matching message is received.

 MPI also provides non-blocking receive functions that can check for a message and
return immediately, whether or not a message is available.

 MPI ensures message order from the same sender:

o If process q sends two messages to r, the first message must arrive before the
second.
 No order is guaranteed for messages from different senders:
o If processes q and t both send to r, their arrival order at r can vary.
 This is because MPI can’t control network delays or speeds, so it avoids enforcing
cross-process message order.

3.1.12 Some potential pitfalls

 MPI_Recv will block if there’s no matching MPI_Send, causing the process to hang.
 Ensure every receive has a corresponding send, with correct tags and ranks.
 Mistakes like mismatched tags or incorrect ranks can cause deadlocks or incorrect
behavior.
 If MPI_Send blocks and there’s no matching receive, the sender may also hang.
 If MPI_Send is buffered and no receive exists, the message may be lost.

3.2 The trapezoidal rule in MPI

Printing messages from processes is fine, but the real goal of learning MPI is to do more
meaningful parallel work than just printing.

3.2.1 The trapezoidal rule

To approximate the area under a curve y=f(x), the trapezoidal rule divides the interval into n
equal parts. (See Fig. 3.3.) Each part forms a trapezoid with base h=x i+1−xi and vertical sides
f(xi) and f(xi+1). (See Fig. 3.4.) The area of each trapezoid is then computed using these values.

Dept. of CSE- GMIT PROF. SMITHA M S 9

 Parallel Computing [BCS702]

Fig 3.3: The trapezoidal rule: (a) area to be estimated and (b) approximate area using
trapezoids.

Fig 3.4: One trapezoid

Since the interval [a,b] is divided into n equal parts, the width of each subinterval is given by

Thus if we call the leftmost endpoint x0, and the rightmost endpoint xn, we have that

Dept. of CSE- GMIT PROF. SMITHA M S 10

 Parallel Computing [BCS702]

and the sum of the areas of the trapezoids—our approximation to the total area—is

Thus, pseudocode for a serial program might look something like this:

3.2.2 Parallelizing the trapezoidal rule

Parallel programmers often use the term “parallelize” to describe converting a serial program
into a parallel one. To design a parallel program, we typically follow four steps:

1. Partition the problem into tasks,

2. Identify communication between tasks,
3. Group tasks into larger units,
4. Assign these units to cores.

For the trapezoidal rule, this means identifying tasks like computing the area of each trapezoid
and summing them. These tasks then communicate as shown in Fig. 3.5.

Fig 3.5: Tasks and communications for the trapezoidal rule.

Dept. of CSE- GMIT PROF. SMITHA M S 11

 Parallel Computing [BCS702]

To map tasks to cores, we note that using more trapezoids gives a more accurate result, so we
typically have many more trapezoids than cores. To manage this, we divide the interval [a,b] into
comm_sz subintervals, where comm_sz is the number of processes. If comm_sz divides n (the
number of trapezoids) evenly, each process computes the area using n/comm_sz trapezoids.
Finally, process 0 collects and adds all the estimates.

Let’s make the simplifying assumption that comm_sz evenly divides n. Then pseudocode for the
program might look something like the following:

For now, we skip user input and directly assign values to a, b, and n, as shown in Program 3.2.
The Trap function, used in Program 3.3, implements the serial trapezoidal rule.

We use different variable names to distinguish between local and global values. For instance,
local_a, local_b, and local_n are local to each process, while a, b, and n are shared values relevant
to all processes. Note that this concept of "local" and "global" differs from standard programming
terminology, but the context clarifies the meaning.

Dept. of CSE- GMIT PROF. SMITHA M S 12

 Parallel Computing [BCS702]

Program 3.2: First version of MPI trapezoidal rule.

Program 3.3: Trap function in MPI trapezoidal rule.

Dept. of CSE- GMIT PROF. SMITHA M S 13

 Parallel Computing [BCS702]

3.3 Dealing with I/O

The current parallel trapezoidal rule program is limited—it only works for the interval
0,30, 3 with 1024 trapezoids. Editing and recompiling the code each time is inefficient
compared to just entering new values. So, we need to add user input handling, and while
doing that, it’s also worth reviewing how output works in parallel programs.

3.3.1 Output

In the “greetings” and trapezoidal rule programs, process 0 writes to stdout using printf, which
usually works as expected.

 Although MPI doesn’t mandate I/O access rules, most implementations allow all
processes in MPI_COMM_WORLD to access stdout and stderr.
 However, there’s no automatic control over the order of output from multiple processes.
 If several processes write to stdout at the same time, their outputs can appear in an
unpredictable or interleaved manner.
 For example, in a simple program where all processes print messages (Program 3.4), the
output might be jumbled when run with multiple processes.

Program 3.4: Each process just prints a message.

Dept. of CSE- GMIT PROF. SMITHA M S 14

 Parallel Computing [BCS702]

However, when we run it with six processes, the order of the output lines is unpredictable:

 MPI processes compete for access to stdout, so the order of their outputs is unpredictable.
 This results in nondeterminism—the output may change from run to run.
 To avoid this, we must structure our program to control the output order.
 One common approach is to let all processes send their messages to process 0, which
then prints the messages in rank order, as done in the “greetings” program.

3.3.2 Input

 Most MPI implementations allow only process 0 to access stdin.

 Allowing all processes access to stdin creates ambiguity in input distribution.
 Therefore, process 0 reads input (e.g., using scanf) and sends it to the other processes.
 This requires checking process rank and using message-passing—similar to the approach
in the “greetings” program.
 Program 3.5 illustrates this with a Get_input function where process 0 reads a, b, and n
and sends them to the rest.

Dept. of CSE- GMIT PROF. SMITHA M S 15

 Parallel Computing [BCS702]

Program 3.5: A function for reading user input.

To use this function, we can simply insert a call to it inside our main function, being careful to
put it after we’ve initialized my_rank and comm_sz:

3.4 Collective communication

Our trapezoidal rule program can be improved, especially in how it handles the global sum. Right
now, only process 0 performs the final addition, while all other processes just send their results
and stop. This approach is inefficient—like having seven workers hand their tools to one worker
and letting that one person do the rest. Instead, we could divide the summing work more evenly
across all processes to better utilize available resources.

Dept. of CSE- GMIT PROF. SMITHA M S 16

 Parallel Computing [BCS702]

3.4.1 Tree-structured communication

We can improve the global sum by using a binary tree structure (see Fig. 3.6). First, processes 1,
3, 5, and 7 send their values to 0, 2, 4, and 6, which add the received values. Then, in the next
steps, 2 and 6 send to 0 and 4, which again add the values. Finally, process 4 sends its value to
0, which adds it to complete the total.

Fig 3.6: A tree-structured global sum.

 In the original method, process 0 does all the work: 7 receives and 7 additions.
 In the binary tree scheme, process 0 does only 3 receives and 3 adds, while others do at
most 2 operations.
 Multiple processes can work in parallel—for example, processes 0, 2, 4, and 6 can
compute simultaneously in the first step.
 This reduces total computation time significantly, as only 3 steps are needed instead of
7.
 With more processes (e.g., 1024), process 0 would need only 10 steps instead of 1023,
improving efficiency by over 100×.

Designing a tree-structured global sum is helpful but can be complex and time-consuming to
code. There are multiple ways to pair processes, like (0,4), (1,5), etc., and each pairing changes
how the tree is built (see Fig. 3.7). Choosing the best structure isn’t straightforward—it may vary
depending on the problem size or the system. So, manually testing and comparing alternatives
across different systems and scales might be necessary.

Dept. of CSE- GMIT PROF. SMITHA M S 17

 Parallel Computing [BCS702]

Fig 3.7: An alternative tree-structured global sum.

3.4.2 MPI_Reduce

 Writing an optimal global-sum manually is impractical due to too many possibilities.

 MPI includes built-in global-sum operations to save developers time and effort.
 Optimization is handled by the MPI library implementers, who understand the system
well.
 Global-sum involves all processes in a communicator (e.g., MPI_COMM_WORLD).
 These are known as collective communications, unlike MPI_Send/MPI_Recv, which
are point-to-point communications.

Global-sum is just one example of many collective communication operations. Sometimes,

instead of a sum, we may need to find the maximum, minimum, product, or other operations
on values distributed across processes. To handle this, MPI provides a generalized function
that can perform various such operations using a single interface.

Dept. of CSE- GMIT PROF. SMITHA M S 18

 Parallel Computing [BCS702]

The generalization relies on the fifth argument, operator, which is of type MPI_Op, similar to
MPI_Datatype and MPI_Comm. MPI provides several predefined operators like MPI_SUM (see
Table 3.2), and users can also define custom ones. To perform a global sum, we simply use
MPI_SUM as the operator and replace Lines 18–28 of Program 3.2 with a single MPI function
call.

Table 3.2 Predefined Reduction Operators in MPI.

An important feature of MPI_Reduce is that it can work on arrays by setting the count argument
to a value greater than 1. This allows you to, for example, sum N-dimensional vectors—one from
each process—using a single MPI_Reduce call.

3.4.3 Collective vs. point-to-point communications

Collective communications differ from point-to-point communication in key ways:

1. All processes in a communicator must call the same collective function; mismatched calls
(e.g., MPI_Reduce with MPI_Recv) will cause errors or hangs.
2. Arguments must be compatible across processes; inconsistent values (like different
destinations) lead to incorrect behavior.

Dept. of CSE- GMIT PROF. SMITHA M S 19

 Parallel Computing [BCS702]

3. Even though only the destination process uses output_data_p, all processes must pass a
valid argument for it.
4. Collective calls are matched by communicator and call order, not tags—so mismatched
order can lead to unexpected results, as shown in Table 3.3.

Table 3.3 Multiple Calls to MPI_Reduce.

Avoid using the same buffer for both input and output in MPI_Reduce, as it is illegal due to
argument aliasing—where two arguments refer to the same memory block. This can lead to
unpredictable behavior, including incorrect results or crashes. MPI disallows this to maintain
consistency with Fortran, which also prohibits aliasing. In certain cases, MPI offers alternative
constructs to bypass this limitation safely.

3.4.4 MPI_Allreduce

 In some programs, like the trapezoidal rule, only one process needs the result of a global
sum, so MPI_Reduce works well.
 But if all processes need the result, manually distributing it can be complex and
inefficient (e.g., using a reversed tree or a butterfly pattern as shown in Figs. 3.8 and 3.9).
 Choosing and implementing an optimal distribution method is difficult and varies by
system.
 MPI simplifies this with MPI_Allreduce, which performs a reduction and stores the result
on all processes—no dest_process argument is needed.

Dept. of CSE- GMIT PROF. SMITHA M S 20

 Parallel Computing [BCS702]

Fig 3.8: A global sum followed by distribution of the result.

Fig 3.9: A butterfly-structured global sum.

3.4.5 Broadcast

Dept. of CSE- GMIT PROF. SMITHA M S 21

 Parallel Computing [BCS702]

 Just as a tree-structured global sum improved performance, we can use a similar approach
to distribute input data efficiently.
 Reversing the communication pattern from Fig. 3.6 gives a tree-structured broadcast (Fig.
3.10).
 This method allows one process to share its data with all others more efficiently than
sending individually.
 MPI supports this pattern using a collective communication called a broadcast, which
sends data from one process to all in the communicator.

Fig 3.10: A tree-structured broadcast.

The process with rank source_proc sends the contents of the memory referenced by data_p to
all the processes in the communicator comm.

Program 3.6 shows how to modify the Get_input function shown in Program 3.5 so that it uses
MPI_Bcast, instead of MPI_Send and MPI_Recv.

Dept. of CSE- GMIT PROF. SMITHA M S 22

 Parallel Computing [BCS702]

Program 3.6: A version of Get_input that uses MPI_Bcast.

In serial programs, an in/out argument is both used and modified by a function. In MPI_Bcast,
the data_p argument acts as input on the source process and as output on the other processes,
meaning its role can differ depending on the process rank.

3.4.6 Data distributions

Suppose we want to write a function that computes a vector sum:

To implement vector addition in MPI, we treat each addition of corresponding elements as an

independent task. Since there's no communication needed between these tasks, we just need to
divide and assign them to processes. If n is the number of vector elements and comm_sz is the
number of processes, we assume n is divisible by comm_sz and let local_n = n / comm_sz. Each
process then handles a block of local_n consecutive elements—a method known as block
partitioning (as shown in Table 3.4).

Dept. of CSE- GMIT PROF. SMITHA M S 23

 Parallel Computing [BCS702]

Program 3.7: A serial implementation of vector addition.

Table 3.4 Different partitions of a 12-component vector among 3 processes.

In a cyclic partition, vector elements are assigned to processes in a round-robin manner (e.g.,
with n = 12 and comm_sz = 3, process 0 gets components 0, 3, 6, 9, and so on). In a block-cyclic
partition, blocks of elements (e.g., 2 elements per block) are distributed cyclically among
processes. Once a partitioning method is chosen, parallel vector addition is straightforward: each
process adds its assigned elements. Each process handles local_n components, stored locally as
an array, and uses a function similar to the serial version to perform the addition (see Program
3.8).

Program 3.8: A parallel implementation of vector addition.

Dept. of CSE- GMIT PROF. SMITHA M S 24

 Parallel Computing [BCS702]

3.4.7 Scatter

To test the vector addition function, process 0 can read the vector size and broadcast it to all
processes. While broadcasting the entire vectors x and y is possible, it's inefficient—each process
would unnecessarily allocate memory for all components. Instead, it's better for process 0 to read
the full vectors and send only the relevant subvector to each process based on the chosen
distribution (e.g., block). This way, each process only stores and works on the components it
actually needs.

 MPI_Scatter sends parts of a vector from a source process (src_proc, typically 0) to all
processes in the communicator comm.
 The entire vector is stored in send_buf_p on src_proc, and it's divided into comm_sz equal
parts (assuming block distribution).
 Each process receives its chunk into recv_buf_p, with size recv_count, which should be
local_n = n / comm_sz.
 send_count should also be local_n, as it specifies how much data each process receives—
not the total size of send_buf_p.
 send_type and recv_type should be set to MPI_DOUBLE for a vector of doubles.
 With this, each process receives only its required portion of the vector, avoiding
unnecessary memory use.

This approach is used in Read_vector (Program 3.9) to efficiently distribute input data.

Program 3.9: A function for reading and distributing a vector.

Dept. of CSE- GMIT PROF. SMITHA M S 25

 Parallel Computing [BCS702]

MPI_Scatter sends consecutive blocks of send_count elements to each process in rank order,
starting with process 0. Therefore, it only works correctly with block distribution when the total
number of elements n is evenly divisible by the number of processes comm_sz.

3.4.8 Gather

 MPI_Gather collects data from all processes in the communicator comm and gathers it
to a single destination process dest_proc.
 Each process sends send_count elements of type send_type from send_buf_p.
 On the destination process, the received data is stored in recv_buf_p, with recv_count
elements expected from each sender, of type recv_type.
 Data from process 0 goes into the first block of recv_buf_p, from process 1 into the
second block, and so on.

 This setup works well for block-distributed vectors, where each process contributes a
fixed block (e.g., local_n components), and recv_count reflects the number of elements
from each process, not the total. We can use this to print a distributed vector by
implementing Program 3.10.
 The restrictions on MPI_Gather are like those on MPI_Scatter: it only works correctly if
the vector uses a block distribution with equal-sized blocks for all processes.

Dept. of CSE- GMIT PROF. SMITHA M S 26

 Parallel Computing [BCS702]

Program 3.10: A function for printing a distributed vector.

3.4.9 Allgather

We want to write an MPI function that multiplies a matrix A by a vector x. If A = (aᵢⱼ) is an m

× n matrix and x have n components, the result y = Ax will have m components. Each component
y[i] is computed as the dot product of the i-th row of A with x, i.e., This logic forms the basis
for the MPI-based matrix-vector multiplication. See Fig 3.11

Fig 3.11: Matrix-vector multiplication.

Dept. of CSE- GMIT PROF. SMITHA M S 27

 Parallel Computing [BCS702]

So we might write pseudocode for serial matrix multiplication as follows:

C allows actual 2D arrays, but due to certain limitations, programmers often simulate them
using 1D arrays. The most common method is to store rows one after another in a single array,
effectively flattening the 2D structure into 1D.

In a one-dimensional representation of a 2D array, the element at row i and column j is stored

at index i × n + j, where n is the number of columns. For example, the element in row 2,
column 1 (value 9) is at index 2 × 4 + 1 = 9.

Using this one-dimensional scheme, we get the C function shown in Program 3.11.

Program 3.11: Serial matrix-vector multiplication.

To parallelize matrix-vector multiplication, we treat each multiplication and addition as a task.

It’s efficient to assign row i of matrix A and the corresponding element y[i] to the same process.
The simplest way in MPI is to use a block row distribution, assuming comm_sz divides the total
number of rows m evenly.

Dept. of CSE- GMIT PROF. SMITHA M S 28

 Parallel Computing [BCS702]

Since each y[i] is computed using row i of matrix A, both should be assigned to the same
process—so we use block distribution for y as well. Although giving every process a full copy
of x would reduce communication, it's common in practice to distribute x like y, especially when
the result vector y becomes the input x in repeated computations.

To ensure that each process has access to all components of x before executing the loop:

we can use the MPI collective communication function MPI_Allgather. This function gathers
each process’s local block of x and distributes the entire vector to all processes. As a result,
each process will have a full copy of x, allowing the loop to compute its assigned portion of y
correctly.

Instead of using both MPI_Gather and MPI_Bcast, which typically involve two tree-structured
communications, MPI offers a more efficient single function. This function performs the same
task using a butterfly communication pattern for better performance.

This function gathers data from all processes by concatenating their send_buf_p contents into
each process’s recv_buf_p. Typically, recv_count equals send_count, representing the data size
from each process. The parallel matrix-vector multiplication (see Program 3.12) can be
optimized by allocating x once in the calling function and reusing it.

Dept. of CSE- GMIT PROF. SMITHA M S 29

 Parallel Computing [BCS702]

Program 3.12: An MPI matrix-vector multiplication function.

3.5 MPI-derived datatypes

In distributed-memory systems, communication is typically far more expensive than local

computation. Sending multiple small messages is usually much slower than sending the same
amount of data in a single message. For instance, multiple send/recv calls in a loop perform
worse than a single combined call.

Dept. of CSE- GMIT PROF. SMITHA M S 30

 Parallel Computing [BCS702]

On some systems, using loops for multiple sends and receives can be 50 to 100 times slower than
a single combined message. Reducing the number of messages greatly improves performance.
MPI helps with this by offering three methods: using the count argument to group data,
creating derived datatypes, and using MPI_Pack/Unpack. This section focuses on one
method of building derived datatypes.

 In MPI, a derived datatype can represent a collection of data items by capturing their
types and relative memory locations.
 This allows send/receive functions to gather and scatter items to the correct locations
automatically.
 For example, instead of three separate MPI_Bcast calls for a, b, and n in the trapezoidal
rule program, we can create one derived datatype for them.
 With this datatype (two doubles and one int), a single MPI_Bcast call can send all three
variables at once.
 This simplifies communication and reduces the number of function calls.

A derived datatype in MPI is a sequence of basic datatypes, each paired with a displacement
indicating its memory location relative to the start. For example, in the trapezoidal rule, if a, b,
and n are stored at different memory addresses on process 0, we can use their types and relative
positions to define a single derived datatype for broadcasting all three together.

Then the following derived datatype could represent these data items:

Each entry in a derived datatype includes the data type and its displacement from the start.
Assuming the structure begins with a, its displacement is 0. Then, b is 16 bytes after a (40 − 24)
beyond the start of a, and n is 24 bytes after a (48 − 24) beyond the start of a.

MPI provides the function MPI_Type_create_struct to create a derived datatype composed of

elements with different basic types and specific memory displacements. This allows grouping
variables of various types (like int, double, etc.) into a single communication unit.

Dept. of CSE- GMIT PROF. SMITHA M S 31

 Parallel Computing [BCS702]

 The count argument specifies how many elements are in the derived datatype (e.g., 3 for
a, b, and n).
 array_of_block_lengths holds the number of items in each block; it's useful if any element
is itself an array.
 For example, if the first element is an array of 5 items, the first entry in
array_of_block_lengths would be 5.

array_of_blocklengths [0] = 5;

However, in our case, none of the elements is an array, so we can simply define

int array_of_blocklengths [3] = {1, 1, 1};

The third argument to MPI_Type_create_struct, array_of_displacements specifies the

displacements in bytes, from the start of the message. So we want

array_of_displacements [] = {0, 16, 24};

To find these values, we can use the function MPI_Get_address:

To get the values for array_of_displacements, use MPI_Get_address, which returns the
memory address of a variable. The MPI_Aint type is used to store these addresses, as it can
hold any memory address value on the system.

Dept. of CSE- GMIT PROF. SMITHA M S 32

 Parallel Computing [BCS702]

The array_of_datatypes should store the MPI datatypes of the elements. So we can just define

MPI_Datatype array_of_types [3] = { MPI_DOUBLE , MPI_DOUBLE , MPI_INT } ;

With these initializations, we can build the new datatype with the call MPI_Datatype

Before we can use input_mpi_t in a communication function, we must first commit it with a
call to

This allows the MPI implementation to optimize its internal representation of the datatype for
use in communication functions. Now, to use input_mpi_t, we make the following call to
MPI_Bcast on each process:

MPI_Bcast (&a, 1, input_mpi_t , 0, comm ) ;

So we can use input_mpi_t, just as we would use one of the basic MPI datatypes.

When a new MPI derived datatype is created, the MPI system may allocate internal storage.
Once the datatype is no longer needed, you should release this storage using a call to
MPI_Type_free.

We followed the outlined steps to create a Build_mpi_type function, which is now used within
the updated Get_input function. Both functions are shown in Program 3.13.

Dept. of CSE- GMIT PROF. SMITHA M S 33

 Parallel Computing [BCS702]

Program 3.13: The Get_input function with a derived datatype.

3.6 Performance evaluation of MPI programs

Let's evaluate the performance of the matrix-vector multiplication program. Parallel programs
are typically written with the expectation of being faster than their serial counterparts. To confirm
this, we need a method of performance comparison—so we'll begin by revisiting key concepts
previously discussed.

3.6.1 Taking timings

When measuring performance, we're not concerned with total program runtime—only the time
spent on key computations like the matrix-vector multiplication. Input/output tasks like typing
or printing are excluded. To measure just the computation time, we can use MPI_Wtime, a
function that returns the number of seconds elapsed since an arbitrary past time, allowing us to
time specific code sections.

Dept. of CSE- GMIT PROF. SMITHA M S 34

 Parallel Computing [BCS702]

double MPI_Wtime ( void ) ;

We can time a block of MPI code as follows:

For timing serial code, there's no need to use MPI libraries. Instead, the POSIX function
gettimeofday can be used, which returns elapsed time in microseconds. A convenient C macro,
GET_TIME, defined in timer.h (available from the book’s website), simplifies this and should
be used with a double variable to store the time.

After executing this macro, now will store the number of seconds since some time in the past.
So we can get the elapsed time of serial code with microsecond resolution by executing

The GET_TIME macro inserts timing code directly into your source via the preprocessor and
operates on a double variable, not a pointer. Since timer.h isn't a system header, you must specify
its location during compilation. For example, if it's in /home/peter/my_include, use:

to compile your program.

Dept. of CSE- GMIT PROF. SMITHA M S 35
 Parallel Computing [BCS702]

Both MPI_Wtime and GET_TIME return wall clock time, which measures total elapsed time,
including idle time. In contrast, timers like the C clock function report only CPU time, excluding
time spent waiting, such as during MPI_Recv. Wall clock time is more useful for evaluating
parallel program performance.

A few timing issues remain in parallel programs. Each process reports its own time, but we
typically want a single time—the time taken by the slowest process. Although we can't guarantee
all processes start simultaneously, we can approximate this using MPI_Barrier, which blocks all
processes until every process has reached the barrier. This helps synchronize timing more
accurately.

So the following code can be used to time a block of MPI code and report a single elapsed time:

The MPI_Reduce call uses the MPI_MAX operator to find the maximum value among all
local_elapsed times.

Timing results can vary across multiple runs of the same program due to unpredictable
interactions with the operating system. Since these interactions rarely speed up execution, it’s
common practice to report the minimum run-time instead of the average or median.

On hybrid systems with multicore nodes, running one MPI process per node can reduce
interconnect contention and improve run-times. It may also lead to more consistent timing
results.

Dept. of CSE- GMIT PROF. SMITHA M S 36

 Parallel Computing [BCS702]

3.6.2 Results

The timing results for the matrix-vector multiplication program (Table 3.5) show run-times in
milliseconds for square matrices. When using one process (comm_sz = 1), the program runs
serially on a single core. As expected, increasing the matrix size (n) increases the run-time. For
a small number of processes, doubling n roughly quadruples the run-time, but this pattern doesn't
hold for larger numbers of processes.

Table 3.5 Run-times of serial and parallel matrix-vector multiplication (times are in
milliseconds).

When we keep the matrix size (n) fixed and increase the number of processes (comm_sz), the
run-time usually goes down. For large n, doubling the processes can almost cut the run-time in
half. But for small n, increasing comm_sz offers little to no benefit—e.g., using 8 or 16 processes
with n = 1024 gives the same run-time. This pattern is common in parallel programs: as problem
size increases, run-time increases, and adding processes helps only up to a point. Beyond that,
more processes can actually make things slower due to overhead.

 Serial run-time is denoted as Tserial(n) since it depends on input size n.

 Parallel run-time is denoted as Tparallel (n, p), as it depends on both input size n and
number of processes p (comm_sz).
 Typically, a parallel program divides the serial work across multiple processes.
 It also introduces extra time for coordination, called Toverhead.
 So, the relation often becomes:

Tparallel(n,p) = Tserial(n)/p + Toverhead

In MPI programs, parallel overhead usually comes from communication and varies with both
problem size and number of processes. This applies to the matrix-vector multiplication program,
where the main computation is done in nested for loops.

Dept. of CSE- GMIT PROF. SMITHA M S 37

 Parallel Computing [BCS702]

If we only consider floating point operations, each inner loop does n multiplications and n
additions, totalling 2n operations. Since this loop runs m times, the total becomes 2mn
operations. So, when m = n, the total is 2n² operations for some constant a.

In the parallel version, each process handles an {n/p} x {n} times n matrix, performing n2/p
operations, which means the work is split among the processes. However, before this local
computation, all processes must run MPI_Allgather to collect the needed vector data, adding
communication overhead.

Based on the timing data, when the number of processes pp is small and the problem size nn is
large, most of the run time comes from the actual computation T serial(n)/p. This is shown by the
fact that doubling pp (like from 2 to 4) nearly halves the total run time, meaning overhead is less
significant in such cases.

Also, if we fix p at a small value (e.g., p = 2, 4), then increasing n seems to have approximately
the same effect as increasing n for the serial program. For example,

Dept. of CSE- GMIT PROF. SMITHA M S 38

 Parallel Computing [BCS702]

When the input size nn is large and the number of processes pp is small, the parallel run-time
behaves like Tserial(n)/p, meaning overhead (like from MPI_Allgather) has little effect.
However, when n is small and pp is large, this pattern no longer holds, and overhead can
significantly impact performance.

So, it appears that for small n and large p, the dominant term in our formula for Tparallel is
Tallgather.

3.6.3 Speedup and efficiency

Speedup is a common way to compare serial and parallel performance. It is calculated by

dividing the serial run-time by the parallel run-time.

The ideal speedup S (n, p) is equal to the number of processes p, meaning the parallel program
runs pp times faster than the serial version. While this "linear speedup" is rarely achieved in
practice, our matrix-vector multiplication program came close for small pp and large n. However,
for large pp and small n, the speedup dropped significantly—e.g., only 2.4 speedup for n=1024
and p = 16.

Also recall that another widely used measure of parallel performance is parallel efficiency. This
is “per process” speedup:

Dept. of CSE- GMIT PROF. SMITHA M S 39

 Parallel Computing [BCS702]

So linear speedup corresponds to a parallel efficiency of p/p = 1.0, and, in general, we expect
that our efficiencies will usually be less than 1.

Table 3.7 Efficiencies of parallel matrix-vector multiplication.

The efficiency results in Table 3.7 show that when the number of processes pp is small and
the problem size n is large, the parallel program runs efficiently. However, when pp is large
and n is small, efficiency drops significantly, meaning the program doesn't make good use of all
available processes.

3.6.4 Scalability

 Our matrix-vector multiplication program doesn’t show linear speedup for small problem
size n and large number of processes p, but that doesn’t mean it’s a bad program.
 Scalability refers to how well a program maintains efficiency as the number of processes
increases.
 Strong scalability means efficiency stays constant even if problem size stays the same.
 Weak scalability means efficiency stays constant if we increase the problem size at the
same rate as the number of processes.
 Example:
o Program A: Always has 0.75 efficiency no matter the size — strongly scalable.
o Program B: Efficiency = n/625p — weakly scalable since efficiency stays the
same if n increases with p.
 Our matrix-vector program behaves more like Program B:
o When both n and p are doubled (for p≥4), efficiency improves or stays the same.
 Thus, the matrix-vector multiplication program is weakly scalable.

Dept. of CSE- GMIT PROF. SMITHA M S 40

 Parallel Computing [BCS702]

3.7 A parallel sorting algorithm

This section describes a sorting algorithm where keys start and end distributed across all
processes. If there are n keys and pp processes (with n divisible by p), each process begins and
ends with n/p keys. Initially, keys are randomly assigned, but after sorting:

 Each process has its keys in increasing order.

 All keys in process q are less than or equal to keys in process r, for q < r.
So, when the keys from all processes are lined up by rank, the full list is sorted. We
assume the keys are standard integers.

3.7.1 Some simple serial sorting algorithms

Bubble sort (see Program 3.14) is a simple serial sorting algorithm where the array a holds the
unsorted values at the start and sorted values at the end. It works by comparing and swapping
adjacent elements if they’re out of order. Each pass pushes the largest unsorted element to its
correct position at the end of the list. As the outer loop progresses, fewer elements need to be
checked, since the largest values settle into place.

Parallelizing bubble sort isn’t useful because its comparisons must follow a strict order to work
correctly. For example, if a[i−1] = 9, a[i] = 5, and a[i+1] = 7, comparing and swapping in the
wrong order leads to incorrect results. The sequence 5, 7, 9 is only achieved if swaps are done
sequentially: first between 9 and 5, then 9 and 7. Doing them out of order can result in 5, 9, 7,
which is incorrect.

Program 3.14: Serial bubble sort.

Dept. of CSE- GMIT PROF. SMITHA M S 41

 Parallel Computing [BCS702]

Odd-even transposition sort is a variation of bubble sort that allows more parallelism. It works
by breaking the compare-swap operations into phases. There are two types of phases: in even
phases, compare-swaps happen between even-indexed pairs, and in odd phases, they happen
between odd-indexed pairs. This separation allows independent swaps to run in parallel.

and during odd phases, compare-swaps are executed on the pairs

Here’s a small example:

Start: 5, 9, 4, 3

Even phase: Compare-swap (5, 9) and (4, 3), getting the list 5, 9, 3, 4.

Odd phase: Compare-swap (9, 3), getting the list 5, 3, 9, 4.

Even phase: Compare-swap (5, 3), and (9, 4) getting the list 3, 5, 4, 9.

Odd phase: Compare-swap (5, 4) getting the list 3, 4, 5, 9.

This example required four phases to sort a four-element list. In general, it may require fewer
phases, but the following theorem guarantees that we can sort a list of n elements in at most n
phases

Theorem. Suppose A is a list with n keys, and A is the input to the odd-even transposition sort
algorithm. Then after n phases, A will be sorted.

Program 3.15: Serial odd-even transposition sort.

Dept. of CSE- GMIT PROF. SMITHA M S 42

 Parallel Computing [BCS702]

3.7.2 Parallel odd-even transposition sort

Odd-even transposition sort allows greater parallelism than bubble sort, as all compare-swaps in
a phase can run at the same time. Using Foster’s methodology, we define each task as computing
the value of a[i] at the end of phase j. Each task communicates with its neighbouring tasks
(a[i−1] or a[i+1]) and must retain its value for the next phase. This task structure is shown in Fig.
3.12.

Fig 3.12: Communications among tasks in odd-even sort. Tasks determining a[k] are labeled
with a[k].

 Each process is initially assigned n/p keys, and aggregation/mapping is partially guided
by this distribution.
 When n = p, Fig. 3.12 helps visualize the algorithm: each process exchanges data with
its neighbor (i−1 or i+1) depending on the phase, and updates its local value accordingly.
 However, sorting when n = p is impractical, as sorting a few thousand elements is easy
for a single processor.
 Even with many processors, the communication overhead from exchanging messages for
each compare-exchange outweighs the benefits, making the program inefficient.
 Thus, communication cost often dominates local computation in parallel sorting. When
each process holds n/p > 1 elements, we start by applying a fast serial sort (e.g., qsort)
locally within each process.
 When each process holds multiple keys (n/p > 1), we first sort the local keys using a
fast serial algorithm like qsort.
 Referring to Table 3.8, with p = 4 and n = 16, each process holds 4 keys.
 In phase 0, processes (0 & 1) and (2 & 3) exchange all their elements. Each pair then
splits the combined keys: the lower half goes to the lower-ranked process, the upper half
to the higher.

Dept. of CSE- GMIT PROF. SMITHA M S 43

 Parallel Computing [BCS702]

 In phase 1, processes 1 and 2 exchange and redistribute similarly, while 0 and 3 are idle.
 After two more alternating phases, each process ends up with keys sorted locally, and
globally ordered so that all keys in process q are ≤ those in process r, for q < r.

Table 3.8 Parallel odd-even transposition sort.

Theorem. If parallel odd-even transposition sort is run with p processes, then after p phases, the
input list will be sorted.

The parallel algorithm is clear to a human computer:

3.7.3 Safety in MPI programs

If a process is not idle, we might try to implement the communication with a call to MPI_Send
and a call to MPI_Recv:

This behavior can cause the program to hang or crash. MPI_Send may either buffer the message
and return or block until MPI_Recv is called. Many MPI implementations switch from buffering

Dept. of CSE- GMIT PROF. SMITHA M S 44

 Parallel Computing [BCS702]

to blocking based on a message size threshold—small messages are buffered, while large ones
cause MPI_Send to block. If all processes call a blocking MPI_Send before any MPI_Recv
begins, none can proceed, resulting in a deadlock. Programs relying on MPI's buffering are
unsafe—they may work for small values of n, but can hang or crash for larger values. So, we
must ask:

1. How do we detect if a program is safe?

2. How can we make parallel odd-even sort communication safe?

To answer to first question, To check if a program is safe, we can replace MPI_Send with
MPI_Ssend, which is a synchronous send and always blocks until the matching MPI_Recv starts.
If the program runs without hanging or crashing using MPI_Ssend, then the original program
using MPI_Send was safe. Both functions use the same arguments.

To answer to the second question, Make a program safe, the communication pattern must be
restructured to avoid all processes sending first and then receiving. This simultaneous sending—
like in our partner exchanges or the “ring pass” where each process q sends to (q + 1) %
comm_sz—can lead to deadlock if the sends block. A safer approach is to have some processes
send while others receive, based on their rank (e.g., even ranks send first, odd ranks receive first).

In both settings, we need to restructure the communications so that some of the processes receive
before sending. For example, the preceding communications could be restructured as follows:

Dept. of CSE- GMIT PROF. SMITHA M S 45

 Parallel Computing [BCS702]

This communication scheme works clearly when comm_sz is even. For example, with comm_sz
= 4, processes 0 and 2 send to 1 and 3, while 1 and 3 receive—then they reverse roles. It’s less
obvious for odd comm_sz, but as shown in Fig. 3.13 (for comm_sz = 5), the alternating send-
receive pattern still allows safe communication by ensuring some processes always receive
before sending.

MPI provides an alternative to scheduling the communications ourselves—we can call the
function MPI_Sendrecv:

Fig 3.13: Safe communication with five processes.

Dept. of CSE- GMIT PROF. SMITHA M S 46

 Parallel Computing [BCS702]

The MPI_Sendrecv function performs a blocking send and receive in one call, allowing the
source and destination to be the same or different. It’s useful because MPI handles the
communication scheduling, preventing deadlocks. This simplifies code by replacing complex
send-receive logic (like odd/even checks) with a single, safe function call.

3.7.4 Final details of parallel odd-even sort

Recall that we had developed the following parallel odd-even transposition sort algorithm:

Dept. of CSE- GMIT PROF. SMITHA M S 47

 Parallel Computing [BCS702]

In light of our discussion of safety in MPI, it probably makes sense to implement the send and
the receive with a single call to MPI_Sendrecv:

 To keep the smallest n/p keys from 2n/p keys, instead of sorting both lists, we can merge
two sorted lists and stop after n/p elements — improving efficiency (see Program 3.16).
 To get the largest n/p keys, we simply reverse the merge—start from the end and move
backward.
 Swapping pointers instead of copying arrays gives one final optimization.
 Table 3.9 shows the run-times using these improvements — for a single process, serial
quicksort is used instead of slower odd-even sort.

Program 3.16: The Merge_low function in parallel odd-even transposition sort.

Table 3.9 Run-times of parallel odd-even sort (times are in milliseconds).

Dept. of CSE- GMIT PROF. SMITHA M S 48