Laboratorium

identyfikacji

procesów

technologicznych

Wykonali :

Kierunek :

Numer wiczenia :

Temat :

Identyfikacja parametryczna

rekurencyjn metod

najmniejszych kwadratów.

Data wykonania :

Data oddania :

Ocena :

Celem wiczenia byo przeprowadzenie procesu identyfikacji parametrycznej danego obiektu dynamicznego rekurencyjn metod najmniejszych kwadratów (RMNK), opierajc si na modelu typu ARX. W tym celu naley wykorzysta odpowiednie funkcje dostpne w bibliotece Identification Toolbox Matlaba.

1. Wstp teoretyczny.

Wiadomoci dotyczce struktury modeli typu ARX oraz ich identyfikacji metod parametryczn, a dokadniej jedn z jej odmian, czyli metod najmniejszych kwadratów - zostay zamieszczone w wiczeniu nr 5.

Rozpatrujemy obiekt SISO. Po wykonaniu eksperymentu pomiarowego otrzymujemy zestaw danych pomiarowych, który przedstawiamy w postaci macierzy :

0x01 graphic

gdzie N oznacza ilo wykonanych pomiarów.

Na podstawie tych danych dymy do znalezienia modelu badanego ukadu. Zakadamy typ poszukiwanego model jako ARX. Celem identyfikacji jest wic wyznaczenie parametrów równania rónicowego :

0x01 graphic

Przyjmujemy take struktur modelu poprzez okrelenie wartoci wspóczynników m i n. Moemy zatem utworzy nastpujcy ukad równa :

0x01 graphic

Stosujc zapis macierzowy powyszego ukadu, otrzymujemy :

0x01 graphic

przy czym :

0x01 graphic
; 0x01 graphic
; 0x01 graphic

Poniewa zazwyczaj ilo pomiarów N >> n+m (w celu zmniejszenia wpywu zakóce), a wic otrzymany ukad równa jest ukadem nadokrelonym. Aby rozwizanie tego ukadu byo moliwe, stosujemy przeksztacenie sprowadzajc go do normalnego ukadu równa, dziki czemu uzyskujemy :

0x01 graphic

Równanie to jest wynikiem zastosowania MNK. Na jego podstawie funkcja arx oszacowuje parametry modelu. Zakadamy jednak, e po otrzymaniu N wyników pomiarowych, dochodzi jeszcze jeden (nowy) pomiar. Pomiar ten uwzgldniamy w naszym procesie identyfikacji, gdy moe on mie wpyw na wynik procesu. Tak wic okrelamy aktualny wektor parametrów modelu N+1. Wtedy :

0x01 graphic
; 0x01 graphic

przy czym : 0x01 graphic

Posugujc si w dalszym cigu `zwyk' MNK, otrzymamy :

0x01 graphic

Wida wic, e naley ponownie wykona operacji odwracania i mnoenia macierzy, których wymiary w tym przypadku ulegy ju powikszeniu o jeden. Pojawianie si nastpnych pomiarów spowoduje dalsze rozbudowywanie macierzy, w efekcie czego czas potrzebny na przeprowadzenie odpowiednich dziaa matematycznych ulega wydueniu (nawet przy wykorzystaniu nowoczesnych rodków obliczeniowych). Dymy zatem do rozwizania tego problemu w sposób pozwalajcy oszacowa nowe parametry modelu po kadorazowym otrzymaniu kolejnego pomiaru, bez koniecznoci ponownego rozwizywania powyszego skomplikowanego równania macierzowego, czyli :

0x01 graphic

Taki zapis rozpatrywanego problemu wskazuje na moliwo wykorzystania oblicze dokonanych wczeniej - na podstawie poprzednich pomiarów - oraz wyeliminowania operacji odwracania macierzy. Aby móc skorzysta z tego zapisu naley przeprowadzi nastpujce transformacje :

oznaczmy 0x01 graphic
; czyli 0x01 graphic
; 0x01 graphic
; 0x01 graphic
;

std 0x01 graphic

aby wydzieli z tego równania macierz PN+1 , naley skorzysta z nastpujcej tosamoci macierzowej algebry :

0x01 graphic

W wyniku tego uzyskujemy :

0x01 graphic

oznaczajc : 0x01 graphic

mamy :

0x01 graphic

Poniewa :

0x01 graphic

std :

0x01 graphic

skd po przeksztaceniach otrzymamy :

0x01 graphic

przy czym :

0x01 graphic
- rzeczywista (zmierzona) warto sygnau wyjciowego w chwili N+1;

0x01 graphic
- warto sygnau wyjciowego w chwili N+1 oszacowana w oparciu o model;

Równanie to stanowi istot rekurencyjnej metody najmniejszych kwadratów, która reprezentowana jest w pakiecie Matlab poprzez funkcj rarx. Wynika z niego take, e RMNK umoliwia cig aktualizacj parametrów modelu w czasie, kiedy dokonywane s pomiary sygnau wejciowego i wyjciowego obiektu. Oznacza to moliwo stosowania tej metody w celach identyfikacji modeli obiektów których parametry s zmienne w czasie, czyli obiektów niestacjonarnych.

2. Przebieg wiczenia.

2.1. Obiekt I .

2.1.1. Zebranie i przygotowanie danych pomiarowych.

0x01 graphic

Rys.1 Schemat blokowy badanego obiektu.

0x01 graphic

Rys.2. Schemat blokowy ukadu pomiarowego.

Zebrane dane umieszczane s w przestrzeni roboczej Matlaba w postaci macierzy z1 o strukturze :

0x01 graphic

gdzie N=300 - ilo wykonanych pomiarów.

Usuwamy z nich warto redni :

z11 = dtrend(z1)

Funkcja :

plot(z11)

daje nam moliwo zaobserwowania danych pomiarowych w postaci graficznej :

0x01 graphic

Rys.3.

2.1.2. Przyjcie struktury modelu.

W naszym wiczeniu przyjmujemy struktur modelu ARX, okrelon jako wektor nn=[2 2 1]. Ogólna posta modelu ARX dla zaoonej struktury :

0x01 graphic

gdzie k = 1..N.

2.1.3. Szacowanie (estymacja) parametrów modelu.

2.1.3.1. Metoda najmniejszych kwadratów (MNK).

W pierwszym kroku tego punktu wiczenia w celu oszacowania parametrów modelu o okrelonej wyej strukturze, zastosujemy `zwyk' metod najmniejszych kwadratów, wywoujc funkcj :

th = arx(z11,nn);

Wynikiem zastosowania tej funkcji jest utworzenie specjalnie skonstruowanej macierzy THETA, z której moemy uzyska potrzebne informacje korzystajc z nastpujcych funkcji :

[num,den] = th2tf(th)

printsys(num,den,'z')

W ten sposób otrzymalimy transmitancj dyskretn :

0x01 graphic

a pena posta modelu:

0x01 graphic

2.1.3.2. Rekurencyjna metoda najmniejszych kwadratów (RMNK).

W celu identyfikacji modelu ARX rekurencyjn metod najmniejszych kwadratów wywoujemy funkcj :

thm = rarx(z11,nn,adm,adg)

przy czym :

adm - oznacza mechanizm adaptacji (u nas jest on oparty na wspóczynniku `zapominania', czyli `ff')

adg - oznacza wspóczynnika `zapominania' (w naszym przypadku ``), którego warto okrela w jakim stopniu pomiary wykonane w stosunkowo odlegych chwilach czasu, wpywaj na warto aktualnie szacowanych parametrów modelu, i tak dla `` = 1 wszystkie wykonane pomiary posiadaj jednakowe znaczenie, a wic w tym przypadku polecenie rarx funkcjonuje tak samo jak polecenie arx, natomiast `` < 1 oznacza zmniejszony wpyw `starych' pomiarów na proces identyfikacji. Taka waciwo jest szczególnie cenna kiedy analizowany obiekt jest niestacjonarny.

Wywoujemy wic funkcj rarx dla `` = 1 :

thm = rarx(z11,nn,'ff',1);

Funkcja ta zwraca macierz thm o wymiarach N×4 (poniewa kada kolumna oznacza jeden z parametrów modelu, a w nasz przypadku model posiada cztery parametry : a1, a2, b1, b2). Wygenerowaa ona zatem zestaw parametrów modelu po kadorazowym pojawieniu si nowego pomiaru, korzystajc z oblicze wykonanych dla wczeniejszych pomiarów. Ostatni wiersz otrzymanej macierzy thm zawiera parametry modelu wyznaczone na podstawie N=300 pomiarów i przedstawia si on w sposób nastpujcy :

a1 a2 b1 b2

-1,4105 0,4412 0,0038 0,0184

Aby móc obserwowa rozkad zmiennoci parametrów modelu w zalenoci od liczby dokonywanych pomiarów, wykonujemy funkcj :

plot(thm)

0x01 graphic

Rys.4.

Podobnie postpujemy dla :

- `` = 0,99, czyli

thm = rarx(z11,nn,'ff',0.99);

Ostatni wiersz otrzymanej macierzy thm :

a1 a2 b1 b2

-1.3493 0.3846 0.0005 0.0242

Wykonujemy funkcj :

plot(thm)

0x01 graphic

Rys.5.

- `` = 0,95, czyli

thm = rarx(z11,nn,'ff',0.95);

Ostatni wiersz otrzymanej macierzy thm :

a1 a2 b1 b2

-1.3455 0.3881 -0.0110 0.0393

Wykonujemy funkcj :

plot(thm) 0x01 graphic

Rys.6.

- `` = 0,93, czyli

thm = rarx(z11,nn,'ff',0.93);

Ostatni wiersz otrzymanej macierzy thm :

a1 a2 b1 b2

-1.3438 0.3854 -0.0132 0.0408

Wykonujemy funkcj :

plot(thm)

0x01 graphic

Rys.7.

2.2. Obiekt II .

2.2.1. Zebranie i przygotowanie danych pomiarowych.

0x01 graphic

Rys.8 Schemat blokowy badanego obiektu.

0x01 graphic

Rys.9. Schemat blokowy ukadu pomiarowego.

Zebrane dane umieszczane s w przestrzeni roboczej Matlaba w postaci macierzy z1 o strukturze :

0x01 graphic

gdzie N=1000 - ilo wykonanych pomiarów.

Usuwamy z nich warto redni :

z11 = dtrend(z1)

Funkcja :

plot(z11)

daje nam moliwo zaobserwowania danych pomiarowych w postaci graficznej :

0x01 graphic

Rys.10.

2.2.2. Przyjcie struktury modelu.

Przyjmujemy - podobnie jak wczeniej - struktur modelu ARX, okrelon jako wektor nn=[2 2 1]. Ogólna posta modelu ARX dla zaoonej struktury :

0x01 graphic

gdzie k = 1..N.

2.2.3. Szacowanie (estymacja) parametrów modelu RMNK.

Wywoujemy wic funkcj rarx (tym razem) dla `` = 0,99 :

thm = rarx(z11,nn,'ff',0.99);

Funkcja ta zwraca macierz thm o wymiarach N×4. Ostatni wiersz otrzymanej macierzy thm zawiera parametry modelu wyznaczone na podstawie N=1000 pomiarów i przedstawia si on w sposób nastpujcy :

a1 a2 b1 b2

-1.6255 0.8089 0.0495 0.1272

Stosujemy funkcj aby otrzyma posta `graficzn' macierzy thm :

plot(thm)

0x01 graphic

Rys.11.

2.3. Obiekt II .

2.3.1. Zebranie i przygotowanie danych pomiarowych.

Schemat ukadu pomiarowego jest w tym przypadku identyczny jak w poprzednim punkcie (Rys.9.)

0x01 graphic

Rys.12. Schemat blokowy badanego obiektu.

Zebrane dane umieszczane s w przestrzeni roboczej Matlaba w postaci macierzy z1 o strukturze :

0x01 graphic

gdzie N=1000 - ilo wykonanych pomiarów.

Usuwamy z nich warto redni :

z11 = dtrend(z1)

Funkcja :

plot(z11)

powoduje :

0x01 graphic

Rys.13.

2.3.2. Przyjcie struktury modelu.

Przyjmujemy - ponownie - struktur modelu ARX, okrelon jako wektor nn=[2 2 1]. Ogólna posta modelu ARX dla zaoonej struktury :

0x01 graphic

gdzie k = 1..N.

2.3.3. Szacowanie (estymacja) parametrów modelu RMNK.

Wywoujemy wic funkcj rarx (tym razem) dla `` = 1 :

thm = rarx(z11,nn,'ff',1);

Funkcja ta zwraca macierz thm o wymiarach N×4. Ostatni wiersz otrzymanej macierzy thm zawiera parametry modelu wyznaczone na podstawie N=1000 pomiarów i przedstawia si on w sposób nastpujcy :

a1 a2 b1 b2

-1.6143 0.8169 0.0582 0.1331

Stosujemy funkcj aby otrzyma posta `graficzn' macierzy thm :

plot(thm)

0x01 graphic

Rys.14.

3. Tre m-pliku.

z11 = dtrend(z1);

nn = [2 2 1];

th = arx(z11,nn);

[num,den] = th2tf(th);

printsys(num,den,'z')

%thm = rarx(z11,nn,'ff',1)

%thm = rarx(z11,nn,'ff',0.99)

%thm = rarx(z11,nn,'ff',0.95)

thm = rarx(z11,nn,'ff',0.93)

plot(thm)

4. Uwagi i wnioski.