1

WOJSKOWA AKADEMIA T E C H N I C Z N A

im. Jarosława Dąbr owskiego

ZAKŁAD AWIONIKI I UZBROJENIA LOTNICZEGO

Przedmiot:

PODSTAWY AUTOMATYKI I AUTOMATYZACJI

ĆWICZENIE LABORATORYJNE Nr 2

MODELOWANIE UKŁADU REGULACJI AUTOMATYCZNEJ Z

WYKORZYSTANIEM PAKIETU MATLAB/SIMULINK

Warszawa 2014

1.1. C

EL ĆWICZENIA

Celem ćwiczenia jest zapoznanie studentów z zapisem modeli

dynamicznych projektowanych układów automatycznego sterowania w
środowisku MATLAB-SIMULINK oraz symulacji i analizy własności tych
układów.

1.2.

PRZEDMIOT ĆWICZENIA

Przedmiotem

ćwiczenia

jest

szczegółowa

analiza

układów

automatycznego sterowania z wykorzystaniem nowoczesnych narzędzi
modelowania i symulacji komputerowej. Studenci zapoznają się ze sposobem
zapisu i preze

ntacji układów automatyki, określą parametry statyczne

i

dynamiczne układu na podstawie otrzymanych charakterystyk czasowych

i

częstotliwościowych.

1.3.

WIADOMOŚCI OGÓLNE

1.3.1. Wprowadzenie

MATLAB

jest interakcyjnym środowiskiem do wykonywania naukowych

i in

żynierskich obliczeń oraz do wizualizacji danych. Zakres zastosowań

obejmuje różne dziedziny nauki i techniki, w tym biologię, medycynę,
ekonomię, meteorologię i wiele innych.

Zasadnicze jego zalety to: możliwość szybkiego uzyskania rezultatów

złożonych obliczeń i przedstawienia ich w postaci wykresów dwu- lub
trójwymiarowych, a także jako mapy wielobarwnej.

MATLAB jest w zasadzie językiem programowania wysokiego poziomu.

Podstawowym typem danych, na których operuje jest macierz rzeczywista lub
zespolona. W

programie tym nie deklaruje się zmiennych, są one

przechowywane w przestrzeni roboczej i są dostępne poprzez nazwę.
Najważniejszymi zaletami pakietu MATLAB są:



Otwarta architektura pakietu, na którą składają się:

­

M-pliki

, umożliwiające definiowanie poleceń i algorytmów

obliczeniowych;

­

MEX-pliki i pliki ASCII

służące do wymiany danych i wyników

obliczeń pomiędzy MATLAB’em, a innymi programami;

­

Grafika służy do wizualizacji danych i wyników obliczeń,
możliwa jest prosta animacja i efekty dźwiękowe;

­

GUI (ang. Graphics User Interface)-

interfejs graficzny dający

możliwość pracy interaktywnej za pomocą okienek
edycyjnych, przycisków, suwaków i menu.

­

Usługi DDE (ang. Dynamic Data Exchange) realizują
statyczną lub dynamiczną wymianę danych tekstowych i
graficznych

pomiędzy programami w środowisku MS

Windows. Za pomocą tej usługi można przesyłać obiekty
graficzne oraz teksty i dane liczbowe pomiędzy MATLAB’em,
a innymi programami Windows.



Toolbox-

jest to ponad 20 wyspecjalizowanych pakietów

oprogramowania, m.in. SIMULINK;

SIMULINK jest to interaktywny pakiet do modelowania i symulacji

układów dynamicznych. Umożliwia on tworzenie wielopoziomowych

3

schematów blokowych. Obiekty biblioteczne są umieszczone w okienkach jako
ikony. Program ten może współpracować z:



Rael Time Toolbox (RTW)- pakietem czasu rzeczywistego z
generatorem

kodu

C

i asemblera;



Z rzeczywistym obiektem regulacji poprzez układy sprzęgające;



SIMULINK accelerator- automatycznie kompiluje i linkuje modele
SIMULINK’a, dając znaczną poprawę szybkości symulowania
modelowanego systemu.

1.3.2.

Modele matematyczne układów dynamicznych

Układem nazywamy umownie wyodrębniony ze środowiska układ

fizyczny lub jego część. Wielkości charakteryzujące oddziaływanie środowiska
na układ nazywać będziemy wymuszeniami lub wielkościami wejściowymi.
Wymuszenia dzielimy na wielkości sterujące i wielkości zakłócające.
Wielkościami sterującymi u(t) nazywamy wymuszenia zmieniane celowo, a
wielkościami

zakłócającymi

z(t)

wielkości

podlegające

zmianom

przypadkowym (losowym). Wielkości charakteryzujące oddziaływanie układu
na środowisko nazywać będziemy odpowiedziami lub wielkościami
wyjściowymi y(t).

 

 

 

 



























t

x

t

x

t

x

t

x

n



2

1





 

t

u

1

 

t

u

2

 

t

u

n

 

t

y

n

 

t

y

2

 

t

y

1

wejście

wyjście

Rys.1. Schemat układu dynamicznego.

Przebieg procesu dynamicznego w czasie zależy nie tylko od wartości

wymusz

eń w danej chwili, ale również od wartości tych wymuszeń w

przeszłości. Można więc powiedzieć, że proces dynamiczny ma pamięć, w
której są gromadzone skutki przeszłych oddziaływań. Dla prezentacji tej
pamięci wprowadza się pojęcie stanu procesu.

Stanem procesu

nazywamy najmniejszy liczebnie zbiór wielkości x

1

(t),

x

2

(t), ..., x

n

(t)

określających w pełni skutki przeszłych oddziaływań na układ,

który jest wystarczający do przewidywania przebiegu procesu w przeszłości.
Wielkości te nazywamy zmiennymi stanu. Znajomość stanu procesu w chwili t

0

oraz wymuszeń w przedziale [t

0

,t

1

)

wystarcza do określenia przebiegów

odpowiedzi i stanu procesu w przedziale [t

0

,t

1

).

1.3.2.1

Opis układów dynamicznych w przestrzeni stanu

Układy, których miarą zmiany procesu w czasie jest pochodna wektora

stanu x(t), a stan procesu dla t>t

0

zależy tylko od stanu w chwili t

0

oraz od

sygnału sterującego u(t) w przedziale [t

0

,t]

można opisać równaniem:

 

   





t

t

u

t

x

f

t

x

,

,





(1)

4

o danym warunku początkowym:

 

0

0

x

t

x



(2)

Dla pełnego opisu układu potrzebne jest równanie wiążące wektor

odpowiedzi y(t) z wektorem stanu x(t) i wektorem wymuszenia u(t) postaci:

 

   





t

t

u

t

x

g

t

y

,

,



(3)

Dla układów liniowych niestacjonarnych (o parametrach zależnych od

czasu) równania (1) i (2) przyjmują postać:

         

t

u

t

B

t

x

t

A

t

x







(4)

 

       

t

u

t

D

t

x

t

C

t

y





(5)

gdzie: A(t)-

macierz układu o wymiarze n x n, B(t)- macierz sterowania o wymiarze n x p,

C(t) - macierz odpowiedzi o wymiarze q x n, D(t)- macierz transmisyjna o wymiarze q x p.

Schemat blokowy takiego układu przedstawia rys.2.

B(t)

D(t)

C(t)

A(t)

+

+

+

+

u(t)

x(t)

x(t)

y(t)

.



Rys.2. Schemat blokowy ciągłego układu liniowego niestacjonarnego.

Dla

układów liniowych stacjonarnych (o parametrach niezależnych od

czasu) elementy macierzy A, B, C, D

są stałe i nie zależą od czasu. W tym

przypadku równania (4) i (5) przyjmują postać:

 

 

 

t

Bu

t

Ax

t

x







(6)

 

 

 

t

Du

t

Cx

t

y





(7)

1.3.2.

Opis

układow

dynamicznych

za

pomocą

transmitancji

operatorowe i widmowej

Transmitancją operatorową G(s) jednowymiarowego układu liniowego

stacjonarnego nazywamy stosunek transformaty odpowiedzi Y(s) do
transformaty wymuszenia U(s)

tego układu przy zerowych warunkach

początkowych.

5

 

 

 

s

X

s

Y

s

G



(8)

Dla

układów

liniowych

stacjonarnych

opisanych

równaniami

różniczkowymi

liniowymi

o

stałych

współczynnikach

transmitancja

operatorowa jest ilorazem d

wóch wielomianów o postaci:

 

 

 

s

M

s

L

s

G



(9)

gdzie:

 







m

i

i

i

s

b

s

L

0

(10)

 

m

n

s

a

s

M

n

j

j

j









,

0

(11)

Pierwiastki równania L(s)=0 nazywamy zerami, a pierwiastki równania

M(s)=0

nazywamy

biegunami

transmitancji

operatorowej

opisanej

równaniem (9).

Transmitancją widmową G(j



)

liniowego układu stacjonarnego nazywać

będziemy wielkość określoną jako stosunek wartości zespolonej składowej
wymuszonej odpowiedzi Y(j



)

wywołanej wymuszeniem sinusoidalnym do

wartości zespolonej tego wymuszenia U(j



):

 

 

 







j

U

j

Y

j

G



(12)

Transmitancja widmowa jest wielkością zespoloną zależną od

parametrów układu i pulsacji wymuszenia



:

   

 







jQ

P

j

G





(13)

gdzie:

 

 









j

G

P

Re



(14)

 

 









j

G

Q

Im



(15)

Transmitancja widmowa G(j



)

z transmitancją operatorową G(s)

związana jest zależnością:

 

 





j

s

s

G

j

G





(16)

Zależność tą otrzymamy z definicji transmitancji operatorowej (8)

i

transmitancji widmowej (12) po uwzględnieniu, że przekształcenie Fouriera

6

jest szczególnym przypadkiem dla s=j



przekształcenia Laplace’a oraz, że

transformata Fouriera prze

biegu sinusoidalnie zmiennego jest równa wartości

zespolonej tego przebiegu.

1.3.3.

Kryteria jakości regulacji

Do najczęściej stosowanych kryteriów jakości regulacji należą:



Kryteria zapasu stabilności;



Kryteria rozkładu pierwiastków równania charakterystycznego;



Kryteria czasowe;



Kryteria częstotliwościowe;



Kryteria całkowe.

Miarą zapasu stabilności jest odległość charakterystyki amplitudowo-

fazowej układu otwartego od punktu (-1,j0) (rys.3.a). Zapas stabilności
amplitudy

(modułu) jest określony następującymi wielkościami:

 

 

2

2

1

1

log

20

,

log

20

k

L

k

L









(17)

gdzie k

1

i k

2

są określone w sposób podany na rys.3.

Zapas amplitudy jest tym większy, im większe są



L

1

i



L

2

. W dobrze

tłumiących układach wielkości te wynoszą od 6 do 12 dB.

Rys.3. Okr

eślenia zapasu stabilności i fazy na podstawie charakterystyki

amplitudowo-

fazowej układu otwartego.

W przypadku szczególnym, gdy charakterystyka amplitudowo-fazowa

układu otwartego ma punkty przecięcia z osią liczb rzeczywistych tylko na
prawo od punktu (-

1,j0), zapas stabilności amplitudy określamy za pomocą

jednej wielkości (rys.3.b):

k

L

log

20





(18)

Zapasem stabilności fazy



nazywamy kąt zawarty między osią liczb

rzeczywistych, a półprostą łączącą początek układu współrzędnych 0 z
punktem B

przecięcia okręgu o promieniu równym 1 i środku w początku

układu współrzędnych z charakterystyką amplitudowo-fazową układu
otwartego (rys.3.), czyli:

7













180

(19)

gdzie:



jest argumentem tr

ansmitancji widmowej układu otwartego K(j



), odp

owiadającym

modułowi równemu 1.

Zapas stabilności fazy jest tym większy, im większy jest kąt



, który

zwykle przyjmuje się w granicach od 30



do 60



. Zadając określoną wartość

zapasu stabilności amplitudy i fazy określamy obszar, przez który nie może
przechodzić charakterystyka amplitudowo-fazowa układu otwartego.

Rys.4. Określenia zapasu stabilności amplitudy i fazy na podstawie logarytmicznych

charakterystyk amplitudowej i fazowej układu otwartego

Dla lo

garytmicznych charakterystyk amplitudowych i fazowych układu

otwartego zapas stabilności zaznaczamy jak na rys.4, a definicje tych
wielkości są następujące:



Zapasem stabilności amplitudy nazywamy ujemne odchylenie
logarytmicznej charakterystyki amplitudowe

j układu otwartego od

wartości 0 dB dla pulsacji odcięcia fazy



0

, przy której

logarytmiczna charakterystyka fazowa przecina prostą -180



.



Zapasem stabilności fazy nazywamy dodatnie odchylenie
logarytmicznej charakterystyki fazowej układu otwartego od
wart

ości -180



dla pulsacji odcięcia modułu



m

, przy której

logarytmiczna charakterystyka amplitudowa przecina oś odciętych
(wartość 0 dB).



Wymagania stawiane układom dotyczące zapasu stabilności i
szybkości działania można sformułować w postaci warunków, aby
pierwiastki równania charakterystycznego układu znajdowały się w
określonej części lewej półpłaszczyzny zmiennej s. W ten sposób
budujemy kryterium położenia biegunów i zer układu,
z którym związane są pojęcia stopnia stabilności



oraz

oscylacyjności



.

8

Stopniem stabilności



nazywamy wartość bezwzględną części

rzeczywistej pierwiastka położonego w lewej półpłaszczyźnie najbliżej osi
urojonej (rys.5). Natomiast

oscylacyjnością



nazywamy największą wartość

bezwzględną stosunku części urojonej



k

do części rzeczywistej



k

pierwiastka

zespolonego

s

k

=-



k



j



k

spośród wszystkich pierwiastków równania

charakterystycznego układu:

k

k

s

k

k

s

k

k

s

s







max

Re

Im

max





(20)

Rys.5. Określenia stopnia stabilności.

Założenie dla układu określonej oscylacyjności



jest równoważne

wymaganiu, aby pierwiastki równania charakterystycznego układu leżały w
części lewej półpłaszczyzny zmiennej zespolonej s ograniczonej z prawej
strony dwiema półprostymi nachylonymi względem osi rzeczywistej pod kątem
nie mniejszym

niż



=arctg(



).

Kryteria czasowe określają pewne parametry lub właściwości przebiegów

wielkości regulowanych w stanie nieustalonym wywołanych skokową zmianą
wielkości zadających lub zakłóceń. Skłonność układu do oscylacji można w
przybliżeniu ocenić za pomocą wielkości zwanej przeregulowaniem, którą
definiujemy jako maksymalną różnicę między wartością chwilową wielkości
regulowanej w stanie nieustalonym, a jej wartością ustaloną, odniesioną do tej
wartości ustalonej:

0

,

max







u

u

u

y

y

y

y



(21)

9

gdzie: y

max

jest maksymalną wartością chwilową wielkości regulowanej, a y

u

jest wartością

ustaloną tej wielkości. (rys.6.). Przyjmuje się, że dopuszczalne jest przeregulowanie w
granicach od 0,1 do 0,3.

Kolejnym parametrem określającym szybkość działania układu jest czas

regulacji t

r

zdefiniowany jako czas (mierzony od chwilowej skokowej zmiany

wymuszenia), po upływie którego wielkość regulowana różni się od wartości
ustalonej trwale mniej niż założona wartość



, przyjmowana zwykle od 2 do

5% w

artości ustalonej.

Między czasem regulacji t

r

, a stopniem stabilności



zachodzi zależność:





1

ln

1



r

t

(22)

Z kryterium częstotliwościowym oceny jakości układu związane są

następujące pojęcia: szczyt rezonansowy M

p,

pulsacja rezonansowa



r

oraz

polsacja odcięcia



o

.

Szczytem

rezonansowym

M

p

nazywamy największą wartość

charakterystyki amplitudowej układu zamkniętego. Im większy jest zapas
stabilności i im słabsze jest tłumienie oscylacji, tym większy jest szczyt
rezonansowy M

p

. Dla układów dobrze tłumiących oscylacje wartość tego

parametru powinna być zawarta w granicach od 1,1 do 1,5.

Pulsacją rezonansową



r

nazywamy pulsację odpowiadającą

największej wartości charakterystyki amplitudowej układu zamkniętego.
Korzystając z charakterystyki amplitudowo-fazowej układu otwartego możemy
w przybliżeniu określić pulsację rezonansową



r

jako pulsację odpowiadającą

punktowi tej charakterystyki położonemu najbliżej punktu (-1,j0). Pulsacja
rezonansowa niewiele różni się od pulsacji drgań występujących w układzie w
stanie nieustalonym. Czas regulacji t

r

z pulsacją rezonansową związany jest

zależnością:

r

r

t





3



(23)

Pulsacją odcięcia



o

nazywać będziemy pulsację odpowiadającą

punktowi przecięcia charakterystyki amplitudowo-fazowej układu otwartego
K(j



)

z okręgiem o promieniu jednostkowym i środku w początku układu

współrzędnych, |K(j



)|=1

. Pulsację odcięcia można określić również jako

pulsację, przy której logarytmiczna charakterystyka amplitudowa układu
otwartego przecina oś odciętych, L(



o

)=0.

Ostatnim kryterium określającym jakość układu regulacji jest kryterium

całkowe. Pozwala ono ocenić jakość regulacji w stanie ustalonym (dokładność
statyczna), jak i w stanie nieustalonym (zapas stabilności i szybkość działania
układu). Za całkowe kryteria jakości regulacji przyjmuje się funkcjonały typu:









k

t

dt

e

e

e

t

f

I

0

,...

,

,

,







(24)

gdzie f jest funkcja, e=e(t)- uchybem regulacji, a t

k

=



.

10

Najczęściej za całkowe kryteria jakości regulacji przyjmuje się

następujące przypadki szczególne funkcjonału (24):





k

t

dt

e

I

0

2

20

(25)

)

(

0

2

2

naturaIna

liczba

dana

p

dt

e

t

I

k

t

p

p







(26)





)

(

0

2

2

2

2

ik

wspólczynn

dany

T

dt

e

T

e

I

k

t











(27)





k

t

dt

e

I

0

0

(28)

)

(

0

2

naturaIna

liczba

dana

p

dt

e

t

I

k

t

p

p







(29)

)

(

0

2

naturaIna

liczba

dana

p

dt

e

t

I

k

t

p

p







(30)

1.4.

ANALIZA UKŁADÓW DYNAMICZNYCH W ŚRODOWISKU MATLAB-
SIMULINK

1.4.1. Postaci modeli dynamicznych wykorzystywanych w pakiecie

Control Toolbox

środowiska MATLAB

Rys.6. Okno programu MATLAB

W pakiecie tym wykorzystywane są następujące postaci liniowych modeli

dynamicznych:



Równania stanu;



Macierze transmitancji dla układów SIMO (jedno wejście wiele
wyjść);

11



Macierze transmitancji (układy SIMO) w postaci iloczynu zer,
biegunów i wzmocnienia.

Transmitancja układu dynamicznego jest funkcją wymierną operatora s

(dla układów ciągłych) lub z (dla układów dyskretnych). W bibliotece Control
Toolbox

transmitancję podaje się w postaci pary wektorów zawierających

współczynników licznika i mianownika, przy czym umieszcza się je tam wg.
malejących potęg operatorów s lub z, np.: L=[1 2]; M=[1 3 2];

Takiemu zapisowi odpowiada transmitancja:

 

2

3

2

2









s

s

s

s

G

. Aby okre

ślić

transmitancję układu SIMO należy:



Podać wektor współczynników mianownika transmitancji;



Podać macierz zawierającą w kolejnych wierszach współczynniki
liczników odpowiadających kolejnym wyjściom układu.

Na przykład, zapis L=[1 2; 3 1]; M=[1 3 2] odpowiada transmitancji:

 



































2

3

1

3

2

3

2

2

2

s

s

s

s

s

s

s

G

.

W programie MATLAB posługujemy się wyłącznie transmitancjami

spełniającymi warunek mówiący, ze stopień mianownika powinien być większy
od stopnia licznika. Aby określić macierz transmitancji przez podanie zer,
biegunów i wzmocnienia należy podać:



Kolumnowy wektor biegunów transmitancji P (ang. poles);



Macierz Z (ang. zeros) zawierającą w kolejnych kolumnach zera
odpowiadające kolejnym wyjściom układu;



Podać wektor kolumnowy G (ang. gains) zawierający wzmocnienia
odpowiadające kolejnym wyjściom układu.

Opisując układ dynamiczny w programie MATLAB należy pamiętać, że
najkorzystniejsze własności ma przedstawienie układu w postaci równań
stanu. Wykorzystujące podczas obliczeń postać transmitancji operaatorowej
układu program staje się bardziej podatny na błędy. Dotyczy to zwłaszcza
transmitancji wysokich rzędów lub o bardzo dużych lub bardzo małych
współczynnikach.
W programie tym przy konstruowaniu modeli dynamicznych wykorzystuje się
następujące polecenia:



Zmiana postaci modelu:

Skrót

Zapis

Opis

ss2tf

[N, D]=ss2tf(A,B,C,D,iu)

Funkcja ta zamienia równania stanu na odpowiadającą mu
transmitancję liczoną względem wejścia o numerze iu.
Macierz N

zawiera w kolejnych wierszach współczynniki

liczników transmitancji odpowiadające kolejnym wyjściom
układu. Wektor D zawiera współczynniki mianownika
transmitancji.

tf2ss

[A, B, C, D]=tf2ss(N, D) Funkcja ta wykonuje operacje odwrotne do funkcji ss2tf.

ss2zp

[Z,P,K]=ss2zp(A,B,C,D,i
u)

Funkcja ta zamienia równania stanu na odpowiadającą im
macierz transmitancji liczoną względem wejścia o numerze
iu

i przedstawia ją w postaci zer Z, biegunów P i

wzmocnień K.

zp2ss

[A, B, C, D]=zp2ss(Z, P,
K)

Funkcja ta wykonuje operacje odwrotne do funkcji ss2zp.

tf2zp

[Z,P,K]=tf2zp(N, D)

Funkcja ta znajduje zera Z, bieguny P i wzmocnienia K dla

12

transmitancji lub macierzy transmitancji układu SIMO
określonej parametrami N- współczynniki licznika, D-
współczynniki mianownika.

zp2tf

[N, D]=zp2tf(Z, P, K)

Funkcja ta wykonuje operacje odwrotne do funkcji zp2tf.



Zmiana struktury równań:

Skrót Zapis

Opis

appen
d

[A,B,C,D]=
append(A1,B1,C1 D1,A2,B2, C2,
D2)

Funkcja ta dołącza do układu opisanego
równaniami stanu określonymi macierzami A1, ...
układ o równaniach stanu wyznaczonymi przez
macierze A2,... .

augst
ate

[Ab, Bb, Cb, Db]=augstate(A,B,C,D) Funkcja ta wykonuje operacje odwrotne do funkcji

ss2tf.

ss2zp [Z,P,K]=ss2zp(A,B,C,D,iu)

Funkcja ta rozszerza wektor wyjść układu o
współrzędne stanu. Zmiana ta wpływa wyłącznie
na równanie wyjścia układu.

ssdel
ete

[A,

B,

C,

D]=

ssdelete(A1,B1,C1,D1,WEJ,WYJ,S
TANY)

Funkcja ta usuwa zadane wejścia WEJ, wyjścia
WYJ i stany STANY

z równań stanu o macierzach

A1,... .

sssele
ct

[A,

B,

C,

D]=

ssdelete(A1,B1,C1,D1,WEJ,WYJ,S
TANY)

Funkcja ta umożliwia utworzenie podukładu na
podstawie układu zadanego przez macierze
A1,....Parametry WEJ, WYJ, STANY

są wektorami

zawierającymi indeksy wejść, itd., które mają
zostać włączone do tworzenia podukładu.



Tworzenie modeli dla różnych struktur układu:

Skrót

Zapis

Opis

cloop

[A,B,C,D]=cloop(A1,B1,C1,D1,ZNAK)
[A,B,C,D]=cloop(A1,B1,C1,D1,WEJ,WYJ)
[L,M]=cloop(L1,M1, ZNAK)

Funkcja ta oblicza macierze stanu układu
zamkniętego sprzężeniem zwrotnym o
wzmocnieniu

1.

ZNAK

określa znak

sprzężenia zwrotnego

feedba
ck

[A,B,C,D]=
feedback(A1,B1,C1,D1,A2,B2,C2,D2,ZN
AK)
[A,B,C,D]=
feedback(A1,B1,C1,D1,A2,B2,C2,D2,WE
J,WYJ)
[L,M]= feedback (L1,M1,L2,M2,ZNAK)

Funkcja ta oblicza macierze stanu układu
zamkniętego

sprzężeniem

zwrotnym

spisanym macierzami A2,.. i macierzami w
torze głównym A1,... ZNAK określa znak
sprzężenia zwrotnego

paralle
l

[A,B,C,D]=parallel(A1,B1,C1,D1,A2,B2,C
2,D2)
[A,B,C,D]=
feedback(A1,B1,C1,D1,A2,B2,C2,D2,
WEJ1,WYJ1,WEJ2,WYJ2)
[L,M]= parallel(L1,M1,L2,M2)

Funkcja ta oblicza macierze stanu dla
obiektów połączonych równolegle.

series

[A,B,C,D]=

series

(A1,B1,C1,D1,A2,B2,C2,D2)
[A,B,C,D]=

series

(A1,B1,C1,D1,A2,B2,C2,D2,
WYJ1,WEJ2)
[L,M]= series (L1,M1,L2,M2)

Funkcja ta oblicza macierze stanu dla
obiektów połączonych szeregowo.

13



Inne funkcje wspomagające konstruowanie modeli dynamicznych:

Skrót

Zapis

Opis

rmode
l
dmod
el

[L,M]=rmodel(n)
[L,M]=rmodel(n,p)
[A,B,C,D]=rmodel(n)
[A,B,C,D]=rmodel(n,p,
m)

Funkcje te generują losowe modele liniowe i dyskretne.
n-

rząd modelu, p- ilość wyjść modelu, m-ilość wejść modelu

Ord2

[A,B,C,D]=ord2(wn,dze
ta)
[L,M]=ord2(wn,dzeta)

Funkcja ta generuje odel układu oscylacyjnego II rzędu.
wn-



n

-

częstotliwośc rezonansowa;

dzeta-



-

współczynnik tłumienia.

Pade

[A,B,C,D]=pade(T,n)
[L,M]= pade(T,n)

Funkcja ta zwraca macierze równań stanu lub transmitancję
układu przybliżającego opóźnienie o dany czas T.
Przybliżenie to jest realizowane w oparciu o rozwinięci funkcji
w szereg Taylora.

1.4.2.

Charakterystyki czasowe i częstotliwościowe układów

W programie MATLAB możemy wygenerować odpowiedź układu w dziedzinie
czasu na następujące wymuszenia:



Odpowiedź na skok jednostkowy;



Odpowiedź impulsowa;



Odpowiedź układu na stan początkowy- przebiegi stanu i wyjścia
jakie powstaną w układzie przy braku pobudzenia z zewnątrz i
danych warunkach początkowych;



W dz

iedzinie częstotliwości możemy wykreślić następujące

charakterystyki układu dynamicznego:



Charakterystykę Bode;



Charakterystykę Nyquista;

Powyższe zadania realizują następujące funkcje:

Skrót

Zapis

Opis

impuls
e
dimpul
se

impulse(A,B,C,D,iu,
t)
impulse(L,M,t)

Funkcje te wyznaczają odpowiedź układu ciągłego (dimpulse-
dyskretnego) na impuls jednostkowy. W przypadku układu
opisanego równaniami stanu kreślona jest odpowiedź wszystkich
wyjść na impuls pojawiający się na wejściu o numerze iu. Dla
układów ciągłych można dołączyć własny wektor czasu symulacji
odpowiedzi. W przypadku pominięcia parametru t wektor czasu
zostanie pominięty.

initial
dinitial

initial(A,B,C,D,X0,T
)

Funkcja ta wyznacza odpowiedź układu opisanego równaniami
stanu na warunek początkowy określony parametrem X0. Czas
symulacji i wektor chwil czasu jest określany automatycznie lub
przez użytkownika (parametr T).

step
dstep

step(A,B,C,D,iu)
step(A,B,C,D,iu,t)
step(L,M)

Funkcja ta wyznacza odpowiedź układu ciągłego na skok
jednostkowy. W przyp

adku układu opisanego równaniami stanu

kreślona jest odpowiedź wszystkich wyjść na skok jednostkowy
pojawiający się na wejściu o numerze iu. Dla układów ciągłych
można dołączyć własny wektor czasu symulacji odpowiedzi. W
przypadku pominięcia parametru t wektor czasu zostanie
pominięty.

bode

bode(A,B,C,D,iu,w)
bode(L,M,w)

Funkcja ta wyznacza charakterystykę częstotliwościową Bode.
Częstotliwości i liczba punktów dobierana jest automatycznie lub
zgodnie z wymaganiami użytkownika (parametr w). W przypadku
uk

ładu opisanego równaniami stanu kreślona jest odpowiedź

wszystkich wyjść na skok jednostkowy pojawiający się na wejściu
o numerze iu.

14

nyquist nyquist(A,B,C,D,iu,

w)
nyquist(L,M,w)

Funkcja ta wyznacza charakterystykę Nyquista. Częstotliwości i
liczba

punktów dla których kreślona jest charakterystyka dobierane

są automatycznie lub zgodnie z wymaganiami użytkownika
(parametr w). W przypadku układu opisanego równaniami stanu
kreślona jest odpowiedź wszystkich wyjść na skok jednostkowy
pojawiający się na wejściu o numerze iu.

1.4.3.

Analiza właściwości dynamicznych obiektu sterowania.

Analiza właściwości obiektu sterowania polega na określeniu stabilności

układu, czy jest on obserwowalny, czy też nie, obliczeniu zapasu amplitudy i
fazy obiektu jak również wykreśleniu dróg pierwiastków układu. Zadania te
realizują następujące funkcje:

Skrót

Zapis

Opis

margin

[Gm,

Pm,

Wcg,

Wcp]

=

margin(A,B,C,D)
[Gm, Pm, Wcg, Wcp] = margin(L,M)

Funkcja ta oblicza zapas amplitudy (Gm) i fazy
(Pm) dla układu opisanego równaniami stanu
lub transmitancją oraz odpowiadające im
częstotliwosci graniczne- odpowiednio Wcg i
Wcp.

imargi
n

[Gm,

Pm,

Wcg,

Wcp]

=

margin(ampl,faza,w)

Funkcja ta działa podobnie jak funkcja margin z
tym, że obliczenia opierają się o zadane
wektory amplitudy i f

azy oraz odpowiadające im

częstotliwości, które można uzyskać funkcją
bode lub dbode, dzięki czemu funkcja ta działa
zarówno

dla

układów

ciągłych,

jak

i

dyskretnych.

rlocus

rlocus(A,B,C,D)
rlocus(L,M)
rlocus (L,M,K)
rlocus(A,B,C,D,K)

Funkcja ta wyznacza d

rógi pierwiastków układu

otwartego.

rlocfin
d

[K,R]=rlocfind(A,B,C,D)
[K,R]=rlocfind (L,M)
[K,R]=rlocfind (L,M,n)
[K,R]=rlocfind (A,B,C,D,n)

Jeżeli w aktywnym oknie znajduje się wykres
zależności położenia biegunów od wzmocnienia
uzyskany funkcją rlocus, to funkcja rlocfind
umożliwia wybranie myszką żądanego bieguna
układu SISO opisanego transmitancją lub
równaniami stanu.

1.5. SIMULINK

Wygodnym i skutecznym narzędziem symulacyjnym oferowanym przez

środowisko MATLAB jest program SIMULINK. Umożliwia on utworzenie w
oknie graficznym struktury układu sterowania zbudowanej z bloków różnych
typów reprezentujących obiekty dynamiczne, źródła sygnału i przyrządy
pomiarowe. Definiując obiekty możemy odwołać się do istniejących w pamięci
zmiennych, dostępnych normalnie z wiersza poleceń MATLAB’a. Po
skonstruowaniu struktury i określeniu parametrów poszczególnych elementów
ustalamy parametry symulacji i wybieramy polecenie Simulation/Start. Wyniki
symulacji zostaną przedstawione w standardowych oknach graficznych
programu MATLAB.

15

Rys.7. Okno programu SIMULINK

W programie tym mamy do dyspozycji nastepujące grupy elementów
dynamicznych (rys.7.):



Źródła (ang. sources)- np. skok jednostkowy, generator sygnałów
zmiennch;



Elementy wizualizacji wyników symulacji (ang. sinks)- np.
oscyloskop, itd.;



Elementy ciągłe (ang. continuous);



Elementy dyskretne (ang. discrete);



Elementy realizujące określone zadania matematyczne (ang.
math);



Funkcje i tabele (ang. functions & tables);



Elementy nieliniowe (ang. nonlinear);



Inne elementy ni

ezbędne do zbudowania układu symulacji. Blok

Signals & Systems-

np. multipleksery i demultipleksery, wejścia i

wyjścia układu, itd.

Rozwinięcie bloków elementów z rys.7. przedstawiono na rys.8. i rys.9.

Rys.8. Biblioteki: Sources, Signal & Systems, Math.

16

Program ten ściśle współpracuje z środowiskiem MATLAB. Dlatego też,

pozwala on na analizę układu zamodelowanego w programie SIMULINK za
pomocą funkcji linmod oraz dlinmod o nastepującym zapisie: [A, B, C,
D]=linmod(‘nazwa’). W miejsce „nazwa” należy wpisać nazwę pliku
zawierającego odpowiednio przygotowany model analizowanego układu.
Właściwe przygotowanie modelu polega na wskazaniu w programie
SIMULINK wejścia i wyjścia układu (rys.9.).

Rys.9. Biblioteki progarmu SIMULINK: Continuous, Discrete, Nonlinear,

Functions & Tables, Sinks

17

1.6.

PRZEBIEG ĆWICZENIA LABORATORYJNEGO

1.6.1. Opis stanowiska laboratoryjnego

Stanowisko laboratoryjne składa się komputera klasy PC wraz z

oprogramowaniem MATLAB- SIMULINK.

1.6.2.

Wykonanie ćwiczenia

1) Badasz układ opisany transmitancją operatorową postaci:

 

1

2

2

2







Ts

s

T

k

s

G



dla danych współczynników:

a)

k=1.3, T=1, ξ =1.9

b)

k=1.2, T=1, ξ =2

c)

k=1.4, T=1, ξ =1.8

d)

k=1.1, T=1, ξ =1.7

e)

k=1.4, T=1, ξ =2.1

f)

k=1.2, T=1, ξ =1.9

g)

k=1.3, T=1, ξ =1.7

h)

k=1.5, T=1, ξ =1.9

2) Korzystając z danych współczynników obliczyć na kartce postać

transmitancji G(s)

3) W programie Matlab zdefiniować transmitancję badanego układu:

a) zdefiniować licznik np.: dla równania

5

3

2



s

zapisać:

L = [3 0 5]

b) zdefiniować mianownik np.: dla równania

7

4

2

2

3





s

s

zapisać:

M = [2 4 0 7]

c) funkcją tf(licznik, mianownik) definiujemy transmitancję np.:

dla powyższych danych:

G = tf (L, M)

4) Wyznaczamy charakterystykę skokową układu opisanego

transmitancją G. Zapisujemy: step(G).
a) Zapisujemy: step(G).
b) Na wykresie klikamy prawym klawiszem myszy →

characteristics i wybieramy wartości charakterystyczne, które
chcemy wyświetlić. np.:
i) peak response – wartość maksymalna
ii) settling time – czas regulacji
iii) rise time – czas narastania
iv) steady state – stan ustalony

c) komentujemy odczytane wartości

5) Wyznaczamy charakterystykę impulsową układu opisanego

transmitancją G.
a) Zapisujemy: impulse(G).

18

b) Na wykresie klikamy prawym klawiszem myszy →

characteristics i wybieramy wartości charakterystyczne, które
chcemy wyświetlić. np.:
i) peak response – wartość maksymalna
ii) settling time – czas regulacji

c) komentujemy odczytane wartości

6) Wyznaczamy charakterystyki amplitudową i fazową układu

opisanego transmitancją G (wykresy Bodego).
a) Zapisujemy: bode(G).
b) Na wykresie klikamy prawym klawiszem myszy →

characteristics i wybieramy wartości charakterystyczne, które
chcemy wyświetlić. np.:
i) peak response – wartość maksymalna
ii) All stability margins – zapas amplitudy i fazy

c) komentujemy odczytane wartości

7) Wyznaczamy charakterystykę amplitudowo fazową układu

opisanego transmitancją G (wykres Nyquista).
a) Zapisujemy: nyquist(G).
b) Na wykresie klikamy prawym klawiszem myszy →

characteristics i wybieramy wartości charakterystyczne, które
chcemy wyświetlić. np.:
i) peak response – wartość maksymalna
ii) All stability margins – zapas amplitudy i fazy

c) komentujemy odczytane wartości

8) Zamykamy układ pętlą ujemnego sprzężenia zwrotnego ze

wzmocnieniem:

wykorzystując funkcję:
[Lf, Mf] = feedback(L, M, L1, M1, znak sprzężenia),
gdzie:
- L – licznik układu w torze głównym
- M – mianownik układu w torze głównym
- L1 – licznik układu w torze sprzężenia zwrotnego
- M1 – mianownik układu w torze sprzężenia zwrotnego
- Lf – licznik układu zamkniętego
- Mf – mianownik układu zamkniętego
np.:
dla wzmocnienia k = 3 i ujemnego sprzężenia zwrotnego:

19

[Lf, Mf] = feedback(L, M, [3], [1], -1)

9) Zdefiniować transmitancję układu zamkniętego zapisując:
10) F = tf (Lf, Mf)
11) Wyznaczyć charakterystyki czasowe układu zamkniętego F wg

punktów 4 i 5

12) Porównać wartości charakterystyczne dla układów otwartych i

zamkniętych.

13) Dobrać wzmocnienie w torze sprzężenia zwrotnego w celu

minimalizacji wybranej wartości charakterystycznej układu
zamkniętego np. czasu regulacji.

14) Ćwiczenie wykonać dla różnych współczynników transmitancji

układu w torze głównym wymienionych w pkt. 1 a) - h).

1.7. WYKONANIE SPRAWOZDANIA

W sprawozdaniu powinno zostać zawarte:



Model badanego układu dynamicznego;



Wydruk M-pliku;



Wydruk charakterystyk dynamicznych układu;



Wnioski do otrzymanych wyników.

1.8. LITERATURA:

1.

B. Mrozek, Z. Mrozek, „MATLAB. Uniwersalne środowisko do obliczeń
naukowo-

technicznych”, Kraków 1995.

2.

T. Kaczorek „Teoria sterowania”, Warszawa 1977

3.

A. Zalewski, R. Cegieła „Matlab- obliczenia numeryczne i ich
zastosowania”, Poznań 1996.