[Date Index]
[Thread Index]
[Author Index]
Re: NDSolve does not work in shooting method
*To*: mathgroup at smc.vnet.net
*Subject*: [mg44159] Re: NDSolve does not work in shooting method
*From*: Rob Knapp <rknapp at wolfram.com>
*Date*: Fri, 24 Oct 2003 04:24:46 -0400 (EDT)
*References*: <bn8du7$nfi$1@smc.vnet.net>
*Sender*: owner-wri-mathgroup at wolfram.com
Alexander Koulidis wrote:
> This comes from a newbie in Mathematica, so have some mercy....
>
> I'm trying to reproduce some examples of the shooting method that have been
> published in this forum, but Mathematica 5.0 shows some strange behavior
> with NDSolve.
>
> In http://forums.wolfram.com/mathgroup/archive/2001/Aug/msg00128.html Brian
> wrote:
>
> The following code implements the shooting technique
>
> system[a_] := {y1'[x] == y2[x], y2'[x] == -Sin[y2[x]] + Cos[5x],
> y1[-1] == 0,
> y2[-1] == a};
> myODEsoln[a_] := NDSolve[system[a], {y1[x], y2[x]}, {x, -1, 1}]
> yend[a_] := (y1[x] /. myODEsoln[a]) /. x -> 1
> bc = FindRoot[First[yend[a]] == 0, {a, -2, 2}];
> Plot[Evaluate[y1[x] /. myODEsoln[a /. bc]], {x, -1, 1},
> AxesLabel -> {"x", "y1(x)"}];
>
>
> While the code was working fine in v. 4.2, now it produces warnings of
> "NDSolve::ndinnt :" or "RuleDelayed::rhs".
> Can someone test this code in v5.0? Does anyone have some idea of what is
> happening and should I implement a typical shooting method?
> Thanks in advance.
> Alexander Koulidis, Athens, Greece
>
>
The strange behavior is not in NDSolve, but arises from the way
expressions are evaluated, in particular in FindRoot. As of Mathematica
5, what FindRoot (also FindMinimum, NIntegrate, NSum) does is give the
variables local definitions, then evaluates the function (or equations
for FindRoot) with the variables having no values: i.e. symbolically.
This means that in your example, the equation FindRoot is working with
is not that you were expecting:
In[1]:= system[a_] := {y1'[x] == y2[x], y2'[x] == -Sin[y2[x]] +
Cos[5x],y1[-1] == 0, y2[-1] == a};
In[2]:= myODEsoln[a_] := NDSolve[system[a], {y1[x], y2[x]}, {x, -1, 1}];
In[3]:= yend[a_] := (y1[x] /. myODEsoln[a]) /. x -> 1;
In[4]:= yend[a]
.... various messages elided ...
Out[4]= y1[1] /.
> NDSolve[{y1'[1] == y2[1], y2'[1] == Cos[5] - Sin[y2[1]], y1[-1] == 0,
> y2[-1] == a}, {y1[1], y2[1]}, {1, -1, 1}]
In[5]:= First[%]
Out[5]= y1[1]
.... so the equation FindRoot tries to work with is
y1[1] == 0, which is not really what you want to solve.
The real problem here is that the function yend does not do what you
intend unless the argument a is a number. If you simply attach this
restriction to the definition to delay evaluation until it is
appropriate to evaluate, then it works as you expect:
In[6]:= Clear[yend]
In[7]:= yend[a_?NumberQ] := First[(y1[x] /. myODEsoln[a]) /. x -> 1]
Now yend does not evaluate with symbolic a
In[8]:= yend[a]
Out[8]= yend[a]
(Note I included the First in the definition of yend, since
First[yend[a]] would evaluate to a with non number values of a)
It does evaluate with numbers.
In[9]:= yend[1]
Out[9]= 0.715203
.... and FindRoot works just fine with this definition.
In[10]:= FindRoot[yend[a] == 0, {a, -2,2}]
Out[10]= {a -> 0.212199}
Prev by Date:
**Re: Re: InterpolatingFunctionAnatomy**
Next by Date:
**Baffling Failure when plotting multiple curves.**
Previous by thread:
**Re: NDSolve does not work in shooting method**
Next by thread:
**[no subject]**
| |