[python-control-discuss] Implementation of lsim, step, initial, impulse

[Eike W. <eik...@gm...>]

Hello!

I have written a patch for the Python Control Systems Library with the 
following contents:

* An implementation of the four simulation functions ``lsim``, ``step``, 
``initial``, and ``impulse`` of the module ``matlab``. 

* It also adds a function ``dcgain`` to the ``matlab`` module, which computes 
the gain of a linear system for steady state and constant input.

* The patch contains a bug fix for class ``StateSpace``, which enables it to 
work properly together with Scipy's ``signal`` module. 

* The simulation functions' return values are changed (back?) to arrays, 
because matrices confuse Matplotlib.

http://sourceforge.net/tracker/?func=detail&aid=3301246&group_id=282523&atid=1198311


Would you include my patch into the library? 

What do you think about the parameter convention below? IMHO a detailed 
section that explains the used convention should be included in the module's 
documentation. The documentation of the individual functions should link to 
it. 

What is the proper communication channel for the library's developers?

If you include my patch, my next step would be to extend ``lsim`` to MIMO 
systems (by copying the code from ``scipy.signal.lsim2``); and to convert the 
module's documentation to proper reStructuredText.


Implementation Notes
-------------------

The implementation is centered around ``lsim``, which checks the parameters 
thoroughly, and generates useful error messages. It then calls 
``scipy.signal.lsim2`` to solve the differential equations. The other 
simulation functions (``step``, ``initial``, ``impulse``) are small and call 
``lsim``.


Parameter Convention
--------------------

The new convention for parameters and return values harmonizes with the way 
how the system equations are written. It can be generalized to MIMO systems, 
and it also works well with Matplotlib:

All parameters can be arrays, matrices, or nested lists. 

The time vector is either 1D, or 2D with shape (1, n)::
 
      T = [[t1,     t2,     t3,     ..., tn    ]]  
      
Input, state, and output all follow the same convention. Columns are different 
points in time, rows are different components. When there is only one row, a 
1D object is accepted or returned, which adds convenience for SISO systems::

      U = [[u1(t1), u1(t2), u1(t3), ..., u1(tn)]
           [u2(t1), u2(t2), u2(t3), ..., u2(tn)]
           ...
           ...
           [ui(t1), ui(t2), ui(t3), ..., ui(tn)]]
      
      Same for X, Y

The initial conditions are either 1D, or 2D with shape (j, 1)::

     X0 = [[u1]
           [u2]
           ...
           ...
           [uj]]

So, U[:,2] is the system's input at the third point in time; and U[1] or 
U[1,:] is the sequence of values for the system's second input. 

As all simulation functions return array plotting is convenient::

    t, y = step(sys)
    plot(t, y)

The output of a MIMO system would be plotted like this::

    t, y, x = lsim(sys, u, t)
    plot(t, y[0], label='y_0')
    plot(t, y[1], label='y_1')


Yours,
Eike.