05 - Kepler-probléma¶
A bolygómozgás potenciálja a következőképpen néz ki:
$$V(x,y)= -\gamma M \frac{m}{\sqrt{x^2+y^2}}$$Legyen most az alternatív naprendszerünkben $M=\gamma=m=1$!
Rajzoljuk fel a bolygó pályáját a sebességekre és helyekre vonatkozó egyenletek numerikus megoldásának segítségével:
\begin{align} m\dot{v_x}&=-\frac{\partial{V}}{\partial{x}} & v_{x,0}&=0\\\\ m\dot{v_y}&=-\frac{\partial{V}}{\partial{y}} & v_{y,0}&=\frac{1}{\sqrt{10}}\\ \dot{x}&=v_x & x_0&=4\\ \dot{y}&=v_y & y_0&=0\\ \end{align}Milyen alakú pályát kaptunk?
Megoldás¶
Először töltsük be a szükséges csomagokat
%pylab inline
from scipy.integrate import *
Ezután megadjuk a szükséges paramétereket, és függvényeket. Mivel a differenciálegyenletben szerepelnek $V$ parciális deriváltjai, szükség van olyan függvényekre is, amelyek kiszámítják ezeket adott $x$ és $y$ helyen.
gamma = 1
M = 1
m = 1
def dVdx(x, y):
return (gamma * M * m * x)/(x**2 + y**2)**(3/2)
def dVdy(x, y):
return (gamma * M * m * y)/(x**2 + y**2)**(3/2)
Ezután definiáljuk a függvényt, ami a kérdéses változó növekményét adja meg. Az előző feladathoz hasonlóan ez itt is egy vektor lesz: $$ \frac{\textrm{d}}{\textrm{d}t} \left(\begin{array}{c} v_x\\ v_y\\ x\\ y \end{array}\right) = \left(\begin{array}{c} - \frac{1}{m} \frac{\partial V}{\partial x}\\ - \frac{1}{m} \frac{\partial V}{\partial y}\\ v_x\\ v_y \end{array}\right) $$ Ez alapján az egyenletet megoljuk:
def f(r, t):
return [-dVdx(r[2], r[3]) / m, -dVdy(r[2], r[3]) / m, r[0], r[1]]
r_0 = [0, 1/sqrt(10), 4, 0]
t = linspace(0, 24.8, 500)
r = odeint(f, r_0, t)
Végül ábrázlhatjuk az így kapott pályát
plot(r[:,2], r[:,3])
xlabel("x-tengely", size=12)
ylabel("y-tengely", size=12)
title("A bolygó pályája", size=18, y=1.05)
A bolygó pályája az ábra alapján ellipszis, és mivel az energiája: $$ E_{teljes} = - 1 \cdot \frac{1}{\sqrt{0 + 16}} + \frac{1}{2} \cdot 1 \cdot \frac{1}{10} = - \frac{1}{5}$$ ami negatív, pontosan ezt is vártuk.