楽しい科学（理論）チャンネル
万有引力のはたらく場にある物体
/* 補足 極座標 */
周期的な運動では、位置座標を\((x,y)\)と表すより、原点\(O\)からどれだけの距離\(r\)離れているか、\(x\)軸と\(r\)を線分とみたときどれだけの角度\(\theta\)があるかで位置を表せる極座標\((r,\theta)\)が便利だ。位置座標を
\[\boldsymbol{x}=\begin{pmatrix}x\\y\end{pmatrix}\]
と表すのであれば、極座標を用いると、
\[\boldsymbol{x}=\begin{pmatrix}r\cos\theta\\r\sin\theta\end{pmatrix}\]
である。速度はそれぞれ、
\[\dot{\boldsymbol{x}}=\begin{pmatrix}\dot x\\\dot y\end{pmatrix}=\begin{pmatrix}\dot r\cos\theta-r\sin\theta\ \dot\theta\\\dot r\sin\theta+r\cos\theta\ \dot\theta\end{pmatrix}\]
これらを使うことで、\(|\boldsymbol{x}|=r,\ |\boldsymbol{x}|=\sqrt{\dot r^2+r^2\dot\theta^2}\)を示すことができる。保存力だった場合エネルギー保存則が使えて、
\[\frac{1}{2}m(\dot r^2+r^2\dot\theta^2)+U(r,\theta)=E\]
向心力だった場合角運動量保存則が使える。運動がうまく\(xy\)平面に収まりように座標系を取れるのであれば、\(\boldsymbol{L}=m\boldsymbol{x}\times\boldsymbol{v}\)は\(z\)成分のみになるので、
\[m(x\dot y-y\dot x)=L\]
\[mr^2\dot\theta=L\]
微分方程式\(\cfrac{dy}{dx}=-\sqrt{1-49(y-3)^2}\)の解答例
\(\cos x\)の微分は
\[\cfrac{d}{dx}\cos x=-\sin x=-\sqrt{1-\cos^2 x}\]
これを応用して、
\[\frac{dy}{dx}=-\sqrt{1-49(y-3)^2}\]
このような微分方程式を解くことができる。\(Y=7(y-3)\)とすると、\(\cfrac{dY}{dx}=7\cfrac{dy}{dx}\)なので、
\[\frac{1}{7}\frac{dY}{dx}=-\sqrt{1-Y^2}\]
\(Y=\cos z(x)\)と置き換えると、
\[\frac{dY}{dx}=\frac{dY}{dz}\frac{dz}{dx}=-\sin z\frac{dz}{dx}\]
であるからもとの微分方程式に代入して、
\[-\frac{1}{7}\sin z\frac{dz}{dx}=-\sqrt{1-\cos^2 z}\]
\[-\frac{1}{7}\sin z\frac{dz}{dx}=-\sin z,\ \frac{dz}{dx}=7\]
両辺を\(x\)で積分すると、
\[z=7x+C\]
(\(C\)は積分定数)以上から、
\[y=\frac{1}{7}Y+3=\frac{1}{7}\cos(7x+C)+3\]
今回万有引力を受ける物体の微分方程式を解くときにこの解法を用いる。
今回は万有引力場の物体の運動の一般解を求めたい。太陽自体も他の惑星に依る万有引力により動いているのだが、太陽は原点に固定して考える。
２次元の話をする。質量\(M\)の物体を原点に貼り付ける。位置\(\boldsymbol{x}=(x,y)\)にある質量\(m\)の物体が、2物体間の距離の2乗に反比例する万有引力
\[\boldsymbol{F}(\boldsymbol{x})=-\frac{GMm}{|\boldsymbol{x}|^2}\hat{\boldsymbol{x}}\tag{1}\]
の力を受けている場合、どのような運動をするだろうか。
\[\frac{d^2}{dt^2}\boldsymbol{x}(t)=-\frac{GMm}{|\boldsymbol{x}|^2}\hat{\boldsymbol{x}}\]
この方程式を時間\(t\)をパラメータとして、\(\boldsymbol{x}\)について解きたい。しかし解を求めることは難しいため、工夫をする。\(-\nabla U=\boldsymbol{F}\)を満たす\(U\)があればエネルギー保存則が使えるので、
\[U(\boldsymbol{x})=-\frac{GMm}{|\boldsymbol{x}|}=-\frac{GMm}{\sqrt{x^2+y^2}}\]
のようなスカラー関数を考える。天下どりになるが\(-\nabla U\)を計算しよう。
\[-\nabla U(\boldsymbol{x}(x,y))=-\frac{\partial U}{\partial x}\boldsymbol{e}_x-\frac{\partial U}{\partial y}\boldsymbol{e}_y\]
\[=-\frac{\partial U(\boldsymbol{x})}{\partial |\boldsymbol{x}|}\frac{\partial|\boldsymbol{x}|}{\partial x}\boldsymbol{e}_x-\frac{\partial U(\boldsymbol{x})}{\partial |\boldsymbol{x}|}\frac{\partial|\boldsymbol{x}|}{\partial y}\boldsymbol{e}_y\]
\[=-\frac{\partial U(\boldsymbol{x})}{\partial |\boldsymbol{x}|}\left(\frac{\partial|\boldsymbol{x}|}{\partial x}\boldsymbol{e}_x+\frac{\partial|\boldsymbol{x}|}{\partial y}\boldsymbol{e}_y\right)\]
\(|\boldsymbol{x}|=\sqrt{x^2+y^2}\)なので、
\[=-\frac{GMm}{|\boldsymbol{x}|^2}\left(\frac{2x}{2\sqrt{x^2+y^2}}\boldsymbol{e}_x+\frac{2y}{2\sqrt{x^2+y^2}}\boldsymbol{e}_y\right)\]
\[=-\frac{GMm}{|\boldsymbol{x}|^2}\frac{x\boldsymbol{e}_x+y\boldsymbol{e}_y}{\sqrt{x^2+y^2}}=-\frac{GMm}{|\boldsymbol{x}|^2}\frac{\boldsymbol{x}}{|\boldsymbol{x}|}=\boldsymbol{F}\]
よし！エネルギー保存則が使える。次に角運動量保存則が使えるか調べよう。
\[\frac{d}{dt}\boldsymbol{L}=\boldsymbol{x}\times\boldsymbol{F}\]
\[=\boldsymbol{x}\times\left(-\frac{GMm}{|\boldsymbol{x}|^2}\hat{\boldsymbol{x}}\right)=\left(-\frac{GMm}{|\boldsymbol{x}|^2}\right)\boldsymbol{x}\times\hat{\boldsymbol{x}}=\boldsymbol{0}\]
\(\boldsymbol{x}\)と\(\hat{\boldsymbol{x}}\)は平行なので外積は\(\boldsymbol{0}\)になるため、角運動量は時間変化しないベクトルとなる。つまり角運動量保存則も使える！
極座標を使うことで、力学的エネルギー保存則は、
\[\frac{1}{2}m(\dot r^2+r^2\dot\theta^2)-\frac{GMm}{r}=E\tag{2}\]
角運動量保存則は、
\[mr^2\dot\theta=L\tag{3}\]
\(L/m=L_m,\ E_m=E/m\)とすると、(2)(3)式は、
\[\frac{1}{2}(\dot r^2+r^2\dot\theta^2)-\frac{GM}{r}=E_m\tag{4}\]
\[\dot\theta=\frac{L_m}{r^2}\tag{5}\]
(5)式を(4)式に代入すると、
\[\frac{1}{2}\left(\dot r^2+\frac{L_m^2}{r^2}\right)-\frac{GM}{r}=E_m\]
\[\dot r=\pm\sqrt{2E_m+\frac{2GM}{r}-\frac{L_m^2}{r^2}}\]
\[\frac{dr}{dt}=\pm L_m\sqrt{\frac{2E_m}{L_m^2}+\frac{2GM}{L_m^2r}-\frac{1}{r^2}}\]
\(r=1/u\)とすると、\(dr=-du/u^2\)であるから、
\[-\frac{1}{u^2}\frac{du}{dt}=\pm L_m\sqrt{\frac{2E_m}{L_m^2}+\frac{2GM}{L_m^2}u-u^2}\]
\[-r^2\frac{du}{dt}=\pm L_m\sqrt{\frac{2E_m}{L_m^2}+\frac{G^2M^2}{L_m^4}-\left(\frac{GM}{L_m^2}-u\right)^2}\]
\(r\)を\(t\)の関数として解を求めるのは難しいため一度\(r\)を\(\theta\)の方程式で表す。(5)式から、\(r^2=L_mdt/d\theta\)であるから、
\[-L_m\frac{dt}{d\theta}\frac{du}{dt}=\pm L_m\sqrt{\frac{2E_m}{L_m^2}+\frac{G^2M^2}{L_m^4}-\left(\frac{GM}{L_m^2}-u\right)^2}\]
\[-\frac{du}{d\theta}=\pm\sqrt{\frac{G^2M^2}{L_m^4}\left(1+\frac{2E_mL_m^2}{G^2M^2}\right)-\left(\frac{GM}{L_m^2}-u\right)^2}\]
\(A=L_m^2/GM,\ e=\sqrt{1+2E_mL^2_m/G^2M^2}\)とすると、
\[-\frac{du}{d\theta}=\pm\sqrt{\frac{e^2}{A^2}-\left(\frac{1}{A}-u\right)^2}\tag{6}\]
解の形を
\[u=\frac{1}{A}(1+e\cos\xi)\tag{7}\]
と仮定する。(6)式に代入すると、
\[\frac{e}{A}\sin\xi\frac{d\xi}{d\theta}=\pm\sqrt{\frac{e^2}{A^2}-(\frac{e}{A}\cos\xi)^2}\]
\[\frac{e}{A}\sin\xi\frac{d\xi}{d\theta}=\pm \frac{e}{A}\sin\xi\]
\[\frac{d\xi}{d\theta}=\pm1\]
両辺を\(\theta\)で積分すると、\(\xi=\pm\theta\)(積分定数は0にする。)これを(7)式に代入して、
\[\frac{1}{r}=u=\frac{1}{A}(1+\cos\theta)\]
\[r=\frac{A}{1+e\cos\theta}\]
定数\(e\)は離心率という値で、\(e=0\)で円、\(e\lt1\)で楕円、\(e=1\)で放物線、\(e\ge1\)で双曲線になる。今回は楕円軌道の範囲で、この極座標の方程式を\((x,y)\)の式に直そう。\(\sqrt{x^2+y^2}=r,\ x=r\cos\theta\)であるから、
\[r(1+e\cos\theta)=A\]
\[r+er\cos\theta=A,\ \sqrt{x^2+y^2}=A-ex\]
両辺を2乗する。
\[x^2+y^2=A^2-2eAx+e^2x^2,\ (1-e^2)x^2+2eAx+y^2=A^2\]
\[x^2+\frac{2eA}{1-e^2}x+\frac{y^2}{1-e^2}=\frac{A^2}{1-e^2}\]
\[(x+\frac{eA}{1-e^2})^2+\frac{y^2}{1-e^2}=\frac{A^2}{1-e^2}+\frac{e^2A^2}{(1-e^2)^2}=\frac{A^2}{(1-e^2)^2}\]
\(\frac{A}{1-e^2}=a,\ \frac{A}{\sqrt{1-e^2}}=b\)とすると、
\[\frac{(x+ea)^2}{a^2}+\frac{y^2}{b^2}=1\tag{8}\]
中心が\((-ea,0)\)、長半径\(a\)、短半径\(b\)の楕円上を物体が運動するのは分かった。新たに(8)式を満たすような次のような媒介変数\(\phi\)を使って\((x,y)=(a\cos\phi-ea,b\sin\phi)\)としてもう一度、角運動量を求める。
\[\dot x=-a\sin\phi\ \dot\phi,\ \dot y=b\cos\phi\ \dot\phi\]
\[\boldsymbol{L}=\boldsymbol{r}\times(m\boldsymbol{v})=\begin{pmatrix}a\cos\phi-ea\\b\sin\phi\\0\end{pmatrix}\times\begin{pmatrix}-ma\sin\phi\ \dot\phi\\mb\cos\phi\ \dot\phi\\0\end{pmatrix}\]
\(x,y\)成分は０になるので、\(z\)成分を計算する。
\[L=mab\cos^2\phi\ \dot\phi-emab\cos\phi\ \dot\phi+mab\sin^2\phi\ \dot\phi\]
\[=mab(1-e\cos\phi)\dot\phi,\ \frac{L_m}{ab}dt=(1-e\cos\phi)d\phi\]
\(L_m,a,b\)は時間に依って変化しないので時間で積分すると
\[\frac{L_m}{ab}(t-t_0)=\phi-e\sin\phi\]
時間とパラメータ\(\phi\)の係数をもっときれいに表したい。\(L_m=\sqrt{AGM},b=\sqrt{1-e^2}a,a=A/(1-e^2)\)であるから、
\[\frac{L_m}{ab}=\frac{\sqrt{AGM}}{\sqrt{1-e^2}a^2}=\sqrt{\frac{GM}{a^3}}\]
万有引力を受ける物体の運動の解は、直接\(\boldsymbol{x}=\boldsymbol{x}(t)\)とするのは難しくてパラメータ\(\phi\)を使って、\(\boldsymbol{x}=\boldsymbol{x}(\phi(t))\)となる。
\[\begin{cases}
\sqrt{GM/a^3}(t-t_0)=\phi-e\sin\phi\\
(x,y)=(a\cos\phi(t)-ea,b\sin\phi(t))\end{cases}\]
となる。