%
% 4 POINTS FOR DECLARING ALL NECESSARY VARIABLES
% first, we declare the constants for the problem
%
%better display values for plotting
fontsize=20;
linewidth=2;
markersize=10;

%default values for plotting
% fontsize=12;
% linewidth=0.5;
% markersize=6;

refine = true;
N_fine = 100;

g = 9.8;
vxo = 0.01;
vyo = 0.01;
to = 0;     %starting time for my ODE solution
tf = 0.0021; %this is ~2*vyo/g, ending time for my ODE solution

% 10 POINTS FOR CORRECT FUNCTION CALL AND CORRECT FUNCTION INTERNALS
%options = odeset('RelTol',1e-3,'AbsTol',[1e-6 1e-6]);   %Matlab default values
options = odeset('RelTol',1e-3,'AbsTol',[1e-6 1e-6]);
[t,v] = ode45(@accelerations,[to tf],[vxo vyo],options);
%the syntax for any ODE integrator is:
%[t, f(t)] = odesolver(@(function that returns [t, df/dt]), [(start time) (end time)], [f(t=0) df/dt(t=0)])

vx = v(:,1);
vy = v(:,2);

% 5 BONUS POINTS IF THEY HAVE ANYTHING LIKE WHAT FOLLOWS, I.E. A TREATMENT
% THAT MAKES THE SOLUTION MORE ACCURATE.
if refine
    t_fine = linspace(t(1),t(end),N_fine*length(t))';
    vx = spline(t,vx,t_fine);
    vy = spline(t,vy,t_fine);
    t = t_fine;
end

% 8 POINTS FOR CORRECT INTEGRATION
x = cumsum(diff(t).*vx(1:end-1));
y = cumsum(diff(t).*vy(1:end-1));

stop = length(y);
for ii=1:length(y)
    if y(ii) < 0
        stop = ii+2;
        break
    end
end

x(stop:end) = [];
y(stop:end) = [];

% 4 POINTS FOR VACUUM SOLUTION
xvac = vxo*t;
yvac = vyo*t - 0.5*g*t.^2;

stop = length(yvac);
for ii=1:length(yvac)
    if yvac(ii) < 0
        stop = ii+2;
        break
    end
end

xvac(stop:end) = [];
yvac(stop:end) = [];

% 4 POINTS FOR ACTUALLY PLOTTING WHAT I ASKED FOR
figure
plot(xvac,yvac,'LineWidth',linewidth)
hold all
plot(x,y,'LineWidth',linewidth)
xlabel('x [m]','FontSize',fontsize)
ylabel('y [m]','FontSize',fontsize)
set(gca,'FontSize',fontsize,'LineWidth',linewidth)