Modelowanie układów dynamiki w pakiecie MATLAB
W pakiecie Matlab modelowanie układów dynamicznych można przeprowadzić na wiele
sposobów. Podstawowymi metodami są:
" Modelowanie z wykorzystaniem równań różniczkowych opisujących badany układ
(np. równań dynamicznych ruchu)
" Modelowanie z wykorzystaniem transmitancji.
" Modelowanie w przestrzeni stanu.
Wszystkie rodzaje modeli są sobie równoważne. Dodatkowo, bazując na metodach podanych
powyżej, możliwe jest modelowanie dynamiki układu, wykorzystując graficzny interfejs
programu MATLAB, tj. SIMULINK a.
1. Podstawy teoretyczne
Układ o jednym stopniu swobody.
Dla układu jak na rysunku poniżej, równania dynamiczne ruchu zapisane zostać mogą w
następującej postaci:
(1)
gdzie:
" M  masa
" C  tłumienie
" K  sztywność
" - przyspieszenie, prędkość przemieszczenie
" f  siła wymuszająca
" t  czas.
K
f(t)
M
x(t)
C
Rys. 1 Przykład układu o jednym stopniu swobody
Wykorzystując transformatę Laplace a, można zalgebraizować powyższe równanie zapisując
je w formie:
(2)
gdzie s jest operatorem Laplace a.
Równanie to można zapisać wykorzystując pojecie sztywności dynamicznej Z:
(3)
Transmitancję układu wyznaczyć możemy korzystając z zależności:
(4)
czyli
(5)
Podstawiając powyższą zależność do równania opisującego dynamikę układu w dziedzinie
zespolonej(2):
(6)
Powyższa zależność definiuje transmitancję układu.
Licznik transmitancji określa się mianem równania charakterystycznego, a jego rozwiązanie
daje informacje o biegunach układu:
(7)
Korzystając z tej zależności, określić można następujące parametry układu:
Nietłumioną częstość drgań własnych (jeśli w układzie nie ma tłumienia):
(8)
Tłumienie krytyczne (minimalna wartość tłumienia przy której ruch jest periodyczny):
(9)
Stosunek tłumienia:
(10)
Rozwiązaniem równanie w dziedzinie czasu jest zależność:
(11)
W zależności od wartość tłumienia układ można sklasyfikować jako układ z tłumieniem
podkrytycznym (ś<1), krytycznym (ś=1) lub nadkrytycznym (ś>1).
Z punktu widzenia identyfikacji interesującym jest przypadek gdy tłumienie jest
podkrytyczne. Dla tego przypadku rozwiązanie równania charakterystycznego składa się z
dwóch sprzężonych pierwiastków:
(12)
gdzie �1 jest współczynnikiem tłumienia, a �1 jest tłumiona częstością drgań własnych.
Inne zależności możliwe do wyznaczenia to:
Układ o dwóch stopniach swobody:
Rys. 2. Układ o dwóch stopniach swobody
Dla układu jak powyżej, równania ruchu można zapisać jako:
(13)
Porządkując współczynniki przy odpowiednich zmiennych, można doprowadzić do zapisu
macierzowego w postaci:
lub inaczej
(14)
gdzie:
" [M]jest macierzą mas
" [C]  jest macierzą tłumień
" [K]  jest macierzą sztywności
" {f(t)}  jest wektorem sił
" {x(t)}  jest wektorem odpowiedzi
Zapisując transformatę Lplace a (zakładając zerowe warunki początkowe):
(15)
Otrzymać możemy macierz sztywności dynamicznej, wyrażoną jako
(16)
Odwracając powyższą zależność otrzymamy macierz transmitancji:
(17)
2. Przykłady.
Sposoby modelowania.
Simulink  graficznie:
Symulacja układu masowo-sprężysto-tłumiącego
Rozpatrzmy układ masowo-sprężysto-tłumiący przedstawiony na rysunku.
Model matematyczny takiego układu opisany jest równaniem:
(2)
gdzie m jest masą układu, c jest współczynnikiem tłumienia oraz k jest
współczynnikiem sztywności. Dodatkowo f(t) jest funkcją wymuszającą układ
w kierunku osi x założonego układu współrzędnych.
W przykładzie pokazane zostanie wykorzystanie Simulinka do symulacji
Rys. E2-1
odpowiedzi układu na wymuszenie sygnałem skoku jednostkowego.
Należy uruchomić Simulinka, utworzyć nowe okno modelu (CTRL+N) oraz za pomocą
myszy przeciągnąć poniższe bloki z okna biblioteki Simulink library:
Bloczki do przeciągnięcia do Lokalizacja bloczka w bibliotece
nowego okna modelu
Step Sources
Gain Math Operation
Krok 1
Sum Math Operation
Integrator Continuous
Scope Sinks
To Workspace Sinks
Równanie (2) opisujące dynamikę układu należy rozwikłać ze względu na
przyspieszenie i przedstawić w postaci opisanej równaniem (3).
(3)
Bazując na równaniu (3), można zestawić układ blokowy pokazany na rys. E2-2.
Dla przypomnienia: skróty klawiszowe CTRL+F oraz CTRL+R pozwalają na odwracanie
(180�) oraz obracanie zaznaczonego bloczka. Można także użyć przycisku CTRL+prawy
klawisz myszy do wykonania rozdzielenia linii łączącej bloczki. Na tym etapie nie
należy przejmować się ustawieniami parametrów poszczególnych bloczków.
Ustawienia dokonane zostaną w następnym kroku. Na tym poziomie konieczne jest
Krok 2
jedynie poprawne i przejrzyste zdefiniowanie układu i wykonanie połączeń pomiędzy
bloczkami.
W tym kroku należy zadać odpowiednie parametry dla każdego z bloków.
W tym przykładzie należy ustawić: m=2.0; c=0.7; k=1. W dalszej części
pracy z przykładem możliwe będzie sprawdzenie jak zmiany parametrów
modelu wpływają na otrzymywane charakterystyki układ. Zmiany można
będzie dokonać identycznie jak w tym momencie.
Sygnały będące wynikiem symulacji mogą być nie tylko prezentowane w
formie wykresów ale możliwe jest także ich dalsze wykorzystanie. W
przykładzie pokazane jest użycie dodatkowego bloczka wyjściowego
 simout . Element może zostać znaleziony w zbiorze Sinks group
biblioteki Simulink Library browser. Wyjście tego elementu
Krok 3
przekazywane jest bezpośrednio do przestrzeni roboczej Matlaba. Dla
pokazania sposobu użycia elementu, należy nazwać parametr wyjściowy
bloczka (przykładowo nadajmy zmiennej wyjściowej taką samą nazwą jak
nazwa bloczka  simout ). Zmianę tę można wykonać w okienku
edycyjnym elementu. Dostęp do niego jest uzyskiwany poprzez
dwukrotne kliknięcie na elemencie w oknie modelu.
Dodatkowo ustawmy typ wynikowej zmiennej jako tablicę wartości
czasowych uzyskiwanych w trakcie symulacji. Dokonamy tego w okienku
Simulation Parameter window (może ono zostać otwarte sekwencją
klawiszy (CTRL+E)) wybierając zakładkę Workspace I/O. Operację
pokazano na rys. E2-3.
Rys. E2-3
W kolejnym kroku można już uruchomić symulację poprzez kliknięcie
przycisku (alternatywnie można użyć skrót klawiaturowy CTRL+T ). W
wyniku symulacji uzyskuje się przebieg zaprezentowany na rys. E2-3.
Rys. E2-3
Krok 4
Uzyskanie przebiegu odpowiedzi układ jest ostatnim elementem
symulacji. Otrzymana odpowiedz
obrazuje zachowanie układu z
tłumieniem podkrytycznym. W celu przetestowania zachowania obiektu
w innych warunkach, konieczne jest dokonanie zmian parametrów masy,
tłumienia i sprężystości (wartość parametrów wpisanych w bloczkach m,
c lub k).
Zobaczmy teraz jak można wykorzystać dodatkowy blok  simout użyty w
trakcie symulacji. W oknie poleceń Matlaba wpiszmy polecenie.
>> who
Powoduje ono wypisanie wszystkich zmiennych zdefiniowanych w
przestrzeni roboczej Matlaba. W wyniku działania polecenia na ekranie
powinny zostać wypisane m.in. zmienne "simout" oraz "time". Zostały
one utworzone przez Simulinka w trakcie wykonywania przez niego
symulacji.
Dane zawarte w tych zmienny mogą zostać użyte w klasyczny sposób
zgodnie z możliwościami jakie daje środowisko Matlab. Przykładowo
wykorzystajmy zawartość uzyskanych zmiennych do pokazania
charakterystyki czasowej używając bezpośrednio poleceń środowiska
Matlab.
>>plot(time,simout);grid
Symulacja układu z przykładu 2 lecz zamodelowanego w postaci funkcji przejścia
W przypadku modelowania układu za pomocą funkcji przejścia pierwszym konieczym do
przeprowadzenia zadaniem jest przekształcenie równań dynamicznych do postaci funkcji
przejścia (nazywanej także funkcją transmitancji). Aby tego dokonać równanie ruchu układu
zapisane w dziedzinie czasu należy przetransformować za pomocą przekształcenia Laplace'a
do dziedziny transformat. W naszym przypadku równanie w dziedzinie czasu ma postać:
Zakładając zerowe warunki początkowe możemy dokonać transformacji:
2
f (0) = a = 0
f(0) = b = 0
Funkcja przejścia przybiera w tym przypadku formę:
Mając równanie układu w postaci funkcji transmitancji można przystąpić do pracy w
Simulinku.
Krok 1 Należy uruchomić Simulinka, utworzyć nowe okno modelu (CTRL+N) oraz za pomocą
myszy przeciągnąć poniższe bloki z okna biblioteki Simulink library:
Należy uruchomić Simulinka, utworzyć nowe okno modelu (CTRL+N) oraz za pomocą
myszy przeciągnąć poniższe bloki z okna biblioteki Simulink library:
Bloczki do przeciągnięcia do Lokalizacja bloczka w bibliotece
nowego okna modelu
Step Sources
Transfer Function Continuous
Scope Sinks
Save File to Workspace Sinks
Krok 2 Bloczki należy poukładać w oknie modelu np. w sposób przedstawiony na poniższym
rysunku E3-1.
Rys. E3-1
Uwaga: Kolory bloczków można ustawić poprzez kliknięcie prawym przyciskiem myszy
na bloczku oraz wybraniu koloru z menu Background Color.
Krok 3 Wprowadzanie wartości parametrów bloczków.
W przykładzie użyte zostaną te same wartości parametrów modelu m, c i k co w
przykładzie 2. Wartości tym razem należy wprowadzić w oknie edycyjnym bloku
Transfer Fcn otwartym poprzez dwukrotne kliknięcie na bloczek. Okno ustawień
pokazane jest na rys. E3-2.
Rys. E3-2
Kliknięcie OK zamknie okno ustawień. W niniejszym przykładzie przedstawiony jest
kolejny blok wyjściowy Save Output To A File. Pozwala on na zapis wyników symulacji
bezpośrednio do pliku dyskowego. Tym razem zestaw parametrów jakie można
ustawić jest podobny do tego przedstawionego w poprzednim przykładzie. Różnica
polega jedynie na dodatkowym polu edycyjnym pozwalającym na wybór nazwy pliku
dyskowego. Plik zapisany na dysk w ten sposób ma typowy dla Matlaba format
binarny. Pliki tego typu często nazywane są mat-file z uwagi na rozszerzenie (.mat).
Użycie zawartości takiego pliku związane jest z wczytaniem go do Matlaba
poleceniem load NazwPlik. Po tej operacji w przestrzeni roboczej programu znajdują
się wszystkie zmienne zawarte w oryginalnym *.mat pliku.
W przypadku tego konkretnego przykładu mat-plik zostanie zachowany na dysku D: w
katalogu "temp". Plik zostanie nazwany "example3out.mat". Oczywiście nazwy i
położenie pliku na dysku może być dowolne. W okienku edycyjnym bloku należy
wpisać wszystkie powyższe informacje tzn: nazwę pliku wraz z lokalizacją oraz nazwę
zmiennej  w tym przypadku także "simout". Wartości pozostałych parametrów
można pozostawić niezmienione.
Rys E3-3
Krok 4
W kolejnym kroku można już uruchomić symulację poprzez kliknięcie przycisku
(alternatywnie można użyć skrót klawiaturowy CTRL+T ). W wyniku symulacji uzyskuje
się przebieg zaprezentowany na rys. E3-4.
Rys E3-4
Z porównania przebiegów rys E2-3 oraz rys. E3-4 wynika iż są one identyczne.
Oczywiście nie ma w tym nic dziwnego, obie symulacje dotyczą przecież tego samego
modelu.
Wykorzystanie danych symulacyjnych zachowanych w pliku dyskowym
Aby wykorzystać wyniki smulacji zapisane w postaci pliku należy najpierw dokonać ich
odczytu. Jak zostało to już powiedziane wcześniej odczytu dokonuje się poleceniem:
>> load example3out
Należy zauważyć iż taka postać powyższego polecenia jest poprawna jedynie wtedy
gdy lokalizacja wczytywanego pliku zawarta jest w predefiniowanym zestawie ścieżek
dostępu Matlaba lub w bieżącym katalogu. Jeśli tak nie jest polecenie powinno mieć
pełną formę zawierającą dodatkowo ścieżkę dostepu do pliku (w tym przypadku load
D:/temp/example3out.mat). Sprawdzenie czy zmienne zostały wczytane poprawnie
wykonane zostanie poleceniem
>> who
Jeśli wszystko przebiegło poprawnie Matlab powinien wyświetlić
Your variables are:
simout
"simout" jest w tym przypadku oczekiwaną przez nas macierzą zawierającą wyniki
wykonanej symulacji. Jej rozmiar w tym przypadku to 2 x m, gdzie m oznacza liczbę
kolumn i jest równe liczbie punktów na danych wygenerowanych w trakcie przebiegu
symulacji. Pierwszy rząd zawiera informacje o punktach na osi czasu natomiast drugi
odpowiadające im wartość. Po to aby wygenerować wykres odpowiedzi układu na
wymuszenie krokowe wystarczające jest użycie funkcji plot w odniesieniu do macierzy
"simout". Dla tego przykładu polecenie to będzie miało postać:
>> plot(simout(1,:),simout(2,:));grid
W wyniku jego działania narysowany zostanie przebieg pokazany na rys. E3-5.
Rys. E3-5
Podobny efekt przyniesie symulacja układu z wykorzystaniem funkcji przejścia w trybie
tekstowym:
Definiujemy licznik i mianownik transmitancji jako wektor wartości odpowiadających
kolejnym współczynnikom przy zmiennej s zaczynając od najwyższej jej potęgi:
L = [1]
M = [2 0.7 1]
Definiujemy funkcje przejścia:
sys = tf(L,M);
badamy podpowiedz
układu na:
Skok jednostkowy  step(sys,t), gdzie t oznacza czas symulacji.
Impuls  impulse(sys, t);
Dowolne wymuszenie  lsim(sys, u, t), gdzie u jest wektorem wymuszenia, a t wektorem
czasu.
Przestrzeń stanu.
Rozważmy układ jak na rysunku (fizyczny model np. elastycznego mocowania silnika).
Równania dynamiczne można zapisać następująco
Zdefiniujmy zmienne stanu jako x1(t) oraz x2(t)
Podstawiając za y(t) i y (t) zastępujemy równanie drugiego rzędu układem równań
pierwszego rzędu.
Równanie wyjść ma postać
W postaci macierzowej równania stanu dla układu możemy zapisać jako
Równanie wyjść:
Schemat blokowy układu przedstawia rysunek
Dla tak wyznaczonego modelu możliwa jest transformacja z opisu za pomocą przestrzeni
stanu do opisu w postaci transmitancji. Służy do tego polecenie:
[num,den] = ss2tf(A,B,C,D)
gdzie num oraz den jest odpowiednio licznikiem i mianownikiem transmitancji. Możliwa jest
również odwrotna transformacja:
[A,B,C,D] = tf2ss(num,den)
Zadania do samodzielnego wykonania.
1. Zamodeluj w środowisku Matlab układ jak na rysunku 2, przyjmując dowolne parametry
układu (m,c, k). Zasymuluj układ, używając jako wymuszenia:
" Funkcji harmonicznej
" Wymuszenia losowego
" Wymuszenia impulsowego
" Skoku jednostkowego.
2. Zamodeluj układ jak na rysunku poniżej, przyjmując dowolne wartości parametrów
mechanicznych. Przeprowadz
symulację jak w punkcie 1.
k1
k2
F
m1 m2
c1
c2
x2 x2