\documentclass[12pt,a4j]{jsarticle}

\title{ITPASS数値計算実習課題その１}
\author{高田 和佳　担当情報実験機:joho07}
\date{\today}

\begin{document}
\maketitle

\section{2体問題}
\subsection{慣性系における中心星と惑星の運動方程式}
中心星に働く力を\manbif{F_{12}}、惑星に働く力を\manbif{F_{21}}とおくと中心星の質量はm_1、惑星の質量はm_2なので、それぞれの運動方程式は

\begin{eqnarray}
\mbox{\boldmath ${F_{12}}$}=m_1a=m_1 \frac{d^2\mbox{\boldmath ${r_1}$}}{dt^2} \\
\mbox{\boldmath ${F_{21}}$}=m_2a=m_2 \frac{d^2\mbox{\boldmath ${r_2}$}}{dt^2} 
\end{eqnarray}

となる。また、 \displaystyle \mbox{\boldmath $\[  r=r_2-r_1\]$}と万有引力の法則を用いると

\begin{eqnarray}
\mbox{\boldmath ${F_{12}}$}=-Gm_1m_2 \frac{\mbox{\boldmath ${r_1-r_2}$}}{|\mbox{\boldmath ${r_1-r_2}$}|^3}\] \\
\mbox{\boldmath ${F_{21}}$}=-Gm_1m_2 \frac{\mbox{\boldmath ${r_2-r_1}$}}{|\mbox{\boldmath ${r_2-r_1}$}|^3}\] 
\end{eqnarray}
よって(1),(2),(3),(4)式より

\begin{eqnarray}
m_1 \frac{d^2\mbox{\boldmath ${r_1}$}}{dt^2}=-Gm_1m_2 \frac{\mbox{\boldmath ${r_1-r_2}$}}{|\mbox{\boldmath ${r_1-r_2}$}|^3}\] \\
m_2 \frac{d^2\mbox{\boldmath ${r_2}$}}{dt^2}=-Gm_1m_2 \frac{\mbox{\boldmath ${r_2-r_1}$}}{|\mbox{\boldmath ${r_2-r_1}$}|^3}\] 
\end{eqnarray}
と書ける。ここから(5),(6)式を整理していく。\\
(5)式を両辺$m_{1}$で、(6)式を両辺$m_{2}$で割ると
\begin{eqnarray}
\frac{d^2\mbox{\boldmath ${r_1}$}}{dt^2}=-Gm_2 \frac{\mbox{\boldmath ${r_1-r_2}$}}{|\mbox{\boldmath ${r_1-r_2}$}|^3}\] \\
\frac{d^2\mbox{\boldmath ${r_2}$}}{dt^2}=-Gm_1 \frac{\mbox{\boldmath ${r_2-r_1}$}}{|\mbox{\boldmath ${r_2-r_1}$}|^3}\] 
\end{eqnarray}
(4)−(3)式に、 \displaystyle \mbox{\boldmath $\[  r=r_2-r_1\]$}　を代入すると

\begin{eqnarray}
\frac{d^2\mbox{\boldmath ${r}$}}{dt^2}=- \frac{G(m_1+m_2)}{r^3}{\mbox{\boldmath ${r}$}}
\end{eqnarray}
が求まる。\ddot{r}=d^2r/dt^2より(9)は

\begin{eqnarray}
\ddot{r}=- \frac{G(m_1+m_2)}{r^3}{\mbox{\boldmath ${r}$}}
\end{eqnarray}
となり、求める式を導出することが出来た。また上記運動方程式は、この2体がケプラー運動をしていることを表している。

\subsection{速度ベクトル}
定義より

\begin{eqnarray}
(\dot{v_x},\dot{v_y})=\frac{d\mbox{\boldmath ${v}$}}{dt}=\frac{d^2\mbox{\boldmath ${r}$}}{dt^2}
\end{eqnarray}
$ \mbox{\boldmath ${r}$}=(x,y) $、(9)を用いて(11)は

\begin{eqnarray}
(\dot{v_x},\dot{v_y})=- \frac{G(m_1+m_2)}{(x^2+y^2)^\frac{3}{2}}{(x,y)}
\end{eqnarray}
と書ける。よって\dot{v_x}と\dot{v_y}は

\begin{eqnarray}
\dot{v_x}=- \frac{G(m_1+m_2)}{(x^2+y^2)^\frac{3}{2}}{x} \\
\dot{v_y}=- \frac{G(m_1+m_2)}{(x^2+y^2)^\frac{3}{2}}{y}
\end{eqnarray}
と書くことが出来る。

\section{3体問題}
\subsection{}
中心星と惑星は$x$軸上にあることから、$y_1$と$y_2$は共に0である。$x$座標について考えると、二つの座標系の原点が中心星と惑星の重心であることから
\begin{eqnarray}
|x_2|:|x_1|=\mu_1:\mu_2
\end{eqnarray}
(15),$r$=1,\mu_1+\mu_2=1より

\begin{eqnarray}
|x_1|=1\times \frac{\mu_2}{\mu_1+\mu_2} =\mu_2 \\
|x_2|=1\times \frac{\mu_1}{\mu_1+\mu_2} =\mu_1
\end{eqnarray}
よって

\begin{eqnarray}
(x_1,y_1)=(-\mu_2,0) \\
(x_2,y_2)=(\mu_1,0)
\end{eqnarray}

\subsection{}
まず、角速度の定義より\theta=\omega t　と書くことが出来る。\\
また、(\xi,\eta) は\theta=\omega t を用いて

\begin{eqnarray}
\xi=x\cos\omega{t}- y\sin\omega{t} \\
\eta=x\sin\omega{t}+ y\cos\omega{t}
\end{eqnarray}
と書くことが出来る。

\subsection{}
(20),(21)を$t$で微分すると

\begin{eqnarray}
\dot{\xi}=(\dot{x}-\omega{y})\cos\omega{t}-(\dot{y}+\omega{x})\sin\omega{t} \\
\dot{\eta}=(\dot{x}-\omega{y})\sin\omega{t}+(\dot{y}+\omega{x})\cos\omega{t}
\end{eqnarray}
さらに(22),(23)を微分すると

\begin{eqnarray}
\ddot{\xi}=(\ddot{x}-2\omega{\dot{y}}-\omega{^2}x)\cos\omega{t}-(\dot{y}+2\omega{\dot{x}}-\omega{^2}y)\sin\omega{t} \\
\ddot{\eta}=(\ddot{y}+2\omega{\dot{x}}-\omega{^2}y)\cos\omega{t}+(\dot{x}-2\omega{\dot{y}}-\omega{^2}x)\sin\omega{t}
\end{eqnarray}
以上より求める式を書くことが出来た。

\subsection{}
準備として(18)(19)(20)(21)を用いて\xi_1,\xi_2,\eta_1,\eta_2 の値をそれぞれ求める。

\begin{eqnarray}
\xi_1&=&x_1\cos\omega{t}-y_1\sin\omega{t}=-\mu_2\cos\omega{t} \\
\xi_2&=&x_2\cos\omega{t}-y_2\sin\omega{t}=\mu_1\cos\omega{t} \\
\eta_1&=&x_1\sin\omega{t}-y_1\cos\omega{t}=-\mu_2\sin\omega{t} \\
\eta_2&=&x_2\sin\omega{t}-y_2\cos\omega{t}=\mu_1\sin\omega{t}
\end{eqnarray}
与式と(24)(25)(26)(27)(28)(29)を用いると

\begin{eqnarray}
\ddot{\xi}&=&-\left(\mu_1\frac{x+\mu_2}{r_1^3}+\mu_2\frac{x-\mu_1}{r_2^3}\right)\cos\omega{t}+\left(\frac{\mu_1}{r_1^3}+\frac{\mu_2}{r_2^3}\right)y\sin\omega{t} \\
\ddot{\eta}&=&-\left(\mu_1\frac{x+\mu_2}{r_1^3}+\mu_2\frac{x-\mu_1}{r_2^3}\right)\sin\omega{t}+\left(\frac{\mu_1}{r_1^3}+\frac{\mu_2}{r_2^3}\right)y\cos\omega{t}
\end{eqnarray}
(30)(31)と(24)(25)より

\begin{eqnarray}
\ddot{x}-2\omega{\dot{y}}-\omega^2{x}&=&-\left[\mu_1\frac{x+\mu_2}{r_1^3}+\mu_2\frac{x-\mu_1}{r_2^3}\right] \\
\ddot{y}+2\omega{\dot{x}}-\omega^2{y}&=&-\left[\frac{\mu_1}{r_1^3}+\frac{\mu_2}{r_2^3}\right]y
\end{eqnarray}
よって回転系における粒子Pについての運動方程式を導出することが出来た。

\subsection{}
$r_1$,$r_2$の長さはそれぞれ

\begin{eqnarray}
r_1=\sqrt{\mathstrut (x+\mu_2)^2+y^2} \\
r_2=\sqrt{\mathstrut (x-\mu_1)^2+y^2}
\end{eqnarray}
と書くことが出来る。$\frac{\mu_1}{r_1}$と$\frac{\mu_2}{r_2}$をそれぞれ$x$で偏微分すると

\begin{eqnarray}
\frac{\partial }{\partial x} \frac{\mu_1}{r_1}&=&-\frac{1}{2}\mu_1\{(x+\mu_2)^2+y^2\}^{-\frac{2}{3}}\times(2x+2\mu_2) \nonumber \\
&=&-\mu_1\frac{(x+\mu_2)}{r_1^3} \\
\frac{\partial }{\partial x} \frac{\mu_2}{r_2}&=&-\frac{1}{2}\mu_2\{(x-\mu_1)^2+y^2\}^{-\frac{2}{3}}\times(2x-2\mu_1) \nonumber \\
&=&-\mu_2\frac{(x-\mu_1)}{r_1^3}
\end{eqnarray}
同様に$\frac{\mu_1}{r_1}$と$\frac{\mu_2}{r_2}$をそれぞれ$y$で偏微分すると

\begin{eqnarray}
\frac{\partial }{\partial y} \frac{\mu_1}{r_1}&=&-\frac{1}{2}\mu_1\{(x+\mu_2)^2+y^2\}^{-\frac{2}{3}}\times2y \nonumber \\
&=&-\frac{\mu_1}{r_1^3}y \\
\frac{\partial }{\partial y} \frac{\mu_2}{r_2}&=&-\frac{1}{2}\mu_2\{(x-\mu_1)^2+y^2\}^{-\frac{2}{3}}\times2y \nonumber \\
&=&-\frac{\mu_2}{r_1^3}y
\end{eqnarray}
$U$を$x$と$y$でそれぞれ偏微分すると(34)(35)(36)(37)より

\begin{eqnarray}
\frac{\partial U}{\partial x}&=&\omega{x^2}-\left[\mu_1\frac{(x+\mu_2)}{r_1^3}+\mu_2\frac{x-\mu_1}{r_2^3}\right] \\
\frac{\partial U}{\partial y}&=&\omega{x^2}-\left[\frac{\mu_1}{r_1^3}+\frac{\mu_2}{r_2^3}\right] y
\end{eqnarray}
よって、(38)(39)より(30)(31)は

\begin{eqnarray}
\ddot{x}-2\omega\dot{y}=\frac{\partial U}{\partial x} \\
\ddot{y}+2\omega\dot{x}=\frac{\partial U}{\partial y}
\end{eqnarray}
と書くことができ、$U$を用いて表すことが出来た。

\subsection{}
(40)(41)の$x$成分に$\dot{x}$を$y$成分に$\dot{y}$をかけて足し合わせる。

\begin{eqnarray}
\dot{x}\ddot{x}-2\omega{\dot{x}}{\dot{y}}=\dot{x}\frac{\partial U}{\partial x} \nonumber \\
\dot{y}\ddot{y}-2\omega{\dot{y}}{\dot{x}}=\dot{y}\frac{\partial U}{\partial y} \nonumber \\
\dot{x}\ddot{x}+\dot{y}\ddot{y}=\dot{x}\frac{\partial U}{\partial x}+\dot{y}\frac{\partial U}{\partial y}
\end{eqnarray}

(44)を時間について積分すると

\begin{eqnarray}
\frac{1}{2}(\dot{x}^2+\dot{y}^2)+C_J=U \nonumber \\
C_J=U-\frac{1}{2}(\dot{x}^2+\dot{y}^2)
\end{eqnarray}
C_Jは積分定数であり、円制限三体問題における保存量のヤコビ定数である。

\end{document}