Wyprowadzenie prawa grawitacji Newtona z praw Keplera  -  K. Siekański, I rok SMP

na podstawie: M. Makowiecki - “Mechanika teoretyczna” (http://pl.wikibooks.org/wiki/Mechanika_teoretyczna/)

Przyspieszenie dośrodkowe

Niech położenie ciała o masie m, krążącego po orbicie eliptycznej wokół ciała o masie M, dane będzie współrzędnymi kartezjańskimi x, y. Wówczas po przejściu na współrzędne biegunowe r, ϕ i obrocie o kąt ϕ, biorąc pierwszą współrzędną przyśpieszenia, otrzymuje się jego składową równoległą do promienia wodzącego r:

x=r[t] Cos[ϕ[t]];
y=r[t] Sin[ϕ[t]];
R={x, y};
A=((∂_({t,2}))#1&)/@R;
OB={{Cos[ϕ[t]],-Sin[ϕ[t]]},{Sin[ϕ[t]],Cos[ϕ[t]]}};
Ar=A.OB;
a=Ar[[1]] // Simplify;
%// TraditionalForm

Druga zasada dynamiki

Wobec powyższego, druga zasada dynamiki dla rozważanego ciała przyjmuje we współrzędnych biegunowych postać:

row3 = F == a*m ;
% // TraditionalForm

Ruch po orbicie eliptycznej

Z kolei zróżniczkowanie zależności promienia wodzącego ciała krążącego po elipsie od kąta obrotu, daje przydatne równanie:

Clear[r,t]
X =r[t]==  k/(1+e Cos[ϕ[t]]);
X = 1/X[[#]] &/@ {1,2};
X = X[[1]] == X[[2]]; 
  XX = D[X,t] // Simplify;
% // TraditionalForm

Moment pędu

Następnie, z definicji momentu pędu wnosimy:

Clear[L,Y]
Y = L==m ((r[t])^2) ϕ'[t];
Y = Y[[#]]/(m r[t]) &/@ {1,2};
 Y = Y[[1]] == Y[[2]];
% // TraditionalForm

Podstawienie powyższego równania do równania poprzedniego, zwraca zależność:

XX = Roots[XX,r'[t]];
rown = Solve[{Y,XX},{ϕ'[t],r'[t]}];
ROW = r'[t]==(r'[t]/.Last[rown]);
% // TraditionalForm

Po powtórnym zróżniczkowaniu i kilku przekształceniach:

ROW = D[ROW,t];
rown7 = Solve[{ROW, Y},{ϕ'[t],r''[t]}];
 row1 =r''[t]==(r''[t] /. Last[rown7]);
row2 = ϕ'[t]==(L/(m ((r[t])^2)));
 row2 = row2[[#]]^2 *r[t] &/@ {1,2};
row2 = row2[[1]] == row2[[2]];
% // TraditionalForm

Siła radialna

Powyższe równania i znajomość geometrii elipsy pozwalają wyprowadzić wzór na siłę radialną:

Clear[F]
 row3 = F == a*m;
 k == b^2 / a1;
rown1 = Solve[{row1,row2,row3,X},{F,ϕ'[t],r''[t],e}] // Simplify;
 row4 = F == (F /.First[rown1]);
% // TraditionalForm

II i III prawo Keplera

Z II prawa Keplera (dla orbity o osiach a1 i b):

row5=π a1 b==(1/2)  ((r[t])^2) ϕ'[t]T;
% // TraditionalForm

co pozwala wyprowadzić okres obiegu :

rown2 = Solve[{row5,row2},{ T,ϕ'[t]}];
row6 = T ==(T /. Last[rown2]);
% // TraditionalForm

Wykorzystanie zaś III prawa Keplera wraz z powyższym wzorem i podstawienie:

((b^2)/a1)->k

daje:

row7 = row6[[#]]^2 / a1^3 &/@ {1,2};
row7 = row7[[1]] == row7[[2]] /. ((b^2)/a1)->k;
% // TraditionalForm

Postać ostatecznego wyniku

Dalsze uproszczenie można uzyskać poprzez wprowadzenie kolejnego oznaczenia:

row8= α==((L^2)/(k (m^2)));
% // TraditionalForm

Wówczas:

row8 = row8[[#]]*(m^2) &/@ {1,2};
row8 = row8[[1]] == row8[[2]];
row4 =row4/. ((L^2)/k)-> α (m^2) // Simplify;
% // TraditionalForm

Należy zauważyć, że jak dotąd nie skorzystano w niniejszym rozumowaniu z prawa grawitacji Newtona - użyto jedynie praw Keplera i zasad dynamiki. Jednak aby uzyskać końcową postać siły grawitacji działającej między dwoma ciałami, konieczne jest zaczerpnięcie od Newtona ważnego stwierdzenia, mówiącego, że źródłem grawitacji jest MASA. Wówczas dopiero można wyciągnąć wniosek, że stała α musi zależeć od masy M. 
	Z kolei z faktu, iż równanie na siłę radialną F powinno zachowywać swą postać niezależnie od tego, z którym z dwóch rozważanych ciał związany jest wybrany na początku układ współrzędnych, wynika, że jest to zależność liniowa w M.
	
Ostatecznie więc, przy podstawieniu:

α->G M

otrzymuje się :

row4/. α->G M // TraditionalForm