2. Given two vectors in x, y ∈ Rn , we want to compute their outer product A = x ⊗ y ∈ Rn×n , where Aij = xi yj .
Give a work-depth algorithm for this problem. Derive
its work, depth, and parallelism as a function of n.
1 def outer (x , y , n ) :
2 A = np . zeros (( n , n ) )
3 for i in range ( n ) :
4 for j in range ( n ) :
5 A [i , j ] = x [ i ] * y [ j ]
6 return A
Line 2: we assume memory allocation costs O(1) work. Line 5: W = O n2 , D = O(1).P = O n2
3. State a non-recursive PRAM algorithm for the reduction algorithm for arbitrary n (problem size) and p (number of processors). Notice that, n and p don’t need to
be powers of two and n ≠ p ). Derive its work and time as a function of n and p. What is the expected speedup and efficiency? Is your algorithm work efficient?
STEP 1 : T1 (n, p) = O(n/p) STEP 2: For every value of k we do O(1) work. T2 (n, p) = O(log p) Thus T (n, p) = O(n/p)+O(log p). The work for the first step is pn/p = n.
Ts (n) p
The work for the second step is O(log p). Thus the total work is O(n) + O(log p), thus the algorithm is work optimal. Speedup: S(n, p) = T (n,p) ≈ O(1)+O(p/n log p)
S 1
Efficiency: E(n, p) = p ≈ O(1)+O(p/n log p)
4. Consider the square matrix-vector multiplication problem y = Ax on a PRAM machine, where A ∈ Rn×n . We have already discussed a PRAM algorithm that uses
row-wise partitioning of A and y. This type of partitioning is called ”one-dimensional”, because partitions across one dimension of the matrix.
Here, we will consider a ”two-dimensional” partitioning. Given a PRAM machine with p2 cores, we partition A to p row blocks and p column blocks and x, y in
p blocks. State a PRAM algorithm that exploits such two-dimensional partitioning. Then, derive its complexity (T (n, p), work W (n, p)) and speedup. Assuming
Wsequential (n) = O n2 , is your algorithm work efficient?
1 def pram_matvec (A ,x ,n ,p , tid ) :
2 # A = input matrix ( n x n ) real , x = input vector ( n )
3 # n = length () , p = number of threads , tid = thread number id
4
5 # STEP 1: figure out partition of matrix A
6 q = sqrt ( p )
7 r = n // q
8 pi = tid // q
9 pj = tid % q
10
11 # Step 2: matrix vector multi plicatio n for each thread
12 xloc = x [ pi * r : ( pi +1) * r ] # global concurrent read
13 Aloc = A [ pi * r : ( pi +1) *r , pj * r :( pj +1) * r ] # global exclusive read
14 y = Aloc . dot ( xloc )
15
16 # this is a generalized reduce in which we have only three threads in the row pi , participating
17 # also , yloc is a vector , but the logic is the same as pram reduce
18 y = p ra m_ ge n _r ed u ce ( yloc , r , pj , q )
19
20 if pj == 0:
21 y [ pi * r :( pi +1) * r ] = y
22
23 return y
√ √
The time of the sequential matvec is O n2 /p and the work is O n2 The time of the reduction for y is O n2 /p + O(n/ p log p) since we have n/ p elements per
√ 2
thread, and every reduction operation costs n/ p. The work is O n .
5. Let A ∈ Rn×n be a lower triangular matrix (i.e., Aij = 0 if j > i ) such that Aii ̸= 0, for 1 ≤ i ≤ n, and let b ∈ Rn . Consider the forward-substitution algorithm for
solving Ax = b for x :
1
x1 = b1 ,
A11
i−1
1 X
xi = bi − Aij xj , i = 2, . . . , n.
Aii j=1
Determine the work and depth of this algorithm (using the work depth model). Then give a PRAM version derive its time and work complexity as a function of the
input size n and the number of processors p. (Optional:) Suggest ways to improve the complexity of this algorithm (in the work-depth model).
1 def forwardsubst (A ,x , n ) :
2 x [0] = b [0]/ A [0 ,0]
3 z = zeros ( n )
4 for i in range (1: n ) : # sequential loop over rows
5 for j in range ( i ) :
6 z [ i ] = A [i , j ] x [ j ]
7 x [ i ] = reduce ( z ) / A [i , i ]
8 return x
P
For every value of i, lines 6 and 7 are executed in parallel in j; the work is O(i) and depth O(log i). But we loop over i sequentially, so the total work is W (n) = i i=
O n2 and the total depth is
P
i log i = O(n log n). The parallelism P = O(n/ log n). An alternative approach is the following
1 def forwardsubst1 (A ,x , n ) :
2 for j in range ( n ) : # loop sequentially over columns
3 x [ j ] = b [ j ]/ A [j , j ]
4 for i in prange ( j +1: n ) :
5 b [ i ] -= A [i , j ] * x [ j ]
6 return x
Now, the complexity of this scheme is W = O n2 and D = O(n) so that the parallelism is improved to P = O(n)
6. Propose a d-dimensional hypercube reduction algorithm on p processors for any array size n ≥ p and for any p; i.e., p does not have to be a power of two. Derive an
expression for the wall-clock time T (n, p) as a function of p and n. Your estimate should include communication costs in terms of latency, bandwidth, and message size.
(Hint: you can either use a single large hypercube or use several smaller hypercubes.)
1 def hy p e r c u b e _ r e d u c e ( A , p , id ) :
2 # id : is this processor id , A : local array owned by processor id , p : number of processors
3 n = len ( A )
4 s =0
5 if n >0: s = sum ( A ) # local sum if array is n o n empty
6 d = ceil ( log2 ( p ) ) # embed p to larger hipercube
7
8 # hypercube loop
9 mask = 0
10 for k in range ( d ) :
11 r =1 < < k
12 if id & mask :
13 partner = id ** r
14 if id & r :
15 send (s , partner )
16 else :
17 if partner <p : # modification of standard hypercube algo
18 recv ( sp , partner )
19 s=s+s p
20 mask = mask ** r
21 return s # processor with id = 0 has the correct answer
The complexity of the sequential part is T (n, p) = O(n/p) The complexity of the reduction is T (n, p) = O(d) = O(log p). The overall computation time is O(n/p) +
O(log p). The overall communication time is O(log p)(l + 1/b), where l is the latency and b is the bandwidth. Notice line 22: that’s the only difference with the standard
hypercube. We run a hypercube with 2d processors where d is defined at line 10 . Since p ≤ 2d we need the check at line 22 .
1. Given x0 and {ai , bi }n−1
i=1 , we wish to compute
xi = ai xi−1 + bi , i = 1, · · · , n − 1
The obvious sequential algorithm for this problem has O(n) work and depth complexity. State a parallel algorithm for this problem using the work-depth language
programming model. Your algorithm should be work optimal and has O(log n) depth.
2. State a PRAM algorithm for merging two sorted arrays. Derive its time and work complexity. You may assume that the number of threads and the array size are
powers of two and that the two arrays have no duplicates.
0: procedure MergeArrays(A[1..n], B[1..n]){A, B are sorted arrays of size n}
0: Let C[1..2n] be a temporary array
0: Distribute the elements of A and B among the processors
0: Each processor Pi compares its two assigned elements: A[i] and B[i], and puts the smaller element in C[i] and the larger in C[n + i]
0: Use log n parallel prefix sums to merge the sorted halves of C into a fully sorted array
0: Copy the sorted array back into A
0: end procedure=0
The time complexity of this algorithm is O(log n), which is the depth of the algorithm. Each iteration of step 2 takes O(n) time since there are n/2 pairs of subarrays
to merge, and each pair requires merging 2(i−1) elements. Therefore, the total time complexity is O(n log n).
The work complexity is O(n log n), which is the product of the number of threads and the depth of the algorithm. Each thread works on a subarray of size 1, so the
total number of elements processed by all threads is n. The work of each iteration of step 2 is also O(n), so the total work is O(n log n).
n(n+1) 2
h i
n(n+1) Pn n(n+1)(2n+1) Pn
•
Pn Pn 2 3
i=1 c = cn i=1 i = 2 i=1 i = 6 i=1 i = 2
• logb (M · N ) = logb M + logb N , logb M
= logb M − logb N , logb M k = k · logb M , logb (1) = 0, logb (b) = 1 logb bk = k, blogb (k) = k
N
Let T be the total wall-clock time,sT be the sequential part (0 ≤ s ≤ 1) and (1 − s)T be the parallelizable part (assume embarrassingly parallel). Assume p cores
• D(n) : Depth, longest chain of dependencies • Speedup S : Best sequential time / time on P • Tp = sT + (1 − s) Tp → S = T
Tp
(equal to number of edges in the longest path Tn,1, best
processors = S = Tn,p
in DAG) This is a lower bound on T (n, p) for • Amdhal’s speed-up law: S = 1
→
all values of p on any parallel machine. s+ 1−s
p
1
• Efficiency: speedup / perfect speedup = E = S
s for large p
• W (n) : Work, the total number of opera- p
tions cumulative across p. (equal to number • Gustafson Fix: Sequential part should be in-
of nodes) • Work Efficiency or Optimality sequential → dependent of the problem size s(p) → 0 as
parallel → ensure work optimality (not always p → ∞. =⇒ Increase problem size, with
W (n)
• Parallelism (higher is better): P (n) = D(n)
possible) increasing parallelism
DAG: directed acyclic graph PRAM Model: PRAM Model:
• Input is a node that has no incoming edge • Programmer writes code for each thread • PRAM (Scheduling Principle): T (n, p) =
O(W (n)/p + D(n))
• Output is a node that has no outgoing edges • Number of threads specified
• Instructions represented by internal nodes in- • Assumes synchronous mode of operation
• Each thread identified by unique ID
coming edges are operands; outgoing edges are
outputs • Memory model • Asynchronous models also exist
• Edge (u, v) indicates instruction v must take • Each thread has local memory along with • Complexity estimates: T (n, p), W (n, p), P =
place after instruction u shared memory W/T
Language: Hypercube: Cache
• Parallel for % data decomposition • d dimensions • Hit” means loading/storing data in cache,
whereas Miss” means data not found in cache
• Parallel do % task parallelism • 2 processes
d
and needing to be retrieved from the main
memory
• Recursion % divide and conquer • d-bit representation
• Cache data is transferred in blocks or lines
Metrics: Depth/Time, Work Opt • High d.... Low d bit
• DAG/Work depth: Small D(n), so that P (n) = • two processes connected iff their ids differ by a • Conceptual categories of misses are as follows:
W (n)/D(n) → ∞ single bit – Compulsory: Occurs during the first ac-
• PRAM: Small T (n, p) so that speed up • diameter = log p cess
T (n)/T (n, p) → p – Capacity: Cache cannot contain all
• Connectivity: log p
• Work optimality: parallel algorithm has blocks needed to execute the program
the same complexity (up to constants) with • Width = p/2
– Conflict: Block replaced but then ac-
the best sequential algorithm i.e T (n, p) ∼
• Total links: p log p/2 cessed soon thereafter
Tseq (n)/p (or as small as possible)
Memory Model • LI CACHE hit, ∼ 4 cycles, L2 CACHE hit,
Scan (or prefix sum)
∼ 12 cycles
• m ≡ No. of words moved from slow to fast
• Sequential W = O(n), D = O(n)
memory • L3 CACHE hit, line unshared ∼ 40 cycles
• Simple parallel(recursive) W =
• f ≡ No. of flops • L3 CACHE hit, shared line in another core ∼ 65
O(n log n), D = O(log n)
cycles
• α ≡ Time per slow memory operation
• Work optimal(recursive) W = O(n), D =
O(log n) • L3 CACHE hit, modified in another core ∼ 75
• τ ≡ Time per flop cycles
Select f
• q≡ m = Flop-to-mop ratio ⇐ Computational • DRAM: ∼ 120 cycles, Off-socket DRAM: ∼ 300
• Sequential W = O(n), D = O(n) intensity cycles
• Scan based W = O(n), D = O(log n) • T =f ·τ +m·α=f ·τ · 1+ α 1
τ · q • HD ∼ 24, 000 cycles