Parallel Sample Sort Foundations and Shared-Memory Solvers

Sorting is a fundamental building block across computational science, databases, and scientific computing. While algorithms like quicksort and mergesort work exceptionally well in sequential computing, parallelizing them efficiently across multi-core CPUs and distributed clusters presents non-trivial challenges in load balancing and communication.

In this lesson, we explore Sample Sort—a parallel generalization of bucket sort designed to handle unknown and non-uniform data distributions—and implement scalable shared-memory versions using OpenMP and Pthreads.

---

7.15 From Bucket Sort to Sample Sort

1. Bucket Sort (Uniform Distribution Assumption)

In classical Bucket Sort, we are given $n$ keys with a known uniform distribution over an interval $[c, d)$. We partition $[c, d)$ into $b$ equal-width buckets of width $h = (d - c) / b$:

$[c, c + h), \quad [c + h, c + 2h), \quad \dots, \quad [c + (b - 1)h, d)$

Keys are routed to buckets using direct arithmetic ($i = \lfloor (x - c) / h \rfloor$), sorted locally, and concatenated. If keys are uniformly distributed, each bucket receives approximately $n / b$ elements, yielding an average $O(n)$ time complexity.

However, if the distribution is non-uniform (e.g., highly skewed or clustered), one bucket may receive the vast majority of keys while others remain empty, destroying parallel load balance and degrading performance to $O(n^2)$.

---

2. Sample Sort (Unknown Distributions)

Sample Sort eliminates the uniform distribution requirement by taking a statistical sample of $s$ keys ($b \le s \le n$) from the original dataset to estimate its actual distribution.

flowchart LR
  subgraph SampleSortFlow["Figure 7.10: Five Phases of Sample Sort"]
      direction LR
      P1["1. Sample Selection<br/>Choose s keys from list"] --> P2["2. Sort Sample and Splitters<br/>Extract b - 1 splitters e_i"]
      P2 --> P3["3. Build Mapping<br/>Count items per bucket<br/>Prefix sums"]
      P3 --> P4["4. Route Keys to Buckets<br/>Push or Pull data"]
      P4 --> P5["5. Local Sort and Concatenate<br/>Final sorted output"]
  end

Splitter Selection Workflow:

1. Sort the $s$ sample keys: $a_0 \le a_1 \le a_2 \le \dots \le a_{s-1}$
2. Select $b - 1$ splitters spaced evenly by intervals of $s / b$ keys: $e_1, e_2, \dots, e_{b-1}$
3. Define the $b$ bucket ranges: $[-\infty, e_1), \quad [e_1, e_2), \quad \dots, \quad [e_{b-1}, \infty)$

Because the splitters reflect the empirical density of the keys, each bucket receives approximately equal numbers of elements, guaranteeing balanced parallel workloads regardless of the underlying distribution!

---

Concrete Splitter Example

Suppose we want $b = 4$ buckets from a sorted sample of $s = 12$ keys:

$\{2, 10, 23, 25, 29, 38, 40, 49, 53, 60, 62, 67\}$

The sample is partitioned into 4 groups of $s / b = 3$ elements each:

• Group 0: $\{2, 10, 23\}$
• Group 1: $\{25, 29, 38\}$
• Group 2: $\{40, 49, 53\}$
• Group 3: $\{60, 62, 67\}$

The splitters separating these groups are:

$e_1 = \dfrac{23 + 25}{2} = 24, \quad e_2 = \dfrac{38 + 40}{2} = 39, \quad e_3 = \dfrac{53 + 60}{2} = 57$

The 4 resulting bucket intervals are:

$[-\infty, 24), \quad [24, 39), \quad [39, 57), \quad [57, \infty)$

---

7.16 Mapping Functions and Prefix Sums

To distribute keys into buckets, we can choose between two fundamental strategies:

Strategy 1: Dynamic Reallocation (realloc)

Each bucket begins with capacity $n / b$. When a bucket overflows, its capacity is doubled using realloc(). While simple to code, repeated calls to realloc() cause memory fragmentation, thread lock contention, and expensive data copies.

Strategy 2: Pre-Allocated Mapping via Exclusive Prefix Sums

Instead of dynamic reallocation, we construct exact bucket offsets in advance using exclusive prefix sums, allowing all keys to be placed into a single pre-allocated contiguous array without reallocation!

---

Mathematical Anatomy of Prefix Sums

Given an array of $n$ numbers $\{x_0, x_1, \dots, x_{n-1}\}$:

• Inclusive Prefix Sum: Each output element includes the corresponding input element: $\text{incl}[i] = \sum_{k=0}^{i} x_k \implies \{x_0, x_0 + x_1, \dots, \sum_{k=0}^{n-1} x_k\}$
• Exclusive Prefix Sum: Each output element sums only the strictly preceding elements, beginning with $0$: $\text{excl}[0] = 0, \quad \text{excl}[i] = \sum_{k=0}^{i-1} x_k \quad (i \ge 1)$

For input array $\{5, 1, 6, 2, 4\}$:

• Inclusive: $\{5, 6, 12, 14, 18\}$
• Exclusive: $\{0, 5, 6, 12, 14\}$

---

Step-by-Step Mapping Matrix Trace

Consider a dataset with $n = 24$ elements divided into $b = 4$ sublists of 6 keys each, with splitters $\{-\infty, 45, 75, 91, \infty\}$:

Sublist 0: {15, 35, 77, 83, 86, 93}
Sublist 1: {21, 27, 49, 62, 86, 92}
Sublist 2: {26, 26, 40, 59, 63, 90}
Sublist 3: {11, 29, 36, 67, 68, 72}

Step 1: Count Matrix (from_to_cts[sublist][bucket])

We count how many elements from each sublist fall into each bucket:

$$\text{from\_to\_cts} = \begin{bmatrix} 2 & 0 & 3 & 1 \\ 2 & 2 & 1 & 1 \\ 3 & 2 & 1 & 0 \\ 3 & 3 & 0 & 0 \end{bmatrix}$$

Step 2: Bucket Offsets within Sublists (bkt_starts_in_slists)

Row-wise exclusive prefix sums of from_to_cts determine the starting index of each bucket within each sublist:

$$\text{bkt\_starts\_in\_slists} = \begin{bmatrix} 0 & 2 & 2 & 5 \\ 0 & 2 & 4 & 5 \\ 0 & 3 & 5 & 6 \\ 0 & 3 & 6 & 6 \end{bmatrix}$$

Step 3: Sublist Offsets within Buckets (slist_starts_in_bkts)

Column-wise exclusive prefix sums of from_to_cts determine where each sublist’s elements begin inside each bucket:

$$\text{slist\_starts\_in\_bkts} = \begin{bmatrix} 0 & 0 & 0 & 0 \\ 2 & 0 & 3 & 1 \\ 4 & 2 & 4 & 2 \\ 7 & 4 & 5 & 2 \end{bmatrix}, \quad \text{Bucket Sizes} = \begin{bmatrix} 10 & 7 & 5 & 2 \end{bmatrix}$$

Step 4: Global Bucket Starting Positions (bkt_starts)

Exclusive prefix sums of the bucket sizes vector $\{10, 7, 5, 2\}$ determine exact starting offsets in the output array:

$$\text{bkt\_starts} = \begin{bmatrix} 0 & 10 & 17 & 22 \end{bmatrix}$$

With these tables, every thread moves its keys directly into non-overlapping memory locations with zero lock contention and zero memory reallocations!

---

7.17 OpenMP Implementations & Benchmarks

Version 1: Master Sampling with Thread-Local Reallocation

In the first OpenMP implementation:

1. Master thread randomly selects and sorts sample keys, computing splitters.
2. Inside #pragma omp parallel, each thread “pulls” all matching keys from the entire input array into its private bucket, calling realloc() as needed.
3. Each thread sorts its bucket with qsort().
4. #pragma omp barrier synchronizes all threads.
5. Threads write their sorted buckets back to the global array.

---

Table 7.10: Version 1 OpenMP Run-times ($n = 2^{24} = 16,777,216$ integers)

Threads	Run-time	Speedup vs. 1 Thread
1	$3.54\,\text{s}$	$1.00\times$
2	$2.29\,\text{s}$	$1.55\times$
4	$1.69\,\text{s}$	$2.09\times$
8	$1.34\,\text{s}$	$2.64\times$
16	$1.40\,\text{s}$	$2.53\times$
32	$1.36\,\text{s}$	$2.60\times$
64	$1.33\,\text{s}$	$2.66\times$
Sequential qsort	$2.59\,\text{s}$	—

Bottleneck Analysis: Performance plateaus beyond 8 threads because every thread scans the entire 16-million-element input list, duplicating memory reads $p$ times and saturating memory bus bandwidth.

---

Version 2: Fully Parallel Deterministic Sampling & Prefix Mapping

To achieve linear scalability, Version 2 eliminates duplicated memory scans:

1. Parallel In-Place Sort: Each thread sorts only its own $n / p$ sublist.
2. Parallel Deterministic Sampling: Each thread takes $s / p$ evenly spaced keys from its sorted sublist.
3. Thread 0 sorts the tiny sample ($s \ll n$) and extracts splitters.
4. Each thread fills its row of from_to_cts using binary search against the splitters ($O(\log b)$ instead of $O(n)$ linear scans).
5. Fast prefix sums allocate contiguous bucket slices.
6. Threads merge and push keys directly to destination buckets without reallocation.

---

Table 7.11: Version 2 OpenMP Run-times ($n = 2^{24} = 16,777,216$ integers)

Threads	Version 1 Run-time	Version 2 Run-time	Speedup vs. Serial qsort ($2.59\,\text{s}$)
1	$3.54\,\text{s}$	$2.65\,\text{s}$	$0.98\times$
2	$2.29\,\text{s}$	$1.33\,\text{s}$	$1.95\times$
4	$1.69\,\text{s}$	$1.01\,\text{s}$	$2.56\times$
8	$1.34\,\text{s}$	$0.35\,\text{s}$	$7.40\times$
16	$1.40\,\text{s}$	$0.21\,\text{s}$	$12.33\times$
32	$1.36\,\text{s}$	$0.15\,\text{s}$	$17.27\times$
64	$1.33\,\text{s}$	$0.13\,\text{s}$	$19.92\times$

By eliminating redundant scans and memory reallocations, Version 2 scales smoothly to 64 cores, sorting 16.7 million integers in just $0.13$ seconds—a $20\times$ speedup over optimized sequential qsort!

---

7.18 Parallelizing Sample Sort with Pthreads

The Pthreads implementation mirrors OpenMP Version 2, replacing pragma directives with explicit thread creation (pthread_create) and synchronization barriers:

Table 7.12 & 7.13: Pthreads Sample Sort Performance ($n = 2^{24}$)

Threads	Version 1 (Pthreads)	Version 2 (Pthreads)	Version 2 (OpenMP)
1	$3.55\,\text{s}$	$2.65\,\text{s}$	$2.65\,\text{s}$
2	$2.31\,\text{s}$	$1.36\,\text{s}$	$1.33\,\text{s}$
4	$1.71\,\text{s}$	$1.08\,\text{s}$	$1.01\,\text{s}$
8	$1.48\,\text{s}$	$0.35\,\text{s}$	$0.35\,\text{s}$
16	$1.45\,\text{s}$	$0.20\,\text{s}$	$0.21\,\text{s}$
32	$1.28\,\text{s}$	$0.14\,\text{s}$	$0.15\,\text{s}$
64	$1.24\,\text{s}$	$0.12\,\text{s}$	$0.13\,\text{s}$

The benchmark results confirm that Pthreads and OpenMP deliver virtually identical performance when using equivalent algorithmic designs. OpenMP offers significantly cleaner and more maintainable code through compiler directives, while Pthreads grants low-level control when non-standard synchronization patterns are required.