Chapter 14: Parallelization — Shared & Distributed Memory, OpenMP, MPI¶

Source slides: V12_Parallelization.pdf. Examples: examples/example01.cpp (MPI hello), examples/example02.cpp (loop split), examples/example03_Cstring.cpp (MPI_Send/MPI_Recv with strings), examples/example03_int.cpp (Send/Recv with ints).

1. Chapter Overview¶

Parallelization is the whole point of HPC. This chapter explains:

Levels of parallelism — bit-level, instruction-level, data-level (SIMD), task-level (multi-thread), node-level (multi-process).
Shared-memory model: many threads see the same RAM. Standard API: OpenMP.
Distributed-memory model: many processes each with private memory; communicate via messages. Standard API: MPI.
Hybrid: MPI between nodes + OpenMP within a node.
Parallel performance: Amdahl's law, strong vs weak scaling, embarrassingly parallel vs inherently serial.
Domain decomposition for CFD: split the mesh into sub-domains; each process owns one; halo (ghost) cells store neighbour data; one MPI exchange per time-step.
Synchronization, race conditions, atomic, critical, barrier.
Compiling and running: mpicc/mpic++, mpirun -n N.

Why it matters in HPC/CFD: every modern solver uses MPI domain decomposition (often + OpenMP inside each rank, or + GPU offload). Without parallel programming, CFD does not scale.

What the examiner asks (very common):

"Difference between shared and distributed memory."
"Compare OpenMP and MPI."
"What is Amdahl's law?"
"Explain MPI_Init, MPI_Comm_rank, MPI_Comm_size, MPI_Finalize."
"Predict the output of mpirun -n 4 ./example01."
"Explain MPI_Send/MPI_Recv arguments."
"What is a halo cell? Why is it needed?"
"Race condition — what is it and how to fix?"

What you must master for top grade:

Memory-model picture: shared vs distributed.
Amdahl's formula.
The MPI 6-call backbone: Init, Comm_rank, Comm_size, Send, Recv, Finalize.
The OpenMP #pragma omp parallel for and reduction(+:sum).
Compile + run: mpic++ -O2 file.cpp -o exe; mpirun -n 4 ./exe.
Domain-decomposition diagram with halo cells.

2. Basics from Zero¶

A modern computer has many CPU cores; a cluster has many such computers ("nodes") connected by a network. To use them all you must split the work.

Two main programming models:

Shared memory (within one node). All threads see the same arrays. You write loops; OpenMP's #pragma omp parallel for distributes iterations to threads.

#pragma omp parallel for reduction(+:sum)
for (int i=0; i<N; ++i) sum += a[i];

বাংলায়: Shared memory মানে একই node-এর সব core একটাই RAM দেখে — তাই thread-গুলো আলাদা করে data পাঠায় না, সরাসরি একই array পড়ে-লেখে। OpenMP-তে শুধু loop-এর উপরে একটা pragma বসালেই iteration-গুলো thread-দের মধ্যে ভাগ হয়ে যায়। কিন্তু সবাই একই memory লেখে বলেই race condition-এর বিপদ আসে — এজন্যই reduction clause লাগে। পরীক্ষায় "OpenMP কোথায় কাজ করে" প্রশ্নের উত্তর: শুধু এক node-এর ভেতরে।

Distributed memory (between nodes). Each process has its own memory; if rank 0 needs rank 1's data, it must receive a message. Standard library: MPI.

int rank;
MPI_Init(&argc,&argv);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
if (rank==0) MPI_Send(&x,1,MPI_INT,1,0,MPI_COMM_WORLD);
if (rank==1) MPI_Recv(&x,1,MPI_INT,0,0,MPI_COMM_WORLD,&status);
MPI_Finalize();

বাংলায়: Distributed memory-তে প্রতিটা process (যাকে MPI বলে rank) সম্পূর্ণ আলাদা memory নিয়ে চলে — rank 0-এর variable rank 1 কখনোই সরাসরি দেখতে পায় না। data দরকার হলে explicit message পাঠাতে হয়: একদিকে MPI_Send, অন্যদিকে মিলিয়ে MPI_Recv। এটাই MPI-র মূল দর্শন, আর এজন্যই MPI হাজার হাজার node-এ scale করে যেখানে OpenMP এক node-এই আটকে থাকে। পরীক্ষায় rank শব্দটার মানে (process-এর ID, 0 থেকে size-1) জিজ্ঞেস করা প্রায় নিশ্চিত।

Real-life analogy. Shared memory = a kitchen with multiple chefs sharing one fridge. Distributed memory = ten kitchens in different cities — chefs phone (Send/Recv) to coordinate.

Real-life HPC use. A CFD solver decomposes a mesh of 1 billion cells into 4 096 sub-domains, one per MPI rank. Each rank computes its sub-domain; once per step they exchange the halo (ghost) layer with neighbours.

What if you misunderstand? You write OpenMP code expecting it to scale across nodes — it can't (only within one node). Or you forget MPI_Init → segfault. Or you let two threads update sum without reduction → wrong sum (race condition).

বাংলায়: তিনটা classic ভুল মুখস্থ রাখো: (১) OpenMP দিয়ে দুই node ব্যবহারের চেষ্টা — অসম্ভব, কারণ দ্বিতীয় node-এর RAM প্রথম node-এর thread দেখতেই পায় না; (২) MPI_Init ভুলে যাওয়া — তখন প্রথম MPI call-এই program ভেঙে পড়ে; (৩) reduction ছাড়া দুই thread-এ একই sum-এ লেখা — উত্তর প্রতিবার আলাদা আসে, এটাই race condition। এই তিনটা trap-ই পরীক্ষায় "what goes wrong here" আকারে আসে।

3. Hard English Made Easy¶

Hard Term	Simple English	বাংলা	Example
Parallelization	Doing many things at once	একসাথে অনেক কাজ	OpenMP/MPI
Speed-up	\(T_s/T_p\)	দ্রুতি বৃদ্ধির অনুপাত	4 cores → ~3.6×
Efficiency	Speed-up / N	দক্ষতা	0.9 = 90%
Strong scaling	Same work, more processors	একই কাজ, বেশি প্রসেসর	speed-up vs N
Weak scaling	Work per processor fixed	প্রতি প্রসেসরে একই কাজ	constant runtime
Amdahl's law	Limit of speed-up	আমডালের সূত্র	\(1/(f+(1-f)/N)\)
Embarrassingly parallel	No communication needed	যোগাযোগহীন সমান্তরাল	Monte Carlo
Shared memory	One RAM, many cores	এক মেমরি, অনেক কোর	OpenMP
Distributed memory	Many RAMs	অনেক মেমরি	MPI
Process	OS-level program	প্রসেস	one MPI rank
Thread	Lightweight execution stream	থ্রেড	OpenMP thread
Rank	MPI process ID	এমপিআই আইডি	0..N-1
Communicator	MPI process group	প্রসেস গ্রুপ	`MPI_COMM_WORLD`
Send / Receive	Inter-rank message	বার্তা পাঠানো / পাওয়া	`MPI_Send/Recv`
Broadcast	One → all	এক থেকে সবাই	`MPI_Bcast`
Scatter	One → distinct each	এক থেকে আলাদা সবাই	`MPI_Scatter`
Gather	All → one	সবাই থেকে এক	`MPI_Gather`
Reduce	Combine across ranks	একত্রিত গণনা	`MPI_Reduce`
Barrier	Wait for all	সবার জন্য অপেক্ষা	`MPI_Barrier`
Race condition	Concurrent unsafe write	কনকারেন্ট ভুল	two threads `++sum`
Critical section	Exclusive region	এক্সক্লুসিভ এলাকা	`#pragma omp critical`
Atomic	Indivisible op	অবিভাজ্য অপারেশন	`#pragma omp atomic`
Deadlock	Mutual wait	পারস্পরিক অপেক্ষা	A waits B, B waits A
Halo / ghost cell	Copy of neighbour cell	প্রতিবেশী কোষের কপি	CFD halo layer
Domain decomposition	Split mesh into pieces	মেশ ভাগ করা	METIS / ParMetis
SIMD	Single Instruction Multiple Data	এক নির্দেশ অনেক ডেটা	AVX-512
Hybrid	MPI + OpenMP	এমপিআই + ওপেনএমপি	per-node threads

4. Deep Theory Explanation¶

4.1 Levels of parallelism (slide 4)¶

Level	Example
Bit	64-bit ALU vs 8-bit
Instruction (ILP)	pipelining, superscalar
Data (SIMD)	AVX-512
Task / thread	OpenMP
Node / process	MPI

4.2 Memory models (slide 5)¶

              SHARED MEMORY  (one node — OpenMP)
┌────────────────────────────────────────────────────────────┐
│  ┌───────┐  ┌───────┐  ┌───────┐  ┌───────┐                │
│  │ Core0 │  │ Core1 │  │ Core2 │  │ Core3 │   threads      │
│  └───┬───┘  └───┬───┘  └───┬───┘  └───┬───┘                │
│      └──────────┴────┬─────┴──────────┘                    │
│                      ▼                                     │
│           ┌─────────────────────┐    every thread reads/   │
│           │      ONE  RAM       │◄── writes the SAME a[i]  │
│           │  (one address space)│    danger: race condition│
│           └─────────────────────┘                          │
└────────────────────────────────────────────────────────────┘

           DISTRIBUTED MEMORY  (many nodes — MPI)
┌───────────────┐   ┌───────────────┐   ┌───────────────┐
│    Node 0     │   │    Node 1     │   │    Node 2     │
│  cores+cores  │   │  cores+cores  │   │  cores+cores  │
│ ┌───────────┐ │   │ ┌───────────┐ │   │ ┌───────────┐ │
│ │PRIVATE RAM│ │   │ │PRIVATE RAM│ │   │ │PRIVATE RAM│ │
│ └─────┬─────┘ │   │ └─────┬─────┘ │   │ └─────┬─────┘ │
└───────┼───────┘   └───────┼───────┘   └───────┼───────┘
        └─────────────┬─────┴───────────────────┘
              ┌───────▼────────┐   no rank can see another
              │  INTERCONNECT  │◄─ rank's RAM; data moves only
              │ (network/IB)   │   as messages: MPI_Send/Recv
              └────────────────┘

Shared memory: thread-level parallelism, atomic / critical / reduction needed for correctness. Distributed memory: explicit message passing, scales to many nodes.

বাংলায়: উপরের ছবিটাই এই chapter-এর মেরুদণ্ড: উপরে এক RAM-এ চারটা core (OpenMP-র জগৎ), নিচে তিনটা node যাদের প্রত্যেকের নিজস্ব private RAM আর মাঝখানে network (MPI-র জগৎ)। Shared memory-তে সমস্যা হলো synchronisation (race), distributed memory-তে সমস্যা হলো data ভাগ করা আর message পাঠানোর খরচ। পরীক্ষায় এই diagram আঁকতে বললে দুটো জিনিস অবশ্যই দেখাবে: এক RAM বনাম অনেক private RAM, আর node-গুলোর মাঝের interconnect।

4.3 Amdahl's law — derived from first principles¶

Definitions (write these first in any exam answer). Let \(T_1\) be the runtime on one processor and \(T_N\) the runtime on \(N\) processors. Then

\[ S(N) = \frac{T_1}{T_N} \quad \text{(speed-up)}, \qquad E(N) = \frac{S(N)}{N} \quad \text{(parallel efficiency, ideal } E=1\text{)}. \]

Derivation. Split the single-processor runtime into a serial fraction \(f\) (cannot be parallelised: I/O, mesh reading, sequential algorithms) and a parallelisable fraction \(1-f\). On \(N\) processors the serial part is unchanged, while the parallel part is ideally divided by \(N\):

\[ T_N = f\,T_1 + \frac{(1-f)\,T_1}{N}. \]

Divide \(T_1\) by this:

\[ S(N) = \frac{T_1}{T_N} = \frac{T_1}{f\,T_1 + (1-f)\,T_1/N} = \frac{1}{f + \dfrac{1-f}{N}}. \]

Taking the limit \(N \to \infty\) makes the second term vanish:

\[ S_{\max} = \lim_{N\to\infty} S(N) = \frac{1}{f}. \]

No matter how many processors you buy, the serial 5 % caps you at 20×. Examples: \(f=0.05\) → max speed-up = 20×; \(f=0.01\) → 100×. CFD codes try to keep \(f<1\%\).

Fully worked example (\(f = 0.05\), \(T_1 = 100\) s). Serial part: \(0.05 \times 100 = 5\) s, parallel part \(95\) s.

For \(N=4\): \(T_4 = 5 + 95/4 = 5 + 23.75 = 28.75\) s, so \(S(4) = 100/28.75 = 3.48\) and \(E(4) = 3.48/4 = 0.87\).

Complete table (all rows computed the same way):

\(N\)	\(T_N = 5 + 95/N\) (s)	\(S(N) = 100/T_N\)	\(E(N) = S/N\)
2	\(5 + 47.50 = 52.50\)	1.90	0.95
4	\(5 + 23.75 = 28.75\)	3.48	0.87
16	\(5 + 5.94 = 10.94\)	9.14	0.57
64	\(5 + 1.48 = 6.48\)	15.42	0.24
256	\(5 + 0.37 = 5.37\)	18.62	0.073
\(\infty\)	\(5 + 0 = 5\)	20.00	\(\to 0\)

Read the last column: at 256 cores you pay for 256 but effectively use 19 of them — efficiency collapses long before the speed-up plateau is reached.

 S(N)  speed-up vs. number of processors (f = 0.05)
  64 ┤                                                ·
     │                                          ·   ideal S = N
  48 ┤                                    ·   (linear)
     │                              ·
  32 ┤                        ·
     │                  ·
  20 ┤┄┄┄┄┄┄┄┄┄┄┄·┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄┄ S_max = 1/f = 20
  16 ┤        ·                  ●━━━━━━━━●━━━━━ Amdahl curve
     │      ·             ●━━━━━┘(64:15.4) (256:18.6)
   8 ┤    ·         ●━━━━┘(16: 9.1)
   4 ┤  · ●━━━━━━━━┘(4: 3.5)
   2 ┤ ●┘(2: 1.9)
   1 ┤●
     └┴───┴────┴───────┴──────────┴──────────────► N (log axis)
      1   2    4      16         64        256
      ● = Amdahl S(N)=1/(0.05+0.95/N)   · = ideal S=N
      The curve PLATEAUS at 20 — more processors buy nothing.

বাংলায়: Amdahl-এর যুক্তিটা খুব সহজ: program-এর serial অংশ f কোনোদিন ভাগ হয় না, শুধু বাকি (1-f) অংশটা N processor-এ ভাগ হয় — তাই total সময় কখনো f×T1-এর নিচে নামে না, আর speed-up আটকে যায় 1/f-এ। f=0.05 মানে মাত্র 5% serial code, তবু হাজার core দিলেও 20×-এর বেশি কিছুই পাবে না। পরীক্ষায় derivation লিখতে বললে আগে T_N = f·T1 + (1-f)·T1/N লাইনটা লেখো, তারপর ভাগ করো — formula মুখস্থ বলার চেয়ে এই দুই-লাইনের derivation-এ বেশি নম্বর। আর efficiency E=S/N-ও সাথে বলো: speed-up বাড়লেও efficiency দ্রুত পড়ে যায়।

4.4 Strong vs weak scaling — and Gustafson's law¶

Formal definitions (state what is held constant!):

Strong scaling: the total problem size \(W\) is held constant while \(N\) varies. Measure \(S_{\text{strong}}(N) = T(W, 1)/T(W, N)\); ideal is \(S = N\) (linear). Governed by Amdahl's law — it eventually saturates.
Weak scaling: the work per processor \(W/N\) is held constant, i.e. the total problem grows proportionally with \(N\) (\(W \propto N\)). Measure the weak-scaling efficiency \(E_{\text{weak}}(N) = T(W,1) / T(N \cdot W, N)\); ideal is a flat runtime curve, \(E_{\text{weak}} = 1\). Governed by Gustafson's law.

Ideal: strong → linear speed-up; weak → constant runtime.

Gustafson's law (derivation). Amdahl assumes the problem stays fixed. Gustafson asks the question the other way round: you run for a fixed wall time on the parallel machine, with normalised runtime \(T_N = f + (1-f)\) (serial part \(f\), parallel part \(1-f\), both measured on the \(N\)-processor run). How long would one processor need for the same scaled job? It would have to do the parallel work \(N\) times sequentially:

\[ T_1 = f + (1-f)\,N. \]

Hence the scaled speed-up

\[ S(N) = \frac{T_1}{T_N} = \frac{f + (1-f)N}{f + (1-f)} = f + (1-f)\,N = N - f\,(N-1). \]

Worked example. \(f = 0.05\), \(N = 64\): \(S = 64 - 0.05 \times 63 = 64 - 3.15 = 60.85\). Compare Amdahl at the same \(f\) and \(N\): \(S = 15.42\) (table in §4.3). The difference is not a contradiction — they answer different questions.

	Amdahl (strong)	Gustafson (weak)
Problem size	fixed	grows with \(N\) (\(W \propto N\))
Question	"same problem faster?"	"bigger problem in same time?"
Formula	\(S = \dfrac{1}{f + (1-f)/N}\)	\(S = N - f(N-1)\)
Limit \(N\to\infty\)	\(1/f\) (plateau)	grows without bound (slope \(1-f\))
\(f=0.05\), \(N=64\)	\(15.42\)	\(60.85\)
CFD reading	refining the same mesh run	finer mesh on more nodes, same walltime

বাংলায়: Strong scaling-এ প্রশ্ন: "একই সমস্যা আরো processor দিলে কত দ্রুত হবে?" — এখানে Amdahl রাজত্ব করে, আর speed-up 1/f-এ আটকে যায়। Weak scaling-এ প্রশ্ন উল্টো: "processor-প্রতি কাজ একই রেখে সমস্যাটাই বড় করলে কী হয়?" — তখন Gustafson-এর সূত্র S = N - f(N-1) প্রায় linear-ভাবে বাড়তে থাকে। পরীক্ষায় সবচেয়ে দামি লাইনটা হলো: কোনটায় কী constant রাখা হয় — strong-এ total problem size, weak-এ per-processor কাজ। দুটো সূত্রে একই f, N বসিয়ে দুটো আলাদা সংখ্যা দেখালে (15.42 বনাম 60.85) পুরো নম্বর নিশ্চিত।

4.5 Domain decomposition (CFD specific)¶

The mesh is split into \(P\) sub-domains, one per MPI rank. Each rank stores a halo (ghost) layer of cells that mirror its neighbours' boundary cells. Each time-step:

Compute fluxes inside the sub-domain.
Exchange halo with neighbours via MPI_Sendrecv or non-blocking MPI_Isend/MPI_Irecv.
Apply boundary conditions.
Update solution.

Tools: METIS / ParMetis partition the mesh.

2D DOMAIN DECOMPOSITION — 4 ranks, each owns a 4×4 block; g = ghost/halo cell

        rank 0                          rank 1
┌───┬───────────────┐          ┌───────────────┬───┐
│ . │ o  o  o  o    │          │    o  o  o  o │ . │
│ . │ o  o  o  o  g │  ◄────►  │ g  o  o  o  o │ . │   vertical halo
│ . │ o  o  o  o  g │ exchange │ g  o  o  o  o │ . │   exchange:
│ . │ o  o  o  o  g │  column  │ g  o  o  o  o │ . │   rank 0's last
└───┴───────┬───────┘          └───────┬───────┴───┘   column is copied
            │ ▲                        │ ▲             into rank 1's g
            ▼ │  horizontal halo       ▼ │             cells, and vice
┌───┬───────┴───────┐          ┌───────┴───────┬───┐   versa — ONE
│ . │ o  o  o  o  g │          │ g  o  o  o  o │ . │   message per
│ . │ o  o  o  o  g │  ◄────►  │ g  o  o  o  o │ . │   neighbour per
│ . │ o  o  o  o  g │          │ g  o  o  o  o │ . │   time-step
│ . │ o  o  o  o    │          │    o  o  o  o │ . │
└───┴───────────────┘          └───────────────┴───┘
        rank 2                          rank 3

 o = cell OWNED and computed by the rank
 g = HALO/GHOST cell: read-only local COPY of the neighbour's
     boundary cells, refreshed once per time-step by MPI exchange
 . = physical boundary (boundary condition, no neighbour rank)

Surface-to-volume communication analysis (why strong scaling dies). Take a cubic mesh of \(n^3\) cells split across \(P\) ranks into sub-cubes of edge \(n/P^{1/3}\). Per rank:

\[ \text{computation} \propto \text{volume} = \frac{n^3}{P}, \qquad \text{communication} \propto \text{halo surface} = 6\left(\frac{n}{P^{1/3}}\right)^{2}. \]

Their ratio is

\[ \frac{\text{comm}}{\text{comp}} = \frac{6\,(n/P^{1/3})^2}{(n/P^{1/3})^3} = \frac{6\,P^{1/3}}{n} \;\propto\; \frac{P^{1/3}}{n}. \]

Adding ranks (\(P\uparrow\)) at fixed mesh (\(n\) fixed) makes the ratio grow — communication overtakes computation.

Worked example (\(n = 512\), \(P = 64\)). \(P^{1/3} = 4\), so each rank owns a sub-cube of edge \(512/4 = 128\):

volume (computed cells): \(128^3 = 2\,097\,152\) cells,
halo surface (communicated cells): \(6 \times 128^2 = 6 \times 16\,384 = 98\,304\) cells,
ratio: \(98\,304 / 2\,097\,152 = 0.047 \approx 4.7\%\) — healthy.

Now grow to \(P = 512\): \(P^{1/3} = 8\), edge \(64\), ratio \(= 6 \times 8 / 512 = 9.4\%\) — multiplying ranks by 8 doubled the relative communication cost (\(P^{1/3}\) doubled). This is the strong-scaling limit in numbers.

1D halo exchange (simplest worked case). \(N = 1000\) cells on \(P = 4\) ranks: each rank owns \(N/P = 250\) cells and allocates \(250 + 2 = 252\), with ghost cells at local index 0 and 251:

rank 1's local array (252 entries):
 index:   0    1    2   ...  250   251
        ┌────┬────┬────┬────┬────┬────┐
        │ G  │ o  │ o  │ .. │ o  │ G  │
        └─┬──┴────┴────┴────┴────┴──┬─┘
          │   owned cells 1..250    │
   copy of rank 0's          copy of rank 2's
   LAST owned cell           FIRST owned cell
   (from MPI_Sendrecv        (from MPI_Sendrecv
    with left neighbour)      with right neighbour)

Per time-step each rank sends 1 boundary cell to each neighbour and receives 1 into each ghost — 2 values communicated against 250 computed; a stencil like \(u_i^{\text{new}} = (u_{i-1} + u_{i+1})/2\) can then be evaluated for all owned cells without further communication.

বাংলায়: Halo (ghost) cell হলো প্রতিবেশী rank-এর boundary cell-এর একটা local কপি — এটা থাকায় প্রতিটা cell-এর হিসাবের সময় বারবার network-এ যেতে হয় না, প্রতি time-step-এ মাত্র একবার exchange করলেই চলে। আর surface-to-volume হিসাবটা মনে রাখো: কাজ বাড়ে volume-এর সাথে (\(n^{3}/P\)), কিন্তু communication বাড়ে surface-এর সাথে — ratio হলো \(6\cdot P^{1/3}/n\)। P বাড়ালে ভাগের subdomain ছোট হয়ে যায়, তখন surface-এর তুলনায় volume কমে যায়, মানে communication-এর ভাগ বেড়ে যায় — এটাই strong scaling মরার আসল কারণ, পরীক্ষায় সংখ্যাসহ (\(512^{3}\), 64 rank → 4.7%) দেখাতে পারলে দারুণ।

4.6 OpenMP (a one-page tour)¶

#include <omp.h>

#pragma omp parallel
{
    int tid = omp_get_thread_num();
    int nt  = omp_get_num_threads();
    #pragma omp single
    std::cout << "Total threads: " << nt << "\n";
}

#pragma omp parallel for reduction(+:sum)
for (int i=0; i<N; ++i) sum += a[i];

#pragma omp parallel for schedule(dynamic, 64)
for (int i=0; i<N; ++i) work(i);

#pragma omp critical
{ shared_log.push_back(x); }

#pragma omp atomic
counter++;

#pragma omp barrier

Compile: g++ -fopenmp prog.cpp -o prog. Run with OMP_NUM_THREADS=8 ./prog.

বাংলায়: OpenMP-র তিনটা সুরক্ষা-অস্ত্রের ক্রম মনে রাখো: reduction সবচেয়ে দ্রুত (প্রতিটা thread নিজের private কপিতে যোগ করে, শেষে একবার মেলায়), atomic মাঝারি (একটা মাত্র update-কে অবিভাজ্য করে), critical সবচেয়ে ধীর (পুরো block-এ এক সময়ে একটাই thread ঢোকে)। sum-জাতীয় কাজে সবসময় reduction লেখো। compile-এ -fopenmp ভুলে গেলে pragma-গুলো নিঃশব্দে উপেক্ষিত হয়ে serial program চলে — কোনো error আসে না, এটা একটা মোক্ষম পরীক্ষা-trap।

4.7 MPI six-call backbone¶

#include <mpi.h>
int main(int argc, char** argv) {
    int rank, size;
    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    MPI_Comm_size(MPI_COMM_WORLD, &size);
    // ... work ...
    MPI_Finalize();
    return 0;
}

Add point-to-point:

if (rank==0) MPI_Send(&x, 1, MPI_INT, 1, 0, MPI_COMM_WORLD);
if (rank==1) MPI_Recv(&x, 1, MPI_INT, 0, 0, MPI_COMM_WORLD, &status);

MPI_Send(buf, count, type, dest, tag, comm) — read it as: send count elements of type starting at buf to rank dest, labelled tag, within communicator comm. MPI_Recv mirrors it with source and a status out-parameter. Send/Recv are blocking: MPI_Send returns when the buffer is reusable; MPI_Recv returns when the message has arrived.

বাংলায়: ছয়টা call-এর গল্পটা এক লাইনে: Init দিয়ে ঢোকা, rank দিয়ে "আমি কে", size দিয়ে "আমরা কয়জন", Send/Recv দিয়ে চিঠি চালাচালি, Finalize দিয়ে বিদায়। Send-এর ছয়টা argument ক্রমে মুখস্থ রাখো — buf, count, type, dest, tag, comm — পরীক্ষায় প্রতিটার মানে আলাদা করে লিখতে বলা হয়।

4.8 Compile & run MPI¶

mpic++ -O2 example01.cpp -o ex01      # wrapper around g++ adding MPI flags
mpirun -n 4 ./ex01                    # launch 4 processes
mpirun --oversubscribe -n 8 ./ex01    # more ranks than cores (testing)

mpic++ --showme reveals the underlying g++ command with -I/-L/-l MPI options.

4.9 Sequential vs causal consistency (slide 22)¶

Sequential consistency: all processes observe all memory/message operations in one single global order. Causal consistency: only causally-related operations must be seen in order; concurrent ones may appear in different orders to different observers. MPI guarantees ordering only per sender-receiver pair and tag — messages from different senders may interleave arbitrarily; that's why mpirun output order varies run to run.

বাংলায়: mpirun-এর output প্রতিবার এলোমেলো কেন — এই প্রশ্নের উত্তর এখানে: MPI শুধু একই sender→receiver জোড়ার বার্তার ক্রম রক্ষা করে; ভিন্ন rank-এর প্রিন্ট কে আগে পৌঁছাবে তার কোনো নিয়ম নেই। পরীক্ষায় "output order may vary" লিখতে ভুলো না — এটাই নম্বরের লাইন।

4.10 The two remaining must-know diagrams¶

(a) Amdahl speed-up curve (f = 0.05) vs ideal:

 S(N)
 64 ┤                                            · ideal S=N
    │                                       ·
 32 ┤                                ·
    │                          ·
 20 ┤  - - - - - - - - - - - - - - - - - - - - - - - -  S_max = 1/f = 20
 16 ┤              ·            ____________------------ Amdahl
    │         ·        ____----
  8 ┤      ·     __---
  4 ┤   · __--
  2 ┤ ·_-
  1 ┼─┬───┬───┬────┬────┬────┬────┬──────────► N
    1 2   4   8   16   32   64  256
    The curve BENDS AWAY from ideal and saturates at 1/f.

(b) MPI Send/Recv timeline — correct pairing vs deadlock:

 CORRECT (0 sends first):                DEADLOCK (both Recv first):
 rank 0          rank 1                  rank 0          rank 1
   │ MPI_Send ─────► MPI_Recv             │ MPI_Recv...    │ MPI_Recv...
   │               │                      │  (waits        │  (waits
   │ MPI_Recv ◄───── MPI_Send             │   forever)     │   forever)
   ▼               ▼                      ▼ nobody ever sends — hang ▼
 time                                    Fix: reorder, or MPI_Sendrecv,
                                         or non-blocking Isend/Irecv+Waitall

5. Command / Syntax / Code Breakdown¶

`#include <mpi.h>` / `#include <omp.h>`¶

Include MPI / OpenMP headers.

`MPI_Init(&argc, &argv)` / `MPI_Finalize()`¶

Begin / end MPI environment. Required.

`MPI_Comm_rank(MPI_COMM_WORLD, &rank)`¶

Get this process's rank (0…size-1).

`MPI_Comm_size(MPI_COMM_WORLD, &size)`¶

Total ranks.

`MPI_Send(buf, count, type, dest, tag, comm)`¶

Blocking send. type is MPI_INT, MPI_DOUBLE, MPI_CHAR, etc.

`MPI_Recv(buf, count, type, source, tag, comm, &status)`¶

Blocking receive. MPI_ANY_SOURCE and MPI_ANY_TAG are wildcards.

`MPI_Bcast(buf, count, type, root, comm)`¶

Root broadcasts to all.

`MPI_Reduce(sendbuf, recvbuf, count, type, op, root, comm)`¶

Combine values across ranks (MPI_SUM, MPI_MAX, MPI_MIN, …).

`mpirun -n N ./prog`¶

Launch N copies of prog connected by MPI.

`#pragma omp parallel for [reduction(+:sum)] [schedule(dynamic,k)]`¶

Distribute loop iterations across threads.

`#pragma omp critical { … }`¶

Mutual exclusion.

`#pragma omp atomic`¶

One-line atomic update.

`#pragma omp barrier`¶

Synchronise all threads.

`OMP_NUM_THREADS=8 ./prog`¶

Set thread count.

6. Mandatory Practical Examples¶

Example 6.1 — MPI hello (example01.cpp)¶

#include <iostream>
#include <mpi.h>
int main(int argc, char** argv) {
    int rank;
    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    std::cout << "Rank number: " << rank << std::endl;
    MPI_Finalize();
    return 0;
}

Compile + run:

mpic++ -O2 example01.cpp -o ex01
mpirun -n 4 ./ex01

Expected output (ORDER MAY VARY):

Rank number: 0
Rank number: 2
Rank number: 1
Rank number: 3

Step-by-Step

MPI_Init initialises the MPI library.
MPI_Comm_rank writes the rank into rank.
4 separate processes run concurrently; each prints its rank.
MPI_Finalize releases resources.

Real-Life Meaning. Foundation of any parallel CFD program.

Written-Exam Relevance. Frequently asked: "Write MPI hello and explain each call."

Example 6.2 — Loop split across ranks (example02.cpp)¶

int numbers = 10;
int per_rank = floor(10.0/n_ranks);
if (numbers % n_ranks > 0) per_rank++;
int my_first = rank*per_rank;
int my_last = my_first + per_rank;
for (int i=my_first; i<my_last; ++i)
    if (i<numbers) std::cout << "Rank "<<rank<<" iter "<<i<<"\n";

mpirun -n 3 distributes 10 iterations across 3 ranks: rank 0 → 0–3, rank 1 → 4–7, rank 2 → 8–9.

Example 6.3 — Send/Recv int (example03_int.cpp)¶

if (n_ranks != 2) { std::cout<<"need 2 ranks\n"; MPI_Finalize(); return 1; }
if (rank==0) {
    int message = 123456789;
    MPI_Send(&message, 1, MPI_INT, 1, 0, MPI_COMM_WORLD);
}
if (rank==1) {
    int message;
    MPI_Status st;
    MPI_Recv(&message, 1, MPI_INT, 0, 0, MPI_COMM_WORLD, &st);
    std::cout << "I received a message: " << message << "\n";
}

(The lecture's literal source uses count=10 for one int — a mistake; the correct count is 1. Mention this in the exam.)

Run with mpirun -n 2 ./ex03i. Output: I received a message: 123456789

Example 6.4 — Send/Recv string (example03_Cstring.cpp)¶

Same structure but with char[16] and MPI_CHAR of count 16.

Example 6.5 — OpenMP parallel sum¶

#include <omp.h>
double sum = 0.0;
#pragma omp parallel for reduction(+:sum)
for (int i=0; i<N; ++i) sum += a[i];

g++ -O3 -fopenmp prog.cpp -o prog; OMP_NUM_THREADS=8 ./prog

Example 6.6 — Race condition fix¶

Wrong:

#pragma omp parallel for
for (int i=0; i<N; ++i) sum += a[i];                // RACE

Right (three options):

#pragma omp parallel for reduction(+:sum)           // best
for (int i=0; i<N; ++i) sum += a[i];

#pragma omp parallel for
for (int i=0; i<N; ++i) {
    #pragma omp atomic
    sum += a[i];
}

#pragma omp parallel for
for (int i=0; i<N; ++i) {
    #pragma omp critical
    sum += a[i];
}

reduction is the fastest of the three; atomic is faster than critical.

বাংলায়: Race fix-এর গতি-ক্রম মুখস্থ রাখো: reduction > atomic > critical। reduction-এ প্রতিটা thread নিজের private sum জমিয়ে শেষে একবার মেলায় — তাই দ্রুততম; critical-এ প্রতি iteration-এ তালা — তাই ধীরতম। পরীক্ষায় "তিনটা উপায় + কোনটা কেন দ্রুত" পুরোটা লিখলেই পূর্ণ নম্বর।

Example 6.7 — Hybrid MPI + OpenMP¶

MPI_Init_thread(&argc,&argv,MPI_THREAD_FUNNELED,&prov);
#pragma omp parallel for
for (int i=local_first; i<local_last; ++i) work(i);
MPI_Allreduce(...);
MPI_Finalize();

Run: mpirun -n 4 -x OMP_NUM_THREADS=8 ./hybrid → 4 ranks × 8 threads = 32 cores total.

7. Real HPC/CFD Workflow¶

# 1. Build the solver
module load gcc openmpi
mpic++ -O3 -fopenmp -DNDEBUG solver.cpp -o solver

# 2. Domain-decompose the mesh
metis_mpiexec -n 256 -p partition mesh.cgns

# 3. Submit
sbatch run.sbatch         # SLURM script (Ch.15)

run.sbatch:

#!/bin/bash
#SBATCH --nodes=8
#SBATCH --ntasks-per-node=8       # 8 MPI ranks per node
#SBATCH --cpus-per-task=16        # 16 OpenMP threads each
#SBATCH --time=24:00:00
export OMP_NUM_THREADS=$SLURM_CPUS_PER_TASK
mpirun -n $SLURM_NTASKS ./solver in.dat

8. Exercises and Solutions¶

The lecture provides only example code, no formal exercise sheet.

Self-made exercise A — MPI sum¶

Write a program where every rank holds the integer equal to its rank, then collect the global sum on rank 0 and print it.

Solution:

#include <iostream>
#include <mpi.h>
int main(int argc,char**argv){
    int rank,size,sum;
    MPI_Init(&argc,&argv);
    MPI_Comm_rank(MPI_COMM_WORLD,&rank);
    MPI_Comm_size(MPI_COMM_WORLD,&size);
    int local = rank;
    MPI_Reduce(&local,&sum,1,MPI_INT,MPI_SUM,0,MPI_COMM_WORLD);
    if (rank==0) std::cout<<"sum = "<<sum<<"\n";   // 0+1+...+(N-1) = N(N-1)/2
    MPI_Finalize();
}

Run with mpirun -n 8 → sum = 28 (check: \(8\cdot7/2 = 28\)).

Self-made exercise B — OpenMP threaded sort¶

Sort a vector with parallel for then merge:

#pragma omp parallel for
for (int b=0; b<NB; ++b) sort_block(b);
merge_blocks();

Marking schemes. MPI hello (8 marks): include + Init + rank + size + body + Finalize + compile mpic++ + run mpirun. OpenMP reduction (5 marks): include omp.h + pragma + reduction clause + compile -fopenmp + correct sum.

Common mistakes.

Forgetting MPI_Init → "Attempting to use MPI before MPI_INIT".
Mismatched send/recv counts → deadlock or wrong size.
Race conditions when forgetting reduction.
Using OpenMP across nodes (impossible).
Mismatching MPI_INT vs MPI_DOUBLE.

Harder version. Replace the loop with non-blocking MPI_Isend/MPI_Irecv + MPI_Waitall.

9. Written Exam Focus¶

9.1 Short Answers¶

Q. Difference between OpenMP and MPI. A. OpenMP is a directive-based shared-memory model: many threads in one process share memory; great within a node. MPI is a message-passing distributed-memory model: many processes (possibly on many nodes) coordinate via explicit Send/Recv. Modern CFD uses hybrid MPI+OpenMP.

Q. State Amdahl's law. A. \(S(N)=\dfrac{1}{f+(1-f)/N}\) with \(f\) the serial fraction; max speed-up \(1/f\) regardless of cores.

Q. What is a halo cell? A. A copy of a neighbouring sub-domain's boundary cells stored locally so each MPI rank can compute its fluxes without per-cell communication.

Q. What does MPI_Comm_rank return? A. The integer ID (0…size-1) of the calling process within the given communicator.

Q. What is a race condition? A. Two or more threads access the same memory concurrently, at least one writes, with nondeterministic ordering → wrong results.

9.2 Medium Answers¶

Q. (8 marks) Compare shared- and distributed-memory programming models.

A. Shared memory: many threads share one address space; communication via shared variables; the hard problem is synchronisation (races, mutual exclusion). API: OpenMP. Distributed memory: each process has its own address space; communication via explicit messages; the hard problem is data distribution and halo exchange. API: MPI. Shared scales only within a node; distributed scales across thousands of nodes. Hybrid combines both.

Q. (5 marks) Explain MPI_Send arguments.

A. MPI_Send(buf, count, type, dest, tag, comm): buf — pointer to data; count — number of elements; type — MPI datatype (MPI_INT, MPI_DOUBLE, MPI_CHAR); dest — destination rank; tag — integer label to match Recvs; comm — communicator (MPI_COMM_WORLD). Returns when the buffer can be reused (blocking).

9.3 Long Answer (12 marks)¶

Q. Discuss MPI in CFD: domain decomposition, halo exchange, scaling.

A.

Introduction. MPI underlies essentially every production CFD code, enabling 1000+-node simulations.

Main concept. The mesh is partitioned into \(P\) sub-domains (one per rank). Each rank stores its own cells plus a halo of neighbour cells. Per time-step: compute → halo exchange → next.

Step-by-step.

MPI_Init, get rank and size.
Read partition; allocate sub-domain arrays + halo.
Time loop: post MPI_Isend/MPI_Irecv for boundary cells; compute interior; MPI_Waitall; compute boundary updates.
Periodically MPI_Allreduce the global residual.
MPI_Finalize.

Scaling. Strong scaling fails when sub-domains get too small (communication > compute, ratio \(\propto P^{1/3}/n\)); weak scaling works when problem size grows with \(N\).

Conclusion. Domain decomposition + MPI is the dominant paradigm for distributed-memory CFD.

9.4 Output Prediction¶

$ mpirun -n 4 ./example01
Rank number: 0
Rank number: 1
Rank number: 2
Rank number: 3

(any order!)

$ mpirun -n 2 ./example03_int
I received a message: 123456789

9.5 Comparison¶

	OpenMP	MPI
Memory	shared	distributed
API	pragmas	function calls
Scope	within node	across nodes
Compile	`-fopenmp`	`mpic++`
Run	`OMP_NUM_THREADS`	`mpirun -n`
Sync	barrier/critical	explicit messages
Difficulty	easy	medium

	Race	Deadlock
Cause	concurrent unsync write	mutual wait
Symptom	wrong result	program hangs
Fix	reduction/atomic/critical	order locks/messages

9.6 Templates¶

OpenMP template: #include <omp.h>; #pragma omp parallel for reduction(+:sum) over the loop.

MPI template: MPI_Init; MPI_Comm_rank; MPI_Comm_size; …; MPI_Finalize.

Domain template: "decompose mesh → assign sub-domain → halo cells → exchange + compute loop."

9.7 Marking Scheme — "MPI Send/Recv" (5 marks)¶

1 each: Init / Comm_rank / Send / Recv / Finalize.

10. Very Hard Questions¶

Beginner

What is OpenMP for? → shared-memory parallelism.
What is MPI for? → distributed-memory.
Compile MPI? → mpic++ -O2 a.cpp -o a.
Run with 4 ranks? → mpirun -n 4 ./a.
Default communicator? → MPI_COMM_WORLD.

Intermediate

State Amdahl's law. → see 9.1.
Difference rank vs size. → ID vs total count.
Why use reduction in OpenMP? → avoid race.
Datatype for double? → MPI_DOUBLE.
Tag in MPI_Send? → message label for matching.

Hard

MPI_Send size doesn't match MPI_Recv? → error or hang (truncation if Recv smaller).
Blocking vs non-blocking. → blocking returns when safe; non-blocking returns immediately, must MPI_Wait.
What is MPI_Allreduce? → reduce + broadcast; every rank gets the result.
Why halo cells? → batch neighbour data once per step instead of per-cell communication.
When is schedule(dynamic) better than static? → when iteration cost varies.

Very Hard

Why does adding ranks beyond a point slow down? → communication > compute (strong-scaling limit).
How to avoid load imbalance in CFD? → ParMETIS partitioning by element count/cost.
MPI_Sendrecv vs separate Send+Recv. → combined exchange; avoids pairwise deadlock.

Deep Integration

MPI + OpenMP + GPU? → MPI between nodes, OpenMP within node, CUDA/OpenACC on GPU.
Why do GPU CFD codes still need MPI? → multi-GPU multi-node runs.

Coding/Command

OpenMP doubling each element:

#pragma omp parallel for
for (int i=0; i<N; ++i) a[i] *= 2;

MPI broadcast of an int from rank 0:

int x;
if (rank==0) x = 42;
MPI_Bcast(&x, 1, MPI_INT, 0, MPI_COMM_WORLD);

Debugging

Program hangs on 2 ranks waiting for Recv. → deadlock; both Recv first; reorder or MPI_Sendrecv.
OpenMP race produces wrong sum. → add reduction(+:sum).

Long Written

(250 words) Explain why hybrid MPI + OpenMP is the dominant CFD paradigm. (Use §4.2, §7.)

11. Debugging and Mistake Analysis¶

Mistake	Why wrong	Correct	Explanation
Compile with g++ (no MPI)	undefined `MPI_Init`	`mpic++` wrapper	brings in MPI flags
Forgot `MPI_Finalize`	resources leak	call at end	etiquette
Mismatched Send/Recv types	undefined behaviour	match exactly	data interpretation
Forgot reduction in OpenMP	wrong sum	add it	race
OpenMP for inter-node	doesn't span nodes	use MPI	model mismatch
Recv before Send (both ranks)	deadlock	`MPI_Sendrecv`	ordering
Wrong tag matching	hang or wrong msg	match tags	identify stream
Uninitialised halo at step 1	wrong fluxes	init from BC	halo bootstrap
MPI call before `MPI_Init`	error	Init always first	order
`mpirun` without `-n`	OS default (often 1)	always specify	clarity

বাংলায়: Parallel bug-এর তিন মহারাজ: race (ফল ভুল, প্রতিবার আলাদা), deadlock (প্রোগ্রাম ঝুলে থাকে), আর mismatch (type/count/tag না মেলা)। উপসর্গ শুনেই রোগ চেনা যায় — "wrong result" = race, "hangs" = deadlock — পরীক্ষার diagnosis প্রশ্নে এই ম্যাপিংটাই উত্তর।

12. Mini Project for Mastery¶

Goal: Compute π via Monte-Carlo with MPI.

#include <iostream>
#include <random>
#include <mpi.h>

int main(int argc,char**argv){
    int rank,size;
    MPI_Init(&argc,&argv);
    MPI_Comm_rank(MPI_COMM_WORLD,&rank);
    MPI_Comm_size(MPI_COMM_WORLD,&size);

    long N = 1000000000L / size;
    std::mt19937 g(rank+12345);
    std::uniform_real_distribution<double> u(0,1);
    long count = 0;
    for (long i=0;i<N;++i) {
        double x=u(g), y=u(g);
        if (x*x+y*y<=1.0) ++count;
    }

    long total;
    MPI_Reduce(&count,&total,1,MPI_LONG,MPI_SUM,0,MPI_COMM_WORLD);
    if (rank==0) std::cout << "pi = " << 4.0*total/( (double)N*size) << "\n";
    MPI_Finalize();
}

mpic++ -O3 mc.cpp -o mc; mpirun -n 8 ./mc → ~3.1416.

Connection to exam: Init/Rank/Size/Reduce/Finalize, no shared state — the canonical embarrassingly parallel example. The math: area of quarter circle/area of square = π/4, so π ≈ 4·(hits/total).

13. Final Chapter Cheat Sheet¶

Item	Memorise
Levels	bit / instr / data / task / node
Shared mem	OpenMP
Dist. mem	MPI
Hybrid	MPI + OpenMP
Speed-up / Efficiency	\(S = T_1/T_N\), \(E = S/N\)
Amdahl	\(S=1/(f+(1-f)/N)\), \(S_{max}=1/f\)
Gustafson	\(S = N - f(N-1)\)
Strong / weak	fixed total work / fixed work-per-processor
Comm-to-compute	\(\propto P^{1/3}/n\) for 3D decomposition
OpenMP compile	`g++ -fopenmp`
OpenMP loop	`#pragma omp parallel for reduction(+:sum)`
OpenMP run	`OMP_NUM_THREADS=8 ./prog`
MPI compile	`mpic++`
MPI 6 calls	Init / Comm_rank / Comm_size / Send / Recv / Finalize
Send args	buf,count,type,dest,tag,comm
Collectives	Bcast, Reduce, Allreduce, Scatter, Gather, Barrier
Race	unsync writes
Deadlock	mutual wait
Halo cell	neighbour copy
Domain decomp	METIS / ParMETIS
Trap	OpenMP across nodes (impossible)
Top phrase	"OpenMP within a node, MPI across nodes — domain decomposition + halo exchange is the CFD canonical pattern."

14. Mock Exam — Four Levels¶

Level 1 — Basic (definitions & syntax)¶

Q1. Write the six MPI backbone calls in the order they appear in a minimal program.

Solution: MPI_Init → MPI_Comm_rank → MPI_Comm_size → MPI_Send / MPI_Recv → MPI_Finalize.

Q2. Give the compile and run commands for an MPI program on 8 processes.

Solution: mpic++ -O2 prog.cpp -o prog and mpirun -n 8 ./prog.

Q3. Define speed-up and efficiency.

Solution: \(S = T_1/T_N\); \(E = S/N\) (fraction of ideal).

Q4. Which OpenMP clause makes a parallel sum correct?

Solution: reduction(+:sum).

Q5. Shared or distributed memory: which model can use all nodes of a cluster?

Solution: Distributed memory (MPI) — shared memory ends at the node boundary.

Level 2 — Intuitive (predict / explain why)¶

Q1. mpirun -n 4 ./ex01 prints ranks in a different order each run. Why, and is it a bug?

Solution: Not a bug: 4 independent processes race to write to the shared terminal; MPI imposes no global ordering between ranks' outputs.

Q2. With \(f = 0.10\), your boss orders 1000 cores expecting ~1000× speed-up. What do you tell them? Compute the bound.

Solution: \(S_{max} = 1/f = 10\). Even with infinite cores, more than 10× is impossible at fixed problem size; at N=1000, \(S = 1/(0.1+0.9/1000) \approx 9.91\). Money would be better spent reducing the serial fraction (or growing the problem — Gustafson).

Q3. Why does the SAME code give 15× on 16 ranks for a \(512^{3}\) mesh but only 3× for a \(64^{3}\) mesh?

Solution: Surface-to-volume: the small mesh's sub-domains are mostly halo — communication dominates compute. Comm/compute ratio \(\propto P^{1/3}/n\): shrinking n by 8 multiplies the ratio by 8.

Q4. Two threads execute counter++ 1000 times each without synchronisation. Why can the result be less than 2000?

Solution: counter++ is read-modify-write; interleaved threads can both read the same old value and write back the same incremented value — increments get lost. That's the race condition.

Q5. Predict: rank 0 calls MPI_Recv from rank 1, and rank 1 calls MPI_Recv from rank 0, then both Send. What happens?

Solution: Deadlock — both block in Recv; no one ever reaches Send. Fix: one rank sends first, or use MPI_Sendrecv/non-blocking calls.

Level 3 — Hard (exam level)¶

Q1. (8 marks) \(f = 0.05\). Compute Amdahl speed-up and efficiency at N = 16 and N = 256, and the limit. Comment on whether buying 256 cores is sensible.

Solution: \(S(16) = 1/(0.05+0.95/16) = 1/0.109375 = 9.14\), \(E = 9.14/16 = 57\%\). \(S(256) = 1/(0.05+0.95/256) = 1/0.05371 = 18.6\), \(E = 18.6/256 = 7.3\%\). Limit \(S_{max} = 20\). Going 16→256 (16× cores) gains only 2× speed-up at 7% efficiency — not sensible for this fixed problem. বাংলা ইঙ্গিত: efficiency-ই আসল বিচারক: \(E = S/N\) যখন এক অঙ্কে নেমে আসে, তখন core বাড়ানো মানে টাকা পোড়ানো।

Q2. (8 marks) Derive Gustafson's law in 4 lines and evaluate at f = 0.05, N = 64. Why does it disagree with Amdahl's 15.4?

Solution: Normalise the parallel run to time \(f + (1-f) = 1\). A single processor doing the same scaled job needs \(f + N(1-f)\). Hence \(S = f + N(1-f) = N - f(N-1)\). At N=64: \(S = 64 - 0.05\cdot63 = 60.85\). No contradiction: Amdahl fixes the problem size (strong scaling); Gustafson grows it with N (weak scaling) — bigger machines run bigger problems, not the same problem faster. বাংলা ইঙ্গিত: দুই সূত্রের পার্থক্য এক বাক্যে: Amdahl — "একই কাজ কত দ্রুত", Gustafson — "একই সময়ে কত বড় কাজ"। পরীক্ষায় এই বাক্যটাই লিখো।

Q3. (10 marks) Write a complete MPI ring-pass: rank r sends its rank to (r+1)%size, receives from (r-1+size)%size, prints what it got. Avoid deadlock for ANY size.

Solution:

#include <iostream>
#include <mpi.h>
int main(int argc,char**argv){
    int rank,size;
    MPI_Init(&argc,&argv);
    MPI_Comm_rank(MPI_COMM_WORLD,&rank);
    MPI_Comm_size(MPI_COMM_WORLD,&size);
    int right = (rank+1)%size, left = (rank-1+size)%size;
    int sendv = rank, recvv = -1;
    MPI_Sendrecv(&sendv, 1, MPI_INT, right, 0,
                 &recvv, 1, MPI_INT, left,  0,
                 MPI_COMM_WORLD, MPI_STATUS_IGNORE);
    std::cout << "Rank " << rank << " received " << recvv << "\n";
    MPI_Finalize();
}

MPI_Sendrecv performs the exchange atomically — no rank ever blocks in a bare Recv, so the ring cannot deadlock (which plain Send→Recv could, if all Sends block). বাংলা ইঙ্গিত: ring/pairwise exchange দেখলেই MPI_Sendrecv — এক call-এ দুই দিক, deadlock-এর প্রশ্নই নেই; modulo-গণিতে (rank-1+size)%size — মাইনাসের negative এড়াতে +size।

Q4. (10 marks) Find ALL the bugs:

#pragma omp parallel for
for (int i = 1; i <= N; ++i) {
    sum += a[i];
    max = (a[i] > max) ? a[i] : max;
}

Solution: (1) race on sum → needs reduction(+:sum); (2) race on max → needs reduction(max:max); (3) bounds: starts at 1 and includes N — if a has N elements indexed 0..N-1, it skips a[0] and reads OOB at a[N]. Fixed:

#pragma omp parallel for reduction(+:sum) reduction(max:max)
for (int i = 0; i < N; ++i) { sum += a[i]; max = a[i]>max ? a[i] : max; }

বাংলা ইঙ্গিত: "find ALL bugs" মানে কমপক্ষে তিনটা — দুটো race আর একটা off-by-one; একটা পেয়ে থেমে গেলেই নম্বর কাটা।

Q5. (10 marks) A 3D mesh of \(512^3\) cells runs on P = 64 ranks. Compute per-rank cells, per-rank halo cells (1-cell-thick, 6 faces), and the halo-to-interior ratio. What happens to the ratio at P = 512?

Solution: Sub-cube edge: \(512/64^{1/3} = 512/4 = 128\) → \(128^3 = 2{,}097{,}152\) cells per rank. Halo: \(6 \times 128^2 = 98{,}304\) cells. Ratio: \(98{,}304/2{,}097{,}152 \approx 4.7\%\). At P=512: edge \(= 512/8 = 64\), halo ratio \(= 6\cdot64^2/64^3 = 6/64 \approx 9.4\%\) — doubles. General rule: ratio \(= 6/(\text{edge})\), growing as \(P^{1/3}\). বাংলা ইঙ্গিত: হিসাবের কাঠামো মুখস্থ: edge = \(n/P^{1/3}\), halo = \(6\cdot\text{edge}^2\), ratio = \(6/\text{edge}\) — তিন লাইনের template।

Level 4 — Beyond the lecture (transfer + coding)¶

Q1. Implement the 1D heat-equation halo exchange skeleton: each rank owns local_n cells in u[1..local_n] with ghosts u[0] and u[local_n+1]. Write ONE timestep's communication (non-blocking) + update loop.

Solution:

MPI_Request req[4];
int L = rank-1, R = rank+1;
int nreq = 0;
if (rank > 0) {
    MPI_Irecv(&u[0],         1, MPI_DOUBLE, L, 0, MPI_COMM_WORLD, &req[nreq++]);
    MPI_Isend(&u[1],         1, MPI_DOUBLE, L, 1, MPI_COMM_WORLD, &req[nreq++]);
}
if (rank < size-1) {
    MPI_Irecv(&u[local_n+1], 1, MPI_DOUBLE, R, 1, MPI_COMM_WORLD, &req[nreq++]);
    MPI_Isend(&u[local_n],   1, MPI_DOUBLE, R, 0, MPI_COMM_WORLD, &req[nreq++]);
}
MPI_Waitall(nreq, req, MPI_STATUSES_IGNORE);
for (int i = 1; i <= local_n; ++i)
    unew[i] = u[i] + alpha*(u[i-1] - 2*u[i] + u[i+1]);

Tags pair each send with the matching recv direction; boundary ranks skip the missing neighbour. বাংলা ইঙ্গিত: ghost-cell index-ই আসল পরীক্ষা: বাঁ ghost u[0] ← বাঁ প্রতিবেশীর শেষ আসল cell; ডান ghost u[local_n+1] ← ডান প্রতিবেশীর প্রথম আসল cell — ছবি এঁকে মেলাও।

Q2. Your strong-scaling study gives: 1 rank 1000 s, 4 ranks 270 s, 16 ranks 90 s, 64 ranks 60 s. Compute the speed-ups, fit the serial fraction f from the N=64 point using Amdahl, and predict T at N=256.

Solution: \(S_4 = 3.70\), \(S_{16} = 11.1\), \(S_{64} = 16.7\). From Amdahl: \(16.7 = 1/(f + (1-f)/64)\) ⇒ \(f + (1-f)/64 = 0.06\) ⇒ \(f(1-1/64) = 0.06-0.015625\) ⇒ \(f = 0.0444/0.9844 \approx 0.045\). Predict: \(S_{256} = 1/(0.045+0.955/256) = 1/(0.04873) \approx 20.5\) → \(T_{256} \approx 1000/20.5 \approx 49\) s. Barely faster than 64 ranks for 4× the cores. বাংলা ইঙ্গিত: মাপা ডেটা থেকে f বের করা = Amdahl উল্টো করে সমাধান — এই reverse-প্রশ্নটাই "out of topic" মার্কা, কিন্তু বীজগণিত সোজা।

Q3. (OpenMP + C++ classes, Ch 11 transfer) Why is this parallel loop subtly broken even though each thread writes a different element?

std::vector<double> v;
#pragma omp parallel for
for (int i = 0; i < N; ++i) v.push_back(f(i));

Solution: push_back mutates shared state (size, capacity, possible reallocation moving the whole buffer) — it is NOT writing "a different element", it's appending to one shared container → data race and corruption. Fix: pre-size, then index:

std::vector<double> v(N);
#pragma omp parallel for
for (int i = 0; i < N; ++i) v[i] = f(i);

বাংলা ইঙ্গিত: "আলাদা element-এ লিখছি" ভ্রমটাই ফাঁদ — push_back মানে shared metadata-তে লেখা; থাম্ব-রুল: parallel লুপে container-এর আকার বদলানো নিষেধ।

Q4. (Bash + MPI, Ch 8 transfer) Write a bash loop that runs the π program on 1, 2, 4, 8, 16 ranks, extracts the runtime printed as time: <sec>, and writes ranks time speedup rows into scaling.dat (speed-up vs the 1-rank time).

Solution:

#!/bin/bash
set -euo pipefail
: > scaling.dat
t1=""
for n in 1 2 4 8 16; do
    t=$(mpirun -n $n ./mc | awk '/^time:/ {print $2}')
    if [ -z "$t1" ]; then t1=$t; fi
    s=$(awk -v a="$t1" -v b="$t" 'BEGIN{printf "%.2f", a/b}')
    echo "$n $t $s" >> scaling.dat
done
column -t scaling.dat

awk does the float division bash cannot; the first iteration pins \(t_1\). বাংলা ইঙ্গিত: bash float পারে না — speed-up-এর ভাগটা awk-এ; আর : > file মানে ফাইল খালি করে শুরু। এক প্রশ্নে তিন chapter — এটাই শেষ বসের চেহারা।

End of Chapter 14.

Chapter 14: Parallelization — Shared & Distributed Memory, OpenMP, MPI¶

1. Chapter Overview¶

2. Basics from Zero¶

3. Hard English Made Easy¶

4. Deep Theory Explanation¶

4.1 Levels of parallelism (slide 4)¶

4.2 Memory models (slide 5)¶

4.3 Amdahl's law — derived from first principles¶

4.4 Strong vs weak scaling — and Gustafson's law¶

4.5 Domain decomposition (CFD specific)¶

4.6 OpenMP (a one-page tour)¶

4.7 MPI six-call backbone¶

4.8 Compile & run MPI¶

4.9 Sequential vs causal consistency (slide 22)¶

4.10 The two remaining must-know diagrams¶

5. Command / Syntax / Code Breakdown¶

#include <mpi.h> / #include <omp.h>¶

MPI_Init(&argc, &argv) / MPI_Finalize()¶

MPI_Comm_rank(MPI_COMM_WORLD, &rank)¶

MPI_Comm_size(MPI_COMM_WORLD, &size)¶

MPI_Send(buf, count, type, dest, tag, comm)¶

MPI_Recv(buf, count, type, source, tag, comm, &status)¶

MPI_Bcast(buf, count, type, root, comm)¶

MPI_Reduce(sendbuf, recvbuf, count, type, op, root, comm)¶

mpirun -n N ./prog¶

#pragma omp parallel for [reduction(+:sum)] [schedule(dynamic,k)]¶

#pragma omp critical { … }¶

#pragma omp atomic¶

#pragma omp barrier¶

OMP_NUM_THREADS=8 ./prog¶

6. Mandatory Practical Examples¶

Example 6.1 — MPI hello (example01.cpp)¶

Example 6.2 — Loop split across ranks (example02.cpp)¶

Example 6.3 — Send/Recv int (example03_int.cpp)¶

Example 6.4 — Send/Recv string (example03_Cstring.cpp)¶

Example 6.5 — OpenMP parallel sum¶

Example 6.6 — Race condition fix¶

Example 6.7 — Hybrid MPI + OpenMP¶

7. Real HPC/CFD Workflow¶

8. Exercises and Solutions¶

Self-made exercise A — MPI sum¶

Self-made exercise B — OpenMP threaded sort¶

9. Written Exam Focus¶

9.1 Short Answers¶

9.2 Medium Answers¶

9.3 Long Answer (12 marks)¶

9.4 Output Prediction¶

9.5 Comparison¶

9.6 Templates¶

9.7 Marking Scheme — "MPI Send/Recv" (5 marks)¶

10. Very Hard Questions¶

11. Debugging and Mistake Analysis¶

12. Mini Project for Mastery¶

13. Final Chapter Cheat Sheet¶

14. Mock Exam — Four Levels¶

Level 1 — Basic (definitions & syntax)¶

Level 2 — Intuitive (predict / explain why)¶

Level 3 — Hard (exam level)¶

Level 4 — Beyond the lecture (transfer + coding)¶

`#include <mpi.h>` / `#include <omp.h>`¶

`MPI_Init(&argc, &argv)` / `MPI_Finalize()`¶

`MPI_Comm_rank(MPI_COMM_WORLD, &rank)`¶

`MPI_Comm_size(MPI_COMM_WORLD, &size)`¶

`MPI_Send(buf, count, type, dest, tag, comm)`¶

`MPI_Recv(buf, count, type, source, tag, comm, &status)`¶

`MPI_Bcast(buf, count, type, root, comm)`¶

`MPI_Reduce(sendbuf, recvbuf, count, type, op, root, comm)`¶

`mpirun -n N ./prog`¶

`#pragma omp parallel for [reduction(+:sum)] [schedule(dynamic,k)]`¶

`#pragma omp critical { … }`¶

`#pragma omp atomic`¶

`#pragma omp barrier`¶

`OMP_NUM_THREADS=8 ./prog`¶