惑星シミュレータ（古典力学の範疇内で）

前フリprologue

宇宙、と聞くとどうしてもプトレマイオスのような同心円の物を浮かべてしまいます．しかし、実際地球が動く軌道は太陽を一つの焦点とする楕円軌道を描いています．これは多体問題を考えなくても古典的な2体問題の解として導くことができます．（基底を(x，y)から(r,θ)に取替えなければ導くのは難しいが…）

(*At this part, I will write the contents of this diary in English.The contents are the same as the one in Japanese.Therefore, you don't need to read this part if you understand this diary in Japanese).
When we think of the planets,we tend to think its orbit's shape as the world of Ptolemaios(,which means not earth-centered theory but the orbit of the planet is concentric circle).But the acutual orbit of the Earth is elliptic one.That result is well-known as the solution of 2 body problems.(It is difficult to lead the solution if we don't use (r,θ) coordinate…).

▲2体問題の解．確かに，原点が楕円軌道の焦点となっている．

The solution of 2 body problems.Acutually, the origin is focal spot of the Earth's elliptic orbit

目的(purpose of this simulation)

もう少し考えてみましょう．宇宙は古典力学では3次元平面と考えられます．僕らは地球と太陽の存在する平面を勝手に設定しましたが、はたしてこれは妥当なのでしょうか？　任意の平面でちゃんと地球は公転するのでしょうか？　例えば，x-y平面では回転するのだけど，そこから45度ほど傾けた平面にすると地球がふっとんでいってしまう，ということはないのでしょうか？　そんなことはありえない．そのことを今回，シミュレーションを使って確かめていきます．

Let's think more detaily.You know we can think our space is 3-dimentional plane space on the classic(Newton's) mechanics.We assumed that the Earth and the Sun exists on x-y plane,but there are some doubts which says the assumption is appropriate? I wonder whether the Earth can revolve on the arbitrary plane? For example, if we put the earth on the plane which leans by 45 degree from x-y plane, doesn't the situation occure where the Earth goes for far? I can say, "the situation is non-sense." I will be sure for it by using the simulation.

問題設定(setting of this simulation)

地球は太陽からしか力を受けない，と仮定します．すると、次の運動方程式が成り立ちます．
$m\vec{a}=\frac{GMm}{r^{3}}\vec{r}$ (M は太陽の質量,m は地球の質量,G は万有引力定数)
$\frac{\vec{r}}{r}$ はr方向の単位ベクトルを表すためのものですので実質的には， $\frac{GMm}{r^{2}}$ を考えてると，これが地球にかかる力の大きさを表していることになります．後は球面座標系を用いて、その力の成分をx,y,z座標ごとに分解してやれば，各時間ごとに地球にかかる力が求められることになります．

We assume that the Earth is forced from only the Sun.Then,we get the following M.O.E..
$m\vec{a}=\frac{GMm}{r^{3}}\vec{r}$ (M is the Sun's mass,m is Earth's mass,G is the gravitational constant)
Where $\frac{\vec{r}}{r}$ is used only for we intend to express the unit vector directed for r, so we can get the summation of force forced to the Earth by using the next term
$\frac{GMm}{r^{2}}$
Next,we can get the force of each direction(x,y,z) to use spherical-coordinate.

使ったアルゴリズム(algorithm)

4次のルンゲ・クッタルーチンを使用しました．詳しい説明はここには載せません．（ググって見てください）
この方式を用いると，方程式系としては次の二つを考えればよくなります．
$m\frac{d\vec{v}}{dt}=\frac{GMm}{r^{3}}\vec{r}$
$\frac{d\vec{r}}{dt}=\vec{v}$

I used the forth degree Runge-Kutta method. I don't explain about this method here.(So, you should google it if you want to get more knowledge on it) .
Then we can get solution only by thinking for following two equations.
$m\frac{d\vec{v}}{dt}=\frac{GMm}{r^{3}}\vec{r}$
$\frac{d\vec{r}}{dt}=\vec{v}$

結果(result)

これらの方程式にルンゲ・クッタルーチンを適用し，公転周期ぐらいの時間をシミュレートしてやります．すると次の軌道を得ました．

I used Runge-Kutta method for the set of equations and got the next orbit.Where I set the time in simulation world is the period of revolution.

▲地球の3次元の公転軌道。

the Earth's orbit of revolution on 3-D world.

原点に太陽があります。太陽からの距離が1.7×10^8kmとなる球面を考え、その1点 $\vec{R}=\biggl(\array{1.202\times10^8\\0\\1.202\times10^8}\biggr)$ に球に接する速度 $\vec{v}=\biggl(\array{0.77\times10^{3}\\0.77\times10^{3}\\-0.77\times10^3}\biggr)$ (km/minute)を置きました。その結果、3次元球面内でも地球は安定軌道を描くことが分かります。

I set the Sun is the origin.I think spherical surface whose distance from the sun is 1.7×10^8km, and set the particle(the Earth) whose velocity is $\vec{v}=\biggl(\array{0.77\times10^{3}\\0.77\times10^{3}\\-0.77\times10^3}\biggr)$ (km/minute) at the one point $\vec{R}=\biggl(\array{1.202\times10^8\\0\\1.202\times10^8}\biggr)$ on that spherical surface.(and the velocity is　tangential from the spherical surface).And I found that the Earth went stable orbit if it is put on not x-y plane.)

おまけ(appendix)
太陽の重さを2倍ほどにしてみたところ，地球は次のような軌道を描きます．