人工衛星が地表の特定地点を正確に撮影したり、通信アンテナを地上局に向けたりするには、衛星の姿勢を精密に制御する必要があります。この精密制御を担うアクチュエータの中でも、最も広く使われているのがリアクションホイール（Reaction Wheel, RW）です。
リアクションホイールの原理は、実は私たちの日常にも潜んでいます。回転椅子に座って両手を円運動させると、身体は逆方向に回転します。両手の角運動量と身体の角運動量の和が保存されるためです。リアクションホイールはこの「角運動量交換」をモーターとフライホイールで実現した装置であり、推進剤を消費せずに衛星の姿勢を変えることができます。
リアクションホイールを理解すると、以下のような実践的な課題に取り組めるようになります。
- 衛星姿勢制御系の設計: RWの台数、配置、容量を合理的に選定できます
- トルク分配問題: 冗長構成において制御トルクを各ホイールに最適に分配するアルゴリズムを設計できます
- 故障耐性設計: ホイール故障時のリコンフィギュレーション（制御の再構成）を検討できます
本記事の内容
- リアクションホイールの動作原理
- トルク生成の数学的定式化
- 3軸制御の最小構成と冗長構成（3+1）
- ピラミッド配置と疑似逆行列による最適トルク分配
- ホイール飽和の問題と対策
- モーメンタムホイール・CMGとの比較
- Pythonによるトルク分配シミュレーション
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
リアクションホイールの動作原理
基本構造
リアクションホイールは、衛星構体に固定されたモーターと、そのモーター軸に取り付けられたフライホイール（慣性質量を持つ回転体）から構成されます。
モーターに電流を流すと、フライホイールにトルクが加わり回転速度が変化します。ニュートンの第3法則（作用・反作用の法則）により、フライホイールに加えたトルクと等大逆向きのトルクが衛星本体に伝わります。
イメージとしては、氷の上に立った人がフリスビーを手で回すと、反作用で身体が逆方向に回転するのと同じです。フリスビーの回転速度を変えるたびに、身体の回転が変化します。
角運動量交換の原理
リアクションホイール1基のスピン軸方向の単位ベクトルを $\hat{\bm{a}}$ とします。ホイールの軸まわり慣性モーメントを $I_w$、ホイールの角速度を $\Omega$ とすると、ホイールの角運動量は
$$ \bm{H}_w = I_w \Omega \hat{\bm{a}} $$
です。モーターがホイールに加えるトルクは $\tau = I_w \dot{\Omega}$ であり、反作用として衛星本体に伝わるトルクは
$$ \bm{N}_{\text{RW}} = -I_w \dot{\Omega} \hat{\bm{a}} $$
です。マイナス符号は反作用（ホイールを加速すると本体には逆方向のトルク）を表しています。
ここで注目すべきは、RWが生成するトルクはホイール回転速度の変化率 $\dot{\Omega}$ に比例するという点です。ホイールが一定速度で回転しているだけではトルクは生じません。トルクを発生させるには、ホイールの回転速度を変化させる必要があります。
制御可能なトルクの方向
1基のRWが生成できるトルクの方向は、スピン軸方向 $\hat{\bm{a}}$ に限定されます（正負両方向）。したがって、衛星の3軸を独立に制御するためには、最低3基のRWが必要であり、そのスピン軸が線形独立でなければなりません。
この制約がRWの配置設計の核心です。どのようにホイールを配置すればよいかを、次のセクションで数学的に定式化します。
トルク生成の数学的定式化
複数ホイールのトルク
$n$ 基のRWを搭載した衛星を考えます。各ホイール $i$ のスピン軸方向を $\hat{\bm{a}}_i$（3次元の単位ベクトル）、ホイールトルクを $\tau_i = I_{w,i}\dot{\Omega}_i$ とすると、衛星本体に作用する合計トルクは
$$ \bm{N}_{\text{ctrl}} = -\sum_{i=1}^{n} \tau_i \hat{\bm{a}}_i $$
これを行列形式で書くと、
$$ \begin{equation} \bm{N}_{\text{ctrl}} = -\bm{A} \bm{\tau} \end{equation} $$
ここで $\bm{A}$ はホイール配置行列（$3 \times n$）、$\bm{\tau} = (\tau_1, \tau_2, \dots, \tau_n)^T$ はホイールトルクベクトルです。
$$ \bm{A} = [\hat{\bm{a}}_1 \quad \hat{\bm{a}}_2 \quad \cdots \quad \hat{\bm{a}}_n] $$
3軸制御の条件
所望のトルク $\bm{N}_{\text{cmd}}$ を実現するには、
$$ \bm{A}\bm{\tau} = -\bm{N}_{\text{cmd}} $$
を解く必要があります。この方程式が任意の $\bm{N}_{\text{cmd}}$ に対して解を持つ条件は、$\bm{A}$ の行ランクが3であること、すなわち $\text{rank}(\bm{A}) = 3$ です。
幾何学的には、$n$ 基のホイールのスピン軸が3次元空間を張ること、つまり全てのスピン軸が1つの平面内に収まらないことが必要です。
ホイール台数と自由度
| ホイール数 $n$ | 方程式の性質 | 解の特性 |
|---|---|---|
| $n = 3$ | 正方行列（$\bm{A}$ が $3 \times 3$） | 解は一意（故障耐性なし） |
| $n = 4$ | 劣決定系（$3 \times 4$） | 解は無限個（1自由度の冗長性） |
| $n > 4$ | 劣決定系 | 解は無限個（$n-3$ 自由度の冗長性） |
| $n < 3$ | 優決定系 | 一般に解なし（3軸制御不可） |
$n = 3$ は最小構成であり、どのホイール1基が故障しても3軸制御を失います。$n = 4$（3+1冗長）は、任意の1基が故障しても残り3基で3軸制御を維持できる構成として広く採用されています。
では、冗長構成（$n > 3$）の場合に解が無限個存在するとき、どの解を選ぶべきでしょうか。これがトルク分配問題です。
ピラミッド配置
配置の基本的な考え方
RWの配置を設計する際の要求は、以下の3点です。
- 3軸制御可能性: $\text{rank}(\bm{A}) = 3$ であること
- 等方性: 各軸方向にバランスよくトルクが生成できること
- 故障耐性: 1基が故障しても3軸制御が維持できること
ピラミッド配置
4輪冗長構成で最も一般的な配置がピラミッド配置です。4基のRWを正四面体の各面の法線方向に配置します。
各ホイールのスピン軸は、衛星の機体座標系（$x, y, z$）に対して以下のように定義されます。
$$ \hat{\bm{a}}_1 = \begin{pmatrix} -\cos\beta \sin\alpha \\ \cos\beta \cos\alpha \\ \sin\beta \end{pmatrix}, \quad \hat{\bm{a}}_2 = \begin{pmatrix} \cos\beta \cos\alpha \\ \cos\beta \sin\alpha \\ \sin\beta \end{pmatrix} $$
$$ \hat{\bm{a}}_3 = \begin{pmatrix} \cos\beta \sin\alpha \\ -\cos\beta \cos\alpha \\ \sin\beta \end{pmatrix}, \quad \hat{\bm{a}}_4 = \begin{pmatrix} -\cos\beta \cos\alpha \\ -\cos\beta \sin\alpha \\ \sin\beta \end{pmatrix} $$
ここで $\beta$ は天頂からの傾斜角（通常 $\beta \approx 35.26°$、つまり $\tan\beta = 1/\sqrt{2}$）、$\alpha = 45°$ です。
$\beta = \arctan(1/\sqrt{2}) \approx 35.26°$ のとき、4基のスピン軸がちょうど正四面体の頂点と中心を結ぶ方向に並び、トルクの等方性が最大になります。この角度では $\sin\beta = 1/\sqrt{3}$, $\cos\beta = \sqrt{2/3}$ です。
配置行列の具体形
$\alpha = 45°$, $\beta = \arctan(1/\sqrt{2})$ のピラミッド配置では、簡略化すると配置行列は
$$ \bm{A} = \frac{1}{\sqrt{3}}\begin{pmatrix} -1 & 1 & 1 & -1 \\ 1 & 1 & -1 & -1 \\ 1 & 1 & 1 & 1 \end{pmatrix} $$
ピラミッド配置の重要な性質として、$\bm{A}\bm{A}^T = \frac{4}{3}\bm{I}_3$ が成り立ちます。つまり各軸に等方的にトルクが分配される理想的な配置です。
また、任意の1列を除いた $3 \times 3$ 部分行列のランクが常に3であることから、1基が故障しても3軸制御が可能です。
次に、この冗長系でトルクを最適に分配する方法を解説します。
疑似逆行列によるトルク分配
トルク分配問題
冗長構成（$n > 3$）では、制御コマンド $\bm{N}_{\text{cmd}}$ を実現するホイールトルク $\bm{\tau}$ は一意ではありません。
$$ \bm{A}\bm{\tau} = -\bm{N}_{\text{cmd}} $$
この劣決定系の解は無限個あります。例えば特解 $\bm{\tau}_0$ に $\bm{A}$ のヌル空間のベクトルを加えた $\bm{\tau}_0 + \bm{A}_{\perp}\bm{z}$（$\bm{z}$ は任意）も解です。
どの解を選ぶかは最適化の問題であり、最も一般的な基準はホイールトルクの2乗和を最小化する（各ホイールの負担を均等にする）ことです。
最小ノルム解と疑似逆行列
$$ \min_{\bm{\tau}} \|\bm{\tau}\|^2 \quad \text{s.t.} \quad \bm{A}\bm{\tau} = -\bm{N}_{\text{cmd}} $$
この問題の解は、ムーア・ペンローズの疑似逆行列（pseudo-inverse）$\bm{A}^+$ を用いて
$$ \begin{equation} \bm{\tau} = -\bm{A}^+ \bm{N}_{\text{cmd}} \end{equation} $$
$$ \bm{A}^+ = \bm{A}^T (\bm{A}\bm{A}^T)^{-1} $$
で与えられます。
ピラミッド配置（$\bm{A}\bm{A}^T = \frac{4}{3}\bm{I}_3$）の場合、疑似逆行列は
$$ \bm{A}^+ = \frac{3}{4}\bm{A}^T $$
と非常にシンプルな形になります。各ホイールの担当トルクは、コマンドトルクの各ホイール方向への射影を均等にスケーリングした値になります。
疑似逆行列の物理的意味
疑似逆行列による分配は、「4基のホイールが均等に負担を分け合う」ことを意味します。これにはいくつかの利点があります。
- 各ホイールのトルクが最小: 1基に過大な負荷がかかることを防ぎます
- 電力消費の均等化: モーターの発熱がバランスされ、熱設計上有利です
- 飽和回避: 特定のホイールだけが先に飽和することを防ぎます
ただし、ホイールが飽和に近い場合は、飽和したホイールの負担を他に再分配する「重み付き疑似逆行列」や最適化ベースの分配法が必要になります。
トルク分配の理論を理解したところで、ホイール飽和の問題について詳しく見ていきましょう。
ホイール飽和の問題
2種類の飽和
リアクションホイールには2種類の「飽和」があります。
1. トルク飽和: モーターが出力できる最大トルクの限界
$$ |\tau_i| \leq \tau_{\max} $$
典型値: $\tau_{\max} = 0.001$ ～ $1$ Nm（ホイールサイズによる）
トルク飽和はマヌーバ中の過渡的な問題であり、マヌーバレートを制限するか、フィードフォワード制御で対処します。
2. 角運動量飽和（速度飽和）: ホイールの最大回転速度の限界
$$ |H_{w,i}| = |I_w \Omega_i| \leq H_{\max} $$
典型値: $H_{\max} = 0.01$ ～ $100$ Nms（ホイールサイズによる）
角運動量飽和は、外乱トルクの蓄積によって生じる長期的な問題です。飽和したホイールは有効なトルクを生成できなくなり、姿勢制御が不能になります。
飽和への対策
角運動量飽和への対策は、前記事で解説したアンローディング（磁気トルカやスラスタによる角運動量排出）です。
トルク飽和への対策としては、
- マヌーバの計画的な制限: 最大角加速度をRWのトルク能力以内に制限する
- フィードフォワード制御: 計画された軌道に沿った理想トルクを事前計算し、フィードバック制御の負担を減らす
- トルク制限付き制御則: 制御トルクが飽和しないようにクリッピングする
$$ \tau_i = \text{clip}(\tau_{i,\text{cmd}}, -\tau_{\max}, \tau_{\max}) $$
冗長構成の故障耐性
4輪ピラミッド配置で1基が故障した場合、残り3基で3軸制御を維持できます。ただし、トルク分配行列を故障前の4輪用から故障後の3輪用に切り替える必要があります。
故障ホイール $j$ を除外した配置行列 $\bm{A}_j$ は $3 \times 3$ の正方行列となり、トルク分配は単純な逆行列で求まります。
$$ \bm{\tau}_j = -\bm{A}_j^{-1} \bm{N}_{\text{cmd}} $$
故障後はトルク等方性が失われるため、特定方向のトルク能力が低下する可能性がありますが、多くのミッション要求は満たせます。
ここまでRWの原理と設計を詳しく見てきましたが、RWとよく比較される他の角運動量交換装置についても整理しておきましょう。
モーメンタムホイール・CMGとの比較
モーメンタムホイール（MW）
モーメンタムホイールは、一定の高い回転速度でフライホイールを回転させ、大きな角運動量バイアスを維持する装置です。
$$ H_{\text{bias}} = I_w \Omega_0 \gg 0 $$
RWとの違いは運用の仕方にあります。MWは名目回転速度が高く（数千rpm）、その速度をわずかに変化させて微小トルクを生成します。一方、RWはゼロ速度を中心に正負両方向に回転します。
MWはモーメンタムバイアス方式の姿勢制御で使用され、ジャイロスコピック安定性を提供します。
コントロールモーメントジャイロ（CMG）
CMGは、一定速度で回転するロータの角運動量ベクトルの方向をジンバルで変化させることでトルクを発生させます。
$$ \bm{N}_{\text{CMG}} = \frac{d\bm{H}_{\text{rotor}}}{dt} = \frac{d(H_0 \hat{\bm{g}})}{dt} = H_0 \dot{\hat{\bm{g}}} $$
ここで $H_0 = I_w \Omega_0$ は一定のロータ角運動量、$\hat{\bm{g}}$ はジンバル角によって変化する角運動量方向です。
CMGの最大の特徴はトルク増幅です。小さなジンバルトルクで大きな姿勢制御トルクを生成できます。
3方式の比較
| 特性 | RW | MW | CMG |
|---|---|---|---|
| 名目回転速度 | ゼロ付近 | 高速（一定） | 高速（一定） |
| トルク生成 | ホイール加速 | ホイール微小速度変化 | ジンバル回転 |
| トルク能力 | 小～中 | 小 | 大 |
| 制御精度 | 高 | 中 | 中～高 |
| 複雑さ | 低 | 低 | 高（ジンバル機構） |
| 特異点問題 | なし | なし | あり |
| 典型的な用途 | 地球観測衛星、科学衛星 | 通信衛星 | 大型宇宙ステーション、アジャイル衛星 |
RWは構造が単純で信頼性が高く、トルク能力は中程度ですが大部分の中小型衛星には十分です。一方、宇宙ステーションのような大型構造物や、短時間で大角度マヌーバを行うアジャイル衛星にはCMGの大きなトルク能力が必要です。
理論的な解説を終えたところで、Pythonでピラミッド配置のトルク分配と姿勢制御シミュレーションを実装しましょう。
Pythonによるトルク分配と姿勢制御シミュレーション
ピラミッド配置の構築
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# ピラミッド配置の定義
beta = np.arctan(1 / np.sqrt(2)) # 傾斜角 ≈ 35.26°
alpha = np.radians(45)
# 4基のスピン軸方向
a1 = np.array([-np.cos(beta)*np.sin(alpha), np.cos(beta)*np.cos(alpha), np.sin(beta)])
a2 = np.array([np.cos(beta)*np.cos(alpha), np.cos(beta)*np.sin(alpha), np.sin(beta)])
a3 = np.array([np.cos(beta)*np.sin(alpha), -np.cos(beta)*np.cos(alpha), np.sin(beta)])
a4 = np.array([-np.cos(beta)*np.cos(alpha), -np.cos(beta)*np.sin(alpha), np.sin(beta)])
A = np.column_stack([a1, a2, a3, a4])
print("配置行列 A:")
print(np.round(A, 4))
print(f"\nA * A^T =")
print(np.round(A @ A.T, 4))
print(f"\nrank(A) = {np.linalg.matrix_rank(A)}")
print(f"\nbeta = {np.degrees(beta):.2f} deg")
# 疑似逆行列
A_pinv = A.T @ np.linalg.inv(A @ A.T)
print(f"\n疑似逆行列 A+:")
print(np.round(A_pinv, 4))
# 3/4 * A^T との比較
print(f"\n3/4 * A^T:")
print(np.round(0.75 * A.T, 4))
print(f"\nA+ と 3/4*A^T の差のノルム: {np.linalg.norm(A_pinv - 0.75*A.T):.2e}")
# ピラミッド配置の3D可視化
fig = plt.figure(figsize=(10, 8))
ax = fig.add_subplot(111, projection='3d')
# 各ホイールのスピン軸を矢印で描画
axes_list = [a1, a2, a3, a4]
colors_rw = ['#e74c3c', '#2ecc71', '#3498db', '#f39c12']
labels_rw = ['RW1', 'RW2', 'RW3', 'RW4']
for a, color, label in zip(axes_list, colors_rw, labels_rw):
ax.quiver(0, 0, 0, a[0], a[1], a[2], color=color, linewidth=2.5,
arrow_length_ratio=0.1, label=label)
# 座標軸
for i, (axis_label, color) in enumerate(zip(['X', 'Y', 'Z'], ['gray', 'gray', 'gray'])):
vec = np.zeros(3)
vec[i] = 1.0
ax.quiver(0, 0, 0, vec[0], vec[1], vec[2], color=color, linewidth=1,
arrow_length_ratio=0.1, linestyle='--', alpha=0.5)
ax.text(vec[0]*1.1, vec[1]*1.1, vec[2]*1.1, axis_label, fontsize=12, color='gray')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('Pyramid configuration of 4 reaction wheels')
ax.legend(loc='upper left')
ax.set_xlim([-1, 1])
ax.set_ylim([-1, 1])
ax.set_zlim([-0.2, 1])
plt.tight_layout()
plt.savefig('pyramid_configuration.png', dpi=150, bbox_inches='tight')
plt.show()
3D可視化から、ピラミッド配置の幾何学的な特徴が確認できます。4基のRWスピン軸は全て上方（+Z方向）に傾いており、$XY$ 平面への射影が90度間隔で対称的に配置されています。この配置により、$X, Y, Z$ の各方向に対して均等なトルク生成能力が得られます。$\bm{A}\bm{A}^T = \frac{4}{3}\bm{I}_3$ が成り立つことが数値的にも確認でき、完全な等方性を持つ配置であることがわかります。
トルク分配の可視化
さまざまな方向のコマンドトルクに対して、4基のホイールにどのようにトルクが分配されるかを可視化します。
# さまざまなコマンドトルクに対するトルク分配
fig, axes = plt.subplots(2, 3, figsize=(15, 8))
test_torques = [
([1, 0, 0], '+X torque'),
([0, 1, 0], '+Y torque'),
([0, 0, 1], '+Z torque'),
([1, 1, 0], 'X+Y torque'),
([1, 1, 1], 'X+Y+Z torque'),
([1, -1, 0.5], 'Arbitrary torque'),
]
for ax, (N_cmd, title) in zip(axes.flat, test_torques):
N_cmd = np.array(N_cmd, dtype=float)
N_cmd = N_cmd / np.linalg.norm(N_cmd) # 正規化
tau = -A_pinv @ N_cmd # トルク分配
bars = ax.bar(labels_rw, tau, color=colors_rw, alpha=0.8, edgecolor='black')
ax.set_title(title, fontsize=11)
ax.set_ylabel('Wheel torque')
ax.grid(True, alpha=0.3, axis='y')
ax.axhline(y=0, color='black', linewidth=0.5)
# 検証: A @ tau が -N_cmd と一致するか
N_achieved = A @ tau
error = np.linalg.norm(N_achieved + N_cmd)
ax.text(0.02, 0.98, f'Error: {error:.2e}', transform=ax.transAxes,
fontsize=8, va='top')
plt.suptitle('Torque distribution for various command torques (pyramid config)', fontsize=13)
plt.tight_layout()
plt.savefig('torque_distribution.png', dpi=150, bbox_inches='tight')
plt.show()
トルク分配の結果から、ピラミッド配置の特性が読み取れます。
- 単軸トルク（+X, +Y, +Z）: 4基全てがトルクを分担していますが、要求トルク方向に対して各ホイールの寄与が異なります。疑似逆行列による分配が各ホイールの負荷を均等化しています。
- 複合トルク（X+Y, X+Y+Z）: 複数軸のトルクが要求される場合でも、各ホイールへの分配は滑らかです。
- 分配の検証: 全てのケースで再構成誤差（Error）が $10^{-16}$ 以下であり、要求トルクが正確に実現されていることが確認できます。
トルクエンベロープの可視化
ホイールのトルク制限を考慮したとき、ピラミッド配置で発生可能な最大トルクの「包絡面（エンベロープ）」を可視化します。
# トルクエンベロープの計算
tau_max = 0.1 # 各ホイールの最大トルク [Nm]
# 球面上の多数の方向に対して最大トルクを計算
n_points = 5000
phi = np.random.uniform(0, 2*np.pi, n_points)
cos_theta = np.random.uniform(-1, 1, n_points)
sin_theta = np.sqrt(1 - cos_theta**2)
directions = np.column_stack([
sin_theta * np.cos(phi),
sin_theta * np.sin(phi),
cos_theta
])
max_torques = np.zeros(n_points)
for k in range(n_points):
d = directions[k]
tau = -A_pinv @ d # 単位方向のトルク分配
# 最もきついホイールのスケーリング
scale = tau_max / np.max(np.abs(tau)) if np.max(np.abs(tau)) > 0 else 0
max_torques[k] = scale
# 最大トルクエンベロープ
envelope = directions * max_torques[:, None]
fig = plt.figure(figsize=(10, 8))
ax = fig.add_subplot(111, projection='3d')
sc = ax.scatter(envelope[:, 0], envelope[:, 1], envelope[:, 2],
c=max_torques, cmap='viridis', s=2, alpha=0.6)
plt.colorbar(sc, label='Max torque magnitude [Nm]', shrink=0.6)
ax.set_xlabel('Nx [Nm]')
ax.set_ylabel('Ny [Nm]')
ax.set_zlabel('Nz [Nm]')
ax.set_title(f'Torque envelope (pyramid config, tau_max = {tau_max} Nm)')
plt.tight_layout()
plt.savefig('torque_envelope.png', dpi=150, bbox_inches='tight')
plt.show()
# 各軸方向の最大トルク
for i, axis in enumerate(['X', 'Y', 'Z']):
d = np.zeros(3)
d[i] = 1.0
tau = -A_pinv @ d
scale = tau_max / np.max(np.abs(tau))
print(f"{axis}軸方向の最大トルク: {scale:.4f} Nm")
トルクエンベロープの可視化から、ピラミッド配置のトルク能力が方向によってほぼ均一であることがわかります。完全な球ではありませんが、ほぼ球に近い形状であり、等方的な配置の効果が表れています。各軸方向の最大トルクもほぼ等しく、特定方向で極端にトルクが弱くなるということがありません。
4輪冗長構成の姿勢制御シミュレーション
最後に、ピラミッド配置の4基RWを使った三軸姿勢制御の完全なシミュレーションを行います。
from scipy.integrate import solve_ivp
# 衛星パラメータ
Ix, Iy, Iz = 20.0, 25.0, 15.0
I_diag = np.array([Ix, Iy, Iz])
Iw = 0.05 # 各ホイールの慣性モーメント [kg m^2]
omega_w_max = 6000 * 2 * np.pi / 60 # 最大ホイール速度 [rad/s]
tau_max_rw = 0.05 # 最大ホイールトルク [Nm]
# PD制御ゲイン
omega_n = np.array([0.05, 0.04, 0.06])
zeta = np.array([0.7, 0.7, 0.7])
Kp = I_diag * omega_n**2
Kd = 2 * zeta * omega_n * I_diag
# 外乱トルク
N_dist = np.array([3e-5, 5e-5, 2e-5])
def satellite_rw_dynamics(t, state):
"""4基RW付き衛星の姿勢ダイナミクス"""
theta = state[:3] # 姿勢角
omega = state[3:6] # 衛星角速度
Omega_w = state[6:10] # 4基のホイール角速度
# PD制御: コマンドトルク
N_cmd = -Kp * theta - Kd * omega
# トルク分配（疑似逆行列）
tau_cmd = -A_pinv @ N_cmd
# トルク飽和
tau = np.clip(tau_cmd, -tau_max_rw, tau_max_rw)
# ホイール速度飽和チェック
for i in range(4):
if abs(Omega_w[i]) >= omega_w_max and np.sign(tau[i]) == np.sign(Omega_w[i]):
tau[i] = 0 # 飽和方向にはトルクを出せない
# 衛星本体に作用するトルク
N_rw = -A @ tau # ホイールからの反作用トルク
N_total = N_rw + N_dist
# 角加速度（線形化モデル）
alpha = N_total / I_diag
# ホイール角加速度
dOmega_w = tau / Iw
return np.concatenate([omega, alpha, dOmega_w])
# 初期条件: 各軸15度のマヌーバ
theta0 = np.radians([15.0, 10.0, 20.0])
state0 = np.concatenate([theta0, np.zeros(3), np.zeros(4)])
t_span = (0, 600)
t_eval = np.linspace(0, 600, 6000)
sol = solve_ivp(satellite_rw_dynamics, t_span, state0,
t_eval=t_eval, rtol=1e-10, atol=1e-12)
fig, axes = plt.subplots(3, 1, figsize=(12, 12))
# 姿勢角
labels_att = ['Roll', 'Pitch', 'Yaw']
colors_att = ['#e74c3c', '#2ecc71', '#3498db']
for i, (label, color) in enumerate(zip(labels_att, colors_att)):
axes[0].plot(sol.t, np.degrees(sol.y[i]),
label=label, color=color, linewidth=1.5)
axes[0].axhline(y=0, color='gray', linestyle='--', alpha=0.5)
axes[0].set_ylabel('Attitude angle [deg]')
axes[0].set_title('Three-axis attitude control with 4 reaction wheels (pyramid config)')
axes[0].legend()
axes[0].grid(True, alpha=0.3)
# 衛星角速度
for i, (label, color) in enumerate(zip(labels_att, colors_att)):
axes[1].plot(sol.t, np.degrees(sol.y[3+i]),
label=label, color=color, linewidth=1.5)
axes[1].set_ylabel('Angular velocity [deg/s]')
axes[1].set_title('Satellite angular velocity')
axes[1].legend()
axes[1].grid(True, alpha=0.3)
# ホイール回転速度
for i, (label, color) in enumerate(zip(labels_rw, colors_rw)):
rpm = sol.y[6+i] * 60 / (2 * np.pi)
axes[2].plot(sol.t, rpm, label=label, color=color, linewidth=1.2)
axes[2].axhline(y=6000, color='gray', linestyle='--', alpha=0.5, label='Max speed')
axes[2].axhline(y=-6000, color='gray', linestyle='--', alpha=0.5)
axes[2].set_xlabel('Time [s]')
axes[2].set_ylabel('Wheel speed [rpm]')
axes[2].set_title('Reaction wheel speeds')
axes[2].legend(loc='upper right')
axes[2].grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('rw_attitude_control.png', dpi=150, bbox_inches='tight')
plt.show()
この完全なシミュレーションから、RWによる三軸姿勢制御の全体像が把握できます。
- 姿勢の収束: 初期姿勢偏差（ロール15度、ピッチ10度、ヨー20度）からゼロ目標に向けて滑らかに収束しています。各軸の固有振動数の違いにより、ヨー（$\omega_n$ 最大）が最も速く、ピッチ（$\omega_n$ 最小）が最も遅く整定しています。
- ホイール速度の変化: マヌーバ中にホイールが加速して衛星にトルクを発生させ、目標に近づくと減速しています。マヌーバ完了後、外乱トルクによりホイール速度が緩やかに蓄積し続けるのが確認できます。
- 4基の均等な負荷分担: ピラミッド配置の疑似逆行列による分配により、4基のホイールがほぼ均等にトルクを分担しています。特定のホイールだけが極端に高速回転することはありません。
- 外乱トルクの蓄積: 姿勢が整定した後も、ホイール速度は徐々に増加しています。このままではいずれ飽和に達するため、定期的なアンローディングが必要であることがシミュレーションからも確認できます。
設計パラメータの典型値
実際の衛星設計で参考になるRWの典型的なパラメータを示します。
| パラメータ | CubeSat (1U-3U) | 小型衛星 (50-200kg) | 中型衛星 (200-1000kg) |
|---|---|---|---|
| ホイール質量 | 0.03 ～ 0.1 kg | 0.5 ～ 2 kg | 2 ～ 10 kg |
| 最大トルク | 0.1 ～ 1 mNm | 1 ～ 50 mNm | 50 ～ 500 mNm |
| 最大角運動量 | 1 ～ 10 mNms | 0.1 ～ 5 Nms | 5 ～ 100 Nms |
| 最大回転速度 | 5000 ～ 8000 rpm | 4000 ～ 6000 rpm | 3000 ～ 6000 rpm |
まとめ
本記事では、リアクションホイールによる姿勢制御の原理、設計、実装について解説しました。
- 角運動量交換の原理: RWはホイールの回転速度を変化させることで反作用トルクを衛星に伝え、推進剤を消費せずに姿勢を制御します
- 配置行列とトルク分配: $n$ 基のRWの配置は行列 $\bm{A}$ で記述され、冗長構成での最適トルク分配は疑似逆行列 $\bm{A}^+$ で求まります
- ピラミッド配置: 4輪ピラミッド配置は $\bm{A}\bm{A}^T = \frac{4}{3}\bm{I}_3$ を満たす等方的な配置であり、1基故障時も3軸制御を維持できます
- 飽和への対処: トルク飽和はマヌーバ制限で、角運動量飽和はアンローディングで対処します
- MW, CMGとの比較: RWは構造がシンプルで信頼性が高く、中小型衛星に最適ですが、大トルクが必要な場合はCMGが選択されます
次のステップとして、以下の記事も参考にしてください。