FIZYKA KOMPUTEROWA SPRAWOZDANIE Z ZADANIA 4. – RUCH W POLU GRAWITACYJNYM GWIAZDY PODWÓJNEJ

FIZYKA KOMPUTEROWA

SPRAWOZDANIE Z ZADANIA 4. – RUCH W POLU GRAWITACYJNYM GWIAZDY PODWÓJNEJ

Krystian Stasiak, nr albumu: 151007 Maciej Malinowski, nr albumu: 154719

Treść zadania:

Dwie gwiazdy o masach M1 i M2 krążą wokół siebie w odległości D. Należy przeanalizować ruch niewielkiej planety w takim układzie (dla uproszczenia przyjąć, że porusza się ona w płaszczyźnie obu gwiazd). W układzie odniesienia obracającym się razem z gwiazdami (prędkość kątową Ω należy samemu wyliczyć) równania ruchu mają postać:

d²x

d x²=−G M₁(x−X₁) [

(

^x−X1

)

^{2+y2]3 /2−}

G M₂(x− X₂)

[

(

^x−X2

)

^{2+y2]3 /2−2 Ω}

dy dt+Ω²x d²y

d y²= −G M₁y [

(

^x−X1

)

^{2+y2]3 /2−}

G M₂y

[

(

^{x− X}2

)

^{2+y2]3/ 2+2 Ω}

dx dt+Ω²y

gdzie G to uniwersalna stała grawitacyjna. Ponadto przyjęto, że w obracającym się układzie odniesienia gwiazdy znajdują się na osi 0X. Ich odległości X1 oraz X2 od ich wspólnego środka masy (który znajduje się w środku układu współrzędnych) należy samemu wyliczyć.

W rozwiązaniu trzeba pokazać trajektorię planety zarówno w układzie obracającym się jak i nieruchomym. Wykreślić zależności czasowe odległości planety od środka układu jak i kąta w stosunku do osi 0X. W przypadku periodycznej orbity określić czas obiegu. Należy

przeanalizować różne jakościowo przypadki – zarówno różne początkowe położenia i prędkości planety, jak i różne proporcje mas M1/ M2. Wśród analizowanych przypadków powinny się znaleźć tzw. punkty Lagrange’a.

Założenia:

Przyjmujemy, że pierwsza gwiazda ma ujemną współrzędną początkową x (X1), a druga – dodatnią. Środek masy znajduje się między gwiazdami, więc równie dobrze można by przyjąć na odwrót, bez jakościowych zmian.

Wykonanie:

Rozważone zostały wszystkie punkty Lagrange’a (L1, L2, L3, L4, L5) oraz 4 inne warte uwagi przypadki. Zauważmy, współrzędne początkowe gwiazd - X1 oraz X2 – spełniają zależności:

X₂−X₁=D M₁X₁=M₂X₂

Zatem:

X₁= −D M₂ (M₁+M₂) X₂= D M₁

(M₁+M₂)

Z warunku równowagi siły odśrodkowej z siłą przyciągania dla dowolnej gwiazdy (aby uniknąć problemów ze znakami, wybieramy drugą), podstawiając prędkość liniową jako iloczyn prędkości obrotowej i promienia, otrzymujemy:

M2(Ω X2)²

X₂ =GM1M2

D²

Co po przekształceniach, oraz rozwinięciu X2 według wyznaczonego powyżej wzoru, daje wynik:

Ω=

√

^{G(MD1+3 M2)}

Przy rozwiązywaniu zadania potrzebne było także zapisanie relacji pomiędzy (traktowanymi z odpowiednimi znakami) siłami przyciągania przez gwiazdy oraz siłą bezwładności planety, które muszą być spełnione w punktach Lagrange’a L1, L2 oraz L3. Współrzędne tych punktów są wyznaczane numerycznie za pomocą procedury fsolve.

Do rozwiązywania podanego układu równań różniczkowych wykorzystano oczywiście procedurę ode. Dla każdego z przypadków rozważono 1000 równomiernie rozłożonych w czasie punktów. Czas mierzony w obrotach układu gwiazd jest elementem każdego z zestawów danych (dobieranym tak, żeby wykresy były czytelne).

Procedura ode pozwala uzyskać informację o współrzędnych x oraz y planety w czasie.

Obróbka tych danych polega na przygotowaniu współrzędnych w układzie obróconym o kąt Ωt, (za pomocą macierzy M=

[

cos ⁡(Ωt ) −sin ⁡(Ωt )

sin ⁡(Ωt ) cos ⁡(Ωt )

]

^przez

^[

^xy

^]

dla każdej chwili czasu t), a także wektora odległości od środka masy układu w czasie (za pomocą procedury norm) oraz kąta nachylenia względem osi OX (za pomocą dwuargumentowego wywołania

procedury atan).

Wskazane wykresy dla każdego z rozważanych przypadków zostały sporządzone za pomocą procedury plot2d, która jest właściwa także do rysowania krzywych na płaszczyźnie, takich

jak trajektoria planety (która porusza się wyłącznie w stałej płaszczyźnie, w której obraca się układ).

W celu przygotowania trajektorii gwiazd w czasie, wystarczyło zapisać ruch po okręgu jako złożenie faz. Dla pierwszej gwiazdy otrzymujemy x (t )= X1cos ⁡(Ωt ) , y (t)= X1sin ⁡(Ωt ) ; dla drugiej – analogicznie, tyle że z X2 zamiast X1. Przeciwne znaki X2 oraz X1 gwarantują, że gwiazdy zawsze będą dokładnie po przeciwnych stronach środka masy układu.

Dobrym pomysłem usprawniającym testowanie różnych zestawów danych było

przygotowanie funkcji generujących całe zestawy – dzięki temu zamiast komentowania wielu linijek kodu w scilabie (co zapewne wprowadziłoby duży bałagan w kodzie), wystarczy w jednym miejscu podmieniać nazwę wykorzystanej funkcji.

Punkty Lagrange’a

Dla następującego zestawu parametrów:

M1 2 * 10³⁰ [kg]

M2 10³⁰ [kg]

D 3 * 10¹¹ [m]

Z których wyniknęło, że:

X1 -1 * 10¹¹ [m]

X2 2 * 10¹¹ [m]

Ω 8,6106 * 10^-8 [Hz]

Rozważone zostały wszystkie punkty Lagrange’a. Dzięki możliwościom programu scilab, nie są potrzebne żadne przybliżenia wynikające z różnicy masy dwóch ciał (takie przybliżenia stosuje się np. do wyznaczania punktów Lagrange’a dla układu Ziemia-Słońce bądź Ziemia-Księżyc) – współrzędne punktów wyznaczane są numerycznie za pomocą fsolve. Potwierdzeniem jakości tych obliczeń są zamieszczone poniżej wykresy – planeta umieszczona w którymkolwiek z wyznaczonych punktów nie poruszy się względem obracającego się układu gwiazd nawet o metr (oczywiście jest to przypadek bardzo wyidealizowany).

Wyznaczony punkt L1 ma współrzędne początkowe [7,1225 * 10¹⁰ m; 0]. Wymienione w treści wykresy dla tego punktu są poniżej. Na pierwszych dwóch wykresach zaznaczono trajektorie oraz położenia początkowe ciał (gwiazd na czerwono, planety na niebiesko).

Gwiazdy można rozróżnić na podstawie położeń początkowych (pierwsza jest po lewej).

Wyznaczony punkt L2 ma współrzędne początkowe [3,7471 * 10¹⁰ m; 0]. Dla tego punktu uzyskano następujące wykresy (wykres odległości od środka masy w czasie pominięto, gdyż jest stałą i wynika to już z pierwszego wykresu):

Wyznaczony punkt L3 ma współrzędne początkowe [-3,4091 * 10¹⁰ m; 0]. Dla tego punktu uzyskano następujące wykresy:

Punkty L4 = [5 * 10¹⁰ m; 2,5981 * 10¹¹ m] oraz L5 = [5 * 10¹⁰ m; -2,5981 * 10¹¹ m] są symetryczne do siebie względem osi OX. Dlatego nie warto zamieszczać niemalże identycznych wykresów dla obu z nich. Dla L4 wykresy mają postać:

Przypadek I

Planeta znajdująca się początkowo wewnątrz układu gwiezdnego może zostać ściągnięta przez pole grawitacyjne jednej z gwiazd i dalej po prostu krążyć wokół tej gwiazdy. Stanie się tak np. w przypadku, gdy masy gwiazd będą różne od siebie, a planeta będzie się znajdywała w środku masy układu. Aby orbita miała atrakcyjniejszy kształt, rozważamy planetę z

prędkością początkową skierowaną w górę. Zacznie ona krążyć wokół cięższej gwiazdy, po orbicie o ciekawym kształcie.

Przyjęte parametry:

M1 2 * 10³⁰ [kg]

M2 10³⁰ [kg]

D 3 * 10¹¹ [m]

x0 0

y0 0

Vx0 0

Vy0 30 000 [m/s]

l. obrotów 0,97 Wielkości wyznaczone:

X1 -1 * 10¹¹ [m]

X2 2 * 10¹¹ [m]

Ω 8,6106 * 10^-8 [Hz]

Wykresy:

W tym przypadku mamy do czynienia z orbitą periodyczną. Jest ona powtarzalna z okresem równym czasowi jednego pełnego obrotu układu gwiezdnego, czyli T =2 π

Ω wynoszącego około 1448 dni. Na taki a nie inny kształt trajektorii istotny wpływ ma oczywiście prędkość początkowa oraz stosunek mas gwiazd – próby dla innych wartości można łatwo

przeprowadzać zmieniając nieznacznie kod programu, bądź dodając nowe funkcje generujące zestawy danych.

Przypadek II

Jeżeli znajdująca się wewnątrz układu planeta będzie miała odpowiednio dużą prędkość, zacznie poruszać się po ciekawym torze wokół obu gwiazd. Wówczas jednak, po dłuższym czasie, doprowadzi to do rozpędzenia się powyżej prędkości ucieczki i opuszczenia układu (oddalania się od niego po spirali).

Przyjęte parametry:

M1 2 * 10³⁰ [kg]

M2 10³⁰ [kg]

D 3 * 10¹¹ [m]

x0 0

y0 0

Vx0 0

Vy0 60 000 [m/s]

Wielkości wyznaczone:

X1 -1 * 10¹¹ [m]

X2 2 * 10¹¹ [m]

Ω 8,6106 * 10^-8 [Hz]

Wykresy dla początkowych 4,5 obrotów:

Dla tego przypadku wykonaliśmy dodatkowe wykresy, przedstawiające przebieg zjawisk dla 13 pełnych obrotów układu – są one mniej czytelne, ale pozwalają stwierdzić ucieczkę planety (po dość długim ruchu w pobliżu gwiazd):

Jak widać, czasem nie wystarczy popatrzeć na trajektorie z początkowych faz ruchu, aby rozstrzygnąć dalsze zachowanie planety (trudno byłoby przewidzieć ucieczkę planety patrząc tylko na pierwszą serię wykresów).

Przypadek III

Planeta znajdująca się nieznacznie bliżej środka masy układu, niż punkt L2 lub L3, zostanie przyciągnięta i zacznie poruszać się po skomplikowanej orbicie wokół obu gwiazd. Nawet niewielka zmiana (w tym przypadku – o jeden metr) robi różnicę.

Przyjęte parametry:

M1 2 * 10³⁰ [kg]

M2 5 * 10²⁹ [kg]

D 3 * 10¹¹ [m]

x0 współrzędna x punktu L2 pomniejszona o 1m

y0 0

Vx0 0

Vy0 0

l. obrotów 5

Wielkości wyznaczone:

X1 -6 * 10¹⁰ [m]

X2 2,4 * 10¹¹ [m]

Ω 7,8604 * 10^-8 [Hz]

Wykresy:

Przez pewien czas planeta niemalże wyłącznie obraca się z układem, ale bardzo powoli zbliża się do środka masy. Jednak w pewnym momencie następuje nagłe, gwałtowne zbliżanie się do środka masy układu, a następnie ruch po skomplikowanej orbicie wokół obu gwiazd.

Przypadek IV

Planeta znajdująca się nieznacznie dalej od środka masy układu, niż punkt L2 lub L3, zostanie zacznie się oddalać od układu. Także w tym przypadku nawet niewielka zmiana (w tym przypadku – o jeden metr) powoduje zajście zjawiska.

Przyjęte parametry:

M1 2 * 10³⁰ [kg]

M2 5 * 10²⁹ [kg]

D 3 * 10¹¹ [m]

x0 współrzędna x punktu L2 zwiększona o 1m

y0 0

Vx0 0

Vy0 0

l. obrotów 4,5

Wielkości wyznaczone:

X1 -6 * 10¹⁰ [m]

X2 2,4 * 10¹¹ [m]

Ω 7,8604 * 10^-8 [Hz]

Wykresy:

Także w tym przypadku, przez pewien czas planeta niemalże wyłącznie obraca się z układem.

Jednakże, w przeciwieństwie do przypadku III, bardzo powoli oddala się od środka masy. W pewnym momencie następuje nagłe, gwałtowne oddalanie się do środka masy układu, które później po prostu równomiernie trwa (przekroczona zostaje prędkość ucieczki).

Wnioski:

Ciało znajdujące się dokładnie w punkcie Lagrange’a postaje nieruchome względem układu gwiazd, jednakże nawet najmniejsze przesunięcie doprowadzi do stopniowej utraty

równowagi, aż w końcu do znacznego oddalenia się od punktu Lagrange’a i szybkiego ruchu względem gwiazd. Jednakże wyposażenie teleskopu bądź stacji kosmicznej która ma być umieszczona w punkcie Lagrange’a nawet w niewielki silnik korygujący położenie mogłoby rozwiązać problem, gdyż początkowa faza wychylania się jest bardzo łagodna. Jednakże, punkty Lagrange’a są punktami równowagi chwiejnej (jeżeli nastąpi wychylenie, to będzie ono samoczynnie narastało).

Ciekawym zjawiskiem jest, że nieruchoma z początku planeta będąca w pobliżu układu dwóch gwiazd może zbliżając się do nich nabrać na tyle dużej prędkości, aby zacząć oddalać się od układu i nigdy do niego nie wrócić (tzn. przekroczyć prędkość ucieczki).

Aby prawidłowo opisać zachowanie planety w układzie dwóch gwiazd, bardzo ważne czasami okazuje się przeprowadzenie symulacji dla dłuższego czasu (tak jak to miało miejsce w przypadku II).

Każdy z dobieranych parametrów ma istotny wpływ na zachowanie planety w układzie. W zależności od doboru parametrów, orbita planety może przybierać ciekawe, często bardzo skomplikowane kształty. Możliwe przypadki można podzielić na przypadki ucieczki, oraz przypadki przyciągnięcia do układu. Jeżeli następuje ucieczka planety, to jej tor poza układem gwiazd jest spiralą.

Bibliografia:

 dr inż. A. Brozi, Wykłady z fizyki komputerowej, Politechnika Łódzka 2009r.

 D. Halliday, R. Resnick, J. Walker, Podstawy Fizyki, tom 2., Wydawnictwo Naukowe PWN, Warszawa 2005

 Informacje z Wikipedii o punktach Lagrange’a (http://en.wikipedia.org/wiki/Lagrangian_point)