Animating the Lorenz attractor with thousands of initial conditions​
​by Ben Bartlett
It’s called the Lorenz attractor because all nearby points are attracted into the set of chaotic orbits, regardless of initial conditions.
In[]:=
SetDirectory[NotebookDirectory[]];
In[]:=
Needs["AnimationTools`"]
In[]:=
σ=10;​​ρ=28;​​β=8/3;​​Tmax=10;
In[]:=
eqn={​​x'[t]σ(-x[t]+y[t]),​​y'[t]ρ*x[t]-y[t]-x[t]*z[t],​​z'[t]x[t]*y[t]-β*z[t]​​};
In[]:=
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}];
In[]:=
nParticles=2500;​​boxRange=30;​​icRange=0.8;​​center={0,0,25};​​initalConditions=center+RandomReal[{-boxRange*icRange,boxRange*icRange},{3,nParticles}];
In[]:=
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200,MaxRecursion6,*)​​ColorFunctionFunction[{x,y,z,t},Directive[ColorData["SolarColors"][(*.05+.85**)
4
((t-tt)/trange)
],
1
Opacity[0.5*(t-tt)/trange]
]],​​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))];(*0tomaxinδt=0.05*)​​numPoints=Clip[numPoints,{1,nParticles}];​​If[ttt+trange<0,pointsAppear[numPoints,imagesize],drawplot[ttt,imagesize]]​​];​​(*renderFrame[1.0]*)
In[]:=
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];​​];
In[]:=
duration=15;​​nFrames=60*duration;​​dt=(0.7-(-0.05))/nFrames;
In[]:=
Do[saveframe[ttt],{ttt,-0.05,0.7,dt}]