%Matlab Simulation of Filter drop: Fnet=-mg+bV^2 (y+ up, y=0 is floor, y=h is starting height)
dt=0.001;
npoints = 2500;
h=1.5;
m=0.0007;%0.7 grams
g=9.8;
b=0.015;
F= zeros(npoints,1);%Preload and fill with 0's  a npoints=rows column vector to store the values.
p= zeros(npoints,1);
v= zeros(npoints,1);
y= zeros(npoints,1);
t= zeros(npoints,1);
t(1)=0;%Initial values corresponds to first row (=1)
v(1)=0;
y(1)=h;
p(1)=m*v(1);
F(1)=-m*g+b*v(1)^2;

for step=1:npoints-1
    t(step+1) = t(step) + dt; %Time Update
    F(step+1) = -m*g+b*v(step)^2 ;%Force Update
    p(step+1) = p(step) + F(step)*dt;%Momentum Update
    v(step+1) = p(step+1)/m;
    y(step+1) = y(step)+v(step+1)*dt;%Position Update

end

% Plot
plot(t,y,'r',t,v,'b');
xlabel('Time(s)')
ylabel('y(red),v(blue)')