% codice oscill armonico (smorzato) come in slide 87 del moto vario. 
% 

s = 2000; %array size (spacial resolution of the simulation)

% create arrays
z = zeros(s,1);  
v = zeros(s,1);  

% parameters
lam  = 2.0;      % coeff J - try es. 1.0
dt   = 0.001;
L    = 0.5;
g    = 9.8;
endTime = s*dt;
t = dt:dt:endTime;

% initial conditions
v(1) = 0.0;
z(1) = 0.5;

% main loop
for k = 2 : 1 : s
  J    = lam*L*v(k-1)*abs(v(k-1))/(2.0*g);
  v(k) = v(k-1) - (g*dt/L)*(2.0*z(k-1) + J);
  z(k) = z(k-1) + v(k-1)*dt;
end
%
% plot 
for k = 2 : 10 : s
    plot(t(1:k),z(1:k),"color","black")
    %hold on
		axis([0 endTime -0.8 0.8])  	
    	xlabel("tempo [secondi]")
        ylabel("altezza [metri]")
		%caxis([-1 1]) %color axis scale
        grid on
		pause ( .00005 )
end
hold on