Rendering Stateless L-System in Mathematica

My interest in this topic was started by a code golf race started by the “Tweet-a-Program” campaign by Wolfram in 2015. Twitter limits the number of characters per tweet, so one should pull some smart tricks to cram in coolness into a tweet. Although I was not the first one nor the only one who posted fractals, I can say that, for this specific topic, my snippet was the most readable, extendible, and concise. Also, I actually learned a lot about Mathematica while doing this - much more than my few years of experience that I already had.

L-System, a.k.a Lindenmayer System, is a string rewriting system, which is often used for generating self-similar fractals. I works by iteratively applying one or more rules to the current state, starting from the initial state. If you have heard about Logo, the turtle drawing language, L-system is pretty much a recursive Logo functions.

The reason why “Stateless” is added to the title is that this post only covers non-state-aware L-systems. Some L-systems depend on save-state feature, which makes sense only in sequential state-aware interpreter. An example of such system is “Fractal Plant”. (It might be possible to simulate this behaviour by tracing back the lines, but it’s too much of work.)

Implementing L-System interpreter

The following is the code from the above tweet:

Graphics@Line@FoldList[Plus,{0,0},
{Sin@#,Cos@#}&/@Accumulate@Flatten[Pi/3
Nest[#/.f->{f,1,f,4,f,1,f}&,{f,4,f,4,f},6]/.f->{}]]

The most important bit here is "/."(ReplaceAll) operator. In fact, Wolfram/Mathematica language itself is a term-rewriting system, and its standard library offers pattern matching functions. That is, implementing L-system, a string rewriting system, in Wolfram language is just a matter of converting the string into a list of terms, which is far much trivial than implementing a full interpreter.

For example, Koch Snowflake can be defined like:

({{Production Rule, :, F → F+F--F+F}, {Initial State, :, F--F--F}})

... and the below code applies the production rule to the initial state twice:

Flatten[Nest[#/.f->{f,1,f,4,f,1,f}&,{f,4,f,4,f},2]]

This is it. This is how we get the state of n-th level. All that’s left is rendering the final state to an actual image, which is also relatively easy using Mathematica’s symbolic graphics library.

However, the code from the tweet is very limited. It assumes that every rotation(+,–) is followed by “Go Forward” symbol, which isn’t right in many cases. The end result looks like this:

Options[RenderStatelessLSystem]={
"ExpandOnlyVars"-> {},
"Level"-> 4,
"Angle"->N[(π/2)]
};
RenderStatelessLSystem[rules_,init_,OptionsPattern[]]:=
Module[{state,steps},
With[{
level=OptionValue["Level"],
angle=OptionValue["Angle"],
ruleRemoveExpOnlyVars=Map[#->Sequence[]&,OptionValue["ExpandOnlyVars"]],
ruleRemoveVars=Map[First[#]->Sequence[]&,rules]
},
state=Flatten[Nest[#/.rules&,Flatten@{0,init},level]];
steps=Reap[Fold[
If[NumberQ[#2],(#1+#2),
Sow[#1];Sow[#2];0]&,
0,state/.ruleRemoveExpOnlyVars
]][[2,1]]/.ruleRemoveVars;
Graphics@Line@FoldList[Plus,{0.,0.},
Map[{Cos@#,Sin@#}&,Accumulate[angle×steps]]]
]];

For the sake of documenting the internals, the function has 3 distinctive steps. First, it computes the target state. Second, it reforms the state into a list of rotation angles b/w lines. Notably, this step is implemented using Sow-Reap pattern. Third, draw lines based on the angles.

Below is just a utility function to simplify export and display.

SetAttributes[runAndSave,HoldFirst];
runAndSave[body_,lvls_,name_]:=
With[{
imgs=ParallelTable[ReleaseHold[body],{lvl,lvls}]
},
Print@GraphicsGrid[
Partition[imgs,2],ImageSize->Full];
ParallelDo[With[{lvl=lvls[[i]],img=imgs[[i]]},
Export[name<>"-"<>ToString[lvl]<>".png",
Show[img,ImageSize->2048],"PNG",Background->None]
],{i,1,Length@imgs}]];

NOTE: you can find the exported image at the bottom of this page.

Koch Snowflake

({{Variables, :, F}, {Constants, :, + -}, {Init String, :, F--F--F}, {Rules, :, F → F+F--F+F}})

runAndSave[Module[{F,G},
RenderStatelessLSystem[
{F->{F,1,F,-2,F,1,F}},
{F,-2,F,-2,F},
"Angle"->N[π/3],
"Level"->lvl
]],{2,3,5,8},"koch-snowflake"]

Sierpinski Triangle

({{Variables, :, F G}, {Constants, :, + -}, {Init String, :, F-G-F}, {Rules, :, F → F-G+F+G-F}, { ,  , G->GG}, {Angle, :, 120°}})

runAndSave[Module[{F,G},
RenderStatelessLSystem[{F->{F,-1,G,1,F,1,G,-1,F},G->{G,G}},
{F,-1,G,-1,F},"Angle"->N[(2/3)π],"Level"->lvl
]],{2,4,6,10},"sierpinski-triangle"]

An approximation of Sierpinski Triangle using Sierpinski arrowhead curve

({{Variables, :, F G}, {Constants, :, + -}, {Init String, :, F-G-F}, {Rules, :, F → F-G+F+G-F}, { ,  , G->GG}, {Angle, :, 60°}})

runAndSave[Module[{F,G},
RenderStatelessLSystem[
{F->{G,-1,F,-1,G},G->{F,1,G,1,F}},
{F},"Angle"->N[(π/3)],"Level"->lvl
]],{2,4,6,10},"sierpinski-arrowhead"]

Koch Curve

runAndSave[Module[{F,G},
RenderStatelessLSystem[
{F->{F,1,F,-1,F,-1,F,1,F}},
{F},"Level"->lvl
]],{2,3,4,8},"koch-curve"]

Dragon Curve

runAndSave[Module[{X,Y,F},
RenderStatelessLSystem[
{X->{X,1,Y,F,1},Y->{-1,F,X,-1,Y},F->F},
{F,X},"Level"->lvl,"ExpandOnlyVars"->{X,Y}
]],{6,8,10,12},"dragon-curve"]

Cesàro Fractal

runAndSave[Module[{X,Y,F},
RenderStatelessLSystem[
{F->{F,1,F,-2,F,1,F}},
{F},"Level"->lvl,"Angle"->N[π/180*85]
]
],{2,4,6,8},"cesaro-fractal"]

Minkowski sausage

runAndSave[Module[{X,Y,F},
RenderStatelessLSystem[
{F->{F,1,F,-1,F,-1,F,F,1,F,1,F,-1,F}},
{F},"Level"->lvl
]
],{1,2,3,5},"minkowski-sausage"]