clf

x=0:0.1:3;
y=0:0.1:3;

y(1)=1;

for i=2:31
  y(i)=y(i-1)-0.1*2*x(i-1)*y(i-1);
end

plot(x,y,'r+');

hold on

plot(x,exp(-x.^2));

legend('approximated solution','exact solution')
xlabel('x');
ylabel('y');
title('y\prime=-2xy');