MATLAB: Ode45 returns NaN, Why

nanode45

Having a system of three odes
I put the odes in a function and the script to solve the given function
I'm getting NaN for output of ode45
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function dx=funF7(t,x)
ft_to_m=0.3048;
lb_to_kg=0.45359273;
S=213.8264*ft_to_m^2;
m=22900*lb_to_kg;
CD=1.34721248450636;
CL=0.387532621392346;
phi1=0;
rho_0=1.226;
beta=0.1387/1000;
R_Earth=6378000;
%State variables
r=x(1);
V=x(2);
gamma=x(3);
h=r-R_Earth;
rho=rho_0*exp(-beta*h);
%Nonlinear equations with recpect to time
q=(V^2*rho)/2;
L=q*S*CL/m;
D=q*S*CD/m;
%Equations of motion
dx=zeros(3,1);
dx(1)=V*sin(gamma);
dx(2)=-D/cos(phi1) - sin(gamma)/r^2;
dx(3)=L/V + V*cos(gamma)/r - cos(gamma)/r^2/V;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Solver
clc; clear all
step_size1=0.0004;
r0=0.3256;
V0=2.4741;
gamma0=-0.1027;
x0=[r0 V0 gamma0];
tspan=[0 0.1];
options = odeset('MaxStep',step_size1);
[tt xx]= ode45(@funF7,tspan,x0,options);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Any help would be much appreciated. Thanks in advance!

Best Answer

rho = 1.226*exp(-1.387e-4*(0.3256-6378000)) ~ 1.226*exp(885) !!!
Best wishes
Torsten.