MATLAB: Slove nonlinear differential equation with varying constants

free fallode45

I wonder how I can solve a nonlinear differential equation with varying constants, eg if I want to solve a free fall from high up in the atmosphere where the gravity acceleration i significantly different from that down on the surface of the Earth.
The equation that I want to solve:
ma=-G*M*m/r^2 + k*v^2
"Newtons 2:nd law" = – "Newtons 3:rd law" + "Air drag depending on the square of the velocity"
I'm trying to solve with ode45 and this is the code that I've produced so far:
function ydot=func(t,y)
G=6.6726e-11; %Gravitational constant
M=5.977e24; %Mass of Earth
r=6.371e6; %Radius of Earth
g=(G*M)./((r+y0).^2); %Here is where I get lost. :/ I want the g to be dependent of the altitude
ydot=[y(2); -g + 0.2*y(2).^2];
And some initial values:
t_span=[0 10000]; %Eg from 0 to 10000 s
initial=[r+5e5 0]; %y0=500 km above surface, v0=0 m/s
[t,y]=ode45(@func,t_span,initial);
plot(t,y(:,1));
xlabel('time (s)')
ylabel('height (m)')
For the moment I'm not including the dimensions of the object, neither proper constants of the air drag (which of course also are changing with altitude).
Thanks in advance!

Best Answer

Something like
G=6.6726e-11; %Gravitational constant
M=5.977e24; %Mass of Earth
r=6.371e6; %Radius of Earth
tspan=[0 10000]; %Eg from 0 to 10000 s
initial=[5e5;0]; %y0=500 km above surface, v0=0 m/s
[t,y]=ode45(@(t,y)func(t,y,G,M,r),tspan,initial);
plot(t,y(:,1));
xlabel('time (s)')
ylabel('height (m)')
function ydot=func(t,y,G,M,r)
g=(G*M)/(r+y(1))^2;
ydot=[y(2); -g + 0.2*y(2)^2];
?
Best wishes
Torsten.