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

• Prev by Date: Smoothing Filters
• Next by Date: init.m question
• Previous by thread: Trouble with NDSolve Function
• Next by thread: Re: XMathematica: can't open display, exiting...