MathGroup Archive 2004

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

Search the Archive

Re: checking accuracy with stepwise ode.

  • To: mathgroup at smc.vnet.net
  • Subject: [mg48503] Re: checking accuracy with stepwise ode.
  • From: rknapp at wolfram.com (Rob Knapp)
  • Date: Wed, 2 Jun 2004 04:22:22 -0400 (EDT)
  • References: <c99dsj$k77$1@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

sean_incali at yahoo.com (sean kim) wrote in message news:<c99dsj$k77$1 at smc.vnet.net>...
> actually I figured it out, I hate it when I answer my own question
> after posting it. I didn't define the functions with [t_]. that's why
> the NDSOlve complained.
> 
> If i use the code below, it works fine, but How do I know the results
> are accurate? Does NDSolve do soemthing different with the stepwise
> functions that it doesn;t do with the normal functions?   how do I
> check my results?  how do I know the answer that i got is correct? I
> ask because the plot looks kinda funky.
> 
> 
> Thanks all in advance for any and all comments. 
> 
> sean 
> 
> k1 = 1/10; k2 = 1/20; 
> a0[t_] := 0 /; t < 0 ;
> a0[t_] := 1/10 /; 0 <= t <= 200 ;
> a0[t_] := 0 /; 200 <= t <= 600 ;
> a0[t_] := 1/10 /; 600 <= t<= 2000;
> 
> ndsolution = NDSolve[{b'[t] == -k2 b[t] y[t], x'[t] == -k1 a0[t]  x[t]
> + k2 b[t] y[t], y'[t] == k1 a0[t]  x[t] - k2 b[t] y[t], b[0] == 1,
> x[0] == 1, y[0] == 0}, {b, x, y}, {t, 0, 2000}][[1]] ;
> 
> Plot[Evaluate[{a0[t], b[t], x[t], y[t]} /. ndsolution], {t, 0, 2000},
> PlotStyle -> {
> {AbsoluteThickness[2], RGBColor[0, 0, 0]}, 
> {AbsoluteThickness[2], RGBColor[.7, 0, 0]}, 
> {AbsoluteThickness[2],RGBColor[0, .7, 0]}, 
> {AbsoluteThickness[2], RGBColor[0, 0, .7]}}, 
> PlotRange -> All, Axes -> False, Frame -> True, 
> PlotLabel -> StyleForm[A StyleForm[" B", FontColor -> RGBColor[.7, 0,
> 0]] StyleForm[" X", FontColor -> RGBColor[0, .7, 0]]StyleForm["Y",
> FontColor -> RGBColor[0, 0, .7]],
> FontFamily -> "Helvetica", FontSize -> 12, FontWeight -> "Bold"]];

The methods that are used in NDSolve will not give you full accuracy
except for smooth functions.  To get full accuracy from a solution to
the problem above, you need to tell NDSolve where the discontinuities
are.  When you do this, NDSolve makes sure to start a new step at the
point of discontinuity.  Further recommended when discontinuities are
present is to use "one-step" methods, such as Runge-Kutta methods,
since there is no cost to starting up again after a discontinuity.

For you example: 

Find the discontinuiity points once you have defined a0 (you could
type them in, but I was too lazy for this)

In[6]:=
discontinuitypoints = 
 Union[Cases[
   Flatten[Cases[
     DownValues[a0],(Less | LessEqual)[x__]\[RuleDelayed]{x}, {0,
      Infinity}]], _?NumberQ]]

Out[6]=
{0,200,600,2000}

The discontinuity points are included in the argument to NDSolve which
gives the range for the independent variable, in this case with
Join[{t,0},discontinuitypoints, {2000}]

In[7]:=
sol = NDSolve[{b'[t] == -k2 b[t] y[t],
   x'[t] == -k1 a0[t] x[t]+k2 b[t] y[t],
   y'[t] == k1 a0[t] x[t]-k2 b[t] y[t],b[0] == 1,x[0] == 1,
   y[0] == 0},{b,x,y},Join[{t,0},discontinuitypoints, {2000}], 
  Method->ExplicitRungeKutta]

Out[7]=
{{b->InterpolatingFunction[{{0.,2000.}},<>],
  x->InterpolatingFunction[{{0.,2000.}},<>],
  y->InterpolatingFunction[{{0.,2000.}},<>]}}

Rob Knapp
Wolfram Research


  • Prev by Date: Re: Extract substrings using metacharacters
  • Next by Date: RE: RE: complex analysis problem in mathematica 3.0
  • Previous by thread: Re: Eigensystem[] bug in Mathematica 5.0 (fixed in 5.0.1)
  • Next by thread: limits on symbol eigenvalues?