Mathematica绘图代码实例6个
标签:
计算机作图mandelbrot集合 |
Mathematica 代码实例6个(可复制粘贴到Mathematica命令代码行)
1)///////////////////////////////
(*runtime:1 minute*)
Mandelbrot[c_] :=
Module[{z = 0, i = 0},
While[i < 100 &&
Abs[z] < 2, z = z^2 + c; i++]; i];
DensityPlot[Mandelbrot[xc + I yc], {xc, -2, 1}, {yc, -1.5, 1.5},
PlotPoints -> 275, Mesh -> False, Frame
-> False,
ColorFunction -> (If[# != 1, Hue[#], Hue[0, 0,
0]] &)]
Mandelbrot[c_] :=
DensityPlot[Mandelbrot[xc + I yc], {xc, -2, 1}, {yc, -1.5, 1.5},
2)//////////////////////////////////////////
(* runtime: 34 seconds *)
order = 12; n = 275; image = Table[0.0, {n}, {n}];
Do[Do[z = N[Root[Sum[(2Mod[Floor[(t - 1)/2^i], 2] - 1) #^(order - i), {i, 0, order}], root]]; {j,i} = Round[n({Re[z], Im[z]}/1.5 + 1)/2]; If[0 < i <= n && 0 < j <= n, image[[i, j]]++], {root, 1, order}], {t, 1, 2^order}];
ListDensityPlot[image, Mesh -> False, Frame -> False, PlotRange -> {0, 4}]
order = 12; n = 275; image = Table[0.0, {n}, {n}];
Do[Do[z = N[Root[Sum[(2Mod[Floor[(t - 1)/2^i], 2] - 1) #^(order - i), {i, 0, order}], root]]; {j,i} = Round[n({Re[z], Im[z]}/1.5 + 1)/2]; If[0 < i <= n && 0 < j <= n, image[[i, j]]++], {root, 1, order}], {t, 1, 2^order}];
ListDensityPlot[image, Mesh -> False, Frame -> False, PlotRange -> {0, 4}]
3)///////////////////////////////////////////////////////
(* runtime: 0.4 second *)
Tree[p1_, theta_, r_, i_] := Module[{p2 = p1 + r{Cos[theta], Sin[theta]}, graphics},graphics = {Thickness[0.07r], RGBColor[0.5r, 1 - r, 0], Line[{p1, p2}]}; If[i > 1, graphics = Join[graphics, {Tree[p2, theta - Pi/6, 0.8r, i - 1], Tree[p2, theta + Pi/4, 0.7r, i - 1]}]]; graphics];
Show[Graphics[Tree[{0, 0}, Pi/2,1, 12], AspectRatio -> 1, Background -> RGBColor[0, 0, 0], PlotRange -> {{-3, 3}, {0, 6}}]]
Tree[p1_, theta_, r_, i_] := Module[{p2 = p1 + r{Cos[theta], Sin[theta]}, graphics},graphics = {Thickness[0.07r], RGBColor[0.5r, 1 - r, 0], Line[{p1, p2}]}; If[i > 1, graphics = Join[graphics, {Tree[p2, theta - Pi/6, 0.8r, i - 1], Tree[p2, theta + Pi/4, 0.7r, i - 1]}]]; graphics];
Show[Graphics[Tree[{0, 0}, Pi/2,1, 12], AspectRatio -> 1, Background -> RGBColor[0, 0, 0], PlotRange -> {{-3, 3}, {0, 6}}]]
4)//////////////////////////////////////////////////
(* runtime: 0.7 second *)
Norm[x_] := x.x; Ry[phi_] := {{Cos[phi], 0, Sin[phi]}, {0,1, 0}, {-Sin[phi], 0, Cos[phi]}}; Rz[theta_] := {{Cos[theta], -Sin[theta], 0}, {Sin[theta], Cos[theta], 0}, {0, 0, 1}};
Tree[p1_, R_, i_] := Module[{p2 = p1 + R.{0, 0, 1}, r, graphics}, r = Sqrt[Norm[p2 - p1]]; graphics = {Thickness[0.07r], RGBColor[0.5r, 1 - r, 0], Line[{p1, p2}]}; If[i > 1, graphics = Join[graphics, Table[Tree[p2, 0.6R.Rz[theta].Ry[-Pi/4], i - 1], {theta, 0, 4 Pi/3, 2Pi/3}]]]; graphics];
Show[Graphics3D[Tree[{0, 0, 0}, IdentityMatrix[3], 8]]]
Norm[x_] := x.x; Ry[phi_] := {{Cos[phi], 0, Sin[phi]}, {0,1, 0}, {-Sin[phi], 0, Cos[phi]}}; Rz[theta_] := {{Cos[theta], -Sin[theta], 0}, {Sin[theta], Cos[theta], 0}, {0, 0, 1}};
Tree[p1_, R_, i_] := Module[{p2 = p1 + R.{0, 0, 1}, r, graphics}, r = Sqrt[Norm[p2 - p1]]; graphics = {Thickness[0.07r], RGBColor[0.5r, 1 - r, 0], Line[{p1, p2}]}; If[i > 1, graphics = Join[graphics, Table[Tree[p2, 0.6R.Rz[theta].Ry[-Pi/4], i - 1], {theta, 0, 4 Pi/3, 2Pi/3}]]]; graphics];
Show[Graphics3D[Tree[{0, 0, 0}, IdentityMatrix[3], 8]]]
5)////////////////////////////////
(* runtime: 65 seconds, increase n for higher resolution
*)
n = 140; r = 1; plight = {0, 0, 1}; Normalize[x_] := x/Sqrt[x.x];
spheres = {{{1, 1.0/Sqrt[3], 0}, {1, 0,0}}, {{-1, 1.0/Sqrt[3], 0}, {0, 1, 0}}, {{0, -2.0/Sqrt[3], 0}, {0, 0, 1}}, {{0, 0, -Sqrt[8.0/3]}, {1, 1, 1}}};
Intersections[p_, dir_] := Module[{a = p.dir, b}, b = r^2 + a^2 - p.p; If[b > 0, Select[Sqrt[b]{-1, 1} - a, # > 10^-10 &], {}]];
RayTrace[p_, dir_, depth_] := Module[{dlist = Map[Min[Intersections[p - #[[1]], dir]] &, spheres], d, i, p2, normal, color, lightdir, intensity}, d = Min[dlist]; If[d < Infinity, p2 = p + d dir; i = Position[dlist, d][[1, 1]]; normal = Normalize[p2 - spheres[[i, 1]]]; color = If[depth > 50, {0, 0, 0}, RayTrace[p2, Normalize[dir - 2(normal.dir)normal], depth + 1]]; lightdir = Normalize[plight - p2]; intensity = normal.lightdir; If[intensity > 0, shadowed = Or @@ Map[Intersections[p2 - #[[1]], lightdir] =!= {} &, spheres]; If[! shadowed,color += intensity spheres[[i, 2]]]]; color, {0, 0, 0}]];
Show[Graphics[RasterArray[Table[RGBColor @@ Map[Min[1, Max[0, #]] &, RayTrace[{0, 0, 1},Normalize[{x,y, -1}], 1]], {y, -0.5, 0.5, 1.0/n}, {x, -0.5, 0.5, 1.0/n}]], AspectRatio -> 1]]
n = 140; r = 1; plight = {0, 0, 1}; Normalize[x_] := x/Sqrt[x.x];
spheres = {{{1, 1.0/Sqrt[3], 0}, {1, 0,0}}, {{-1, 1.0/Sqrt[3], 0}, {0, 1, 0}}, {{0, -2.0/Sqrt[3], 0}, {0, 0, 1}}, {{0, 0, -Sqrt[8.0/3]}, {1, 1, 1}}};
Intersections[p_, dir_] := Module[{a = p.dir, b}, b = r^2 + a^2 - p.p; If[b > 0, Select[Sqrt[b]{-1, 1} - a, # > 10^-10 &], {}]];
RayTrace[p_, dir_, depth_] := Module[{dlist = Map[Min[Intersections[p - #[[1]], dir]] &, spheres], d, i, p2, normal, color, lightdir, intensity}, d = Min[dlist]; If[d < Infinity, p2 = p + d dir; i = Position[dlist, d][[1, 1]]; normal = Normalize[p2 - spheres[[i, 1]]]; color = If[depth > 50, {0, 0, 0}, RayTrace[p2, Normalize[dir - 2(normal.dir)normal], depth + 1]]; lightdir = Normalize[plight - p2]; intensity = normal.lightdir; If[intensity > 0, shadowed = Or @@ Map[Intersections[p2 - #[[1]], lightdir] =!= {} &, spheres]; If[! shadowed,color += intensity spheres[[i, 2]]]]; color, {0, 0, 0}]];
Show[Graphics[RasterArray[Table[RGBColor @@ Map[Min[1, Max[0, #]] &, RayTrace[{0, 0, 1},Normalize[{x,y, -1}], 1]], {y, -0.5, 0.5, 1.0/n}, {x, -0.5, 0.5, 1.0/n}]], AspectRatio -> 1]]
6)///////////////////////////////////////////////
(* runtime: 25 seconds *)
f[z_] := z^36 - 1;
Show[Graphics[RasterArray[Table[data = FixedPointList[# - f[#]/f'[#] &, x + I y, 29]; z = 1 - Length[data]/30; Hue[Arg[Last[data]]/(2Pi), Min[1,2(1 - z)], Min[1, 2z]], {y, -1.5, 1.5, 3.0/275}, {x, -1.5, 1.5, 3.0/275}]]], AspectRatio -> 1]
(* runtime: 25 seconds *)
f[z_] := z^36 - 1;
Show[Graphics[RasterArray[Table[data = FixedPointList[# - f[#]/f'[#] &, x + I y, 29]; z = 1 - Length[data]/30; Hue[Arg[Last[data]]/(2Pi), Min[1,2(1 - z)], Min[1, 2z]], {y, -1.5, 1.5, 3.0/275}, {x, -1.5, 1.5, 3.0/275}]]], AspectRatio -> 1]
/////////////////////////////////////////////////附录
机翼理论示意图
(* runtime: 13 seconds *)
Clear[z]; U = rho = 1; chord = 4; thk = 0.5; alpha = Pi/9; y0 = 0.2; x0 = -thk/5.2; L = chord/4; a = Sqrt[y0^2 + L^2]; gamma = 4Pi a U Sin[alpha + ArcCos[L/a]];
w[z_, sign_] := Module[{zeta = (z + sign Sqrt[z^2 - 4 L^2])/2}, zeta = (zeta - x0 - I y0)Exp[-I alpha]/Sqrt[(1 - x0)^2 + y0^2]; U(zeta + a^2/zeta) + I gamma Log[zeta]/(2Pi)];
sign[z_] := Sign[Re[z]]If[Abs[Re[z]] < chord/2 && 0 < Im[z] < 2y0(1 - (2Re[z]/chord)^2), -1, 1];w[z_] := w[z, sign[z]]; V[z_] = D[w[z, sign], z] /. sign -> sign[z];
ContourPlot[Im[w[(x + I y)Exp[I alpha]]], {x, -3, 3}, {y, -3, 3}, Contours -> Table[x^3 + 0.0208, {x, -2, 2, 0.1}], PlotPoints -> 100, ContourShading -> False, AspectRatio -> Automatic]
DensityPlot[-0.5rho Abs[V[(x + I y)Exp[I alpha]]]^2, {x, -3, 3}, {y, -3, 3}, PlotPoints -> 275, Mesh -> False, Frame -> False, ColorFunction -> (If[# == 1, Hue[0, 0, 0], Hue[(5# - 1)/6]] &), AspectRatio -> Automatic]
///////////////////////////////////////////////////////////
Clear[z]; U = rho = 1; chord = 4; thk = 0.5; alpha = Pi/9; y0 = 0.2; x0 = -thk/5.2; L = chord/4; a = Sqrt[y0^2 + L^2]; gamma = 4Pi a U Sin[alpha + ArcCos[L/a]];
w[z_, sign_] := Module[{zeta = (z + sign Sqrt[z^2 - 4 L^2])/2}, zeta = (zeta - x0 - I y0)Exp[-I alpha]/Sqrt[(1 - x0)^2 + y0^2]; U(zeta + a^2/zeta) + I gamma Log[zeta]/(2Pi)];
sign[z_] := Sign[Re[z]]If[Abs[Re[z]] < chord/2 && 0 < Im[z] < 2y0(1 - (2Re[z]/chord)^2), -1, 1];w[z_] := w[z, sign[z]]; V[z_] = D[w[z, sign], z] /. sign -> sign[z];
ContourPlot[Im[w[(x + I y)Exp[I alpha]]], {x, -3, 3}, {y, -3, 3}, Contours -> Table[x^3 + 0.0208, {x, -2, 2, 0.1}], PlotPoints -> 100, ContourShading -> False, AspectRatio -> Automatic]
DensityPlot[-0.5rho Abs[V[(x + I y)Exp[I alpha]]]^2, {x, -3, 3}, {y, -3, 3}, PlotPoints -> 275, Mesh -> False, Frame -> False, ColorFunction -> (If[# == 1, Hue[0, 0, 0], Hue[(5# - 1)/6]] &), AspectRatio -> Automatic]
///////////////////////////////////////////////////////////
下载曼氏集合绘图软件:
前一篇:数学家Mandelbrot生平
后一篇:Mandelbrot集合浅说

加载中…