Session 17¶

Topic 2: The Go language¶

Today¶

Case studies with parallel algorithms¶

parallel matrix-matrix multiplication
parallel histogram

Parallel matrix multiplication in Go¶

The problem¶

compute \(C = A B\) for \(n \times n\) (square) matrices:

\[ C_{i,j} = \sum_{k=0}^{n-1} A_{i,k} B_{k,j} \]

Strategy:

preallocate \(C\), split rows of \(C\) among k goroutines
each goroutine owns a non-overlapping row block
→ no synchronization needed on writes

Sequential baseline¶

func SeqMatMul(A, B [][]float64) [][]float64 {
  n := len(A)
  C := make([][]float64, n)
  for i := range C {
    C[i] = make([]float64, n)
  }
  // loop order: i-k-j
  for i := 0; i < n; i++ {
    for k := 0; k < n; k++ {
      for j := 0; j < n; j++ {
        C[i][j] += A[i][k] * B[k][j]
      }
    }
  }
  return C
}

Loop order: i-k-j keeps B[k][:] in cache across the inner loop → maximises cache reuse on B

SeqMatMul takes two square matrices A and B as slices of row slices and returns their product C. It first allocates C as n rows, each an n-element float64 slice initialised to zero. The three nested loops then compute the standard matrix product \(C_{i,j} = \sum_k A_{i,k} B_{k,j}\), but with a deliberate i-k-j loop order rather than the more intuitive i-j-k. The reason is cache efficiency: the outermost k loop fixes a single value A[i][k], which is then reused across the entire inner j sweep without reloading; more importantly, the inner j loop strides through B[k][j] and C[i][j] sequentially in memory (since Go stores each row contiguously), so both accesses are cache-friendly. The naive i-j-k order, by contrast, would stride through B column-by-column in the inner loop, jumping between rows on every step and causing frequent cache misses for large n. The accumulated products land directly in C[i][j] via +=, and the completed matrix is returned.

Parallel implementation¶

Distributing the work¶

func ParallelMatMul(A, B [][]float64, numWorkers int) [][]float64 {
  n := len(A)
  C := allocMatrix(n)

  var wg sync.WaitGroup
  rowsPerWorker := (n + numWorkers - 1) / numWorkers

  for w := 0; w < numWorkers; w++ {
    lo := w * rowsPerWorker
    hi := min(lo+rowsPerWorker, n)
    if lo >= n { break }

    // ... launch goroutines on a subset of rows (shown next)
  }
  wg.Wait()
  return C
}

Launching the goroutines¶

  for w := 0; w < numWorkers; w++ {
    // ...

    wg.Add(1)
    go func(lo, hi int) {
      defer wg.Done()
      for i := lo; i < hi; i++ {
        for k := 0; k < n; k++ {
          for j := 0; j < n; j++ {
            C[i][j] += A[i][k] * B[k][j]
          }
        }
      }
    }(lo, hi)
  }

Full parallel implementation:

func ParallelMatMul(A, B [][]float64, numWorkers int) [][]float64 {
  n := len(A)
  C := allocMatrix(n)

  var wg sync.WaitGroup
  rowsPerWorker := (n + numWorkers - 1) / numWorkers

  for w := 0; w < numWorkers; w++ {
    lo := w * rowsPerWorker
    hi := min(lo+rowsPerWorker, n)
    if lo >= n { break }

    wg.Add(1)
    go func(lo, hi int) {
      defer wg.Done()
      for i := lo; i < hi; i++ {
        for k := 0; k < n; k++ {
          for j := 0; j < n; j++ {
            C[i][j] += A[i][k] * B[k][j]
          }
        }
      }
    }(lo, hi)
  }
  wg.Wait()
  return C
}

ParallelMatMul computes the same matrix product as the sequential version but distributes the work across numWorkers goroutines, each responsible for a contiguous horizontal slice of rows in the result matrix C. The output matrix is pre-allocated once with allocMatrix before any goroutine is launched, which is safe because each goroutine writes to a strictly non-overlapping row range — no two goroutines ever touch the same element of C, so no mutex or channel is needed for the writes. The row range for each worker is calculated with ceiling division (n + numWorkers - 1) / numWorkers, which ensures that rows are spread as evenly as possible even when n is not exactly divisible by numWorkers; the min clamp on hi and the lo >= n guard prevent the last worker from running off the end of the matrix when numWorkers exceeds n. Each goroutine receives its lo and hi bounds as explicit function arguments rather than closing over the loop variable w, avoiding the classic loop-variable capture bug. Inside each goroutine the loop order is again i-k-j for the same cache-locality reasons. A sync.WaitGroup is used to coordinate completion: wg.Add(1) is called before each go statement, each goroutine calls defer wg.Done() on exit, and wg.Wait() blocks the caller until every worker has finished before the completed matrix is returned.

Correctness check¶

func matricesEqual(A, B [][]float64, tol float64) bool {
  for i := range A {
    for j := range A[i] {
      if math.Abs(A[i][j]-B[i][j]) > tol {
        return false
      }
    }
  }
  return true
}

// test:
for _, n := range []int{100, 500, 1000} {
  A, B := randMatrix(n), randMatrix(n)
  Cgot  := ParallelMatMul(A, B, runtime.NumCPU())
  Cwant := SeqMatMul(A, B)
  if !matricesEqual(Cgot, Cwant, 1e-9) {
    log.Fatalf("FAIL n=%d", n)
  }
}

Benchmarking & plotting¶

for _, numWorkers := range []int{1, 2, 4, 8, 16, 32} {
  start := time.Now()
  ParallelMatMul(A, B, numWorkers)
  fmt.Printf("numWorkers=%2d  %.3fs\n",
             numWorkers, time.Since(start).Seconds())
}

plot wall-clock time vs. numWorkers for each \(n\)

Expected behaviour¶

`numWorkers`	Behaviour
1 → `NumCPU`	time drops
`numWorkers` = `NumCPU`	near-optimal
`numWorkers` > `NumCPU`	flat or slightly worse

Hands-on tasks¶

Run correctness tests for \(n = 100, 500, 1000\).
Benchmark numWorkers = 1, 2, 4, 8, 16, 32 and print wall-clock time.
Add GOMAXPROCS(1) — what changes and why?
(stretch) Transpose \(B\) before multiplying. Does it help? Why?

What happens when numWorkers > NumCPU?¶

the OS scheduler time-slices goroutines onto physical cores
no additional parallelism — each core still does the same total work
overhead from goroutine creation, scheduling, and context switching increases
for CPU-bound work: speedup plateaus, then regresses

Rule of thumb: set numWorkers = runtime.NumCPU() for pure compute tasks

Amdahl’s law¶

\[ S_p = \frac{1}{\alpha + \frac{1 - \alpha}{p}} \]

for matrix multiply, the serial fraction \(\alpha \approx 0\) (only allocation & loop setup)

\(\alpha\)	Theoretical max speedup
\(0\)%	\(p\times\) (ideal)
\(1\)%	\(50\times\)
\(5\)%	\(20\times\)

In practice:

memory bandwidth and cache contention hit first
large \(n\) → rows of \(B\) exceed cache → bandwidth-bound, not compute-bound

Key takeaways¶

row partitioning: contention-free parallel writes
no channels needed, only sync.WaitGroup
optimal goroutine count ≈ runtime.NumCPU(); more goroutines add overhead
Amdahl’s law predicts near-linear speedup — memory bandwidth is the real ceiling
cache-friendly loop order (i-k-j) matters as much as parallelism

Parallel Histogram in Go¶

The problem¶

Given \(n\) observations \(x_1, \dots, x_n \in \mathbb{R}\), compute a histogram with \(B\) bins over \([\min, \max]\).

*Figure: An example
histogram.*

bin width:

\[ w = \frac{\max - \min}{B} \]

bin index for observation \(x_i\):

\[ b_i = \left\lfloor \frac{x_i - \min}{w} \right\rfloor, \quad b_i \in \{0, \dots, B-1\} \]

Goal: parallelize using \(p\) goroutines, each processing a chunk of size \(\lceil n/p \rceil\)

Algorithm overview¶

data  ──┬──► goroutine 1 ──► hist_1 ──┐
        ├──► goroutine 2 ──► hist_2 ──┤──► merge ──► global hist
        ├──► goroutine 3 ──► hist_3 ──┤
        └──► goroutine p ──► hist_4 ──┘

Distribution of data¶

each goroutine \(g\) owns this slice of the data:

\[ \text{chunk}_g = \text{data}\!\Bigl[\, g \cdot \lceil n/p \rceil \;:\; \min\bigl((g{+}1)\lceil n/p \rceil,\, n\bigr) \Bigr] \]

Local histograms and the global histogram¶

local histogram update (no synchronization needed within a goroutine):

\[ H_g[b] \mathrel{+}= \mathbf{1}\!\left[b_i = b\right] \]

merge (sequential, after all goroutines finish):

\[ H[b] = \sum_{g=0}^{p-1} H_g[b], \quad b = 0, \dots, B-1 \]

Core implementation¶

func ParallelHistogram(data []float64, bins int, numWorkers int) []int {
  n := len(data)
  lo, hi := slices.Min(data), slices.Max(data)
  w := (hi - lo) / float64(bins)

  binOf := func(x float64) int {
    b := int((x - lo) / w)
    if b >= bins { b = bins - 1 } // clamp max value
    return b
  }

  chunkSize := (n + numWorkers - 1) / numWorkers
  ch := make(chan []int, numWorkers)

  // compute local histograms ...
  // compute the global histogram ...

Compute local histograms¶

  // compute local histograms
  for g := 0; g < numWorkers; g++ {
    start, end := g*chunkSize, min((g+1)*chunkSize, n)
    go func(chunk []float64) {
      local := make([]int, bins)
      for _, x := range chunk {
        local[binOf(x)]++
      }
      ch <- local
    }(data[start:end])
  }

Compute the global histogram¶

  // compute the global histogram
  global := make([]int, bins)
  for range numWorkers {
    for b, c := range <-ch {
      global[b] += c
    }
  }
  return global

Full parallel implementation:

func ParallelHistogram(data []float64, bins int, numWorkers int) []int {
  n := len(data)
  lo, hi := slices.Min(data), slices.Max(data)
  w := (hi - lo) / float64(bins)

  binOf := func(x float64) int {
    b := int((x - lo) / w)
    if b >= bins { b = bins - 1 } // clamp max value
    return b
  }

  chunkSize := (n + numWorkers - 1) / numWorkers
  ch := make(chan []int, numWorkers)

  // compute local histograms
  for g := 0; g < numWorkers; g++ {
    start, end := g*chunkSize, min((g+1)*chunkSize, n)
    go func(chunk []float64) {
      local := make([]int, bins)
      for _, x := range chunk {
        local[binOf(x)]++
      }
      ch <- local
    }(data[start:end])
  }

  // compute the global histogram
  global := make([]int, bins)
  for range numWorkers {
    for b, c := range <-ch {
      global[b] += c
    }
  }
  return global
}

Strategy A: channel (message passing)¶

Each goroutine sends its complete local histogram over a buffered channel.

ch := make(chan []int, numWorkers) // buffered: goroutines never block

go func(chunk []float64) {
  local := make([]int, bins)
  for _, x := range chunk { local[binOf(x)]++ }
  ch <- local                      // send full slice once
}(data[start:end])

// main goroutine merges
for range numWorkers {
  for b, c := range <-ch { global[b] += c }
}

Properties:

zero contention during computation
one allocation per goroutine (bins ints)
communication cost: \(O(p \cdot B)\) total values sent through channel

Strategy B: shared histogram + mutex¶

All goroutines write directly to a single shared histogram under a mutex.

var mu sync.Mutex
global := make([]int, bins)

go func(chunk []float64) {
  local := make([]int, bins)
  for _, x := range chunk { local[binOf(x)]++ }

  mu.Lock()
  for b, c := range local { global[b] += c }
  mu.Unlock()
}(data[start:end])

Properties:

lock held only during the \(O(B)\) merge, not during computation
contention is low: at most \(p\) goroutines competing for one short critical section
no channel overhead; merge is in-place

Anti-pattern to avoid: locking inside the per-element loop → \(O(n)\) lock acquisitions can result in poor performance

Benchmarking

func BenchmarkChannel(b *testing.B) {
  data := generateData(1_000_000)
  b.ResetTimer()
  for b.Loop() {
    ParallelHistogramChannel(data, 100, runtime.NumCPU())
  }
}

func BenchmarkMutex(b *testing.B) {
  data := generateData(1_000_000)
  b.ResetTimer()
  for b.Loop() {
    ParallelHistogramMutex(data, 100, runtime.NumCPU())
  }
}

go test ./histogram/ -bench=. -benchmem -count=5 | tee bench.txt
benchstat bench.txt

key metrics:

ns/op
allocs/op
speedup \(S_p = T(1) / T(p)\)
efficiency \(E_p = S_p/p\)

Expected results & complexity¶

	Channel	Mutex
sync overhead	\(O(p \cdot B)\) send	\(O(p \cdot B)\) lock+write
memory	\(p+1\) histograms	\(p+1\) histograms
contention	none (goroutines)	low (\(p\) brief locks)
typical winner	small \(B\), many goroutines	large \(B\) (channel copy cost rises)

Amdahl’s law: the fraction \(\alpha\) that is sequential includes min/max & merge

\[ S_p \leq \frac{1}{\alpha + \frac{1-\alpha}{p}} \]

→ don’t expect ideal speedup

Visualization¶

func PrintHistogram(hist []int, lo, hi float64, barWidth int) {
    bins := len(hist)
    w := (hi - lo) / float64(bins)
    maxCount := slices.Max(hist)
    fmt.Println()
    for i, c := range hist {
        label := fmt.Sprintf("[%5.1f,%5.1f)", lo+float64(i)*w, lo+float64(i+1)*w)
        filled := 0
        if maxCount > 0 {
            filled = barWidth * c / maxCount
        }
        bar := strings.Repeat("█", filled) + strings.Repeat("░", barWidth-filled)
        fmt.Printf("  %s │%s│ %d\n", label, bar, c)
    }
    fmt.Println()
}

Hands-on tasks¶

Correctness: Verify sum(hist) == len(data) and compare parallel vs. serial output.
Sweep \(p\): Benchmark for \(p = 1, 2, 4, 8, 16\); print/plot \(S_p\) and \(E_p\).
Sweep \(B\): How does crossover point between channel and mutex shift with bin count?
Eliminate local buffer in Mutex: Lock per element; measure the penalty.
(stretch) sync/atomic: Replace mutex with atomic.AddInt64; compare.

Key takeaways¶

partition → compute locally → merge is the universal pattern for parallel reductions
both strategies avoid contention during computation — synchronization is \(O(B)\), not \(O(n)\)
channel vs. mutex is a tradeoff in communication cost, not correctness; profile before choosing
speedup is bounded by Amdahl’s Law — the sequential merge sets the ceiling

Looking ahead¶

practice your skills by reviewing the previous sessions
practice concurrent implementions of math problems in Go

Prepare for the next session 🌱¶

sessions/18-todo