Please help! Modify the animate.m script from Lecture 4 to create a movie of the
ID: 3582018 • Letter: P
Question
Please help!
Modify the animate.m script from Lecture 4 to create a movie of the Lorenz attractor similar to the movie embedded on slide 11 of the Lecture 26 notes. To do this, look up the documentation on the MATLAB VideoWriter built-in utility. You should create a movie in either the y1-y2, y2-y3, or y3-y1 planes (your choice, see the embedded movie in the Lecture 26 Notes for inspiration, but note that the example movie labels the axes x,y,z instead of y1, y2, y3 - you do not need to label or show the axes in your movie). Solve the system of equations using a suitable integrator (Lecture 25). For every movie frame, render the past history of the system as a point trace, and the present state marked by a larger dot. Save the movie as an .avi file and upload it on BB. However, to facilitate grading, also include some snapshots of the movie as still images in the Word or PDF document.
Hints: You have several choices and complete creative freedom how to implement this. One way would be to pre-compute the solutions as is done e.g. on slide 12 of Lecture 26 and then iterate through the solutions in a loop, capturing video frames from the saved solutions. An alternative way would be to compute the results on the fly. This would require to write an on-the fly updating function. You could modify rk4sys.m of Lecture 25 to perform only first order Euler integration (since we don’t care so much about accuracy in the movie). Then perform successive integration steps and capture the associated movie frames in real time, terminating the run at will. This is more elegant because it allows the user to end the run interactively, but it requires more sophisticated coding (see also http://stackoverflow.com/questions/23788722/stop-a-infinite-while-loop-pressing-a-key-in-matlab).
clc,clf,clear
g=9.81; theta0=45*pi/180; v0=5;
t(1)=0;x=0;y=0;
plot(x,y,'o','MarkerFaceColor','b','MarkerSize',8)
axis([0 3 0 0.8])
M(1)=getframe;
dt=1/128;
for j = 2:1000
t(j)=t(j-1)+dt;
x=v0*cos(theta0)*t(j);
y=v0*sin(theta0)*t(j)-0.5*g*t(j)^2;
plot(x,y,'o','MarkerFaceColor','b','MarkerSize',8)
axis([0 3 0 0.8])
M(j)=getframe;
if y<=0, break, end
end
pause
movie(M,1)
Explanation / Answer
clc,clf,clear
maxit=1000;
g=9.81; theta0=50*pi/180; v0=5; CR=0.83;
j=1;t(j)=0;x=0;y=0;
xx=x;yy=y;
plot(x,y,'o','MarkerFaceColor','b','MarkerSize',8)
xmax=8; axis([0 xmax 0 0.8])
M(1)=getframe;
dt=1/128; j=1; xxx=0; iter=0;
while(1)
tt=0;
timpact=2*v0*sin(theta0)/g;
ximpact=v0*cos(theta0)*timpact;
while(1)
j=j+1;
h=dt;
if tt+h>timpact,h=timpact-tt;end
t(j)=t(j-1)+h;
tt=tt+h;
x=xxx+v0*cos(theta0)*tt;
y=v0*sin(theta0)*tt-0.5*g*tt^2;
xx=[xx x];
yy=[yy y];
plot(xx,yy,':',x,y,'o','MarkerFaceColor','b','MarkerSize',8)
axis([0 xmax 0 0.8])
M(j)=getframe; iter=iter+1;
if tt>=timpact, break, end
end
v0=CR*v0;
xxx=x;
if x>=xmax|iter>=maxit,break,end
end
pause
clf
axis([0 xmax 0 0.8])
movie(M,1,36)