スペースデブリ除去のロボティクス — 非協力物体の捕獲戦略

地球の周りを秒速7km以上で周回する「宇宙ゴミ」が、いま宇宙開発最大のリスクになりつつあります。国際宇宙ステーション（ISS）は年に数回、デブリ衝突を回避するための軌道変更を行っています。2009年にはイリジウム33号とコスモス2251号が衝突し、数千個もの破片が撒き散らされました。もしデブリ同士の衝突が連鎖的に起こり始めたら — これがケスラーシンドロームと呼ばれる最悪のシナリオです。特定の軌道高度が使い物にならなくなり、衛星通信やGPS、地球観測といった私たちの日常を支えるインフラが丸ごと失われかねません。

この問題に立ち向かうのがADR（Active Debris Removal: 能動的デブリ除去）のロボティクス技術です。しかし、デブリは「非協力物体」です。制御不能で回転しており、捕獲のためのハンドルもなければ、こちらに合わせて動いてもくれません。宇宙でゴミを拾うことが、なぜこれほど難しいのかを理解するには、回転力学、接触ダイナミクス、ロボット工学の知識が交差する地点に立つ必要があります。

スペースデブリ除去のロボティクスを学ぶと、以下のような応用と理解が広がります。

• 軌道上サービス（OOS）: 故障衛星の修理・燃料補給など、協力的なターゲットへのサービスにも同じ捕獲技術が活きる
• 惑星探査: 小惑星サンプルリターン（はやぶさ2やOSIRIS-REx）における非協力天体への接近・接触戦略
• 宇宙持続可能性: 将来の大規模コンステレーションの安全な運用と廃棄計画の設計

本記事の内容

• スペースデブリ問題の現状とケスラーシンドローム
• ADR（能動的デブリ除去）の概念と必要性
• 非協力物体の定義と捕獲の難しさ
• タンブリングデブリの回転運動の数理
• 5つの捕獲手法の原理と比較（ロボットアーム・ネット・ハープーン・接着パッド・テザー）
• 捕獲後のデオービット戦略
• 主要ミッション（RemoveDEBRIS, ELSA-d, ClearSpace-1）の概要
• Pythonによるタンブリングデブリの回転運動シミュレーション

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

• 微小重力下の接触ダイナミクス — 宇宙空間での物体同士の接触・把持に関する力学的基礎

スペースデブリ問題の現状

軌道上は「過密」になりつつある

スペースデブリ（宇宙ゴミ）とは、役目を終えた人工衛星、ロケットの上段、衝突や爆発で生じた破片など、軌道上で機能していない人工物体の総称です。2024年時点で、地上レーダーで追跡可能な10cm以上のデブリは約36,500個、1cm以上まで含めると約100万個、1mm以上では1億個を超えると推定されています。

なぜ小さな破片がこれほど問題になるのでしょうか。低軌道（LEO）における相対速度は最大で秒速15kmにも達します。この速度では、わずか1cmの破片でも運動エネルギーは手榴弾に匹敵し、衛星の機能を破壊するのに十分です。10cmの破片なら衛星を完全に粉砕し、さらに大量の二次デブリを生成します。

ケスラーシンドローム — 連鎖的衝突の悪夢

1978年、NASAのドナルド・ケスラーは衝撃的な予測を発表しました。軌道上のデブリ密度がある臨界値を超えると、デブリ同士の衝突が新たなデブリを生み、その新しいデブリがさらに衝突を引き起こす連鎖反応が発生するというものです。この現象はケスラーシンドロームと呼ばれています。

ケスラーシンドロームのメカニズムを整理すると次のようになります。

1. 既存デブリの衝突: 大型デブリ同士が衝突し、数千個の破片が生成される
2. 破片の拡散: 生成された破片がそれぞれ異なる軌道に乗り、衝突確率を増大させる
3. 連鎖的増殖: 新たな破片が別の物体と衝突し、さらに破片が増える
4. 不可逆な環境劣化: 特定の軌道高度帯が実質的に使用不能になる

恐ろしいのは、この連鎖反応は人間が新しい物体を一切打ち上げなくても進行するという点です。NASAやESAのシミュレーションによれば、現在の軌道環境はすでに臨界に近く、年間5〜10個の大型デブリを能動的に除去しなければ、ケスラーシンドロームの進行を食い止められないとされています。

デブリ問題を数値で見る

デブリの増加傾向を把握するために、いくつかの重要な指標を見ておきましょう。

特に高度750〜1,000kmのLEO帯は、太陽同期軌道の需要が集中する一方で大気抵抗による自然落下が極めて遅く、デブリが何十年〜何百年も滞留します。この軌道帯こそ、ケスラーシンドロームのリスクが最も高い領域です。

ここまでで、デブリ問題がいかに深刻かがわかりました。では、この問題に対して具体的に何ができるのでしょうか。受動的な対策（衛星の防護シールド、衝突回避マヌーバ）だけでは不十分で、能動的にデブリを「取り除く」技術が求められています。

ADR（Active Debris Removal）の概念と必要性

受動的対策の限界

デブリ問題への対処は大きく3段階に分けられます。

1. デブリ緩和（Mitigation）: 新たなデブリの発生を防ぐ。パッシベーション（残留燃料の排出）、25年ルール（運用終了後25年以内に大気圏再突入）などのガイドライン
2. デブリ防護（Protection）: 既存デブリからの被害を防ぐ。ウィップルシールド、衝突回避マヌーバ
3. デブリ除去（Remediation）: 既存デブリを能動的に除去する。これがADR

緩和策は「これ以上増やさない」ための策であり、すでに軌道上にある約11,000トンのデブリには効果がありません。防護策は小型デブリには有効ですが、大型デブリの衝撃には耐えられません。したがって、ケスラーシンドロームを防ぐには、軌道上のデブリを「実際に取り除く」ADRが不可欠です。

ADRの基本アーキテクチャ

ADRミッションの典型的なシーケンスは以下のとおりです。

1. ターゲット選定: 衝突リスクの高いデブリを選定（質量 $\times$ 衝突確率で優先順位を付ける）
2. ランデブー: チェイサー衛星がデブリに接近し、相対距離を数十mまで縮める
3. 近傍運動: 数十mから数mまで接近し、ターゲットの状態（姿勢、回転率）を詳細に観測する
4. 捕獲: いずれかの手段でデブリを物理的に把持・拘束する
5. 安定化: 捕獲後、デブリの回転を抑制し、チェイサーとの結合体を安定させる
6. デオービット: 結合体の軌道を変更し、大気圏再突入させる（あるいはグレーブヤード軌道へ移動）

このシーケンスの中で、技術的に最も困難なステップが「捕獲」です。なぜなら、ターゲットとなるデブリは非協力物体（non-cooperative object）だからです。

ターゲット選定の戦略

ADRの効果を最大化するには、「どのデブリを優先的に除去するか」が重要です。主要な選定基準は次のとおりです。

• 質量: 大型デブリほど衝突時に大量の二次デブリを生成する。ロケット上段や大型衛星が優先
• 衝突確率: 混雑軌道にあるデブリほど衝突リスクが高い
• 複合指標: 質量 $M$ と衝突確率 $P_c$ の積 $M \times P_c$ で優先順位を付けるのが一般的

ESAのモデルでは、質量が1トン以上で高度750〜1,000kmの太陽同期軌道にある使用済みロケット上段が、最優先ターゲットとされています。

それでは、ADRの最大の技術的ハードルである「非協力物体の捕獲」とは、具体的に何が難しいのかを詳しく見ていきましょう。

非協力物体の特性 — なぜ捕獲が難しいのか

「非協力」の意味

ADRの文脈で非協力物体とは、以下の条件を持つターゲットを指します。

地上の産業用ロボットが対象物を把持する場合、対象物は通常静止しており、形状は既知で、専用の把持点が設計されています。ADRではこれらの前提がすべて崩れます。まるで暗闘の中で、回転しながら飛んでくる不定形の物体を、素手で捕まえるようなものです。

タンブリングの問題

非協力物体の捕獲を最も困難にしている要因がタンブリング（tumbling）です。姿勢制御を失った衛星やロケット上段は、外部トルク（残留磁気トルク、重力傾斜トルク、太陽輻射圧トルク）の累積効果により、複雑な回転運動を行っています。

観測データによると、典型的なデブリの回転率は以下のとおりです。

• 低回転: 0.5 deg/s 以下（比較的穏やかなケース）
• 中回転: 0.5〜5 deg/s（多くのデブリがこの範囲）
• 高回転: 5 deg/s 以上（捕獲が極めて困難）

特に問題なのは、回転が単一軸まわりとは限らない点です。3軸すべてに回転成分がある場合、ターゲット表面上の各点は複雑な軌跡を描きます。ロボットアームのエンドエフェクタをターゲットの特定点に合わせ続けるには、高精度・高速な追従制御が必要です。

形状不確かさの問題

デブリの多くは、運用中の衛星としての形状データが残っていますが、長年の微小デブリ衝突や熱サイクルにより、表面状態が変化している可能性があります。太陽電池パネルが折れ曲がっている、断熱材（MLI）が剥がれてひらひらしている — このような予期しない形状変化は、把持計画を狂わせます。

さらに、ターゲットに接近する際の照明条件も問題です。太陽光が強く当たる面は白飛びし、影の面は完全に暗い。こうした極端なコントラストの下で、3次元形状を推定しながらリアルタイムで接近戦略を修正する必要があります。

ここまで、非協力物体の捕獲が「通常のロボティクスの常識が通じない」極めて困難なタスクであることを説明しました。次に、この困難さの核心であるタンブリング運動を、数理的に記述します。回転の物理を理解することが、適切な捕獲戦略を設計するための出発点です。

タンブリングデブリの回転運動解析

オイラーの回転方程式

タンブリングデブリの回転運動を記述する基本方程式は、オイラーの回転方程式です。デブリの重心に固定されたボディ座標系（主軸座標系）で記述すると、外部トルクがない自由回転の場合は次のようになります。

ここで重要なのは、ボディ座標系は物体とともに回転する座標系であるという点です。慣性系（空間に固定された座標系）で書くと慣性テンソルが時間とともに変化して扱いにくいのですが、ボディ座標系で書けば慣性テンソルは定数になります。この座標系の選び方が、オイラーの方程式を簡潔にするカギです。

デブリの主慣性モーメントを $I_x, I_y, I_z$、ボディ座標系での角速度成分を $\omega_x, \omega_y, \omega_z$ とすると、トルクフリー条件下でのオイラーの方程式は次のとおりです。

$$ \begin{equation} I_x \dot{\omega}_x – (I_y – I_z)\omega_y \omega_z = 0 \end{equation} $$

$$ \begin{equation} I_y \dot{\omega}_y – (I_z – I_x)\omega_z \omega_x = 0 \end{equation} $$

$$ \begin{equation} I_z \dot{\omega}_z – (I_x – I_y)\omega_x \omega_y = 0 \end{equation} $$

これら3つの連立微分方程式が、デブリの自由回転の全貌を決定します。注目すべきは、各式の右辺（ここでは左辺の第2項）に他の2軸の角速度の積が含まれている点です。これは回転運動の本質的な非線形性を表しており、一般には解析的に解くことができません。

保存量 — 角運動量とエネルギー

トルクフリーの回転には、2つの保存量が存在します。これらは運動の大域的な性質を理解する上で欠かせません。

角運動量の大きさ $L$ は、慣性系で見た角運動量ベクトルの大きさであり、外部トルクがゼロなら一定に保たれます。

$$ \begin{equation} L^2 = I_x^2 \omega_x^2 + I_y^2 \omega_y^2 + I_z^2 \omega_z^2 = \text{const.} \end{equation} $$

回転の運動エネルギー $T$ もまた保存されます。

$$ \begin{equation} T = \frac{1}{2}\left(I_x \omega_x^2 + I_y \omega_y^2 + I_z \omega_z^2\right) = \text{const.} \end{equation} $$

この2つの保存量は、$(\omega_x, \omega_y, \omega_z)$ の3次元空間において、それぞれ楕円体の表面を定義します。角運動量楕円体とエネルギー楕円体の交線が、角速度ベクトルの取りうる軌跡（ポインソーの軌跡）です。

軸対称体の解析解

一般の非対称剛体では解析解が楕円関数を用いた複雑な形になりますが、デブリのモデルとして実用的な軸対称体（$I_x = I_y \equiv I_\perp, \, I_z \equiv I_s$）の場合は、きれいな解析解が得られます。

まず、$I_x = I_y = I_\perp$ を代入すると、オイラーの第3式は次のように簡略化されます。

$$ I_s \dot{\omega}_z – (I_\perp – I_\perp)\omega_x \omega_y = I_s \dot{\omega}_z = 0 $$

したがって $\omega_z = \text{const.}$ であり、対称軸まわりの回転成分は一定です。

次に、第1式と第2式に $I_x = I_y = I_\perp$ を代入すると

$$ I_\perp \dot{\omega}_x = (I_\perp – I_s)\omega_y \omega_z $$

$$ I_\perp \dot{\omega}_y = (I_s – I_\perp)\omega_x \omega_z $$

ここでボディノーテーションレート $\Omega$ を次のように定義します。

$$ \begin{equation} \Omega = \frac{I_s – I_\perp}{I_\perp}\omega_z \end{equation} $$

すると、第1式と第2式は次のように書き換えられます。

$$ \dot{\omega}_x = -\Omega \omega_y, \quad \dot{\omega}_y = \Omega \omega_x $$

これは $\omega_x, \omega_y$ が角振動数 $\Omega$ で振動する単振動の連立方程式です。初期条件 $\omega_x(0) = \omega_\perp, \, \omega_y(0) = 0$ とすると、解は次のようになります。

$$ \begin{equation} \omega_x(t) = \omega_\perp \cos(\Omega t), \quad \omega_y(t) = \omega_\perp \sin(\Omega t), \quad \omega_z(t) = \omega_z(0) \end{equation} $$

この解の物理的意味は明確です。角速度ベクトル $\bm{\omega}$ のボディ座標系での先端は、$z$軸まわりを角振動数 $\Omega$ で円錐運動します。これは慣性系から見ると、物体の対称軸が空間に固定された角運動量ベクトルの周りを歳差運動する、いわゆる自由歳差運動（free precession）に対応します。

非対称体の運動 — 中間軸定理

実際のデブリは完全な軸対称ではありません。3つの主慣性モーメントがすべて異なる一般の非対称剛体では、有名な中間軸定理（tennis racket theorem）が成り立ちます。

慣性モーメントを $I_x < I_y < I_z$ と並べたとき、

• $I_x$（最小）の軸まわりの回転 → 安定
• $I_z$（最大）の軸まわりの回転 → 安定
• $I_y$（中間）の軸まわりの回転 → 不安定

テニスラケットを空中に放り上げて回転させると、面に垂直な軸（中間軸）まわりに回すと予測不能な反転が起きますが、グリップ軸（最小慣性モーメント軸）まわりなら安定して回ります。これが中間軸定理の日常的な現れです。

デブリ除去の観点では、この定理は重要な意味を持ちます。時間とともに内部のエネルギー散逸（構造のたわみ、液体のスロッシングなど）により、非対称体は最終的に最大慣性モーメント軸まわりの純回転に近づきます（フラットスピン）。これは回転の運動エネルギーが最小の状態であり、角運動量は保存されつつエネルギーだけが減少する結果です。

大型のロケット上段などでは、制御喪失後数年でフラットスピンに移行していることが多く、捕獲戦略を設計する際にはこの回転モードを考慮する必要があります。

ここまでで、タンブリングデブリの回転運動の数理的な基礎が整いました。回転の性質がわかれば、「どこをどうやって掴むか」の戦略が見えてきます。次に、各種捕獲手法の原理と特徴を詳しく比較していきましょう。

捕獲手法の分類と比較

デブリ捕獲の手法は大きく分けて接触式（contact-based）と非接触式（contactless）に分類されます。非接触式にはイオンビーム・レーザーアブレーション・電磁誘導などがありますが、現在技術的に成熟度が高いのは接触式です。本記事では、主要な接触式捕獲手法を5つ取り上げ、それぞれの原理・長所・短所を解説します。

手法1: ロボットアーム（剛体接触型）

原理

チェイサー衛星に搭載された多関節ロボットアームを伸ばし、ターゲットデブリの突起部分（ランチアダプタリング、ノズル、構造部材など）を直接把持する手法です。ISS上のカナダアーム2（SSRMS）やJEMリモートマニピュレータ（JEMRMS）と基本原理は同じですが、ターゲットが非協力物体である点が根本的に異なります。

制御の要件

タンブリングターゲットをロボットアームで捕獲するには、以下の制御が必要です。

視覚追跡（Visual Tracking）: ターゲットの3次元姿勢と回転率をリアルタイムで推定します。LiDARやステレオカメラを組み合わせ、ターゲット表面の特徴点を追跡します。ターゲットが既知の形状（CADモデルが存在するロケット上段など）であればモデルベース追跡が可能ですが、形状が未知の場合はモデルフリーの点群マッチングが必要になります。

同期運動（Synchronization）: エンドエフェクタの速度と角速度をターゲットの把持点の運動に同期させます。把持点の速度 $\bm{v}_g$ は、ターゲットの重心速度 $\bm{v}_c$ と回転による寄与の和で表されます。

$$ \begin{equation} \bm{v}_g = \bm{v}_c + \bm{\omega} \times \bm{r}_g \end{equation} $$

ここで $\bm{\omega}$ はターゲットの角速度ベクトル、$\bm{r}_g$ は重心から把持点までの位置ベクトルです。回転率が高いほど $\bm{v}_g$ は大きくなり、エンドエフェクタに要求される速度追従性能が厳しくなります。

インピーダンス制御: 把持の瞬間には衝撃力が発生します。この衝撃がチェイサー衛星の姿勢を乱したり、ターゲットを弾き飛ばしたりしないよう、エンドエフェクタのインピーダンス（見かけの質量・ダンパー・バネ特性）を適切に調整する必要があります。目標インピーダンスモデルは

$$ \begin{equation} \bm{M}_d \ddot{\bm{x}} + \bm{D}_d \dot{\bm{x}} + \bm{K}_d \bm{x} = \bm{F}_{ext} \end{equation} $$

と記述されます。$\bm{M}_d, \bm{D}_d, \bm{K}_d$ はそれぞれ目標慣性、目標ダンピング、目標剛性行列です。微小重力下では、地上のように重力が安定化力として働かないため、接触時のダイナミクスが地上とは本質的に異なります。

長所と短所

• 長所: 把持後の拘束が確実。ターゲットの特定部位を精密に掴める。再利用可能
• 短所: ターゲットへの極近接が必要（数m以内）。タンブリング率に上限あり（概ね2〜3 deg/s以下）。制御系の複雑さが高い。把持に失敗した場合のリカバリが難しい

手法2: ネットキャプチャ

原理

チェイサーからネット（網）を射出し、ターゲットを包み込んで拘束する手法です。漁師が投網で魚を獲るのと同じ原理ですが、微小重力環境では網の展開ダイナミクスが地上とは全く異なります。

ネットは通常、四隅に質量体（コーナーマス）を取り付けた正方形または円形で、射出時にコーナーマスの遠心力で展開します。ネットの展開面積はターゲットよりも十分に大きくし、ターゲットの形状やタンブリング状態によらず包絡できるように設計します。

力学的特徴

ネットがターゲットに到達し包み込む過程は、柔軟体の非線形ダイナミクスで記述されます。ネットの各節点に作用する力は

$$ \begin{equation} m_i \ddot{\bm{r}}_i = \sum_{j \in \mathcal{N}(i)} \bm{T}_{ij} + \bm{F}_{contact,i} \end{equation} $$

ここで $m_i$ は節点 $i$ の質量、$\bm{r}_i$ はその位置、$\bm{T}_{ij}$ は隣接節点 $j$ からの張力、$\bm{F}_{contact,i}$ はターゲット表面との接触力です。張力は、ネットの糸が自然長よりも伸びた場合にのみ発生する片側拘束です。

$$ \bm{T}_{ij} = \begin{cases} k(|\bm{r}_j – \bm{r}_i| – L_0)\frac{\bm{r}_j – \bm{r}_i}{|\bm{r}_j – \bm{r}_i|} & (|\bm{r}_j – \bm{r}_i| > L_0) \\ \bm{0} & (|\bm{r}_j – \bm{r}_i| \leq L_0) \end{cases} $$

ここで $k$ は糸の剛性、$L_0$ は自然長です。この片側拘束条件がネットの挙動を複雑にしています。糸は引っ張りには耐えますが、圧縮力は伝達しません。

長所と短所

• 長所: 比較的遠距離（10〜50m程度）から射出可能で、チェイサーの安全性が高い。タンブリング率の制約が緩い。形状不確かさに対してロバスト
• 短所: 一度きりの使い捨て（再利用困難）。包絡後の安定化が難しい。ネットとデブリの結合体の姿勢制御が複雑。ネットの展開失敗リスクがある

手法3: ハープーン

原理

チェイサーからハープーン（銛）を射出し、ターゲットの構造体に突き刺して固定する手法です。原理は単純で、捕鯨のモリと同じです。テザー（紐）で連結されたハープーンをターゲットに打ち込み、テザーを通じて引き寄せます。

設計の要点

ハープーンの設計には以下の考慮が必要です。

貫入力学: ハープーンの先端がターゲットの外壁（アルミ合金、ハニカムパネル等）を貫通し、返しが効くだけの十分な運動エネルギーが必要です。貫入深さ $d$ は射出速度 $v_0$ と材料特性で決まります。

$$ \begin{equation} \frac{1}{2}m_h v_0^2 = \int_0^d F_{resist}(x) \, dx \end{equation} $$

ここで $m_h$ はハープーンの質量、$F_{resist}(x)$ は貫入抵抗力です。射出速度は通常20〜30 m/s程度で、火薬やバネで加速します。

二次デブリの防止: ハープーンが貫通した際に、破片（二次デブリ）が生成されないことが求められます。これはデブリ除去ミッションにとって致命的な矛盾 — デブリを除去しようとして新たなデブリを作ってしまう — を避けるために極めて重要な設計制約です。

テザー管理: 射出後のテザーが絡まないよう、テザーの巻き取りリールと張力制御が必要です。

長所と短所

• 長所: 機構がシンプルで軽量。比較的遠距離（10〜30m）から射出可能。確実な物理的結合。タンブリング率の許容範囲が広い
• 短所: 非可逆（やり直し不可）。二次デブリ生成のリスク。ターゲットの構造強度が不明な場合、貫通しすぎたり弾かれたりする可能性。デブリ内部に推進剤タンクがある場合は爆発リスク

手法4: 接着パッド / ゲッコーグリッパ

原理

ターゲットの表面に接着力で張り付いて捕獲する手法です。2つのアプローチがあります。

ドライ接着（ゲッコーグリッパ）: ヤモリ（gecko）の足裏にヒントを得た技術です。ヤモリの足裏には微細な毛（セタ）が数百万本生えており、各毛の先端がファンデルワールス力でガラスや壁に付着します。この原理を工学的に再現した微細構造（マイクロウェッジアレイ）を用いると、宇宙の真空環境でも劣化しない接着力が得られます。

接着力 $F_{adhesion}$ は、接触面積 $A$ と単位面積あたりの接着強度 $\sigma_{ad}$ の積で表されます。

$$ \begin{equation} F_{adhesion} = \sigma_{ad} \cdot A \end{equation} $$

ゲッコーグリッパの特徴は、剪断方向に引くと接着し、法線方向に引くと容易に剥がれるという異方性にあります。これにより、「接着 → 保持 → 解放」のサイクルが可能になります。

粘着パッド: 粘着性のある材料（エラストマー等）でターゲット表面に貼り付きます。粘着パッドはゲッコーグリッパよりも接着力が大きいですが、繰り返し使用で性能が劣化し、宇宙環境（紫外線、熱サイクル、アウトガス）での長期耐久性に課題があります。

長所と短所

• 長所: ターゲットの形状・材質によらず接着可能。二次デブリの生成がゼロ。構造にダメージを与えない（特にゲッコーグリッパ）。繰り返し使用可能（ゲッコー型）
• 短所: 接触面が汚染されている場合（アウトガス付着物、微小デブリ衝突痕など）に接着力が低下。法線方向の荷重に弱い。大きなトルクには耐えられない。ターゲットへの極近接が必要

手法5: テンタクル / テザー

原理

テンタクル型は、複数の柔軟な腕（テンタクル）でターゲットを包み込むように把持する手法です。タコの腕が獲物を巻き付けるイメージに近く、ターゲットの形状が不規則でも柔軟に適応できます。

テザー型は、ターゲットに接続したテザー（ワイヤー）を通じて力学的に結合する手法です。テザーの長さを数m〜数十mとることで、チェイサーとデブリの間に安全な距離を保てます。

テザーの動力学

テザーで結合された2体系の運動は興味深い力学問題です。チェイサーの質量を $m_c$、デブリの質量を $m_d$、テザーの長さを $l$ とすると、軌道面内でのライブレーション（振り子運動）の角振動数は

$$ \begin{equation} \omega_{lib} = \sqrt{3}\, n \end{equation} $$

ここで $n$ は軌道の平均運動（角速度）です。この結果は、テザー長さや質量比によらないという注目すべき特徴を持っています。これは重力傾斜安定化の基本的な性質であり、テザー結合体のデオービット戦略を設計する際に重要です。

長所と短所

• 長所: 柔軟な適応性（テンタクル型）。安全距離の確保（テザー型）。形状不確かさに対するロバスト性
• 短所: 柔軟体の制御が複雑。テザーの絡まりリスク。テンタクルの把持力に限界。デブリの質量が大きいとテザーに要求される強度が高い

手法の比較表

以上の5つの手法を横断的に比較しましょう。

手法の選定は、ターゲットの特性（質量、回転率、形状）、ミッション要件（コスト、リスク許容度、再利用性）、技術成熟度を総合的に勘案して行います。複数の手法を組み合わせるハイブリッドアプローチ（例: ネットで大まかに捕獲してからロボットアームで安定化）も研究されています。

捕獲手法の全体像が見えたところで、次は「捕獲した後どうするか」— デオービットの戦略について考えましょう。

捕獲後のデオービット戦略

デオービットとは

デブリを捕獲しただけでは問題は解決しません。捕獲したデブリを軌道から除去する — すなわちデオービット（de-orbit）させる必要があります。デオービットの方法は、軌道高度と最終的な処分方法によって異なります。

制御された大気圏再突入

LEOデブリの最も一般的なデオービット方法は、制御された大気圏再突入です。チェイサーの推進系を用いて逆噴射し、近地点高度を大気圏（約120km以下）まで下げます。必要な $\Delta v$ は軌道力学から計算できます。

円軌道からの離脱に必要な $\Delta v$ は、ホーマン遷移の考え方で近似できます。高度 $h$ の円軌道（半径 $r = R_E + h$、$R_E$ は地球半径）から、近地点が $r_p = R_E + h_p$（$h_p \approx 80$ km）の楕円軌道に遷移する場合、

$$ \begin{equation} \Delta v = v_{circ} – v_{transfer} = \sqrt{\frac{\mu}{r}} – \sqrt{\mu\left(\frac{2}{r} – \frac{2}{r + r_p}\right)} \end{equation} $$

ここで $\mu = GM_E \approx 3.986 \times 10^{14} \, \text{m}^3/\text{s}^2$ は地球の重力パラメータです。括弧の中を整理すると

$$ \begin{equation} \Delta v = \sqrt{\frac{\mu}{r}}\left(1 – \sqrt{\frac{2r_p}{r + r_p}}\right) \end{equation} $$

例えば、高度800kmの円軌道（$r = 7,178$ km）からデオービットする場合、$r_p = 6,458$ km（高度80km）として計算すると、$\Delta v \approx 190$ m/s となります。1トンのデブリを結合体として運ぶ場合、チェイサーにはこの $\Delta v$ を出せるだけの推進剤が必要です。

ドラッグ増加デバイス

推進剤を節約するアプローチとして、ドラッグセイルやバルーンなどの大気抵抗増加デバイスがあります。これらは断面積 $A$ を増大させることで、大気抵抗による減速を加速させます。高度 $h$ における大気抵抗の加速度は

$$ \begin{equation} a_{drag} = -\frac{1}{2}\frac{C_D A}{m}\rho(h) v^2 \end{equation} $$

ここで $C_D$ は抗力係数、$\rho(h)$ は高度 $h$ での大気密度、$v$ は軌道速度です。$C_D A / m$ を弾道係数の逆数（$B^{-1}$）と呼びます。ドラッグセイルで断面積を10〜100倍に増大させれば、高度800kmからの自然落下に数百年かかるところを、数年〜数十年に短縮できます。

ただし、ドラッグ増加デバイスは制御された再突入ではないため、地上への落下地点が予測しにくいという問題があります。大型デブリの場合は、再突入時に燃え残る可能性があるため、制御された再突入の方が望ましいです。

グレーブヤード軌道への移動

GEO（静止軌道）近傍のデブリに対しては、大気圏再突入は $\Delta v$ が非常に大きく非現実的です。代わりに、GEOより約300km高いグレーブヤード軌道（墓場軌道）へ移動させます。必要な $\Delta v$ は約11 m/s で、再突入に比べてはるかに小さくなります。

推力方式の選択

デオービットの推力源としては、以下が候補です。

• 化学推進: 高推力で即座に $\Delta v$ を出せるが、推進剤質量が大きい
• 電気推進: 高比推力で推進剤を節約できるが、推力が小さく時間がかかる
• テザー推力（EDT）: 導電性テザーが地球磁場中を運動する際のローレンツ力を利用。推進剤不要だが、テザー長が数km必要

電気推進の場合、比推力 $I_{sp}$ が大きいため、同じ $\Delta v$ をはるかに少ない推進剤で達成できます。ツィオルコフスキーの式

$$ \begin{equation} \Delta v = I_{sp} \, g_0 \, \ln\frac{m_0}{m_f} \end{equation} $$

から、$I_{sp} = 300$ s（化学）と $I_{sp} = 3000$ s（電気）では、同じ $\Delta v$ に対して必要な推進剤比率 $m_0/m_f$ が大きく異なります。数トンクラスのデブリを扱うADRミッションでは、この推進剤の差がミッション全体の質量バジェットを左右します。

ここまでで、デブリの「捕獲」から「除去」までの技術的な流れを一通り見てきました。次に、これらの技術が実際に宇宙で実証されたミッションを紹介します。

主要ミッションの紹介

RemoveDEBRIS（2018年）

RemoveDEBRISは、サリー大学が主導しEUが資金提供した技術実証ミッションで、2018年にISSから放出されました。世界初の軌道上デブリ除去技術実証であり、以下の4つの実験を計画しました。

1. ネット捕獲実験: CubeSatターゲットを放出し、約7m離れた位置からネットを射出して捕獲。成功 — ネットがターゲットを包み込む過程が初めて軌道上で実証された
2. ハープーン実験: ブームの先端に取り付けたターゲットパネルに向けてハープーンを射出。成功 — ハープーンがパネルを貫通・固定
3. ビジョンベースナビゲーション実験: LiDARとカメラによるターゲット追跡技術の検証。成功
4. ドラッグセイル実験: ミッション終了後にドラッグセイルを展開しデオービットを加速。部分的成功 — セイルは展開したが意図どおりの姿勢を維持できなかった

RemoveDEBRISの最大の貢献は、ネットとハープーンという2つの捕獲手法が軌道上で実際に機能することを示した点です。ただし、ターゲットは自ら放出した協力的なCubeSatであり、本物のタンブリングデブリの捕獲とは条件が大きく異なります。

ELSA-d（2021年）

ELSA-d（End-of-Life Services by Astroscale – demonstration）は、日本のAstroscale社が開発したランデブー・近傍運動・ドッキング技術の実証ミッションです。2021年3月にカザフスタンから打ち上げられました。

ELSA-dはチェイサー（約175kg）とクライアント（約17kg）の2機で構成されます。クライアントには磁気ドッキングプレートが搭載されており、チェイサーはこのプレートに磁力で結合します。

ミッションのシーケンスは以下のとおりです。

1. チェイサーがクライアントを放出
2. チェイサーがランデブー・近傍運動を行い、クライアントに接近
3. 磁気ドッキングプレートに結合して捕獲
4. クライアントをタンブリング状態で放出し、タンブリングターゲットの捕獲に挑戦

ELSA-dは、複数回のランデブー・ドッキングと姿勢推定に成功し、軌道上サービスの基本技術を実証しました。ただし、磁気ドッキングプレートを使用しているため、完全な非協力物体とは言えません。Astroscaleは後継ミッションADRAS-Jで実際のロケット上段デブリに接近する実験を行い、2024年に世界で初めて実デブリの詳細な近傍観測に成功しています。

ClearSpace-1（計画中）

ClearSpace-1は、ESAが主導し、スイスのClearSpace社が開発を進めているミッションです。ターゲットは、2013年にVegaロケットの打ち上げで使用されたVESPA（Vega Secondary Payload Adapter）上段アダプタ（約112kg、高度約660km）です。

ClearSpace-1の特徴は、4本のロボットアームでターゲットを包み込むように把持する手法を採用している点です。これは前述のテンタクル型に近い設計思想で、ターゲットの正確な把持点を必要とせず、全体を抱え込むことで形状不確かさに対するロバスト性を確保しています。

捕獲後は、結合体として大気圏再突入しデブリを処分します。打ち上げは2026年以降に予定されており、実現すればESA初の本格的ADRミッションとなります。

ミッションの比較

これらのミッションを見ると、ADR技術が「技術実証」から「実運用」に向かって段階的に進化していることがわかります。初期のRemoveDEBRISでは協力的なターゲットに対する基本的な手法を検証し、ELSA-dではランデブー・ドッキングの自律性を高め、ClearSpace-1ではいよいよ実デブリの捕獲・除去に挑みます。

ここまでの理論的な議論を、Pythonシミュレーションで視覚的に確認してみましょう。タンブリングデブリの回転運動を数値積分し、角速度ベクトルの時間発展とポインソーの軌跡を可視化します。

Pythonによるタンブリングデブリの回転運動シミュレーション

シミュレーション1: 軸対称体の自由歳差運動

まず、軸対称体（$I_x = I_y$）のトルクフリー回転をシミュレーションし、解析解と比較します。解析解の節で導いたように、角速度ベクトルのボディ座標系での先端は、対称軸まわりを円運動するはずです。

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# 軸対称体の慣性モーメント（ロケット上段を想定）
I_perp = 800.0   # kg m^2 (横軸: Ix = Iy)
I_s = 200.0       # kg m^2 (対称軸: Iz) — 細長い円柱

# 初期角速度 [rad/s]
omega_x0 = 0.05   # 約2.9 deg/s
omega_y0 = 0.0
omega_z0 = 0.10   # 約5.7 deg/s

# オイラーの回転方程式（トルクフリー）
def euler_equations(t, omega, Ix, Iy, Iz):
    wx, wy, wz = omega
    dwx = (Iy - Iz) * wy * wz / Ix
    dwy = (Iz - Ix) * wz * wx / Iy
    dwz = (Ix - Iy) * wx * wy / Iz
    return [dwx, dwy, dwz]

# 数値積分
t_span = (0, 120)  # 120秒
t_eval = np.linspace(0, 120, 3000)
sol = solve_ivp(euler_equations, t_span, [omega_x0, omega_y0, omega_z0],
                args=(I_perp, I_perp, I_s), t_eval=t_eval, rtol=1e-12, atol=1e-14)

# 解析解
Omega_body = (I_s - I_perp) / I_perp * omega_z0  # ボディノーテーションレート
omega_perp = omega_x0  # 初期横軸角速度
wx_analytic = omega_perp * np.cos(Omega_body * t_eval)
wy_analytic = omega_perp * np.sin(Omega_body * t_eval)

# 保存量の確認
L2 = I_perp**2 * sol.y[0]**2 + I_perp**2 * sol.y[1]**2 + I_s**2 * sol.y[2]**2
T  = 0.5 * (I_perp * sol.y[0]**2 + I_perp * sol.y[1]**2 + I_s * sol.y[2]**2)

fig, axes = plt.subplots(2, 2, figsize=(12, 9))

# (a) 角速度の時間発展
axes[0, 0].plot(sol.t, np.degrees(sol.y[0]), 'b-', label=r'$\omega_x$ (numerical)', linewidth=1.5)
axes[0, 0].plot(sol.t, np.degrees(sol.y[1]), 'r-', label=r'$\omega_y$ (numerical)', linewidth=1.5)
axes[0, 0].plot(sol.t, np.degrees(sol.y[2]), 'g-', label=r'$\omega_z$ (numerical)', linewidth=1.5)
axes[0, 0].plot(t_eval, np.degrees(wx_analytic), 'b--', alpha=0.5, label=r'$\omega_x$ (analytic)')
axes[0, 0].plot(t_eval, np.degrees(wy_analytic), 'r--', alpha=0.5, label=r'$\omega_y$ (analytic)')
axes[0, 0].set_xlabel('Time [s]')
axes[0, 0].set_ylabel('Angular velocity [deg/s]')
axes[0, 0].set_title('(a) Angular velocity components')
axes[0, 0].legend(fontsize=8)
axes[0, 0].grid(True, alpha=0.3)

# (b) ボディ座標系での角速度ベクトルの軌跡
axes[0, 1].plot(np.degrees(sol.y[0]), np.degrees(sol.y[1]), 'b-', linewidth=0.8)
axes[0, 1].plot(np.degrees(sol.y[0][0]), np.degrees(sol.y[1][0]), 'go', markersize=8, label='Start')
axes[0, 1].set_xlabel(r'$\omega_x$ [deg/s]')
axes[0, 1].set_ylabel(r'$\omega_y$ [deg/s]')
axes[0, 1].set_title(r'(b) $\omega$ trajectory in body frame')
axes[0, 1].set_aspect('equal')
axes[0, 1].legend()
axes[0, 1].grid(True, alpha=0.3)

# (c) 角運動量の保存
axes[1, 0].plot(sol.t, L2 / L2[0] - 1, 'k-', linewidth=1.0)
axes[1, 0].set_xlabel('Time [s]')
axes[1, 0].set_ylabel(r'$(L^2 - L_0^2)/L_0^2$')
axes[1, 0].set_title('(c) Angular momentum conservation')
axes[1, 0].ticklabel_format(style='scientific', axis='y', scilimits=(-3, 3))
axes[1, 0].grid(True, alpha=0.3)

# (d) 回転エネルギーの保存
axes[1, 1].plot(sol.t, T / T[0] - 1, 'k-', linewidth=1.0)
axes[1, 1].set_xlabel('Time [s]')
axes[1, 1].set_ylabel(r'$(T - T_0)/T_0$')
axes[1, 1].set_title('(d) Rotational kinetic energy conservation')
axes[1, 1].ticklabel_format(style='scientific', axis='y', scilimits=(-3, 3))
axes[1, 1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('axisymmetric_tumbling.png', dpi=150, bbox_inches='tight')
plt.show()

上の4つのグラフから、以下のことが読み取れます。

1. (a) 角速度の時間発展: $\omega_x$ と $\omega_y$ が一定振幅で振動し、$\omega_z$ は完全に一定に保たれています。数値解（実線）と解析解（破線）が完全に重なっており、オイラーの方程式の数値積分が正確に行われていることが確認できます。ボディノーテーションレート $\Omega = (I_s – I_\perp)/I_\perp \cdot \omega_z = (200 – 800)/800 \times 0.1 = -0.075$ rad/s に対応する周期は $2\pi / |\Omega| \approx 83.8$ 秒であり、グラフの振動周期と一致しています。
2. (b) 角速度ベクトルの軌跡: ボディ座標系での角速度ベクトルの $(\omega_x, \omega_y)$ 平面への射影が正確な円を描いています。これは解析解で示した「角速度ベクトルが対称軸まわりを円錐運動する」ことの視覚的な確認です。
3. (c)(d) 保存量の検証: 角運動量 $L^2$ と回転エネルギー $T$ の相対誤差が $10^{-12}$ オーダー以下であり、数値計算の精度が十分に高いことがわかります。

シミュレーション2: 非対称体のタンブリングと中間軸不安定性

次に、3つの慣性モーメントがすべて異なる非対称体のタンブリングをシミュレーションします。中間軸定理を視覚的に確認するため、3つの主軸それぞれの近傍に初期角速度を設定した場合を比較します。

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# 非対称体の慣性モーメント（使用済み衛星を想定）
Ix = 400.0    # kg m^2（最小）
Iy = 700.0    # kg m^2（中間）
Iz = 1000.0   # kg m^2（最大）

# オイラーの回転方程式
def euler_eqs(t, omega):
    wx, wy, wz = omega
    dwx = (Iy - Iz) * wy * wz / Ix
    dwy = (Iz - Ix) * wz * wx / Iy
    dwz = (Ix - Iy) * wx * wy / Iz
    return [dwx, dwy, dwz]

omega0 = 0.08  # 主回転 ~4.6 deg/s
eps = 0.005     # 微小摂動 ~0.3 deg/s

# 3つの初期条件: 各軸近傍に回転を設定
cases = {
    r'Near min axis ($I_x$)': [omega0, eps, eps],
    r'Near mid axis ($I_y$)': [eps, omega0, eps],
    r'Near max axis ($I_z$)': [eps, eps, omega0],
}

t_span = (0, 300)
t_eval = np.linspace(0, 300, 5000)

fig, axes = plt.subplots(1, 3, figsize=(15, 5))

for idx, (title, ic) in enumerate(cases.items()):
    sol = solve_ivp(euler_eqs, t_span, ic, t_eval=t_eval, rtol=1e-12, atol=1e-14)

    axes[idx].plot(sol.t, np.degrees(sol.y[0]), 'b-', label=r'$\omega_x$', linewidth=1.0)
    axes[idx].plot(sol.t, np.degrees(sol.y[1]), 'r-', label=r'$\omega_y$', linewidth=1.0)
    axes[idx].plot(sol.t, np.degrees(sol.y[2]), 'g-', label=r'$\omega_z$', linewidth=1.0)
    axes[idx].set_xlabel('Time [s]')
    axes[idx].set_ylabel('Angular velocity [deg/s]')
    axes[idx].set_title(title)
    axes[idx].legend(fontsize=9)
    axes[idx].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('asymmetric_tumbling.png', dpi=150, bbox_inches='tight')
plt.show()

3つのグラフを比較すると、中間軸定理が鮮やかに確認できます。

1. 最小慣性モーメント軸（$I_x$）近傍の回転（左）: $\omega_x$ がほぼ一定の大きな値を維持し、$\omega_y, \omega_z$ は小さな振幅で安定に振動しています。微小摂動が成長せず、回転は安定です。
2. 中間慣性モーメント軸（$I_y$）近傍の回転（中央）: 初期にはほぼ$y$軸まわりの回転だったものが、微小摂動が急激に成長し、3軸すべてに大きな角速度成分が現れます。これが中間軸の不安定性であり、テニスラケット効果です。デブリが「予測不能な反転」を起こしているように見えます。
3. 最大慣性モーメント軸（$I_z$）近傍の回転（右）: 最小軸と同様に安定です。$\omega_z$ が一定を維持し、摂動は成長しません。

この結果はデブリ除去にとって実践的な示唆を与えます。エネルギー散逸により多くのデブリが最大慣性モーメント軸まわりの回転（フラットスピン）に落ち着いている一方で、中間軸近傍の回転状態にあるデブリは姿勢が大きく乱れ、捕獲タイミングの予測が困難になります。

シミュレーション3: ポインソーの楕円体と角速度の軌跡

角速度ベクトルの3次元的な振る舞いを、保存量が定義する楕円体とともに可視化します。

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# 慣性モーメント
Ix, Iy, Iz = 400.0, 700.0, 1000.0

def euler_eqs(t, omega):
    wx, wy, wz = omega
    dwx = (Iy - Iz) * wy * wz / Ix
    dwy = (Iz - Ix) * wz * wx / Iy
    dwz = (Ix - Iy) * wx * wy / Iz
    return [dwx, dwy, dwz]

# 中間軸近傍の初期条件（不安定ケース）
omega_init = [0.005, 0.08, 0.005]

t_span = (0, 600)
t_eval = np.linspace(0, 600, 10000)
sol = solve_ivp(euler_eqs, t_span, omega_init, t_eval=t_eval, rtol=1e-12, atol=1e-14)

# 保存量の計算
L2_val = Ix**2 * sol.y[0][0]**2 + Iy**2 * sol.y[1][0]**2 + Iz**2 * sol.y[2][0]**2
T_val = 0.5 * (Ix * sol.y[0][0]**2 + Iy * sol.y[1][0]**2 + Iz * sol.y[2][0]**2)

# エネルギー楕円体の表面を生成
u = np.linspace(0, 2 * np.pi, 100)
v = np.linspace(0, np.pi, 50)
# T = 0.5 * (Ix*wx^2 + Iy*wy^2 + Iz*wz^2) → 楕円体
a_T = np.sqrt(2 * T_val / Ix)
b_T = np.sqrt(2 * T_val / Iy)
c_T = np.sqrt(2 * T_val / Iz)
x_ell = a_T * np.outer(np.cos(u), np.sin(v))
y_ell = b_T * np.outer(np.sin(u), np.sin(v))
z_ell = c_T * np.outer(np.ones_like(u), np.cos(v))

fig = plt.figure(figsize=(10, 8))
ax = fig.add_subplot(111, projection='3d')

# エネルギー楕円体（半透明）
ax.plot_surface(x_ell, y_ell, z_ell, alpha=0.15, color='cyan', edgecolor='none')

# 角速度の軌跡
ax.plot(sol.y[0], sol.y[1], sol.y[2], 'r-', linewidth=0.8, label=r'$\omega(t)$ trajectory')
ax.plot([sol.y[0][0]], [sol.y[1][0]], [sol.y[2][0]], 'go', markersize=8, label='Start')

ax.set_xlabel(r'$\omega_x$ [rad/s]')
ax.set_ylabel(r'$\omega_y$ [rad/s]')
ax.set_zlabel(r'$\omega_z$ [rad/s]')
ax.set_title("Poinsot's construction: angular velocity on energy ellipsoid")
ax.legend()

plt.tight_layout()
plt.savefig('poinsot_ellipsoid.png', dpi=150, bbox_inches='tight')
plt.show()

この3次元プロットから、角速度ベクトルがエネルギー楕円体の表面上を動いている様子がわかります。角運動量楕円体との交線に沿って角速度の先端が軌道を描き、中間軸（$\omega_y$ 軸）近傍では大きく揺れ動く不安定な挙動を示しています。この可視化は、ポインソーが19世紀に導いた幾何学的解釈を数値的に裏付けるものです。

デブリの捕獲計画を立てる際には、この角速度ベクトルの軌跡から「デブリの表面上の特定点がどのような速度で動くか」を予測し、ロボットアームやネットの射出タイミングを決定することになります。

シミュレーション4: デブリ表面の把持点の速度プロファイル

最後に、タンブリングデブリの表面にある把持候補点の速度を計算し、ロボットアームの追従要件を定量的に評価します。

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# 慣性モーメント（ロケット上段、非対称）
Ix, Iy, Iz = 500.0, 750.0, 900.0

def euler_eqs(t, omega):
    wx, wy, wz = omega
    dwx = (Iy - Iz) * wy * wz / Ix
    dwy = (Iz - Ix) * wz * wx / Iy
    dwz = (Ix - Iy) * wx * wy / Iz
    return [dwx, dwy, dwz]

# 典型的なタンブリングデブリの初期角速度（約3 deg/s）
omega_init = [0.03, 0.04, 0.02]  # rad/s

t_span = (0, 200)
t_eval = np.linspace(0, 200, 5000)
sol = solve_ivp(euler_eqs, t_span, omega_init, t_eval=t_eval, rtol=1e-12, atol=1e-14)

# 把持候補点のボディ座標系での位置（重心からの距離）
# ロケット上段のランチアダプタリング（端部）を想定
r_grip = np.array([0.5, 0.0, 2.5])  # [m] — 半径0.5m, 長さ方向2.5m

# 各時刻での把持点の速度（ω × r）
v_grip = np.zeros((3, len(t_eval)))
for i in range(len(t_eval)):
    omega_vec = sol.y[:, i]
    v_grip[:, i] = np.cross(omega_vec, r_grip)

# 把持点の速度の大きさ
v_mag = np.linalg.norm(v_grip, axis=0)

fig, axes = plt.subplots(2, 1, figsize=(12, 8))

# (a) 把持点の速度成分
axes[0].plot(sol.t, v_grip[0] * 100, 'b-', label=r'$v_x$', linewidth=1.0)
axes[0].plot(sol.t, v_grip[1] * 100, 'r-', label=r'$v_y$', linewidth=1.0)
axes[0].plot(sol.t, v_grip[2] * 100, 'g-', label=r'$v_z$', linewidth=1.0)
axes[0].set_xlabel('Time [s]')
axes[0].set_ylabel('Grip point velocity [cm/s]')
axes[0].set_title('(a) Velocity components of grip point on tumbling debris')
axes[0].legend()
axes[0].grid(True, alpha=0.3)

# (b) 速度の大きさとロボットアームの追従限界
axes[1].plot(sol.t, v_mag * 100, 'k-', linewidth=1.5, label='Grip point speed')
axes[1].axhline(y=5.0, color='orange', linestyle='--', linewidth=1.5, label='Arm tracking limit (5 cm/s)')
axes[1].axhline(y=10.0, color='red', linestyle='--', linewidth=1.5, label='Arm tracking limit (10 cm/s)')
axes[1].fill_between(sol.t, 0, 5.0, alpha=0.1, color='green', label='Safe capture zone')
axes[1].set_xlabel('Time [s]')
axes[1].set_ylabel('Grip point speed [cm/s]')
axes[1].set_title('(b) Grip point speed vs. robot arm tracking capability')
axes[1].legend(fontsize=9)
axes[1].grid(True, alpha=0.3)
axes[1].set_ylim(0, max(v_mag * 100) * 1.3)

plt.tight_layout()
plt.savefig('grip_point_velocity.png', dpi=150, bbox_inches='tight')
plt.show()

# 数値情報の出力
print(f"把持点の平均速度: {np.mean(v_mag)*100:.2f} cm/s")
print(f"把持点の最大速度: {np.max(v_mag)*100:.2f} cm/s")
print(f"把持点の最小速度: {np.min(v_mag)*100:.2f} cm/s")
print(f"初期角速度の大きさ: {np.degrees(np.linalg.norm(omega_init)):.2f} deg/s")

上のグラフから、デブリ捕獲の実践的な課題が定量的に見えてきます。

1. (a) 速度成分: 把持候補点の速度は各成分が複雑に振動しており、単純な周期運動ではありません。これは非対称体のタンブリングが3軸間で結合していることの直接的な反映です。ロボットアームの追従制御には、この非周期的な変動に対応できる高帯域の制御系が必要です。
2. (b) 速度の大きさとアーム追従限界: ロボットアームのエンドエフェクタが追従可能な速度には技術的な上限があります。現在の宇宙用マニピュレータの追従能力は5〜10 cm/s程度とされています。グラフで把持点の速度がこの限界を超える時間帯があれば、その時間帯での捕獲は不可能です。逆に、速度が限界以下になる「捕獲ウィンドウ」を特定し、そのタイミングで把持動作を実行する戦略が有効です。
3. 重心から離れた把持点ほど、$|\bm{\omega} \times \bm{r}|$ が大きくなるため速度が増します。把持点をなるべく重心に近い位置に選ぶことが、追従要件を緩和するための基本的な設計指針です。

まとめ

本記事では、スペースデブリ除去のロボティクス技術について、問題の背景から具体的な捕獲戦略、デオービットまでを体系的に解説しました。

• ケスラーシンドロームの脅威により、年間5〜10個の大型デブリの能動的除去（ADR）が不可欠である
• 非協力物体は「通信不能・制御不能・把持点なし・タンブリング」という4重の困難を捕獲ロボットに突きつける
• タンブリングデブリの回転運動はオイラーの方程式で記述され、中間軸定理により回転の安定性が慣性モーメントの大小関係で決まる
• 捕獲手法はロボットアーム（精密だが近接が必要）、ネット（ロバストだが使い捨て）、ハープーン（確実だが二次デブリリスク）、接着パッド（非破壊だが接着力に限界）、テザー（安全距離を確保できるが制御が複雑）と、それぞれ一長一短がある
• 捕獲後のデオービットでは、軌道力学に基づく $\Delta v$ の計算と推進方式の選定が重要である
• RemoveDEBRIS、ELSA-d、ClearSpace-1と、技術実証は着実に進展している
• Pythonシミュレーションにより、タンブリングの非線形性、中間軸不安定性、把持点速度の定量的評価を視覚的に確認した

スペースデブリ除去は、宇宙ロボティクスの中でも最も挑戦的な応用分野の一つです。本記事で扱ったオイラーの方程式やインピーダンス制御は、ロボット工学と軌道力学が交差する領域の基礎をなしています。

次のステップとして、以下の記事も参考にしてください。

• 微小重力下の接触ダイナミクス — 本記事の前提となる接触力学の詳細
• 遠隔操作と時間遅延補償 — 地上からデブリ捕獲ロボットを操作する際の通信遅延問題