第1回 軌道解析（はやぶさカプセル）[python]
小惑星からサンプルを採取し地球に持って帰ってくるとても技術力の必要なミッション（はやぶさミッション）を日本は成功させています。
目的の惑星に行くまでも大変ですが、さらに地球にサンプルを持って帰ってくる必要があります。
カプセルが地球に着地するまでに、地球の大気圏を突破しなければなりません。大気圏を飛行中に高温に包まれると聞きますが、どののような環境に晒されるかシミュレーションしてみたいと思います．
使用したコードはGithubにあげてあります．
次回はこちら。
1 支配方程式
今回は、重力と抗力、揚力を考慮した２次元の運動方程式を解きます。
$$ \frac{d V}{dt} = \frac{-\rho V^2}{2} \frac{C_{d} A}{m} – G sin \gamma $$
$$ \frac{d \gamma}{dt} = \frac{\rho V^2}{2} \frac{C_{l} A}{m} – \left( \frac{G}{V} – \frac{V}{r} \right) cos \gamma $$
$$ \frac{d h}{dt} = V sin \gamma $$
$$ \frac{d \theta}{dt} = \frac{V cos \gamma }{r}$$
ここで、\(V\)は速度、\(\gamma\)は経路角、\(\rho\)は空気の密度、\(C_{d}\)は抗力係数、\(A\)はカプセルの投影面積、\(m\)はカプセルの質量、\(G\)は重力加速度、\(C_{l}\)はカプセルの揚力係数、\(r\)は地球の中心からカプセル中心までの距離，\(\theta\)は初期位置に対し地球中心とカプセルのなす角ですです。
ここで、上の図において下記の式が成り立ちます。
$$r = h+R_{E}$$
また、重力加速度Gは高度によって変わるため、下記の式に置き換えます。
$$ G = g \frac{R_{E}^2}{r^2}, g = 9.81 \ m/s^2 $$
空気の密度も高度によって変わりますが、補間式を作成することは難しいので、大気モデルとしてNRLMSIS Atmosphere Modelを使用します。
下記のページに日付やほしい高度等を入力することで、データを取得することが来ます。
atmospheremodel : https://kauai.ccmc.gsfc.nasa.gov/instantrun/msis
線形補間により、各高度における密度を取得します。また，空気抵抗が大事になる今回のような解析では，日付の影響はほとんどないと思います．
以上より支配方程式は下記のようになります。
$$ \frac{d V}{dt} = \frac{-\rho V^2}{2} \frac{C_{d} A}{m} – g \frac{R_{E}^2}{r^2} sin \gamma $$
$$ \frac{d \gamma}{dt} =\frac{\rho V^2}{2} \frac{C_{l} A}{m} – \left( \frac{g}{V} \frac{R_{E}^2}{r^2} – \frac{V}{r} \right) cos \gamma $$
$$ \frac{d h}{dt} = V sin \gamma $$
$$ \frac{d \theta}{dt} = \frac{V cos \gamma }{r}$$
$$r = h+R_{E}$$
ただし\(\theta\)は独立しているため，解かなくても問題ありません．
2 解析条件
今回は、下記の解析条件[1]を用いてはやぶさカプセルの解析を実施します。
はやぶさカプセルは抗力型カプセルと呼ばれており、重心位置がカプセルの中心線上にあり左右対象形状のため、カプセルがまっすぐ向いている限りは揚力は発生しません。
以上を踏まえて、下記のような条件で解析を実施しました。
また、時間積分として、1次後退差分と4次ルンゲクッタの2ケース実施しました。
次回は、結果と考察を述べます。
hide
参考文献
[1] Takahashi, Yusuke, and Kazuhiko Yamada, “Aerodynamic-heating analysis of sample-return capsule in future Trojan-asteroid exploration,” Journal of thermophysics and heat transfer 32.3 (2018): 547-559.
[…] 前回はこちら。 […]