Ghostscript wrapper for D:\Digitalizacja\MTS87_t25_z1_4_PDF_artykuly\02mts87_t25_zeszyt_3.pdf M E C H AN I K A TEORETYCZNA I STOSOWANA 3, 25, (1987) IDEN TYFIKACJA PARAMETRYCZNA M OD ELU MATEMATYCZN EG O SAM OLOTU WŁAD YSŁAW JAROM IN EK Polska Akademia N auk, W arszawa TAD EU SZ STEFAŃ SKI Politechnika Ś wię tokrzyska, Kielce 1. Wprowadzenie Jednym z podstawowych problemów syntezy ukł adu stabilizacji samolotu jest problem identyfikacji, którego zadaniem jest otrzymanie niezbę dnych informacji o charakterysty- kach statycznych i dynamicznych samolotu. Problem ten posiada szczególnie istotne znaczenie w przypadku samolotu naddź wię kowego, gdzie informacje o aktualnych zmia- nach jego parametrów muszą być cią gle w czasie lotu korygowane. Identyfikacja obiektów niestacjonarnych — a takim jest samolot — sprowadza się gł ównie do zał oż enia quasi- stacjonarnoś ci w pewnym dopuszczalnym przedziale czasu AT i na znalezieniu takiego algorytmu identyfikacji, którego czas procesu przejś ciowego jest mniejszy od czasu AT . W niniejszym artykule do syntezy algorytmów identyfikacji parametrów modelu ma- tematycznego samolotu w kanale podł uż nym zastosowano metodę najmniejszych kwa- dratów. M etoda ta wykazuje w tym przypadku wiele istotnych zalet, w tym : otrzymanie algorytmu identyfikacji w postaci ukł adu równań normalnych lub w postaci rekurencyj- nej — prostej i wygodnej do obliczeń numerycznych, moż liwość uwzglę dnienia szumów pomiarowych lub zakł óceń dział ają cych na samolot, wysoka zbież ność i dokł adność identyfikacji parametrów modelu matematycznego. Ponieważ przedstawione algorytmy pozwalają identyfikować parametry dyskretnego modelu matematycznego samolotu, w pracy przedstawiono także problem okreś lania parametrów cią gł ego modelu matema- tycznego na podstawie modelu dyskretnego. 2. Model matematyczny samolotu W wielu praktycznych zastosowaniach, np. w przypadku stabilizacji ukł adu pilotażu rę cznego samolotu w kanale podł uż nym, istotne znaczenie ma ruch krótkookresowy. Wynika to z faktu, że pierwszą reakcją na wychylenie steru wysokoś ci jest zmiana ką ta natarcia, podczas gdy zmiana prę dkoś ci lotu zachodzi znacznie póź niej. 430 W. JAROMIN EK, T. STEFAŃ SKI Liniowy m odel matematyczny krótkookresowego ruchu podł uż nego samolotu, dla ustalonego lotu poziomego, ma postać [2]: 4 = óc + P a a + P ó d + P w (1) gdzie: a — ką t n atarcia; # — ką t pochylenia; <5— wychylenie steru; n — przyś pieszenie n o rm aln e; V— prę dkość lo t u ; g — przyś pieszenie ziemskie; P a , P s , Mi, M a , M- a , M a ~ współ czynniki zależ ne od geometrii, kinematyki i parametrów lotu; P w , M w — zakł ó- cenia dział ają ce na sam olot. Z równań (1) otrzymano transmitancje wią ż ą ce sygnał y wyjś ciowe • & i n z wejś ciowym sygnał em sterują cym d dla P„ = 0 oraz M w = 0, które mają postać: d(s) "" ' d(s) gdzie: H + MJ = J l f P " M sPa- MaPt M 5 P a - M a P ó V A = A Współ czynniki transmitancji (2) i (3) są zmienne i zależ ne od prę dkoś ci V oraz wysokoś ci H lo t u ; ogólnie stwierdza się , że są funkcją czasu. U kł ad równ ań (1) m oż na zapisać w postaci macierzowej x(t)=A c x(t)+Bj(t)+W c (4) przy czym: Również w tym przypadku param etry modelu (4), tj. macierze A c , B c i W c są funkcjami czasu. 3. Algorytmy identyfikacji Sam olot jest obiektem, którego wł aś ciwoś ci statyczne i dynamiczne zmieniają się w czasie lotu, n a ogół jedn ak zmiany są powolne. Z tych też wzglę dów dopuszczalne jest zał oż enie, że istnieje przedział czasu AT , w którym parametry modelu matematycznego I D E N T YF I K AC JA P AR AM E TR YC Z N A SAM OLOTU 431 (2) lub (4) nie ulegają istotnym zmianom. D o identyfikacji parametrów modelu m ate- matycznego należy wówczas zastosować metodę , której proces przejś ciowy nie bę dzie dłuż szy niż AT . Tak postawiony problem rozwią zany zostanie przy pomocy metody najmniejszych kwadratów, której algorytmy cechują się wysoką zbież noś cią i dokł adnoś cią, identyfikacji. W warunkach ą uasistacjonarnoś ci, dyskretna postać równania (4) okreś lona jest: nastę pują co x(k+1) - Ax(k)+Bó(k) + W (5> gdzie: A = cxv(A c T ) B = A^ (exp(A c T )- T )B e W = Ać1 (expO4c T ) - T )W e W przypadku mał ej wartoś ci okresu impulsowania T moż na ograniczyć się do przybli- ż enia liniowego rozwinię cia macierzy QXp(A c T ) w szereg potę gowy i wóczas A = I+A C T W = W C T Identyfikację parametrów modelu matematycznego (5) moż na przeprowadzić w opar- ciu o pomiar ó(k) i &(k) lub w oparciu o pomiar d{k) oraz x(k). W metodzie bazują cej na pomiarze skalarnych sygnał ów (tj. d(k) i &(k)) konieczna jest przedstawienie modelu (5) w postaci kanonicznej F robeniusa. D la kwadratowego- wskaź nika jakoś ci algorytm identyfikacji ma postać: k- \ - )] k)P(k- l)G T (k)]- l G(k)P(k- l) l ) gdzie:
(k- l)+P(k)G
r
(k)Q(k)t&(k)~G(k)ę (k~l)]
P(k) = R(k)- R(k)GT(k)[Q- i(k) + G(k)R(k)GT(k)]- 1G(k)R(k) ( 9 )
przy czym najczę ś ciej przyjmuje się R(k) = yl, gdzie / jest macierzą jednostkową a y
współ czynnikiem o duż ej wartoś ci (rzę du 103- ^10s). Zbież ność tego algorytmu moż na
dodatkowo poprawić dobierając eksperymentalnie ciąg wartoś ci współ czynników wagi
Q{k).
W przypadku, gdy niestacjonarność identyfikowanych parametrów jest silna, algorytmy
rekurencyjne najczę ś ciej są rozbież ne. Lepsze wyniki osią ga się stosując tzw. ukł ad równań
normalnych. N iech A T bę dzie przedział em czasu, w którym parametry obiektu nie ulegają
istotnym zmianom. D la tego przedział u estymator identyfikowanych parametrów, uzyskany
metodą najmniejszych kwadratów ma postać
£ = [GT(N )G(N )]~ 1GT(N )y(N ) (10)
gdzie:
G(N ) = [G '(l)|(?r(2)! ... \ GT(N )]T
y(N )= $ ( 1) , *( 2) , . . . , *( # ) ]
Idea tej metody polega na tym, że w każ dym przedziale czasu A T należy dokonać N pomia-
rów sygnał u wejś ciowego d(k) i wyjś ciowego &{k) w chwilach dyskretnych k, przy czym
N > 4, a nastę pnie przetwarzać je zgodnie z zależ noś cią (10). Najczę ś ciej równania (10)
rozwią zuje się w oparciu o dekompozycję macierzy, n p. metodą ortogonalizacji H ousehol-
dera lub metodą SVD. Zwię kszenie iloś ci pomiarów N ma duże znaczenie w przypadku
dział ania n a ukł ad zakł óceń i szumów pomiarowych. Wpł yw wymienionych zjawisk na
dokł adność identyfikacji moż na zmniejszyć stosując odpowiednią filtrację.
Identyfikacji parametrów modelu matematycznego samolotu moż na także dokonać
poprzez pomiar sygnał u wejś ciowego d{k) oraz wektora stanu x{k). Wówczas algorytm
identyfikacji ma postać:
$r(fc) = y)T(k- \ ) + P(k)S(k)Q{k)[xT{k)- ST(k)y>T{k- \ )]
P(k) =
przy czym:
S(k) =
\ x(k-
r- i)-J
Również w tym przypadku macierz kowariancji P(k — 1) moż na zastą pić macierzą R(k),
analogicznie jak w algorytmie (9). Zaletą algorytmu (11) jest wyż sza, niż algorytmu (9),
zbież ność I dokł adność identyfikacji, wadą —• konieczność pomiaru ką ta pochylenia
• &(k) i ką ta natarcia a.(k). G dy zmiany parametrów są szybkie, zamiast algorytmu reku-
rencyjnego (11) należy zastosować ukł ad równań normalnych, wynikają cych z metody
.najmniejszych kwadratów.
W wielu przypadkach zachodzi konieczność okreś lania parametrów cią gł ego modelu
matematycznego. Cią gły model matematyczny obiektu moż na uzyskać jako przypadek
I D E N T YF I K AC JA P AR AM E TR YC Z N A SAM OLOTU 433
graniczny m odelu dyskretn ego, przyjmują c, że okres im pulsowan ia T dą ży do zera. Wów-
czas macierze A
c
i B
c
ró wn an ia (4) m oż na okreś lić z zależ noś ci (6)
A
e
= (Ar- JOT -
1
B
c
= BT '
1 . ( 1 2 )
, , - , • , . t • • •
. • - . . , . • • • - • . . . - • - . . - - , • .
Jeż eli okres im pulsowan ia T ma. dużą wartoś ć, czł ony wyż szego rzę du rozwin ię cia m acierzy
exp(A
c
T ) mają duży wpł yw n a A (nie są jeszcze bliskie zeru) i zależ n oś ci (12) n ie m o ż na
stosować, gdyż bł ą d obliczenia macierzy A
c
i B
c
jest zbyt duż y. R ówn ież ba r d zo m ał a
wartość T n ie jest wskazan a, przy stosowan iu tych zależ noś ci, gdyż bł ą d iden tyfikacji
param etrów m o d elu dyskretn ego silnie przen osi się n a wartoś ci p aram et ró w m o d elu
cią gł ego.
M acierze A
c
i B
c
okreś lić m o ż na ze stosun kowo m ał ym bł ę dem, korzystają c z ró wn ań
(5) i t a k:
A
c
= - jrlnA
(13)
W - Ij- UeB
Zasadniczą trudn ość sprawia t u obliczenie m acierzy A
c
. Jeś li wartoś ci wł asn e Aj m acierzy
A speł niają waru n ek
\ Xj- \ \ < 1 y - I , 2 (14)
to sł uszna jest zależ n ość
kT
fc- i
D la ukł adów stabiln ych i lj > 0 warun ek (14) jest zawsze speł n ion y. M acierz B
c
wyznaczyć
m oż na z zależ noś ci (13), uprzedn io okreś lając m acierz exp(A
c
T ), lu b też z zależ noś ci
przybliż onej
co
[ 2 ] e B (16)
4. Badania symulacyjne algorytmów identyfikacji
P rzedstawion e algorytm y identyfikacji przetestowan o n a E M C w celu ocen y ich zbież-
noś ci i dokł adn oś ci identyfikacji. Analizy d o ko n an o dla stacjon arn ego o raz n iestacjo-
n arn ego m odelu m atem atyczn ego. W pierwszym przypadku ba d a n o p o d st awo we wł aś-
ciwoś ci algorytm ów identyfikacji, w drugim — przydatn ość tych algorytm ów d o iden ty-
fikacji m odeli n iestacjon arn ych .
Symulowany obiekt identyfikacji okreś lony został przez m odel m at em at yczn y ( 4) ;
param etry tego m odelu dla zakresu prę dkoś ci lot u M a = 0,5- r- 1,8 i zakr esu wysokoś ci
lot u H — 0- 7- 13500 m zaczerpn ię to z pracy [4].
8 Mech. Teoret. i Stos. 3/ 87
434 W. JAROMIN EK, T. STEFAŃ SKI
Przyję to nastę pują ce wartoś ci parametrów identyfikowanego obiektu stacjonarnego
[0,95 - 0,2J _ ro,281
^ - [ 0 ,1 0,95 J ^ ~ i0,26j
W wyniku symulacji stwierdzono (rys. 1), że algorytm (11) wykazuje dużo lepszą zbież-
ność niż algorytm (7). Również w przypadku uwzglę dnienia szumów pomiarowych (rys. 2)
a)
T
0.5 1.0
kTts]
R ys. 1. Zbież ność procesu identyfikacji dla obserwacji:
a) — skalarnego sygnał u wyjś ciowego,
b) — peł nego wektora stanu.
1.5
kTls]
R ys. 2, Z bież ność procesu identyfikacji w przypadku uwzglę dnienia szumów pomiarowych dla obserwacji:
a) — skalarnego sygnał u wyjś ciowego,
b) — peł nego wektora stanu.
I D E N T YF I K AC JA P AR AM E TR YC Z N A SAM O LO T U 435
.— wniosek jest identyczny. Podczas symulacji przyję to, że szumy pomiarowe są sygnał ami
losowymi o zerowych wartoś ciach oczekiwanych i maksymalnych amplitudach nie prze-
kraczają cych 10% amplitud mierzonych sygnał ów.
Zbież ność algorytmu identyfikacji dla obiektu niestacjonarnego przedstawiono na rys. 3.
0 15 20 25 30 35 40 45 50 55 60 65 70
Rys. 3. Identyfikacja obiektu niestacjonarnego przy obserwacji:
a) — skalarnego sygnał u wyjś ciowego,
b) — peł nego wektora stan u.
Podczas symulacji procesu identyfikacji zał oż ono, że samolot nabiera prę dkoś ci do war-
toś ci 2Ma i wysokoś ci do 15000 m w cią gu 75s. Linią cią głą oznaczono zmiany wartoś ci
parametrów obiektu, natomiast linią przerywaną — zmiany wartoś ci parametrów samolotu
otrzymanych z procesu identyfikacji. Również w tym przypadku dużo wyż szą zbież ność
i dokł adność identyfikacji wykazuje algorytm (11). Algorytm (7) przy duż ej szybkoś ci
zmian parametrów obiektu utracił zbież noś ć. Oczywiś cie lepsze wyniki uzyska się stosują c
algorytm (8) lub (10).
5. Podsumowanie
Przedstawione algorytmy identyfikacji moż na z powodzeniem wykorzystać do iden-
tyfikacji parametrów stacjonarnego i niestacjonarnego modelu matematycznego. Jeś li
jest moż liwość obserwacji ką ta natarcia