Control and polynomial functions

[Federico Zenith]

Hi all,
I am a PhD student at NTNU Trondheim (Norway). My tutoring professor, Sigurd 
Skogestad, recently asked me to present him Scilab and Octave as alternative 
tools to Matlab in process control. I had already used Scilab in that context 
before, and knew it was working. Given the quantity of control-related 
functions that I had already looked into while working on my documentation 
patch, I thought that it would behave just as well.

However, Octave is sadly using an out-of-date way of representing linear 
systems. It is painstakingly difficult only to make a feedback loop, which 
involves multiple operations like connecting outputs and inputs together. 
Instead of using overloaded +, -, *, / and \, a number of ad hoc functions 
are present, like sysadd(), sysmult() and syssub(). I had no choice but to 
tell him that Scilab and Matlab outperform Octave in control functions. I 
thought this was sad because Octave is my favourite, since it's the only one 
of the three which is FOSS (the Scilab licence is technically neither free 
nor open-source - it's free of charge for non-commercial applications).

This out-of-date representation of systems stems partly from the out-of-date 
representation of polynomials; I am convinced that polynomials should be 
fixed first so that control functions can later be revamped, since using 
polynomial matrices in defining systems is a very handy approach used in 
Scilab. Here's my idea:

1- Imitate Scilab polynomials. They are handy. They define a polynomial 
object, which is used to define a "rational" object (ratio of two 
polynomials, a list with structure ( [info], num,den)) and "rational" 
matrices. (This is at least the idea I got, I still have to dig in their 
code)
The net result is that, instead of calling functions to sum or multiply 
polynomials, the user simply uses + and *, which is much more intuitive.
(I'm not sure how Matlab 7.0 does this. Worth checking out?)
For example, defining a polynomial in Scilab is as follows:

s = poly(0,"s")  // s is a polynomial whose variable is "s", with a single
    // root in 0; therefore the result is "s".
    // This is the most common definition of a polynomial,
    // and is often found in startup scripts.

y = poly([0 2], "x") // y is a polynomial whose variable is "x", with two
    // roots in 0 and 2. The result is x^2 - 2 * x

y = poly([0 2], "x", "coeff")
    // This time, instead of the roots, the coefficients are
    // defined, and the result is 2 * x.

2- Once this is implemented, such polynomial matrices may be used to define 
linear systems. It is much easier to define a feedback loop with 
T=inv(1+GK)*GK instead of hacking around with connections. This is at least 
how control people like to do it today. In general, it would be wise to 
imitate the current Matlab 7.0 implementation, which in turn took many ideas 
from Scilab 2.6 and later ("make people work with systems, not jungles of 
matrices").

Operating plan:
1- Define a structure for (rational) polynomials;
2- Implement a user-defined data structure (maybe writing the missing 
documentation in octave.dvi while we're at it) for rational polynomials 
matrices;
3- Rewrite the polynomial functions to accomodate for the new standard. They 
are not incredibly many, just 1237 lines of code with comments and all;
4- Make a patch and submit to sources mailing list.
At some later point... - start working on the control functions.

I'd like to know whether there are others interested in this approach. I'm not 
the definitive 1337 h4x0r, and someone more experienced may be of help.
I'm currently engaged in the Italian translation of the KDE project, which 
takes a lot of my spare time, but I could evaluate shifting to Octave 
development (it would be about time I move to something more challenging than 
translations, even though it improved my Italian orthography..., ehm!), 
especially since it's definitely more relevant for my job.

Is there someone who would like to start a small work group for the polynomial 
issue?

Cheers,
-Federico

pgpUQtNGEC5Q4.pgp
Description: PGP signature