数百キログラムの地球観測衛星を、ある地点から別の地点へ数秒で振り向けたい——こういう「機敏な姿勢変更（アジャイル・スルー）」を実現する装置がコントロールモーメントジャイロ（CMG: Control Moment Gyroscope）です。国際宇宙ステーション（ISS）の姿勢維持、ハッブル宇宙望遠鏡の視線指向、そして高解像度の光学偵察衛星の高速な首振りは、いずれもCMGなしには成り立ちません。
なぜリアクションホイールではなくCMGなのでしょうか。リアクションホイールは「ホイールを加減速させた反作用」でトルクを作りますが、大きなトルクを出そうとするとモーターも大きく重くなります。一方CMGは、高速で回り続けるホイールの向きをほんの少し傾けるだけで、桁違いに大きなトルクを取り出せます。これをトルク増幅と呼びます。その代わり、CMG特有の厄介な問題——ある姿勢方向にまったくトルクを出せなくなる特異点（シンギュラリティ）——が現れます。
この記事では、CMGが生み出すジャイロトルクを角運動量の交換として定式化し、出力トルク $\bm{\tau} = \bm{\Omega} \times \bm{h}$ を導出します。リアクションホイールと対比してトルク増幅の原理を明らかにし、複数CMGを組んだピラミッド配置のヤコビアンから特異点がどう生じるか、そして擬似逆行列によるステアリング則で指令トルクをどう各ジンバルに配分するかまで、順を追って解説します。最後にPythonで、単一CMGのトルクと剛体の姿勢応答、ピラミッドクラスタの特異点近傍でトルク発生能力が落ちる様子を数値計算して可視化します。
本記事の内容
- CMGの動作原理をジャイロ効果の直感から理解する
- 出力トルク $\bm{\tau} = \bm{\Omega} \times \bm{h}$ を角運動量交換から導出する
- リアクションホイールと対比してトルク増幅を定量化する
- ピラミッド配置のヤコビアンと特異点、擬似逆行列ステアリング則を導く
- Pythonで単一CMGの姿勢応答とクラスタの特異点を可視化する
CMGは大きく3つの部品からなります。高速で回り続けるフライホイール（青の角運動量 $h$ は大きさ一定）、それを傾けるジンバル（オレンジの回転 $\dot{\delta}$）、そしてその結果として胴体に現れる出力トルク（赤の $\tau=\Omega\times h$）です。ホイールの回転そのものではなく、その「向きを傾ける」操作がトルクを生む点が、この図の要になります。
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
CMGとは
自転車に乗っているとき、なぜ倒れずに進めるのか不思議に思ったことはないでしょうか。速く回転している前輪は「向きを変えられることに抵抗する」性質を持っています。回っているコマがなかなか倒れないのも同じ理由です。この「回転している物体は自分の回転軸の向きを保とうとする」という性質がジャイロ効果であり、CMGはこれを姿勢制御に積極的に利用します。
もう少し具体的に想像してみましょう。手のひらに、勢いよく回っている自転車のホイールを軸で持ったとします。このホイールの回転軸を無理やり左へ傾けようとすると、ホイールは素直に左へは行かず、なぜか上や下へぐいっと反発してきます。この「傾けた方向とは直交する向きに現れる力」こそがジャイロトルクです。CMGは、この反発トルクを衛星の胴体に伝えて姿勢を変える装置です。
CMGの構造はシンプルです。まず、高速（毎分数千〜数万回転）で回り続けるホイール（フライホイール、ロータ）があります。このホイールは大きな角運動量 $\bm{h}$ を持ちます。ホイールはジンバルという枠に取り付けられており、ジンバルを回すとホイールの回転軸ごと向きが変わります。ホイールを回すモーター（スピンモーター）は一定速度を保つだけで、姿勢制御の「操作」はもっぱらジンバルを回すこと（ジンバル角 $\delta$ を変えること）で行います。
ここで、リアクションホイールとの本質的な違いを押さえておきましょう。リアクションホイールは「ホイールの回転速度（角運動量の大きさ）」を変えてトルクを作ります。対してCMGは「ホイールの回転速度は一定に保ったまま、その向きを変える」ことでトルクを作ります。角運動量というベクトルの、大きさをいじるのがリアクションホイール、向きをいじるのがCMG——この対比が両者を理解する鍵になります。
左のリアクションホイールは角運動量ベクトル $h$ を同じ向きのまま伸縮させ（加減速）、その大きさの変化 $\Delta h$ がトルクになります。右のCMGはベクトルの長さを保ったまま向きだけを回し、その回転がトルクを生みます。同じ「角運動量を変化させる」でも、操作する自由度が大きさか向きかで装置の性格が正反対になることが見て取れます。
では、向きを変えるとなぜトルクが出るのか。それを角運動量の言葉で厳密に定式化していきましょう。
角運動量の交換とジャイロトルク
CMGの出力トルクを導くための出発点は、ニュートン・オイラーの回転運動の法則です。ある物体に働くトルク $\bm{\tau}$ は、その物体の角運動量 $\bm{h}$ の時間変化率に等しい、という関係です。
$$ \begin{equation} \bm{\tau} = \frac{d\bm{h}}{dt} \end{equation} $$
この式は「トルクとは角運動量を変える働きである」と述べています。ここで重要なのは、角運動量 $\bm{h}$ はベクトルだということです。ベクトルが時間変化する原因は2つあります——大きさが変わるか、向きが変わるかです。リアクションホイールは前者を、CMGは後者を使います。
CMGのホイールが持つ角運動量を $\bm{h}$ とし、その大きさ $h_0 = |\bm{h}|$ は一定に保たれているとします。すると $\bm{h}$ の変化は向きの変化だけになります。単位ベクトル $\hat{\bm{h}}$ を使って $\bm{h} = h_0 \hat{\bm{h}}$ と書き、ジンバルの回転を考えます。
ジンバルの回転を、ジンバル軸まわりの角速度ベクトル $\bm{\Omega}$ で表します。剛体上に固定されたベクトルが角速度 $\bm{\Omega}$ で回転するとき、そのベクトルの時間微分は外積で書ける——これはベクトルの回転に関する基本公式（ポアソンの公式）です。ホイールの角運動量ベクトル $\bm{h}$ はジンバルとともに向きを変えるので、この公式を $\bm{h}$ に適用すると
$$ \begin{equation} \bm{\tau} = \frac{d\bm{h}}{dt} = \bm{\Omega} \times \bm{h} \end{equation} $$
が得られます。これがCMGの出力トルクの基本式です。念のため、なぜ回転するベクトルの微分が外積になるのかを確認しておきましょう。
出力トルク式の導出
ジンバル軸方向の単位ベクトルを $\hat{\bm{g}}$、ジンバル角を $\delta$、ジンバル角速度を $\dot{\delta}$ とします。ホイールの角運動量ベクトルは、ジンバル軸に垂直な平面内でジンバル角 $\delta$ とともに回転します。この平面を張る2つの直交する単位ベクトルを $\hat{\bm{u}}, \hat{\bm{v}}$（ともに $\hat{\bm{g}}$ に垂直）とすると、角運動量の単位ベクトルは
$$ \hat{\bm{h}}(\delta) = \cos\delta \, \hat{\bm{u}} + \sin\delta \, \hat{\bm{v}} $$
と表せます。両辺に大きさ $h_0$ を掛けて $\bm{h}(\delta) = h_0 \hat{\bm{h}}(\delta)$ とし、時間で微分します。$\delta$ が時間の関数であることに注意して連鎖律を使うと
$$ \begin{align} \frac{d\bm{h}}{dt} &= h_0 \frac{d\hat{\bm{h}}}{d\delta}\dot{\delta} \\ &= h_0 \dot{\delta}\left(-\sin\delta \, \hat{\bm{u}} + \cos\delta \, \hat{\bm{v}}\right) \end{align} $$
となります。ここで括弧の中身 $-\sin\delta\,\hat{\bm{u}} + \cos\delta\,\hat{\bm{v}}$ が何者かを見極めます。$\hat{\bm{g}}, \hat{\bm{u}}, \hat{\bm{v}}$ が右手系（$\hat{\bm{g}} \times \hat{\bm{u}} = \hat{\bm{v}}$、$\hat{\bm{g}} \times \hat{\bm{v}} = -\hat{\bm{u}}$）をなすとして、$\hat{\bm{g}} \times \hat{\bm{h}}$ を計算すると
$$ \begin{align} \hat{\bm{g}} \times \hat{\bm{h}} &= \hat{\bm{g}} \times (\cos\delta\,\hat{\bm{u}} + \sin\delta\,\hat{\bm{v}}) \\ &= \cos\delta\,(\hat{\bm{g}}\times\hat{\bm{u}}) + \sin\delta\,(\hat{\bm{g}}\times\hat{\bm{v}}) \\ &= \cos\delta\,\hat{\bm{v}} – \sin\delta\,\hat{\bm{u}} \end{align} $$
となり、これはまさに先ほどの括弧の中身と一致します。したがって
$$ \begin{equation} \bm{\tau} = \frac{d\bm{h}}{dt} = h_0 \dot{\delta}\,(\hat{\bm{g}} \times \hat{\bm{h}}) = (\dot{\delta}\,\hat{\bm{g}}) \times (h_0 \hat{\bm{h}}) = \bm{\Omega} \times \bm{h} \end{equation} $$
と、$\bm{\Omega} = \dot{\delta}\,\hat{\bm{g}}$（ジンバル角速度ベクトル）を使って外積の形にまとまります。導出のポイントは、角運動量の向きの変化率がジンバル軸との外積として現れるという一点です。
ジンバル角 $\delta$ が進むと、大きさ一定の角運動量ベクトル $h(\delta)$ の先端は円周上を移動します。その移動の向き、すなわち $dh/dt$ は円の接線方向、つまり $h$ に対して直角の方向を向きます。トルクは $\tau=dh/dt$ ですから、この図は「トルクが常に角運動量と直交して現れる」という外積の幾何学的な意味そのものを表しています。
この式から、出力トルクの向きと大きさが読み取れます。まず向きは $\bm{\Omega} \times \bm{h}$、すなわちジンバル軸とホイール角運動量の両方に垂直な方向です。単一ジンバルのCMGでは $\bm{\Omega}$（ジンバル軸方向）と $\bm{h}$（それに垂直な平面内）が常に直交するので、大きさは
$$ |\bm{\tau}| = |\bm{\Omega}||\bm{h}|\sin 90^\circ = \dot{\delta}\, h_0 $$
となります。トルクの大きさは「ジンバル角速度 $\dot{\delta}$」と「ホイール角運動量 $h_0$」の積です。ここに、次に述べるトルク増幅のカラクリが隠れています。
3本のベクトルの関係を3次元で見ておきます。ジンバル軸方向の $\Omega$（オレンジ）とホイールの角運動量 $h$（シアン）は直交し、出力トルク $\tau=\Omega\times h$（赤）はその両方に垂直な残る1方向を向きます。右手系の外積そのままの配置で、単一CMGでは $\Omega$ と $h$ が常に直角なので、大きさは $|\tau|=\dot{\delta}\,h_0$ に簡単化されます。
トルク増幅の原理
出力トルクの大きさ $|\bm{\tau}| = \dot{\delta}\,h_0$ をよく見てください。この式は、ホイールの角運動量 $h_0$ が大きければ大きいほど、ほんの小さなジンバル角速度 $\dot{\delta}$ でも巨大なトルクが出ることを意味します。これがトルク増幅です。
リアクションホイールと数値で比べてみましょう。リアクションホイールは、慣性モーメント $I_w$ のホイールを角加速度 $\dot{\omega}_w$ で回してトルクを作ります。
$$ \tau_{\text{RW}} = I_w \dot{\omega}_w $$
大きなトルクを出すには、モーターが大きな角加速度 $\dot{\omega}_w$ を生み出さなければなりません。しかもホイールの回転数には上限（飽和速度）があるので、大トルクを長く出し続けることはできません。トルクの大きさは、そのままモーターが供給すべきトルクの大きさです。
一方CMGでは、スピンモーターはホイールを一定速度で回し続けるだけです。姿勢制御のトルクを担うジンバルモーターが供給すべきトルクは、ジンバルと支持部の慣性モーメント $J_g$ に対して $J_g \ddot{\delta}$ 程度でしかありません。ところが衛星胴体が受け取る出力トルクは $\dot{\delta}\,h_0$ です。両者の比を取ると
$$ \frac{|\bm{\tau}_{\text{out}}|}{|\bm{\tau}_{\text{gimbal}}|} \sim \frac{\dot{\delta}\,h_0}{J_g \ddot{\delta}} $$
となり、ホイール角運動量 $h_0$ が大きいほど、ジンバルに与えた小さな仕事から大きな出力トルクが取り出せることがわかります。具体的な数字を入れてみましょう。ホイール角運動量 $h_0 = 50\,\text{N·m·s}$ のCMGを、ジンバル角速度 $\dot{\delta} = 1\,\text{rad/s}$ で回すと、出力トルクは
$$ |\bm{\tau}| = \dot{\delta}\, h_0 = 1 \times 50 = 50\,\text{N·m} $$
にもなります。同じ $50\,\text{N·m}$ をリアクションホイールで出そうとすると、慣性モーメント $I_w = 0.05\,\text{kg·m}^2$ のホイールを $\dot{\omega}_w = 1000\,\text{rad/s}^2$ もの猛烈な角加速度で回す必要があり、モーターにそのまま $50\,\text{N·m}$ の負荷がかかります。CMGならジンバルをゆっくり回すだけで済むのです。
対数軸の棒グラフで3つのトルクを比べています。衛星胴体が受け取る出力トルク（赤、$50\,\text{N·m}$）に対し、ジンバルモーターが実際に供給する負荷（緑）は桁違いに小さく、これがCMGの増幅効果です。一方リアクションホイール（オレンジ）は同じ出力を得るためにモーターがそのまま $50\,\text{N·m}$ を負担しなければならず、CMGとの負担の差が一目でわかります。
このトルク増幅こそが、大型で機敏な衛星がCMGを選ぶ理由です。同じ質量・電力のアクチュエータでも、CMGはリアクションホイールの数十倍のトルクを叩き出せます。ただし、この強力さには代償があります。角運動量の「向き」を操作するという仕組み自体が、次に述べる特異点問題を生むのです。まずは、複数のCMGを束ねたときの定式化を用意しましょう。
複数CMGクラスタとヤコビアン
単一のCMGは、出力トルクの向きが $\bm{\Omega}\times\bm{h}$ という1つの方向にほぼ縛られてしまいます。3軸すべての姿勢を自由に制御するには、複数のCMGを異なる向きに組み合わせる必要があります。実用衛星では、4基のCMGを四角錐（ピラミッド）の側面に配置するピラミッド配置が定番です。冗長性（1基故障しても3軸制御を維持）と、なるべく等方的な角運動量エンベロープを得られるのが理由です。
$n$ 基のCMGクラスタを考えます。$i$ 番目のCMGはジンバル軸 $\hat{\bm{g}}_i$ を持ち、ジンバル角 $\delta_i$ に応じた角運動量 $\bm{h}_i(\delta_i)$ を出します（大きさは $h_0$ で共通とします）。クラスタ全体の合成角運動量は各CMGの和
$$ \bm{H}(\bm{\delta}) = \sum_{i=1}^{n} \bm{h}_i(\delta_i) $$
です。ここで $\bm{\delta} = (\delta_1, \dots, \delta_n)^\top$ はジンバル角ベクトルです。クラスタが衛星胴体に及ぼす出力トルクは、角運動量保存則から合成角運動量の変化率の符号を反転したもの
$$ \bm{\tau} = -\frac{d\bm{H}}{dt} $$
になります（クラスタが角運動量 $+d\bm{H}$ を得れば、衛星胴体は反作用として $-d\bm{H}$ を受け取る、という反作用の関係です）。合成角運動量を時間で微分し、連鎖律を各ジンバル角に適用すると
$$ \begin{align} \frac{d\bm{H}}{dt} &= \sum_{i=1}^{n} \frac{\partial \bm{h}_i}{\partial \delta_i}\dot{\delta}_i \\ &= \bm{A}(\bm{\delta})\,\dot{\bm{\delta}} \end{align} $$
とまとまります。ここで $\bm{A}(\bm{\delta})$ は $3 \times n$ のヤコビアン行列で、その第 $i$ 列は
$$ \bm{a}_i = \frac{\partial \bm{h}_i}{\partial \delta_i} = \hat{\bm{g}}_i \times \bm{h}_i(\delta_i) $$
です。各列 $\bm{a}_i$ は「$i$ 番目のジンバルを単位速度で回したときにクラスタが出せるトルク方向」を表しています。したがって出力トルクとジンバル角速度の関係は
$$ \begin{equation} \bm{\tau} = -\bm{A}(\bm{\delta})\,\dot{\bm{\delta}} \end{equation} $$
という線形写像になります。ヤコビアン $\bm{A}$ は $3 \times n$ 行列で、入力（$n$ 次元のジンバル角速度）を出力（3次元のトルク）へ写します。姿勢制御では逆に「欲しいトルク $\bm{\tau}$ を出すには、どのジンバルをどれだけ回せばよいか（$\dot{\bm{\delta}}$）」を解きたい——これがステアリングの問題です。
ピラミッド配置のヤコビアンを具体的に書く
$n=4$ のピラミッド配置を具体的に書き下します。スキュー角（ピラミッドの側面が底面となす角）を $\beta$ とし、4基のジンバル軸を四角錐の側面法線に沿って配置します。標準的な取り方では、各CMGの角運動量ベクトルは次のようになります。
$$ \bm{h}_1 = h_0\begin{bmatrix} -\cos\beta\sin\delta_1 \\ \cos\delta_1 \\ \sin\beta\sin\delta_1 \end{bmatrix},\quad \bm{h}_2 = h_0\begin{bmatrix} -\cos\delta_2 \\ -\cos\beta\sin\delta_2 \\ \sin\beta\sin\delta_2 \end{bmatrix} $$
$$ \bm{h}_3 = h_0\begin{bmatrix} \cos\beta\sin\delta_3 \\ -\cos\delta_3 \\ \sin\beta\sin\delta_3 \end{bmatrix},\quad \bm{h}_4 = h_0\begin{bmatrix} \cos\delta_4 \\ \cos\beta\sin\delta_4 \\ \sin\beta\sin\delta_4 \end{bmatrix} $$
各列 $\bm{a}_i = \partial\bm{h}_i/\partial\delta_i$ は、上を $\delta_i$ で微分すればすぐ得られます。たとえば第1列は
$$ \bm{a}_1 = h_0\begin{bmatrix} -\cos\beta\cos\delta_1 \\ -\sin\delta_1 \\ \sin\beta\cos\delta_1 \end{bmatrix} $$
です。スキュー角には $\beta = 53.13^\circ$（$\cos\beta = 0.6$, $\sin\beta = 0.8$）がよく使われます。この角度だと角運動量エンベロープがほぼ球に近くなり、どの方向にも均等に角運動量を蓄えられるためです。
4基のジンバル軸 $g_1,\dots,g_4$ を四角錐（ピラミッド）の側面に沿って外向きに配置した様子です。各軸が互いに異なる方向を向くことで、4基の出力トルク方向が3次元空間を張り、どの姿勢軸にもトルクを回せるようになります。スキュー角 $\beta=53.13^\circ$ はこの配置が最も等方的（球に近い角運動量エンベロープ）になる定番値です。
次節では、このヤコビアンがある姿勢で「潰れる」——ランクが落ちてトルクを出せなくなる——特異点を調べます。
特異点問題
CMGクラスタの最大の弱点が特異点（シンギュラリティ）です。特異点とは、すべてのジンバルをどう回しても、ある特定の方向にはトルクをまったく出せなくなる状態を指します。これは数学的には、ヤコビアン $\bm{A}$ のランクが3未満に落ちることに対応します。
なぜランク落ちが起きるのか、直感的に理解しましょう。ヤコビアンの各列 $\bm{a}_i = \hat{\bm{g}}_i \times \bm{h}_i$ は、$i$ 番目のCMGが出せるトルク方向でした。もしすべてのCMGのトルク方向 $\bm{a}_i$ がたまたま同一平面上に並んでしまうと、その平面に垂直な方向にはどのCMGもトルクを出せません。この状態が特異点です。角運動量の向きだけをいじるという仕組みの必然として、ジンバル角の組み合わせによってはこういう「打つ手なし」の配置が生じてしまうのです。
左は特異点から遠い配置で、4本のトルク方向ベクトル $a_i$ が3次元空間に広がり可操作度 $m=1.27$ と余裕があります。右は特異配置で、ベクトルがほぼ同一平面に潰れて $m=0$ となり、その平面に垂直な方向にはどのCMGもトルクを出せません。ランク落ち＝特異点とは、まさにこの「トルク方向が張れなくなる」幾何を指しています。
特異点への近さを測る指標が可操作度（マニピュラビリティ） $m$ です。$3 \times n$ のヤコビアンに対して
$$ \begin{equation} m(\bm{\delta}) = \sqrt{\det\!\left(\bm{A}\bm{A}^\top\right)} \end{equation} $$
と定義します。$\bm{A}\bm{A}^\top$ は $3\times3$ の対称行列で、その行列式が正のときクラスタは3方向すべてにトルクを出せます。$m$ が大きいほど「特異点から遠い（余裕がある）」、$m \to 0$ で特異点に到達します。行列式がゼロになるのは $\bm{A}\bm{A}^\top$ が非正則、すなわち $\bm{A}$ のランクが3未満になるときですから、$m=0 \Leftrightarrow \text{特異点}$ が成り立ちます。
$\bm{A}\bm{A}^\top$ の固有値を $\sigma_1^2 \ge \sigma_2^2 \ge \sigma_3^2 \ge 0$（$\sigma_i$ は $\bm{A}$ の特異値）とすれば $m = \sigma_1\sigma_2\sigma_3$ です。ひとつでも特異値がゼロに近づくと $m$ は急落し、その特異値に対応する方向のトルク発生能力が失われます。可操作度は、CMGクラスタが「今どれだけ元気にトルクを出せるか」を1つの数字で表す健康診断のようなものだと思ってください。
特異点には性質の異なるものがあります。飽和特異点は角運動量が最大まで蓄積されてもうこれ以上ためられない外側の限界（エンベロープ境界）で起き、これは物理的に避けられません。一方内部特異点はエンベロープの内側で起きるもので、ジンバル配置の巡り合わせによって生じます。厄介なのは内部特異点で、通り抜けようとするとジンバル角速度が発散してしまいます。この特異点をどう避けながら指令トルクを配分するか——それがステアリング則の役目です。
4基のジンバル角をランダムに振って合成角運動量 $H$ を打った点群です。到達可能な $H$ の範囲がほぼ球状の塊をなしており、これがピラミッド配置の角運動量エンベロープです。飽和特異点はこの塊の外殻（境界）で起こり物理的に避けられない一方、内部特異点は塊の内側で巡り合わせによって生じる——という2種類の特異点の住み分けが、この分布から直感できます。
擬似逆行列によるステアリング則
姿勢制御則（PD制御やLQRなど）が「この瞬間これだけのトルク $\bm{\tau}_{\text{cmd}}$ が欲しい」と要求してきたとき、CMGクラスタはそれを各ジンバルの回転速度 $\dot{\bm{\delta}}$ へ翻訳しなければなりません。解くべき方程式は、符号を吸収して
$$ \bm{A}(\bm{\delta})\,\dot{\bm{\delta}} = \bm{\tau}_{\text{cmd}} $$
です。ここで $\bm{A}$ は $3 \times n$ 行列。冗長構成（$n>3$、ピラミッドなら $n=4$）では未知数 $\dot{\bm{\delta}}$ が4個に対して方程式が3本しかなく、解が無数に存在します。この不定性を「もっとも小さいジンバル速度で済ませたい」という基準で1つに絞るのが最小ノルム解で、ムーア・ペンローズの擬似逆行列を使って書けます。
最小ノルム解を導きましょう。制約 $\bm{A}\dot{\bm{\delta}} = \bm{\tau}_{\text{cmd}}$ のもとで $\|\dot{\bm{\delta}}\|^2$ を最小化する問題を、ラグランジュの未定乗数 $\bm{\lambda}$ で解きます。ラグランジアンは
$$ L = \tfrac{1}{2}\dot{\bm{\delta}}^\top\dot{\bm{\delta}} + \bm{\lambda}^\top\left(\bm{\tau}_{\text{cmd}} – \bm{A}\dot{\bm{\delta}}\right) $$
です。$\dot{\bm{\delta}}$ で偏微分してゼロと置くと $\dot{\bm{\delta}} = \bm{A}^\top\bm{\lambda}$ を得ます。これを制約に代入すると $\bm{A}\bm{A}^\top\bm{\lambda} = \bm{\tau}_{\text{cmd}}$ となり、$\bm{A}\bm{A}^\top$ が正則なら $\bm{\lambda} = (\bm{A}\bm{A}^\top)^{-1}\bm{\tau}_{\text{cmd}}$ です。これを $\dot{\bm{\delta}} = \bm{A}^\top\bm{\lambda}$ に戻すと、ステアリング則
$$ \begin{equation} \dot{\bm{\delta}} = \bm{A}^\top\left(\bm{A}\bm{A}^\top\right)^{-1}\bm{\tau}_{\text{cmd}} = \bm{A}^{+}\bm{\tau}_{\text{cmd}} \end{equation} $$
が導かれます。$\bm{A}^{+} = \bm{A}^\top(\bm{A}\bm{A}^\top)^{-1}$ が擬似逆行列です。導出のポイントは、最小ノルムという素直な要求からラグランジュ乗数法だけで擬似逆行列が自然に出てくる、という点です。
しかしこのステアリング則には落とし穴があります。特異点近傍では $\bm{A}\bm{A}^\top$ の行列式（＝可操作度 $m^2$）がゼロに近づくため、逆行列 $(\bm{A}\bm{A}^\top)^{-1}$ の成分が発散します。すると要求トルクが小さくても、計算されるジンバル速度 $\dot{\bm{\delta}}$ が跳ね上がり、物理的に実現不可能な指令になってしまいます。
この発散を抑える実用的な手法が特異点ロバスト逆行列（SR-inverse）です。$\bm{A}\bm{A}^\top$ に微小な正の項 $\lambda\bm{I}$ を足して
$$ \begin{equation} \dot{\bm{\delta}} = \bm{A}^\top\left(\bm{A}\bm{A}^\top + \lambda\bm{I}\right)^{-1}\bm{\tau}_{\text{cmd}} \end{equation} $$
とします。$\lambda > 0$ を入れることで、たとえ $\bm{A}\bm{A}^\top$ が特異になっても $\bm{A}\bm{A}^\top + \lambda\bm{I}$ は正則を保ち、ジンバル速度の発散を防げます。代償として特異点方向のトルク精度は落ちますが（指令通りのトルクが出ない）、少なくともアクチュエータが破綻せずに特異点を「通り抜ける」ことができます。$\lambda$ を可操作度 $m$ に応じて可変にする（特異点から遠いときは $\lambda \approx 0$、近づくと $\lambda$ を増やす）のが定石です。
理屈が出そろったので、これらをPythonで確かめましょう。単一CMGの姿勢応答、ピラミッドクラスタの可操作度マップ、そして特異点近傍でジンバル速度が発散する様子を数値計算します。
Pythonでの実装
単一CMGによる姿勢スルー
まず、1軸まわりに回る衛星を1基のCMGで駆動する最小モデルを作ります。CMGの出力トルクは $\tau(t) = h_0\,\dot{\delta}(t)$。台形のジンバル角速度プロファイル（加速→一定→減速）を与えて、剛体の角加速度・角速度・姿勢角を積分し、「小さなジンバル運動で衛星が大きく首を振る」ことを確認します。
import numpy as np
import matplotlib, matplotlib.pyplot as plt
# 日本語フォント設定
for cand in ["Hiragino Sans", "Yu Gothic", "Noto Sans CJK JP", "IPAexGothic", "Meiryo"]:
if any(cand == f.name for f in matplotlib.font_manager.fontManager.ttflist):
plt.rcParams["font.family"] = cand
break
plt.rcParams["axes.unicode_minus"] = False
# --- パラメータ ---
h0 = 50.0 # ホイール角運動量 [N·m·s]
I_sc = 500.0 # 衛星の慣性モーメント [kg·m^2]
dt = 0.001 # 時間刻み [s]
T = 6.0 # シミュレーション時間 [s]
t = np.arange(0, T, dt)
# --- 台形ジンバル角速度プロファイル ---
ddelta = np.zeros_like(t) # ジンバル角速度 [rad/s]
peak = 0.30 # 最大ジンバル角速度 [rad/s]
ddelta[(t >= 0.5) & (t < 1.5)] = np.linspace(0, peak, np.sum((t >= 0.5) & (t < 1.5)))
ddelta[(t >= 1.5) & (t < 4.0)] = peak
ddelta[(t >= 4.0) & (t < 5.0)] = np.linspace(peak, 0, np.sum((t >= 4.0) & (t < 5.0)))
delta = np.cumsum(ddelta) * dt # ジンバル角 [rad]（積分）
# --- CMG出力トルクと剛体応答（オイラー積分）---
tau = h0 * ddelta # 出力トルク τ = h0 * δ_dot
omega = np.cumsum(tau / I_sc) * dt # 衛星角速度 [rad/s]
theta = np.cumsum(omega) * dt # 衛星姿勢角 [rad]
fig, ax = plt.subplots(2, 2, figsize=(11, 7))
ax[0, 0].plot(t, np.degrees(delta), color="tab:green")
ax[0, 0].set_title("ジンバル角 δ（入力：小さい）"); ax[0, 0].set_ylabel("δ [deg]")
ax[0, 1].plot(t, tau, color="tab:red")
ax[0, 1].set_title("CMG出力トルク τ = h₀·δ̇"); ax[0, 1].set_ylabel("τ [N·m]")
ax[1, 0].plot(t, np.degrees(omega), color="tab:blue")
ax[1, 0].set_title("衛星の角速度 ω"); ax[1, 0].set_ylabel("ω [deg/s]"); ax[1, 0].set_xlabel("時刻 [s]")
ax[1, 1].plot(t, np.degrees(theta), color="tab:purple")
ax[1, 1].set_title("衛星の姿勢角 θ（出力：大きい）"); ax[1, 1].set_ylabel("θ [deg]"); ax[1, 1].set_xlabel("時刻 [s]")
for a in ax.ravel(): a.grid(True, alpha=0.3)
plt.tight_layout(); plt.savefig("cmg_single.png", dpi=140); plt.show()
print(f"最大ジンバル角: {np.degrees(delta).max():.1f} deg")
print(f"最大出力トルク: {tau.max():.1f} N·m")
print(f"最終姿勢角: {np.degrees(theta)[-1]:.1f} deg")
このシミュレーションでは、ジンバル角の変化はたかだか数十度に収まっているのに、衛星の姿勢角は大きく変化します。トルクのグラフ（右上）が台形になっているのは $\tau = h_0\dot{\delta}$ のためで、ジンバル角速度の形をそのまま反映しています。加速局面でプラスのトルク、減速局面でマイナスのトルクが出て、最終的に衛星は静止（角速度ゼロ）した状態で新しい姿勢角に到達します。「小さなジンバル運動が大きな姿勢変化を生む」というトルク増幅の効果が、姿勢角のグラフ（右下）にはっきり現れています。
ピラミッドクラスタのヤコビアンと可操作度
次に、4基ピラミッド配置のヤコビアン $\bm{A}(\bm{\delta})$ を実装し、可操作度 $m = \sqrt{\det(\bm{A}\bm{A}^\top)}$ を計算します。まずはヤコビアンを組み立てる関数を用意します。
import numpy as np
beta = np.deg2rad(53.13) # スキュー角
cb, sb = np.cos(beta), np.sin(beta)
h0 = 1.0 # 正規化（h0=1 とする）
def jacobian(delta):
"""ジンバル角ベクトル delta(4,) からヤコビアン A(3x4) を返す"""
d1, d2, d3, d4 = delta
a1 = h0 * np.array([-cb*np.cos(d1), -np.sin(d1), sb*np.cos(d1)])
a2 = h0 * np.array([ np.sin(d2), -cb*np.cos(d2), sb*np.cos(d2)])
a3 = h0 * np.array([ cb*np.cos(d3), np.sin(d3), sb*np.cos(d3)])
a4 = h0 * np.array([-np.sin(d4), cb*np.cos(d4), sb*np.cos(d4)])
return np.column_stack([a1, a2, a3, a4]) # 3x4
def manipulability(delta):
A = jacobian(delta)
return np.sqrt(max(np.linalg.det(A @ A.T), 0.0))
# 全ジンバル角を同時に動かしたときの可操作度をスキャン
gs = np.linspace(-np.pi, np.pi, 400)
m_scan = np.array([manipulability(np.array([g, g, g, g])) for g in gs])
import matplotlib.pyplot as plt
plt.figure(figsize=(9, 4.5))
plt.plot(np.degrees(gs), m_scan, color="tab:orange")
plt.axhline(0, color="k", lw=0.8)
plt.xlabel("共通ジンバル角 δ [deg]"); plt.ylabel("可操作度 m = √det(AAᵀ)")
plt.title("ジンバル角に沿った可操作度（0に触れる点＝特異点）")
plt.grid(True, alpha=0.3); plt.tight_layout()
plt.savefig("cmg_manip.png", dpi=140); plt.show()
print(f"可操作度の最小値: {m_scan.min():.4f}")
print(f"可操作度の最大値: {m_scan.max():.4f}")
このグラフは、4基のジンバルを共通角 $\delta$ で一斉に動かしたときの可操作度の変化を示します。曲線が谷になってゼロ付近まで落ち込む角度が存在し、そこがまさに特異点です。可操作度が高い（山の）領域ではクラスタは3軸すべてに元気にトルクを出せますが、谷に落ちると特定方向のトルク発生能力を失います。CMGを運用するとは、この谷を避けながらジンバル角を動かし続けることに他なりません。
特異点近傍でジンバル速度が発散する様子
最後に、擬似逆行列ステアリング則 $\dot{\bm{\delta}} = \bm{A}^{+}\bm{\tau}_{\text{cmd}}$ が特異点近傍でどうなるかを見ます。可操作度 $m$ を横軸に、一定の指令トルクに対して要求されるジンバル速度のノルム $\|\dot{\bm{\delta}}\|$ を、通常の擬似逆行列とSR-inverse（$\lambda>0$）で比較します。
import numpy as np
tau_cmd = np.array([1.0, 0.0, 0.0]) # 指令トルク（特異配置で到達不能になる方向）
def steer_pinv(delta, tau, lam=0.0):
A = jacobian(delta)
M = A @ A.T + lam * np.eye(3)
return A.T @ np.linalg.solve(M, tau) # δ_dot (4,)
# 特異点へ向かってジンバル角を連続的に動かし、m と ||δ_dot|| を記録
alphas = np.linspace(0.0, 1.0, 300)
# 特異配置の一例（内部特異点に近づく経路）
d_sing = np.array([np.pi/2, 0.0, -np.pi/2, 0.0])
d_start = np.array([0.3, 0.3, 0.3, 0.3])
m_list, n_pinv, n_sr = [], [], []
for a in alphas:
d = (1 - a) * d_start + a * d_sing
m_list.append(manipulability(d))
n_pinv.append(np.linalg.norm(steer_pinv(d, tau_cmd, lam=0.0)))
n_sr.append(np.linalg.norm(steer_pinv(d, tau_cmd, lam=1e-2)))
import matplotlib.pyplot as plt
fig, ax1 = plt.subplots(figsize=(9, 4.8))
ax1.plot(alphas, n_pinv, color="tab:red", label="擬似逆行列 (λ=0)")
ax1.plot(alphas, n_sr, color="tab:blue", label="SR-inverse (λ=0.01)")
ax1.set_yscale("log"); ax1.set_xlabel("特異配置への接近パラメータ α")
ax1.set_ylabel("要求ジンバル速度 ‖δ̇‖ (対数)"); ax1.legend(loc="upper left")
ax2 = ax1.twinx()
ax2.plot(alphas, m_list, color="tab:green", ls="--", label="可操作度 m")
ax2.set_ylabel("可操作度 m", color="tab:green")
plt.title("特異点近傍：擬似逆行列は発散、SR-inverse は有界")
plt.tight_layout(); plt.savefig("cmg_singularity.png", dpi=140); plt.show()
print(f"接近端での可操作度: {m_list[-1]:.4f}")
print(f"接近端の ||δ̇|| (擬似逆): {n_pinv[-1]:.1f}")
print(f"接近端の ||δ̇|| (SR-inv): {n_sr[-1]:.1f}")
このグラフには2つの重要な事実が現れます。第一に、緑の破線（可操作度 $m$）が特異配置へ近づくにつれてゼロへ落ちること。第二に、そのとき赤い曲線（通常の擬似逆行列）の要求ジンバル速度が対数軸で天井知らずに発散すること——つまり、わずかなトルクを出すのに無限大に近いジンバル速度が要求されてしまいます。一方、青いSR-inverse は微小項 $\lambda$ のおかげで要求ジンバル速度が有界にとどまります。特異点方向のトルク精度は犠牲になりますが、アクチュエータを破綻させずに済むという実用上の折り合いが、この一枚から読み取れます。
まとめ
本記事では、コントロールモーメントジャイロ（CMG）の理論と実装を解説しました。
- 動作原理: CMGは高速回転ホイールの角運動量 $\bm{h}$ の「向き」をジンバルで変え、その変化率としてトルクを取り出す。出力トルクは $\bm{\tau} = \bm{\Omega}\times\bm{h}$、大きさは $\dot{\delta}\,h_0$
- トルク増幅: ホイール角運動量 $h_0$ が大きいほど、小さなジンバル角速度で巨大なトルクが出る。角運動量の大きさを変えるリアクションホイールに対し、CMGは向きを変えるだけで数十倍のトルクを実現する
- ヤコビアンと特異点: $n$ 基クラスタの出力トルクは $\bm{\tau} = -\bm{A}(\bm{\delta})\dot{\bm{\delta}}$。ヤコビアン $\bm{A}$ のランクが3未満に落ちる配置が特異点で、可操作度 $m = \sqrt{\det(\bm{A}\bm{A}^\top)}$ が特異点への近さを測る
- ステアリング則: 冗長構成では最小ノルム解 $\dot{\bm{\delta}} = \bm{A}^{+}\bm{\tau}_{\text{cmd}}$ を使うが、特異点近傍でジンバル速度が発散する。SR-inverse（$\lambda\bm{I}$ 付加）で発散を抑える
CMGは、リアクションホイールの延長線上にありながら「角運動量ベクトルの向きを操作する」という発想の転換によって圧倒的なトルクを生み、その代償として特異点という幾何学的な難題を抱えます。特異点回避ステアリング則の設計や、CMGが飽和したときに角運動量を逃がすアンローディング（磁気トルカやスラスタとの併用）は、実運用でさらに重要になるテーマです。
次のステップとして、以下の記事も参考にしてください。