On 9-Aug-2008, Torquil Macdonald Sørensen wrote:
| Marc Normandin wrote:
| > Torquil Macdonald Sørensen wrote:
| >> Hi, I'm getting very different results when solving the following
| >> initial value ODE problem in Matlab and Octave:
| >>
| >> dy/dt=1/sqrt(y^2 + 1)+y-y^2 on t \in [0,10] with y(0) = 0
| >>
| >> From looking at the equation, I believe that the Matlab solution is the
| >> correct one, so I'm wondering if I have not converted the Matlab-file
| >> correctly to Octave? Or does the Octave algorithm not work in this
| >> situation? It does not change when I lower the values for the absolute
| >> and relative tolerances with lsode_options, so maybe I have done
| >> something wrong and I'm solving a different ODE in Octave.
| >>
| >> Here are the Octave and Matlab programs that I thought were solving the
| >> same equation:
| >>
| >>
| >> *** Octave version ***
| >>
| >> function nonlinear_de
| >>
| >> % Time span
| >> t = linspace(0,10,100);
| >>
| >> % Initial condition
| >> y0=[0];
| >>
| >> % Solve the DE
| >> lsode_options('absolute tolerance', 0.0001);
| >> [y, istate, msg] = lsode("ode",y0,t);
| >>
| >> istate
| >> msg
| >>
| >> % Plot the solution
| >> plot(t,y(:,1),'-')
| >>
| >> % Differential equation
| >> function dydt = ode(t,y)
| >> dydt = [ 1/sqrt(y(1)^2 + 1)+y(1)-y(1)^2];
| >>
| >>
| >> *** Matlab version ***
| >>
| >> function nonlinear_de
| >>
| >> % Time span
| >> tspan=[0 10];
| >>
| >> % Initial condition
| >> y0=[0];
| >>
| >> % Solve the DE
| >> [t,y] = ode45(@ode,tspan,y0);
| >>
| >> % Plot the solution
| >> plot(t,y(:,1),'-')
| >>
| >> % Differential equation
| >> function dydt = ode(t,y)
| >> dydt = [ 1/sqrt(y(1)^2 + 1)+y(1)-y(1)^2];
| >> _______________________________________________
| >> Help-octave mailing list
| >> address@hidden
| >> https://www.cae.wisc.edu/mailman/listinfo/help-octave
| >
| > The arguments of your differential function should be swapped for the
| > code using lsode. From "help lsode":
| >
| > The first argument, FCN, is a string, or cell array of strings,
| > inline or function handles, that names the function to call to
| > compute the vector of right hand sides for the set of equations.
| > The function must have the form
| >
| > XDOT = f (X, T)
| >
| > Changing your Octave code to read
| > function dydt = ode(y,t)
| > instead of
| > function dydt = ode(t,y)
| > gives results in agreement with ode45 using the other program.
| >
| > You might also want to look at odepkg in Octave-Forge, it includes ode45
| > and other handy solvers.
| >
|
| Thanks Marc, for some reason I didn't spot that one. Yes, I had found
| that there was an octave-odepkg package in Debian which contained ode45,
| among others, and gave the same results as Matlab after also specifying
| some sensible tolerance values.
Octave's lsode function should also give good results if you use it
properly. Does it?
jwe