Ghostscript wrapper for D:\BBB-ARCH\ARCHIWUM-lata-78-71\MTS76_t14z1_4\mts76_t14z1.pdf M E C H A N I K A T E O R E T Y C Z N A I S T O S O W A N A 1. 14 (1976) O PEWNEJ NOWEJ METODZIE ANALIZY STATECZNOŚ CI ROZWIĄ ZAŃ UKŁADÓW NIELINIOWYCH O JEDNYM STOPNIU SWOBODY A L I C J A P I E N I Ą Ż E K, WIESŁAW P I E N I Ą Ż EK ( K R A K Ó W ) 1. W s t ę p W pracy [1] została przedstawiona pewna nowa metoda analizy statecznoś ci nielinio wych u k ł a d ó w dynamicznych opisanych równaniami róż niczkowymi zwyczajnymi dru giego rzę du. Autorzy oparli się na nastę pują cym rozumowaniu: każ dy układ dynamiczny m o ż na uważ ać jako pewne pole dynamiczne, które działa z pewną siłą na znajdują cy się w nim punkt materialny. Stateczność tego pola dynamicznego zależy od tego, czy energia zgromadzona przez punkt materialny roś nie czy też maleje pod wpływem sił pola. Rozumowanie to doprowadziło do opracowania pewnego, wygodnego w stosowaniu, algorytmu, na podstawie którego m o ż na wnioskować o statecznoś ci u k ł a d u . Algorytm ten umoż liwia także badanie statecznoś ci cykli granicznych. W niniejszej pracy zostanie przedstawiona idea metody, zgodnie z [1], jej p o r ó w n a n i e z innymi, istnieją cymi metodami, oraz jej zastosowanie do badania statecznoś ci r ó w n a n i a typu Rayleigha z nieliniową charakterystyką sprę ż ystoś ci. 2. Opis metody [1) i jej sens fizyczny Niech układ dynamiczny bę dzie opisany równaniem róż niczkowym zwyczajnym rzę du drugiego w postaci (2.1) x=f(x,k) lub równoważ nym mu układem równań pierwszego rzę du (2.2) k = y, y=f(x,y), gdzie f(x,y) jest funkcją nieparzystą, ze wzglę du na x, w ogólnym przypadku nieliniową. Jeż eli układ ten bę dziemy uważ ali za pewne pole dynamiczne, to wówczas siła tego pola wyniesie (2.3) F = mx = my = mf(x, y), gdzie m jest masą punktu materialnego znajdują cego się w polu. Wobec tego, równanie (2.1) m o ż na zapisać w postaci (2.4) mx = mf(x,y). 138 A . P I E N I Ą Ż E K, W . P I E N I Ą Ż EK Przekształcimy powyż sze równanie w ten sposób, że dodamy i odejmiemy wyraż enie mx (o wymiarze siły) do prawej strony. Wyraż enie to zinterpretujemy jako siłę zachowa wczą. Otrzymujemy (2.5) mx = mf(x, x) + mx — mx. Jeż eli w powyż szym r ó w n a n i u wprowadzimy oznaczenie (2.6) Fy = mf(x, x) + mx, to wówczas ruch punktu materialnego bę dzie się odbywał pod wpływem siły zachowa wczej i siły Fx, co m o ż na ująć zależ noś cią (2.7) mx = Fx—mx. Jeś li obliczymy pracę siły Ft podczas jednego, pełnego okresu ruchu T, to w przypadku gdy jest ona dodatnia punkt materialny znajdują cy się w polu dynamicznym powię ksza swoją energię, co oznacza, że układ (pole) jest niestateczny. W przeciwnym przypadku, u k ł a d posiada cechy statecznoś ci (gdy praca siły F , jest ujemna, to wówczas układ wydaje p r a c ę i energia punktu materialnego maleje). Z r ó w n o w a ż my siłę Ft siłą F0 równą co do wartoś ci, lecz przeciwnie skierowaną, tak aby punkt materialny wykonywał ostatecznie ruch zachowawczy opisany r ó w n a n i e m (2.8) x = x lub r ó w n o w a ż n ym mu u k ł a d e m r ó w n a ń pierwszego rzę du (2.8a) x = у , у = —x. Znak pracy wykonywanej przez siłę F0 bę dzie przeciwny do znaku pracy siły Fx, a więc gdy praca LF0 siły F0 bę dzie ujemna, to układ bę dzie niestateczny, zaś w przeciwnym przypadku u k ł a d bę dzie stateczny. Obliczymy obecnie pracę siły F0 w jednym, pełnym okresie ruchu T. Otrzymujemy T (2.9) LF0 = f F0xdt. o Ponieważ jest (2.10) F0 = Ft = m[f(x,x)+x], to praca ta wyniesie T (2.9a) LF0 = m j [f(x, x) + x]xdt. o Uwzglę dniając (2.2) otrzymamy г (2.9b) LF0 = mf(xX+yy)dt. o Jak powiedzieliś my wcześ niej, ruch bę dzie odbywał się po okrę gu x = r c o s 0 , (2.11) . _ r = const > 0. у = rsmG, O M E T O D Z I E A N A L I Z Y S T A T E C Z N O Ś CI R O Z W I Ą Z AŃ 139 Róż niczkując którekolwiek z r ó w n a ń (2.11) wzglę dem czasu i wykorzystując odpowied nio zależ ność (2.8a) otrzymujemy (2.12) ^ = 1 , d8= dt. Wprowadzimy obecnie, zgodnie z [1], funkcję xx+yy xy+yf(x,y) (2.13) S'(x,y) = ]/x2+y2 \/x2+y2 Jak widać, jest to składowa radialna prę dkoś ci fazowej: vr = r. Z a u w a ż my teraz, że wyraż enie podcałkowe w (2.9b) m o ż na zapisać, wykorzystując (2.13) , (2.14) xx+yy = / ^ T ^ j ) . Biorąc pod uwagę (2.14), a także (2.11) i (2.12), po zmianie granic całkowania wzór (2.9b) przyjmie p o s t a ć (2.15) LF0 = mrj S(r, G)dG, m > 0, r > 0. o Wprowadzimy obecnie funkcję 2л (216) g(r) = f S(r, &)d9. o Wobec tego, praca LF0 posiada znak przeciwny do znaku wyraż enia (2.16), k t ó r e łat wo obliczymy znając funckję S(r, &). Funkcję g(r) wykorzystamy zatem do okreś lania charakteru statecznoś ci u k ł a d u (2.1); mianowicie; — jeż eli g(r) = 0, V/ > 0 to układ jest zachowawczy, (2.17) —jeż eli g(r) > 0, V/ > 0 (LF0 < 0), to układ jest niestateczny, — jeż eli g(r) < 0, V/ > 0 (LF0 > 0), to układ jest stateczny. Przy pomocy tej metody m o ż na b a d a ć stateczność cykli granicznych. Posługujemy się tutaj także funkcją g(r). D l a statecznego cyklu granicznego mamy: g(r) > 0 dla 0 < r < r0, (218) g(r0) = 0, g(r) < 0 dla r > r0. D l a niestatecznego cyklu granicznego jest: g(r) < 0 dla 0 < /• < r0 (2.19) g(r0) = 0, g(r) > 0 dla r > r0. 140 A . P I E N I Ą Ż E K, W . P I E N I Ą Ż EK Wprowadzimy jeszcze funkcję: R'(x, y) okreś loną wzorem (2.20) R'(x,y) Х ' У ~ к у \/x2+y2 \/x2+y2 Jest to składowa transwersalna v@ — rŚ prę dkoś ci fazowej. W dalszych rozważ aniach bę dziemy stosować funkcje R'(x,y), S'(x,y) lub R(r,0), S(r,0) otrzymane po podsta wieniu w (2.13) i (2.20):x = rcosG, у = rsin<9. O statecznoś ci u k ł a d u m o ż na są dzić na podstawie analizy funkcji S(r, 0) i R(r, 0). Mianowicie, zgodnie z [1], jeż eli w punktach, gdzie R(r, 0) = 0 funkcja S(r, 0) < 0, to wówczas układ jest stateczny. G d y w punktach, w których R(r,0) = 0 funkcja S(r, 0) > 0, to układ jest niestateczny. Interpretację fizyczną powyż szych w a r u n k ó w m o ż na p o d a ć opierając się na pracy [2]. Podano w niej metodę dwóch funkcji, pozwalają cą na ocenę statecznoś ci u k ł a d u na podstawie przebiegu trajektorii fazowych, które m o ż na przewi dzieć na podstawie wykresów tych funkcji na płaszczyź nie fazowej. W pracy [2] układ (2.2) zinterpretowano jako u k ł a d okreś lają cy prę dkość punktu ma terialnego na płaszczyź nie fazowej. Wspomniane dwie funkcje są okreś lone w z o r a m i 0 (2.21) Ф (х ,у )~х у +у /(х ,у ), (222) W(x,y) = xf(x,y)y2. Przedstawiają one odpowiednio: iloczyn skalarny i współrzę dną iloczynu wektorowe go odległoś ci т (х , y) punktu na trajektorii fazowej od począ tku u k ł a d u współrzę dnych i wektora prę dkoś ci fazowej v(x,y) tego punktu. Jak łatwo stwierdzić, licznik (2.13) jest to funkcja Ф (х , у ), zaś licznik (2.20) jest to funkcja W(x, y). Znając wykresy tych funkcji na płaszczyź nie fazowej moż emy są dzić o przebiegu trajektorii fazowych u k ł a d u (2.2). W cytowanej pracy znajduje się szereg przykładów takiej analizy dla róż nych typów u k ł a dów. Porównajmy analizę statecznoś ci na podstawie metody dwóch funkcji w [2] z analizą na podstawie funkcji R'(x,y) i S'(x,y) w nowej metodzie. Warunek R'(x, y) = 0 oznacza to samo со \P(x, у ) = 0, czyli że wektory г i v pokry wają się, a ich iloczyn wektorowy r x v = 0. Jeż eli równocześ nie przy tym S'(x,y) < 0, a więc Ф (х ,у ) < 0, to układ jest stateczny. Istotnie, iloczyn skalarny r x v < 0 , a to ś wiadczy o tym, że wektor v posiada zwrot przeciwny do wektora r i punkt porusza się w kierunku począ tku u k ł a d u współrzę dnych. Przeciwny znak S'(x,y), przy równoczes nym R'(x, y) = 0, wskazuje na zgodne zwroty r i v i punkt oddala się wówczas od po czą tku u k ł a d u współrzę dnych, a to ś wiadczy o niestatecznoś ci u k ł a d u dynamicznego. Dalsze p o r ó w n a n i a z innymi metodami podamy póź niej, a obecnie przedstawimy wspomniany algorytm z [1] w jego pełnej, usystematyzowanej postaci. A . l . Tworzymy funkcje S'(x, y), R'(x, y) S ' ( X > Y) = ХУ +УЛХ>У \ \/x2+y2 " W pracy [2] funkcja f(x, y) ma znak przeciwny, p o n i e w a ż u k ł a d (2.2) jest zapisany w postaci: x = у , у = —Д .х, у ), co nie ma w p ł y w u na nasze dalsze r o z w a ż a n i a. O M E T O D Z I E A N A L I Z Y S T A T E C Z N O Ś CI R O Z W I Ą Z AŃ 141 R'Qc, У ) = » \/x2+y2 lub funkcje R(r, 0), S(r, 0), po podstawieniu w powyż szych wzorach A = rcos0, у = = r s i n 0 . A.2. Badamy powyż sze funkcje; jeż eli w punktach, w k t ó r y c h : a) R(r,0) = 0 funkcja S(r,0) < 0, to układ jest stateczny, b) R(r, 0) = 0 funkcja S(r, 0) > 0, to u k ł a d jest niestateczny. A.3. Tworzymy funkcję g(r) g(r) = f S(r, 0)d0, r = const, o A.4. N a podstawie charakterystyki funkcji g(r) oceniamy stateczność u k ł a d u . Jeż eli: a) g(r) < 0 V> > 0 i g(r) > — oo, to układ jest stateczny, r>oo b) g(r) = 0 W > 0, to układ jest zachowawczy, c ) Ł ( r ) > 0 V r > 0, to układ jest niestateczny. d) układ posiada stateczny cykl graniczny, jeż eli: g(r) > 0 dla 0 < r < r0, g(r0) = 0, g(r) < 0 dla r > r0; e) układ posiada niestateczny cykl graniczny, jeż eli: g(r) < 0 dla 0 < r < r0, g(r0) 0, g(r) > 0 dla r > r0. Jeż eli funkcja R(r, 0) jest ujemna wszę dzie, to wówczas punkt A . 4 powyż szego al gorytmu daje pełną odpowiedź o statecznoś ci u k ł a d u . Autorzy niniejszej pracy zbadali przydatność przedstawionej metody do okreś lania statecznoś ci układów o róż nych typach p u n k t ó w osobliwych. Jak się okazało , dla ukła d ó w o siodłowym punkcie osobliwym metoda ta nie daje odpowiedzi o charakterze sta tecznoś ci. W pozostałych przypadkach metoda daje wyniki zgodne z otrzymanymi przy zastosowaniu innych metod. N a zakoń czenie tego rozdziału podamy dalsze p o r ó w n a n i a i moż liwoś ci zastosowania metody. Napiszmy równanie trajektorii dla u k ł a d u (2.2). Po wyrugowaniu czasu otrzymujemy (2.23) 4 = • а х у W metodzie krzywych stykowych Poincare'go, (patrz np. [3]), bierze się rodzinę okrę gów koncentrycznych, ze ś rodkiem w punkcie osobliwym, o równaniu (2.24) x2+y2 = с , с = const > 0. 142 A . P I E N I Ą Ż E K, W . P I E N I Ą Ż EK Po zróż niczkowaniu powyż szego r ó w n a n i a otrzymujemy dy x dx у Wobec powyż szego, krzywa stykowa jest opisana r ó w n a n i e m (2.25) £*Ż L = * У У lub po przekształceniu (2.26) xy+yf(x,y) = 0. Jak łatwo stwierdzić warunek ten m o ż na otrzymać przyrównując funkcję S'(x, y) do zera. Utwórzmy jeszcze stosunek [3] rd@ R(r, в ) (2.27) _ _ = _ ^ _ j . = ^ = c o n 8 t . Otrzymaliś my równanie izokliny wzglę dem wektora wodzą cego. K ą t
0, /?o > 0, y0 > 0.
Wprowadzimy transformację czasu i oznaczenia według w z o r ó w :
(3.2) x = co01, — = a, pow0 = P, ~ = y.
a>0 Wg
P o uwzglę dnieniu powyż szych zależ noś ci otrzymujemy równanie w postaci bezwymia
rowej
(3.3) х +х + у хъ = (