Scilab Recipe 2: with Control Engineering Flavor

A couple of years ago, I wrote Scilab Recipe 1: Introduction to Scilab/Scicos. By the time people begin to find it on Google, some information has become outdated. For example, the file extension for an xcos model now changes from .xcos to .zcos. The Scilab user interface has also been improved noticeably. The basics and syntax for Scilab commands have not changed much though, so the article could still serve as an introduction to this open-source software.

As a dude whose main interests are in control systems, I now decide to write Recipe 2, tailored more to this particular field of study. My first encounter with Scilab was dated several years back when I read some online course materials from Prof. Alex Megretski website. He introduced this software in place of Matlab. At that time it was quite immature and not user-friendly. I recall one time I was so clumsy connecting a signal line to a port in a simulation model in front of class that a student asked why I bothered to use this at all.

Frankly, I began to feel comfortable with Scilab and Xcos after version 5.3. Now it has become my main computing software for class and research. Indeed, if the commands and usage are not quite alike, I could forget how to use Matlab/Simulink entirely.

So, let’s discuss some common task for control system analysis and design. It is suitable for an introductory document to limit the scope to only contiuous-time Linear Time-Invariant (LTF) systems, such that the plant and controller can be described as transfer functions in the Laplace domain. We also consider only Single-Input-Single-Output (SISO) systems.

Transfer Function Representation

Note: For a more comprehensive treatment of achieving a transfer function of a dynamical system, read Module 1 of our Scilab Control Engineering Basics.

An LTI plant can be modeled from physics using differential equation mathematics, then converted to a transfer function by Laplace transform. This is where we start our discussion. Given a plant represented by

(1)   \begin{equation*}  P(s) = \frac{2}{s(s+0.5)} \end{equation*}

which is a simplified example of permanent magnet DC motor. To construct this plant in Scilab workspace, start with a seed s representing Laplace variable.

-->s=poly(0,'s')
 s  =
    s

Then we can build on s any polynomial or rational function. For (1),

-->P = 2/(s*(s+0.5))
 P  = 
       2       
    -------    
            2  
    0.5s + s

Those who are familiar with MATLAB may notice that the power of s is displayed in reverse fashion in Scilab; i.e., starting from lowest order term. You somehow have to get used to this, especially when representing the numerator and denominator in vector form.

At this point, our plant appears like a rational function , but is not yet a linear transfer function that Scilab could interpret. We need to do one more step by using syslin command

-->P = syslin('c',P)
 P  =
       2       
    -------    
            2  
    0.5s + s

We do this in two steps to make it clear. Next time you can save time by generating a transfer function in one step.

Now we can observe Bode frequency response of P

-->bode(P)

which yields Figure 1.



Figure 1 Bode plot of P(s)

Some other useful commands that you can try by your own

  • nyquist(P) : generate a Nyquist plot
  • plzr(P): plot pole/zero mapping in complex plane

Now, we can do standard computation with transfer functions in Scilab. Say, suppose we have a controller

(2)   \begin{equation*}  C(s) = \frac{s+1}{s+10} \end{equation*}

Construct it in Scilab

-->C=syslin('c',(s+1)/(s+10))
 C  =
    1 + s    
    -----    
    10 + s

To compute a loop transfer function L(s) = C(s)P(s)

-->L = C*P
 L  = 
        2 + 2s       
    --------------   
              2   3  
    5s + 10.5s + s

Likewise, for the sensitivity S(s) = 1/(1+L(s)) and complementary sensitivity T(s) = L(s)/(1+L(s))

-->S = 1/(1+L)
 S  =
                2   3    
      5s + 10.5s + s     
    ------------------   
                  2   3  
    2 + 7s + 10.5s + s   
 
-->T = L/(1+L)
 T  =
          2 + 2s         
    ------------------   
                  2   3  
    2 + 7s + 10.5s + s

Stability of closed-loop system can be determined from the poles of S(s) or T(s). That raises a question. How to compute the pole of a transfer function? MATLAB users would think of a convenient command poles. Surprisingly enough, there is no such command in Scilab! I wonder how the developers forget a function so useful. Or does it come with other name I still cannot find?

Anyway, the way I used to get around is to convert the data to state-space format, say, for T(s)

-->Tss = tf2ss(T);

Then we can check the eigenvalue of matrix A in Tss by command spec

-->spec(Tss.A)
 ans  = 
  - 9.8070205               
  - 0.3464898 + 0.2896211i  
  - 0.3464898 - 0.2896211i

Since all poles have negative real parts, this feedback system is stable

Xcos Simulation

With this plant and controller, we now want to simulate how the closed-loop system behave in time-domain. The most common to observe is a response to step command. So we have to create an Xcos simulation model. Type xcos at Scilab command prompt to open a blank untitled window. Then open the Palette browser from View menu, find and drag the following palettes to the blank Xcos model.

  • STEP_FUNCTION and CLOCK_c from Sources Palettes
  • CSCOPE from Sinks Palettes
  • MUX from Signal Routing Palettes
  • BIGSOM_f from Mathematical Operations Palettes
  • CLR from Continuous time systems Palettes
  • TEXT_f from Annotations Palettes (optional)

A finished zcos file is available for download below, but I recommend that you try building one by your own first to clearly understand the process.

Connect the diagram as shown in Figure 2. The black and red signal routes cannot be mixed. Here we have red line only to the scope representing simulation clock. It is a good practice to label the blocks and signals with TEXT_f so that your diagram is understandable to others, but not required for the simulation to work.



Figure 2 feedback structure for simulation

Save the diagram to myfeedback.zcos, or whatever name you like.

Note that the CLR block has defaulted data 1/(1+s). We must put in data for our controller and plant. Figure 2 shows a CLR block after we double-click to input data



Figure 3 parameter window of CLR block

This poses another problem. The CLR block needs data in the form of numerator and denominator polynomials, but we have in Scilab workspace transfer functions C and P. MATLAB users would recall once again a function tfdata to extract the num and den information. They would be disappointed to find none. However, the data is already hidden in the transfer function C and P. For example,

-->P.num
 ans  =
    2   
 
-->P.den
 ans  =
            2  
    0.5s + s

So we simply put C.num, C.den and P.num, P.den into the controller and plant CLR blocks, respectively, as shown in Figure 4.



Figure 4 putting plant and controller data into CLR blocks

We’re almost done. The rest is just trivial task like setting simulation and scope parameters. You might not get good values at start, especially when this controller does not come from a design spec. We have no idea on its performance. Try adjusting the Y-axis range of scope to Ymin = 0, Ymax = 2. Set final integration time = 30 (and set the scope Refresh period to match). Run the simulation and see what happens then readjust them.

The resulting response is unstable. Can you figure out what goes wrong?

It turns out that we forget to set minus sign for the feedback path. Doule-click on the summing junction block to open its parameter window. Change the vector value in Inputs ports signs/gain field from [1; 1] to [1, -1]. The block now shows minus sign for the feedback.

Run the simulation again. You should now see the step response in Figure 5.



Figure 5 step response from the simulation

Bravo! Now you have some basics of Scilab and Xcos to get started on using them on a more challenged problem. To learn more, the Scilab Control Engineering Basics on our website might be helpful.

Exercise

1. Change the controller block in myfeedback.zcos to a PID block (in Continuous time systems palette). Double-click the block to adjust PID gains and run the simulation. Could you get a better response with some choice of PID gains?

Supplement

myfeedback.zcos: Xcos model file used in this article

Comments

comments

Comments

Scilab Recipe 2: with Control Engineering Flavor — 2 Comments

  1. Pingback: Scilab Recipe 1: Introduction to Scilab/Scicos (or Xcos) | Scilab Ninja

  2. Pingback: Scilab Recipe 3: Xcos Blocks Seasoning | Scilab Ninja