[pypy-issue] [issue731] Problem handling set operations
Jose Antonio Martin H.
tracker at bugs.pypy.org
Fri May 27 19:53:53 CEST 2011
Jose Antonio Martin H. <jamartinh at fdi.ucm.es> added the comment:
Hi Alex, try the attached script in python27
python GraphLib.py runs without errors
However
pypy GraphLib.py after about 10 iterations gives an error.
Thanks for your quick hands on!.
Best,
Jose
El 27/05/2011 19:28, Alex Gaynor escribió:
>
> Alex Gaynor<alex.gaynor at gmail.com> added the comment:
>
> Can we get a slightly more complete example? Preferable a test case we can
> actually execute.
>
> ----------
> nosy: +agaynor
> status: unread -> chatting
>
> ________________________________________
> PyPy bug tracker<tracker at bugs.pypy.org>
> <https://bugs.pypy.org/issue731>
> ________________________________________
--
======================================================================
José Antonio MartÃn H. PhD. E-Mail: jamartinh at fdi.ucm.es
Computer Science Faculty Phone: (+34) 91 3947650
Complutense University of Madrid Fax: (+34) 91 3947527
C/ Prof. José GarcÃa Santesmases,s/n http://www.dacya.ucm.es/jam/
28040 Madrid, Spain
El orden es el recurso no renovable Order is the truly nonrenewable
más importante resource.
======================================================================
________________________________________
PyPy bug tracker <tracker at bugs.pypy.org>
<https://bugs.pypy.org/issue731>
________________________________________
-------------- next part --------------
#-------------------------------------------------------------------------------
# Name: GraphLib.py
# Description: small graph library
#
# Author: Jose Antonio Martin H. (JAMH)
# Contact: <jamartinh at fdi.ucm.es>
#
# Created: 17/05/2011
# Copyright: (c) Jose Antonio Martin H. 2011
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
#
#-------------------------------------------------------------------------------
#!/usr/bin/env python
from copy import deepcopy
import random
from itertools import combinations,ifilter,ifilterfalse,tee,izip
from math import sqrt
class Stack:
""" Simple Stack class implemented as a Python list """
def __init__(self):
self.contents = []
def Push(self, v):
self.contents.append(v)
def Pop(self):
v = self.contents[-1]
self.contents = self.contents[:-1]
return v
def Clear(self):
self.contents = []
def IsEmpty(self):
return (len(self.contents) == 0)
def IsNotEmpty(self):
return (len(self.contents) > 0)
def Contains(self, v):
return v in self.contents
################################################################################
#
# Graph
#
################################################################################
class Graph:
def __init__(self):
self.vertices = list()
self.edges = set()
self.adjLists = {}
self.highVertexID = 0 # INTERNAL
def AddVertex(self):
""" Add an isolated vertex. Returns the id of the new vertex """
x = self.GetNextVertexID()
self.vertices.append(x)
self.adjLists[x] = set()
return x
def DeleteVertex(self, v):
""" Delete the vertex v and its incident edges """
for u in self.adjLists[v]:
self.adjLists[u].discard(v)
self.edges.discard( {u,v} )
del self.adjLists[v]
def QVertex(self, v):
""" Check whether v is a vertex """
return v in self.vertices
def AddEdge(self, x, y):
""" Add an edge (x,y). Returns nothing
"""
if x == y: raise KeyError
self.adjLists[x].add(y)
self.adjLists[y].add(x)
self.edges.add( frozenset((x,y)) )
def DeleteEdge(self, x, y):
""" Deletes edge (tail,head).
"""
self.edges.discard( {x,y} )
self.adjLists[x].discard(y)
self.adjLists[y].discard(x)
def QEdge(self, x, y):
""" Returns 1 if (x,yd) is an edge in G"""
return {x,y} in self.edges
def Neighborhood(self, x):
""" Returns the vertices which are connected to v. """
return self.adjLists[x]
def Edges(self):
return self.edges
def Vertices(self):
""" Returns all edges """
return self.vertices
def GetNextVertexID(self):
""" *Internal* returns next free vertex id """
self.highVertexID += 1
return self.highVertexID
def Order(self):
""" Returns order i.e., the number of vertices """
return len(self.vertices)
def Size(self):
""" Returns size i.e., the number of edge """
return len(self.edges)
def Degree(self, v):
""" Returns the degree of the vertex v """
return len(self.adjLists[v])
def QIsolated(self, v):
""" Returns 1 if the vertex v is isolated"""
return len(self.adjLists[v]) == 0
#-------------------------------------------------------------------------------
# GraphCreator
#-------------------------------------------------------------------------------
def MaximalPlanarEdges(G, n, direction):
Edges=[] #6*n
AdjEdges={}
for v in G.Vertices():
AdjEdges[v]=[]
index=0
a=G.vertices[index]
index+=1
b=G.vertices[index]
index+=1
Edges.append((a,b))
AdjEdges[a].append((a,b))
Edges.append((b,a))
AdjEdges[b].append((b,a))
m=2
while index < n:
e = Edges[random.randint(0,m-1)]
v = G.vertices[index]
index = index+1
while e[1]!=v:
x=(v,e[0])
Edges.append(x)
m=m+1
AdjEdges[v].append(x)
y=(e[0],v)
Edges.append(y)
m=m+1
AdjEdges[e[0]].insert(AdjEdges[e[0]].index(e)+1,y)
index2=AdjEdges[e[1]].index((e[1],e[0]))
if index2==0:
e = AdjEdges[e[1]][-1]
else:
e = AdjEdges[e[1]][index2-1]
if direction==0: # undirected
m=m-1
while m>0:
del Edges[m]
m=m-2
return Edges
#----------------------------------------------------------------------
#-------------------------------------------------------------------------------
# PlanarityTest
#-------------------------------------------------------------------------------
#=============================================================================#
class List:
def __init__(self,el=[]):
elc=deepcopy(el)
self.elements=elc
# a) Access Operations
def length(self):
return len(self.elements)
def empty(self):
if self.length()==0:
return 1
else:
return 0
def head(self):
return self.elements[0]
def tail(self):
return self.elements[-1]
# b)Update Operations
def push(self,x):
self.elements.insert(0,x)
return x
def Push(self,x):
self.elements.append(x)
return x
def append(self,x):
self.Push(x)
def pop(self):
x=self.elements[0]
self.elements=self.elements[1:]
return x
def Pop(self):
x=self.elements[-1]
self.elements=self.elements[:-1]
return x
def clear(self):
self.elements=[]
def conc(self,A):
self.elements=self.elements+A.elements
A.elements=[]
#=============================================================================#
#=============================================================================#
class pt_graph:
def __init__(self):
self.V = []
self.E = []
self.adjEdges = {}
# a) Access operations
def source(self,e):
return e[0]
def target(self,e):
return e[1]
def number_of_nodes(self):
return len(self.V)
def number_of_edges(self):
return len(self.E)
def all_nodes(self):
return self.V
def all_edges(self):
return self.E
def adj_edges(self, v):
return self.adjEdges[v]
def adj_nodes(self, v):
nodelist=[]
for e in self.adj_edges(v):
nodelist.append(e[1])
return nodelist
def first_node(self):
return self.V[0]
def last_node(self):
return self.V[-1]
def first_edge(self):
return self.E[0]
def last_edge(self):
return self.E[-1]
def first_adj_edge(self,v):
if len(self.adj_edges(v))>0:
return self.adj_edges(v)[0]
else:
return None
def last_adj_edge(self,v):
if len(self.adj_edges(v))>0:
return self.adj_edges(v)[-1]
else:
return None
# b) Update operations
def new_node(self, v):
self.V.append(v)
self.adjEdges[v]=[]
return v
def new_edge(self, v, w):
if v==w: # Loop
raise KeyError
if (v,w) in self.E: # Multiple edge
raise KeyError
self.E.append((v,w))
self.adjEdges[v].append((v,w))
return (v,w)
def del_node(self, v):
try:
for k in self.V:
for e in self.adj_edges(k):
if source(e)==v or target(e)==v:
self.adjEdges[k].remove(e)
self.V.remove(v)
for e in self.E:
if source(e)==v or target(e)==v:
self.E.remove(e)
except KeyError:
raise KeyError
def del_edge(self,e):
try:
self.E.remove(e)
self.adjEdges[source(e)].remove((source(e),target(e)))
except KeyError:
raise KeyError
def del_nodes(self, node_list): # deletes all nodes in list L from self
L=deepcopy(node_list)
for l in L:
self.del_node(l)
def del_edges(self, edge_list): # deletes all edges in list L from self
L=deepcopy(edge_list)
for l in L:
self.del_edge(l)
def del_all_nodes(self): # deletes all nodes from self
self.del_nodes(self.all_nodes())
def del_all_edges(self): # deletes all edges from self
self.del_edges(self.all_edges())
def sort_edges(self, cost):
sorted_list = cost.items()
sorted_list.sort(key = lambda x: x[1])
self.del_all_edges()
for i in sorted_list:
self.new_edge(source(i[0]),target(i[0]))
def source(e):
return e[0]
def target(e):
return e[1]
def reversal(e):
return (e[1],e[0])
#=============================================================================#
#=============================================================================#
class block:
# The constructor takes an edge and a list of attachments and creates
# a block having the edge as the only segment in its left side.
#
# |flip| interchanges the two sides of a block.
#
# |head_of_Latt| and |head_of_Ratt| return the first elements
# on |Latt| and |Ratt| respectively
# and |Latt_empty| and |Ratt_empty| check these lists for emptyness.
#
# |left_interlace| checks whether the block interlaces with the left
# side of the topmost block of stack |S|.
# |right_interlace| does the same for the right side.
#
# |combine| combines the block with another block |Bprime| by simply
# concatenating all lists.
#
# |clean| removes the attachment |w| from the block |B| (it is
# guaranteed to be the first attachment of |B|).
# If the block becomes empty then it records the placement of all
# segments in the block in the array |alpha| and returns true.
# Otherwise it returns false.
#
# |add_to_Att| first makes sure that the right side has no attachment
# above |w0| (by flipping); when |add_to_Att| is called at least one
# side has no attachment above |w0|.
# |add_to_Att| then adds the lists |Ratt| and |Latt| to the output list
# |Att| and records the placement of all segments in the block in |alpha|.
def __init__(self,e,A):
self.Latt=List(); self.Ratt=List() # list of attachments "ints"
self.Lseg=List(); self.Rseg=List() # list of segments represented by
# their defining "edges"
self.Lseg.append(e)
self.Latt.conc(A) # the other two lists are empty
def flip(self):
ha=List() # "ints"
he=List() # "edges"
# we first interchange |Latt| and |Ratt| and then |Lseg| and |Rseg|
ha.conc(self.Ratt); self.Ratt.conc(self.Latt); self.Latt.conc(ha);
he.conc(self.Rseg); self.Rseg.conc(self.Lseg); self.Lseg.conc(he);
def head_of_Latt(self):
return self.Latt.head()
def empty_Latt(self):
return self.Latt.empty()
def head_of_Ratt(self):
return self.Ratt.head()
def empty_Ratt(self):
return self.Ratt.empty()
def left_interlace(self,S):
# check for interlacing with the left side of the
# topmost block of |S|
if (S.IsNotEmpty() and not((S.contents[-1]).empty_Latt()) and
self.Latt.tail()<(S.contents[-1]).head_of_Latt()):
return 1
else:
return 0
def right_interlace(self,S):
# check for interlacing with the right side of the
# topmost block of |S|
if (S.IsNotEmpty() and not((S.contents[-1]).empty_Ratt()) and
self.Latt.tail()<(S.contents[-1]).head_of_Ratt()):
return 1
else:
return 0
def combine(self,Bprime):
# add block Bprime to the rear of |this| block
self.Latt.conc(Bprime.Latt)
self.Ratt.conc(Bprime.Ratt)
self.Lseg.conc(Bprime.Lseg)
self.Rseg.conc(Bprime.Rseg)
Bprime = None
def clean(self,dfsnum_w,alpha,dfsnum):
# remove all attachments to |w|; there may be several
while not(self.Latt.empty()) and self.Latt.head()==dfsnum_w:
self.Latt.pop()
while not(self.Ratt.empty()) and self.Ratt.head()==dfsnum_w:
self.Ratt.pop()
if not(self.Latt.empty()) or not(self.Ratt.empty()):
return 0
# |Latt| and |Ratt| are empty;
# we record the placement of the subsegments in |alpha|.
for e in self.Lseg.elements:
alpha[e]=left
for e in self.Rseg.elements:
alpha[e]=right
return 1
def add_to_Att(self,Att,dfsnum_w0,alpha,dfsnum):
# add the block to the rear of |Att|. Flip if necessary
if not(self.Ratt.empty()) and self.head_of_Ratt()>dfsnum_w0:
self.flip()
Att.conc(self.Latt)
Att.conc(self.Ratt)
# This needs some explanation.
# Note that |Ratt| is either empty or {w0}.
# Also if |Ratt| is non-empty then all subsequent
# sets are contained in {w0}.
# So we indeed compute an ordered set of attachments.
for e in self.Lseg.elements:
alpha[e]=left
for e in self.Rseg.elements:
alpha[e]=right
#=============================================================================#
#=============================================================================#
# GLOBALS:
left=1
right=2
G=pt_graph()
reached={}
dfsnum={}
parent={}
dfs_count=0
lowpt={}
Del=[]
lowpt1={}
lowpt2={}
alpha={}
Att=List()
cur_nr=0
sort_num={}
tree_edge_into={}
#=============================================================================#
#=============================================================================#
def planarity_test(Gin):
# planarity_test decides whether the InputGraph is planar.
# it also order the adjecentLists in counterclockwise.
n = Gin.Order() # number of nodes
if n<3: return 1
if Gin.Size()>3*n-6: return 0 # number of edges
if Gin.Size()>6*n-12: return 0
#--------------------------------------------------------------
# make G a copy of Gin and make G bidirected
global G,cur_nr
G = pt_graph()
for v in Gin.vertices:
G.new_node(v)
for e in Gin.Edges():
e = list(e)
G.new_edge(source(e),target(e))
cur_nr=0
nr={}
cost={}
n=G.number_of_nodes()
for v in G.all_nodes():
nr[v]=cur_nr
cur_nr=cur_nr+1
for e in G.all_edges():
if nr[source(e)] < nr[target(e)]:
cost[e]=n*nr[source(e)] + nr[target(e)]
else:
cost[e]=n*nr[target(e)] + nr[source(e)]
G.sort_edges(cost)
L=List(G.all_edges())
while not(L.empty()):
e=L.pop()
if (not(L.empty()) and source(e)==target(L.head())
and source(L.head())==target(e)):
L.pop()
else:
G.new_edge(target(e),source(e))
#--------------------------------------------------------------
#--------------------------------------------------------------
# make G biconnected
Make_biconnected_graph()
#--------------------------------------------------------------
#--------------------------------------------------------------
# make H a copy of G
#
# We need the biconnected version of G (G will be further modified
# during the planarity test) in order to construct the planar embedding.
# So we store it as a graph H.
H=deepcopy(G)
#--------------------------------------------------------------
#--------------------------------------------------------------
# test planarity
global dfsnum,parent,alpha,Att
dfsnum={}
parent={}
for v in G.all_nodes():
parent[v]=None
reorder()
alpha={}
for e in G.all_edges():
alpha[e]=0
Att=List()
alpha[G.first_adj_edge(G.first_node())] = left
if not(strongly_planar(G.first_adj_edge(G.first_node()),Att)):
return 0
#--------------------------------------------------------------
#--------------------------------------------------------------
# construct embedding
global sort_num,tree_edge_into
T=List()
A=List()
cur_nr=0
sort_num={}
tree_edge_into={}
embedding(G.first_adj_edge(G.first_node()),left,T,A)
# |T| contains all edges incident to the first node except the
# cycle edge into it.
# That edge comprises |A|.
T.conc(A)
for e in T.elements:
sort_num[e]=cur_nr
cur_nr=cur_nr+1
H.sort_edges(sort_num)
#--------------------------------------------------------------
return 1 if H.all_edges() else 0 # ccwOrderedEges
#=============================================================================#
#=============================================================================#
def pt_DFS(v):
global G,reached
S=Stack()
if reached[v]==0:
reached[v]=1
S.Push(v)
while S.IsNotEmpty():
v=S.Pop()
for w in G.adj_nodes(v):
if reached[w]==0:
reached[w]=1
S.Push(w)
#=============================================================================#
#=============================================================================#
def Make_biconnected_graph():
# We first make it connected by linking all roots of a DFS-forest.
# Assume now that G is connected.
# Let a be any articulation point and let u and v be neighbors
# of a belonging to different biconnected components.
# Then there are embeddings of the two components with the edges
# {u,a} and {v,a} on the boundary of the unbounded face.
# Hence we may add the edge {u,v} without destroying planarity.
# Proceeding in this way we make G biconnected.
global G,reached,dfsnum,parent,dfs_count,lowpt
#--------------------------------------------------------------
# We first make G connected by linking all roots of the DFS-forest.
reached={}
for v in G.all_nodes():
reached[v]=0
u=G.first_node()
for v in G.all_nodes():
if not(reached[v]):
# explore the connected component with root v
pt_DFS(v)
if u!=v:
# link v's component to the first component
G.new_edge(u,v)
G.new_edge(v,u)
#--------------------------------------------------------------
#--------------------------------------------------------------
# We next make G biconnected.
for v in G.all_nodes():
reached[v]=0
dfsnum={}
parent={}
for v in G.all_nodes():
parent[v]=None
dfs_count=0
lowpt={}
dfs_in_make_biconnected_graph(G.first_node())
#--------------------------------------------------------------
#=============================================================================#
#=============================================================================#
def dfs_in_make_biconnected_graph(v):
# This procedure determines articulation points and adds appropriate
# edges whenever it discovers one.
global G,reached,dfsnum,parent,dfs_count,lowpt
dfsnum[v]=dfs_count
dfs_count=dfs_count+1
lowpt[v]=dfsnum[v]
reached[v]=1
if not(G.first_adj_edge(v)): return # no children
u=target(G.first_adj_edge(v)) # first child
for e in G.adj_edges(v):
w=target(e)
if not(reached[w]):
# e is a tree edge
parent[w]=v
dfs_in_make_biconnected_graph(w)
if lowpt[w]==dfsnum[v]:
# v is an articulation point. We now add an edge.
# If w is the first child and v has a parent then we
# connect w and parent[v], if w is a first child and v
# has no parent then we do nothing.
# If w is not the first child then we connect w to the
# first child.
# The net effect of all of this is to link all children
# of an articulation point to the first child and the
# first child to the parent (if it exists).
if w==u and parent[v]:
G.new_edge(w,parent[v])
G.new_edge(parent[v],w)
if w!=u:
G.new_edge(u,w)
G.new_edge(w,u)
lowpt[v]=min(lowpt[v],lowpt[w])
else:
lowpt[v]=min(lowpt[v],dfsnum[w]) # non tree edge
#=============================================================================#
#=============================================================================#
def reorder():
# The procedure reorder first performs DFS to compute dfsnum, parent
# lowpt1 and lowpt2, and the list Del of all forward edges and all
# reversals of tree edges.
# It then deletes the edges in Del and finally reorders the edges.
global G,dfsnum,parent,reached,dfs_count,Del,lowpt1,lowpt2
reached={}
for v in G.all_nodes():
reached[v]=0
dfs_count = 0
Del=[]
lowpt1={}
lowpt2={}
dfs_in_reorder(G.first_node())
#--------------------------------------------------------------
# remove forward and reversals of tree edges
for e in Del:
G.del_edge(e)
#--------------------------------------------------------------
#--------------------------------------------------------------
# we now reorder adjacency lists
cost={}
for e in G.all_edges():
v = source(e)
w = target(e)
if dfsnum[w]<dfsnum[v]:
cost[e]=2*dfsnum[w]
elif lowpt2[w]>=dfsnum[v]:
cost[e]=2*lowpt1[w]
else:
cost[e]=2*lowpt1[w]+1
G.sort_edges(cost)
#--------------------------------------------------------------
#=============================================================================#
#=============================================================================#
def dfs_in_reorder( v):
global G,dfsnum,parent,reached,dfs_count,Del,lowpt1,lowpt2
#--------------------------------------------------------------
dfsnum[v]=dfs_count
dfs_count=dfs_count+1
lowpt1[v]=lowpt2[v]=dfsnum[v]
reached[v]=1
for e in G.adj_edges(v):
w = target(e);
if not(reached[w]):
# e is a tree edge
parent[w]=v
dfs_in_reorder(w)
lowpt1[v]=min(lowpt1[v],lowpt1[w])
else:
lowpt1[v]=min(lowpt1[v],dfsnum[w]) # no effect for forward edges
if dfsnum[w]>=dfsnum[v] or w==parent[v]:
# forward edge or reversal of tree edge
Del.append(e)
#--------------------------------------------------------------
#--------------------------------------------------------------
# we know |lowpt1[v]| at this point and now make a second pass over all
# adjacent edges of |v| to compute |lowpt2|
for e in G.adj_edges(v):
w = target(e)
if parent[w]==v:
# tree edge
if lowpt1[w]!=lowpt1[v]:
lowpt2[v]=min(lowpt2[v],lowpt1[w])
lowpt2[v]=min(lowpt2[v],lowpt2[w])
else:
# all other edges
if lowpt1[v]!=dfsnum[w]:
lowpt2[v]=min(lowpt2[v],dfsnum[w])
#--------------------------------------------------------------
#=============================================================================#
#=============================================================================#
def strongly_planar(e0,Att):
# We now come to the heart of the planarity test: procedure strongly_planar.
# It takes a tree edge e0=(x,y) and tests whether the segment S(e0) is
# strongly planar.
# If successful it returns (in Att) the ordered list of attachments of S(e0)
# (excluding x); high DFS-numbers are at the front of the list.
# In alpha it records the placement of the subsegments.
#
# strongly_planar operates in three phases.
# It first constructs the cycle C(e0) underlying the segment S(e0).
# It then constructs the interlacing graph for the segments emanating >from the
# spine of the cycle.
# If this graph is non-bipartite then the segment S(e0) is non-planar.
# If it is bipartite then the segment is planar.
# In this case the third phase checks whether the segment is strongly planar
# and, if so, computes its list of attachments.
global G,alpha,dfsnum,parent
#--------------------------------------------------------------
# DETERMINE THE CYCLE C(e0)
# We determine the cycle "C(e0)" by following first edges until a back
# edge is encountered.
# |wk| will be the last node on the tree path and |w0|
# is the destination of the back edge.
x=source(e0)
y=target(e0)
e=G.first_adj_edge(y)
wk=y
while dfsnum[target(e)]>dfsnum[wk]: # e is a tree edge
wk=target(e)
e=G.first_adj_edge(wk)
w0=target(e)
#--------------------------------------------------------------
#--------------------------------------------------------------
# PROCESS ALL EDGES LEAVING THE SPINE
# The second phase of |strongly_planar| constructs the connected
# components of the interlacing graph of the segments emananating
# from the the spine of the cycle "C(e0)".
# We call a connected component a "block".
# For each block we store the segments comprising its left and
# right side (lists |Lseg| and |Rseg| contain the edges defining
# these segments) and the ordered list of attachments of the segments
# in the block;
# lists |Latt| and |Ratt| contain the DFS-numbers of the attachments;
# high DFS-numbers are at the front of the list.
#
# We process the edges leaving the spine of "S(e0)" starting at
# node |wk| and working backwards.
# The interlacing graph of the segments emanating from
# the cycle is represented as a stack |S| of blocks.
w=wk
S=Stack()
while w!=x:
count=0
for e in G.adj_edges(w):
count=count+1
if count!=1: # no action for first edge
# TEST RECURSIVELY
# Let "e" be any edge leaving the spine.
# We need to test whether "S(e)" is strongly planar
# and if so compute its list |A| of attachments.
# If "e" is a tree edge we call our procedure recursively
# and if "e" is a back edge then "S(e)" is certainly strongly
# planar and |target(e)| is the only attachment.
# If we detect non-planarity we return false and free
# the storage allocated for the blocks of stack |S|.
A=List()
if dfsnum[w]<dfsnum[target(e)]:
# tree edge
if not(strongly_planar(e,A)):
while S.IsNotEmpty(): S.Pop()
return 0
else:
A.append(dfsnum[target(e)]) # a back edge
# UPDATE STACK |S| OF ATTACHMENTS
# The list |A| contains the ordered list of attachments
# of segment "S(e)".
# We create an new block consisting only of segment "S(e)"
# (in its L-part) and then combine this block with the
# topmost block of stack |S| as long as there is interlacing.
# We check for interlacing with the L-part.
# If there is interlacing then we flip the two sides of the
# topmost block.
# If there is still interlacing with the left side then the
# interlacing graph is non-bipartite and we declare the graph
# non-planar (and also free the storage allocated for the
# blocks).
# Otherwise we check for interlacing with the R-part.
# If there is interlacing then we combine |B| with the topmost
# block and repeat the process with the new topmost block.
# If there is no interlacing then we push block |B| onto |S|.
B=block(e,A)
while 1:
if B.left_interlace(S): (S.contents[-1]).flip()
if B.left_interlace(S):
B = None
while S.IsNotEmpty(): S.Pop()
return 0
if B.right_interlace(S): B.combine(S.Pop())
else: break
S.Push(B)
# PREPARE FOR NEXT ITERATION
# We have now processed all edges emanating from vertex |w|.
# Before starting to process edges emanating from vertex
# |parent[w]| we remove |parent[w]| from the list of attachments
# of the topmost
# block of stack |S|.
# If this block becomes empty then we pop it from the stack and
# record the placement for all segments in the block in array
# |alpha|.
while (S.IsNotEmpty() and
(S.contents[-1]).clean(dfsnum[parent[w]],alpha,dfsnum)):
S.Pop()
w=parent[w]
#--------------------------------------------------------------
#--------------------------------------------------------------
# TEST STRONG PLANARITY AND COMPUTE Att
# We test the strong planarity of the segment "S(e0)".
# We know at this point that the interlacing graph is bipartite.
# Also for each of its connected components the corresponding block
# on stack |S| contains the list of attachments below |x|.
# Let |B| be the topmost block of |S|.
# If both sides of |B| have an attachment above |w0| then
# "S(e0)" is not strongly planar.
# We free the storage allocated for the blocks and return false.
# Otherwise (cf. procedure |add_to_Att|) we first make sure that
# the right side of |B| attaches only to |w0| (if at all) and then
# add the two sides of |B| to the output list |Att|.
# We also record the placements of the subsegments in |alpha|.
Att.clear()
while S.IsNotEmpty():
B = S.Pop()
if (not(B.empty_Latt()) and not(B.empty_Ratt()) and
B.head_of_Latt()>dfsnum[w0] and B.head_of_Ratt()>dfsnum[w0]):
B = None
while S.IsNotEmpty(): S.Pop()
return 0
B.add_to_Att(Att,dfsnum[w0],alpha,dfsnum)
B = None
# Let's not forget that "w0" is an attachment
# of "S(e0)" except if w0 = x.
if w0!=x: Att.append(dfsnum[w0])
return 1
#--------------------------------------------------------------
#=============================================================================#
#=============================================================================#
def embedding(e0,t,T,A):
# embed: determine the cycle "C(e0)"
#
# We start by determining the spine cycle.
# This is precisley as in |strongly_planar|.
# We also record for the vertices w_r+1, ...,w_k, and w_0 the
# incoming cycle edge either in |tree_edge_into| or in the local
# variable |back_edge_into_w0|.
global G,dfsnum,cur_nr,sort_num,tree_edge_into,parent
x=source(e0)
y=target(e0)
tree_edge_into[y]=e0
e=G.first_adj_edge(y)
wk=y
while (dfsnum[target(e)]>dfsnum[wk]): # e is a tree edge
wk=target(e)
tree_edge_into[wk]=e
e=G.first_adj_edge(wk)
w0=target(e)
back_edge_into_w0=e
# process the subsegments
w=wk
Al=List()
Ar=List()
Tprime=List()
Aprime=List()
T.clear()
T.append(e) # |e=(wk,w0)| at this point
while w!=x:
count=0
for e in G.adj_edges(w):
count=count+1
if count!=1: # no action for first edge
# embed recursively
if dfsnum[w]<dfsnum[target(e)]:
# tree edge
if t==alpha[e]:
tprime=left
else:
tprime=right
embedding(e,tprime,Tprime,Aprime)
else:
# back edge
Tprime.append(e)
Aprime.append(reversal(e))
# update lists |T|, |Al|, and |Ar|
if t==alpha[e]:
Tprime.conc(T)
T.conc(Tprime) # T = Tprime conc T
Al.conc(Aprime) # Al = Al conc Aprime
else:
T.conc(Tprime) # T = T conc Tprime
Aprime.conc(Ar)
Ar.conc(Aprime) # Ar = Aprime conc Ar
# compute |w|'s adjacency list and prepare for next iteration
T.append(reversal(tree_edge_into[w])) # (w_j-1,w_j)
for e in T.elements:
sort_num[e]=cur_nr
cur_nr=cur_nr+1
# |w|'s adjacency list is now computed; we clear |T| and
# prepare for the next iteration by moving all darts incident
# to |parent[w]| from |Al| and |Ar| to |T|.
# These darts are at the rear end of |Al| and at the front end
# of |Ar|.
T.clear()
while not(Al.empty()) and source(Al.tail())==parent[w]:
# |parent[w]| is in |G|, |Al.tail| in |H|
T.push(Al.Pop()) # Pop removes from the rear
T.append(tree_edge_into[w]) # push would be equivalent
while not(Ar.empty()) and source(Ar.head())==parent[w]:
T.append(Ar.pop()); # pop removes from the front
w=parent[w]
# prepare the output
A.clear()
A.conc(Ar)
A.append(reversal(back_edge_into_w0))
A.conc(Al)
#=============================================================================#
#-------------------------------------------------------------------------------
# PlanarGraph
#-------------------------------------------------------------------------------
class PlanarGraph(Graph):
def __init__(self):
Graph.__init__(self)
self.identities = dict() #a dictionary of identified vertices
def AddVertex(self):
""" Add an isolated vertex. Returns the id of the new vertex """
x = Graph.AddVertex(self)
self.identities[x] = {x}
return x
def DeleteVertex(self, v):
Graph.DeleteVertex(self,v)
del self.identities[v]
def QPlanar(self,e=None):
""" test if graph is planar
"""
if e:
x,y = e
self.AddEdge(x,y)
Q = planarity_test(self)
self.DeleteEdge(x,y)
return Q
if not planarity_test(self): return 0
return 1
def QTriangleFree(self):
""" test if graph is triangle free
"""
for t in self.Triangles(): return 0
return 1
def Contract(self, x, y):
""" contract to vertices,e.g. identify it if there is no edge between x,y
or contract edge x,y if x,y is an edge of G.
This routine does not insert a new vertex, just copy the neighbors of y in x
and then deletes vertex y
This can be optimized by making x the vertex with higher degree
"""
for z in ( self.Neighborhood(y) - self.Neighborhood(x) - {x} ):
self.AddEdge(x,z)
self.identities[x] |= self.identities[y]
self.DeleteVertex(y)
def Subdivide(self, x, y):
"""
Subdivide edge xy and return new vertex z between x and y
"""
z = self.AddVertex()
self.AddEdge(x,z)
self.AddEdge(y,z)
self.DeleteEdge(x,y)
return z
def QComplete(self, v_list):
for x,y in combinations(v_list,2):
if not self.QEdge(x,y): return False
return True
def QRegular(self, r=0):
r = r if r else self.Degre(self.vertices[0])
for v in self.Vertices():
if self.Degree(v) != r: return False
return True
def Triangles(self):
for x,y in self.Edges():
for z in ( self.Neighborhood(x) & self.Neighborhood(y) ):
yield (x,y,z)
def FourClique(self):
for x,y,z in self.Triangles():
for w in ( self.Neighborhood(x) & self.Neighborhood(y) & self.Neighborhood(z) ):
return 1,[x,y,z,w]
return 0,[None,None,None,None]
def FindK112(self):
""" Finds a 3-partite complete subgraph $K_{1,1,2}$
"""
for x,y,z in self.Triangles():
for w in ( self.Neighborhood(x) & self.Neighborhood(y) - {z} ):
return 1,[x,y,z,w]
return 0,[]
def FindT31(self):
""" Finds a tadpole subgraph $T{3,1}$
"""
for x,y,z in self.Triangles():
for w in ( self.Neighborhood(z) - {x,y} ):
if self.QEdge(x,w) or self.QEdge(y,w): continue
yield (x,y,z,w)
def FindTC4(self):
""" Finds a square $C_4$
"""
for x,y in self.Edges():
for z in ( self.Neighborhood(y) - self.Neighborhood(x) - {x} ):
for w in ( (set(self.Neighborhood(x)) & set(self.Neighborhood(z))) - {y} ):
yield (x,y,z,w)
#-------------------------------------------------------------------------------------------------------------------------------------------
def Triangles(G):
for x,y in G.Edges():
for z in ( G.Neighborhood(x) & G.Neighborhood(y) ):
T.append([x,y,z])
return T
def FourClique(G):
for x,y,z in G.Triangles():
for w in ( G.Neighborhood(x) & G.Neighborhood(y) & G.Neighborhood(z) ):
return 1,[x,y,z,w]
return 0,[None,None,None,None]
def PlanarDensity( N, M):
return round(100.0*(M/(3.0*N-6.0)),2)
def EdgesFromDensity( N, D):
return round((D*(3.0*N-6.0))/100)
def FindPlanarPreservingEdge(G):
for x,y in combinations(G.Vertices(),2):
if G.QEdge(x,y): continue # assure not an edge
if G.QPlanar((x,y)): return (x,y) # assure planarity preserving operation
return () # the planar graph should be maximal
#-------------------------------------------------------------------------------
def randomPlanarGraphX( n, m):
G = PlanarGraph()
for v in range(n): G.AddVertex()
Edges = MaximalPlanarEdges(G,n,0)
for i in range(m):
pos = random.randint(0,len(Edges)-1)
G.AddEdge(Edges[pos][0],Edges[pos][1])
del Edges[pos]
return G
def randomPlanarG( N, M):
# ensure a pseudo-random planar graph with low probability of K4
G = randomPlanarGraphX(N,M)
nmax = 100
counter = 1
Q,(x,y,z,w) = G.FourClique()
while Q and counter < nmax:
G = randomPlanarGraphX(N,M)
counter += 1
Q,(x,y,z,w) = G.FourClique()
return G
#-------------------------------------------------------------------------------
# GraphLib
#-------------------------------------------------------------------------------
def verify_colorig( C, G):
for i in range(3):
for x,y in combinations(C[i],2):
if G.QEdge(x,y): return False
return True
def BT( n, i, v_list, C, G):
# 3^N brute force algorithm for 3 coloring with backtracking
if verify_colorig(C,G):
if i>n: return 1,C
else: return 0,C
for j in range(3):
x = deepcopy(C)
x[j].append(v_list[i])
Q,x = BT(n,i+1,v_list,x,G)
if Q: return 1,x
return 0,[[],[],[]]
def AlgBF(v_list, G):
n = len(v_list)
Q,(w,x,y,z) = FourClique(G)
if Q: return 0,[[],[],[]]
return BT(n-1,0,v_list,[[],[],[]],G)
def Is3ColorableBF(G):
"""
brute force backtracking coloring algorithm
"""
v_list = sorted( G.Vertices(), key = lambda v: G.Degree(v), reverse = True )
Q, C = AlgBF(v_list,G)
return Q,C,None
#-------------------------------------------------------------------------------
# 4-regular planar graphs
#-------------------------------------------------------------------------------
def random_permutation(iterable, r = 0):
"Random selection from itertools.permutations(iterable)"
pool = tuple(iterable)
r = (len(pool) if not r else r)
return tuple(random.sample(pool, r))
def pairwise(iterable):
"s -> (s0,s1), (s1,s2), (s2, s3), ..."
a, b = tee(iterable)
next(b, None)
return izip(a, b)
def Find_phi_A(G):
"""
two nonadjacent edges
"""
adbc = ifilterfalse(lambda ((a,d),(b,c)): a==b or a==c or b==d or d==c or
G.QEdge(a,b) or G.QEdge(c,d) or
not G.QPlanar((a,b)) or not G.QPlanar((d,c))
, combinations(G.Edges(),2) )
try: (a,d),(b,c) = random_permutation(adbc,1)
except ValueError: return ()
return (a,b,c,d)
def phi_A( G, CASE_A):
"""
N increases by 1
"""
a,b,c,d = CASE_A
G.DeleteEdge(a,d)
G.DeleteEdge(b,c)
x = G.AddVertex()
G.AddEdge(x,a)
G.AddEdge(x,b)
G.AddEdge(x,c)
G.AddEdge(x,d)
def Find_phi_B(G):
"""
for every triangular face
"""
try: a,b,u = random_permutation(G.Triangles(),1)
except ValueError: return ()
Nu = G.Neighborhood(u)
c,d = Nu-{a,b}
return (a,b,c,d,u)
def phi_B( G, CASE_B):
"""
N increases by 2
"""
v_list = []
a,b,c,d,u = CASE_B
G.DeleteEdge(a,b)
G.DeleteVertex(u)
for i in range(3):
v_list.append( G.AddVertex() )
x,y,x = v_list
G.AddEdge(x,y)
G.AddEdge(x,z)
G.AddEdge(y,z)
G.AddEdge(x,a)
G.AddEdge(x,b)
G.AddEdge(y,b)
G.AddEdge(y,c)
G.AddEdge(z,a)
G.AddEdge(z,d)
def Find_phi_C(G):
"""
for every vertex
"""
u = random.choice(G.Vertices())
a,b,c,d = G.Neighborhood(u)
return (a,b,c,d,u)
def phi_C( G, CASE_C):
"""
N increase by 4
"""
v_list = []
a,b,c,d,u = CASE_C
for i in (a,b,c,d):
v_list.append( G.Subdivide(u,i) )
v,w,x,y = v_list
c_f_permutations = ( (v,w,x,y,v), (v,w,y,x,v), (v,x,w,y,v) ) #free circular permutations
for v_tuple in c_f_permutations:
for i,j in pairwise(v_tuple):
G.AddEdge(i,j)
if G.QPlanar(): return
for i,j in pairwise(v_tuple): G.DeleteEdge(i,j)
def Find_phi_F(G):
"""
two square faces with a common edge
"""
for a,x in random.sample(G.Edges(),G.Size()):
for b in random_permutation( G.Neighborhood(a) - G.Neighborhood(x) - {x} ):
for y in random_permutation( (G.Neighborhood(b) & G.Neighborhood(x)) - G.Neighborhood(a) - {a} ):
if not G.QPlanar((a,y)) or not G.QPlanar((b,x)): continue
for c in random_permutation( G.Neighborhood(b) - G.Neighborhood(y) - {a,x,y} ):
for z in random_permutation( (G.Neighborhood(y) & G.Neighborhood(c)) - G.Neighborhood(b) - {a,x,b} ):
if not G.QPlanar((c,y)) or not G.QPlanar((b,z)): continue
return (a,b,c,x,y,z)
return ()
def phi_F( G, CASE_F):
"""
N increase by 2
"""
a,b,c,x,y,z = CASE_F
G.DeleteEdge(a,x)
G.DeleteEdge(c,z)
v = G.Subdivide(b,y)
w = G.Subdivide(v,y)
G.AddEdge(a,v)
G.AddEdge(c,v)
G.AddEdge(x,w)
G.AddEdge(z,w)
def Octahedron():
G = PlanarGraph()
v_list = []
for i in range(6):
v_list.append( G.AddVertex() )
v,w,x,y,z1,z2 = v_list
G.AddEdge(v,w)
G.AddEdge(w,x)
G.AddEdge(x,y)
G.AddEdge(y,v)
G.AddEdge(v,z1)
G.AddEdge(w,z1)
G.AddEdge(x,z1)
G.AddEdge(y,z1)
G.AddEdge(v,z2)
G.AddEdge(w,z2)
G.AddEdge(x,z2)
G.AddEdge(y,z2)
return G
def Random4RegularPlanarGraph( N):
G = Octahedron()
while( G.Order() < N ):
r = random.randint(0,3)
if r==0:
conf = Find_phi_A(G)
if not conf: continue
phi_A(G,conf)
## if not G.QRegular(4) :
## print "A"
## return 0
elif r==1:
conf = Find_phi_B(G)
if not conf: continue
phi_B(G,conf)
## if not G.QRegular(4) :
## print "B"
## return 0
elif r==2:
conf = Find_phi_C(G)
if not conf: continue
phi_C(G,conf)
## if not G.QRegular(4) :
## print "C"
## return 0
elif r==3:
conf = Find_phi_F(G)
if not conf: continue
phi_F(G,conf)
## if not G.QRegular(4) :
## print "F"
## return 0
if not G.QPlanar() or not G.QRegular(4): return False
return G
def main():
import time
for i in range(1000):
t1 =time.clock()
G = Random4RegularPlanarGraph(20)
t2 = time.clock()
if G:
print "Four regular planar graph of order ",G.Order(),'size',G.Size(),' in time',t2-t1
else:
print "not a planar or regular graph found"
break
if __name__ == '__main__':
main()
More information about the pypy-issue
mailing list