function hw2()

tmin = 0;
tmax = 1;
h = 0.1;
v0 = 8;

[v1,t1] = ForEuler(v0,tmin,tmax,h);
[v2,t2] = Taylor2(v0,tmin,tmax,h);
[V,T] = ode45(@FuncVel,[tmin tmax],v0);
plot(t1,v1,'g',t2,v2,'y',V,T,'r');

function [v,t] = ForEuler(v0,tmin,tmax,h)

t = tmin:h:tmax;
n = length(t);
v = zeros(n,1);
v(1) = v0;

for i = 2:n;
   v(i) = v(i-1)+h*FuncVel(t(i-1),v(i-1));
end



function [v,t] = Taylor2(v0,tmin,tmax,h)

t = tmin:h:tmax;
n = length(t);
v = zeros(n,1);
v(1) = v0;

for i = 2:n;
    if (v(i-1) >= 0)
        v(i) = v(i-1)+h*(FuncVel(t(i-1),v(i-1))+h/2*PosVel(t(i-1),v(i-1)));
    else
        v(i) = v(i-1)+h*(FuncVel(t(i-1),v(i-1))+h/2*NegVel(t(i-1),v(i-1)));
    end
end


function v = FuncVel(t,v)
v = -9.8 -0.002*v*abs(v)/0.11; 

function vder = PosVel(t,v)
vder = -2*0.002*v*FuncVel(t,v)/0.11;

function vder = NegVel(t,v)
vder = 2*0.002*v*FuncVel(t,v)/0.11;

