2 Bryła w 3D: tensor bezwładności, równania Eulera i oś pośrednia

0 Prolog: satelita, który zgubił Słońce

W sierpniu 1997 roku NASA wyniosła na orbitę małego satelitę obserwacyjnego Lewis. Po kilku dniach był stracony. Markley i Crassidis opisują przyczynę w pierwszym rozdziale swojej książki. W trybie awaryjnym satelita miał trzymać jedną swoją oś skierowaną na Słońce. Była to oś o pośrednim momencie bezwładności. Niewielka asymetria silniczków rozkręciła satelitę wokół tej osi, a obrót zaczął przelewać się na oś, którą sterowanie próbowało utrzymać. Silniczki pracowały coraz więcej, aż program je wyłączył. Wtedy obrót przeszedł na inną oś główną, satelita obrócił się o $90^\circ$ , panele ustawiły się prawie bokiem do Słońca i wyczerpały się baterie.

Ten sam efekt można zobaczyć w domu. Podrzuć książkę (zabezpieczoną gumką) tak, żeby wirowała wokół osi przechodzącej przez grzbiet, potem wokół osi prostopadłej do okładki, a na koniec wokół trzeciej osi, równoległej do okładki i prostopadłej do grzbietu. Dwa pierwsze obroty są spokojne. Trzeci prawie zawsze kończy się przewrotem. W 1985 roku kosmonauta Władimir Dżanibekow zauważył na stacji Salut 7, że odkręcona nakrętka motylkowa co kilka sekund sama przewraca się o $180^\circ$ i wraca. Nie działała na nią żadna siła ani moment.

Pytanie na start. Moment pędu swobodnej bryły jest stały, energia kinetyczna też. Jak to możliwe, że oś obrotu sama się przewraca? Odpowiedź w części 3.

1 Moment pędu w 3D: tensor bezwładności

Moment pędu bryły względem środka masy to suma $\mathbf r\times(\boldsymbol\omega\times\mathbf r)$ po wszystkich jej kawałkach, pomnożonych przez masę. Rozwijając podwójny iloczyn wektorowy, $\mathbf r\times(\boldsymbol\omega\times\mathbf r) = r^2\boldsymbol\omega - (\mathbf r\cdot\boldsymbol\omega)\mathbf r$ , i to jest wyrażenie liniowe w $\boldsymbol\omega$ . Każdą zależność liniową zapisujemy macierzą (lekcja 2.2 Przekształcenia liniowe: obroty, zmiana bazy, transformacja Lorentza):

Tensor bezwładności

\mathbf L = I\boldsymbol\omega,\qquad I = \int\left(r^2\,\mathbb 1 - \mathbf r\,\mathbf r^{\mathsf T}\right)dm = \begin{bmatrix} \int(y^2 + z^2)\,dm & -\int xy\,dm & -\int xz\,dm\\ -\int xy\,dm & \int(x^2 + z^2)\,dm & -\int yz\,dm\\ -\int xz\,dm & -\int yz\,dm & \int(x^2 + y^2)\,dm\end{bmatrix} .

Energia kinetyczna obrotu: $T = \tfrac12\boldsymbol\omega^{\mathsf T}I\boldsymbol\omega$ . Na przekątnej są znane momenty bezwładności względem osi, poza nią iloczyny dewiacji.

Macierz $I$ jest symetryczna, więc (lekcja 2.3 Wartości i wektory własne) ma trzy rzeczywiste wartości własne i trzy wzajemnie prostopadłe wektory własne. To są osie główne bryły i główne momenty bezwładności $I_1, I_2, I_3$ . Tylko przy obrocie wokół osi głównej $\mathbf L\parallel\boldsymbol\omega$ . Przy każdym innym $\mathbf L$ odchyla się od $\boldsymbol\omega$ w stronę osi o większym momencie.

Jedna rzecz decyduje o tym, w jakich osiach piszemy równania. W osiach związanych z bryłą macierz $I$ jest stała (masa nie przesuwa się względem bryły), a w osiach głównych diagonalna. W osiach inercjalnych zmienia się przy każdym obrocie: $I^n = R\,I^b R^{\mathsf T}$ . Dlatego równania ruchu obrotowego zawsze zapisuje się w osiach ciała, choć prawo Newtona (moment siły równa się pochodnej momentu pędu) obowiązuje w układzie inercjalnym.

Zadania · część 1

4 zadania · ok. 35 min

Hantla: dwie masy po 1 kg na końcach nieważkiego pręta o długości 1 m. Pręt leży w płaszczyźnie $xz$ pod kątem $30^\circ$ do osi $z$ , a hantla obraca się wokół osi $z$ z $\omega = 2$ rad/s. Policz $\mathbf L$ i kąt między $\mathbf L$ a $\boldsymbol\omega$ .

Odpowiedź
Masy w $\pm\mathbf r$ , $\mathbf r = 0{,}5\,[\sin 30^\circ;\ 0;\ \cos 30^\circ] = [0{,}25;\ 0;\ 0{,}433]$ . $\;\mathbf L = 2m\,\mathbf r\times(\boldsymbol\omega\times\mathbf r) = 2m\left(r^2\boldsymbol\omega - (\mathbf r\cdot\boldsymbol\omega)\mathbf r\right) = 2\left(0{,}25\cdot[0;0;2] - 0{,}866\,[0{,}25;\ 0;\ 0{,}433]\right) = [-0{,}433;\ 0;\ 0{,}25]$ kg·m²/s. Kąt z osią $z$ : $\arctan(0{,}433/0{,}25) = 60^\circ$ . $\mathbf L$ jest prostopadły do pręta i obraca się razem z nim, więc łożyska muszą dawać zmienny moment: to „bicie” niewyważonego koła.
Książka $24\times 17\times 3$ cm o masie 1 kg (jednorodny prostopadłościan, $I_x = \frac{m}{12}(b^2 + c^2)$ itd.). Policz główne momenty bezwładności. Która oś jest pośrednia?

Odpowiedź
$I_x = \frac{1}{12}(0{,}17^2 + 0{,}03^2) = 2{,}48\cdot 10^{-3}$ kg·m² (oś wzdłuż 24 cm), $I_y = \frac{1}{12}(0{,}24^2 + 0{,}03^2) = 4{,}87\cdot 10^{-3}$ (oś wzdłuż 17 cm), $I_z = \frac{1}{12}(0{,}24^2 + 0{,}17^2) = 7{,}21\cdot 10^{-3}$ (oś prostopadła do okładki). Pośrednia: oś równoległa do krótszego boku okładki.
Pokaż, że dla cienkiej płytki w płaszczyźnie $xy$ zachodzi $I_z = I_x + I_y$ . Sprawdź na książce, dlaczego równość nie jest dokładna.

Odpowiedź
Dla płytki $z \approx 0$ : $I_x = \int y^2dm$ , $I_y = \int x^2dm$ , $I_z = \int(x^2 + y^2)dm$ . Książka: $2{,}48 + 4{,}87 = 7{,}35$ wobec $7{,}21$ (w $10^{-3}$ kg·m²): różnica $2\cdot\frac{m}{12}c^2 = 0{,}15\cdot 10^{-3}$ pochodzi od grubości 3 cm.
Pokaż, że $T = \frac12\boldsymbol\omega^{\mathsf T}I\boldsymbol\omega \ge 0$ dla każdego $\boldsymbol\omega$ i że wszystkie główne momenty bezwładności są nieujemne. Dlaczego każdy z nich jest nie większy od sumy dwóch pozostałych?

Odpowiedź
$\boldsymbol\omega^{\mathsf T}I\boldsymbol\omega = \int\left(r^2\omega^2 - (\mathbf r\cdot\boldsymbol\omega)^2\right)dm = \int|\boldsymbol\omega\times\mathbf r|^2dm \ge 0$ . Dla wektora własnego daje to $I_k \ge 0$ . Nierówność trójkąta: $I_1 + I_2 - I_3 = 2\int z^2dm \ge 0$ w osiach głównych (i analogicznie dla pozostałych). Równość tylko dla płaskiej bryły.

2 Równania Eulera

Prawo dynamiki ruchu obrotowego z lekcji Dział 4: Ruch obrotowy i bryła sztywna: pochodna momentu pędu w układzie inercjalnym równa się momentowi sił zewnętrznych. Chcemy je zapisać w osiach ciała, gdzie $I$ jest stała. Wektor zapisany w osiach obracających się z prędkością $\boldsymbol\omega$ zmienia się z dwóch powodów: bo zmieniają się jego składowe i bo obracają się same osie. To ten sam rachunek co $\dot R = R[\boldsymbol\omega]_\times$ z lekcji 1 Orientacja w 3D: macierz obrotu, kąty Eulera, oś i kąt, kwaterniony: $\dot{\mathbf L}\big|_\text{inercjalnie} = \dot{\mathbf L}\big|_\text{ciało} + \boldsymbol\omega\times\mathbf L$ .

Równania Eulera (osie ciała, środek masy)

I\dot{\boldsymbol\omega} + \boldsymbol\omega\times I\boldsymbol\omega = \boldsymbol\tau .

W osiach głównych:

I_1\dot\omega_1 = (I_2 - I_3)\,\omega_2\omega_3 + \tau_1,\qquad I_2\dot\omega_2 = (I_3 - I_1)\,\omega_3\omega_1 + \tau_2,\qquad I_3\dot\omega_3 = (I_1 - I_2)\,\omega_1\omega_2 + \tau_3 .

Razem z kinematyką $\dot q = \frac12 q\otimes(0, \boldsymbol\omega)$ to pełny opis ruchu obrotowego bryły.

Człon $\boldsymbol\omega\times I\boldsymbol\omega$ jest nieliniowy i nazywa się żyroskopowym. Znika, gdy bryła obraca się wokół osi głównej ( $I\boldsymbol\omega\parallel\boldsymbol\omega$ ) albo gdy jest kulisto symetryczna. W pozostałych przypadkach obrót wokół jednej osi wywołuje przyspieszenie kątowe wokół innej, bez żadnego momentu siły. Nie jest to nowa siła: to skutek tego, że piszemy równania w osiach, które same się obracają.

Dla samolotu oznaczamy $\boldsymbol\omega = (p, q, r)$ , a osie ciała zwykle nie są osiami głównymi: przy symetrii lewo-prawo zostaje jeden iloczyn dewiacji $I_{xz}$ , bo płaszczyzna symetrii zawiera osie $x$ i $z$ . Równania Eulera z $I_{xz}$ zobaczymy w lekcji 5 Lot w płaszczyźnie symetrii: model podłużny i pochodne stateczności i w pełnej symulacji (lekcja 10 Symulacja 6 stopni swobody: składamy samolot).

Zadania · część 2

3 zadania · ok. 30 min

Bryła osiowo symetryczna, $I_1 = I_2$ , bez momentów. Co wynika z trzeciego równania Eulera?

Odpowiedź
$I_3\dot\omega_3 = (I_1 - I_2)\omega_1\omega_2 = 0$ , więc $\omega_3$ jest stałe. Prędkość obrotu wokół osi symetrii nie zmienia się, nawet gdy bryła się chwieje.
Dla tej samej bryły pokaż, że $(\omega_1, \omega_2)$ obraca się w osiach ciała ze stałą prędkością $\Omega = \frac{I_3 - I_1}{I_1}\,\omega_3$ .

Odpowiedź
$\dot\omega_1 = -\Omega\,\omega_2$ , $\;\dot\omega_2 = \Omega\,\omega_1$ , z $\Omega = (I_3 - I_1)\omega_3/I_1$ . To równania oscylatora: $\omega_1 + i\omega_2 = Ae^{i\Omega t}$ (lekcja 1.7 Liczby zespolone i fazory). Dla bryły spłaszczonej ( $I_3 > I_1$ ) wektor $\boldsymbol\omega$ obiega oś symetrii w tę samą stronę co obrót, dla wydłużonej w przeciwną.
Pokaż na równaniach Eulera z $\boldsymbol\tau = 0$ , że $|\mathbf L|^2 = \sum I_k^2\omega_k^2$ i $2T = \sum I_k\omega_k^2$ są stałe.

Odpowiedź
$\frac{d}{dt}\sum I_k\omega_k^2 = 2\sum\omega_k I_k\dot\omega_k = 2\omega_1\omega_2\omega_3\left[(I_2 - I_3) + (I_3 - I_1) + (I_1 - I_2)\right] = 0$ . $\;\frac{d}{dt}\sum I_k^2\omega_k^2 = 2\omega_1\omega_2\omega_3\left[I_1(I_2 - I_3) + I_2(I_3 - I_1) + I_3(I_1 - I_2)\right] = 0$ . Ogólniej: $\dot{\mathbf L}|_\text{ciało} = -\boldsymbol\omega\times\mathbf L \perp\mathbf L$ i $\dot T = \boldsymbol\omega\cdot I\dot{\boldsymbol\omega} = -\boldsymbol\omega\cdot(\boldsymbol\omega\times I\boldsymbol\omega) = 0$ .

3 Ruch swobodny i oś pośrednia

Ułóżmy osie tak, żeby $I_1 < I_2 < I_3$ . Czy obrót wokół osi $k$ ze stałą prędkością $\Omega$ jest stabilny? Linearyzujemy jak w lekcji 4.6 Układy nieliniowe, płaszczyzna fazowa i stabilność: $\boldsymbol\omega = \Omega\mathbf e_k + \delta\boldsymbol\omega$ i pomijamy iloczyny małych zaburzeń. Dla obrotu wokół osi 2:

I_1\delta\dot\omega_1 = (I_2 - I_3)\,\Omega\,\delta\omega_3,\qquad I_3\delta\dot\omega_3 = (I_1 - I_2)\,\Omega\,\delta\omega_1\quad\Longrightarrow\quad \delta\ddot\omega_1 = \underbrace{\frac{(I_3 - I_2)(I_2 - I_1)}{I_1I_3}}_{>0}\,\Omega^2\,\delta\omega_1 .

Stabilność obrotu wokół osi głównej

Zaburzenie obrotu wokół osi $k$ spełnia $\delta\ddot\omega = -\dfrac{(I_k - I_i)(I_k - I_j)}{I_iI_j}\,\Omega^2\,\delta\omega$ , gdzie $i, j$ to dwie pozostałe osie.

Oś najmniejszego i największego momentu: współczynnik dodatni, chybotanie z częstością $\Omega\sqrt{\frac{(I_k - I_i)(I_k - I_j)}{I_iI_j}}$ , stabilne.

Oś pośrednia: współczynnik ujemny, zaburzenie rośnie jak $e^{\lambda t}$ , $\;\lambda = \Omega\sqrt{\dfrac{(I_3 - I_2)(I_2 - I_1)}{I_1I_3}}$ . Niestabilne.

Odpowiedź na pytanie z prologu daje obraz geometryczny. W osiach ciała koniec wektora $\mathbf L$ musi leżeć na sferze o promieniu $|\mathbf L|$ (stały moment pędu) i jednocześnie na elipsoidzie $L_1^2/I_1 + L_2^2/I_2 + L_3^2/I_3 = 2T$ (stała energia). Porusza się więc po krzywej ich przecięcia, zwanej polhodią. Wokół osi 1 i 3 przecięcia są zamkniętymi pętlami. Przez oś 2 przechodzi krzywa szczególna, separatrysa: dwa wielkie koła przecinające się na osi 2. Ruch rozpoczęty minimalnie obok osi 2 biegnie wzdłuż separatrysy aż na drugą stronę sfery, czyli $\mathbf L$ w osiach ciała odwraca się o $180^\circ$ . W przestrzeni $\mathbf L$ stoi w miejscu, więc to bryła odwraca się względem $\mathbf L$ . Po czasie wraca: przewroty powtarzają się okresowo. To dokładnie nakrętka Dżanibekowa. Oś 2 jest punktem siodłowym w języku lekcji 4.6 Układy nieliniowe, płaszczyzna fazowa i stabilność.

Zadania · część 3

3 zadania · ok. 30 min

Książka z zadania 1.2 wiruje wokół osi pośredniej z $\Omega = 10$ rad/s. Ile wynosi $\lambda$ ? Po jakim czasie zaburzenie rośnie $e$ razy?

Odpowiedź
$\lambda = 10\sqrt{\frac{(7{,}21 - 4{,}87)(4{,}87 - 2{,}48)}{2{,}48\cdot 7{,}21}} = 10\sqrt{0{,}312} = 5{,}58$ 1/s. Czas $1/\lambda = 0{,}18$ s, czyli mniej niż jedna trzecia obrotu.
Wyprowadź warunek dla obrotu wokół osi 3 (największej) i pokaż, że chybotanie jest stabilne.

Odpowiedź
$I_1\delta\dot\omega_1 = (I_2 - I_3)\Omega\,\delta\omega_2$ , $I_2\delta\dot\omega_2 = (I_3 - I_1)\Omega\,\delta\omega_1$ , więc $\delta\ddot\omega_1 = \frac{(I_2 - I_3)(I_3 - I_1)}{I_1I_2}\Omega^2\delta\omega_1$ . Pierwszy czynnik ujemny, drugi dodatni: współczynnik ujemny, czyli oscylator. Książka: częstość $10\sqrt{0{,}911} = 9{,}5$ rad/s.
Rakieta tenisowa: która oś jest najmniejsza, która największa, a która pośrednia? Którym rzutem rakieta się przewraca?

Odpowiedź
Najmniejszy moment: wokół osi rączki (masa blisko niej). Największy: wokół osi prostopadłej do naciągu (cała masa daleko). Pośredni: wokół osi leżącej w płaszczyźnie naciągu, prostopadłej do rączki. Rzut „koziołkiem” wokół tej osi kończy się dodatkowym półobrotem wokół rączki: tzw. twierdzenie o rakiecie tenisowej.

4 Rozpraszanie energii i reguła osi największej

Bryła doskonale sztywna stabilnie wiruje zarówno wokół osi największego, jak i najmniejszego momentu. Prawdziwy satelita nie jest sztywny: anteny drgają, paliwo się przelewa, konstrukcja się ugina. Siły wewnętrzne nie zmieniają momentu pędu (działają parami), ale zamieniają energię obrotu w ciepło. Przy stałym $|\mathbf L|$ energia czystego obrotu wokół osi $k$ wynosi

T = \frac{L^2}{2I_k},

więc jest największa dla osi najmniejszego momentu, a najmniejsza dla osi największego. Bryła, która traci energię przy stałym $\mathbf L$ , musi więc w końcu wirować wokół osi największego momentu. To reguła osi największej (ang. major axis rule).

Explorer 1, pierwszy amerykański satelita (1958), był wydłużonym walcem z czterema giętkimi antenami prętowymi. Po wejściu na orbitę wirował z prędkością 750 obr/min wokół osi podłużnej, czyli osi najmniejszego momentu. Drgania anten rozpraszały energię, kąt chybotania (nutacji) rósł, aż satelita koziołkował w „płaskim obrocie” wokół osi poprzecznej. Moment pędu przez cały czas był ten sam. Satelity wirujące projektuje się od tamtej pory jako spłaszczone (oś wirowania jest osią największego momentu) albo dwuspinowe (zadanie 4.2), a spłaszczonym celowo dodaje się tłumiki nutacji: rurki z cieczą lub masą na sprężynie, które rozpraszają energię chybotania i stabilizują obrót.

Zadania · część 4

2 zadania · ok. 20 min

Explorer 1 (Kluever, przykład 12.4): $I_\text{osiowy} = 0{,}17$ kg·m², $I_\text{poprzeczny} = 5$ kg·m², początkowo 750 obr/min wokół osi podłużnej. Policz energię początkową i końcową oraz prędkość końcowego płaskiego obrotu.

Odpowiedź
$\omega_0 = 78{,}5$ rad/s, $L = 0{,}17\cdot 78{,}5 = 13{,}35$ kg·m²/s. $T_0 = L^2/(2\cdot 0{,}17) = 524$ J, $T_k = L^2/(2\cdot 5) = 17{,}8$ J, $\omega_k = L/5 = 2{,}67$ rad/s $= 25{,}5$ obr/min. Energia spadła 29 razy, tyle co stosunek momentów.
Dlaczego tłumik nutacji pomaga satelicie spłaszczonemu, a szkodziłby wydłużonemu? Co zrobić, jeśli satelita musi być wydłużony (np. ze względu na owiewkę rakiety)?

Odpowiedź
Tłumik rozprasza energię; przy stałym $L$ to zawsze przybliża bryłę do obrotu wokół osi największego momentu. Dla satelity spłaszczonego to jest oś wirowania (chybotanie maleje), dla wydłużonego oś poprzeczna (chybotanie rośnie). Satelita wydłużony potrzebuje aktywnego sterowania albo konstrukcji dwuspinowej: wirnik w środku daje moment pędu, a reszta wiruje wolno; rozwiązanie z kursu astronautyki (lekcja 11 Satelita bez sterowania: wirowanie, dwuspin i gradient grawitacyjny).

5 Sprzężenie żyroskopowe: śmigło i wirniki

Samolot ze śmigłem to bryła z dodatkowym wirującym elementem o momencie pędu $\mathbf h$ (stałym w osiach ciała, jeśli obroty silnika są stałe). Całkowity moment pędu to $I\boldsymbol\omega + \mathbf h$ i równania Eulera dostają jeden człon więcej:

Bryła z wirnikiem

I\dot{\boldsymbol\omega} + \boldsymbol\omega\times(I\boldsymbol\omega + \mathbf h) = \boldsymbol\tau .

Człon $\boldsymbol\omega\times\mathbf h$ działa jak dodatkowy moment $-\boldsymbol\omega\times\mathbf h$ , prostopadły zarówno do osi wirnika, jak i do osi obrotu kadłuba.

Śmigło obracające się zgodnie z ruchem wskazówek zegara, patrząc z kabiny, ma $\mathbf h = (h, 0, 0)$ wzdłuż osi $x$ (do przodu). Gdy samolot pochyla się z prędkością $q$ , $\;-\boldsymbol\omega\times\mathbf h = -(0, q, 0)\times(h, 0, 0) = (0, 0, qh)$ : moment odchylający. Podnoszenie nosa odchyla samolot w prawo, opuszczanie w lewo. Tę drugą sytuację zna każdy pilot samolotu z kółkiem ogonowym: przy podnoszeniu ogona na rozbiegu samolot ucieka w lewo (obok innych przyczyn, jak strumień zaśmigłowy).

Ten sam człon tłumaczy, dlaczego wielowirnikowce mają wirniki obracające się parami w przeciwne strony: suma $\mathbf h$ wynosi zero, więc sprzężenie żyroskopowe znika (a przy okazji znosi się moment oporu wirników, który inaczej obracałby kadłubem). Satelity robią odwrotnie: celowo montują wirnik o dużym $\mathbf h$ (koło zamachowe), żeby usztywnić orientację (lekcja 13 Sterowanie orientacją: koła reakcyjne, magnetyki i silniczki).

Zadania · część 5

3 zadania · ok. 25 min

Przyjmij śmigło o momencie bezwładności 3 kg·m² przy 2400 obr/min (liczby przyjęte dla przykładu). Samolot pochyla się z $q = 0{,}5$ rad/s. Jaki moment odchylający daje sprzężenie żyroskopowe i jakie przyspieszenie kątowe odchylania przy $I_z = 2500$ kg·m²?

Odpowiedź
$h = 3\cdot 2400\cdot 2\pi/60 = 754$ N·m·s; moment $qh = 377$ N·m; $\dot r = 377/2500 = 0{,}15$ rad/s². Po sekundzie gwałtownego pochylania samolot odchyla się z prędkością rzędu $0{,}15$ rad/s, czyli $9^\circ$ /s: wyraźnie odczuwalne, łatwe do skorygowania sterem kierunku.
Samolot z tym samym śmigłem wykonuje zakręt w prawo z prędkością odchylania $r = 0{,}2$ rad/s (bez pochylania). Jaki moment żyroskopowy powstaje i wokół której osi?

Odpowiedź
$-\boldsymbol\omega\times\mathbf h = -(0, 0, r)\times(h, 0, 0) = (0, -rh, 0)$ : moment pochylający $-rh = -151$ N·m, czyli w dół (nos w dół). Przy zakręcie w lewo nos idzie w górę. Pochylanie sprzęga się z odchylaniem w obie strony.
Dlaczego w modelu samolotu z lekcji 5 Lot w płaszczyźnie symetrii: model podłużny i pochodne stateczności zwykle pomija się $\mathbf h$ , a w modelu śmigłowca nie można?

Odpowiedź
W samolocie lekkim moment żyroskopowy jest mały w porównaniu z momentami aerodynamicznymi od usterzenia (zadanie 5.1: setki N·m wobec tysięcy), więc liniowe modele go pomijają lub traktują jako zakłócenie. W śmigłowcu wirnik nośny ma moment pędu porównywalny z momentem pędu całego kadłuba i właśnie przez niego sterujemy (przechylanie płaszczyzny wirnika), więc sprzężenie jest częścią dynamiki.

6 Metoda i przykłady rozwiązane

Metoda: dynamika obrotowa bryły

Tensor: policz $I$ w osiach ciała (z tablic i twierdzenia Steinera); znajdź osie główne jako wektory własne.
Równania: $I\dot{\boldsymbol\omega} + \boldsymbol\omega\times(I\boldsymbol\omega + \mathbf h) = \boldsymbol\tau$ w osiach ciała, kinematyka na kwaternionie.
Stabilność: linearyzacja wokół obrotu ustalonego; oś pośrednia jest siodłem.
Bryła giętka: stałe $\mathbf L$ , malejące $T$ ; dąży do osi największego momentu.
Sprawdzenie: bez momentów $|\mathbf L|$ i $T$ muszą być stałe w symulacji (dobry test integratora).

Przykład 6.1 Ile obrotów do przewrotu książki

Książka z zadania 1.2 wiruje wokół osi pośredniej z $\Omega = 10$ rad/s, a początkowe zaburzenie na każdej z pozostałych osi to 1% prędkości (wektor $\boldsymbol\omega$ odchylony od osi o ok. $0{,}8^\circ$ ). Kiedy nastąpi przewrót?

$\lambda = 5{,}58$ 1/s (zadanie 3.1).
Zaburzenie musi urosnąć z 0,01 do rzędu 1, czyli $e^{\lambda t} \approx 100$ : $t \approx \ln 100/5{,}58 = 0{,}83$ s.
W tym czasie książka wykonuje $10\cdot 0{,}83/(2\pi) \approx 1{,}3$ obrotu.
Symulacja nieliniowa (symulacja 1): pierwsza zmiana znaku $\omega_y$ przy $t = 1{,}14$ s, kolejne co $2{,}03$ s.

Wynik: przewrót po około sekundzie, zanim książka spadnie na ziemię z ręki.

Sprawdzenie: zmniejszenie zaburzenia dziesięciokrotnie opóźnia przewrót tylko o $\ln 10/\lambda = 0{,}41$ s: dlatego „idealny” rzut nie pomaga ✓.

Przykład 6.2 Explorer 1: dokąd prowadzi rozpraszanie

Dane z zadania 4.1. Pokaż, że obrót wokół osi podłużnej przy stałym $L$ ma więcej energii niż każdy inny ruch i że stan końcowy to płaski obrót.

Dla bryły osiowo symetrycznej z kątem nutacji $\theta$ (między osią symetrii a $\mathbf L$ ): $T = \frac{L^2}{2}\left(\frac{\sin^2\theta}{I_t} + \frac{\cos^2\theta}{I_a}\right)$ , gdzie $I_a = 0{,}17$ , $I_t = 5$ .
$\frac{dT}{d\theta} = \frac{L^2}{2}\sin 2\theta\left(\frac{1}{I_t} - \frac{1}{I_a}\right) < 0$ dla $0 < \theta < 90^\circ$ , bo $I_a < I_t$ .
Malejąca energia wymusza więc rosnące $\theta$ , aż do $90^\circ$ , gdzie $T = L^2/(2I_t) = 17{,}8$ J.

Wynik: nutacja rośnie monotonicznie do $90^\circ$ ; płaski obrót z 25,5 obr/min.

Sprawdzenie: dla bryły spłaszczonej ( $I_a > I_t$ ) znak pochodnej się odwraca i rozpraszanie zmniejsza $\theta$ (symulacja 3) ✓.

Przykład 6.3 Równania ruchu obrotowego samolotu

Zapisz równania Eulera dla samolotu z płaszczyzną symetrii $xz$ (tensor z $I_{xz} \ne 0$ ) i śmigłem $\mathbf h = (h, 0, 0)$ , z momentami aerodynamicznymi $(\mathcal L, \mathcal M, \mathcal N)$ .

$I = \begin{bmatrix} I_x & 0 & -I_{xz}\\ 0 & I_y & 0\\ -I_{xz} & 0 & I_z\end{bmatrix}$ , $\;I\boldsymbol\omega + \mathbf h = (I_xp - I_{xz}r + h,\ I_yq,\ I_zr - I_{xz}p)$ .
$\boldsymbol\omega\times(I\boldsymbol\omega + \mathbf h)$ , składowa $y$ : $r(I_xp - I_{xz}r + h) - p(I_zr - I_{xz}p) = (I_x - I_z)pr + I_{xz}(p^2 - r^2) + rh$ .
Pochylanie: $I_y\dot q = \mathcal M - (I_x - I_z)pr - I_{xz}(p^2 - r^2) - rh$ . Analogicznie dla przechylania i odchylania, które są sprzężone przez $I_{xz}$ : $I_x\dot p - I_{xz}\dot r = \mathcal L + \dots$ , $\;I_z\dot r - I_{xz}\dot p = \mathcal N + \dots$

Wynik: trzy równania z członami żyroskopowymi kadłuba ( $pr$ , $p^2 - r^2$ , $pq$ , $qr$ ) i śmigła ( $rh$ , $qh$ ). W locie prostym ( $p, q, r$ małe) człony kwadratowe znikają po linearyzacji, a przy pominięciu $h$ ruch podłużny ( $q$ ) oddziela się od bocznego ( $p, r$ ). Człony śmigła $rh$ i $qh$ są liniowe i sprzęgają oba ruchy.

Sprawdzenie: dla $h = 0$ , $I_{xz} = 0$ wracamy do równań Eulera w osiach głównych; człon $-rh$ zgadza się z zadaniem 5.2 ✓.

7 Błędne intuicje, których trzeba się oduczyć

„Moment pędu zawsze jest równoległy do prędkości kątowej.”: Tylko przy obrocie wokół osi głównej. W ogólności $\mathbf L = I\boldsymbol\omega$ i oba wektory mogą tworzyć duży kąt (hantla: $60^\circ$ ).
„Bez momentu siły prędkość kątowa jest stała.”: Stały jest moment pędu w przestrzeni. Prędkość kątowa w osiach ciała zmienia się zgodnie z $I\dot{\boldsymbol\omega} = -\boldsymbol\omega\times I\boldsymbol\omega$ , czasem gwałtownie (nakrętka Dżanibekowa).
„Przewrót przy osi pośredniej to efekt niedokładnego rzutu.”: Niedokładność tylko go uruchamia. Zaburzenie rośnie wykładniczo, więc dziesięć razy staranniejszy rzut opóźnia przewrót o niecałe pół sekundy.
„Obrót wokół osi najmniejszego momentu jest stabilny, więc satelita może tak wirować.”: Dla bryły sztywnej tak. Każda giętkość rozprasza energię i przy stałym $\mathbf L$ przenosi obrót na oś największego momentu.
„Siła żyroskopowa to nowa siła.”: To skutek zapisu równań w osiach obracających się razem z bryłą. W układzie inercjalnym jest tylko $\dot{\mathbf L} = \boldsymbol\tau$ .

8 Zadania do przemyślenia

Te zadania wymagają więcej myślenia niż liczenia.

Zadania do przemyślenia

5 zadań · ok. 1,5 godziny

Chybotanie Ziemi. Ziemia jest spłaszczona: $(C - A)/A \approx 1/305$ , gdzie $C$ to moment względem osi biegunowej, $A$ względem osi w płaszczyźnie równika. Jeśli oś obrotu nie pokrywa się dokładnie z osią figury, jaki okres chybotania przewiduje zadanie 2.2? Obserwowany okres (chybotanie Chandlera) to ok. 433 dni. Skąd różnica?
Rozwiązanie
1. $\Omega = \frac{C - A}{A}\omega_3$ , więc okres w osiach Ziemi to $305$ dób: przewidywanie Eulera.
2. Ziemia nie jest sztywna: płaszcz i oceany odkształcają się pod wpływem zmieniającej się siły odśrodkowej, co częściowo „dopasowuje” wybrzuszenie do chwilowej osi i wydłuża okres. Chandler zmierzył go w 1891 roku. To przykład, w którym model bryły sztywnej daje poprawny rząd wielkości, ale nie dokładną liczbę.
Lewis. Wyjaśnij na polhodiach, dlaczego obrót satelity Lewis wokół osi pośredniej „przelał się” na inną oś i dlaczego silniczki sterujące jedną osią nie mogły tego powstrzymać.
Rozwiązanie
1. Obrót wokół osi pośredniej leży na separatrysie: najmniejsze zaburzenie przenosi $\mathbf L$ w osiach ciała daleko od osi 2, w stronę osi o skrajnych momentach. Energia i moment pędu się nie zmieniają, zmienia się tylko oś, wokół której satelita wiruje.
2. Według raportu komisji badającej awarię obrotu wokół tej osi nikt nie mierzył: jedyny czujnik prędkości kątowej w trybie awaryjnym był dwuosiowy i tej osi nie obejmował, a sterowanie próbowało utrzymać inną oś. Silniczki walczyły z przepływem momentu pędu między osiami, zużywały paliwo, aż zabezpieczenie je wyłączyło. Właściwe rozwiązanie: wybrać do trybu awaryjnego oś największego momentu.
Bąk i żyroskop. W lekcji Dział 4: Ruch obrotowy i bryła sztywna żyroskop z ciężarkiem nie przewraca się, tylko krąży (precesja). Wyprowadź prędkość precesji z $\dot{\mathbf L} = \boldsymbol\tau$ dla szybkiego wirnika i porównaj z członem $-\boldsymbol\omega\times\mathbf h$ z części 5.
Rozwiązanie
1. Dla szybkiego wirnika $\mathbf L\approx\mathbf h$ , poziomo. Moment ciężaru $\tau = mgd$ jest poziomy i prostopadły do $\mathbf h$ , więc obraca $\mathbf h$ w płaszczyźnie poziomej z prędkością $\Omega_p = mgd/h$ .
2. W osiach związanych z ramą żyroskopu: obrót ramy z $\boldsymbol\Omega_p$ (pionowo) daje człon $\boldsymbol\Omega_p\times\mathbf h$ , który musi być równy momentowi ciężaru: $\Omega_ph = mgd$ . Ten sam rachunek co moment żyroskopowy śmigła.
Przy równych momentach. Co się dzieje z niestabilnością osi pośredniej, gdy $I_2 \to I_1$ albo $I_2 \to I_3$ ? Co to znaczy dla brył prawie symetrycznych?
Rozwiązanie
1. $\lambda\to 0$ : niestabilność słabnie, a separatrysa przylega do okręgu wokół osi o równych momentach. Przy $I_1 = I_2$ oś 2 przestaje być wyróżniona, a każda oś w płaszczyźnie 1–2 jest osią główną.
2. Bryła prawie symetryczna ma bardzo wolny przewrót, łatwo zasłonięty przez rozpraszanie energii i momenty zakłócające. Klucz T i książka mają momenty wyraźnie różne, dlatego efekt jest szybki i widoczny.
Test integratora. W symulacji 1 zmiana energii po 8 s jest rzędu $10^{-6}$ %. Zaproponuj, jak to wykorzystać do sprawdzania kodu symulacji samolotu, w którym momenty sił nie są zerowe.
Rozwiązanie
1. Wyłącz siły i momenty zewnętrzne (lub ustaw je na zero na chwilę) i sprawdź, czy $|\mathbf L|$ (w układzie inercjalnym: wektor $R\,I\boldsymbol\omega$ ) i $T$ pozostają stałe. Dryf wskazuje błąd w równaniach (np. zły znak członu żyroskopowego) albo za duży krok.
2. Z momentami: sprawdź bilans $\Delta T = \int\boldsymbol\tau\cdot\boldsymbol\omega\,dt$ i $\Delta\mathbf L = \int\boldsymbol\tau^n\,dt$ . To tani test, który wyłapuje większość pomyłek w lekcji 10 Symulacja 6 stopni swobody: składamy samolot.

9 Sprawdź się

Odpowiedz w myślach, zanim rozwiniesz odpowiedź.

1 Kiedy

\mathbf L\parallel\boldsymbol\omega

Gdy bryła obraca się wokół osi głównej, czyli wektora własnego tensora bezwładności.

2 Dlaczego równania obrotu pisze się w osiach ciała?

Tam tensor bezwładności jest stały (w osiach głównych diagonalny). Ceną jest człon $\boldsymbol\omega\times I\boldsymbol\omega$ .

3 Jak brzmią równania Eulera?

$I\dot{\boldsymbol\omega} + \boldsymbol\omega\times I\boldsymbol\omega = \boldsymbol\tau$ ; w osiach głównych $I_1\dot\omega_1 = (I_2 - I_3)\omega_2\omega_3 + \tau_1$ i cyklicznie.

4 Wokół której osi obrót jest niestabilny i jak szybko rośnie zaburzenie?

Wokół osi pośredniej; $\lambda = \Omega\sqrt{(I_3 - I_2)(I_2 - I_1)/(I_1I_3)}$ .

5 Co mówi reguła osi największej?

Bryła rozpraszająca energię przy stałym $\mathbf L$ kończy w obrocie wokół osi największego momentu, bo tam $T = L^2/(2I)$ jest najmniejsza.

6 Jaki moment daje wirnik o momencie pędu

\mathbf h

w obracającym się kadłubie?

$-\boldsymbol\omega\times\mathbf h$ ; przy śmigle pochylanie daje odchylanie i odwrotnie.

10 Karta wzorów i co dalej

Tensor i energia

$\mathbf L = I\boldsymbol\omega$ , $I = \int(r^2\mathbb 1 - \mathbf r\mathbf r^{\mathsf T})dm$

$T = \frac12\boldsymbol\omega^{\mathsf T}I\boldsymbol\omega$ ; osie główne = wektory własne $I$

$I^n = RI^bR^{\mathsf T}$

prostopadłościan: $I_x = \frac m{12}(b^2 + c^2)$

Dynamika

$I\dot{\boldsymbol\omega} + \boldsymbol\omega\times(I\boldsymbol\omega + \mathbf h) = \boldsymbol\tau$

oś pośrednia: $\lambda = \Omega\sqrt{\frac{(I_3 - I_2)(I_2 - I_1)}{I_1I_3}}$

symetryczna: $\Omega_n = \frac{I_3 - I_1}{I_1}\omega_3$

stałe $L$ : $T = L^2/(2I_k)$ , minimum dla $I_\text{max}$

Do przeczytania

Obowiązkowo: C. A. Kluever, Space Flight Dynamics, rozdz. 12.2–12.4 (tensor, równania Eulera, stabilność obrotu, Explorer 1 w przykładzie 12.4).
Dla chętnych: H. Schaub, J. L. Junkins, Analytical Mechanics of Aerospace Systems, rozdz. 4 (dynamika bryły, polhodie, bryła z wirnikiem); F. L. Markley, J. L. Crassidis, rozdz. 1 (historia satelity Lewis) i 3.3.

Mamy już orientację i dynamikę bryły. Zanim zajmiemy się siłami aerodynamicznymi, wrócimy do problemu z automatyki: skąd kontroler lotu wie, jak samolot jest zorientowany? W lekcji 4.5 automatyki łączyliśmy żyroskop i akcelerometr w jednej osi. Następna lekcja robi to w trzech osiach, na kwaternionach, z magnetometrem do kursu.