(* You should be able to copy and paste this text file directly into Mathematica. These notes will be interpreted as commented code. *) (* Bugs in the code have now been fixed -- 9:35am, Aug 16th -- E.D. *) (* For problem 1 *) f1a[z_] := z + z^2 (z-I)(z+1) f1b[z_] := z^2 - 3z + 5 f1c[z_] := I z^2 + z + I/4 f1d[z_] := 1 + (4I) / (z + 2 - 3I) (* For problem 2 *) f2a[z_] := 1 - I + I z f2b[z_] := z^2 + I f2c[z_] := 1 - (3/2) z^2 - (1/2) z^3 f2d[z_] := 3z + 4/z f2e[z_] := z^2 (* For problem 3 *) f3a[z_] := z^2 - 0.4 - 0.1 I f3b[z_] := z^2 + 0.2 - I f3c[z_] := z^2 + 1. I f3d[z_] := z^2 - 0.53 + 0.6 I (* For problem 4 *) f4a[z_] := z^3 + 0.4z - I f4b[z_] := (z^3 - 1.) / (z + I) f4c[z_] := z - (z^3 - 1.) / (3 z^2) (* For problem 5 *) c5a := -0.1 + 0.7 I c5b := 0.3 + 0.55 I c5c := -0.5 + 0.55 I c5d := -0.625 + 0.425 I c5e := 0.05 + 0.63 I c5f := 0.375 + 0.27 I (* For problem 6 *) c6a := -0.21 + 0.8 I c6b := 0.13 - 0.63 I c6c := 0.390 + 0.232 I c6d := -0.547 - 0.557 I (* Here is code for invoking Mathematica's Julia set plotting algorithms, for problems 3-5 *) f3a[z_] := z^2 - 0.4 - 0.1 I JuliaSetPlot[f3a[z], z, PlotLegends -> Automatic, ImageResolution -> 200] (* This will plot the filled Julia set for a function with default escape-time colorings, and also display the legend for comparing colorings to the number of iterates required for escape *) JuliaSetPlot[f3a[z], z, ColorFunction -> None, ImageResolution -> 200] (* This will plot the Julia set for a function, with default blue points *) JuliaSetPlot[f3a[z], z, ColorFunction -> (If[#3 > 9/10, RGBColor[0.5, 0.19, 0.85], White] &), ImageResolution -> 200] (* This will plot the filled Julia set for a function, with points colored in purple *) (* This method will not be effective for plotting Julia sets that are very sparse *) ListPlot[JuliaSetPoints[f3a[z], z, "ClosenessTolerance" -> 0.005], PlotStyle -> Directive[RGBColor[0.5, 0.19, 0.85], PointSize[0.005]], Axes -> False, Frame -> True] (* This will plot the Julia set for a function with purple points. To increase the resolution of the plot, decrease the value of "ClosenessTolerance" -- note that this method may take a long time if the Julia set is complicated, so it is recommended you first try plotting the set with ClosenessTolerance 0.005 and then lower the value as needed *) (* Here is code for plotting a Julia set for a quadratic map using backwards iteration (the chaos game procedure) for problems 3 and 5 *) h1[c_][z_] := N[Sqrt[z - c]] h2[c_][z_] := N[-Sqrt[z - c]] hch[c_][z_] := RandomChoice[{1, 1} -> {h1[c], h2[c]}][z] (* This defines the two components of the inverse map and a function which randomly chooses one of the maps to apply *) julpoints[c_, n_] := Drop[NestList[hch[c], 0.01, n + 1000], 1000] (* This will compute n+1000 backwards iterates then delete the first 1000 terms, to clean up any slow convergence at the beginning *) complexlisttopointlist[list_] := Transpose[{Re[list], Im[list]}] julplot[c_, n_, pointsize_] := ListPlot[complexlisttopointlist[julpoints[c,n]], Axes -> True, AxesLabel -> {Re, Im}, PlotStyle -> Directive[RGBColor[0.64, 0.21, 0.76], PointSize[pointsize]], AspectRatio -> 1, PlotLabel -> StringJoin["Julia set for q(z)=z^2+(", ToString[c], ") plotted with ", ToString[n], " total points"]] (* This will plot n points on the Julia set for q(z)=z^2+c, labeled with the parameter c and the number of points. It may be advisable to change the point size if more points are plotted to create a higher-resolution picture. *) julplot[-0.5, 20000, 0.0010] (* A plot of the Julia set for z^2 - 0.5 with 20000 points and point size of 0.0010. *) (* Here is code for looking at the critical orbit for a quadratic map, for problems 3, 5, 6 *) q[c_][z_] := z^2 + c c3a := 0.4 - 0.1 I (* Defining functions for ease of use *) NestList[q[c3a], 0.0, 10] (* This will compute the first ten terms on the critical orbit, which you can use to get an idea of the behavior *) Drop[NestList[q[c3a], 0.0, 1030],1000] (* This will compute the 1001st through 1030th terms on the orbit, which will generally give a very clear picture of the asymptotic behavior *) (* Here is code for plotting the Mandelbrot set, for problems 5 and 6 *) MandelbrotSetPlot[{-0.4 + 0.5 I, 0.1 + 1.2 I}, ImageResolution->500] (* This plots the portion of the Mandelbrot set ranging from the lower left corner -0.4+0.5i to upper right corner 0.1+1.2i *) mandelbrotbox[center_, width_] := MandelbrotSetPlot[{center - (1 + I) width/2, center + (1 + I) width/2}, ImageResolution -> 500] (* This plots the Mandelbrot set centered at the indicated point with the given frame size *) (* To increase the plot detail, increase ImageResolution *) mandelbrotbox[-0.1+0.7I, 0.9] (* Plots the Mandelbrot set on the box centered at -0.1 + 0.7I with width 0.9 *) (* Here is code for plotting the multibrot set, for problem 9 *) multibrotpoints[n_] := With[{multibrot = Block[{z = #, c = #}, Catch@Do[If[Abs[z] > 2, Throw@i]; (* To change the escape criterion, change 2 to the desired value *) z = z^n + c, (* To plot a different family, simply change the function here *) {i, 100}]] &}, (* To change the number of iterations used, change 100 to the desired value *) Compile[{}, Table[multibrot[x - y I], (* The minus sign is needed because otherwise the plot will be upside-down *) {y, -1.5, 1.5, 0.005}, (* To change the y-range or step size, adjust these values *) {x, -1.5, 1.5, 0.005}]]]; (* To change the x-range or step size, adjust these values *) (* This computes an array of integers, indicating the number of iterates required to escape under the given map *) multibrotset[n_] := ArrayPlot[multibrotpoints[n][], FrameTicks -> {{{0, 1.5}, {100, 1}, {200, 0.5}, {300, 0}, {400, -0.5}, {500, -1}, {600, -1.5}}, {{0, -1.5}, {100, -1}, {200, -0.5}, {300, 0}, {400, 0.5}, {500, 1}, {600, 1.5}}}, FrameLabel -> {Re, Im}, RotateLabel -> False, PlotLabel ->StringJoin[ "Multibrot set for p(z) = z^",ToString[n],"+c"]] (* This plots the array and shows the multibrot set M_3, with manual labeling of the frame to indicate the coordinate values *) multibrotset[3] (* This plots the cubic multibrot set *)