(*^ ::[ Information = "This is a Mathematica Notebook file. It contains ASCII text, and can be transferred by email, ftp, or other text-file transfer utility. It should be read or edited using a copy of Mathematica or MathReader. If you received this as email, use your mail application or copy/paste to save everything from the line containing (*^ down to the line containing ^*) into a plain text file. On some systems you may have to give the file a name ending with ".ma" to allow Mathematica to recognize it as a Notebook. The line below identifies what version of Mathematica created this file, but it can be opened using any other version as well."; FrontEndVersion = "NeXT Mathematica Notebook Front End Version 2.2"; NeXTStandardFontEncoding; fontset = title, inactive, noPageBreakBelow, noPageBreakInGroup, nohscroll, preserveAspect, groupLikeTitle, center, M7, bold, L3, e8, 24, "Times"; ; fontset = subtitle, inactive, noPageBreakBelow, noPageBreakInGroup, nohscroll, preserveAspect, groupLikeTitle, center, M7, bold, L1, e6, 18, "Times"; ; fontset = subsubtitle, inactive, noPageBreakBelow, noPageBreakInGroup, nohscroll, preserveAspect, groupLikeTitle, center, M7, bold, L1, e6, 14, "Times"; ; fontset = section, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, grayBox, M22, bold, L1, a20, 14, "Times"; ; fontset = subsection, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, blackBox, M19, bold, L1, a15, 12, "Times"; ; fontset = subsubsection, inactive, noPageBreakBelow, nohscroll, preserveAspect, groupLikeSection, whiteBox, M18, bold, L1, a12, 10, "Times"; ; fontset = text, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 14, "Times"; ; fontset = smalltext, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 10, "Times"; ; fontset = input, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeInput, M42, N23, bold, L1, 12, "Courier"; ; fontset = output, output, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L-5, 12, "Courier"; ; fontset = message, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, R65535, L1, 12, "Courier"; ; fontset = print, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L1, 12, "Courier"; ; fontset = info, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeOutput, M42, N23, L1, 12, "Courier"; ; fontset = postscript, PostScript, formatAsPostScript, output, inactive, noPageBreakInGroup, nowordwrap, preserveAspect, groupLikeGraphics, M7, l34, w415, h422, L1, 10, "Courier"; ; fontset = name, inactive, noPageBreakInGroup, nowordwrap, nohscroll, noKeepOnOnePage, preserveAspect, M7, italic, B65535, L1, 10, "Times"; ; fontset = header, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 10, "Times"; ; fontset = leftheader, 12; fontset = footer, inactive, nohscroll, noKeepOnOnePage, preserveAspect, center, M7, L1, 12; fontset = leftfooter, 12; fontset = help, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 10, "Times"; ; fontset = clipboard, inactive, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12; fontset = completions, inactive, nowordwrap, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12; fontset = special1, inactive, nowordwrap, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12; fontset = special2, inactive, nowordwrap, nohscroll, noKeepOnOnePage, preserveAspect, center, M7, L1, 12; fontset = special3, inactive, nowordwrap, nohscroll, noKeepOnOnePage, preserveAspect, right, M7, L1, 12; fontset = special4, inactive, nowordwrap, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12; fontset = special5, inactive, nowordwrap, nohscroll, noKeepOnOnePage, preserveAspect, M7, L1, 12; paletteColors = 128; currentKernel; ] :[font = title; inactive; dontPreserveAspect] Appendix: Supplementary Programs :[font = section; inactive; Cclosed; dontPreserveAspect; startGroup] Chapter 1 :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Wilson's Theorem Primality Test :[font = input; dontPreserveAspect; endGroup] WilsonPrimeTest[n_] := Block[{i = 1}, Nest[Mod[# i++, n]&, 1, n-1] == n-1] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] primeQ :[font = input; dontPreserveAspect] BeginPackage["primeQ`"]; :[font = input; wordwrap; dontPreserveAspect] primeQ::usage = "primeQ[n] modifies PrimeQ by first testing for divisibility by any of the first 1,000 primes; it should be used instead of PrimeQ when numbers being tested are large"; :[font = input; dontPreserveAspect] Begin["`Private`"]; :[font = input; dontPreserveAspect; endGroup] smallprimes = Prime[Range[1000]]; testforsmalldivisors[x_] := (n = 1; While[n < 1000 && Mod[x, smallprimes[[n]]] != 0, n++]; n == 1000) primeQ[n_] := (n <= 7919 || testforsmalldivisors[n]) && PrimeQ[n] End[]; EndPackage[]; :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] encode & decode (for version 1.2) :[font = input; dontPreserveAspect; endGroup] encode[message_] := ToExpression[ Apply[StringJoin, Map[ToString[ToASCII[#] - 22]&, Characters[message] /. "z" -> "["]]] decode[codedmessage_Integer] := Block[{out}, out = Characters[codedmessage]; out = Table[out[[{2i - 1, 2i}]], {i, Length[out]/2}]; out = Map[Apply[StringJoin, #]&, out]; out = Map[FromASCII[ToExpression[#] + 22]&, out]; Apply[StringJoin, (out /. "[" -> "z" )]] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Version 2.0 versions of encode and decode :[font = input; dontPreserveAspect; endGroup] encode[message_] := ToExpression[ StringJoin @@ ToString /@ (ToCharacterCode[message]-22 /. 100 -> 71) ] decode[codedmessage_] := StringJoin @@ FromCharacterCode /@ (22 + (ToExpression[StringJoin @@ #] & /@ Partition[Characters[ToString @ codedmessage], 2, 2] /. 71 -> 100)) :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] PrimePi (due to Ilan Vardi) :[font = input; dontPreserveAspect; endGroup] PrimePi::usage = "PrimePi[x] returns the number of primes less than or equal to x"; PrimePi[x_] := PrimePi[Floor[x]] /; Not[IntegerQ[x]] PrimePi[n_Integer] := 0 /; n <= 1 PrimePi[4] = 2; PrimePi[n_Integer] := Block[{t, diff, m}, t = Round[N[Log[n]]]; m = Round[N[LogIntegral[n]]]; diff = n - Prime[m]; While[ Abs[diff] > t^2 / (1.16)^2, m += Round[ diff / t]; diff = n - Prime[m];] If[diff < 0, While[n < Prime[m], m--;], While[n >= Prime[m], m++;]; m--]; Return[m]] /; n > 1 && n != 4 Attributes[PrimePi] = Listable; :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] RiemannR :[font = input; dontPreserveAspect; endGroup] RiemannR::usage = "Riemann[x] returns the value of Riemann's prime-approximating function R(x)."; rlist = Map[1/(#! # N[Zeta[# + 1]]) &, Range[100]] RiemannR[x_] := Round[1 + (Floor[N[Log[x]] ^ Range[100]]) . rlist] RiemannR[0] = 0; Attributes[RiemannR] = {Listable}; :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] PrimePiModn :[font = input; dontPreserveAspect; endGroup] PrimePiModn[x_, n_] := Block[ {i, temp, modlist = Mod[Prime[Range[x]], n]}, Table[(temp = Flatten[Position[modlist, i]]; Transpose[{Prime[temp], Range[Length[temp]]}]), {i, n - 1}]] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] PrimePiMod4 exercise :[font = input; dontPreserveAspect; endGroup; endGroup] ListPlot[Accumulate[Plus, Mod[Prime[Range[3000]], 4] - 2], PlotJoined->True]; :[font = section; inactive; Cclosed; dontPreserveAspect; startGroup] Chapter 2 :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Trochoid Animation Generator :[font = input; dontPreserveAspect; endGroup] r = Input["This program traces the path of a point on a spoke of a rolling wheel. The radius of the rolling wheel is 1. Enter the radius of the spoke."] rr = Max[r - 1, 0]; pi = N[Pi]; cycloid[t_] := {t - r Sin[t], 1 - r Cos[t]}; rotate = {{Cos[theta],Sin[theta]}, {-Sin[theta],Cos[theta]}} spokepoints := Table[{t,1} + rotate . (cycloid[t] - {t,1}), {theta, 0, 2 pi, pi/4}] spokes := {Thickness[.004], Map[Line[{{t,1}, #}] &, spokepoints]} baseline = Line[{{-2,0}, {10.5,0}}]; dots = {PointSize[.006]}; Do[ AppendTo[dots, Point[cycloid[t]]]; Show[Graphics[{{GrayLevel[.5], Disk[{t,1}, 1]}, Line[{{t,1}, cycloid[t]}], (* spoke *) spokes, baseline, Circle[{t,1}, 1], dots}], AspectRatio->(3+2 rr)/(12+2 rr), PlotRange->{{-1.5-rr, 10.5+rr}, {-.5-rr, 2.5+rr}}], {t, 0., 2 pi 59/40, 2 pi/40}] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Epicycloid and Hypocycloid Animation Generator :[font = input; dontPreserveAspect; endGroup] Clear[rotate] r = Input["This program traces the path of a point on the circumference of a circle rolling around (positive r) or inside (negative r) a circle of radius 1. Enter the radius of the rolling circle."] rr = Max[1, 1 + 2 r]; pi = N[Pi] rotate[point_, theta_] := {{Cos[theta], Sin[theta]}, {-Sin[theta], Cos[theta]}} . point fixedcircle = {Thickness[.006], Circle[{0, 0}, 1]} rollingdisk := {GrayLevel@.5, Disk[center, Abs[r]]} rollingcircle := Circle[center, Abs[r]] dots = {PointSize[.006]} Do[ center = rotate[{-1-r, 0}, theta]; spoke1 = rotate[{1,0}, theta (1+1/r)]; spokeline1 = {center + r spoke1, center - r spoke1}; spoke2 = rotate[{0,1}, theta (1+1/r)]; spokeline2 = {center + r spoke2, center - r spoke2}; spokes = {Thickness[.007], Line[spokeline1], GrayLevel[1], Line[spokeline2]}; AppendTo[dots, Point[center + r spoke1]]; Show[Graphics[ {rollingdisk, spokes, fixedcircle, rollingcircle, dots}], AspectRatio->1, PlotRange->{{-rr - .03, rr + .03}, {-rr - .03, rr + .03}}], {theta, 0., 2 pi , 2 pi/60}] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Rolling Polygons :[font = input; dontPreserveAspect; endGroup] n = Input["How many sides does the rolling polygon have?"]; pi = Pi//N; cr = Csc[pi/n]; (* circumradius of polygon *) ir = Cot[pi/n]; (* inradius of polygon *) a = ir ArcSinh[1/ir]; dots = {PointSize[.0018]}; rotate[p_, t_, center_] := center + {{Cos[t], Sin[t]}, {-Sin[t], Cos[t]}} . (p - center) plot = Table[Plot[cr - ir Cosh[(x - 2 a i)/ir], {x, -a + 2 a i, a + 2 a i}, DisplayFunction->Identity, PlotStyle->Thickness[.001]], {i, 0, 5}]; startpolygon = Map[# + {0, cr} &, Map[{Re[#], Im[#]} &, cr N[Exp[I Range[-Pi/2 + Pi/n, 3 Pi/2 + Pi/n, 2 Pi/n]]]]]; polygon[t_] := Block[{loop = N[Floor[(t + a)/(2a)]]}, (* loop tells which loop the square is on *) Map[{t,0} + rotate[#, loop 2pi/n + ArcTan[Sinh[(t-loop 2 a)/ir]], {0,cr}]&, startpolygon]] Do[ newpolygon = polygon[t]; AppendTo[dots, Point[newpolygon[[n]]]]; Show[ Graphics[{PointSize[.008], Thickness[.001], {GrayLevel[.75], Polygon[newpolygon]}, Line[newpolygon], Map[Line[{{t, cr},#}] &, newpolygon], Line[{{-cr, cr}, {11.5 a,cr}}], Point[{t, cr}], Point[newpolygon[[n]]], dots}], plot, PlotRange->{{-cr, 11 a + cr}, {-.2, 2 cr}}, AspectRatio->(2 cr +.2) / (11 a + 2 cr), DisplayFunction->$DisplayFunction, Ticks->None, Axes->None], {t, 0., 11 a, a/10}]; :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Reuleaux Triangle Rotating inside a Square :[font = input; dontPreserveAspect; endGroup] Clear[t, rotate] centroid = {1 - 2 Sqrt[3.]/3, 0} choose[t_] := Floor[.5 (Floor[t / (Pi/6)] + 3)] rotate[point_, t_, center_] := center + {{Cos[t], -Sin[t]}, {Sin[t], Cos[t]}} . (point - center) arc[{center_, startangle_}] := Block[{t}, Table[center + 2 {Cos[t], Sin[t]}, {t, startangle, N[startangle + Pi/3], 0.08}]] verts = {{1 - Sqrt[3.], -1}, {1, 0}, {1 - Sqrt[3.], 1}} angles = {Pi/6, 5 Pi/6, 9 Pi/6}; dots = {}; Do[tempverts = Map[rotate[#, t, centroid] &, verts]; translate = 2 Mod[{choose[t], choose[t-Pi/2]}, 2] - 1 - {tempverts[[Mod[choose[t], 3] + 1, 1]], tempverts[[Mod[choose[t - Pi/2], 3] + 1, 2]]}; newverts = Map[ # + translate &, tempverts]; reuleaux = Flatten[Map[arc, Transpose@{newverts, t + angles}], 1]; AppendTo[reuleaux, First@reuleaux]; dots = Join[dots, {newverts[[3]], centroid + translate, rotate[{0,0}, t, centroid] + translate}]; Show[Graphics[{{GrayLevel[.75], Polygon[reuleaux]}, Line[reuleaux], Line[Append[newverts, First[newverts]]], PointSize[.02], Point /@ { newverts[[3]], centroid + translate, rotate[{0, 0}, t, centroid] + translate}, PointSize[.007], Point /@ dots, {Thickness[.001], Line[{{1,1},{1,-1},{-1,-1},{-1,1},{1, 1}}]}, Map[Line[{centroid + translate, #}] &, newverts]}], PlotRange->{{-1.05, 1.05}, {-1.05, 1.05}}, AspectRatio->1], {t, 0, 1.99 Pi, Pi/12}] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Cycloid as tautochrone :[font = input; Cclosed; dontPreserveAspect; startGroup] pi = Pi//N; box = {Thickness[.003], Line[{{2 pi + .2, .2}, {2 pi + .2, -2.2}, {-.2, -2.2}, {-.2, .2}, {2 pi + .2, .2}}]}; cycloid[t_] := {t - Sin[t], -1 + Cos[t]}; cycloidplot = ParametricPlot[cycloid[s], {s, 0, 2 pi}, Axes->None, PlotStyle->Thickness[.0005], DisplayFunction->Identity]; options = {PlotRange->{{-.2, 2 pi + .2}, {-2.2, .2}}, AspectRatio->2.4/(2 pi + .4), DisplayFunction->$DisplayFunction}; stopline = {Thickness[.002], Line[{{pi, -2.2}, {pi, -1.8}}]}; showframe := Show[Graphics[{Thickness[.0004], Disk[{2 pi-posn[[1]], posn[[2]]},radius], Disk[posn1, radius], stopline}], cycloidplot, options]; showstill := Show[cycloidplot, Graphics[{box, Thickness[.0004], list, stopline}], options]; g = 980 2 pi /15; initial = Input["Enter initial t-value for start time on auxiliary cycloid: 0 < t < 3"] //N; choice = InputString["Do you want the whole movie [movie] or just a single image [still]?"] inc = Input["Enter the time increment, which controls the number of frames generated. Values between 1/100 and 1/400 are best."] //N t = time = 0.; radius = .07; list = {}; (* MAIN LOOP BEGINS HERE. The variables posn and t refer to the case of the bead starting at the top of the cycloid; posn1 and t1 refer to the userÕs starting position. *) While[t < pi - radius/2, t = Sqrt[g]*time; t1 = 2 ArcCos[Cos[initial/2.] * Cos[Sqrt[g] time / 2]]; posn = cycloid[Min[t, pi - radius/2]]; posn1 = cycloid[Min[t1, pi - radius/2]]; time += inc; If[choice == "movie", showframe, list = Join[list, {Disk[{2 pi - posn[[1]], posn[[2]]}, radius], Disk[posn1, radius]}]]]; If[choice != "movie", showstill] :[font = postscript; PostScript; formatAsPostScript; output; inactive; dontPreserveAspect; pictureLeft = 2; pictureWidth = 456; pictureHeight = 163; endGroup; endGroup] %! %%Creator: Mathematica %%AspectRatio: 0.35911 MathPictureStart % Scaling calculations 0.02993 0.14963 0.32918 0.14963 [ [ 0 0 0 0 ] [ 1 0.35911 0 0 ] ] MathScale % Start of Graphics 1 setlinecap 1 setlinejoin newpath %%Object: Graphics [ ] 0 setdash 0 setgray gsave gsave grestore grestore 0 0 moveto 1 0 lineto 1 0.35911 lineto 0 0.35911 lineto closepath clip newpath 0 setgray gsave gsave gsave gsave 0.0005 setlinewidth 0.02993 0.32918 moveto 0.02993 0.32918 lineto 0.02993 0.32916 lineto 0.02993 0.32914 lineto 0.02993 0.3291 lineto 0.02993 0.329 lineto 0.02993 0.32886 lineto 0.02994 0.32868 lineto 0.02995 0.32846 lineto 0.02998 0.3279 lineto 0.03003 0.32719 lineto 0.03011 0.32631 lineto 0.03037 0.32409 lineto 0.0308 0.32124 lineto 0.03142 0.31779 lineto 0.03346 0.30914 lineto 0.03677 0.29826 lineto 0.04164 0.28536 lineto 0.05703 0.25437 lineto 0.08126 0.21828 lineto 0.11533 0.17956 lineto 0.15961 0.14083 lineto 0.21373 0.10474 lineto 0.24419 0.08847 lineto 0.27668 0.07375 lineto 0.31098 0.06085 lineto 0.34684 0.04997 lineto 0.38398 0.04132 lineto 0.40294 0.03787 lineto 0.4221 0.03502 lineto 0.44143 0.0328 lineto 0.45114 0.03192 lineto 0.46088 0.03121 lineto 0.47064 0.03065 lineto 0.47553 0.03043 lineto 0.48042 0.03025 lineto 0.48531 0.03011 lineto 0.48776 0.03005 lineto 0.49021 0.03001 lineto 0.49266 0.02997 lineto 0.4951 0.02995 lineto 0.49755 0.02993 lineto 0.5 0.02993 lineto 0.50245 0.02993 lineto 0.5049 0.02995 lineto 0.50734 0.02997 lineto 0.50979 0.03001 lineto 0.51469 0.03011 lineto 0.51958 0.03025 lineto 0.52447 0.03043 lineto Mistroke 0.52936 0.03065 lineto 0.53912 0.03121 lineto 0.54886 0.03192 lineto 0.55857 0.0328 lineto 0.5779 0.03502 lineto 0.59706 0.03787 lineto 0.61602 0.04132 lineto 0.65316 0.04997 lineto 0.68902 0.06085 lineto 0.72332 0.07375 lineto 0.78627 0.10474 lineto 0.84039 0.14083 lineto 0.88467 0.17956 lineto 0.91874 0.21828 lineto 0.94297 0.25437 lineto 0.95168 0.27064 lineto 0.95836 0.28536 lineto 0.96323 0.29826 lineto 0.96654 0.30914 lineto 0.96858 0.31779 lineto 0.9692 0.32124 lineto 0.96963 0.32409 lineto 0.96989 0.32631 lineto 0.96997 0.32719 lineto 0.97002 0.3279 lineto 0.97005 0.32846 lineto 0.97006 0.32868 lineto 0.97007 0.32886 lineto 0.97007 0.329 lineto 0.97007 0.32906 lineto 0.97007 0.3291 lineto 0.97007 0.32914 lineto 0.97007 0.32916 lineto 0.97007 0.32918 lineto 0.97007 0.32918 lineto Mfstroke grestore grestore grestore gsave gsave gsave 0.003 setlinewidth 1 0.35911 moveto 1 0 lineto 0 0 lineto 0 0.35911 lineto 1 0.35911 lineto stroke grestore gsave 0.0004 setlinewidth 0.97007 0.32918 moveto 0.97007 0.32918 0.01047 0 365.73 arc fill 0.19313 0.11729 moveto 0.19313 0.11729 0.01047 0 365.73 arc fill 0.51047 0.03002 moveto 0.51047 0.03002 0.01047 0 365.73 arc fill 0.48953 0.03002 moveto 0.48953 0.03002 0.01047 0 365.73 arc fill grestore gsave 0.002 setlinewidth 0.5 0 moveto 0.5 0.05985 lineto stroke grestore grestore grestore grestore % End of Graphics MathPictureEnd :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Cycloid as brachistochrone :[font = input; dontPreserveAspect; endGroup; endGroup] Clear[i, t, tt] pi = Pi//N; norm[p_] := Sqrt[p.p]; cycloid[t_] := {t - Sin[t], -1 + Cos[t]}; cycloidplot = ParametricPlot[cycloid[tt], {tt, pi, 2 pi}, {Axes->None, PlotStyle->Thickness[.0005], DisplayFunction->Identity}]; degree = Input["Enter the degree (between 1 and 15) of the polynomial curve connecting the two points."]; n = -1 + Input["Enter the number of frames desired. Values between 10 and 30 are best."]//N; g = 980 2 pi /15; inc = pi/Sqrt[g]/(3 n); f[t_] := {t, -2 (t/pi) ^ (1/degree)}; directionf[t_] := Block[{v}, If[(degree > 1 && t != 0) || degree == 1, (v = f'[t])/norm[v], {0,-1}]]; choice = InputString["Do you want the whole movie [movie] or just a single image [still]?"]; options = {PlotRange->{{-.2, 2 pi + .2}, {-2.2, .2}}, AspectRatio->2.4/(2 pi + .4), DisplayFunction->$DisplayFunction}; box = {Thickness[.003], Line[{{2 pi + .2, .2}, {2 pi + .2, -2.2}, {-.2, -2.2}, {-.2, .2}, {2 pi + .2, .2}}]}; stopline = {Thickness[.002], Line[{{pi, -2.1}, {pi, -1.9}}]}; showstill := Show[cycloidplot, plotf, Graphics[{box, list, stopline}], options]; showframe := Show[cycloidplot, plotf, Graphics[{Thickness[.0004], Disk[{2 pi-posn[[1]], posn[[2]]}, radius], Disk[posn1, radius], stopline}], options]; newposn[{x_,y_}, t1_] := FindRoot[ f[tt][[2]] - y == -(tt - x) Divide @@ d, {tt, Max[t1, 10.^-60]}, MaxIterations->100]; plotf = ParametricPlot[f[t], {t, 0, pi}, Axes->None, PlotStyle->Thickness[.0005], DisplayFunction->Identity]; radius = .07; posn = posn1 = {0,0}; t1 = 0; d = directionf[0]; list = {Thickness[.0004], Disk[{2 pi, 0}, radius], Disk[posn1, radius]}; If[choice == "movie", showframe]; Do[ t = pi i / (3 n); {x, y} = posn1 + d*inc*(Sqrt[-2 g posn1[[2]]] - g d[[2]] inc/2); (* Preceding line computes the point on the tangent line that the bead would roll to; the following line sends this position to newposn to get the corresponding point on the curve *) t1 = If[degree == 1, x, tt /. newposn[{x,y}, t1]]; d = directionf[t1]; posn = cycloid[Min[t, pi - radius/2]]; posn1 = f[Min[t1, pi - radius]]; If[Mod[i, 3] == 0, If[choice == "movie", showframe, list = Join[list, {Disk[{2 pi - posn[[1]], posn[[2]]}, radius], Disk[posn1, radius]}]]], {i, 3 n}]; If[choice != "movie", showstill]; :[font = section; inactive; Cclosed; dontPreserveAspect; startGroup] Chapter 4 :[font = subsection; inactive; dontPreserveAspect] Cantor Set via Intervals :[font = input; dontPreserveAspect] spawn1[{a_, b_}] := Block[{w = (b - a)/3}, {{a, a + w}, {b - w, b}}] /; (Head[a] == Head[b] == Rational) || a==0 || b == 1; spawn1[intervals_List] := Flatten[Map[spawn1, intervals],1] /; Length[intervals[[1]]] > 1; lines[n_Integer, l_List] := Map[Line[{{#[[1]], n + .5}, {#[[2]], n + .5}}] &, l]; level = 2; Show[Graphics[{Thickness[.0005], Map[lines[level--, #] &, NestList[spawn1, {{0,1}}, level]]}], Ticks->{N[Range[0,1,1/9],3], None}, Axes->Automatic] :[font = subsection; inactive; dontPreserveAspect] Animation of Cantor functions :[font = input; dontPreserveAspect] Clear[spawn, f]; 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]; connect[{a_, b_}] := Block[{y = f[a]}, Line[{{a, y}, {b, y}}]]; f[0] = 0; f[x_] := f[3x]/(r + 1) /; x <= 1/3; f[x_] := 1 - r f[1 - x] /; 2/3 <= x; CantorFunction[parameter_] := (r = parameter; Show[Graphics[{Thickness[.001], Map[connect, removedintervals[6]]}], Axes->{0, 0}, PlotRange->{{0,1}, {0,1}}, Ticks->{N[Range[0., 1, 1/9], 3], Automatic}]); Attributes[CantorFunction] = {Listable} CantorFunction[{.01, .1, .3, .5, 1, 4, 10, 100}]; :[font = subsection; inactive; dontPreserveAspect] The Sierpinski Triangle :[font = input; dontPreserveAspect; endGroup] Clear[spawn] spawn[{a_, b_, c_, a_}] := Block[{ab = (a+b)/2, bc = (b+c)/2, ac = (a+c)/2}, {{a, ab, ac, a}, {b, bc, ab, b}, {c, ac, bc, c}}]; spawn[l_List] := Flatten[Map[spawn, l], 1] /; Length[l[[1]]] > 1; Show[Graphics[{Thickness[.001], Line /@ Nest[spawn, {{0, 0}, {.5, Sqrt[3.]/2}, {1, 0}, {0, 0}}, 7]}], AspectRatio->Sqrt[3]/2] :[font = section; inactive; Cclosed; dontPreserveAspect; startGroup] Chapter 5 :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Chaos Game Using a Universal Sequence :[font = input; dontPreserveAspect; endGroup] top = {.5, Sqrt[3.]/2} f1[x_] := .5 x; f2[x_] := .5(x + {1,0}); f3[x_] := .5(x + top) list = Append[Table[0, {8-1}], 1]; Show[Graphics[Map[Point, NestList[( list[[1]] = Mod[list[[1]] - list[[6]], 3]; list = RotateLeft[list]; {f1[#], f2[#], f3[#]}[[Last[list]+1]] ) &, {0,0}, 3^8-1]]], PlotRange->{{0, 1}, {0, .87}}, Axes->Automatic, AspectRatio->.87, Ticks->{{.25, .5, .75, 1}, {.25, .5, .75, 1}}] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Filled Julia Sets for the Real Quadratic Map :[font = input; dontPreserveAspect; endGroup] FilledJuliaSetReal::usage = "FilledJuliaSetReal[r, meshx, meshy, x, y, iters] shows the filled-in Julia set for r*z*(1-z), by using the Escape Time Algorithm on a meshx-by-meshy grid over the rectangle where x varies from 0 to x and y from 0 to y. Fourfold symmetry is then used to generate an image where x ranges from 0 to 1 and y from -y to y. Thus x should be set to 0.5 to get the entire set. The number of iterations used is given by iters. A Print statement in the subroutine orbitcheckreal acts as a progress checker." orbitcheckreal[z_, r_, iters_] := (s = z; i = 0; If[Re[z] != xold, xold=Re[z]; Print[xold]]; While[++i <= iters && Abs[s *= r (1-s)] < r0]; If[i != iters + 1, {}, a = Re[z]; b = Im[z]; Flatten[Outer[List, {a, 1-a}, {b, -b}], 1]]) FilledJuliaSetReal[r_, meshx_, meshy_, x_, y_, iters_] := (r0 = 1 + 1/Abs[r]; xold = 0; Show[Graphics[{PointSize[.002], Point /@ Flatten[ Outer[orbitcheckreal[#1 + I #2, N[r], iters]&, Range[0, x, x/meshx], Range[0, y, y/meshy]], 2]}], PlotRange->All, AspectRatio->y/x]) :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Tree Animations :[font = input; dontPreserveAspect] Tree::usage = "Tree[depth, length, angle, ratio] uses the recursive subroutine, branch, to generate and display a binary tree with small round leaves at the end of each terminal branch. " depth::usage = "The depth of the recursion of the tree function. The tree will have 2^depth terminal branches and leaves." length::usage = "The original length of the trunk of the tree in the tree function. This argument determines the overall size of the tree." angle::usage = "The whole angle, in degrees, between paired subbranches of a single parent branch. This argument determines how much the tree spreads out." ratio::usage = "The reduction ratio which determines how much smaller a subbranch is than its parent branch."; Attributes[Tree] = Listable branch[depth_Integer, length_, angle_, ratio_] := ( forward[length]; left[angle/2]; If[depth > 0, branch[depth - 1, length*ratio, angle, ratio], AppendTo[points, X]]; right[angle]; If[depth > 0, branch[depth - 1, length*ratio, angle, ratio]] ; left[angle/2]; back[length];) Tree[depth_Integer, length_, angle_, ratio_] := ( initialize[]; points = {}; (* to be a list of dots at the branch-ends *) left[90]; If[depth > 0, branch[depth, length, angle, ratio]]; finished; Show[Graphics[{Line[path], PointSize[0.03], Point /@ points}], PlotRange->All]) :[font = input; dontPreserveAspect] Tree[5, 1, Range[45, 90, 2.5], .8]; (* Tree pushups *) :[font = text; inactive; dontPreserveAspect] The following command generates a nice animation if the PlotRange setting in the definition of Tree is changed to {{-2.4, 2.4}, {0, 3.5}}. The sequence of frames should be animated in the forward direction ;[s] 7:0,0;56,1;65,2;95,3;99,4;115,5;138,6;207,-1; 7:1,13,10,Times,0,14,0,0,0;1,11,9,Courier,1,14,0,0,0;1,13,10,Times,0,14,0,0,0;1,11,9,Courier,1,14,0,0,0;1,13,10,Times,0,14,0,0,0;1,11,9,Courier,1,14,0,0,0;1,13,10,Times,0,14,0,0,0; :[font = input; dontPreserveAspect; endGroup; endGroup] Tree[5, 1, 90, Range[.05, .8, .05]]; (* Tree growth *) :[font = section; inactive; Cclosed; dontPreserveAspect; startGroup] Chapter 6 ÑÑ UltraTurtle and RecursiveTurtle3D placed in Chapter 6 file. :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Sierpinski and original Peano space-fillers :[font = input; dontPreserveAspect; endGroup] SierpinskiCurve[n_] := Show[RecursiveTurtle[ {"X"-> "XF-F+F-XF+F+XF-F+F-X"}, "F+XF+F+XF", n, 90, side = Sqrt[2.] 2^-(n+2), ({temp = ((2^(n-1)6 + 2^(n) -2) side /Sqrt[2])}; (*stepsize*) {(temp+1)/2 - side/Sqrt[2], (1-temp)/2}), (*start direction*) Sqrt[{2,2}]/2, (*path length*) (4^(n+2) - 1)/3 ], Axes->{0, 0}, PlotRange->{{0, 1}, {0, 1}}, Ticks->{{.25, .5, .75, 1}, {.25, .5, .75, 1}}] OriginalPeanoCurve[n_] := Show[RecursiveTurtle[ {"X"-> "XFYFX+F+YFXFY-F-XFYFX", "Y"->"YFXFY-F-XFYFX+F+YFXFY"}, "X", n, 90, 3^-n, {3, 3}^-n / 2, 9^n], Axes->{0,0}, PlotRange->{{0, 1}, {0, 1}}, Ticks->{{.33, .67, 1}, {.33, .67, 1}}] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Illustration of a Function from the Interval Onto the Square :[font = input; dontPreserveAspect] CantorBijection::usage = "CantorBijection[t] gives the value of a bijection from the unit interval to the unit square based on even- and odd-indexed binary digits. The input and outputs are base-10 numbers in the unit interval and unit square, respectively."; CantorImage::usage = "CantorImage[n] gives the path in the unit square obtained by using straight lines to connect the values of CantorBijection on {0, 1/n, 2/n, . . . , 1}."; SetAttributes[{CantorBijection, CantorImage}, Listable]; base2digits[t_] := Block[{tt = N[t]}, Map[If[2. ^ -# <= tt, tt -= 2. ^ -#; 1, 0] &, Range[30]]] CantorBijection[t_] := ( temp = Flatten@Position[base2digits[t], 1]; temp = {Select[temp, !EvenQ[#]&], Select[temp, EvenQ]}; {Plus @@ (2^((-temp[[1]]-1)/2)), Plus @@ (2^(-temp[[2]]/2)) } //N ) CantorBijection[1] = {0,1}; CantorImage[n_] := Show[Graphics[Line[CantorBijection[Range[0, 1, 1/n]]]], PlotRange->{{0, 1}, {0, 1}}, AspectRatio->1, Axes->Automatic, AxesStyle->GrayLevel[1], Ticks->Range[{0, 0}, {1, 1}, {.25, .25}]] :[font = input; dontPreserveAspect; endGroup; endGroup] CantorImage[13] (* n should be a power of 2 *) :[font = section; inactive; Cclosed; dontPreserveAspect; startGroup] Chapter 8 :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Test for Whether Integers Are Pairwise Relatively Prime :[font = input; dontPreserveAspect; endGroup] PairwiseCoprime::usage = "PairwiseCoprime[mlist] returns True if the integers in mlist are pairwise relatively prime, False otherwise." PairwiseCoprime[mlist_List] := Length[mlist] == Length[Union[mlist]] && And @@ Flatten[Outer[GCD[#1,#2] == 1 || #1 == #2 &, mlist, mlist]] :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] The Splitting Algorithm for Egyptian Fractions :[font = input; dontPreserveAspect] EgyptianSplitting::usage = "EgyptionSplitting[q] uses the splitting method to write q as a sum of distinct unit fractions" EgyptianSplitting[q_] := ( list = Table[Denominator[q], {Numerator[q]}]; While[Length[list] != Length[Union[list]], list = Sort[Flatten[ Map[(t = list[[#]]; If[{#} != First[Position[list, t]], {t+1, t(t+1)}, t]) &, Range[Length[list]]]]]]; list) :[font = input; Cclosed; dontPreserveAspect; startGroup] EgyptianSplitting[4/7] :[font = output; output; inactive; dontPreserveAspect; endGroup; endGroup] {7, 8, 9, 10, 56, 57, 58, 72, 73, 90, 3192, 3193, 3306, 5256, 10192056} ;[o] {7, 8, 9, 10, 56, 57, 58, 72, 73, 90, 3192, 3193, 3306, 5256, 10192056} :[font = subsection; inactive; Cclosed; dontPreserveAspect; startGroup] Code to Generate Table 8.2 :[font = input; dontPreserveAspect; endGroup; endGroup] Map[Length[Digits[Last[#]]]&, EgyptianFraction[(10 ^ Range[9] - 1) / 10 ^ Range[9]]] :[font = section; inactive; dontPreserveAspect] Chapter 9 routines placed in Chapter 9 file ^*)