WYDZIAŁ ELEKTROTECHNIKI, AUTOMATYKI I INFORMATYKI

INSTYTUT AUTOMATYKI I INFORMATYKI

KIERUNEK AUTOMATYKA I ROBOTYKA

STUDIA STACJONARNE I STOPNIA

PRZEDMIOT : LABORATORIUM PODSTAW AUTOMATYKI

ĆW nr 2

TEMAT: REPREZENTACJA UKŁADÓW LTI W MATLABIE

NAZWISKO: DĄBEK IMIĘ: DOMINIKA

TERMIN WYKONANIA: 17-03-2011 TERMIN ODDANIA: : 24-03-2011

Prowadzący:

Dr inż. Grzegorz Bialic

W tej części ćwiczeń zajmowaliśmy się sposobami reprezentacji układów LTI za pomocą transmitancji; zer, biegunów, wzmocnienia układu oraz w przestrzeni stanów.

1. Transmitancja.

Transmitancją nazywamy stosunek transformaty Laplace’a sygnału wyjściowego Y(s) do transformaty Laplace’a sygnału wejściowego U(s) przy zerowych warunkach początkowych. Oznaczając ten stosunek jako G(s) mamy:

$$G\left( s \right) = \frac{Y(s)}{U(s)} = \frac{b_{1}s^{n - 1} + \ldots + b_{n - 1}s + b_{n}}{a_{1}s^{m - 1} + \ldots + a_{m - 1}s + a_{m}}$$

Aby zapisać transmitancję w Matlabie należy podać współczynniki licznika i mianownika, zaczynając od najwyższej potęgi s:

>> licz = [0 0 1]

licz =

0 0 1

>> mian = [1000 500 400]

mian =

1000 500 400

Licznik i mianownik transmitancji w pełni charakteryzują obiekt, można ich używać jako argumenty np. funkcji step (odpowiedź skokowa) i impulse (odpowiedź impulsowa):

>> step(licz,mian)

>> impulse(licz,mian)

Obiekt jest strukturą przechowującą wszystkie informacje o obiekcie.

>> obiekt = tf(licz,mian)

Transfer function:

1

----------------------

1000 s^2 + 500 s + 400

Aby wyświetlić pola tej struktury należy zastosować polecenie:

>> get(obiekt)

num: {[0 0 1]}

den: {[1e+003 500 400]}

ioDelay: 0

Variable: 's'

Ts: 0

InputDelay: 0

OutputDelay: 0

InputName: {''}

OutputName: {''}

InputGroup: [1x1 struct]

OutputGroup: [1x1 struct]

Name: ''

Notes: {}

UserData: []

W wersji sfaktoryzowanej transmitancję można wyrazić jako:

$$G\left( s \right) = k\frac{\left( s - z_{1} \right)\left( s - z_{2} \right)\ldots(s - z_{n})}{\left( s - p_{1} \right)\left( s - p_{2} \right)\ldots(s - p_{m})}$$

gdzie: z1...zn – zera, p1...pm – bieguny, k – wzmocnienie.

W celu znalezienia zer (pierwiastków licznika transmitancji), biegunów (pierwiastków mianownika transmitancji) oraz wzmocnienia układu, można zastosować funkcję tf2zp.

>> [z, p, k] = tf2zp(licz, mian)

z =

Empty matrix: 0-by-1

p =

-0.2500 + 0.5809i

-0.2500 - 0.5809i

k =

1.0000e-003

Zera i bieguny można przedstawić graficznie na płaszczyźnie zespolonej za pomocą funkcji:

>> pzmap(p,z)

albo

>> pzmap(licz,mian)

albo

>>pzmap(obiekt)

Zad. 1

1. Czy bieguny są rzeczywiste ?

Na podstawie polecenia [z, p, k] = tf2zp(licz, mian) użytego powyżej, widać, że:

p =

-0.2500 + 0.5809i

-0.2500 - 0.5809i

Czyli bieguny nie są rzeczywiste.

Czy układ jest stabilny ?

Patrząc na zamieszczony wyżej wykres można wywnioskować, że układ jest stabilny (części rzeczywiste p są ujemne).

Czy układ jest minimalnofazowy ?

Jak widać z powyższego wykresu, układ nie posiada zer w prawej półpłaszczyźnie. Układ jest zatem nieminimalnofazowy.

[z, p, k] = tf2zp(licz, mian)

z =

Empty matrix: 0-by-1

p =

-0.2500 + 0.5809i

-0.2500 - 0.5809i

k =

1.0000e-003

Postać sfaktoryzowana:

$$G\left( s \right) = \frac{0,001}{\left( s + \left( - 2500 + 0,5809i \right) \right)(s + \left( - 2500 - 0,5809i \right))}$$

1. Przypadek oscylacyjny:

>> M=1000; c=0.1; a=1;

>> epsilon=a/(2*sqrt(M*c));

>> obiekt=tf([0 0 1],[M a c]);

>> impulse(obiekt);

Przypadek tłumiony:

>> M=10000; c=1; a=200;

>> epsilon=a/(2*sqrt(M*c));

>> obiekt=tf([0 0 1],[M a c]);

>> step(obiekt);

1. Zera, bieguny, wzmocnienie.

Dany jest układ o transmitancji:

$$G\left( s \right) = \frac{3s + 1}{s\left( s + 1 \right)(s + 3)} = \frac{3(s + \frac{1}{3})}{s\left( s + 1 \right)(s + 3)}$$

Jak widać jest jedno zero z =$- \frac{1}{3}$ oraz trzy bieguny p1 = 0, p2 = –1, p3 = –3. Wzmocnienie wynosi k = 3. Licznik i mianownik transmitancji w postaci wielomianowej łatwo jest znaleźć stosując funkcję zp2tf konwersji reprezentacji zera/bieguny/wzmocnienie do transmitancji:

>> [licz,mian] = zp2tf(-1/3,[0 -1 -3],3);

Aby wyświetlić transmitancję należy wpisać:

>> printsys(licz, mian)

num/den =

3 s + 1

-----------------

s^3 + 4 s^2 + 3 s

Innym sposobem jest zastosowanie funkcji:

>> obiekt = zpk(-1/3,[0 -1 -3],3)

Zero/pole/gain:

3 (s+0.3333)

-------------

s (s+1) (s+3)

zad.2

>> obiekt = zpk (-1/4, [0 -0.1 -5], 4)

Zero/pole/gain:

4 (s+0.25)

---------------

s (s+0.1) (s+5)

1. Przestrzeń stanów.

Obiekt w tzw. „przestrzeni stanów” jest opisany za pomocą układu równań:

$\left\{ \begin{matrix} x = A\dot{x} + Bx \\ y = Cx + Du \\ \end{matrix} \right.\ $

Gdzie:

u – wejściem układu

y – wyjściem układu

x – stanem układu

$$A = \begin{bmatrix} 0 & 1 \\ - \frac{c}{M} & - \frac{\alpha}{M} \\ \end{bmatrix}\ \ \ \ \ \ \ \ \ B = \begin{bmatrix} 0 \\ \frac{1}{M} \\ \end{bmatrix}\ \ \ \ \ \ \ \ \ C = \begin{bmatrix} 1 & 0 \\ \end{bmatrix}\ \ \ \ \ \ \ \ \ D = 0$$

Aby zamienić reprezentację zera/bieguny/wzmocnienie lub transmitancję na reprezentację w przestrzeni stanów, należy zastosować funkcję:

>> [A,B,C,D] = zp2ss(z,p,k)

A =

-0.5000 -0.6325

0.6325 0

B =

1

0

C =

0 0.0016

D =

0

albo:

>> [A,B,C,D] = tf2ss(licz,mian)

A =

-4 -3 0

1 0 0

0 1 0

B =

1

0

0

C =

0 3 1

D =

0

Macierze A,B,C,D w pełni charakteryzują obiekt, można ich używać jako argumenty np. funkcji step i impulse:

>> step(A,B,C,D)

>> impulse(A,B,C,D)

Aby obliczyć wzmocnienie w stanie ustalonym można zastosować funkcję:

>> k = dcgain(A,B,C,D)

k =

Inf

Zad. 3

>> [A1,B1,C1,D1] = zp2ss(z,p,k);

>> [A1,B1,C1,D1] = zp2ss(z,p,k)

A1 =

-0.5000 -0.6325

0.6325 0

B1 =

1

0

C1 =

0 0.0016

D1 =

0

>> [A2,B2,C2,D2] = tf2ss(licz,mian);

>> [A2,B2,C2,D2] = tf2ss(licz,mian)

A2 =

-0.5000 -0.4000

1.0000 0

B2 =

1

0

C2 =

1.0e-003 *

0 1.0000

D2 =

0

Macierze A,B,C,D nie są takie same. A i C różnią się.

>> subplot(1,2,1);

>> step(A1,B1,C1,D1);

>> subplot(1,2,2);

>> step(A2,B2,C2,D2);

Odpowiedzi skokowe są takie same.

1. Dokładność obliczeń.

Funkcje MATLAB’a dają zazwyczaj wiarygodne wyniki. Należy jednak uważać w przypadku obliczeń zer i biegunów wielokrotnych oraz w przypadku układów wysokiego rzędu. Mogą pojawić się wtedy błędy numeryczne. Jako przykład weźmy układ bez zer i z 10 biegunami wielokrotnymi w s = -1, oraz ze wzmocnieniem k = 1:

>> z = [];

>> p = [-1 -1 -1 -1 -1 -1 -1 -1 -1 -1];

>> k = 1’

Przechodzimy do transmitancji:

>>[licz, mian] = zp2tf(z, p, k)

licz =

Columns 1 through 10

0 0 0 0 0 0 0 0 0 0

Column 11

1

mian =

Columns 1 through 10

1 10 45 120 210 252 210 120 45 10

Column 11

1

A następnie z powrotem do reprezentacji z/p/k :

>> [z1, p1, k1] = tf2zp(licz, mian)

z1 =

Empty matrix: 0-by-1

p1 =

-1.0461

-1.0377 + 0.0269i

-1.0377 - 0.0269i

-1.0148 + 0.0443i

-1.0148 - 0.0443i

-0.9857 + 0.0446i

-0.9857 - 0.0446i

-0.9621 + 0.0276i

-0.9621 - 0.0276i

-0.9532

k1 =

1

Jak widać występują już błędy numeryczne.

1. Wpływ zer i biegunów na kształt odpowiedzi skokowej.

Ten przykład pokazuje wpływ położenia biegunów układu II rzędu na jego odpowiedź skokową. Poniższy kod generuje trójwymiarowy wykres 12 odpowiedzi skokowych dla układu z biegunami położonymi w s = -$\frac{n}{4}$±3i gdzie n zmienia się od 1 do 12.

>> t=0:0.05:5;

>> dl=length(t);

>> LiczbaWykresow=12;

>> y=zeros(dl,LiczbaWykresow);

>> n=1;

>> while(n<=LiczbaWykresow)

[licz,mian]=zp2tf([],[-n/4+3*i -n/4-3*i], (n/4)^2+9);

[y(1:dl,n),x,tt]=step(licz,mian,t);

n=n+1;

end

>> mesh(t,1:12,y');

Wykres trójwymiarowy:

Część urojona biegunów jest stałą i wynosi ±3i, natomiast część rzeczywista zmienia się od –0.25 do –3 z krokiem 0.25. Ostatni argument funkcji zp2tf (wzmocnienie) zapewnia normalizację stanu ustalonego do 1 dla wszystkich odpowiedzi. Jak widać na rysunku, w miarę jak bieguny przesuwają się w lewo, układ staje się coraz mniej oscylacyjny (wzrasta współczynnik tłumienia) i układ staje się wolniejszy (wzrasta czas narastania).

1. Łączenie modeli.

Łączenie szeregowe

Dwa obiekty sys1 i sys2, połączone szeregowo, można połączyć w jeden obiekt sys za pomocą funkcji series:

sys = series(sys1,sys2)

Transmitancja obiektu sys jest nazywana transmitancją zastępczą dla całego układu.

Łączenie równoległe

Dwa obiekty sys1 i sys2, połączone równolegle, można połączyć w jeden obiekt sys za pomocą funkcji parallel:

sys = parallel(sys1,sys2)

Sprzężenie zwrotne

Dwa obiekty sys1 i sys2, połączone za pomocą ujemnego sprzężenia zwrotnego, można sprowadzić do jednego obiektu sys za pomocą funkcji feedback:

sys = feedback(sys1,sys2)

W przypadku dodatniego sprzężenia zwrotnego należy zastosować:

sys = feedback(sys1,sys2,+1)

zad. 6

Wprowadzenie danych transmitancji obiektów:

>> licz1 = [0 1 1];

>> mian1 = [1 5 1];

>> sys1=tf(licz1, mian1)

Transfer function:

s + 1

-------------

s^2 + 5 s + 1

>> licz2= [0 0 0 1];

>> mian2 = [ 1 1 -2 1];

>> sys2 = tf(licz2,mian2)

Transfer function:

1

-------------------

s^3 + s^2 - 2 s + 1

Obliczenie transmitancji dla połączenia szeregowego:

>> sys=series(sys1,sys2)

Transfer function:

s + 1

-------------------------------------

s^5 + 6 s^4 + 4 s^3 - 8 s^2 + 3 s + 1

Obliczenie transmitancji dla połączenia równoległego:

>> sys = parallel(sys1,sys2)

Transfer function:

s^4 + 2 s^3 + 4 s + 2

-------------------------------------

s^5 + 6 s^4 + 4 s^3 - 8 s^2 + 3 s + 1

Obliczenie transmitancji dla ujemnego sprzężenia zwrotnego:

>> sys=feedback(sys1,sys2)

Transfer function:

s^4 + 2 s^3 - s^2 - s + 1

-------------------------------------

s^5 + 6 s^4 + 4 s^3 - 8 s^2 + 4 s + 2