MathGroup Archive 2009

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Solving Ordinary differential equations by NDSolve


It makes no sense to enter Pi to two decimal places. In general, enter all constants exactly and let the subsequent processes define the overall precision.

d = 2*10^-2;
A = d^2 Pi/4;
Po = 5*101325;
Pa = 1*101325;
rho = 1000;
V0 = 5*10^-3;
gamma = 114/100;

sol = V /. NDSolve[{
      V'[t] == A Sqrt[2 (Po (V0/V[t])^gamma - Pa)/rho],
      V[0] == V0}, {V}, {t, 0, 9},
     MaxSteps -> 1000000,
     AccuracyGoal -> 10,
     PrecisionGoal -> 10][[1]];

Use Chop to eliminate the numeric noise (imaginary values smaller than your accuracy and precision).

v[t_] := Chop[sol[t]]

Plot[v[t]*1000, {t, 0, 9},
 PlotRange -> All]

Grid[Table[{t, v[t]}, {t, 0, 9, 1}]]


Bob Hanlon

---- Allamarein <matteo.diplomacy at gmail.com> wrote: 

=============
I'd solve this ODE:

V'[t] == A Sqrt[2 (Po (V0/V[t])^gamma - Pa)/rho
IC: V[0] == V0

I wrote this code:

d = 2*10^-2 ;
A = d^2  3.14/4;
Po = 5 *101325;
Pa = 1*101325 ;
rho = 1000 ;
V0 = 5*10^-3 ;
gamma = 1.14;
sol = NDSolve[{
    V'[t] == A Sqrt[2 (Po (V0/V[t])^gamma - Pa)/ rho],
    V[0] == V0},
   {V}, {t, 0, 9},
   MaxSteps -> 1000000, AccuracyGoal -> 10, PrecisionGoal -> 10];
v[t_] := V[t] /. sol[[1]];
Plot[Evaluate[V[t] /. sol]*1000, {t, 0, 9}, PlotRange -> All]
Grid[Table[{t, v[t]}, {t, 0, 9, 1}]]

If it can be useful, I can underline units of these variables:
d [m]
P0, Pa [Pa]
rho [kg/m^3]
V [m^3]
gamma [--]
t [sec]

Running this code, V has got comlex part. This is impossible, because
it's a volume.
I should re-write my ODE in order NDSolve can digest better or I can
set an option in my code where I suggest " V must be positive and
real"?



  • Prev by Date: Re: Mathematica skill level snippet(s)
  • Next by Date: Re: Solving Ordinary differential equations by NDSolve
  • Previous by thread: Re: Solving Ordinary differential equations by NDSolve
  • Next by thread: Re: Solving Ordinary differential equations by NDSolve