Re: Solving systems of differential equations
- To: mathgroup at smc.vnet.net
- Subject: [mg105248] Re: Solving systems of differential equations
- From: dh <dh at metrohm.com>
- Date: Wed, 25 Nov 2009 02:35:08 -0500 (EST)
- References: <hegdl4$kk$1@smc.vnet.net>
Yun Zhao wrote:
> I have a question as to an apparent inconsistency I noticed while using
> Mathematica to solve differential equations, I have 2 ODE like this, and I
> solve them symbolically using Mathematica 7.0:
>
> system = {p'[t] == -(rho + d)*p[t],
> n'[t] == rho*p[t] + (rho - d)*n[t]};
> sol = DSolve[system, {p, n}, t][[1]]
>
> I get the following solution:
>
> {p -> Function[{t}, E^((-d - rho) t) C[1]],
> n -> Function[{t},
> 1/2 (-E^((-d - rho) t) + E^((-d + rho) t)) C[1] +
> E^((-d + rho) t) C[2]]}
>
> It looks to be correct.
>
> Then, manually, I know p(0)=C1, n(0)=0, so C1=p0, C2=0.
>
> I know the values for rho and d. So for t from 0 to 120, I can calculate and
> plot the values for p(t) and n(t).
>
> Up to now, this is all good.
>
> Now I want to solve the system of ODEs numerically (even though I have
> already solved it, I have good reasons to want to learn to solve it
> numerically, as I am going to add more complexity to this model), so I did
> the following:
>
> sol4 = NDSolve[{p'[t] == -(0.0258 + 0.0123)*p[t],
> n'[t] == 0.0258*p[t] + (0.0258 - 0.0123)*n[t], p[0] == 30000,
> n[0] == 0}, {p, n}, {t, 0, 120}]
>
> Now I get the following solution:
>
> {{p-> InterpolatingFunction[{{0.`,120.`}},"<>"],
> n->InterpolatingFunction[{{0.`,120.`}},"<>"]}}
>
> It may be correct, but I then plotted it using
>
> Plot[Evaluate[{n[t]} /. First[sol4]], {t, 0, 120},
> PlotRange -> All]
>
> Now my question is:
>
> I generated plot #1 of n(t) by solving the system of ODEs first
> symbolically, then manually inputting in values (C1, C2, rho, d).
> I generated plot #2 of n(t) by solving the system of ODEs numerically, then
> had mathematica plot the solutions of and n(t) for me.
>
> Why does plot #1 and #2 look so different? By that I mean, n(t) in plot 2 is
> much lower than n(t) in plot 1.
>
> Could you please explain what I did wrong? Thank you.
>
> I know mathgroup doesn't like attachment, but is there anyway I can show
> these plots?
>
>
Hi,
I can not see a difference:
rho0 = 0.0258; d0 = 0.0123; p0 = 30000; n0 = 0;
system = {p'[t] == -(rho + d)*p[t],
n'[t] == rho*p[t] + (rho - d)*n[t]};
sol = DSolve[system, {p, n}, t][[1]] /. {rho -> rho0, d -> d0,
C[1] -> p0, C[2] -> n0}
Plot[Evaluate[{p[t], n[t]} /. sol], {t, 0, 120}]
sol4 = NDSolve[{p'[t] == -(rho + d)*p[t],
n'[t] == rho*p[t] + (rho - d)*n[t], p[0] == p0,
n[0] == n0} /. {rho -> rho0, d -> d0}, {p, n}, {t, 0, 120}]
Plot[Evaluate[{p[t], n[t]} /. sol4], {t, 0, 120}]
Daniel