matrix exponential used by the control package

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

matrix exponential used by the control package

Torsten Lilge
Hi everyone,

the control package uses the SLICOT function mb05nd [1]
(__sl_mb05nd__()) instead of the octave function expm() for the matrix
exponential. The problem with mb05nd is that when discretizing a static
gain in state space results in

>> [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
a = [](1x0)
b = 0
c = [](1x0)
d = 1

Instead, using expm() leads to the correct result

>> [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
a = [](0x0)
b = [](0x1)
c =
[](1x0)
d = 1

Are there any reasons against using expm()?

Torsten


[1] http://slicot.org/objects/software/shared/doc/MB05ND.html


Reply | Threaded
Open this post in threaded view
|

Re: matrix exponential used by the control package

Doug Stewart-4


On Mon, Jul 20, 2020 at 4:38 PM Torsten Lilge1960 Blue Heron Drive, <[hidden email]> wrote:
Hi everyone,

the control package uses the SLICOT function mb05nd [1]
(__sl_mb05nd__()) instead of the octave function expm() for the matrix
exponential. The problem with mb05nd is that when discretizing a static
gain in state space results in

>> [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
a = [](1x0)
b = 0
c = [](1x0)
d = 1

Instead, using expm() leads to the correct result

>> [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
a = [](0x0)
b = [](0x1)
c =
[](1x0)
d = 1

Are there any reasons against using expm()?

Torsten


[1] http://slicot.org/objects/software/shared/doc/MB05ND.html


Is there a difference in speed?
If it is slower then we should look and see where the '0" is coming from.

--
DASCertificate for 206392

Reply | Threaded
Open this post in threaded view
|

Re: matrix exponential used by the control package

Torsten Lilge
On Mon, 2020-07-20 at 17:22 -0400, Doug Stewart wrote:

>
>
> On Mon, Jul 20, 2020 at 4:38 PM Torsten Lilge1960 Blue Heron Drive, <
> [hidden email]> wrote:
> > Hi everyone,
> >
> > the control package uses the SLICOT function mb05nd [1]
> > (__sl_mb05nd__()) instead of the octave function expm() for the
> > matrix
> > exponential. The problem with mb05nd is that when discretizing a
> > static
> > gain in state space results in
> >
> > >> [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
> > a = [](1x0)
> > b = 0
> > c = [](1x0)
> > d = 1
> >
> > Instead, using expm() leads to the correct result
> >
> > >> [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
> > a = [](0x0)
> > b = [](0x1)
> > c =
> > [](1x0)
> > d = 1
> >
> > Are there any reasons against using expm()?
> >
> > Torsten
> >
> >
> > [1] http://slicot.org/objects/software/shared/doc/MB05ND.html
> >
> >
>
> Is there a difference in speed?
> If it is slower then we should look and see where the '0" is coming
> from.


__sl_mb05nd__ is a compiled function and - at least on my system -
approx. 20 times faster.

Torsten



Reply | Threaded
Open this post in threaded view
|

Re: matrix exponential used by the control package

Torsten Lilge
On Tue, 2020-07-21 at 07:22 +0200, Torsten Lilge wrote:

> On Mon, 2020-07-20 at 17:22 -0400, Doug Stewart wrote:
> >
> > On Mon, Jul 20, 2020 at 4:38 PM Torsten Lilge1960 Blue Heron Drive,
> > <
> > [hidden email]> wrote:
> > > Hi everyone,
> > >
> > > the control package uses the SLICOT function mb05nd [1]
> > > (__sl_mb05nd__()) instead of the octave function expm() for the
> > > matrix
> > > exponential. The problem with mb05nd is that when discretizing a
> > > static
> > > gain in state space results in
> > >
> > > > > [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
> > > a = [](1x0)
> > > b = 0
> > > c = [](1x0)
> > > d = 1
> > >
> > > Instead, using expm() leads to the correct result
> > >
> > > > > [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
> > > a = [](0x0)
> > > b = [](0x1)
> > > c =
> > > [](1x0)
> > > d = 1
> > >
> > > Are there any reasons against using expm()?
> > >
> > > Torsten
> > >
> > >
> > > [1] http://slicot.org/objects/software/shared/doc/MB05ND.html
> > >
> > >
> >
> > Is there a difference in speed?
> > If it is slower then we should look and see where the '0" is coming
> > from.
>
> __sl_mb05nd__ is a compiled function and - at least on my system -
> approx. 20 times faster.
>

Hello Doug,

The routine mb05nd does not accept input dimensions 0x0 for an empty
matrix but it accepts (1x0) which is forced by the wrapper sl_mb05nd for
empty matrices. I have therefore  pushed changeset 7a8204f62fe1 to the
control package which adds the following test to sl_mb05nd:

+        if (a.numel () == 0)
+          {
+            // given matrix is empty, just return the empty matrix,
+            // otherise mb05d would make a [](1x0) matrix from it
+            retval(0) = a;
+            retval(1) = a;
+            return retval;
+          }

I have also prepared a patch for adding the first order-hold method
'foh' to c2d. This enables us to use 'foh' in lsim() or ramp(). Should I
push this to the default branch or should I create a feature branch for
this?

Torsten







Reply | Threaded
Open this post in threaded view
|

Re: matrix exponential used by the control package

Doug Stewart-4


On Sat, Jul 25, 2020 at 10:28 AM Torsten Lilge <[hidden email]> wrote:
On Tue, 2020-07-21 at 07:22 +0200, Torsten Lilge wrote:
> On Mon, 2020-07-20 at 17:22 -0400, Doug Stewart wrote:
> >
> > On Mon, Jul 20, 2020 at 4:38 PM Torsten Lilge1960 Blue Heron Drive,
> > <
> > [hidden email]> wrote:
> > > Hi everyone,
> > >
> > > the control package uses the SLICOT function mb05nd [1]
> > > (__sl_mb05nd__()) instead of the octave function expm() for the
> > > matrix
> > > exponential. The problem with mb05nd is that when discretizing a
> > > static
> > > gain in state space results in
> > >
> > > > > [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
> > > a = [](1x0)
> > > b = 0
> > > c = [](1x0)
> > > d = 1
> > >
> > > Instead, using expm() leads to the correct result
> > >
> > > > > [a,b,c,d] = ssdata (c2d (ss (tf (1,1)),1))
> > > a = [](0x0)
> > > b = [](0x1)
> > > c =
> > > [](1x0)
> > > d = 1
> > >
> > > Are there any reasons against using expm()?
> > >
> > > Torsten
> > >
> > >
> > > [1] http://slicot.org/objects/software/shared/doc/MB05ND.html
> > >
> > >
> >
> > Is there a difference in speed?
> > If it is slower then we should look and see where the '0" is coming
> > from.
>
> __sl_mb05nd__ is a compiled function and - at least on my system -
> approx. 20 times faster.
>

Hello Doug,

The routine mb05nd does not accept input dimensions 0x0 for an empty
matrix but it accepts (1x0) which is forced by the wrapper sl_mb05nd for
empty matrices. I have therefore  pushed changeset 7a8204f62fe1 to the
control package which adds the following test to sl_mb05nd:

+        if (a.numel () == 0)
+          {
+            // given matrix is empty, just return the empty matrix,
+            // otherise mb05d would make a [](1x0) matrix from it
+            retval(0) = a;
+            retval(1) = a;
+            return retval;
+          }

I have also prepared a patch for adding the first order-hold method
'foh' to c2d. This enables us to use 'foh' in lsim() or ramp(). Should I
push this to the default branch or should I create a feature branch for
this?

Torsten

I am sorry but I ave been sick (nothing too bad) and did not look at this yet :-(
I am better now!

I would push this to the default branch .
After you push this I hope to do more testing using 6.1
Doug


--
DASCertificate for 206392