Re: Control Function With NDsolve

• To: mathgroup at smc.vnet.net
• Subject: [mg19386] Re: Control Function With NDsolve
• From: "Allan Hayes" <hay at haystack.demon.co.uk>
• Date: Mon, 23 Aug 1999 13:57:06 -0400
• References: <7p5dm5\$127@smc.vnet.net> <7pl53k\$c74@smc.vnet.net>
• Sender: owner-wri-mathgroup at wolfram.com

```Ekhardt,

A version ot your solution. I have added two print commands to show where
the switching occurs. One problem is that the accuracy of the switching is
dependent on the step sizes.

f[y_Real] := Which[
y >= ymax && fval == 1, Print[y]; fval = 0,
y <= ymin && fval == 0, Print[y]; fval = 1,
True, fval
]

ymax = 0.9; ymin = 0.1; a = 1; b = 1;

eqs = {y'[t] == a*f[y[t]] - b*y[t], y[0] == 0};

fval = 1; NDSolve[eqs, y[t], {t, 0, 15}];

0.903906
0.0810061
0.913014
0.0807329
0.912961
0.0807266
0.90911

y1[t_] = y[t] /. First[%];

Plot[y1[t], {t, 0, 15}, PlotRange -> All]

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

Eckhard Hennig <hennig at itwm.uni-kl.de> wrote in message
news:7pl53k\$c74 at smc.vnet.net...
>
> Don Paddleford schrieb in Nachricht <7p5dm5\$127 at smc.vnet.net>...
> >In solving a control type dif eq with NDSolve I have the following
> >question. Suppose the eq is of the following simplified form
> >
> > y'[t]==a*f[y[t]]-b*y[t]
> > y[0]==0
> >
> >How to define f so that it starts at f=1, and changes to f=0 when y
> >reaches ymax, and then changes back to f=1 when y reaches ymin, and so
> >on in oscilatory fashion?
> >
>
> Don,
>
> you can define such a function as follows. Note that it is important to
> define the pattern for f such that it applies only to numeric arguments.
> Otherwise, f[y[t]] would be evaluated prematurely in In[3].
>
> In[1]:= f[y_Real] :=
>           If[(y > ymax && fval == 1) || (y < ymin && fval == 0),
>              fval = 1 - fval,
>              fval]
>
> In[2]:= ymax = 0.9; ymin = 0.1; a = 1; b = 1;
>
> In[3]:= eqs = {y'[t] == a*f[y[t]] - b*y[t], y[0] == 0};
>
> In[4]:= fval = 1; NDSolve[eqs, y[t], {t, 0, 10}];
>
> In[5]:= y1[t_] = y[t] /. First[%];
>
> In[6]:= Plot[y1[t], {t, 0, 10}, PlotRange->All]
>
>
> -- Eckhard
>
> -----------------------------------------------------------
> Dipl.-Ing. Eckhard Hennig      mailto:hennig at itwm.uni-kl.de
> Institut fuer Techno- und Wirtschaftsmathematik e.V. (ITWM)
> Erwin-Schroedinger-Strasse,  67663 Kaiserslautern,  Germany
>   Voice: +49-(0)631-205-3126    Fax: +49-(0)631-205-4139
>     http://www.itwm.uni-kl.de/as/employees/hennig.html
>
>         http://www.itwm.uni-kl.de/as/products/ai
> -----------------------------------------------------------
>
>
>
>
>

```

• Prev by Date: Re: ContourShading with ListContourPlot
• Next by Date: Re: Simulation
• Previous by thread: Re: Control Function With NDsolve
• Next by thread: Re: Control Function With NDsolve