# MATH 3600 - HW24

## Graphs

A graph is an ordered pair, G = (V,E), where V is a set of vertices and  E = {{u,v} | u ε V and v ε V} is the set of Edges.

More precisely, E is a collection of subsets of V, each of size 2, and the elements of E are called edges.

We will start be defining the class of Graph objects. This class will store V and E as sets.  The vertices can be any kind of immutable object.  The edges in E will be ordered pairs (tuples) of vertices.  Tuples have the advantage of being immutable objects.  This use of immutable objects will allow us to implement the several of the methods of the Graph class using Python's extremely efficient Dictionary class.

A Graph will own the sets V and E which define it. By using the set class for V and E we eliminate the concern about multiple copies of a vertex.  It is also that case that a tuple (u,v) cannot be accidentally stored more than once, but we will need to avoid accidentally storing both (u,v) and (v,u). A Graph will also own a Dictionary called Adj, which is a particularly efficient data structure for algorithms that would naturally exploit the relationships stored in the Adjacency matrix.  It also provides an efficient record of the degree for each vertex. V, E, and Adj are not treated as private objects so that we can directly explore possible algorithms which use them.

In HW13 we introduced the initial version of the Graph class which was very limited while we began to explore methods to be added to the class.  Methods which are both relevant to the effective use of Graphs and capable of being efficiently implemented. In that homework we developed a code for plotting a graph and two codes for determining if a graph is connected. One approach was a recursive code that used a depth-first search. The second approach used a queue and a breadth first search. The Graph class in this homework incorporates what we learned by including methods for plotting the graph, for determining whether the graph is connected, and for creating a list of the components of the graph. These last two methods are based on the breadth-first search algorithm, which requires that we include the code for the Queue class.

Be sure to check out this new Graph class by doing tab completion on one or more of the Graph objects that you create.

A Queue Class

The following block defines a queue class with five methods - enqueue, dequeue, peek, clear, and display. However, instead of following the Stack class approach and having this own a private deque object, this queue is a subclass of the deque class. So this queue object is a deque. One consequence is that this queue object will inherit the methods of the deque class, which includes the clear method. This can be confirmed with the tab completion in the block following the definition of the Queue class.

import collections ''' This block defines a Queue class as a subclass of the collections.deque class This class inherits all the methods of the deque class and adds the methods enqueue, dequeue, peek, and display. Daniel D. Warner February 26, 2012 ''' class Queue(collections.deque): def enqueue(self,x): ''' enqueue adds an item to the end of the queue.''' self.appendleft(x) def dequeue(self): ''' dequeue removes an item from the front of the queue.''' if 0 < len(self): return self.pop() else: return None def peek(self): if 0 < len(self): p = self.pop() self.append(p) return p else: return None def display(self): print self
# Build a simple example V = ['A', 'B', 'C', 'D', 'E', 'F'] Vcoord = {'A':(1,10),'B':(5,5),'C':(1,-10),'D':(5,-5),'E':(3,1),'F':(-5,1)} E2 = [('A','E'),('C','F'),('B','D'),('A','F'),('C','E'),('A','B'),('C','D')] G2 = Graph(V,E2) print 'G2.V',G2.V print 'G2.E',G2.E print 'G2.Adj',G2.Adj g2 = G2.plot(Vcoord) show(g2,figsize=(3,3),axes=False) print 'Degrees' for v in V: print v,G2.Degree(v) print G2.build_Components() if G2.is_connected(): print 'G2 is connected' else: print 'G2 is NOT connected'

Minimum Weighted Spanning Tree

If we include a weight for each edge, then we can make a simple modification to the build_spanning_tree algorithm to obtain Prim's algorithm for finding a Minimum Weighted Spanning.  This is a greedy algorithm that will find a spanning tree where the sum of the weights on the edges is minimized.  That is, the sum of the weights on the edges of the tree found by Prim's algorithm will be less than or equal to the sum of the weights on the edges of any other spanning tree.

The essential step is to replace the queue, which stores the vertices that have been added to the component, with a priority queue, which contains the edges (prioritizedby their weights) that cross the boundary of the component.  The code for Prim's algorithm contains three print statements to provide some insight into how the algorithm behaves.

def Prim(G,W): # simplify edge lookup WW = copy(W) for edge in W: s,t = edge WW[(t,s)] = W[edge] nV = len(G.V) Visited = {} for v in G.V: Visited[v] = False spanning_tree = [] PQ = PriorityQueue() weight = 0 for u in G.V: Visited[u] = True for v in G.Adj[u]: w = WW[(u,v)] PQ.push(w,(u,v)) while PQ.is_not_empty() and len(spanning_tree)<(nV-1): # print PQ._pq w,edge = PQ.pop() s,t = edge # print (s,t),WW[(s,t)] if not Visited[t]: Visited[t] = True spanning_tree.append((s,t)) weight += WW[(s,t)] # print spanning_tree,weight for v in G.Adj[t]: if not Visited[v]: w = WW[(t,v)] PQ.push(w,(t,v)) return spanning_tree,weight

## A Priority Queue Class

A priority queue is a container where each item is inserted with a certain priority.  When an item is removed, the item with the highest priority is the item that is removed.  Priority queues can be implemented using a heap.  A heap is a Complete Binary Tree with a partial ordering.  A complete binary tree is a binary tree where the tree is filled out  from top to bottom, left to right.  To be a heap, the complete binary tree must satisfy the following partial ordering criteria.

The priority of the parent must always be equal to or higher than the priority of its children.

This is often referred to as the heap property.  The heap property is an invariant which must be true before and after pushing or popping an item.

There is a natural embedding of a complete binary tree into an array using the functions

• parent(n) = (n-1)//2
• leftchild(n) = 2*n + 1
• rightchild(n) = 2*n + 2

The following version of the Priority Queue makes the bubble-up and sift-down routines explicit, so that it is possible to extend the class to allow for changes in the priority and other features that are important for discrete event simulation.

This class has been extended to include the routines priority_of(object)  and adjust(new_priority, object). These extensions are necessary for the efficient implementation of Dijkstra's algorithm.

class PriorityQueue(): ''' The arguments passed to a PriorityQueue must consist of a priority and an object. It must be possible to compare priorities using the < operator. By adding the dictionary indx, it is possible to change the priority of an object in the priority queue. ''' def __init__(self): self._pq = [] self._indx = {} def _parent(self,n): return (n-1)//2 def _leftchild(self,n): return 2*n + 1 def _rightchild(self,n): return 2*n + 2 def push(self, pri, obj): self._pq.append((pri,obj)) n = len(self._pq) self._indx[obj] = n-1 self._bubble_up(n-1) # print 'push pq',self._pq # print 'push index',self._indx def _bubble_up(self,c): while 0<c: c_item = self._pq[c] p = self._parent(c) p_item = self._pq[p] if c_item < p_item: self._pq[p] = c_item self._pq[c] = p_item c_obj = c_item p_obj = p_item self._indx[c_obj] = p self._indx[p_obj] = c c = p else: break def pop(self): n = len(self._pq) if 0==n: obj = None elif 1==n: pri,obj = self._pq.pop() self._indx[obj] = None else: pri,obj = self._pq self._indx[obj] = None self._pq = self._pq.pop() pri1,obj1 = self._pq self._indx[obj1] = 0 self._sift_down(0) # print 'pop pq',self._pq # print 'pop index',self._indx return pri,obj def _sift_down(self,p): n = len(self._pq) while p<n: p_item = self._pq[p] lc = self._leftchild(p) if n <= lc: break c_item = self._pq[lc] c = lc rc = self._rightchild(p) if rc < n: r_item = self._pq[rc] if r_item < c_item: c_item = r_item c = rc if p_item <= c_item: break self._pq[p] = c_item self._pq[c] = p_item c_obj = c_item p_obj = p_item self._indx[c_obj] = p self._indx[p_obj] = c p = c def priority_of(self,obj): if obj in self._indx.keys(): i = self._indx[obj] if not i == None: pri,obj = self._pq[i] return pri return None def adjust(self, new_pri, obj): if obj in self._indx: i = self._indx[obj] item = self._pq[i] old_pri = item self._pq[i] = (new_pri,obj) if new_pri < old_pri: self._bubble_up(i) else: self._sift_down(i) else: self.push(new_pri,obj) def is_not_empty(self): if 0 < len(self._pq): return True else: return False

Create a Network by adding weights to the Graph, G2.  Next find the Minimum Weighted Spanning Tree using Prim's algorithm. Then display the Network and the MST.

W = {('A','B'):7, ('A','E'):18, ('C','F'):21, ('D','C'):13, ('B','D'):14, ('A','F'):19, ('C','E'):16} G2 = Graph(V,W.keys()) G2plot = G2.plot(Vcoord) WW = copy(W) for edge in W.keys(): w = W[edge] s,t = edge WW[(t,s)] = w tree,weight = Prim(G2,WW) if len(tree)==(len(G2.V)-1): print 'Minimum Spanning Tree:',tree print 'Weight:',weight for edge in G2.E: s,t = edge w = WW[(s,t)] xs,ys = Vcoord[s] xt,yt = Vcoord[t] xw,yw = ((xs+xt)/2,1/2+(ys+yt)/2) G2plot += text(w,(1.1*xw,1.1*yw)) for edge in tree: s,t = edge xs,ys = Vcoord[s] xt,yt = Vcoord[t] G2plot += line2d([(xs,ys),(xt,yt)],rgbcolor=(1,0,0)) show(G2plot,figsize=(3,3),axes=False) else: print 'Disconnected Subtree:',tree

The Prime Maze Graph

The following blocks define and display the Prime Maze Graph.

# The dictionary of edges and weights for the Prime Maze WPM = {('A', 'B'): 2, ('A', 'C'): 8, ('A', 'D'): 10, ('A', 'E'): 5,\ ('B', 'D'): 14, ('B', 'J'): 17, ('C', 'D'): 12, ('C', 'G'): 3,\ ('C', 'H'): 10, ('D', 'H'): 4, ('E', 'F'): 2, ('E', 'N'): 6,\ ('F', 'G'): 6, ('F', 'O'): 6, ('G', 'O'): 4, ('H', 'T'): 8,\ ('I', 'J'): 8, ('I', 'K'): 4, ('I', 'L'): 6, ('J', 'M'): 10,\ ('J', 'S'): 6, ('K', 'P'): 12, ('K', 'Q'): 12, ('L', 'M'): 12,\ ('L', 'Q'): 8, ('M', 'R'): 8, ('N', 'O'): 2, ('N', 'V'): 4,\ ('O', 'V'): 12, ('P', 'T'): 6, ('P', 'U'): 4, ('Q', 'R'): 10,\ ('R', 'S'): 12, ('S', 'W'): 13, ('T', 'U'): 14, ('U', 'W'): 10,\ ('V', 'W'): 2 }
PMV = ['A','B','C','D','E','F','G','H','I','J','K','L','M','N','P','Q','R','S','T','U','V','W'] PMcoord = {'B':(-3,3),'A':(0,3),'E':(3,3)} PMcoord['J']=(-3,2); PMcoord['D']=(-1,2); PMcoord['C']=(1,2); PMcoord['F']=(2,2) PMcoord['I']=(-1,1); PMcoord['H']=(0,1); PMcoord['G']=(1,1) PMcoord['M']=(-2,0); PMcoord['L']=(-1,0); PMcoord['K']=(0,0); PMcoord['T']=(1,0); PMcoord['O']=(2,0); PMcoord['N']=(3,0) PMcoord['R']=(-2,-1); PMcoord['Q']=(-1,-1); PMcoord['P']=(0,-1); PMcoord['U']=(2,-1) PMcoord['S']=(-3,-2); PMcoord['W']=(0,-2); PMcoord['V']=(3,-2) # print PMcoord

Using the Prime Maze Graph, PMG, and its weights, WPM, find the minimum weighted spanning tree using the Prim function.  Using a copy of Gpm, the plot for the Prime Maze Graph, superimpose a plot of the resulting tree with the edges drawn in red.

# Construct the Prime Maze Graph PMG = Graph([],WPM.keys()) # print 'Graph of the Prime Maze' Gpm = PMG.plot(PMcoord) for edge in PMG.E: s,t = edge xs,ys = PMcoord[s] xt,yt = PMcoord[t] w = WPM[edge] Gpm += text(w,((xs+xt)/2, 0.2+(ys+yt)/2)) # show(Gpm,axes=False) WWW = copy(WPM) for edge in PMG.E: w = WWW[edge] s,t = edge WWW[(t,s)] = w tree,weight = Prim(PMG,WWW) if len(tree)==(len(PMG.V)-1): print '\nMinimum Spanning Tree:',tree print 'Weight:',weight for edge in tree: s,t = edge xs,ys = PMcoord[s] xt,yt = PMcoord[t] Gpm += line2d([(xs,ys),(xt,yt)],rgbcolor=(1,0,0)) show(Gpm,axes=False) else: print 'Disconnected Subtree:',tree

Shortest Paths

We will start by looking at 3 algorihtms for finding the shortest path between two vertices in a network.  The two networks will be G2 and PMG from the previous examples.

For G2 we will find the shortest paths from A to B, A to C, and A to D.

For PMG we will find the shortest paths from A to B, A to R, A to O, and A to W.

The first algorithm is a depth-first, recursive, back-tracking search.  Once a vertex has been added to the path it is no longer available.  We will use a dictionary labeled, available, rather than a linear search to see if the vertex in already in the path.  After we have found one path to the last vertex, we can prune all future branches of the search tree.

def Shortest_Path_df(G,W,path,available,last_v,distance,best,count,max_path): # print path,distance,best,count count += 1 # Only used to obtain efficiency statistics max_path = max(max_path,len(path)) current = path[-1] if current == last_v: # Found a path print path,distance best = distance else: for v in G.Adj[current]: if available[v]: w = W[(current,v)] if distance+w < best: path.append(v) available[v] = False best,count,max_path = Shortest_Path_df(G,W,path,available,last_v,distance+w,best,count,max_path) path.pop() available[v] = True return best,count,max_path def Shortest_Path_DF(G,W,first_v,last_v): available = {} for v in G.V: available[v] = True available[first_v] = False # Start best with a gross overestimate of the distance best = 0 for w in W: best += W[w] best,count,max_path = Shortest_Path_df(G,W,[first_v],available,last_v,0,best,0,0) print best,count,max_path
Shortest_Path_DF(G2,WW,'A','B') timeit("Shortest_Path_DF(G2,WW,'A','B')")
Shortest_Path_DF(G2,WW,'A','C') timeit("Shortest_Path_DF(G2,WW,'A','C')")
Shortest_Path_DF(G2,WW,'A','D') timeit("Shortest_Path_DF(G2,WW,'A','D')")
Shortest_Path_DF(PMG,WWW,'A','B') timeit("Shortest_Path_DF(PMG,WWW,'A','B')")
Shortest_Path_DF(PMG,WWW,'A','R') timeit("Shortest_Path_DF(PMG,WWW,'A','R')")
Shortest_Path_DF(PMG,WWW,'A','O') timeit("Shortest_Path_DF(PMG,WWW,'A','O')")
Shortest_Path_DF(PMG,WWW,'A','W') timeit("Shortest_Path_DF(PMG,WWW,'A','W')")

The Breadth-First Search uses a queue that stores tuples consisting of a path and the distance.

Note that if the distance between every adjacent vertex was the same then the breadth-first search could stop as soon as it found the last vertex.  However, since this isn't the case we must keep searching as long as there are candidates. It is also the case that since there are many different paths in the queue, checking for availability requires a linear search to see if the vertex is in the path.

def Shortest_Path_BF(G,W,first_v,last_v): # Start best with a gross overestimate of the distance best = 0 for w in W: best += W[w] Q = Queue() Q.enqueue(([first_v],0)) q_size = 1 max_q = q_size # Each item in the queue is a (path,distance) tuple while 0<len(Q): path,distance = Q.dequeue() q_size -= len(path) if distance < best: # This path is still a contender current = path[-1] for v in G.Adj[current]: if not v in path: new_distance = distance + W[(current,v)] if new_distance < best: new_path = copy(path) new_path.append(v) if v == last_v: # Found a path best = new_distance print new_path,'\tdistance:',best,'\tmax queue',max_q else: Q.enqueue((new_path,new_distance)) q_size += len(new_path) max_q = max(max_q,q_size)
Shortest_Path_BF(G2,WW,'A','B') timeit("Shortest_Path_BF(G2,WW,'A','B')")
Shortest_Path_BF(G2,WW,'A','C') timeit("Shortest_Path_BF(G2,WW,'A','C')")
Shortest_Path_BF(G2,WW,'A','D') timeit("Shortest_Path_BF(G2,WW,'A','D')")
Shortest_Path_BF(PMG,WWW,'A','B') timeit("Shortest_Path_BF(PMG,WWW,'A','B')")
Shortest_Path_BF(PMG,WWW,'A','R') timeit("Shortest_Path_BF(PMG,WWW,'A','R')")
Shortest_Path_BF(PMG,WWW,'A','O') timeit("Shortest_Path_BF(PMG,WWW,'A','O')")
Shortest_Path_BF(PMG,WWW,'A','W') timeit("Shortest_Path_BF(PMG,WWW,'A','W')")

Exercise 1:  In the block below discuss the advantages and disadvantages of shortest path algorithms based on depth-first search versus breadth-first search.

Dijkstra's Algorithm replaces the queue with a priority queue and builds a tree of the shortest paths from the source to all the other vertices in the network. Each node in the tree records the current parent of the vertex and the distance from the source to the vertex along the shortest path to the parent. As better paths are found the parent and the distance are updated. This necessitates that the priority queue keep track of the index of each object so that the priority of that object can be adjusted.

Once the tree is constructed, finding a specific path is extremely efficient.

def Dijkstra(G,W,source): # Dijkstra's algorithm actually produces a tree of the # shortest distances to all the vertices from the root # This implementation will stop when it reaches last_v # The entries in the tree are the parent and # the minimum distance from the root tree = {} tree[source] = (None,0) # Start with a gross overestimate of the distance too_big = 0 for w in W: too_big += W[w] PQ = PriorityQueue() for v in G.V: if not v == source: PQ.push(too_big,v) tree[v] = (source,too_big) u0 = source while PQ.is_not_empty(): d0 = tree[u0] for v in G.Adj[u0]: d1 = PQ.priority_of(v) #print 'Adj',u0,v,d1 if d1: d2 = d0 + W[(u0,v)] #print u0,d0,d1,d2 if d2 < d1: PQ.adjust(d2,v) tree[v] = (u0,d2) #print PQ._pq #print PQ._indx #print tree #print v,d1,d2 #print tree d1,u1 = PQ.pop() u0 = u1 return tree def Short_Path(tree, target): # This routine finds the shortest path from a pre-specified # source to the requested target path = [target] v,d0 = tree[target] while v: path.insert(0,v) v,d = tree[v] print path,d0
Tree = Dijkstra(G2,WW,'A') print Tree Short_Path(Tree,'B') Short_Path(Tree,'C') Short_Path(Tree,'D') timeit("Tree = Dijkstra(G2,WW,'A')") timeit("Short_Path(Tree,'B');Short_Path(Tree,'C');Short_Path(Tree,'D')")
Tree = Dijkstra(PMG,WWW,'A') print Tree Short_Path(Tree,'B') Short_Path(Tree,'R') Short_Path(Tree,'O') Short_Path(Tree,'W') timeit("Dijkstra(PMG,WWW,'A')") timeit("Short_Path(Tree,'B');Short_Path(Tree,'R');Short_Path(Tree,'O');Short_Path(Tree,'W')")

Exercise 2:  In the block below discuss the advantages and disadvantages of the overhead in computing Dijkstra's Shortest Path Tree versus the two previous Shortest Path algorithms.

Exercise 3:  In the block below discuss the similarities and differences between the Minimum Weighted Spanning Tree and Dijkstra's Shortest Path Tree.