Re: Matlab-ode45 vs Octave-lsode for a nonlinear ODE

[John W. Eaton]

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