Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
SetDirectory[NotebookDirectory[]];
\[Sigma] = 10;
\[Rho] = 28;
\[Beta] = 8/3;
Tmax = 10;
eqn = {
x'[t] == \[Sigma] (-x[t] + y[t]),
y'[t] == \[Rho]* x[t] - y[t] - x[t]*z[t],
z'[t] == x[t]*y[t] - \[Beta]*z[t]
};
trange = 0.1;
soln = ParametricNDSolveValue[
Join[eqn, {x[0] == x0, y[0] == y0, z[0] == z0}], {x[t], y[t],
z[t]}, {t, -1.1*trange, Tmax}, {x0, y0, z0}];
nParticles = 2500;
boxRange = 30;
icRange = 0.8;
center = {0, 0, 25};
initalConditions =
center + RandomReal[{-boxRange*icRange, boxRange*icRange}, {3,
nParticles}];
allPoints = Evaluate[MapThread[soln, initalConditions]];
drawplot[tt_, imagesize_ : 1080] := Module[{t0, t1},
t0 = Max[0, 0 + tt];
t1 = Min[Tmax, trange + tt];
Show[
Graphics3D[{
PointSize -> .004,
Opacity[0.9],
White,
Glow[White],
Point[allPoints /. t -> tt + trange]
}],
ParametricPlot3D[
allPoints, {t, t0, t1},
(*PlotPoints\[Rule]200,MaxRecursion\[Rule]6,*)
ColorFunction ->
Function[{x, y, z, t},
Directive[
ColorData["SolarColors"][(*.05+.85**)((t - tt)/trange)^4],
Opacity[0.5*(t - tt)/trange]^1]],
ColorFunctionScaling -> False],
PlotRange ->
center + {{-boxRange, boxRange}, {-boxRange,
boxRange}, {-boxRange, boxRange}},
BoxRatios -> {1, 1, 1},
ImagePadding -> 50,
PlotRangeClipping -> False,
Axes -> False,
Background -> Black,
Boxed -> True, BoxStyle -> White,
ImageSize -> imagesize,
ViewPoint -> 5*{1.3, -2.4, 1.},
ViewVertical -> {0, 0, 1}
]
];
pointsAppear[numPoints_, imagesize_ : 1080] := Module[{},
Graphics3D[{
PointSize -> .004,
White,
Opacity[0.9],
Point[allPoints[[;; numPoints]] /. t -> 1*^-6],
Opacity[0.0],
Point[allPoints[[numPoints + 1 ;;]] /. t -> 1*^-6]
},
PlotRange ->
center + {{-boxRange, boxRange}, {-boxRange,
boxRange}, {-boxRange, boxRange}},
BoxRatios -> {1, 1, 1},
ImagePadding -> 50,
Axes -> False,
Background -> Black,
Boxed -> True, BoxStyle -> White,
ImageSize -> imagesize,
ViewPoint -> 5*{1.3, -2.4, 1.},
ViewVertical -> {0, 0, 1}
]
];
renderFrame[ttt_, imagesize_ : 1080] := Module[{numPoints},
numPoints =
Ceiling[nParticles*(1 + ((ttt +
trange)/.05))]; (*0 to max in \[Delta]t=0.05*)
numPoints = Clip[numPoints, {1, nParticles}];
If[ttt + trange < 0, pointsAppear[numPoints, imagesize],
drawplot[ttt, imagesize]]
];
saveframe[ttt_] := Module[{frame, title},
frame = renderFrame[ttt - trange + 1*^-6];
title = IntegerString[Floor[(ttt + 1)*10000], 10, 9] <> ".png";
Export["frames/" <> title, frame];
Print[title];
];
duration = 15;
nFrames = 60 * duration;
dt = (0.7 - (-0.05))/nFrames;
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment