F=@(x) [16e4*x+20e4*(x-0.01).*(x>0.01)]/4;
t=linspace(0,0.02);dt=t(5)-t(4);
z(1)=0;z(2)=4*dt;
for i=2:99;
z(i+1)=2*z(i)-z(i-1)-F(z(i))*dt^2;
end;
for i=1:100;
figure(1);set(h,'ydata',[0.04 0.04]-z(i));
set(h1,'ydata',y+linspace(0,-z(i)));
pause(0.4);
end;