Chapter Four:

Iterative Graphics

4.1 Roses

rect[theta_, r_] := {r Cos[theta], r Sin[theta]} //N
ParametricPlot[rect[theta, Sin[2 theta]], {theta, 0, 2 Pi},
			Ticks->None, AspectRatio->1]

rose[i_] := rect[i Degree, Sin[n i Degree]]
	Attributes[rose] = {Listable}

n = 2;
d = 71;
walk = Line[rose[Range[0, 360 d, d]]];
Show[Graphics[{Thickness[.0005], walk}], AspectRatio->1]

n = 2; d = 30;

cosets = Table[Line[rose[Range[m, m+LCM[d, 360], d]]],
			{m, 0, 0}];

Show[Graphics[{Thickness[.0005], cosets}], AspectRatio->1]

n = 2; d = 180;

cosets = Table[Line[rose[Range[m, m + LCM[d, 360], d]]],
		{m, 0, GCD[d, 360] - 1}];

Show[Graphics[{Thickness[.0005], cosets}], AspectRatio->1]

MaurerRose::usage = "MaurerRose[n, d, z:360] shows the
walk along an n- (or 2n-) leafed rose in steps of d degrees 
(including cosets). Varying z changes the size of a 'degree'.";

MaurerRose[n_Integer, d_Integer, z_Integer:360] :=
   Block[{zz = If[!EvenQ[n] && z == 360, 180, z], cosets},

   rect[theta_, r_] := {r Cos[theta], r Sin[theta]} //N;

   rose[i_] := rect[i 2 Pi / z, Sin[n i 2 Pi / z]];
		Attributes[rose] = {Listable};

   cosets = Table[Line[rose[Range[m, m + LCM[d, zz], d]]],
								{m, 0, GCD[d, zz] - 1}];

   Show[Graphics[{Thickness[.0005], cosets}], AspectRatio->1]]

MaurerRose[4,120];

MaurerRose[6,72];

MaurerRose[2,150]

MaurerRose[90, 90, 359]

MaurerRose[181, 90, 359]

4.2 The Cantor Function

spawn[{a_Rational, b_Rational}] :=
	Block[{w = (b - a)/3}, {{a - 2w, a - w}, {b + w, b + 2w}}]

spawn[{1/3, 2/3}]

spawn[{1/9,2/9}]

spawn[intervals_List] := Flatten[Map[spawn, intervals], 1] /; 
	Length[intervals[[1]]] > 1

spawn[spawn[{1/3, 2/3}]]

removedintervals[n_] := 
	Flatten[NestList[spawn, {{1/3, 2/3}}, n - 1], 1]
removedintervals[4]

BaseForm[removedintervals[3]//N, 3]

spawn[{a_Rational, b_Rational}] :=
	Block[{w = (b - a)/3}, {{a - 2w, a - w}, {b + w, b + 2w}}]
spawn[intervals_List] := Flatten[Map[spawn, intervals], 1] /; 
Length[intervals[[1]]] > 1 
removedintervals[n_] := 
	Flatten[NestList[spawn, {{1/3, 2/3}}, n - 1], 1]

f[0]  =  0;
f[x_] := f[3x]/2       /;  x   <= 1/3;
f[x_] := 1 - f[1 - x]  /;  2/3 <= x;

connect[{a_, b_}] := Block[{y = f[a]}, Line[{{a, y}, {b, y}}]]

CantorFunction[levels_] := Show[Graphics[
 		{Thickness[.001], 
 		 Map[connect, removedintervals[levels]]}],
		 Axes->{0, 0},
		 Ticks->{N[Range[0., 1, 1/9], 3], Automatic}]

CantorFunction[7]

4.4 Penrose Tiles

phi = GoldenRatio //N;
toppoint = {1/2, Sin[72 Degree] phi} //N;
starttriangle = {toppoint, {1, 0.}, {0., 0.}, 1};

dissect[{p_, q_, r_, 1}] :=
			(newpoint = (phi q + p) * (2 - phi);
			{{r, newpoint, q, 1}, {r, newpoint, p, -1}})

dissect[{p_, q_, r_, -1}] := (newpoint = (phi r + p)*(2 - phi);
		{{r, newpoint, q, -1}, {p, q, newpoint, 1}})

dissect[list_] := (type *= -1;
    Select[list, Last[#] != type &] 
			~Join~
    Flatten[Map[dissect, Select[list, Last[#] == type &]],1] /; 
									TensorRank[list] != 1)

TriangleDissection[n_] := (type = 1;
 tiles = Nest[dissect, starttriangle, n] /.
	      {p_, q_, r_, x_} -> 
		     {Line[{q, r, p}], Thickness[.015], Line[{p, q}]};
 Show[Graphics[{Thickness[.0005], tiles}],
		AspectRatio->toppoint[[2]],
		PlotRange->All])

Attributes[TriangleDissection] = Listable;

TriangleDissection[5]

KitesAndDarts[n_] := (type = 1;
	tiles = Nest[dissect, starttriangle, 2 n] /.
  	   {{p_, q_, r_, -1} -> Line[{q, r, p}],
	   {p_ ,q_, r_,  1} -> {GrayLevel[.5], Polygon[{q, r, p}], 
	    			      GrayLevel[0],  Line[{q, r, p}]}};
Show[Graphics[{Thickness[.0005], tiles}],
	AspectRatio->toppoint[[2]], PlotRange->All])

KitesAndDarts[5]

First[x/.N[Solve[Tan[72 Degree] x == 1-2x, x]]]

Rhombs[n_] := (type = 1;
			tiles = Nest[dissect, starttriangle, 2 n + 1] /.
      {  {p_, q_, r_, 1} -> Line[{r, p, q}],
         {p_, q_, r_, -1} ->
		  {{GrayLevel[.5], Polygon[{p,q,r}]}, Line[{p,q,r}]}};
	Show[Graphics[{Thickness[.003], tiles}],
	AspectRatio->toppoint[[2]], PlotRange->All])

Rhombs[5]

4.5 Dynamics of the Quadratic Map

NestList[Cos, 1.5, 20]

NestList[Cos, Nest[Cos, 1.5, 100], 10]

Table[Nest[Cos, Random[Real, {-100, 100}], 100], {20}]

Plot[Nest[Cos, x, 4], {x, 0, 2 Pi}, PlotRange->{-.5, 1}];

Trajectory::usage = "Trajectory[f, x0, initial, length, xmin:0,
 xmax:1] returns a graphical trace of the iterations of f. 
 Starting value is x0; orbit consists of orbitlength many 
 iterations starting from the initial'th; the graph is 
 displayed between xmin and xmax.";

Trajectory[f_, x0_, initial_, orbitlength_, xmin_:0, xmax_:1] :=
	Block[{start, orbit, plot, lines},
	  start = Nest[f, N[x0], initial];
	  orbit = NestList[f, start, orbitlength];
	  plot = Plot[f[x], {x, xmin, xmax}, 
						DisplayFunction->Identity];
	  lines = Line[Rest[Partition[
		Flatten[Transpose[{orbit, orbit}]], 2, 1]]];

      Show[plot,
      	Graphics[{{Thickness[.0001], PointSize[.02],
      	lines,
		Point[{start, f[start]}],
		Line[{{xmin, xmin}, {xmax, xmax}}]}}],
	   Axes->{xmin, xmin},
	   DisplayFunction->$DisplayFunction,
	   PlotRange->{{xmin, xmax}, {xmin, xmax}}]]

Attributes[Trajectory] = Listable;

Trajectory[Cos, 1.5, 0, 20, 0, 1.6]

f[x_] := r x (1 - x)
r = 2
Trajectory[f, .1, 0, 10, 0, 1]

r = 3.45
Trajectory[f, .1, 0, 100, 0, .9]
Trajectory[f, .1, 300, 50, .3, .9]

r = 3.839
Trajectory[f, .1, 150, 50, 0, 1]

r = 4
Trajectory[f, .1, 200, 75, 0, 1]
Trajectory[f, .1, 1000, 50, 0, 1]

r = 3.45
Table[NestList[f, Nest[f, Random[Real, {0, 1}], 400], 3], {10}]

r = 345/100;
f[49/69]

r = 4; n = 80;
traj  = Transpose[{Range[0, n], NestList[f, 0.1, n]}];
traj1 = Transpose[{Range[0, n], NestList[f, 0.1 + 10^-10, n]}];
	
Show[ListPlot[traj, PlotJoined->True,
					PlotStyle->Thickness[.0005]],
	ListPlot[traj1, PlotJoined->True,
					PlotStyle->Thickness[.003]]]

OrbitDiagram::usage = "OrbitDiagram[start, end, n, init, final]
 produces an orbit diagram of the quadratic map as r varies 
 from start to end in n steps. The number of initial iterations
 is init; the number of displayed iterations is final.";

f = r # (1-#) &;

(* following is the version with plotstart and plotend *)

OrbitDiagram[
	start_, end_, plotstart_, plotend_, n_, init_, final_] :=  

Show[Graphics[{PointSize[.001], Table[
   Map[Point[{r, #}] &, NestList[f, Nest[f, .5, init], final]],
	   					{r, end, start, (start - end)/n}]}],
	PlotRange->{{plotstart,plotend},{0, 1}}, 
	Axes->{plotstart, 0}]

Following two cells take a while

OrbitDiagram[2.9,  3.45, 2.9, 4, 300, 120, 4];

OrbitDiagram[3.45, 4,    2.9, 4, 300, 120, 150];

Clear[r]
Do[Trajectory[f, 0.5, 500, 5, 0, 1], {r, 2.9, 3.5, .025}]

f[f[x]]

g[x_] := r^2*(1 - x)*x*(1 - r*(1 - x)*x)
r = 3.2;
Do[Trajectory[g, x0, 0, 10], {x0, .1, .9, .15}]

r = .; Do[Trajectory[g, .1, 0, 10], {r, 2.9, 3.2, .05}]

g[r_, x_] := r^2*(1 - x)*x*(1 - r*(1 - x)*x)

(* this plot takes a while *)

Plot3D[{z = g[r,x], GrayLevel[If[Abs[z -x] <=.0028, 0, 1]]}, 
							{r, 2.8, 3.2}, {x, 0, 1},
PlotRange->{0,.9}, Shading->False, Lighting-> False,
AxesLabel->{"r", "x", None},
ClipFill->None, Boxed-> False, PlotPoints->{30,90},
ViewPoint-> {2, .8, 1.5}];

f = r # (1-#) &; limitr = 3.569945557391440;
Show[Graphics[{PointSize[.001],
	Table[(r = limitr - 10 ^ logr;
         Map[Point[{-logr, #}] &,
         NestList[f, Nest[f, .5, 150], 32]]),
    {logr, start = Log[10, limitr-2.5], -3, (-3-start)/200}]}], 
  PlotRange->{0.3, .9}, Axes->{0, 0.3}];

4.6 Variations on Circular Motion

ParametricPlot[
	(2 + Sin[8 t]/2) {Cos[tt = (t + Sin[16 t]/4)], Sin[tt]},
		{t, 0, 2 Pi},
	MaxBend->5, PlotPoints->50, PlotDivision->30,
	AspectRatio->1]

f[t_, a_, b_, c_] :=
    (2 + Sin[a t]/2) * {Cos[tt = (t + Sin[b t]/c)], Sin[tt]};

Attributes[f] = Listable;

Starr[a_, b_, c_, n_] := Block[{r}, 
	 values = f[Range[0., 2 Pi, 2 Pi/n //N], a, b, c];
	 Show[Graphics[{Thickness[.0005],
	       Map[Line[{{0, 0}, #}] &, values],
 	       Table[Line[values r ^ 0.2 ], {r, 0., 1, .05}]}],
      PlotRange->{{-2.5, 2.5}, {-2.5, 2.5}},
      AspectRatio->1]]

Starr[8, 16,  4, 720]

Starr[9,  6,  6, 720]

Starr[8, 16, 16, 720]

Starr[6, 18, 18, 720]