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.