Re: Solving Ordinary differential equations by NDSolve
- To: mathgroup at smc.vnet.net
- Subject: [mg104976] Re: [mg104962] Solving Ordinary differential equations by NDSolve
- From: "David Park" <djmpark at comcast.net>
- Date: Sun, 15 Nov 2009 20:47:38 -0500 (EST)
- References: <27585103.1258283977493.JavaMail.root@n11>
One possible solution is to construct a second FunctionInterpolation, taking
only the Re part, since the imaginary part is just round off errors. Also
use Set instead of SetDelayed when calculating v[t].
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_] = FunctionInterpolation[Re[V[t]] /. sol[[1]], {t, 0, 9}][t]
Plot[v[t]*1000, {t, 0, 9}, PlotRange -> All] Grid[
Table[{t, v[t]}, {t, 0, 9, 1}]]
David Park
djmpark at comcast.net
http://home.comcast.net/~djmpark/
From: Allamarein [mailto:matteo.diplomacy at gmail.com]
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"?