Problem using ode "rkf"
Bouchra BENSIALI
bouchra.bensiali at yahoo.fr
Wed Apr 28 15:42:33 CEST 2010
Hello world !
I have a problem using ode "rkf". I defined a function like that :
---
function Xdot = F(t,X)
[v1,v2,v3] = fonc(t,X(1),X(2));
xdot = v1;
ydot = v2;
Adot = v3;
Xdot = [xdot; ydot; Adot];
endfunction
---
then I used :
---
Sol = ode("rkf",X0,t0,t0:pas_temps:tf,F);
---
In this case, xdot and ydot don't depend on A=X(3), I obtain however a different result for x=X(1) and y=X(2) compared with the result obtained by integrating only x and y (same initial condition, etc.).
I noticed that it depended on how I define Adot. For example, doing this :
---
function Xdot = F(t,X)
printf("t=%f \n",t);
[v1,v2,v3] = fonc(t,X(1),X(2));
xdot = v1;
ydot = v2;
Adot = 1;
Xdot = [xdot; ydot; Adot];
endfunction
Sol = ode("rkf",X0,t0,t0:pas_temps:tf,F);
---
with t0=0.17, pas_temps=0.001 et tf=t0+pas_temps
leads to :
---
t=0.170000
t=0.170250
t=0.170375
t=0.170923
t=0.171000
t=0.170500
t=0.171000
---
I verified that these times correspond to the intermediate times used by the "rkf" method : t0, t0+1/4*pas_temps, t0+3/8*pas_temps, t0+12/13*pas_temps, t0+pas_temps, t0+1/2*pas_temps, as said here http://math.fullerton.edu/mathews/n2003/RungeKuttaFehlbergMod.html.
Now if I use for example :
---
function Xdot = F(t,X)
printf("t=%f \n",t);
[v1,v2,v3] = fonc(t,X(1),X(2));
xdot = v1;
ydot = v2;
Adot = v3;
Xdot = [xdot; ydot; Adot];
endfunction
Sol = ode("rkf",X0,t0,t0:pas_temps:tf,F);
---
then Scilab returns that :
---
t=0.170000
t=0.170250
t=0.170375
t=0.170923
t=0.171000
t=0.170500
t=0.170218
t=0.170326
t=0.170803
t=0.170870
t=0.170435
t=0.170870
t=0.170903
t=0.170919
t=0.170990
t=0.171000
t=0.170935
t=0.171000
---
and as a consequence that changes the results for x and y (they are supposed not to change because they don't depend on A). I checked the definition of v3 and don't see errors. For me, the number of intermediate times is supposed to be the same, am I right ?
An explanation of that mystery ?
Thank you in advance.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.scilab.org/pipermail/users/attachments/20100428/3d6623b1/attachment.htm>
More information about the users
mailing list