Chapter 1 Graphical modeling using Lsystems Lindenmayer systems — or Lsystems for short — were conceived as a mathematical theory of plant development [82]. Originally, they did not include enough detail to allow for comprehensive modeling of higher plants. The emphasis was on plant topology, that is, the neighborhood relations between cells or larger plant modules. Their geometric aspects were beyond the scope of the theory. Subsequently, several geometric interpretations of Lsystems were proposed with a view to turning them into a versatile tool for plant modeling. Throughout this book, an interpretation based on turtle geometry is used [109]. Basic notions related to Lsystem theory and their turtle interpretation are presented below.
1.1
Rewriting systems
The central concept of Lsystems is that of rewriting. In general, rewriting is a technique for deﬁning complex objects by successively replacing parts of a simple initial object using a set of rewriting rules or productions. The classic example of a graphical object deﬁned in terms of rewriting rules is the snowﬂake curve (Figure 1.1), proposed in 1905 by von Koch [155]. Mandelbrot [95, page 39] restates this construction as follows: One begins with two shapes, an initiator and a generator. The latter is an oriented broken line made up of N equal sides of length r. Thus each stage of the construction begins with a broken line and consists in replacing each straight interval with a copy of the generator, reduced and displaced so as to have the same end points as those of the interval being replaced.
Koch construction
2
Chapter 1. Graphical modeling using Lsystems
initiator
generator Figure 1.1: Construction of the snowﬂake curve
Grammars
Lsystems
While the Koch construction recursively replaces open polygons, rewriting systems that operate on other objects have also been investigated. For example, Wolfram [160, 161] studied patterns generated by rewriting elements of rectangular arrays. A similar arrayrewriting mechanism is the cornerstone of Conway’s popular game of life [49, 50]. An important body of research has been devoted to various graphrewriting systems [14, 33, 34]. The most extensively studied and the best understood rewriting systems operate on character strings. The ﬁrst formal deﬁnition of such a system was given at the beginning of this century by Thue [128], but a wide interest in string rewriting was spawned in the late 1950s by Chomsky’s work on formal grammars [13]. He applied the concept of rewriting to describe the syntactic features of natural languages. A few years later Backus and Naur introduced a rewritingbased notation in order to provide a formal deﬁnition of the programming language ALGOL60 [5, 103]. The equivalence of the BackusNaur form (BNF) and the contextfree class of Chomsky grammars was soon recognized [52], and a period of fascination with syntax, grammars and their application to computer science began. At the center of attention were sets of strings — called formal languages — and the methods for generating, recognizing and transforming them. In 1968 a biologist, Aristid Lindenmayer, introduced a new type of stringrewriting mechanism, subsequently termed Lsystems [82]. The essential diﬀerence between Chomsky grammars and Lsystems lies in
1.2. DOLsystems
3
Figure 1.2: Relations between Chomsky classes of languages and language classes generated by Lsystems. The symbols OL and IL denote language classes generated by contextfree and contextsensitive Lsystems, respectively.
the method of applying productions. In Chomsky grammars productions are applied sequentially, whereas in Lsystems they are applied in parallel and simultaneously replace all letters in a given word. This diﬀerence reﬂects the biological motivation of Lsystems. Productions are intended to capture cell divisions in multicellular organisms, where many divisions may occur at the same time. Parallel production application has an essential impact on the formal properties of rewriting systems. For example, there are languages which can be generated by contextfree Lsystems (called OLsystems) but not by contextfree Chomsky grammars [62, 128] (Figure 1.2).
1.2
DOLsystems
This section presents the simplest class of Lsystems, those which are deterministic and contextfree, called DOLsystems. The discussion starts with an example that introduces the main idea in intuitive terms. Consider strings (words) built of two letters a and b, which may occur many times in a string. Each letter is associated with a rewriting rule. The rule a → ab means that the letter a is to be replaced by the string ab, and the rule b → a means that the letter b is to be replaced by a. The rewriting process starts from a distinguished string called the axiom. Assume that it consists of a single letter b. In the ﬁrst derivation step (the ﬁrst step of rewriting) the axiom b is replaced
Example
4
Chapter 1. Graphical modeling using Lsystems
Figure 1.3: Example of a derivation in a DOLsystem
Lsystem
Derivation
by a using production b → a. In the second step a is replaced by ab using production a → ab. The word ab consists of two letters, both of which are simultaneously replaced in the next derivation step. Thus, a is replaced by ab, b is replaced by a, and the string aba results. In a similar way, the string aba yields abaab which in turn yields abaababa, then abaababaabaab, and so on (Figure 1.3). Formal deﬁnitions describing DOLsystems and their operation are given below. For more details see [62, 127]. Let V denote an alphabet, V ∗ the set of all words over V , and V + the set of all nonempty words over V . A string OLsystem is an ordered triplet G = V, ω, P where V is the alphabet of the system, ω ∈ V + is a nonempty word called the axiom and P ⊂ V × V ∗ is a ﬁnite set of productions. A production (a, χ) ∈ P is written as a → χ. The letter a and the word χ are called the predecessor and the successor of this production, respectively. It is assumed that for any letter a ∈ V , there is at least one word χ ∈ V ∗ such that a → χ. If no production is explicitly speciﬁed for a given predecessor a ∈ V , the identity production a → a is assumed to belong to the set of productions P . An OLsystem is deterministic (noted DOLsystem) if and only if for each a ∈ V there is exactly one χ ∈ V ∗ such that a → χ. Let µ = a1 . . . am be an arbitrary word over V . The word ν = χ1 . . . χm ∈ V ∗ is directly derived from (or generated by) µ, noted µ ⇒ ν, if and only if ai → χi for all i = 1, . . . , m. A word ν is generated by G in a derivation of length n if there exists a developmental sequence of words µ0 , µ1 , . . . , µn such that µ0 = ω, µn = ν and µ0 ⇒ µ1 ⇒ . . . ⇒ µn .
1.2. DOLsystems
5
Figure 1.4: Development of a ﬁlament (Anabaena catenula) simulated using a DOLsystem
The following example provides another illustration of the operation of DOLsystems. The formalism is used to simulate the development of a fragment of a multicellular ﬁlament such as that found in the bluegreen bacteria Anabaena catenula and various algae [25, 84, 99]. The symbols a and b represent cytological states of the cells (their size and readiness to divide). The subscripts l and r indicate cell polarity, specifying the positions in which daughter cells of type a and b will be produced. The development is described by the following Lsystem: ω p1 p2 p3 p4
: : : : :
ar ar → al b r al → b l ar b r → ar b l → al
(1.1)
Starting from a single cell ar (the axiom), the following sequence of words is generated: ar al b r b l ar ar al al b r al b r b l ar b l ar ar b l ar ar ···
Anabaena
6
Chapter 1. Graphical modeling using Lsystems
Under a microscope, the ﬁlaments appear as a sequence of cylinders of various lengths, with atype cells longer than btype cells. The corresponding schematic image of ﬁlament development is shown in Figure 1.4. Note that due to the discrete nature of Lsystems, the continuous growth of cells between subdivisions is not captured by this model.
1.3
Previous methods
Turtle
Turtle interpretation of strings
The geometric interpretation of strings applied to generate schematic images of Anabaena catenula is a very simple one. Letters of the Lsystem alphabet are represented graphically as shorter or longer rectangles with rounded corners. The generated structures are onedimensional chains of rectangles, reﬂecting the sequence of symbols in the corresponding strings. In order to model higher plants, a more sophisticated graphical interpretation of Lsystems is needed. The ﬁrst results in this direction were published in 1974 by Frijters and Lindenmayer [46], and Hogeweg and Hesper [64]. In both cases, Lsystems were used primarily to determine the branching topology of the modeled plants. The geometric aspects, such as the lengths of line segments and the angle values, were added in a postprocessing phase. The results of Hogeweg and Hesper were subsequently extended by Smith [136, 137], who demonstrated the potential of Lsystems for realistic image synthesis. Szilard and Quinton [141] proposed a diﬀerent approach to Lsystem interpretation in 1979. They concentrated on image representations with rigorously deﬁned geometry, such as chain coding [43], and showed that strikingly simple DOLsystems could generate the intriguing, convoluted curves known today as fractals [95]. These results were subsequently extended in several directions. Siromoney and Subramanian [135] speciﬁed Lsystems which generate classic spaceﬁlling curves. Dekking investigated the limit properties of curves generated by Lsystems [32] and concentrated on the problem of determining the fractal (Hausdorﬀ) dimension of the limit set [31]. Prusinkiewicz focused on an interpretation based on a LOGOstyle turtle [1] and presented more examples of fractals and plantlike structures modeled using Lsystems [109, 111]. Further applications of Lsystems with turtle interpretation include realistic modeling of herbaceous plants [117], description of kolam patterns (an art form from Southern India) [112, 115, 133, 134], synthesis of musical scores [110] and automatic generation of spaceﬁlling curves [116]. The basic idea of turtle interpretation is given below. A state of the turtle is deﬁned as a triplet (x, y, α), where the Cartesian coordinates (x, y) represent the turtle’s position, and the angle α, called the heading, is interpreted as the direction in which the turtle is facing. Given the step size d and the angle increment δ, the turtle can respond to
1.3. Turtle interpretation of strings
7
Figure 1.5: (a) Turtle interpretation of string symbols F , +, −. (b) Interpretation of a string. The angle increment δ is equal to 90◦ . Initially the turtle faces up.
commands represented by the following symbols (Figure 1.5a): F
Move forward a step of length d. The state of the turtle changes to (x , y , α), where x = x + d cos α and y = y + d sin α. A line segment between points (x, y) and (x , y ) is drawn.
f
Move forward a step of length d without drawing a line.
+
Turn left by angle δ. The next state of the turtle is (x, y, α +δ). The positive orientation of angles is counterclockwise.
−
Turn right by angle δ. The next state of the turtle is (x, y, α − δ).
Given a string ν, the initial state of the turtle (x0 , y0 , α0 ) and ﬁxed parameters d and δ, the turtle interpretation of ν is the ﬁgure (set of lines) drawn by the turtle in response to the string ν (Figure 1.5b). Speciﬁcally, this method can be applied to interpret strings which are generated by Lsystems. For example, Figure 1.6 presents four approximations of the quadratic Koch island taken from Mandelbrot’s book [95, page 51]. These ﬁgures were obtained by interpreting strings generated by the following Lsystem: ω : F −F −F −F p : F → F − F + F + FF − F − F + F The images correspond to the strings obtained in derivations of length 0 to 3. The angle increment δ is equal to 90◦ . The step length d is decreased four times between subsequent images, making the distance
Interpretation
8
Chapter 1. Graphical modeling using Lsystems
Figure 1.6: Generating a quadratic Koch island
between the endpoints of the successor polygon equal to the length of the predecessor segment. The above example reveals a close relationship between Koch constructions and Lsystems. The initiator corresponds to the axiom and the generator corresponds to the production successor. The predecessor F represents a single edge. Lsystems speciﬁed in this way can be perceived as codings for Koch constructions. Figure 1.7 presents further examples of Koch curves generated using Lsystems. A slight complication occurs if the curve is not connected; a second production (with the predecessor f ) is then required to keep components the proper distance from each other (Figure 1.8). The ease of modifying Lsystems makes them suitable for developing new Koch curves. For example, one can start from a particular Lsystem and observe the results of inserting, deleting or replacing some symbols. A variety of curves obtained this way are shown in Figure 1.9.
Koch constructions vs. Lsystems
1.3. Turtle interpretation of strings
9
Figure 1.7: Examples of Koch curves generated using Lsystems: (a) Quadratic Koch island [95, page 52], (b) A quadratic modiﬁcation of the snowﬂake curve [95, page 139]
Figure 1.8: Combination of islands and lakes [95, page 121]
10
Chapter 1. Graphical modeling using Lsystems
Figure 1.9: A sequence of Koch curves obtained by successive modiﬁcation of the production successor
1.4. Synthesis of DOLsystems
a n=10, δ=90◦ Fl Fl →Fl +Fr + Fr →Fl Fr
11
b n=6, δ=60◦ Fr Fl →Fr +Fl +Fr Fr →Fl Fr Fl
Figure 1.10: Examples of curves generated by edgerewriting Lsystems: (a) the dragon curve [48], (b) the Sierpi´ nski gasket [132]
1.4
Synthesis of DOLsystems
Random modiﬁcation of productions gives little insight into the relationship between Lsystems and the ﬁgures they generate. However, we often wish to construct an Lsystem which captures a given structure or sequence of structures representing a developmental process. This is called the inference problem in the theory of Lsystems. Although some algorithms for solving it were reported in the literature [79, 88, 89], they are still too limited to be of practical value in the modeling of higher plants. Consequently, the methods introduced below are more intuitive in nature. They exploit two modes of operation for Lsystems with turtle interpretation, called edge rewriting and node rewriting using terminology borrowed from graph grammars [56, 57, 87]. In the case of edge rewriting, productions substitute ﬁgures for polygon edges, while in node rewriting, productions operate on polygon vertices. Both approaches rely on capturing the recursive structure of ﬁgures and relating it to a tiling of a plane. Although the concepts are illustrated using abstract curves, they apply to branching structures found in plants as well.
1.4.1
Edge rewriting
Edge rewriting can be viewed as an extension of Koch constructions. For example, Figure 1.10a shows the dragon curve [21, 48, 95] and the Lsystem that generated it. Both the Fl and Fr symbols represent edges created by the turtle executing the “move forward” command. The productions substitute Fl or Fr edges by pairs of lines forming
12
Chapter 1. Graphical modeling using Lsystems
b a
n=4, δ=60◦ Fl Fl →Fl +Fr ++Fr Fl Fl Fl Fr + Fr →Fl +Fr Fr ++Fr +Fl Fl Fr
n=2, δ=90◦ Fr Fl →Fl Fl Fr Fr +Fl +Fl Fr Fr Fl + Fr +Fl Fl Fr Fl +Fr +Fl Fl + Fr Fl Fr Fr Fl +Fl +Fr Fr Fr →+Fl Fl Fr Fr +Fl +Fl Fr +Fl Fr Fr Fl Fr +Fl Fr Fr Fl Fr Fl +Fl +Fr Fr Fl +Fl +Fr Fr
Figure 1.11: Examples of FASS curves generated by edgerewriting Lsystems: (a) hexagonal Gosper curve [51], (b) quadratic Gosper curve [32] or Ecurve [96]
FASS curve construction
left or right turns. Many interesting curves can be obtained assuming two types of edges, “left” and “right.” Figures 1.10b and 1.11 present additional examples. The curves included in Figure 1.11 belong to the class of FASS curves (an acronym for spacefilling, selfavoiding, simple and selfsimilar) [116], which can be thought of as ﬁnite, selfavoiding approximations of curves that pass through all points of a square (spaceﬁlling curves [106]). McKenna [96] presented an algorithm for constructing FASS curves using edge replacement. It exploits the relationship between such a curve and a recursive subdivision of a square into tiles. For example, Figure 1.12 shows the tiling that corresponds to the Ecurve of Figure 1.11b. The polygon replacing an edge Fl (Figure 1.12a) approximately ﬁlls the square on the left side of Fl (b). Similarly, the polygon replacing an edge Fr (c) approximately ﬁlls the square on the right side of that edge (d). Consequently, in the next derivation step, each of the 25 tiles associated with the curves (b) or (d) will be covered by their reduced copies (Figure 1.11b). A recursive application of this argument indicates that the whole curve is approximately spaceﬁlling. It is also selfavoiding due to the following two properties:
1.4. Synthesis of DOLsystems
13
Fl →Fl Fl +Fr +Fr Fl Fl +Fr +Fr Fl Fr Fl Fl Fr + Fl Fr Fl Fl Fr +Fl Fr +Fr +Fl Fl Fr Fr + Fr →Fl Fl +Fr +Fr Fl Fl Fr Fl +Fr Fr +Fl +Fr Fl Fr Fr +Fl +Fr Fl Fl Fr +Fr +Fl Fl Fr Fr Figure 1.12: Construction of the Ecurve on the square grid. Left and right edges are distinguished by the direction of ticks.
• the generating polygon is selfavoiding, and • no matter what the relative orientation of the polygons lying on two adjacent tiles, their union is a selfavoiding curve. The ﬁrst property is obvious, while the second can be veriﬁed by considering all possible relative positions of a pair of adjacent tiles. Using a computer program to search the space of generating polygons, McKenna found that the Ecurve is the simplest FASS curve obtained by edge replacement in a square grid. Other curves require generators with more edges (Figure 1.13). The relationship between edge rewriting and tiling of the plane extends to branching structures, providing a method for constructing and analyzing Lsystems which operate according to the edgerewriting paradigm (see Section 1.10.3).
1.4.2
Node rewriting
The idea of node rewriting is to substitute new polygons for nodes of the predecessor curve. In order to make this possible, turtle interpretation is extended by symbols which represent arbitrary subﬁgures. As shown in Figure 1.14, each subﬁgure A from a set of subﬁgures A is represented by: • two contact points, called the entry point PA and the exit point QA , and • two direction vectors, called the entry vector pA and the exit vector qA . During turtle interpretation of a string ν, a symbol A ∈ A incorporates the corresponding subﬁgure into a picture. To this end, A is translated
Subﬁgures
14
Chapter 1. Graphical modeling using Lsystems
Figure 1.13: Examples of FASS curves generated on the square grid using edge replacement: (a) a SquaRecurve (grid size 7 × 7), (b) an Etour (grid size 9 × 9). Both curves are from [96].
Figure 1.14: Description of a subﬁgure A
1.4. Synthesis of DOLsystems
Ln
L n+1
L n+2
15
Rn
R n+1
R n+2
Figure 1.15: Recursive construction of the Hilbert curve [63] in terms of node replacement
and rotated in order to align its entry point PA and direction pA with the current position and orientation of the turtle. Having placed A, the turtle is assigned the resulting exit position QA and direction qA . For example, assuming that the contact points and directions of subﬁgures Ln and Rn are as in Figure 1.15, the ﬁgures Ln+1 and Rn+1 are captured by the following formulas: Ln+1 = +Rn F − Ln F Ln − F Rn + Rn+1 = −Ln F + Rn F Rn + F Ln − Suppose that curves L0 and R0 are given. One way of evaluating the string Ln (or Rn ) for n > 0 is to generate successive strings recursively, in the order of decreasing value of index n. For example, the computation of L2 would proceed as follows: L2 = +R1 F − L1 F L1 − F R1 + = +(−L0 F + R0 F R0 + F L0 −)F − (+R0 F − L0 F L0 − F R0 +) F (+R0 F − L0 F L0 − F R0 +) − F (−L0 F + R0 F R0 + F L0 −)+
Recursive formulas
16
Chapter 1. Graphical modeling using Lsystems
Thus, the generation of string Ln can be considered as a stringrewriting mechanism, where the symbols on the left side of the recursive formulas are substituted by corresponding strings on the right side. The substitution proceeds in a parallel way with F, + and − replacing themselves. Since all indices in any given string have the same value, they can be dropped, provided that a global count of derivation steps is kept. Consequently, string Ln can be obtained in a derivation of length n using the following Lsystem: ω : L p1 : L → +RF − LF L − F R+ p2 : R → −LF + RF R + F L− Pure curves
FASS curve construction
In order to complete the curve deﬁnition, it is necessary to specify the subﬁgures represented by symbols L and R. In the special case of pure curves [116], these subﬁgures are reduced to single points. Thus, one can assume that symbols L and R are erased (replaced by the empty string) at the end of the derivation. Alternatively, they can be left in the string and ignored by the turtle during string interpretation. This second approach is consistent with previous deﬁnitions of turtle interpretation [109, 112]. A general discussion of the relationship between recurrent formulas and Lsystems is presented in [61, 62]. Construction of the Lsystem generating the Hilbert curve can be extended to other FASS curves [116]. Consider an array of m×m square tiles, each including a smaller square, called a frame. The edges of the frame run at some distance from the tile’s edges. Each frame bounds an open selfavoiding polygon. The endpoints of this polygon coincide with the two contact vertices of the frame. Suppose that a singlestroke line running through all tiles can be constructed by connecting the contact vertices of neighboring frames using short horizontal or vertical line segments. A FASS curve can be constructed by the recursive repetition of this connecting pattern. To this end, the array of m × m connected tiles is considered a macrotile which contains an open polygon inscribed into a macroframe. An array of m × m macrotiles is formed, and the polygons inscribed into the macroframes are connected together. This construction is carried out recursively, with m × m macrotiles at level n yielding one macrotile at level n + 1. Tile arrangements suitable for the generation of FASS curves can be found algorithmically, by searching the space of all possible arrangements on a grid of a given size. Examples of curves synthesized this way are given in Figures 1.16 and 1.17. As in the case of edge rewriting, the relationship between node rewriting and tilings of the plane extends to branching structures. It oﬀers a method for synthesizing Lsystems that generate objects with a given recursive structure, and links methods for plant generation based on Lsystems with those using iterated function systems [7] (see Chapter 8).
1.4. Synthesis of DOLsystems
a
17
b
◦
n=3,δ=90 L L→LF+RFR+FLFLFLFLFRFR+ R→LFLF+RFRFR+F+RFLFLFR
n=2, δ=90◦ L L→LFLF+RFR+FLFLFRFLFLFR+F+RFLFLFRFRFR+ R→LFLFLF+RFR+FLFLF+RFR+ FLF+RFRFLFLFRFR
Figure 1.16: Sample FASS curves constructed using tiles with contact points positioned along a tile edge: (a) 3 × 3 tiles form a macrotile, (b) 4 × 4 tiles form a macrotile
a n=2, δ=90◦ L L→LFRFLFRFLFR+F+LFRFL R→RFLFR+F+LFRFLFRFLFR
b n=2, δ=45◦ L L→L+F+RFL+F+RFLFR+F+LFRFL+F+RFLFRFL+F+R+F+L+F+RFL+F+R+F+LRF+F+L+F+RFL+F+RFL R→RFL+F+RFL+F+R+F+LFR+F+L+F+RFL+F+R+F+L+F+ RFLFRFL+F+RFLFR+F+LFRFL+F+RFL+F+R Figure 1.17: Sample FASS curves constructed using tiles with contact points positioned diagonally: (a) 3 × 3 tiles form a macrotile (Peano curve [106]), (b) 5 × 5 tiles form a macrotile
18
1.4.3
Chapter 1. Graphical modeling using Lsystems
Relationship between edge and node rewriting
The classes of curves that can be generated using the edgerewriting and noderewriting techniques are not disjoint. For example, reconsider the Lsystem that generates the dragon curve using edge replacement: ω : Fl p1 : Fl → Fl + Fr + p2 : Fr → −Fl − Fr Assume temporarily that a production predecessor can contain more than one letter; thus an entire subword can be replaced by the successor of a single production (a formalization of this concept is termed pseudoLsystems [109]). The dragongenerating Lsystem can be rewritten as: ω : Fl p1 : F l → F l + rF + p2 : rF → −F l − rF where the symbols l and r are not interpreted by the turtle. Production p1 replaces the letter l by the string l + rF − while the leading letter F is left intact. In a similar way, production p2 replaces the letter r by the string −F l − r and leaves the trailing F intact. Thus, the Lsystem can be transformed into noderewriting form as follows: ω : Fl p1 : l → l + rF + p2 : r → −F l − r In practice, the choice between edge rewriting and node rewriting is often a matter of convenience. Neither approach oﬀers an automatic, general method for constructing Lsystems that capture given structures. However, the distinction between edge and node rewriting makes it easier to understand the intricacies of Lsystem operation, and in this sense assists in the modeling task. Speciﬁcally, the problem of ﬁlling a region by a selfavoiding curve is biologically relevant, since some plant structures, such as leaves, may tend to ﬁll a plane without overlapping [38, 66, 67, 94].
1.5
Modeling in three dimensions
Turtle interpretation of Lsystems can be extended to three dimensions following the ideas of Abelson and diSessa [1]. The key concept is to represent the current orientation of the turtle in space by three vectors H ,L, U , indicating the turtle’s heading, the direction to the left, and the direction up. These vectors have unit length, are perpendicular to each
1.5. Modeling in three dimensions
19
Figure 1.18: Controlling the turtle in three dimensions
other, and satisfy the equation H × L = U . Rotations of the turtle are then expressed by the equation
L U H
=
L U H
R,
where R is a 3 × 3 rotation matrix [40]. Speciﬁcally, rotations by angle α about vectors U ,L and H are represented by the matrices:
cos α sin α 0 RU (α) = − sin α cos α 0 0 0 1 cos α 0 − sin α 0 1 0 RL (α) = sin α 0 cos α
1 0 0 RH (α) = 0 cos α − sin α 0 sin α cos α The following symbols control turtle orientation in space (Figure 1.18): +
Turn left by angle δ, using rotation matrix RU (δ).
−
Turn right by angle δ, using rotation matrix RU (−δ).
&
Pitch down by angle δ, using rotation matrix RL (δ).
∧
Pitch up by angle δ, using rotation matrix RL (−δ).
\
Roll left by angle δ, using rotation matrix RH (δ).
/
Roll right by angle δ, using rotation matrix RH (−δ).

Turn around, using rotation matrix RU (180◦ ).
20
Chapter 1. Graphical modeling using Lsystems
n=2, A A → B → C → D →
δ=90◦ BF+CFC+FD&F∧DF+&&CFC+F+B// A&F∧CFB∧F∧D∧∧FD∧F∧BFC∧F∧A// D∧F∧BF+C∧F∧A&&FA&F∧C+F+B∧F∧D// CFBF+BFA&F∧A&&FBF+BFC//
Figure 1.19: A threedimensional extension of the Hilbert curve [139]. Colors represent threedimensional “frames” associated with symbols A (red), B (blue), C (green) and D (yellow).
1.6. Branching structures
21
As an example of a threedimensional object created using an Lsystem, consider the extension of the Hilbert curve shown in Figure 1.19. The Lsystem was constructed with the nodereplacement technique discussed in the previous section, using cubes and “macrocubes” instead of tiles and macrotiles.
1.6
Branching structures
According to the rules presented so far, the turtle interprets a character string as a sequence of line segments. Depending on the segment lengths and the angles between them, the resulting line is selfintersecting or not, can be more or less convoluted, and may have some segments drawn many times and others made invisible, but it always remains just a single line. However, the plant kingdom is dominated by branching structures; thus a mathematical description of treelike shapes and methods for generating them are needed for modeling purposes. An axial tree [89, 117] complements the graphtheoretic notion of a rooted tree [108] with the botanically motivated notion of branch axis.
1.6.1
Axial trees
A rooted tree has edges that are labeled and directed. The edge sequences form paths from a distinguished node, called the root or base, to the terminal nodes. In the biological context, these edges are referred to as branch segments. A segment followed by at least one more segment in some path is called an internode. A terminal segment (with no succeeding edges) is called an apex. An axial tree is a special type of rooted tree (Figure 1.20). At each of its nodes, at most one outgoing straight segment is distinguished. All remaining edges are called lateral or side segments. A sequence of segments is called an axis if: • the ﬁrst segment in the sequence originates at the root of the tree or as a lateral segment at some node, • each subsequent segment is a straight segment, and • the last segment is not followed by any straight segment in the tree. Together with all its descendants, an axis constitutes a branch. A branch is itself an axial (sub)tree. Axes and branches are ordered. The axis originating at the root of the entire plant has order zero. An axis originating as a lateral segment of an norder parent axis has order n+1. The order of a branch is equal to the order of its lowestorder or main axis.
22
Chapter 1. Graphical modeling using Lsystems
Figure 1.20: An axial tree
Figure 1.21: Sample tree generated using a method based on Horton– Strahler analysis of branching patterns
1.6. Branching structures
23
Figure 1.22: A tree production p and its application to the edge S in a tree T1
Axial trees are purely topological objects. The geometric connotation of such terms as straight segment, lateral segment and axis should be viewed at this point as an intuitive link between the graphtheoretic formalism and real plant structures. The proposed scheme for ordering branches in axial trees was introduced originally by Gravelius [53]. MacDonald [94, pages 110–121] surveys this and other methods applicable to biological and geographical data such as stream networks. Of special interest are methods proposed by Horton [70, 71] and Strahler, which served as a basis for synthesizing botanical trees [37, 152] (Figure 1.21).
1.6.2
Tree OLsystems
In order to model development of branching structures, a rewriting mechanism can be used that operates directly on axial trees. A rewriting rule, or tree production, replaces a predecessor edge by a successor axial tree in such a way that the starting node of the predecessor is identiﬁed with the successor’s base and the ending node is identiﬁed with the successor’s top (Figure 1.22). A tree OLsystem G is speciﬁed by three components: a set of edge labels V , an initial tree ω with labels from V , and a set of tree productions P . Given the Lsystem G, an axial tree T2 is directly derived from a tree T1 , noted T1 ⇒ T2 , if T2 is obtained from T1 by simultaneously replacing each edge in T1 by its successor according to the production set P . A tree T is generated by G in a derivation of length n if there exists a sequence of trees T0 , T1 , . . . , Tn such that T0 = ω, Tn = T and T0 ⇒ T1 ⇒ . . . ⇒ Tn .
24
Chapter 1. Graphical modeling using Lsystems
Figure 1.23: Bracketed string representation of an axial tree
1.6.3
Bracketed OLsystems
The deﬁnition of tree Lsystems does not specify the data structure for representing axial trees. One possibility is to use a list representation with a tree topology. Alternatively, axial trees can be represented using strings with brackets [82]. A similar distinction can be observed in Koch constructions, which can be implemented either by rewriting edges and polygons or their string representations. An extension of turtle interpretation to strings with brackets and the operation of bracketed Lsystems [109, 111] are described below. Two new symbols are introduced to delimit a branch. They are interpreted by the turtle as follows: Stack operations
2D structures
Bushlike structure
[
Push the current state of the turtle onto a pushdown stack. The information saved on the stack contains the turtle’s position and orientation, and possibly other attributes such as the color and width of lines being drawn.
]
Pop a state from the stack and make it the current state of the turtle. No line is drawn, although in general the position of the turtle changes.
An example of an axial tree and its string representation are shown in Figure 1.23. Derivations in bracketed OLsystems proceed as in OLsystems without brackets. The brackets replace themselves. Examples of twodimensional branching structures generated by bracketed OLsystems are shown in Figure 1.24. Figure 1.25 is an example of a threedimensional bushlike structure generated by a bracketed Lsystem. Production p1 creates three new branches from an apex of the old branch. A branch consists of an edge F forming the initial internode, a leaf L and an apex A (which will subsequently create three new branches). Productions p2 and p3
1.6. Branching structures
a
b
25
c
n=5,δ=25.7 F F →F[+F]F[F]F
n=5,δ=20 F F →F[+F]F[F][F]
n=4,δ=22.5◦ F F→FF[F+F+F]+ [+FFF]
d
e
f
◦
n=7,δ=20◦ X X →F[+X]F[X]+X F →FF
◦
n=7,δ=25.7◦ X X →F[+X][X]FX F →FF
n=5,δ=22.5◦ X X→F[[X]+X]+F[+FX]X F→FF
Figure 1.24: Examples of plantlike structures generated by bracketed OLsystems. Lsystems (a), (b) and (c) are edgerewriting, while (d), (e) and (f) are noderewriting.
26
Chapter 1. Graphical modeling using Lsystems
n=7, δ=22.5◦ ω p1 p2 p3 p4
: : : : :
A A F S L
→ → → →
[&FL!A]/////’[&FL!A]///////’[&FL!A] S ///// F F L [’’’∧∧{f+f+ff+f+f}]
Figure 1.25: A threedimensional bushlike structure
Plant with ﬂowers
specify internode growth. In subsequent derivation steps the internode gets longer and acquires new leaves. This violates a biological rule of subapical growth (discussed in detail in Chapter 3), but produces an acceptable visual eﬀect in a still picture. Production p4 speciﬁes the leaf as a ﬁlled polygon with six edges. Its boundary is formed from the edges f enclosed between the braces { and } (see Chapter 5 for further discussion). The symbols ! and are used to decrement the diameter of segments and increment the current index to the color table, respectively. Another example of a threedimensional plant is shown in Figure 1.26. The Lsystem can be described and analyzed in a way similar to the previous one.
1.6. Branching structures
n=5, δ=18◦ ω : p1 : p2 p3 p4 p5
: : : :
p6 : p7 :
plant plant → internode + [ plant + ﬂower] − − // [ − − leaf ] internode [ + + leaf ] − [ plant ﬂower ] + + plant ﬂower internode → F seg [// & & leaf ] [// ∧ ∧ leaf ] F seg seg → seg F seg leaf → [’ { +f−ﬀ−f+  +f−ﬀ−f } ] ﬂower → [ & & & pedicel ‘ / wedge //// wedge //// wedge //// wedge //// wedge ] pedicel → FF wedge → [‘ ∧ F ] [ { & & & & −f+f  −f+f } ] Figure 1.26: A plant generated by an Lsystem
27
28
1.7
Lsystem
Derivation
Example
Chapter 1. Graphical modeling using Lsystems
Stochastic Lsystems
All plants generated by the same deterministic Lsystem are identical. An attempt to combine them in the same picture would produce a striking, artiﬁcial regularity. In order to prevent this eﬀect, it is necessary to introduce specimentospecimen variations that will preserve the general aspects of a plant but will modify its details. Variation can be achieved by randomizing the turtle interpretation, the Lsystem, or both. Randomization of the interpretation alone has a limited eﬀect. While the geometric aspects of a plant — such as the stem lengths and branching angles — are modiﬁed, the underlying topology remains unchanged. In contrast, stochastic application of productions may aﬀect both the topology and the geometry of the plant. The following deﬁnition is similar to that of Yokomori [162] and Eichhorst and Savitch [35]. A stochastic OLsystem is an ordered quadruplet Gπ = V, ω, P, π. The alphabet V , the axiom ω and the set of productions P are deﬁned as in an OLsystem (page 4). Function π : P → (0, 1], called the probability distribution, maps the set of productions into the set of production probabilities. It is assumed that for any letter a ∈ V , the sum of probabilities of all productions with the predecessor a is equal to 1. We will call the derivation µ ⇒ ν a stochastic derivation in Gπ if for each occurrence of the letter a in the word µ the probability of applying production p with the predecessor a is equal to π(p). Thus, diﬀerent productions with the same predecessor can be applied to various occurrences of the same letter in one derivation step. A simple example of a stochastic Lsystem is given below. ω p1 p2 p3
Flower ﬁeld
: : : :
F .33 F → F [+F ]F [−F ]F .33 F → F [+F ]F .34 F → F [−F ]F
The production probabilities are listed above the derivation symbol →. Each production can be selected with approximately the same probability of 1/3. Examples of branching structures generated by this Lsystem with derivations of length 5 are shown in Figure 1.27. Note that these structures look like diﬀerent specimens of the same (albeit ﬁctitious) plant species. A more complex example is shown in Figure 1.28. The ﬁeld consists of four rows and four columns of plants. All plants are generated by a stochastic modiﬁcation of the Lsystem used to generate Figure 1.26.
1.7. Stochastic Lsystems
Figure 1.27: Stochastic branching structures
Figure 1.28: Flower ﬁeld
29
30
Chapter 1. Graphical modeling using Lsystems
The essence of this modiﬁcation is to replace the original production p3 by the following three productions: p3 : p3 : p 3 :
.33
seg → seg [ // & & leaf ] [// ∧∧ leaf ] F seg .33 seg → seg F seg .34 seg → seg
Thus, in any step of the derivation, the stem segment seg may either grow and produce new leaves (production p3 ), grow without producing new leaves (production p3 ), or not grow at all (production p 3 ). All three events occur with approximately the same probability. The resulting ﬁeld appears to consist of various specimens of the same plant species. If the same Lsystem was used again (with diﬀerent seed values for the random number generator), a variation of this image would be obtained.
1.8 Context in string Lsystems
Contextsensitive Lsystems
Productions in OLsystems are contextfree; i.e. applicable regardless of the context in which the predecessor appears. However, production application may also depend on the predecessor’s context. This eﬀect is useful in simulating interactions between plant parts, due for example to the ﬂow of nutrients or hormones. Various contextsensitive extensions of Lsystems have been proposed and studied thoroughly in the past [62, 90, 128]. 2Lsystems use productions of the form al < a > ar → χ, where the letter a (called the strict predecessor) can produce word χ if and only if a is preceded by letter al and followed by ar . Thus, letters al and ar form the left and the right context of a in this production. Productions in 1Lsystems have onesided context only; consequently, they are either of the form al < a → χ or a > ar → χ. OLsystems, 1Lsystems and 2Lsystems belong to a wider class of ILsystems, also called (k,l)systems. In a (k,l)system, the left context is a word of length k and the right context is a word of length l. In order to keep speciﬁcations of Lsystems short, the usual notion of ILsystems has been modiﬁed here by allowing productions with diﬀerent context lengths to coexist within a single system. Furthermore, contextsensitive productions are assumed to have precedence over contextfree productions with the same strict predecessor. Consequently, if a contextfree and a contextsensitive production both apply to a given letter, the contextsensitive one should be selected. If no production applies, this letter is replaced by itself as previously assumed for OLsystems.
1.8. Contextsensitive Lsystems
31
Figure 1.29: The predecessor of a contextsensitive tree production (a) matches edge S in a tree T (b)
The following sample 1Lsystem makes use of context to simulate signal propagation throughout a string of symbols:
Signal propagation
ω : baaaaaaaa p1 : b < a → b b →a p2 : The ﬁrst few words generated by this Lsystem are given below: baaaaaaaa abaaaaaaa aabaaaaaa aaabaaaaa aaaabaaaa ··· The letter b moves from the left side to the right side of the string. A contextsensitive extension of tree Lsystems requires neighbor edges of the replaced edge to be tested for context matching. A predecessor of a contextsensitive production p consists of three components: a path l forming the left context, an edge S called the strict predecessor, and an axial tree r constituting the right context (Figure 1.29). The asymmetry between the left context and the right context reﬂects the fact that there is only one path from the root of a tree to a given edge, while there can be many paths from this edge to various terminal nodes. Production p matches a given occurrence of the edge S in a tree T if l is a path in T terminating at the starting node of S, and r is a subtree
Context in tree Lsystems
32
Context in bracketed Lsystems
Chapter 1. Graphical modeling using Lsystems
of T originating at the ending node of S. The production can then be applied by replacing S with the axial tree speciﬁed as the production successor. The introduction of context to bracketed Lsystems is more diﬃcult than in Lsystems without brackets, because the bracketed string representation of axial trees does not preserve segment neighborhood. Consequently, the context matching procedure may need to skip over symbols representing branches or branch portions. For example, Figure 1.29 indicates that a production with the predecessor BC < S > G[H]M can be applied to symbol S in the string ABC[DE][SG[HI[JK]L]M N O], which involves skipping over symbols [DE] in the search for left context, and I[JK]L in the search for right context. Within the formalism of bracketed Lsystems, the left context can be used to simulate control signals that propagate acropetally, i.e., from the root or basal leaves towards the apices of the modeled plant, while the right context represents signals that propagate basipetally, i.e., from the apices towards the root. For example, the following 1Lsystem simulates propagation of an acropetal signal in a branching structure that does not grow: #ignore : +− ω : Fb [+Fa ]Fa [−Fa ]Fa [+Fa ]Fa p1 : Fb < Fa → Fb Symbol Fb represents a segment already reached by the signal, while Fa represents a segment that has not yet been reached. The #ignore statement indicates that the geometric symbols + and − should be considered as nonexistent while context matching. Images representing consecutive stages of signal propagation (corresponding to consecutive words generated by the Lsystem under consideration) are shown in Figure 1.30a. The propagation of a basipetal signal can be simulated in a similar way (Figure 1.30b): #ignore : +− ω : Fa [+Fa ]Fa [−Fa ]Fa [+Fa ]Fb p1 : Fa > Fb → Fb
1.8. Contextsensitive Lsystems
33
Figure 1.30: Signal propagation in a branching structure: (a) acropetal, (b) basipetal
The operation of contextsensitive Lsystems is examined further using examples obtained by Hogeweg and Hesper [64]. In 1974, they published the results of an exhaustive study of 3,584 patterns generated by a class of bracketed 2Lsystems deﬁned over the alphabet {0,1}. Some of these patterns had plantlike shapes. Subsequently, Smith signiﬁcantly improved the quality of the generated images using stateoftheart computer imagery techniques [136, 137]. Sample structures generated by Lsystems similar to those proposed by Hogeweg and Hesper are shown in Figure 1.31. The diﬀerences are related to the geometric interpretation of the resulting strings. According to the original interpretation, consecutive branches are issued alternately to the left and right, whereas turtle interpretation requires explicit speciﬁcation of branching angles within the Lsystem.
Lsystems of Hogeweg, Hesper and Smith
34
Chapter 1. Graphical modeling using Lsystems
Figure 1.31: Examples of branching structures generated using Lsystems based on the results of Hogeweg and Hesper [64]
1.8. Contextsensitive Lsystems
35
a n=30,δ=22.5◦
b n=30,δ=22.5◦
#ignore: F1F1F1 0 < 0 > 0 < 0 > 0 < 1 > 0 < 1 > 1 < 0 > 1 < 0 > 1 < 1 > 1 < 1 > * < + > * <  >
#ignore: F1F1F1 0 < 0 > 0 < 0 > 0 < 1 > 0 < 1 > 1 < 0 > 1 < 0 > 1 < 1 > 1 < 1 > * < + > * <  >
+F 0 1 0 1 0 1 0 1 * *
→ → → → → → → → → →
0 1[+F1F1] 1 1 0 1F1 0 0 +
c n=26, δ=25.75◦ #ignore: F1F1F1 0 < 0 > 0 < 0 > 0 < 1 > 0 < 1 > 1 < 0 > 1 < 0 > 1 < 1 > 1 < 1 > * <  > * < + >
→ → → → → → → → → →
0 1 0 1 0 1 0 1 * *
→ → → → → → → → → →
1 1[F1F1] 1 1 0 1F1 1 0 +
d n=24, δ=25.75◦
+F 0 1 0 1 0 1 0 1 * *
+F
#ignore: F0F1F1 0 < 0 > 0 < 0 > 0 < 1 > 0 < 1 > 1 < 0 > 1 < 0 > 1 < 1 > 1 < 1 > * < + > * <  >
0 1 0 1[+F1F1] 0 1F1 0 0 + 
+F 0 1 0 1 0 1 0 1 * *
→ → → → → → → → → →
1 0 0 1F1 1 1[+F1F1] 1 0 +
e n=26, δ=22.5◦ #ignore: F1F1F1 0 < 0 > 0 < 0 > 0 < 1 > 0 < 1 > 1 < 0 > 1 < 0 > 1 < 1 > 1 < 1 > * < + > * <  >
+F 0 1 0 1 0 1 0 1 * *
→ → → → → → → → → →
0 1[F1F1] 1 1 0 1F1 1 0 +
Figure 1.31 (continued): Lsystems of Hogeweg and Hesper
36
Chapter 1. Graphical modeling using Lsystems
1.9 Exponential growth
Basic properties
Growth functions
During the synthesis of a plant model it is often convenient to distinguish productions that specify the branching pattern from those that describe elongation of plant segments. This separation can be observed in some of the Lsystems considered so far. For example, in Lsystems (d), (e) and (f) from Figure 1.24 the ﬁrst productions capture the branching patterns, while the remaining productions, equal in all cases to F → F F , describe elongation of segments represented by sequences of symbols F . The number of letters F in a string χn originating from a single letter F is doubled in each derivation step, thus the elongation is exponential, with length(χn ) = 2n . A function that describes the number of symbols in a word in terms of its derivation length is called a growth function. The theory of Lsystems contains an extensive body of results on growth functions [62, 127]. The central observation is that the growth functions of DOLsystems are independent of the letter ordering in the productions and derived words. Consequently, the relation between the number of letter occurrences in a pair of words µ and ν, such that µ ⇒ ν, can be conveniently expressed using matrix notation. Let G = V, ω, P be a DOLsystem and assume that letters of the alphabet V have been ordered, V = {a1 , a2 , . . . , am }. Construct a square matrix Qm×m , where entry qij is equal to the number of occurrences of letter aj in the successor of the production with predecessor ai . Let aki denote the number of occurrences of letter ai in the word x generated by G in a derivation of length k. The deﬁnition of direct derivation in a DOLsystem implies that
ak1 ak2 · · · akm
q11 q12 · · · q1m q21 q22 · · · q2m .. . qm1 qm2 · · · qmm
= ak+1 1
ak+1 · · · ak+1 2 m
.
This matrix notation is useful in the analysis of growth functions. For example, consider the following Lsystem: ω : a p1 : a → ab p2 : b → a
(1.2)
The relationship between the number of occurrences of letters a and b in two consecutively derived words is
ak b k
1 1
1 0
=
ak+1 bk+1
1.9. Growth functions
or
37
ak+1 = ak + bk = ak + ak−1
for k = 1, 2, 3, . . . . From the axiom it follows that a0 = 1 and a1 = b0 = 0. Thus, the number of letters a in the strings generated by the Lsystem speciﬁed in equation (1.2) grows according to the Fibonacci series: 1, 1, 2, 3, 5, 8, . . . . This growth function was implemented by productions p2 and p3 in the Lsystem generating the bush in Figure 1.25 (page 26) to describe the elongation of its internodes. Polynomial growth functions of arbitrary degree can be obtained using Lsystems of the following form: ω p1 p2 p3 p4
: : : : :
a0 a0 a1 a2 a3 .. .
1 0 0 0
1 1 0 0
→ a0 a1 → a1 a2 → a2 a3 → a3 a4
The matrix Q is given below: Q=
0 1 1 0 .. .
0 0 1 1
··· ··· ··· ···
Thus, for any i, k ≥ 1, the number aki of occurrences of symbol ai in the string generated in a derivation of length k satisﬁes the equality aki + aki+1 = ak+1 i+1 . Taking into consideration the axiom, the distribution of letters ai as a function of the derivation length is captured by the following table (only nonzero terms are shown): k
ak0
ak1
ak2
0 1 2 3 4 5 6 7
1 1 1 1 1 1 1 1
1 2 3 4 5 6 7
1 3 6 10 15 21
ak3
ak4
ak5
ak6
ak7
1 4 10 20 35 .. .
1 5 15 35
1 6 21
1 7
1
Polynomial growth
38
Chapter 1. Graphical modeling using Lsystems
This table represents the Pascal triangle, thus for any k ≥ i ≥ 1 its terms satisfy the following equality:
aki
Characterization
=
k i
=
k(k − 1) · · · (k − i + 1) 1 · 2···i
Consequently, the number of occurrences of letter ai as a function of the derivation length k is expressed by a polynomial of degree i. By identifying letter ai with the turtle symbol F , it is possible to model internode elongation expressed by polynomials of arbitrary degree i ≥ 0. This observation was generalized by Szilard [140], who developed an algorithm for constructing a DOLsystem with growth functions speciﬁed by any positive, nondecreasing polynomials with integer coeﬃcients [62, page 276]. The examples of growth functions considered so far include exponential and polynomial functions. Rozenberg and Salomaa [127, pages 30–38] show that, in general, the growth function fG (n) of any DOLsystem G = V, ω, P is a combination of polynomial and exponential functions: fG (n) =
s
Pi (n)ρni for n ≥ n0 ,
(1.3)
i=1
Sigmoidal growth
where Pi (n) denotes a polynomial with integer coeﬃcients, ρi is a nonnegative integer, and n0 is the total number of letters in the alphabet V . Unfortunately, many growth processes observed in nature cannot be described by equation (1.3). Two approaches are then possible within the framework of the theory of Lsystems. The ﬁrst is to extend the size n0 of the alphabet V , so that the growth process of interest will be captured by the initial derivation steps, ω = µ0 ⇒ µ1 ⇒ · · · ⇒ µn0 , before equation (1.3) starts to apply. For example, the Lsystem ω : a0 pi : ai → ai+1 b0 for i < k pk+j : bj → bj+1 F for j < l
Squareroot growth
(1.4)
over the alphabet V = {a0 , a1 , ..., ak } ∪ {b0 , b1 , ..., bl } ∪ {F } can be used to approximate a sigmoidal elongation of a segment represented by a sequence of symbols F (Figure 1.32). The term sigmoidal refers to a function with a plot in the shape of the letter S. Such functions are commonly found in biological processes [143], with the initial part of the curve representing the growth of a young organism, and the latter part corresponding to the organism close to its ﬁnal size. The second approach to the synthesis of growth functions outside the class captured by equation (1.3) is to use contextsensitive Lsystems. For example, √ the following 2Lsystem has a growth function given by fG (n) = n + 4, where x is the ﬂoor function.
1.9. Growth functions
39
Figure 1.32: A sigmoidal growth function implemented using the Lsystem in equation (1.4), for k = l = 20
ω p1 p2 p3 p4 p5 p6
: : : : : : :
XFu Fa X Fu < Fa > Fa Fu < Fa > X Fa < Fa > Fd X < Fa > Fd Fu Fd
→ → → → → →
Fu Fd Fa Fd Fu Fa Fa
(1.5)
The operation of this Lsystem is illustrated in Figure 1.33. Productions p1 and p3 , together with p5 and p6 , propagate symbols Fu and Fd up and down the string of symbols µ. Productions p2 and p4 change the propagation direction, after symbol X marking a string end has been reached by Fu or Fd , respectively. In addition, p2 extends the string with a symbol Fa . Thus, the number of derivation steps increases by two between consecutive applications of production p2 . As a result, string extension occurs at derivation steps n expressed√by the square of the string length, which yields the growth function n + 4. In practice it is often diﬃcult, if not impossible, to ﬁnd Lsystems with the required growth functions. Vit´anyi [153] illustrates this by referring to sigmoidal curves: If we want to obtain sigmoidal growth curves with the original Lsystems then not even the introduction of cell interaction can help us out. In the ﬁrst place, we end up constructing quite unlikely ﬂows of messages through the organism, which are more suitable to electronic computers, and in fact give the organism universal computing power. Secondly, and this is more fundamental, we can not obtain
Limitations
40
Chapter 1. Graphical modeling using Lsystems
Figure 1.33: Squareroot growth implemented using the Lsystem speciﬁed in equation (1.5)
growth which, always increasing the size of the organism, tends towards stability in the limit. The slowest increasing growth we can obtain by allowing cell interaction is logarithmic and thus can not at all account for the asymptotic behavior of sigmoidal growth functions. In the next section we present an extension of Lsystems that makes it possible to avoid this problem by allowing for explicit inclusion of growth functions into Lsystem speciﬁcations.
1.10 Motivation
@ √ @ 2 1 @ @ @
1
Parametric Lsystems
Although Lsystems with turtle interpretation make it possible to generate a variety of interesting objects, from abstract fractals to plantlike branching structures, their modeling power is quite limited. A major problem can be traced to the reduction of all lines to integer multiples of the unit segment. As a result, even such a simple ﬁgure as an isosceles rightangled triangle cannot be traced exactly, since the ratio of its hypotenuse √ length to the length of a side is expressed by the irrational number 2. Rational approximation of line length provides only a limited solution, because the unit step must be the smallest common denominator of all line lengths in the modeled structure. Consequently, the representation of a simple plant module, such as an internode, may require a large number of symbols. The same argument applies to angles. Problems become even more pronounced while simulating changes to the modeled structure over time, since some growth functions cannot be expressed conveniently using Lsystems. Generally, it is diﬃcult
1.10. Parametric Lsystems
41
to capture continuous phenomena, since the obvious technique of discretizing continuous values may require a large number of quantization levels, yielding Lsystems with hundreds of symbols and productions. Consequently, model speciﬁcation becomes diﬃcult, and the mathematical beauty of Lsystems is lost. In order to solve similar problems, Lindenmayer proposed that numerical parameters be associated with Lsystem symbols [83]. He illustrated this idea by referring to the continuous development of branching structures and diﬀusion of chemical compounds in a nonbranching ﬁlament of Anabaena catenula. Both problems were revisited in later papers [25, 77]. A deﬁnition of parametric Lsystems was formulated by Prusinkiewicz and Hanan [113] and is presented below.
1.10.1
Parametric OLsystems
Parametric Lsystems operate on parametric words, which are strings of modules consisting of letters with associated parameters. The letters belong to an alphabet V , and the parameters belong to the set of real numbers . A module with letter A ∈ V and parameters a1 , a2 , ..., an ∈ is denoted by A(a1 , a2 , ..., an ). Every module belongs to the set M = V × ∗ , where ∗ is the set of all ﬁnite sequences of parameters. The set of all strings of modules and the set of all nonempty strings are denoted by M ∗ = (V × ∗ )∗ and M + = (V × ∗ )+ , respectively. The realvalued actual parameters appearing in the words correspond with formal parameters used in the speciﬁcation of Lsystem productions. If Σ is a set of formal parameters, then C(Σ) denotes a logical expression with parameters from Σ, and E(Σ) is an arithmetic expression with parameters from the same set. Both types of expressions consist of formal parameters and numeric constants, combined using the arithmetic operators +, −, ∗, /; the exponentiation operator ∧, the relational operators <, >, =; the logical operators !, &,  (not, and, or); and parentheses (). Standard rules for constructing syntactically correct expressions and for operator precedence are observed. Relational and logical expressions evaluate to zero for false and one for true. A logical statement speciﬁed as the empty string is assumed to have value one. The sets of all correctly constructed logical and arithmetic expressions with parameters from Σ are noted C(Σ) and E(Σ). A parametric OLsystem is deﬁned as an ordered quadruplet G = V, Σ, ω, P , where • V is the alphabet of the system, • Σ is the set of formal parameters, • ω ∈ (V × ∗ )+ is a nonempty parametric word called the axiom, • P ⊂ (V × Σ∗ ) × C(Σ) × (V × E(Σ))∗ is a ﬁnite set of productions.
Parametric words
Expressions
Parametric OLsystem
42
Chapter 1. Graphical modeling using Lsystems
The symbols : and → are used to separate the three components of a production: the predecessor, the condition and the successor. For example, a production with predecessor A(t), condition t > 5 and successor B(t + 1)CD(t ∧ 0.5, t − 2) is written as A(t) : t > 5 → B(t + 1)CD(t ∧ 0.5, t − 2). Derivation
(1.6)
A production matches a module in a parametric word if the following conditions are met: • the letter in the module and the letter in the production predecessor are the same, • the number of actual parameters in the module is equal to the number of formal parameters in the production predecessor, and • the condition evaluates to true if the actual parameter values are substituted for the formal parameters in the production.
Example
A matching production can be applied to the module, creating a string of modules speciﬁed by the production successor. The actual parameter values are substituted for the formal parameters according to their position. For example, production (1.6) above matches a module A(9), since the letter A in the module is the same as in the production predecessor, there is one actual parameter in the module A(9) and one formal parameter in the predecessor A(t), and the logical expression t > 5 is true for t = 9. The result of the application of this production is a parametric word B(10)CD(3, 7). If a module a produces a parametric word χ as the result of a production application in an Lsystem G, we write a → χ. Given a parametric word µ = a1 a2 ...am , we say that the word ν = χ1 χ2 ...χm is directly derived from (or generated by) µ and write µ =⇒ ν if and only if ai → χi for all i = 1, 2, ..., m. A parametric word ν is generated by G in a derivation of length n if there exists a sequence of words µ0 , µ1 , ..., µn such that µ0 = ω, µn = ν and µ0 =⇒ µ1 =⇒ ... =⇒ µn . An example of a parametric Lsystem is given below. ω p1 p2 p3 p4
: : : : :
B(2)A(4, 4) A(x, y) : y <= 3 A(x, y) : y > 3 B(x) :x<1 B(x) : x >= 1
→ → → →
A(x ∗ 2, x + y) B(x)A(x/y, 0) C B(x − 1)
(1.7)
As in the case of nonparametric Lsystems, it is assumed that a module replaces itself if no matching production is found in the set P . The words obtained in the ﬁrst few derivation steps are shown in Figure 1.34.
1.10. Parametric Lsystems
43
Figure 1.34: The initial sequence of strings generated by the parametric Lsystem speciﬁed in equation (1.7)
1.10.2
Parametric 2Lsystems
Productions in parametric OLsystems are contextfree, i.e., applicable regardless of the context in which the predecessor appears. A contextsensitive extension is necessary to model information exchange between neighboring modules. In the parametric case, each component of the production predecessor (the left context, the strict predecessor and the right context) is a parametric word with letters from the alphabet V and formal parameters from the set Σ. Any formal parameters may appear in the condition and the production successor. A sample contextsensitive production is given below:
Example
A(x) < B(y) > C(z) : x + y + z > 10 → E((x + y)/2)F ((y + z)/2) It can be applied to the module B(5) that appears in a parametric word · · · A(4)B(5)C(6) · · ·
(1.8)
since the sequence of letters A, B, C in the production and in parametric word (1.8) are the same, the numbers of formal parameters and actual parameters coincide, and the condition 4 + 5 + 6 > 10 is true. As a result of the production application, the module B(5) will be replaced by a pair of modules E(4.5)F (5.5). Naturally, the modules A(4) and C(6) will be replaced by other productions in the same derivation step. Parametric 2Lsystems provide a convenient tool for expressing developmental models that involve diﬀusion of substances throughout an organism. A good example is provided by an extended model of the pattern of cells observed in Anabaena catenula and other bluegreen bacteria [99]. This model was proposed by de Koster and Lindenmayer [25].
Anabaena with heterocysts
44
Chapter 1. Graphical modeling using Lsystems
#define CH 900 #define CT 0.4 #define ST 3.9 #include H #ignore f ∼ H ω : p1 : p2 : p3 : p4 : p5 :
/* /* /* /*
high concentration */ concentration threshold */ segment size threshold */ heterocyst shape specification */
(90)F(0,0,CH)F(4,1,CH)F(0,0,CH) F(s,t,c) : t=1 & s>=6 → F(s/3*2,2,c)f(1)F(s/3,1,c) F(s,t,c) : t=2 & s>=6 → F(s/3,2,c)f(1)F(s/3*2,1,c) F(h,i,k) < F(s,t,c) > F(o,p,r) : s>STc>CT → F(s+.1,t,c+0.25*(k+r3*c)) F(h,i,k) < F(s,t,c) > F(o,p,r) : !(s>STc>CT) → F(0,0,CH) ∼ H(1) H(s) : s<3 → H(s*1.1) Lsystem 1.1: Anabaena catenula
Generally, the bacteria under consideration form a nonbranching ﬁlament consisting of two classes of cells: vegetative cells and heterocysts. Usually, the vegetative cells divide and produce two daughter vegetative cells. This mechanism is captured by the Lsystem speciﬁed in equation (1.1) and Figure 1.4 (page 5). However, in some cases the vegetative cells diﬀerentiate into heterocysts. Their distribution forms a welldeﬁned pattern, characterized by a relatively constant number of vegetative cells separating consecutive heterocysts. How does the organism maintain constant spacing of heterocysts while growing? The model explains this phenomenon using a biologically wellmotivated hypothesis that heterocyst distribution is regulated by nitrogen compounds produced by the heterocysts, transported from cell to cell across the ﬁlament, and decayed in the vegetative cells. If the compound’s concentration in a young vegetative cell falls below a speciﬁc level, this cell diﬀerentiates into a heterocyst (Lsystem 1.1). The #deﬁne statements assign values to numerical constants used in the Lsystem. The #include statement speciﬁes the shape of a heterocyst (a disk) by referring to a library of predeﬁned shapes (see Section 5.1). Cells are represented by modules F (s, t, c), where s stands for cell length, t is cell type (0  heterocyst, 1 and 2  vegetative types1 ), 1
The model of Anabaena introduced in Section 1.2 distinguished between four types of cells: ar , br , al and bl . Cells b do not divide and can be considered as young forms of the corresponding cells a. Thus, the vegetative type 1 considered here embraces cells ar and br , while type 2 embraces al and bl . The formal relationship between the fourcell and twocell models is further discussed in Chapter 6.
1.10. Parametric Lsystems
45
Figure 1.35: Development of Anabaena catenula with heterocysts, simulated using parametric Lsystem 1.1
and c represents the concentration of nitrogen compounds. Productions p1 and p2 describe division of the vegetative cells. They each create two daughter cells of unequal length. The diﬀerence between cells of type 1 and 2 lies in the ordering of the long and short daughter cells. Production p3 captures the processes of transportation and decay of the nitrogen compounds. Their concentration c is related to the concentration in the neighboring cells k and r, and changes in each derivation step according to the formula c = c + 0.25(k + r − 3 ∗ c). Production p4 describes diﬀerentiation of a vegetative cell into a heterocyst. The condition speciﬁes that this process occurs when the cell length does not exceed the threshold value ST = 3.9 (which means that the cell is young enough), and the concentration of the nitrogen compounds falls below the threshold value CT = 0.4. Production p5 describes the subsequent growth of the heterocyst. Snapshots of the developmental sequence of Anabaena are given in Figure 1.35. The vegetative cells are shown as rectangles, colored according to the concentration of the nitrogen compounds (white means low concentration). The heterocysts are represented as red disks. The values of parameters CH, CT and ST were selected to provide the
46
Chapter 1. Graphical modeling using Lsystems
correct distribution of the heterocysts, and correspond closely to the values reported in [25]. The mathematical model made it possible to estimate these parameters, although they are not directly observable.
1.10.3
Turtle interpretation of parametric words
If one or more parameters are associated with a symbol interpreted by the turtle, the value of the ﬁrst parameter controls the turtle’s state. If the symbol is not followed by any parameter, default values speciﬁed outside the Lsystem are used as in the nonparametric case. The basic set of symbols aﬀected by the introduction of parameters is listed below.
Row of trees
F (a)
Move forward a step of length a > 0. The position of the turtle changes to (x , y , z ), where x = x + aH x , y = y + aH y and z = z + aH z . A line segment is drawn between points (x, y, z) and (x , y , z ).
f (a)
Move forward a step of length a without drawing a line.
+(a)
Rotate around U by an angle of a degrees. If a is positive, the turtle is turned to the left and if a is negative, the turn is to the right.
&(a)
Rotate around L by an angle of a degrees. If a is positive, the turtle is pitched down and if a is negative, the turtle is pitched up.
/(a)
Rotate around H by an angle of a degrees. If a is positive, the turtle is rolled to the right and if a is negative, it is rolled to the left.
It should be noted that symbols +, &, ∧, and / are used both as letters of the alphabet V and as operators in logical and arithmetic expressions. Their meaning depends on the context. The following examples illustrate the operation of parametric Lsystems with turtle interpretation. The ﬁrst Lsystem is a coding of a Koch construction generating a variant of the snowﬂake curve (Figure 1.1 on page 2). The initiator (production predecessor) is the hypotenuse AB of a right triangle ABC (Figure 1.36). The ﬁrst and the fourth edge of the generator subdivide AB into segments AD and DB, while the remaining two edges traverse the altitude CD in opposite directions. From elementary geometry it follows that the lengths of these segments satisfy the equations √ q =c−p and h = pq. The edges of the generator can be associated with four triangles that are similar to ABC and tile it without gaps. According to the relationship between curve construction by edge rewriting and planar tilings
1.10. Parametric Lsystems
47
Figure 1.36: Construction of the generator for a “row of trees.” The edges are associated with triangles indicated by ticks.
(Section 1.4.1), the generated curve will approximately ﬁll the triangle ABC. The corresponding Lsystem is given below: #deﬁne #deﬁne #deﬁne #deﬁne
c p q h
1 0.3 c−p (p ∗ q) ∧ 0.5
ω : F (1) p1 : F (x) → F (x ∗ p) + F (x ∗ h) − −F (x ∗ h) + F (x ∗ q) The resulting curve is shown in Figure 1.37a. In order to better visualize its structure, the angle increment has been set to 86◦ instead of 90◦ . The curve ﬁlls diﬀerent areas with unequal density. This results from the fact that all edges, whether long or short, are replaced by the generator in every derivation step. A modiﬁed curve that ﬁlls the underlying triangle in a more uniform way is shown in Figure 1.37b. It was obtained by delaying the rewriting of shorter segments with respect to the longer ones, as speciﬁed by the following Lsystem. ω : F (1, 0) p1 : F (x, t) : t = 0 → F (x ∗ p, 2) + F (x ∗ h, 1)− −F (x ∗ h, 1) + F (x ∗ q, 0) p2 : F (x, t) : t > 0 → F (x, t − 1) The next example makes use of node rewriting (Section 1.4.2). The construction recursively subdivides a rectangular tile ABCD into two tiles, AEF D and BCF E, similar to ABCD (Figure 1.38). The lengths of the edges form the proportion a b = 1 , b a 2
Branching structure
48
Chapter 1. Graphical modeling using Lsystems
Figure 1.37: Two curves suggesting a “row of trees.” Curve (b) is from [95, page 57].
√ which implies that b = a/ 2. Each tile is associated with a singlepoint frame lying in the tile center. The tiles are connected by a branching line speciﬁed by the following Lsystem: #deﬁne R 1.456 ω : A(1) p1 : A(s) → F (s)[+A(s/R)][−A(s/R)]
(1.9)
The √ ratio of branch sizes R slightly exceeds the theoretical value of 2. As a result, the branching structure shown in Figure 1.39 is selfavoiding. The angle increment was set arbitrarily to δ = 85◦ . The Lsystem in equation (1.9) operates by appending segments of decreasing length to the structures obtained in previous derivation steps. Once a segment has been incorporated, its length does not change. A structure with identical proportions can be obtained by
1.10. Parametric Lsystems
49
Figure 1.38: Tiling associated with a spaceﬁlling branching pattern
Figure 1.39: A branching pattern generated by the Lsystem speciﬁed in equation (1.9)
50
Chapter 1. Graphical modeling using Lsystems
Figure 1.40: Initial sequences of ﬁgures generated by the Lsystems speciﬁed in equations (1.9) and (1.10)
appending segments of constant length and increasing the lengths of previously created segments by constant R in each derivation step. The corresponding Lsystem is given below. ω : A p1 : A → F (1)[+A][−A] p2 : F (s) → F (s ∗ R)
(1.10)
The initial sequence of structures obtained by both Lsystems are compared in Figure 1.40. Sequence (a) emphasizes the fractal character of the resulting structure. Sequence (b) suggests the growth of a tree. The next two chapters show that this is not a mere coincidence, and the Lsystem speciﬁed in equation (1.10) is a simple, but in principle correct, developmental model of a sympodial branching pattern found in many herbaceous plants and trees.