From 41f917a3cbf5d774a9fdf9774e5ef2fd17b237ed Mon Sep 17 00:00:00 2001 From: nojhan Date: Thu, 20 Jan 2011 23:22:12 +0100 Subject: [PATCH] first commit --- lindenmayer.py | 196 ++++++++++++++++++++++++++++++++++++++ tsplib.py | 249 +++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 445 insertions(+) create mode 100644 lindenmayer.py create mode 100644 tsplib.py diff --git a/lindenmayer.py b/lindenmayer.py new file mode 100644 index 0000000..a210c41 --- /dev/null +++ b/lindenmayer.py @@ -0,0 +1,196 @@ +from collections import deque + +class IndexedGenerator(object): + """Add a way to get a generator item by its index""" + def __init__(self, generator): + self.generator = generator + self.cache = [] + + def __getitem__(self, index): + for i in xrange( index - len(self.cache) + 1 ): + self.cache.append( self.generator.next() ) + return self.cache[index] + + +class LindenmayerSystem(IndexedGenerator): + """Base virtual class for a Lindenmayer system""" + def __init__(self, axiom, rules, angle, heading=0): + self.angle = angle + self.heading = heading + self.states = deque() + self.actions = { + 'F': self.forward, + '+': self.right, + '-': self.left, + '[': self.save, + ']': self.restore, + } + super(LindenmayerSystem, self).__init__(self.lindenmayer(axiom, rules)) + + def lindenmayer(self, axiom, rules): + rules = rules.items() + + # def inside the lindenmayer function, so as to use "axiom" at instanciation + def apply(axiom, (symbol, replacement)): + return axiom.replace(symbol, replacement.lower()) + + while True: + yield axiom + axiom = reduce(apply, rules, axiom).upper() + + def forward(self): + raise NotImplementedError + + def right(self): + raise NotImplementedError + + def left(self): + raise NotImplementedError + + def save(self): + raise NotImplementedError + + def restore(self): + raise NotImplementedError + + +class TurtleLSystem(LindenmayerSystem): + """Draw a L-System using the Turtle module""" + def __init__(self, turtle, axiom, rules, angle, heading=0, size=1): + self.turtle = turtle + self.size = size + super(TurtleLSystem, self).__init__( axiom, rules, angle, heading ) + + def draw(self, depth): + self.turtle.setheading(self.heading) + + for char in self[depth]: + if char in self.actions: + self.actions[char]() + + def forward(self): + self.turtle.forward(self.size) + + def left(self): + self.turtle.left(self.angle) + + def right(self): + self.turtle.right(self.angle) + + def save(self): + x = self.turtle.xcor() + y = self.turtle.ycor() + h = self.turtle.heading() + self.states.append( (x, y, h) ) + + def restore(self): + turtle.up() + x, y, h = self.states.pop() + turtle.setx(x) + turtle.sety(y) + turtle.setheading(h) + turtle.down() + + +class DumpTurtleLSystem(TurtleLSystem): + """Keep the set of uniques L-System segments drawn by the Turtle""" + def __init__(self, turtle, axiom, rules, angle, heading=0, size=1, rounding=10): + # using a set avoid duplicate segments + self.segments = set() + # nb of significant digits for rounding + self.rounding=10 + super(DumpTurtleLSystem, self).__init__( turtle, axiom, rules, angle, heading, size ) + + def forward(self): + """Store segment coordinates and do a forward movement""" + # without rounding, there may be the same node with different coordinates, + # because of error propagation + x1 = round( self.turtle.xcor(), self.rounding ) + y1 = round( self.turtle.ycor(), self.rounding ) + start = ( x1, y1 ) + super(DumpTurtleLSystem, self).forward() + x2 = round( self.turtle.xcor(), self.rounding ) + y2 = round( self.turtle.ycor(), self.rounding ) + end = ( x2, y2 ) + self.segments.add( (start,end) ) + + def draw(self, depth): + """Call the draw function, then clean the data""" + super(DumpTurtleLSystem, self).draw(depth) + self.clean() + + def clean(self): + """Remove segments that have duplicated clones in the reverse direction + (the segments is a set, that guarantees that no other clone exists)""" + for segment in self.segments: + for start,end in segment: + # FIXME surely faster to catch the exception than to do two search? + if (end,start) in self.segments: + self.segments.remove( (end,start) ) + + def __str__(self): + dump = "" + for segment in self.segments: + for coords in segment: + for x in coords: + dump += str(x)+" " + dump += "\n" + return dump + + +def plot_segments( segments ): + + import matplotlib.pyplot as plot + from matplotlib.path import Path + import matplotlib.patches as patches + + fig = plot.figure() + ax = fig.add_subplot(111) + + for segment in segments: + start,end = segment + verts = [start,end,(0,0)] + codes = [Path.MOVETO,Path.LINETO,Path.STOP] + path = Path(verts, codes) + patch = patches.PathPatch(path, facecolor='none', lw=1) + ax.add_patch(patch) + + ax.set_xlim(-50,50) + ax.set_ylim(-50,50) + + plot.show() + + + +if __name__=="__main__": + import sys + + depth = 1 + if len(sys.argv) > 1: + depth = int( sys.argv[1] ) + + segment_size = 10 + float_rounding = 10 + + from turtle import Turtle + turtle = Turtle() + turtle.speed('fastest') + penrose = DumpTurtleLSystem(turtle, + axiom="[X]++[X]++[X]++[X]++[X]", + rules={ + 'F': "", + 'W': "YF++ZF----XF[-YF----WF]++", + 'X': "+YF--ZF[---WF--XF]+", + 'Y': "-WF++XF[+++YF++ZF]-", + 'Z': "--YF++++WF[+ZF++++XF]--XF" + }, + angle=36, heading=0, size=segment_size, rounding=float_rounding ) + + penrose.draw( depth ) + #print penrose + + #plot_segments( penrose.segments ) + + import tsplib + + tsplib.write_segments( penrose.segments, segment_size, depth, float_rounding, fd=sys.stdout ) diff --git a/tsplib.py b/tsplib.py new file mode 100644 index 0000000..47d2cc8 --- /dev/null +++ b/tsplib.py @@ -0,0 +1,249 @@ +import sys +import scipy + +def write_segments( segments, size, depth, rounding, fd = sys.stdout, + node_coord_section=False, + edge_data_section=False, + edge_weight_section=True, + display_data_section=True ): + + # construct a {coords:id} dictionary + nodes = {} + nb = 0 + for segment in segments: + for coords in segment: + if not nodes.has_key(coords): + nodes[coords] = nb + nb += 1 + + fd.write( "NAME : penrose3_%i\n" % depth) + fd.write("COMMENT : Rhombus Penrose tiling (type P3) as generated by a L-system, at depth %i\n" % depth) + fd.write("TYPE : TSP\n") + fd.write("DIMENSION : %i\n" % nb ) + if edge_weight_section: + fd.write("EDGE_WEIGHT_TYPE : EXPLICIT\n") + fd.write("EDGE_WEIGHT_FORMAT : FULL_MATRIX\n") + if edge_data_section: + fd.write("EDGE_DATA_FORMAT : ADJ_LIST\n") # via the weight matrix? + if node_coord_section: + fd.write("NODE_COORD_TYPE : TWOD_COORDS\n") # do not work with concord + if display_data_section: + fd.write("DISPLAY_DATA_TYPE : TWOD_DISPLAY\n") + + if node_coord_section: + fd.write("NODE_COORD_SECTION\n") + fmt = "%"+str(len(str(nb)))+"i %"+str(rounding)+"f %"+str(rounding)+"f\n" + for x,y in nodes: + fd.write(fmt % (nodes[(x,y)],x,y)) + + if edge_data_section: + fd.write("EDGE_DATA_SECTION\n") + for segment in segments: + start,end = segment + fd.write( str(nodes[start])+" "+str(nodes[end])+"\n" ) + + if edge_weight_section: + fd.write("EDGE_WEIGHT_SECTION\n") + # fill the weights matrix with size where necessary + weights = scipy.zeros((nb,nb), type(size)) + for segment in segments: + start,end = segment + weights[nodes[start],nodes[end]] = size + weights[nodes[end],nodes[start]] = size + #fd.write(nodes[start],nodes[end] + + fmt = "%"+str(len(str(size)))+"i " + for i in xrange(weights.shape[0]): + # full matrix + for j in xrange(weights.shape[1]): + fd.write(fmt % weights[i,j]) + fd.write('\n') + + if display_data_section: + fd.write("DISPLAY_DATA_SECTION\n") + fmt = "%"+str(len(str(nb)))+"i %"+str(rounding)+"f %"+str(rounding)+"f\n" + for x,y in nodes: + fd.write(fmt % (nodes[(x,y)],x,y)) + + fd.write("EOF\n") + + +def read_tour_index( fd ): + tour = [] + nb = int(fd.readline().strip()) + for line in fd: + tour += line.split() + + return map(int, tour) + + +def read_nodes( fd ): + """Parse a .tsp file and returns a dictionary of nodes, of the form {id:(x,y)}""" + nodes = {} + data_section = False + for line in fd: + if line.strip() == "DISPLAY_DATA_SECTION": + data_section = True + continue + + data = line.strip().split() + if len(data) != 3 or ( len(data) > 1 and data[0].isalpha() ): + data_section = False + + if data_section == True: + nodes[ int(data[0]) ] = ( float(data[1]),float(data[2]) ) + + return nodes + + +def read_vertices( fd ): + vertices = set() + + data_section = False + i=-1 + for line in fd: + if line.strip() == "EDGE_WEIGHT_SECTION": + data_section = True + i=0 + continue + + data = line.strip().split() + if len(data)==0 or ( len(data) >= 1 and data[0][0].isalpha() ): + data_section = False + + if data_section == True: + for j in xrange(len(data)): + if float(data[j]) != 0 and (j,i) not in vertices: + vertices.add( (i,j) ) + i += 1 + + return vertices + + +def plot_segments( segments ): + + import matplotlib.pyplot as plot + from matplotlib.path import Path + import matplotlib.patches as patches + + fig = plot.figure() + ax = fig.add_subplot(111) + + for segment in segments: + start,end = segment + verts = [start,end,(0,0)] + codes = [Path.MOVETO,Path.LINETO,Path.STOP] + path = Path(verts, codes) + patch = patches.PathPatch(path, facecolor='none', lw=1) + ax.add_patch(patch) + + ax.set_xlim(-50,50) + ax.set_ylim(-50,50) + + plot.show() + + +def plot_segments_tour( segments_1, segments_2 ): + + import matplotlib.pyplot as plot + from matplotlib.path import Path + import matplotlib.patches as patches + + fig = plot.figure() + ax = fig.add_subplot(111) + + for segment in segments_1: + start,end = segment + verts = [start,end,(0,0)] + codes = [Path.MOVETO,Path.LINETO,Path.STOP] + path = Path(verts, codes) + patch = patches.PathPatch(path, facecolor='0.5', lw=1) + ax.add_patch(patch) + + for segment in segments_2: + start,end = segment + verts = [start,end,(0,0)] + codes = [Path.MOVETO,Path.LINETO,Path.STOP] + path = Path(verts, codes) + patch = patches.PathPatch(path, facecolor='1.0', lw=3) + ax.add_patch(patch) + + + ax.set_xlim(-50,50) + ax.set_ylim(-50,50) + + plot.show() + + + + +if __name__=="__main__": + import sys + +# segments = [ +# ( (0,0),(0,2) ), +# ( (0,2),(2,2) ), +# ( (2,2),(2,0) ), +# ( (2,0),(0,0) ) +# ] +# +# filename = "test.tsp" +# with open(filename,"w") as fd: +# write_segments( segments, fd=fd, size=1, depth=0, rounding=10 ) +# write_segments( segments, fd=sys.stdout, size=1, depth=0, rounding=10 ) +# +# with open(filename,"r") as fd: +# nodes = read_nodes( fd ) +# +# print "Nodes: id (x, y)" +# for idx,node in nodes.items(): +# print idx,node +# +# with open(filename,"r") as fd: +# vertices = read_vertices( fd ) +# +# print "Segments: (x1,y1) (x2,y2)" +# segments = [] +# for i1,i2 in vertices: +# print nodes[i1],nodes[i2] +# segments.append( (nodes[i1],nodes[i2]) ) +# +# plot_segments( segments ) + + finstance = sys.argv[1] + ftour = sys.argv[2] + + print "Read nodes" + with open(finstance,"r") as fd: + nodes = read_nodes( fd ) + + print "Read vertices" + with open(finstance,"r") as fd: + vertices = read_vertices( fd ) + + print "Build segments" + segments = [] + for i1,i2 in vertices: + #print nodes[i1],nodes[i2] + segments.append( (nodes[i1],nodes[i2]) ) + +# print "Plot segments" +# plot_segments( segments ) + + + print "Read tour" + with open(ftour,"r") as fd: + tour = read_tour_index( fd ) + + print "Build tour segments" + tour_segments = [] + for i in xrange(0,len(tour)-1): + tour_segments.append( ( nodes[tour[i]],nodes[tour[i+1]] ) ) + print tour_segments[-1] + + tour_segments.append( ( nodes[tour[i+1]], nodes[tour[0]] ) ) + print tour_segments[-1] + + print "Plot tour segments" + plot_segments_tour( segments, tour_segments ) +