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.
>
>
>
> 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?