A Lindenmayer system or L-system is a mathematical system that can be used to describe growth process such as the growth of plants. Formally, it is a symbol expansion system whereby rewrite rules are applies iteratively to generate a longer string of symbols starting from a simple initial state. In this notebook, we will see how various types of fractal, including plant-like ones can be generated with L-systems and visualized with HoloViews.
In [ ]:
import holoviews as hv
import numpy as np
hv.extension('bokeh')
This notebook makes extensive use of the Path element and we will want to keep equal aspects and suppress the axes:
In [ ]:
%opts Path {+framewise +axiswise} [xaxis=None, yaxis=None show_title=False] (color='black')
In this notebook, we will be drawing paths relative to an agent, in the spirit of turtle graphics. For this we define a simple agent class that has a path property to show us the path travelled from the point of initialization:
In [ ]:
class SimpleAgent(object):
def __init__(self, x=0,y=0, heading=0):
self.x, self.y = x,y
self.heading = heading
self.trace = [(self.x, self.y)]
def forward(self, distance):
self.x += np.cos(2*np.pi * self.heading/360.0)
self.y += np.sin(2*np.pi * self.heading/360.0)
self.trace.append((self.x,self.y))
def rotate(self, angle):
self.heading += angle
def back(self, distance):
self.heading += 180
self.forward(distance)
self.heading += 180
@property
def path(self):
return hv.Path([self.trace])
We can now test our SimpleAgent by drawing some spirographs:
In [ ]:
def pattern(angle= 5):
agent = SimpleAgent()
for i in range(360//angle):
for i in range(4):
agent.forward(1)
agent.rotate(90)
agent.rotate(angle)
return agent
(pattern(20).path + pattern(10).path + pattern(5).path
+ pattern(5).path * pattern(10).path * pattern(20).path).cols(2)
We can also draw some pretty rose patterns, adapted from these equations:
In [ ]:
def roses(l,n,k):
agent = SimpleAgent()
n * 10
x = (2.0 * k -n) / (2.0 * n)
for i in range(360*n):
agent.forward(l)
agent.rotate(i + x)
return agent
roses(5, 7, 3).path + roses(5, 12, 5).path
We now want to the capabilites of our agent with the ability to read instructions, telling it which path to follow. Let's define the meaning of the following symbols:
F: Move forward by a pre-specified distance.
B: Move backwards by a pre-specified distance.
+: Rotate anti-clockwise by a pre-specified angle.
-: Rotate clockwise by a pre-specified angle.
Here is an agent class that can read strings of such symbols to draw the corresponding pattern:
In [ ]:
class Agent(SimpleAgent):
"An upgraded agent that can follow some rules"
default_rules = {'F': lambda t,d,a: t.forward(d),
'B': lambda t,d,a: t.back(d),
'+': lambda t,d,a: t.rotate(-a),
'-': lambda t,d,a: t.rotate(a)}
def __init__(self, x=0,y=0, instructions=None, heading=0,
distance=5, angle=60, rules=default_rules):
super(Agent,self).__init__(x,y, heading)
self.distance = distance
self.angle = angle
self.rules = rules
if instructions: self.process(instructions, self.distance, self.angle)
def process(self, instructions, distance, angle):
for i in instructions:
self.rules[i](self, distance, angle)
L-systems are defined with a rewrite system, making use of a set of production rules). What this means is that L-systems can generate instructions for our agent to follow, and therefore generate paths.
Now we define the expand_rules function which can process some expansion rules to repeatedly substitute an initial set of symbols with new symbols:
In [ ]:
def expand_rules(initial, iterations, productions):
"Expand an initial symbol with the given production rules"
expansion = initial
for i in range(iterations):
intermediate = ""
for ch in expansion:
intermediate = intermediate + productions.get(ch,ch)
expansion = intermediate
return expansion
To demonstrate expand_rules, let's define two different rules:
In [ ]:
koch_curve = {'F':'F+F-F-F+F'} # Replace 'F' with 'F+F-F-F+F'
koch_snowflake = {'F':'F-F++F-F'} # Replace 'F' with 'F-F++F-F'
Here are the first three steps using the first rule:
In [ ]:
for i in range(3):
print('%d: %s' % (i, expand_rules('F', i, koch_curve)))
Note that these are instructions our agent can follow!
In [ ]:
%%opts Path {+axiswise} (color=Cycle())
k1 = Agent(-200, 0, expand_rules('F', 4, koch_curve), angle=90).path
k2 = Agent(-200, 0, expand_rules('F', 4, koch_snowflake)).path
k1 + k2 + (k1 * k2)
This shows two variants of the Koch snowflake where koch_curve is a variant that uses right angles.
The following example introduces a mutual relationship between two symbols, 'A' and 'B', instead of just the single symbol 'F' used above:
In [ ]:
sierpinski_triangle = {'A':'B-A-B', 'B':'A+B+A'}
for i in range(3):
print('%d: %s' % (i, expand_rules('A', i,sierpinski_triangle)))
Once again we can use these instructions to draw an interesting shape although we also need to define what these symbols mean to our agent:
In [ ]:
%%opts Path (color='green')
sierpinski_rules = {'A': lambda t,d,a: t.forward(d),
'B': lambda t,d,a: t.forward(d),
'+': lambda t,d,a: t.rotate(-a),
'-': lambda t,d,a: t.rotate(a)}
instructions = expand_rules('A', 9,sierpinski_triangle)
Agent(x=-200, y=0, rules=sierpinski_rules, instructions=instructions, angle=60).path
We see that with our L-system expansion in terms of 'A' and 'B', we have defined the famous Sierpinski_triangle fractal.
Now for another famous fractal:
In [ ]:
dragon_curve = {'X':'X+YF+', 'Y':'-FX-Y'}
We now have two new symbols 'X' and 'Y' which we need to define in addition to 'F', '+' and '-' which we used before:
In [ ]:
dragon_rules = dict(Agent.default_rules, X=lambda t,d,a: None, Y=lambda t,d,a: None)
Note that 'X' and 'Y' don't actual do anything directly! These symbols are important in the expansion process but have no meaning to the agent. This time, let's use a HoloMap to view the expansion:
In [ ]:
%%opts Path {+framewise}
def pad_extents(path):
"Add 5% padding around the path"
minx, maxx = path.range('x')
miny, maxy = path.range('y')
xpadding = ((maxx-minx) * 0.1)/2
ypadding = ((maxy-miny) * 0.1)/2
path.extents = (minx-xpadding, miny-ypadding, maxx+xpadding, maxy+ypadding)
return path
hmap = hv.HoloMap(kdims='Iteration')
for i in range(7,17):
path = Agent(-200, 0, expand_rules('FX', i, dragon_curve), rules=dragon_rules, angle=90).path
hmap[i] = pad_extents(path)
hmap
This fractal is known as the Dragon Curve.
We have seen how to generate various fractals with L-systems, but we have not yet seen the plant-like fractals that L-systems are most famous for. This is because we can't draw a realistic plant with a single unbroken line: we need to be able to draw some part of the plant then jump back to an earlier state.
This can be achieved by adding two new actions to our agent: push to record the current state of the agent and pop to pop back to the state of the last push:
In [ ]:
class AgentWithState(Agent):
"Stateful agent that can follow instructions"
def __init__(self, x,y, instructions, **kwargs):
super(AgentWithState, self).__init__(x=x,y=y, instructions=None, **kwargs)
self.traces = []
self.state = []
self.process(instructions, self.distance, self.angle)
def push(self):
self.traces.append(self.trace[:])
self.state.append((self.heading, self.x, self.y))
def pop(self):
self.traces.append(self.trace[:])
[self.heading, self.x, self.y] = self.state.pop()
self.trace = [(self.x, self.y)]
@property
def path(self):
traces = self.traces + [self.trace]
return hv.Path(traces)
Let's look at the first three expansions of a new ruleset we will use to generate a plant-like fractal:
In [ ]:
plant_fractal = {'X':'F-[[X]+X]+F[+FX]-X', 'F':'FF'}
for i in range(3):
print('%d: %s' % (i, expand_rules('X', i, plant_fractal)))
The new symbols '[' and ']' correspond to the new push and pop state actions:
In [ ]:
plant_rules = dict(Agent.default_rules, X=lambda t,d,a: None,
**{'[': lambda t,d,a: t.push(), ']': lambda t,d,a: t.pop()})
We can now generate a nice plant-like fractal:
In [ ]:
%%opts Path {+framewise} (color='g' line_width=1)
hmap = hv.HoloMap(kdims='Iteration')
for i in range(7):
instructions = expand_rules('X', i, plant_fractal)
if i > 2:
hmap[i] = AgentWithState(-200, 0, instructions, heading=90, rules=plant_rules, angle=25).path
hmap