Translating ODE from Matlab to Octave

classic Classic list List threaded Threaded
10 messages Options
Reply | Threaded
Open this post in threaded view
|

Translating ODE from Matlab to Octave

RICHARD Dominique
Hello everybody,
  I am switching from Matlab to Octave and I got some code that I am not sure how to deal with.
For example these lines:

options = odeset('Mass',M);
f = @(t,y) daesystem(t,y,params); 
[t,y] = ode15s(f,t,y0,options);

which I transformed to :

lsode_options ('Mass',M);

f = @(t,y) daesystem(t,y,params);

[y, t]= lsode(@(y,t) f, y0, t)


which does not seem to run properly.

I had a look at the Octave notice but besides the definitions, I did not found clear enough examples for me to use.
I would be very grateful for any help to get me to the free side of software.


-- 
Dominique RICHARD
Laboratoire de Génie des Procédés Catalytiques - CNRS-CPE Lyon
3, rue Victor Grignard F-69616 Villeurbanne cedex 
Tél. (33) 4 72 43 17 52 - Fax (33) 4 72 43 16 73

_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

c.-2
Hi,

On 26 Nov 2013, at 15:19, RICHARD Dominique <[hidden email]> wrote:

> Hello everybody,
>   I am switching from Matlab to Octave and I got some code that I am not sure how to deal with.
> For example these lines:
>
> options = odeset('Mass',M);
> f = @(t,y) daesystem(t,y,params);  
> [t,y] = ode15s(f,t,y0,options);
>
> which I transformed to :
>
> lsode_options
>       ('Mass',M);

I don't think "lsode_options" has an option called 'Mass',
have you looked at the help for lsode [1] or at the help for lsode in
the manual [2] before writing that command?
furthermore "ode15s" is often (although not exclusively) used to solve
Differential Algebraic Equations (DAE), i.e. ODEs with a singular mass matrix.
From that and from the fact that your function is named "daesystem" I infer that
probably your problem is a system of DAEs.

If this is the case you should not be using lsode but rather daspk[5, 6]
to solve your problem.

Otherwise, if your mass matrix is not singular
you can just redefine

  f = @(t,y) M \ daesystem(t,y,params);  

to take the mass matrix into account.

> f = @(t,y) daesystem(t,y,params);
> [y,
>       t]= lsode(@(y,t) f, y0, t)

that is not the correct syntax for calling lsode, [1] and [2] will explain how to use it.

>
> which does not seem to run properly.
>
> I had a look at the Octave notice but besides the definitions, I did not found clear enough examples for me to use.
do these examples [3, 4] help?

> I would be very grateful for any help to get me to the free side of software.

HTH,
c.


[1] type "help lsode" at the octave prompt to see the help
[2] http://www.gnu.org/software/octave/doc/interpreter/Ordinary-Differential-Equations.html#Ordinary-Differential-Equations
[3] http://wiki.octave.org/Cookbook#Parametrized_Functions
[4] http://wiki.octave.org/Bim_package#3D_Time_dependent_problem
[5] type "help daspk" at the octave prompt to see the help
[6] http://www.gnu.org/software/octave/doc/interpreter/Differential_002dAlgebraic-Equations.html#Differential_002dAlgebraic-Equations



_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

Bård Skaflestad
In reply to this post by RICHARD Dominique
On Tue, 2013-11-26 at 15:19 +0100, RICHARD Dominique wrote:

> f = @(t,y) daesystem(t,y,params);  
> [t,y] = ode15s(f,t,y0,options);
>
> which I transformed to :

[...]

> f = @(t,y) daesystem(t,y,params);
> [y, t]= lsode(@(y,t) f, y0, t)

Try f = @(y, t) daesystem(t, y, params)

The call

        lsode(@(y,t) f, ...)

makes no sense because 'f' is a function handle so

        @(y,t) f

is a function that, for all parameter pairs (y,t) returns the original
function handle ('f') rather than the value of the function (daesystem)
at (y,t).

Alternatively, you *could* use

        [y, t] = lsode(@(y, t) f(t, y), y0, t)

if you don't want to change the definition of 'f'.


Sincerely,
--
Bård Skaflestad <[hidden email]>
SINTEF ICT, Applied mathematics

_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

c.-2

On 26 Nov 2013, at 16:54, Bård Skaflestad <[hidden email]> wrote:

> Alternatively, you *could* use
>
>        [y, t] = lsode(@(y, t) f(t, y), y0, t)
>
> if you don't want to change the definition of 'f'.

no, you can't call 'lsode' like that.

'lsode' does not return the time vector 't', the calling sequence is different from ode15s.
the set of times at which you want your solution evaluate should be passed as an input parameter.

c.
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

Bård Skaflestad
On Tue, 2013-11-26 at 17:05 +0100, c. wrote:

> On 26 Nov 2013, at 16:54, Bård Skaflestad <[hidden email]> wrote:
>
> > Alternatively, you *could* use
> >
> >        [y, t] = lsode(@(y, t) f(t, y), y0, t)
> >
> > if you don't want to change the definition of 'f'.
>
> no, you can't call 'lsode' like that.
>
> 'lsode' does not return the time vector 't',

Thanks a lot for alerting me to that fact.  Now I just have to ask
myself why I didn't derive that conclusion from the documentation...

> the calling sequence is different from ode15s. the set of times at
>  which you want your solution evaluate should be passed as an input
>  parameter.

Are you saying that LSODE does not support automatic step size selection
or does the solver "simply" not expose the selected time steps to the
caller?


Sincerely,
--
Bård Skaflestad <[hidden email]>
SINTEF ICT, Applied mathematics

_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

c.-2

On 26 Nov 2013, at 17:18, Bård Skaflestad <[hidden email]> wrote:

>> no, you can't call 'lsode' like that.
>>
>> 'lsode' does not return the time vector 't',
>
> Thanks a lot for alerting me to that fact.  Now I just have to ask
> myself why I didn't derive that conclusion from the documentation...

The documentation [1] gives two possibilities for calling lsode:

 — Loadable Function: [x, istate, msg] = lsode (fcn, x_0, t)
 — Loadable Function: [x, istate, msg] = lsode (fcn, x_0, t, t_crit)

where 'istate' and æmsg' are flags describing the succes / failure of the integration. a

As you see there is no 't' in the list of outputs, maybe you just assumed it would work that way because matlab ode solvers do ...
'lsode' does not have a matlab compatible interface because there is no 'lsode' function in matlab, you can find other solvers with
more matlab-like interface in the odepkg forge package [2].

a matlab compatible ode15s is not implemented yet and it won't be in 3.8 either, but that is a planned feature that might be in 4.0.

> Are you saying that LSODE does not support automatic step size selection
> or does the solver "simply" not expose the selected time steps to the
> caller?

lsode does use internally adaptive time-stepping but only returns the solution evaluated at a set of requested time steps.

c.


[1] http://www.gnu.org/software/octave/doc/interpreter/Ordinary-Differential-Equations.html#Ordinary-Differential-Equations
[2] http://octave.sourceforge.net/odepkg/


_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

Torsten
On 26.11.2013 18:09, c. wrote:
> As you see there is no 't' in the list of outputs, maybe you just assumed it would work that way because matlab ode solvers do ...
> 'lsode' does not have a matlab compatible interface because there is no 'lsode' function in matlab, you can find other solvers with
> more matlab-like interface in the odepkg forge package [2].

And odepkg provides the mass matrix parameter "Mass" in case there might
be numerical issues with the inverse of M. I use ode45 on a regular
basis on octave and matlab systems (without code modifications).

Torsten

_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

c.-2

On 26 Nov 2013, at 18:53, Torsten <[hidden email]> wrote:

>
> And odepkg provides the mass matrix parameter "Mass" in case there might
> be numerical issues with the inverse of M. I use ode45 on a regular
> basis on octave and matlab systems (without code modifications).

ode45 uses the one-step explicit RK45  method which is not a viable choice for integrating stiff systems,
and even more so if the system to be integrated is a set of DAEs which may be roughly considered as
"infinitely stiff" ODEs.
Indeed I see in the source code of ode45 [3] that an inversion of the mass matrix is used, which cannot work for DAEs.

There are other matlab compatible solvers in odepkg that are intended for solving stiff problems, e.g. ode23s [1],
but I am quite sure ode23s won't work for DAEs, as I see in the code that an inversion of the mass matrix [2]  is used in the algorithm.

There are other solvers in odepkg that do work for DAEs [4] but I'm not sure whether they are matlab compatible.

The recomended solver to use for DAEs in matlab is ode15s [4], which is not (yet?) implemented in Octave or Octave Forge, but
the algorithm used by daspk is exactly the same as that of ode15s so I'd recommennd using daspk for solving DAEs in Octave.

> Torsten

c.


[1] http://octave.sourceforge.net/odepkg/function/ode23s.html
[2] http://sourceforge.net/p/octave/odepkg/ci/default/tree/inst/ode23s.m#l179
[3] http://sourceforge.net/p/octave/odepkg/ci/default/tree/inst/ode45.m#l342
[4] see the 5th table in the page http://www.mathworks.it/it/help/matlab/ref/ode15s.html
_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

Torsten
On 26.11.2013 19:33, c. wrote:

>
> On 26 Nov 2013, at 18:53, Torsten <[hidden email]> wrote:
>
>>
>> And odepkg provides the mass matrix parameter "Mass" in case there might
>> be numerical issues with the inverse of M. I use ode45 on a regular
>> basis on octave and matlab systems (without code modifications).
>
> ode45 uses the one-step explicit RK45  method which is not a viable choice for integrating stiff systems,
> and even more so if the system to be integrated is a set of DAEs which may be roughly considered as
> "infinitely stiff" ODEs.
> Indeed I see in the source code of ode45 [3] that an inversion of the mass matrix is used, which cannot work for DAEs.
>
> There are other matlab compatible solvers in odepkg that are intended for solving stiff problems, e.g. ode23s [1],
> but I am quite sure ode23s won't work for DAEs, as I see in the code that an inversion of the mass matrix [2]  is used in the algorithm.
>
> There are other solvers in odepkg that do work for DAEs [4] but I'm not sure whether they are matlab compatible.
>
> The recomended solver to use for DAEs in matlab is ode15s [4], which is not (yet?) implemented in Octave or Octave Forge, but
> the algorithm used by daspk is exactly the same as that of ode15s so I'd recommennd using daspk for solving DAEs in Octave.
>
>> Torsten
>
> c.
>
>
> [1] http://octave.sourceforge.net/odepkg/function/ode23s.html
> [2] http://sourceforge.net/p/octave/odepkg/ci/default/tree/inst/ode23s.m#l179
> [3] http://sourceforge.net/p/octave/odepkg/ci/default/tree/inst/ode45.m#l342
> [4] see the 5th table in the page http://www.mathworks.it/it/help/matlab/ref/ode15s.html
>

You are right.

I did not want to recommend ode45 for this specific simulation case. I
just wanted to confirm your statement that odepkg provides a matlab
compatible interface. Sorry if my post was misleading.

Torsten

_______________________________________________
Help-octave mailing list
[hidden email]
https://mailman.cae.wisc.edu/listinfo/help-octave
Reply | Threaded
Open this post in threaded view
|

Re: Translating ODE from Matlab to Octave

RICHARD_D
Thank you everybody for your help.
I was away from the lab at the end of last week and only had a look at the thread this morning.
A special thank to C;-2 for the detailed answer.
As a matter of fact the problem I was dealing with was ODE, but the Matlab code I received was ambiguous (ill written) and could be seen as a DAE problem, which was not the case.
By browsing through the examples (ref.[3] and [4]) I succeeded in getting the syntax right, which was what I was looking for..

Thanks again !

Dominique
Octave 3.2.4
Kubuntu 12.04 (LTS)