Re: Trouble with NDSolve Function
- To: mathgroup at smc.vnet.net
- Subject: [mg5420] Re: Trouble with NDSolve Function
- From: saarinen (Sirpa Saarinen)
- Date: Sat, 7 Dec 1996 00:25:49 -0500
- Organization: Wolfram Research, Inc.
- Sender: owner-wri-mathgroup at wolfram.com
wdavis at csranet.com (Wayne Davis) writes:
>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.}, <>]}}
Here we need to re-write part of your equations. If you are using V 2.2
you can use:
imax=3;ya=100.;k=.5;c=.1;tmax=2.;
SetAttributes[y,NProtectedFirst];
ClearAttributes[If, HoldAll];
q[1][t_]:=If[t<1,10.,0];
q[imax+1][t_]:=0.;
q[i_][t_]:=If[y[i-1][t]>ya,q[i-1][t],k (y[i-1][t]-y[i][t])];
difeqs=Table[c y[i]'[t]==q[i][t]-q[i+1][t],{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}]
and if you are using V 3.0:
imax=3;ya=100.;k=.5;c=.1;tmax=2.;
ClearAttributes[If, HoldRest];
q[1][t_]:=If[t<1,10.,0];
q[imax+1][t_]:=0.;
q[i_][t_]:=If[y[i-1][t]>ya,q[i-1][t],k (y[i-1][t]-y[i][t])];
difeqs=Table[c y[i]'[t]==q[i][t]-q[i+1][t],{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}]
Part of the problem was the If-statement and part of the problem lay in the fact
that the q-function used no argument (t) in its definition.
I hope this helps,
Sirpa
saarinen at wolfram.com