はじめに
アルテミスⅡのオリオン宇宙船が無事帰還しました。
調べてみると、オリオン宇宙船は、地球周回軌道上から１回のエンジン噴射 （TLI：Trans-Lunar Injection）で月の裏側を通過する「自由帰還軌道」に載せ、月の重力によるスイングバイで地球に帰還しました。
以前、地球と月の軌道を状態方程式にしたシミュレーション 地球の運動について を書いたので、オリオン宇宙船の有人月周回ミッションの軌道もシミュレーションできるのではないかと考えました。 が、完全な軌道のシミュレーションは難しいので、面白そうなところだけ部分的にやってみます。
面白そうなところとは、地球周回軌道上からのエンジン噴射 （TLI：Trans-Lunar Injection）と、月の重力によるスイングバイです。
地球と月の軌道シミュレーションでは、制御入力はありませんでしたが、オリオン宇宙船は、エンジン噴射の推力により、地球周回から月周回に軌道を変えていきます。また、月の重力によるスイングバイも軌道の変更方法として大変興味深いです。
シミュレーションについて
上で書いたように、【地球周回軌道上からのエンジン噴射】と【月の重力によるスイングバイ】について、アルテミスⅡの計画を参考に軌道をシミュレーションします。
ただし、アルテミスⅡのオリオン宇宙船の軌道を正確に再現するわけではありません。厳密な条件などが分からないところもあるので、設定条件を簡略化してシミュレーションします。
宇宙船モデルの追加
地球と月の軌道を状態方程式にした 地球の運動について に宇宙船のモデルを追加します。
宇宙船は太陽と地球と月の引力の影響を受けますが、質量の差から太陽と地球と月への影響は無いとします。また、太陽と地球と月の公転軌道は同一平面としており、宇宙船の軌道も同じ2次元平面 \((X_{r},Y_{r})\) で考えていきます。（宇宙船のパラメータの添え字は r です。）
宇宙船の質量を\(M_{r}\)、太陽、地球、月の質量を\(M_{s}\), \(M_{e}\), \(M_{m}\)
宇宙船と各星の距離を\(L_{sr}\), \(L_{er}\), \(L_{mr}\)、ロケットエンジンの推力を\(Fi\) とすると
宇宙船に作用する力は、
\[F_{r}= -G \frac{M_{s} M_{r}}{{L_{sr}}^2} -G \frac{M_{e} M_{r}}{{L_{er}}^2} -G \frac{M_{m} M_{r}}{{L_{mr}}^2}+Fi\hspace{1cm}[N]\]さらに、宇宙船の移動速度 \(V_{r}\) とすると、
\[F_{r}=M_{r} \frac{dV_{r}}{dt}\]より\[ \frac{dV_{r}}{dt}= -G \frac{M_{s}}{{L_{sr}}^2} -G \frac{M_{e}}{{L_{er}}^2} -G \frac{M_{m}}{{L_{mr}}^2} + \frac{1}{M_{r}} Fi \]
2次元平面の座標 \(X_{r}, Y_{r}\) 、移動速度 \(V_{rx}, V_{ry}\) とすると
\[\frac{dX_{r}}{dt}=V_{rx} \hspace{1cm} \frac{dY_{r}}{dt}=V_{ry} \tag{1} \]\[ \frac{dV_{rx}}{dt}= -G \frac{M_{s}}{{L_{sr}}^2}\frac{X_{r}}{L_{sr}} -G \frac{M_{e}}{{L_{er}}^2}\frac{(X_{r} -X_{e})}{L_{er}} -G \frac{M_{m}}{{L_{mr}}^2}\frac{(X_{r} -X_{m})}{L_{mr}} + \frac{1}{M_{r}} \frac{V_{erx}}{V_{er}} Fi \]\[ \hspace{1cm} = G \frac{M_{e}}{{L_{er}}^3}X_{e} +G \frac{M_{m}}{{L_{mr}}^3}X_{m} -G(\frac{M_{s}}{{L_{sr}}^3} +\frac{M_{e}}{{L_{er}}^3} +\frac{M_{m}}{{L_{mr}}^3}) X_{r} +\frac{1}{M_{r}} \frac{V_{erx}}{V_{er}} Fi \tag{2} \]
\[ \frac{dV_{ry}}{dt}= -G \frac{M_{s}}{{L_{sr}}^2}\frac{Y_{r}}{L_{sr}} -G \frac{M_{e}}{{L_{er}}^2}\frac{(Y_{r} -Y_{e})}{L_{er}} -G \frac{M_{m}}{{L_{mr}}^2}\frac{(Y_{r} -Y_{m})}{L_{mr}} + \frac{1}{M_{r}} \frac{V_{ery}}{V_{er}} Fi \]\[ \hspace{1cm} = G \frac{M_{e}}{{L_{er}}^3}Y_{e} +G \frac{M_{m}}{{L_{mr}}^3}Y_{m} -G(\frac{M_{s}}{{L_{sr}}^3} +\frac{M_{e}}{{L_{er}}^3} +\frac{M_{m}}{{L_{mr}}^3}) Y_{r} +\frac{1}{M_{r}} \frac{V_{ery}}{V_{er}} Fi \tag{3} \]
宇宙船と太陽の距離 \(L_{sr}=\sqrt{ {X_{r}}^2+{Y_{r}}^2 }\)
宇宙船と地球の距離 \(L_{er}=\sqrt{ {(X_{e}-X_{r})}^2 + {(Y_{e}-Y_{r})}^2 }\)
宇宙船と月の距離 \(L_{mr}=\sqrt{ {(X_{m}-X_{r})}^2 + {(Y_{m}-Y_{r})}^2 }\)
ロケットエンジンの推力\(Fi\) のベクトル成分には、下記の宇宙船と地球の相対速度を使っています。その理由は 地球周回軌道上からのエンジン噴射（TLI：Trans-Lunar Injection） で説明します。
宇宙船と地球の相対速度
\(V_{erx}=V_{ex} -V_{rx}\), \(V_{ery}=V_{ey} -V_{ry}\)
\(V_{er}=\sqrt{ {(V_{ex}-V_{rx})}^2 + {(V_{ey}-V_{ry})}^2 }\)
以上の宇宙船のモデルを 地球の運動について の地球と月の軌道の状態方程式に追加します。
Excelでは、シート【状態方程式】の下記の状態方程式[A][B]のピンクの行と列が宇宙船モデルの追加部分です。
地球低軌道(LEO)での確認
組み込んだ宇宙船のモデルについて地球の低軌道(LEO)を周回するか確認します。
Excelファイル ⑧-3 非線形モデル SpaceCraft LEO(低軌道).xlsm は、こちら→リンク からダウンロードできます。
Excelファイルのワークシート【状態方程式】の上部の黄色のセルが宇宙船モデルの初期条件の情報になります。
宇宙船の重量は 26.5t 、地球を低軌道(LEO)で周回する場合の高度は 190km（地球の半径 6371km）、速度は 7.8km/sec、50secで離散化します。
Excelファイルのワークシート【シミュレーション】の上部にある【軌道計算】をクリックすると計算を開始します。計算を止める場合は【STOP】で止まります。
非線形なモデルであるため、離散時間毎に逐次式のPとQを再計算しています。詳しくは、こちら→ 地球の運動について でご確認ください。
青のラインが地球の軌道で、黄色が月、紫が宇宙船の軌道です。宇宙船が地球を周回する軌道になっています。
宇宙船のモデルを地球と月の軌道シミュレーションの中に組み込みました。アルテミスⅡでは、ロケットの打ち上げから帰還まで10日程の計画です。上記の軌道は下図の 〇 の部分になります。
地球周回軌道上からのエンジン噴射 （TLI：Trans-Lunar Injection）
実はオリオン宇宙船の地球周回軌道は、遠地点（地球から最も離れる点）が約74,000kmに達する大きな楕円軌道です。月に向かう軌道に切り替えるエンジン噴射（TLI：Trans-Lunar Injection）は、楕円軌道の近地点（高度約600〜2,000kmの範囲）で行われるようですが、このブログではアルテミスⅡの計画も参考にしながら、より条件を簡略化してシミュレーションします。
設定条件
オリオン宇宙船の重量 26.5 [t] （結構重いですね）
メインエンジンの推力 26700 [N]
エンジン推力印加時間 350 [sec]
TLI開始時の速度 10.5 [km/sec]
TLI修了時の速度 10.9 [km/sec]
TLI開始時の高度 900 [km]
高度 900 km の地球周回軌道から、350秒間推力を周回軌道の接線方向に印加します。TLI終了後の速度は 10.9 km/sec ですが、月に向かう過程で速度は低下していきます。
単位系について
通常単位系はSI単位系を使いますが、扱う数値が大きいので、単位系を変更しています。
単位時間：10000秒、単位距離：10000km、質量：kg、速度：km/秒
メインエンジンの推力 26700 [N] は、この単位系での分解能は 0.1 [N] になっているので、入力値は 267000 になります。
周回軌道の接線方向への推力印加
重力場の周回軌道上にある宇宙船は、速度ベクトルと同じ方向（接線方向）に推力を印加すると効率よく加速できるようです。このため、状態方程式の入力になる宇宙船の推力には、地球との相対速度のベクトル \( V_{er},V_{erx},V_{ery} \) を使っています。
シミュレーション
Excelファイル ⑧-3 非線形モデル SpaceCraft SwingBy.xlsm は、こちら→リンク からダウンロードできます。
この地球周回軌道上からのTLI（Trans-Lunar Injection）のシミュレーションでは、月への軌道が地球の周回軌道の接線方向となるように、宇宙船の初期位置を地球の6時の方向に設定しています。
その初期位置から、350秒間エンジン噴射による推力を軌道の接線方向に印加します。50秒で離散化しているので、推力の入力は7回になります。
エンジン噴射の有無を比較すると若干の違いには見えますが、月に向かう軌道になっています。
月の重力によるスイングバイ
月の軌道上で月の初期位置を調整すると、宇宙船を地球への帰還方向にスイングバイできるポイントがあります。
宇宙船は月の左側から月の軌道に接近。月の後方を横切り、月の右側（裏側）に出る。
さらに、宇宙船は月の重力によるスイングバイで月の前方を横切り、地球への帰還方向に軌道を変更しています。
Excelファイルのワークシート【シミュレーション】のセルF4の設定値 1475 は月の軌道上の初期位置を指定しています。ワークシート【シミュレーション】のセルY23から月の初期位置テーブルがあり、設定値 0～2000 の間で初期位置を指定します。
ただし、スイングバイできる月の初期位置の範囲はかなり狭く、そこを外すと月の重力による軌道変更が発生しなかったり、逆に下図のように月の重力に捕捉されてしまったりが発生します。
取り組んでみた感想
月の重力による宇宙船のスイングバイの物理事象を状態方程式にしてシミュレーションで確認できたので、当初の目的は果たせたと思います。
ただし、スイングバイできたと書いていますが、実は宇宙船と月の距離を確認すると、月の半径より内側を通過していました。残念。さらに、地球への帰還軌道方向に宇宙船をスイングバイしますが、完全な地球への帰還軌道には載っていません。
試行錯誤の途中で気が付いたのは、地球への帰還軌道側から遡って軌道設計していく必要があるということ。地球への帰還軌道に載せるための月の重力によるスイングバイのタイミングと月の位置、そのタイミングに合わせたエンジン噴射 （TLI：Trans-Lunar Injection）の位置、さらに打ち上げのタイミングや地球周回軌道 と決めていく必要があります。それに軌道設計するための基本的な軌道工学の知識も必要です。軌道工学を学ぶには、こちら→リンク1 , リンク2 が参考になるかもしれません。
最後までお読みいただき、ありがとうございました。
しかし、アポロ8号が最初に有人で月を周回したのは 1968年です。どうやって軌道計算していたのかなぁ。当時使っていた計算機は IBM System/360。調べてみると最上位機種でも演算性能は 1.7MIPS 程です。
さらに、宇宙船に搭載されていたアポロ誘導コンピュータ（AGC）は、クロック周波数 2MHz、RAM 2KB、ROM 36KB。Arduino UNO R3 と同等スペックのコンピュータ（Arduino UNO R3は16MHzだけど）に宇宙船の誘導、制御ロジック、カルマンフィルタを実装して月に行きました！すごいねぇ。
0 件のコメント:
コメントを投稿