Re: Solving Ordinary differential equations by NDSolve

*To*: mathgroup at smc.vnet.net*Subject*: [mg104982] Re: [mg104962] Solving Ordinary differential equations by NDSolve*From*: Bob Hanlon <hanlonr at cox.net>*Date*: Sun, 15 Nov 2009 20:48:47 -0500 (EST)*Reply-to*: hanlonr at cox.net

The code that I sent you runs on my system. Your latest code appears to have an error in the definition of V0, i.e., V0 = V0 = 5*10 - 3 rather than V0 = 5*10^- 3 Bob Hanlon ---- Matteo <matteo.diplomacy at gmail.com> wrote: ============= So..you suggest to modify in this way my code: 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]]; v[t_] := Chop[sol[t]] Plot[v[t]*1000, {t, 0, 9}, PlotRange -> All] Grid[Table[{t, v[t]}, {t, 0, 9, 1}]] Does it run on you machine? I get this error message: DSolve::mxst: Maximum number of 1000000 steps reached at the point t == 0.0789357392769894`. I tried to set up MaxStep -> 10^7 but the new error is DSolve::mxst: Maximum number of 1000000 steps reached at the point t == 0.11680804227781108`. I had the problem to have imaginary part for variables that I know it must be real. I would solve my trouble definitively by this example-problem. Bob Hanlon ha scritto: > 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"? > > > -- Bob Hanlon