function explicit_euler
T = 10; y_0 = 0; v_0 = 1;
s = 5; tau = 2^(-s); K = floor(T/tau);
y = zeros(K+1,2); y(1,:) = [y_0,v_0];
for k = 1:K
    y(k+1,:) = y(k,:)+tau*Phi((k-1)*tau,y(k,:),tau);
    plot(tau*(0:k),y(1:k+1,1),'r');
    axis([0,T,-5,5]); drawnow;
end
D = 1; r = 1/10; omega = sqrt(D-r^2/4); t = K*tau; 
function val = Phi(t,y,tau)
val = f(t,y);
function vec = f(t,y)
r = 1/10; D = 1; ell = 0;
vec = [y(2),-r*y(2)-D*(y(1)-ell)];
