はじめに
前回の記事では、太陽系を 円軌道近似でシミュレーションするところまでを実装しました。
👉 前回の記事
https://qiita.com/tsuruoka91/items/6a046db3d5ee3ae82bb8
しかし実際の惑星軌道は円ではなく楕円軌道です。
そこで今回は ケプラー方程式を用いて、より現実に近い軌道シミュレーションを行うように改良しました。
成果物はこちらで公開しています。
👉 GitHub
https://github.com/tsuruoka91/solar-system
コペルニクスとケプラー
ニコラウス・コペルニクス (1473-1543)
ポーランドの天文学者で、地球が宇宙の中心という天動説に対し、太陽が中心で地球を含む惑星がその周りを回るという地動説（太陽中心説）を提唱しました。
晩年の著書『天体の回転について』でその理論を発表し、天文学に革命をもたらしました。
しかし、彼のモデルはまだ円軌道を前提としており、予測精度は高くありませんでした
ヨハネス・ケプラー (1571-1630)
ドイツの天文学者で、コペルニクス地動説の熱烈な支持者でした。
偉大な観測家ティコ・ブラーエの膨大な観測データを受け継ぎ、特に火星の軌道研究に没頭しました。
その結果、惑星は真円ではなく楕円軌道を描くこと、そして「ケプラーの3法則」を発見し、コペルニクス理論に数学的な精密さと実証性をもたらしました
前回の記事で円軌道を使って太陽系のシミュレートを行いましたが、これはコペルニクスのモデルだったんですね。今回はより制度を上げるためにケプラーのモデルに挑戦してみます。
なぜケプラー方程式が必要か
円軌道の場合、角速度は一定なので
θ = ωt
で位置が簡単に求まります。
しかし楕円軌道では、
- 惑星は一定速度では動かない
- 太陽に近いほど速く、遠いほど遅くなる（面積速度一定）
という性質があります。
この非線形な運動を扱うために登場するのが ケプラー方程式です。
ケプラー方程式の概要
楕円軌道では、以下の3つの角度が登場します。
- 平均近点角 M
- 離心近点角 E
- 真近点角 ν
その中核となるのがケプラー方程式：
M = E - e\sin E
- e : 離心率
- M : 時刻から一意に決まる
- E : この式を 数値的に解く必要がある
ここが今回の実装の一番のポイントです。
離心近点角 E を数値的に解く
ケプラー方程式は解析的に解けないため、今回は ニュートン法を使って解いています。
擬似コードで書くと以下のようになります。
let E = M;
for (let i = 0; i < 10; i++) {
E = E - (E - e * Math.sin(E) - M) / (1 - e * Math.cos(E));
}
- 初期値は E = M
- 数回の反復で十分に収束します
- 惑星数が少ないため計算コストも問題になりません
座標への変換
離心近点角 E が求まれば、軌道上の位置 (x, y) は次の式で計算できます。
\begin{aligned}
x &= a (\cos E - e) \\
y &= a \sqrt{1 - e^2} \sin E
\end{aligned}
- a : 軌道長半径
- 楕円軌道が自然に表現できるようになります
この座標をそのまま描画に使っています。
冥王星の追加
実際に描画させてみると、各惑星が楕円になっているのが確認できるのですが、もっと分かりやすい冥王星を追加してみます。
現在の定義では冥王星は準惑星ですが、
- 軌道の離心率が非常に大きい
- 他の惑星と明確に異なる挙動をする
- ケプラー方程式による楕円軌道の効果が分かりやすい
という理由から、シミュレーション題材としては非常に分かりやすい存在です。
| 名前 | 離心率 | 公転周期 |
|---|---|---|
| 冥王星(Pluto) | 0.2488 | 90560日 |
冥王星の軌道を見てみると一時的に冥王星の方が海王星よりも太陽に近くなる期間が存在することを確認できました。
ハレー彗星の追加
折角なので今度はハレー彗星を追加してみましょう。
ハレー彗星は約76年周期で公転する短周期彗星です。肉眼でも観測可能であり、私が子供の頃地球へ接近したことがあります。次に地球に接近するのは2061年とのことなので、今どのあたりを回っているのか確認してみましょう。
| 名前 | 離心率 | 公転周期 |
|---|---|---|
| ハレー彗星(Halley) | 0.967 | 27500日 |
問題はハレー彗星のユリウス日(J2000)の経度が調べても分からなかったことです。仕方がないのでシミュレーターの日時を前回ハレー彗星がもっとも太陽に近づいた1986年2月9日に設定して、数値を何度も調整しながら経度を予測してみました。その結果45という値でハレー彗星のおおよその位置を設定しています。
現在2025年はちょうど折り返し地点を過ぎた辺りにいることが確認できます。まだまだ旅は長いようです。
また、時間を進めてアニメーションしていくと、ゆっくり動いているハレー彗星が太陽に近づくとスピードが上がるという、ケプラーの第二法則(面積速度一定の法則)を確認できます
おわりに
今回は ケプラー方程式を使った楕円軌道シミュレーションを実装しました。
数式としては有名ですが、実際に「動くもの」として見ると理解が一段深まります。
物理や天文をコードで触ってみたい方の参考になれば幸いです。