CS 267: Applications of Parallel Computers Graph Partitioning James Demmel and Horst Simon www.cs.berkeley.edu/~demmel/cs267_Spr10 02/11/2010 CS267 Lecture 8 1 Outline of Graph Partitioning Lecture Review definition of Graph Partitioning problem Overview of heuristics Partitioning with Nodal Coordinates Ex: In finite element models, node at point in (x,y) or (x,y,z) space Partitioning without Nodal Coordinates Ex: In model of WWW, nodes are web pages Multilevel Acceleration BIG IDEA, appears often in scientific computing
Comparison of Methods and Applications 02/11/2010 CS267 Lecture 8 Definition of Graph Partitioning Given a graph G = (N, E, WN, WE) N = nodes (or vertices), WN = node weights E = edges WE = edge weights 1 (2) 2 (2) 4 2 5 (1) 3 5
6 (2) 1 3 (1) 2 4 (3) 1 2 1 8 (1) 6 7 (3) Ex: N = {tasks}, WN = {task costs}, edge (j,k) in E means task j sends WE(j,k) words to task k Choose a partition N = N1 U N2 U U NP such that The sum of the node weights in each Nj is about the same The sum of all edge weights of edges connecting all different pairs Nj and Nk is minimized
Ex: balance the work load, while minimizing communication Special case of N = N1 U N2: Graph Bisection 02/11/2010 CS267 Lecture 8 3 Definition of Graph Partitioning Given a graph G = (N, E, WN, WE) N = nodes (or vertices), WN = node weights E = edges WE = edge weights 1 (2) 2 (2) 4 2 5 (1)
3 5 6 (2) 1 3 (1) 2 4 (3) 1 2 1 8 (1) 6 7 (3) Ex: N = {tasks}, WN = {task costs}, edge (j,k) in E means task j sends WE(j,k) words to task k Choose a partition N = N1 U N2 U U NP such that
The sum of the node weights in each Nj is about the same The sum of all edge weights of edges connecting all different pairs Nj and Nk is minimized (shown in black) Ex: balance the work load, while minimizing communication Special case of N = N1 U N2: Graph Bisection 02/11/2010 CS267 Lecture 8 4 Some Applications Telephone network design Original application, algorithm due to Kernighan Load Balancing while Minimizing Communication Sparse Matrix times Vector Multiplication Solving PDEs N = {1,,n}, (j,k) in E if A(j,k) nonzero, WN(j) = #nonzeros in row j, WE(j,k) = 1 VLSI Layout N = {units on chip}, E = {wires}, WE(j,k) = wire length
Sparse Gaussian Elimination Used to reorder rows and columns to increase parallelism, and to decrease fill-in Data mining and clustering Physical Mapping of DNA Image Segmentation 02/11/2010 CS267 Lecture 8 5 Sparse Matrix Vector Multiplication y = y +A*x declare A_local, A_remote(1:num_procs), x_local, x_remote, y_local y_local = y_local + A_local * x_local for all procs P that need part of x_local send(needed part of x_local, P) for all procs P owning needed part of x_remote receive(x_remote, P) y_local = y_local + A_remote(P)*x_remote 02/11/2010
CS267 Lecture 8 6 Cost of Graph Partitioning Many possible partitionings to search Just to divide in 2 parts there are: n choose n/2 = n!/((n/2)!)2 ~ sqrt(2/(n)*2n possibilities Choosing optimal partitioning is NP-complete (NP-complete = we can prove it is a hard as other well-known hard problems in a class Nondeterministic Polynomial time) Only known exact algorithms have cost = exponential(n) We need good heuristics 02/11/2010 CS267 Lecture 8 7
Outline of Graph Partitioning Lectures Review definition of Graph Partitioning problem Overview of heuristics Partitioning with Nodal Coordinates Ex: In finite element models, node at point in (x,y) or (x,y,z) space Partitioning without Nodal Coordinates Ex: In model of WWW, nodes are web pages Multilevel Acceleration BIG IDEA, appears often in scientific computing Comparison of Methods and Applications 02/11/2010 CS267 Lecture 8 First Heuristic: Repeated Graph Bisection To partition N into 2k parts bisect graph recursively k times Henceforth discuss mostly graph bisection 02/11/2010
CS267 Lecture 8 9 Edge Separators vs. Vertex Separators Edge Separator: Es (subset of E) separates G if removing Es from E leaves two ~equal-sized, disconnected components of N: N 1 and N2 Vertex Separator: Ns (subset of N) separates G if removing Ns and all incident edges leaves two ~equal-sized, disconnected components of N: N1 and N2 G = (N, E), Nodes N and Edges E Es = green edges or blue edges Ns = red vertices Making an Ns from an Es: pick one endpoint of each edge in Es |Ns| |Es| Making an Es from an Ns: pick all edges incident on Ns 02/11/2010 CS267 8 |Es| d * |Ns| where d is theLecture
maximum degree of the graph 10 Overview of Bisection Heuristics Partitioning with Nodal Coordinates Each node has x,y,z coordinates partition space Partitioning without Nodal Coordinates E.g., Sparse matrix of Web documents A(j,k) = # times keyword j appears in URL k Multilevel acceleration (BIG IDEA) Approximate problem by coarse graph, do so recursively 02/11/2010 CS267 Lecture 8 11 Outline of Graph Partitioning
Lectures Review definition of Graph Partitioning problem Overview of heuristics Partitioning with Nodal Coordinates Ex: In finite element models, node at point in (x,y) or (x,y,z) space Partitioning without Nodal Coordinates Ex: In model of WWW, nodes are web pages Multilevel Acceleration BIG IDEA, appears often in scientific computing Comparison of Methods and Applications 02/11/2010 CS267 Lecture 8 Nodal Coordinates: How Well Can We Do? A planar graph can be drawn in plane without edge crossings Ex: m x m grid of m2 nodes: vertex separatorNs with | Ns| = m = sqrt(|N|) (see earlier slide for m=5 ) Theorem (Tarjan, Lipton, 1979): If G is planar, Ns such that N = N1 U Ns U N2 is a partition,
|N1| <= 2/3 |N| and |N2| <= 2/3 |N| |Ns| <= sqrt(8 * |N|) Theorem motivates intuition of following algorithms 02/11/2010 CS267 Lecture 8 13 Nodal Coordinates: Inertial Partitioning For a graph in 2D, choose line with half the nodes on one side and half on the other In 3D, choose a plane, but consider 2D for simplicity Choose a line L, and then choose a line L perpendicular to it, with half the nodes on either side L 1. Choose a line L through the points L given by a*(x-xbar)+b*(y-ybar)=0, with a2+b2=1; (a,b) is unit vector to L 2.
3. L Project each point to the line For each nj = (xj,yj), compute coordinate Sj = -b*(xj-xbar) + a*(yj-ybar) along L Compute the median (xbar,ybar) Let Sbar = median(S1,,Sn) 4. Use median to partition the nodes Let nodes with Sj < Sbar be in N1, rest in N2 02/11/2010 CS267 Lecture 8 (a,b) 14 Inertial Partitioning: Choosing L
Clearly prefer L, L on left below N1 N1 N2 L L L N2 L Mathematically, choose L to be a total least squares fit of the nodes Minimize sum of squares of distances to L (green lines on last slide) Equivalent to choosing L as axis of rotation that minimizes the moment of inertia of nodes (unit weights) - source of name 02/11/2010
CS267 Lecture 8 15 Inertial Partitioning: choosing L (continued) L (xj , yj ) (a,b) is unit vector perpendicular to L (xbar,ybar) (a,b) j (length of j-th green line)2 = j [ (xj - xbar)2 + (yj - ybar)2 - (-b*(xj - xbar) + a*(yj - ybar))2 ] Pythagorean Theorem = a2 * j (xj - xbar)2 + 2*a*b* j (xj - xbar)*(xj - ybar) + b2 j (yj - ybar)2 = a2 * X1 + 2*a*b* X2 + b 2 * X3 = [a b] * X1 X2 * a X2 X3 b Minimized by choosing
(xbar , ybar) = (j xj , j yj) / n = center of mass (a,b) = eigenvector of smallest eigenvalue of X1 X2 02/11/2010 CS267 Lecture 8 X2 X3 16 Nodal Coordinates: Random Spheres Generalize nearest neighbor idea of a planar graph to higher dimensions Any graph can fit in 3D without edge crossings Capture intuition of planar graphs of being connected to nearest neighbors but in higher than 2 dimensions For intuition, consider graph defined by a regular 3D mesh An n by n by n mesh of |N| = n3 nodes Edges to 6 nearest neighbors Partition by taking plane parallel to 2 axes Cuts n2 =|N|2/3 = O(|E|2/3) edges For the general graphs Need a notion of well-shaped like mesh
02/11/2010 CS267 Lecture 8 17 Random Spheres: Well Shaped Graphs Approach due to Miller, Teng, Thurston, Vavasis Def: A k-ply neighborhood system in d dimensions is a set {D1,,Dn} of closed disks in Rd such that no point in Rd is strictly interior to more than k disks Def: An (,k) overlap graph is a graph defined in terms of 1 and a k-ply neighborhood system {D1,,Dn}: There is a node for each Dj, and an edge from j to i if expanding the radius of the smaller of Dj and Di by > causes the two disks to overlap Ex: n-by-n mesh is a (1,1) overlap graph Ex: Any planar graph is (,k) overlap for some ,k 02/11/2010 CS267 Lecture 8
2D Mesh is (1,1) overlap graph 18 Generalizing Lipton/Tarjan to Higher Dimensions Theorem (Miller, Teng, Thurston, Vavasis, 1993): Let G=(N,E) be an (,k) overlap graph in d dimensions with n=|N|. Then there is a vertex separator Ns such that N = N1 U Ns U N2 and N1 and N2 each has at most n*(d+1)/(d+2) nodes Ns has at most O( * k1/d * n(d-1)/d ) nodes When d=2, same as Lipton/Tarjan Algorithm: Choose a sphere S in Rd Edges that S cuts form edge separator Es Build Ns from Es Choose S randomly, so that it satisfies Theorem with high probability 02/11/2010 CS267 Lecture 8
19 Stereographic Projection Stereographic projection from plane to sphere In d=2, draw line from p to North Pole, projection p of p is where the line and sphere intersect p p p = (x,y) p = (2x,2y,x2 + y2 1) / (x2 + y2 + 1) Similar in higher dimensions 02/11/2010 CS267 Lecture 8 20 Choosing a Random Sphere
Do stereographic projection from Rd to sphere S in Rd+1 Find centerpoint of projected points Any plane through centerpoint divides points ~evenly There is a linear programming algorithm, cheaper heuristics Conformally map points on sphere Rotate points around origin so centerpoint at (0,0,r) for some r Dilate points (unproject, multiply by sqrt((1-r)/(1+r)), project) this maps centerpoint to origin (0,,0), spreads points around S Pick a random plane through origin Intersection of plane and sphere S is circle Unproject circle yields desired circle C in Rd Create Ns: j belongs to Ns if *Dj intersects C 02/11/2010 CS267 Lecture 8 21
Random Sphere Algorithm (Gilbert) 04/11/2007 CS267 Lecture 23 22 Random Sphere Algorithm (Gilbert) 04/11/2007 CS267 Lecture 23 23 Random Sphere Algorithm (Gilbert) 02/11/2010 CS267 Lecture 23 24
Random Sphere Algorithm (Gilbert) 02/11/2010 CS267 Lecture 8 25 Random Sphere Algorithm (Gilbert) 02/11/2010 CS267 Lecture 23 26 Random Sphere Algorithm (Gilbert) 02/11/2010 CS267 Lecture 23 27
Nodal Coordinates: Summary Other variations on these algorithms Algorithms are efficient Rely on graphs having nodes connected (mostly) to nearest neighbors in space algorithm does not depend on where actual edges are! Common when graph arises from physical model Ignores edges, but can be used as good starting guess for subsequent partitioners that do examine edges Can do poorly if graph connection is not spatial: Details at www.cs.berkeley.edu/~demmel/cs267/lecture18/lecture18.html www.cs.ucsb.edu/~gilbert www.cs.bu.edu/~steng 02/11/2010 CS267 Lecture 8 28 Outline of Graph Partitioning
Lectures Review definition of Graph Partitioning problem Overview of heuristics Partitioning with Nodal Coordinates Ex: In finite element models, node at point in (x,y) or (x,y,z) space Partitioning without Nodal Coordinates Ex: In model of WWW, nodes are web pages Multilevel Acceleration BIG IDEA, appears often in scientific computing Comparison of Methods and Applications 02/11/2010 CS267 Lecture 8 Coordinate-Free: Breadth First Search (BFS) Given G(N,E) and a root node r in N, BFS produces A subgraph T of G (same nodes, subset of edges) T is a tree rooted at r Each node assigned a level = distance from r root Level 0
Level 1 N1 Level 2 Level 3 N2 Level 4 Tree edges Horizontal edges Inter-level edges 02/11/2010 CS267 Lecture 8 30 Breadth First Search Queue (First In First Out, or FIFO)
root Enqueue(x,Q) adds x to back of Q x = Dequeue(Q) removes x from front of Q Compute Tree T(NT,ET) NT = {(r,0)}, ET = empty set Initially T = root r, which is at level 0 Enqueue((r,0),Q) Put root on initially empty Queue Q Mark r Mark root as having been processed While Q not empty While nodes remain to be processed (n,level) = Dequeue(Q) Get a node to process For all unmarked children c of n NT = NT U (c,level+1) Add child c to NT ET = ET U (n,c) Add edge (n,c) to ET Enqueue((c,level+1),Q)) Add child c to Q for processing Mark c Mark c as processed
Endfor Endwhile 02/11/2010 CS267 Lecture 8 31 Partitioning via Breadth First Search root BFS identifies 3 kinds of edges Tree Edges - part of T Horizontal Edges - connect nodes at same level Interlevel Edges - connect nodes at adjacent levels No edges connect nodes in levels differing by more than 1 (why?) BFS partioning heuristic N = N1 U N2, where N1 = {nodes at level <= L},
N2 = {nodes at level > L} Choose L so |N1| close to |N2| BFS partition of a 2D Mesh using center as root: N1 = levels 0, 1, 2, 3 N2 = levels 4, 5, 6 02/11/2010 CS267 Lecture 8 32 Coordinate-Free: Kernighan/Lin Take a initial partition and iteratively improve it Kernighan/Lin (1970), cost = O(|N|3) but easy to understand Fiduccia/Mattheyses (1982), cost = O(|E|), much better, but more complicated Given G = (N,E,WE) and a partitioning N = A U B, where |A| = |B|
T = cost(A,B) = {W(e) where e connects nodes in A and B} Find subsets X of A and Y of B with |X| = |Y| Consider swapping X and Y if it decreases cost: newA = (A X) U Y newT = cost(newA , newB) < T = cost(A,B) and newB = (B Y) U X Need to compute newT efficiently for many possible X and Y, choose smallest (best) 02/11/2010 CS267 Lecture 8 33 Kernighan/Lin: Preliminary Definitions
T = cost(A, B), newT = cost(newA, newB) Need an efficient formula for newT; will use E(a) = external cost of a in A = {W(a,b) for b in B} I(a) = internal cost of a in A = {W(a,a) for other a in A} D(a) = cost of a in A = E(a) - I(a) E(b), I(b) and D(b) defined analogously for b in B Consider swapping X = {a} and Y = {b} newA = (A - {a}) U {b}, newB = (B - {b}) U {a} newT = T - ( D(a) + D(b) - 2*w(a,b) ) T - gain(a,b) gain(a,b) measures improvement gotten by swapping a and b Update formulas newD(a) = D(a) + 2*w(a,a) - 2*w(a,b) for a in A, a a newD(b) = D(b) + 2*w(b,b) - 2*w(b,a) for b in B, b b 02/11/2010 CS267 Lecture 8 34 Kernighan/Lin Algorithm
Compute T = cost(A,B) for initial A, B cost = O(|N|2) Repeat One pass greedily computes |N|/2 possible X,Y to swap, picks best Compute costs D(n) for all n in N cost = O(|N|2) Unmark all nodes in N cost = O(|N|) While there are unmarked nodes |N|/2 iterations Find an unmarked pair (a,b) maximizing gain(a,b) cost = O(|N|2) Mark a and b (but do not swap them) cost = O(1) Update D(n) for all unmarked n, as though a and b had been swapped cost = O(|N|) Endwhile At this point we have computed a sequence of pairs (a1,b1), , (ak,bk) and gains gain(1),., gain(k) where k = |N|/2, numbered in the order in which we marked them Pick m maximizing Gain = k=1 to m gain(k) cost = O(|N|) Gain is reduction in cost from swapping (a1,b1) through (am,bm)
If Gain > 0 then it is worth swapping Update newA = A - { a1,,am } U { b1,,bm } cost = O(|N|) Update newB = B - { b1,,bm } U { a1,,am } cost = O(|N|) Update T = T - Gain cost = O(1) endif Until Gain <= 0 02/11/2010 CS267 Lecture 8 35 Comments on Kernighan/Lin Algorithm Most expensive line shown in red, O(n3) Some gain(k) may be negative, but if later gains are large, then final Gain may be positive can escape local minima where switching no pair helps How many times do we Repeat? K/L tested on very small graphs (|N|<=360) and got
convergence after 2-4 sweeps For random graphs (of theoretical interest) the probability of convergence in one step appears to drop like 2-|N|/30 02/11/2010 CS267 Lecture 8 36 Coordinate-Free: Spectral Bisection Based on theory of Fiedler (1970s), popularized by Pothen, Simon, Liou (1990) Motivation, by analogy to a vibrating string Basic definitions Vibrating string, revisited Implementation via the Lanczos Algorithm To optimize sparse-matrix-vector multiply, we graph partition To graph partition, we find an eigenvector of a matrix associated with the graph To find an eigenvector, we do sparse-matrix vector multiply No free lunch ... 02/11/2010
CS267 Lecture 8 37 Motivation for Spectral Bisection Vibrating string Think of G = 1D mesh as masses (nodes) connected by springs (edges), i.e. a string that can vibrate Vibrating string has modes of vibration, or harmonics Label nodes by whether mode - or + to partition into N- and N+ Same idea for other graphs (eg planar graph ~ trampoline) 02/11/2010 CS267 Lecture 8 38 Basic Definitions Definition: The incidence matrix In(G) of a graph G(N,E) is an |N| by |E| matrix, with one row for each node and one column for each edge. If edge e=(i,j) then column e of In(G) is zero except for the i-th and j-th entries, which
are +1 and -1, respectively. Slightly ambiguous definition because multiplying column e of In(G) by -1 still satisfies the definition, but this wont matter... Definition: The Laplacian matrix L(G) of a graph G(N,E) is an |N| by |N| symmetric matrix, with one row and column for each node. It is defined by L(G) (i,i) = degree of node i (number of incident edges) L(G) (i,j) = -1 if i j and there is an edge (i,j) L(G) (i,j) = 0 otherwise 02/11/2010 CS267 Lecture 8 39 Example of In(G) and L(G) for Simple Meshes 02/11/2010 CS267 Lecture 8 40
Properties of Laplacian Matrix Theorem 1: Given G, L(G) has the following properties (proof on 1996 CS267 web page) L(G) is symmetric. This means the eigenvalues of L(G) are real and its eigenvectors are real and orthogonal. In(G) * (In(G))T = L(G) The eigenvalues of L(G) are nonnegative: 0 = 1 2 n The number of connected components of G is equal to the number of i equal to 0. Definition: 2(L(G)) is the algebraic connectivity of G The magnitude of 2 measures connectivity
In particular, 2 0 if and only if G is connected. 02/11/2010 CS267 Lecture 8 42 Spectral Bisection Algorithm Spectral Bisection Algorithm: Compute eigenvector v2 corresponding to 2(L(G)) For each node n of G if v2(n) < 0 put node n in partition N- else put node n in partition N+ Why does this make sense? First reasons... Theorem 2 (Fiedler, 1975): Let G be connected, and N- and N+ defined as above. Then N- is connected. If no v2(n) = 0, then
N+ is also connected. (proof on 1996 CS267 web page) Recall2(L(G)) is the algebraic connectivity of G Theorem 3 (Fiedler): Let G1(N,E1) be a subgraph of G(N,E), so that G1 is less connected than G. Then 2(L(G1)) 2(L(G)) , i.e. the algebraic connectivity of G1 is less than or equal to the algebraic connectivity of G. (proof on 1996 CS267 web page) 02/11/2010 CS267 Lecture 8 43 Spectral Bisection Algorithm Spectral Bisection Algorithm: Compute eigenvector v2 corresponding to 2(L(G)) For each node n of G if v2(n) < 0 put node n in partition N- else put node n in partition N+
Why does this make sense? More reasons... Theorem 4 (Fiedler, 1975): Let G be connected, and N1 and N2 be any partition into part of equal size |N|/2. Then the number of edges connecting N1 and N2 is at least .25 * |N| * 2(L(G)). (proof on 1996 CS267 web page) 02/11/2010 CS267 Lecture 8 44 Motivation for Spectral Bisection (recap) Vibrating string has modes of vibration, or harmonics Modes computable as follows Model string as masses connected by springs (a 1D mesh) Write down F=ma for coupled system, get matrix A Eigenvalues and eigenvectors of A are frequencies and shapes of modes Label nodes by whether mode - or + to get N- and N+ Same idea for other graphs (eg planar graph ~ trampoline) 02/11/2010
CS267 Lecture 23 45 Details for Vibrating String Analogy Force on mass j = k*[x(j-1) - x(j)] + k*[x(j+1) - x(j)] = -k*[-x(j-1) + 2*x(j) - x(j+1)] F=ma yields m*x(j) = -k*[-x(j-1) + 2*x(j) - x(j+1)] Writing (*) for j=1,2,,n yields m * d2 dx2 (*) x(1) 2*x(1) - x(2) 2 -1 x(1) x(1) x(2) -x(1) + 2*x(2) - x(3) -1 2 -1 x(2)
x(2) =-k* =-k* * =-k*L* x(j) -x(j-1) + 2*x(j) - x(j+1) -1 2 -1 x(j) x(j) x(n) 2*x(n-1) - x(n) -1 2 x(n) x(n) (-m/k) x = L*x
02/11/2010 CS267 Lecture 23 46 Details for Vibrating String (continued) -(m/k) x = L*x, where x = [x1,x2,,xn ]T Seek solution of form x(t) = sin(*t) * x0 L*x0 = (m/k)*2 * x0 = * x0 For each integer i, get = 2*(1-cos(i*/(n+1)), x0 = sin(1*i*/(n+1)) sin(2*i*/(n+1)) sin(n*i*/(n+1)) Thus x0 is a sine curve with frequency proportional to i Thus 2 = 2*k/m *(1-cos(i*/(n+1)) or ~ sqrt(k/m)**i/(n+1) L = 2 -1 not quite Laplacian of 1D mesh, -1 2 -1 but we can fix that ...
. -1 02/11/2010 2 CS267 Lecture 8 47 Motivation for Spectral Bisection Vibrating string has modes of vibration, or harmonics Modes computable as follows Model string as masses connected by springs (a 1D mesh) Write down F=ma for coupled system, get matrix A Eigenvalues and eigenvectors of A are frequencies and shapes of modes Label nodes by whether mode - or + to get N- and N+ Same idea for other graphs (eg planar graph ~ trampoline) 02/11/2010 CS267 Lecture 8
48 Eigenvectors of L(1D mesh) Eigenvector 1 (all ones) Eigenvector 2 Eigenvector 3 02/11/2010 CS267 Lecture 8 49 2nd eigenvector of L(planar mesh) 02/11/2010 CS267 Lecture 8 50
4th eigenvector of L(planar mesh) 02/11/2010 CS267 Lecture 8 51 Computing v2 and2 of L(G) using Lanczos Given any n-by-n symmetric matrix A (such as L(G)) Lanczos computes a k-by-k approximation T by doing k matrix-vector products, k << n Choose an arbitrary starting vector r b(0) = ||r|| j=0 repeat j=j+1 q(j) = r/b(j-1) scale a vector (BLAS1) r = A*q(j) matrix vector multiplication, the most expensive step r = r - b(j-1)*v(j-1) axpy, or scalar*vector + vector (BLAS1)
a(j) = v(j)T * r dot product (BLAS1) r = r - a(j)*v(j) axpy (BLAS1) b(j) = ||r|| compute vector norm (BLAS1) until convergence details omitted T = a(1) b(1) b(1) a(2) b(2) b(2) a(3) b(3) b(k-2) a(k-1) b(k-1) b(k-1) a(k) Approximate As eigenvalues/vectors using Ts 02/11/2010 CS267 Lecture 8 52
Spectral Bisection: Summary Laplacian matrix represents graph connectivity Second eigenvector gives a graph bisection Roughly equal weights in two parts Weak connection in the graph will be separator Implementation via the Lanczos Algorithm To optimize sparse-matrix-vector multiply, we graph partition To graph partition, we find an eigenvector of a matrix associated with the graph To find an eigenvector, we do sparse-matrix vector multiply Have we made progress? 02/11/2010 The first matrix-vector multiplies are slow, but use them to learn how to make the rest faster CS267 Lecture 8 53 Outline of Graph Partitioning
Lectures Review definition of Graph Partitioning problem Overview of heuristics Partitioning with Nodal Coordinates Ex: In finite element models, node at point in (x,y) or (x,y,z) space Partitioning without Nodal Coordinates Ex: In model of WWW, nodes are web pages Multilevel Acceleration BIG IDEA, appears often in scientific computing Comparison of Methods and Applications 02/11/2010 CS267 Lecture 8 Introduction to Multilevel Partitioning If we want to partition G(N,E), but it is too big to do efficiently, what can we do? 1) Replace G(N,E) by a coarse approximation Gc(Nc,Ec), and partition Gc instead
2) Use partition of Gc to get a rough partitioning of G, and then iteratively improve it What if Gc still too big? Apply same idea recursively 02/11/2010 CS267 Lecture 8 55 Multilevel Partitioning - High Level Algorithm (N+,N- ) = Multilevel_Partition( N, E ) recursive partitioning routine returns N+ and N- where N = N+ U Nif |N| is small (1) Partition G = (N,E) directly to get N = N+ U NReturn (N+, N- ) else (2) Coarsen G to get an approximation Gc = (Nc, Ec) (3)
(Nc+ , Nc- ) = Multilevel_Partition( Nc, Ec ) (4) (5) Expand (Nc+ , Nc- ) to a partition (N+ , N- ) of N Improve the partition ( N+ , N- ) Return ( N+ , N- ) endif V - cycle: How do we Coarsen? Expand? Improve? (2,3) (4) (5) (2,3) (4) (5)
(2,3) 02/11/2010 (5) (4) (1) 56 Multilevel Kernighan-Lin Coarsen graph and expand partition using maximal matchings Improve partition using Kernighan-Lin 02/11/2010 CS267 Lecture 8 57 Maximal Matching Definition: A matching of a graph G(N,E) is a subset Em of
E such that no two edges in Em share an endpoint Definition: A maximal matching of a graph G(N,E) is a matching Em to which no more edges can be added and remain a matching A simple greedy algorithm computes a maximal matching: 02/11/2010 let Em be empty mark all nodes in N as unmatched for i = 1 to |N| visit the nodes in any order if i has not been matched mark i as matched if there is an edge e=(i,j) where j is also unmatched, add e to Em mark j as matched endif endif endfor CS267 Lecture 8 58 Maximal Matching: Example
02/11/2010 CS267 Lecture 8 59 Example of Coarsening 02/11/2010 CS267 Lecture 8 60 Coarsening using a maximal matching 1) Construct a maximal matching Em of G(N,E) for all edges e=(j,k) in Em 2) collapse matched nodes into a single one Put node n(e) in Nc W(n(e)) = W(j) + W(k)
gray statements update node/edge weights for all nodes n in N not incident on an edge in Em 3) add unmatched nodes Put n in Nc do not change W(n) Now each node r in N is inside a unique node n(r) in N c 4) Connect two nodes in Nc if nodes inside them are connected in E for all edges e=(j,k) in Em for each other edge e=(j,r) or (k,r) in E Put edge ee = (n(e),n(r)) in Ec W(ee) = W(e) If there are multiple edges connecting two nodes in Nc, collapse them, adding edge weights CS267 Lecture 8 02/11/2010 61 Expanding a partition of Gc to a partition of G 02/11/2010 CS267 Lecture 8
62 Multilevel Spectral Bisection Coarsen graph and expand partition using maximal independent sets Improve partition using Rayleigh Quotient Iteration 02/11/2010 CS267 Lecture 8 63 Maximal Independent Sets Definition: An independent set of a graph G(N,E) is a subset Ni of N such that no two nodes in Ni are connected by an edge Definition: A maximal independent set of a graph G(N,E) is an independent set Ni to which no more nodes can be added and remain an independent set A simple greedy algorithm computes a maximal independent set: let Ni be empty for k = 1 to |N| visit the nodes in any order
if node k is not adjacent to any node already in Ni add k to Ni endif endfor 02/11/2010 CS267 Lecture 8 64 Example of Coarsening - encloses domain Dk = node of Nc 02/11/2010 CS267 Lecture 8 65 Coarsening using Maximal Independent Sets Build domains D(k) around each node k in Ni to get nodes in Nc Add an edge to Ec whenever it would connect two such domains
Ec = empty set for all nodes k in Ni D(k) = ( {k}, empty set ) first set contains nodes in D(k), second set contains edges in D(k) unmark all edges in E repeat choose an unmarked edge e = (k,j) from E if exactly one of k and j (say k) is in some D(m) mark e add j and e to D(m) else if k and j are in two different D(m)s (say D(mk) and D(mj)) mark e add edge (mk, mj) to Ec else if both k and j are in the same D(m) mark e add e to D(m) else leave e unmarked endif until no unmarked edges 02/11/2010 CS267 Lecture 8 66
Expanding a partition of Gc to a partition of G Need to convert an eigenvector vc of L(Gc) to an approximate eigenvector v of L(G) Use interpolation: For each node j in N if j is also a node in Nc, then v(j) = vc(j) else use same eigenvector component v(j) = average of vc(k) for all neighbors k of j in Nc end if endif 02/11/2010 CS267 Lecture 8 67 Example: 1D mesh of 9 nodes
02/11/2010 CS267 Lecture 8 68 Improve eigenvector: Rayleigh Quotient Iteration j=0 pick starting vector v(0) from expanding vc repeat j=j+1 r(j) = vT(j-1) * L(G) * v(j-1) r(j) = Rayleigh Quotient of v(j-1) = good approximate eigenvalue v(j) = (L(G) - r(j)*I)-1 * v(j-1) expensive to do exactly, so solve approximately using an iteration called SYMMLQ, which uses matrix-vector multiply (no surprise) v(j) = v(j) / || v(j) || normalize v(j) until v(j) converges Convergence is very fast: cubic 02/11/2010
CS267 Lecture 8 69 Example of convergence for 1D mesh 02/11/2010 70 Outline of Graph Partitioning Lectures Review definition of Graph Partitioning problem Overview of heuristics Partitioning with Nodal Coordinates Ex: In finite element models, node at point in (x,y) or (x,y,z) space Partitioning without Nodal Coordinates Ex: In model of WWW, nodes are web pages Multilevel Acceleration BIG IDEA, appears often in scientific computing
Comparison of Methods and Applications 02/11/2010 CS267 Lecture 8 Available Implementations Multilevel Kernighan/Lin METIS (www.cs.umn.edu/~metis) ParMETIS - parallel version Multilevel Spectral Bisection S. Barnard and H. Simon, A fast multilevel implementation of recursive spectral bisection , Proc. 6th SIAM Conf. On Parallel Processing, 1993 Chaco (www.cs.sandia.gov/CRF/papers_chaco.html) Hybrids possible Ex: Using Kernighan/Lin to improve a partition from spectral bisection Recent package, collection of techniques Zoltan (www.cs.sandia.gov/Zoltan) 02/11/2010
CS267 Lecture 8 72 Comparison of methods Compare only methods that use edges, not nodal coordinates CS267 webpage and KK95a (see below) have other comparisons Metrics Speed of partitioning Number of edge cuts Other application dependent metrics Summary No one method best Multi-level Kernighan/Lin fastest by far, comparable to Spectral in the number of edge cuts www-users.cs.umn.edu/~karypis/metis/publications/main.html
see publications KK95a and KK95b Spectral give much better cuts for some applications Ex: image segmentation See Normalized Cuts and Image Segmentation by J. Malik, J. Shi 02/11/2010 CS267 Lecture 8 73 Number of edges cut for a 64-way partition For Multilevel Kernighan/Lin, as implemented in METIS (see KK95a) # of Nodes Graph
144 4ELT ADD32 AUTO BBMAT FINAN512 LHR10 MAP1 MEMPLUS SHYY161 TORSO 144649 15606 4960 448695 38744 74752 10672 267241 17758 76480 201142
# of Edges 1074393 45878 9462 3314611 993481 261120 209093 334931 54196 152002 1479989 # Edges cut Expected Expected for 64-way # cuts for # cuts for Description 2D mesh 3D mesh partition 6427 31805 3D FE Mesh 88806 2111 7208
2D FE Mesh 2965 1190 3357 32 bit adder 675 11320 67647 3D FE Mesh 194436 3326 13215 2D Stiffness M. 55753 4620 20481 Lin. Prog. 11388 1746 5595 Chem. Eng. 58784 8736 47887 Highway Net. 1388 2252 7856 Memory circuit 17894 4674
20796 Navier-Stokes 4365 7579 39623 3D FE Mesh 117997 Expected # cuts for 64-way partition of 2D mesh of n nodes n1/2 + 2*(n/2)1/2 + 4*(n/4)1/2 + + 32*(n/32)1/2 ~ 17 * n1/2 Expected # cuts for 64-way partition of 3D mesh of n nodes = n2/3 + 2*(n/2)2/3 + 4*(n/4)2/3 + + 32*(n/32)2/3 ~ 11.5 * n2/3 02/11/2010 74 Speed of 256-way partitioning (from KK95a) # of Nodes Graph 144 4ELT ADD32 AUTO
BBMAT FINAN512 LHR10 MAP1 MEMPLUS SHYY161 TORSO 144649 15606 4960 448695 38744 74752 10672 267241 17758 76480 201142 Partitioning time in seconds # of Multilevel Multilevel
Edges Spectral Kernighan/ Bisection Lin 1074393 607.3 48.1 45878 25.0 3.1 9462 18.7 1.6 3314611 2214.2 179.2 993481 474.2 25.5 261120 311.0 18.0 209093
142.6 8.1 334931 850.2 44.8 54196 117.9 4.3 152002 130.0 10.1 1479989 1053.4 63.9 Description 3D FE Mesh 2D FE Mesh 32 bit adder 3D FE Mesh 2D Stiffness M. Lin. Prog. Chem. Eng. Highway Net.
Memory circuit Navier-Stokes 3D FE Mesh Kernighan/Lin much faster than Spectral Bisection! 02/11/2010 CS267 Lecture 8 75 Beyond Simple Graph Partitioning Undirected graphs model symmetric matrices, not unsymmetric ones More general graph models include: Hypergraph: nodes are computation, edges are communication, but connected to a set (>= 2) of nodes HMETIS package Bipartite model: use bipartite graph for directed graph Multi-object, Multi-Constraint model: use when single structure
may involve multiple computations with differing costs For more see Bruce Hendricksons web page www.cs.sandia.gov/~bahendr/partitioning.html Load Balancing Myths, Fictions & Legends 02/11/2010 CS267 Lecture 8 76 Extra Slides 02/11/2010 CS267 Lecture 8 77 Coordinate-Free Partitioning: Summary Several techniques for partitioning without coordinates Breadth-First Search simple, but not great partition Kernighan-Lin good corrector given reasonable partition
Spectral Method good partitions, but slow Multilevel methods Used to speed up problems that are too large/slow Coarsen, partition, expand, improve Can be used with K-L and Spectral methods and others Speed/quality For load balancing of grids, multi-level K-L probably best For other partitioning problems (vision, clustering, etc.) spectral may be better Good software available 03/09/2009 CS267 Lecture 13 78 Is Graph Partitioning a Solved Problem? Myths of partitioning due to Bruce Hendrickson 1. Edge cut = communication cost 2. Simple graphs are sufficient
3. Edge cut is the right metric 4. Existing tools solve the problem 5. Key is finding the right partition 6. Graph partitioning is a solved problem Slides and myths based on Bruce Hendricksons: Load Balancing Myths, Fictions & Legends 03/09/2009 CS267 Lecture 13 79 Myth 1: Edge Cut = Communication Cost Myth1: The edge-cut deceit edge-cut = communication cost Not quite true: #vertices on boundary is actual communication volume Do not communicate same node value twice
Cost of communication depends on # of messages too ( term) Congestion may also affect communication cost Why is this OK for most applications? Mesh-based problems match the model: cost is ~ edge cuts Other problems (data mining, etc.) do not 03/09/2009 CS267 Lecture 13 80 Myth 2: Simple Graphs are Sufficient Graphs often used to encode data dependencies Do X before doing Y Graph partitioning determines data partitioning Assumes graph nodes can be evaluated in parallel Communication on edges can also be done in parallel Only dependence is between sweeps over the graph More general graph models include:
Hypergraph: nodes are computation, edges are communication, but connected to a set (>= 2) of nodes Bipartite model: use bipartite graph for directed graph Multi-object, Multi-Constraint model: use when single structure may involve multiple computations with differing costs 03/09/2009 CS267 Lecture 13 81 Myth 3: Partition Quality is Paramount When structure are changing dynamically during a simulation, need to partition dynamically Speed may be more important than quality Partitioner must run fast in parallel Partition should be incremental Change minimally relative to prior one Must not use too much memory Example from Touheed, Selwood, Jimack and Bersins
1 M elements with adaptive refinement on SGI Origin Timing data for different partitioning algorithms: Repartition time from 3.0 to 15.2 secs Migration time : 17.8 to 37.8 secs Solve time: 2.54 to 3.11 secs 03/09/2009 CS267 Lecture 13 82 References Details of all proofs on Jim Demmels 267 web page A. Pothen, H. Simon, K.-P. Liou, Partitioning sparse matrices with eigenvectors of graphs, SIAM J. Mat.
Anal. Appl. 11:430-452 (1990) M. Fiedler, Algebraic Connectivity of Graphs, Czech. Math. J., 23:298-305 (1973) M. Fiedler, Czech. Math. J., 25:619-637 (1975) B. Parlett, The Symmetric Eigenproblem, Prentice-Hall, 1980 www.cs.berkeley.edu/~ruhe/lantplht/lantplht.html www.netlib.org/laso 03/09/2009 CS267 Lecture 13 83 Summary Partitioning with nodal coordinates: Inertial method Projection onto a sphere Algorithms are efficient Rely on graphs having nodes connected (mostly) to nearest neighbors in space Partitioning without nodal coordinates:
Breadth-First Search simple, but not great partition Kernighan-Lin good corrector given reasonable partition Spectral Method good partitions, but slow Today: Spectral methods revisited Multilevel methods 03/109/2009 CS267 Lecture 13 84 Another Example Definition: The Laplacian matrix L(G) of a graph G(N,E) is an |N| by |N| symmetric matrix, with one row and column for each node. It is defined by L(G) (i,i) = degree of node I (number of incident edges) L(G) (i,j) = -1 if i != j and there is an edge (i,j) L(G) (i,j) = 0 otherwise 1 4
G= 2 L(G) = 3 5 2 -1 -1 -1 2 -1 -1 -1 4 0 0 -1 0 0 -1 0 0 0 0 -1 -1 2 -1 -1 2 Hidden slide 03/09/2009
CS267 Lecture 13 85
Walking School Bus and Its Sufficiency to Overcome Barriers ...
Walking School Bus potential "big bang" for virtually no "bucks"! Uses existing school personnel or volunteers. Simple strategy. Does not interfere with instruction time. Summary. Walking School Bus sufficient to overcome many barriers to active commute to school.