Re: Trouble with NDSolve Function
- To: mathgroup at smc.vnet.net
- Subject: [mg51089] Re: Trouble with NDSolve Function
- From: pein <petsie at arcor.de>
- Date: Mon, 4 Oct 2004 06:17:58 -0400 (EDT)
- References: <57tm0o$1di@dragonfly.wolfram.com> <cjojb7$aq6$1@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Terry schrieb:
> I was interested to see your 1996 Mathematica query because I have a
> similar problem in 2004. I could not find any responses to your query
> on the interent so I'm taking the liberty of sending this email to ask
> if you managed to solve the problem and, if so, how?
>
> Any feedback would be most welcome.
>
> Terry Fairclough
>
>
> On 2 Dec 1996 04:30:16 GMT, Wayne Davis wrote:
>
>>I have had some trouble getting the Mathematica NDSolve function to
>
> accept
>
>>and solve problems containing certain types of conditionals in the
>
> dif eq,
>
>>even though the problem seems properly stated to me. I am attaching
>
> an
>
>>example and would like to know what I'm doing wrong. (It appears that
>
> the
>
>>assignments don't get substituted into the conditional terms so that
>>NDSolve can handle them as numerical values).
>>
>>imax=3;ya=100.;k=.5;c=.1;tmax=2.;
>>q[1]=If[t<1,10.,0];
>>q[imax+1]=0.;
>>q[i_]=If[y[i-1][t]>ya,q[i-1],k (y[i-1][t]-y[i][t])];
>>difeqs=Table[c y[i]'[t]==q[i]-q[i+1],{i,1,imax}];
>>ics=Table[y[i][0]==25,{i,1,imax}];
>>eqs=Flatten[Join[difeqs,ics]];
>>vars=Table[y[i],{i,1,imax}];
>>soln=NDSolve[eqs,vars,{t,0,tmax}]
>>NDSolve::ndnum:
>> Differential equation does not evaluate to a number at
>> t = 0..
>>{{y[1] -> InterpolatingFunction[{0., 0.}, <>],
>>
>> y[2] -> InterpolatingFunction[{0., 0.}, <>],
>>
>> y[3] -> InterpolatingFunction[{0., 0.}, <>]}}
>>
>>Thanks in advance for any help.
>>
>>Don Paddleford
>>Westinghouse Savannah River Company
>>
>>--
>>Wayne Davis Q: What has four legs and one arm?
>>wdavis at csranet.com A: A happy Pit Bull
>
>
If you rewrite the conditional parts of the definition of q with the
help of UnitStep, there's no problem (in vers. 4.00):
imax = 3; ya = 100; k = 1/2; c = 1/10; tmax = 2;
Clear[q, y];
q[1, t_] := If[t < 1, 10, 0];
q[imax + 1, t_] = 0;
q[i_, t_] :=
#q[i - 1, t] + (1 - #)k (y[i - 1][t] - y[i][t]) &@
UnitStep[y[i - 1][t] - ya];
difeqs = Table[ y[i]'[t] == (q[i, t] - q[i + 1, t])/c, {i, 1, imax}];
ics = y[#][0] == 25 & /@ Range[imax];
eqs = Join[difeqs, ics];
vars = y /@ Range[imax];
soln = NDSolve[eqs, vars, {t, 0, tmax}]