Mathematica 9 is now available
Services & Resources / Wolfram Forums
-----
 /
MathGroup Archive
2000
*January
*February
*March
*April
*May
*June
*July
*August
*September
*October
*November
*December
*Archive Index
*Ask about this page
*Print this page
*Give us feedback
*Sign up for the Wolfram Insider

MathGroup Archive 2000

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

Search the Archive

RE: DSolve problems with system of ODEs. Solutions and Thanks

  • To: mathgroup at smc.vnet.net
  • Subject: [mg21786] RE: DSolve problems with system of ODEs. Solutions and Thanks
  • From: "SANCHEZ DE LEON, Guillermo" <gsl at fab.enusa.es>
  • Date: Thu, 27 Jan 2000 22:57:09 -0500 (EST)
  • Sender: owner-wri-mathgroup at wolfram.com

[My usual E-mail is guillerm at gugu.usal.es]

I posted a few days ago the following question:

> > If I evaluate, in order to obtain the analitycal solution, the System of
> ODEs
> >
> > DSolve[{x'[t] == -0.036 x[t] + 0.0124 y[t] + 0.000035 z[t] + 49.3,
> >     y'[t] == 0.0111 x[t] - 0.0286 y[t], z'[t] == 0.0039 x[t] - 0.000035
> > z[t],
> >     x[0] == 0, y[0] == 0, z[0] == 0}, {x[t], y[t], z[t]}, t]
> >
> >In version 4.0.0  you receive a massage like this
> > (RootSum::"pfn" :
> >     "(DSolve`DSolveDump`dysfunction$6[273147 + <<1>>...
> > is not a pure function.) 
In version 2.2 you will receive the correct solutions, and in version 3 you
will receive a error message

I have received several ideas for menber of Mathgroup (Thanks every body, it
a pleasure to have the opportunity of share this group with excelents
Mathematica gurus and belong). 

Solution (in Math v4.0.0) forcing the numerical solutions accoding with
Jens-Peer Kuska suggestion

> eqns = 
> 
> {x'[t] == -0.036 x[t] + 0.0124 y[t] + 0.000035 z[t] + 49.3,
> y'[t] == 0.0111 x[t] - 0.0286 y[t],
> z'[t] == 0.0039 x[t] - 0.000035 z[t], x[0] == 0, y[0] == 0, z[0] == 0};
> 
Sol = Chop[Simplify[N[DSolve[eqns], {x[t], y[t], z[t]}, t]]], 10^(-7)]

out[2]:=

{{x[t] -> 1800.097 - 718.4499/E^(0.04466875*t) -  854.97/E^(0.0200*t) - 
     224.897/E^(0.000030632*t), 
   y[t] -> 698.6390 + 496.2920/E^(0.04466875163472879*t) - 
     1108.0984/E^(0.0200356*t) - 87.379/E^(0.000030632*t), 
   z[t] -> 200582.240249792 + 62.7765894/E^(0.04466875*t) + 
     166.7142446/E^(0.0200356*t) - 200811.9232/E^(0.0000306322*t)}}

	Other solution

> I have deleted the graphic display and most of the output.
> 
> 
> eqns = {x'[t] == -0.036 x[t] + 0.0124 y[t] + 0.000035 z[t] + 49.3,
>       y'[t] == 0.0111 x[t] - 0.0286 y[t],
>       z'[t] == 0.0039 x[t] - 0.000035 z[t], x[0] == 0, y[0] == 0, z[0] ==
> 0};
> 
> soln = DSolve[ eqns, {x[t], y[t], z[t]}, t];
> 
> I am not sure what caused the error message.
> We can get what they indicate by Help Browser > GoTo [RootSum::pfn] with
> Other Information chosen from the category buttons. We get
> "Generated when the first argument in RootSum is not a Function
> expression."
> 
> However it looks as if it is serious, since the solution is incorrect:
> 
> soln /. t -> 0
> 
> {{x[0] -> 1.77331, y[0] -> -0.546333, z[0] -> -0.192132}}
> 
> But, in reply to your question about the form of the answer:
> 
> RootSum[f, form]represents the sum of form[x]for allx that satisfy the
> polynomial equationf[x] == 0.
> This is very useful notation - for mathematics in general, not just for
> Mathematica.
> In this case we can convert to standard notation (but it takes a while)
> 
> ToRadicals[soln]
> 
> And, with another wait, to traditional notation
> 
> TraditionalForm[%]
> 
> Since your input includes inexact numbers, you might accept the muchs
> shorter form.
> 
> N[soln]
> 
> Or maybe a numerical solution over some range  of t (I have used -1 to 1)
> 
> Numsoln = NDSolve[eqns, {x[t], y[t], z[t]}, {t, -1, 1}]
> 
> This is OK at t=0:
> 
> Numsoln /. t -> 0
> 
> {{x[0] -> 0., y[0] -> 0., z[0] -> 0.}}
> 
> To check on the values throughout the t-range it is convenient to convert
> the solution to one for the functions x,y and z:
> 
> FunNumsoln = Nsoln /. f_[t] -> f
> 
> and, in eqns  change Equal to Subtract:
> 
> diffs = eqns /. Equal -> Subtract
> 
> Now plot diffs with the functions substitued for x, y, z  (with a perfect
> solution each would give constant zero function)
> 
> Plot[Evaluate[diffs /. FunNumsoln], {t, -1, 1}, PlotRange -> All];
> 
> The maximum absolute error looks to be  less than  0.000035.
> 
> 
> Allan
> ---------------------
> Allan Hayes
> Mathematica Training and Consulting
> Leicester UK
> www.haystack.demon.co.uk
> hay at haystack.demon.co.uk
> Voice: +44 (0)116 271 4198
> Fax: +44 (0)870 164 0565
> 
> "


  • Prev by Date: Re: Could this be improved?
  • Next by Date: Re: cellular automata & OOP stupid?
  • Previous by thread: Sturm-Liouville Problem (differential equation eigenvalue problem)
  • Next by thread: Re: cellular automata & OOP stupid?