Limitations of NDSolve

• To: mathgroup at smc.vnet.net
• Subject: [mg2515] Limitations of NDSolve
• From: David Garza <dgarza at utsi.edu>
• Date: Wed, 15 Nov 1995 02:01:42 -0500
• Organization: University of Tennessee, Knoxville

```I have tried using Mathematica to numerically solve the equations of
motion for a vehicle in flight about a round Earth.  The flight is in
two dimensions and the equations consist of a set of five nonlinear
ordinary differential equations.  I have tried the same thing without
problems using computer programs and Excel.  Even the Runge-Kutta
package supplied with Mathematica does fine.  Unfortunately, NDSolve
will not work, instead reporting an infinite expression.  The initial
conditions normally provide a small, nonzero value for the velocity
to eliminate the problem of velocity in the denominator in the fourth
equation.  This has been something I've tried off and on a couple of
times in the past, and in talking with others I have found that NDSolve
has had trouble with other equations when those equations have not
presented problems to other means of computing numerical solutions,
like writing a program in Fortran or C.  What sort of limitations are
there on NDSolve with regard to the problems it can handle?  The system
of equations I use are included below, along with initial conditions and
the parameters necessary for solution.

ODE's:
{
th'[t]==-(v[t] Cos[gam[t]]/r[t]),
r'[t]==v[t] Sin[gam[t]],
v'[t]==(-drag+thrust Cos[al+eps]-g m[t]*
Sin[gam[t]])/m[t],
gam'[t]==(lift-m[t] g Cos[gam[t]]+ m[t]*
v[t]^2 Cos[gam[t]]/r[t]+thrust Sin[al+eps])/
(m[t]*v[t]),
m'[t]==-thrust/(Isp g)
}

Initial Conditions:
{
r[0]==r0,
th[0]==.1,
v[0]==.1,
gam[0]==Pi*89.995/180,
m[0]==10000
}

Other important functions:
thrust=1.3*m0*g;
drag=.5*rho*v[t]^2*area*cd;
lift=0;      (* This is a rocket, so there is little lift - neglect it now *)
g=9.81;      (* surface gravitational acceleration *)
m0=10000;    (* Initial mass *)
r0=6376000;  (* Radius of the Earth *)
area=10;     (* Rocket frontal area *)
cd=.1;       (* Drag coefficient *)
Isp=400;     (* Engine specific impulse *)
eps=0;       (* Engine gimbal angle, try 0 for now to see if this works *)
al=0;        (* Angle of attack, set to zero to see if this works *)
rho=.00276;  (* Air density, try a constant for now *)

Thanks,
--
David Garza
dgarza at utsi.edu

```

• Prev by Date: Re: turn off symbols in MultipleListPlot?
• Next by Date: Re: 2-d graphs to 3-d
• Previous by thread: Limitations of NDSolve
• Next by thread: Re: Limitations of NDSolve