剛体の回転運動はニュートン力学の中でも特に豊かな現象を生みます。オイラーの運動方程式は、剛体の回転を記述する基本方程式であり、コマの運動から衛星の姿勢制御まで幅広く応用されます。
本記事の内容
- 剛体の角運動量と慣性テンソル
- 主慣性モーメント
- オイラーの運動方程式の導出
- トルクフリー運動（自由回転）
- 地球ゴマの運動
- Pythonでの数値シミュレーション
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
剛体の角運動量と慣性テンソル
剛体の角運動量は角速度ベクトル $\bm{\omega}$ と慣性テンソル $\bm{I}$ を用いて表されます。
$$ \bm{L} = \bm{I}\bm{\omega} $$
慣性テンソルは $3 \times 3$ の対称行列です。
$$ \bm{I} = \begin{pmatrix} I_{xx} & -I_{xy} & -I_{xz} \\ -I_{yx} & I_{yy} & -I_{yz} \\ -I_{zx} & -I_{zy} & I_{zz} \end{pmatrix} $$
各成分は:
$$ I_{xx} = \int (y^2 + z^2) \, dm, \quad I_{xy} = \int xy \, dm $$
主慣性モーメント
慣性テンソルは実対称行列なので、適切な座標系（主軸）を選べば対角化できます。
$$ \bm{I} = \begin{pmatrix} I_1 & 0 & 0 \\ 0 & I_2 & 0 \\ 0 & 0 & I_3 \end{pmatrix} $$
$I_1, I_2, I_3$ を主慣性モーメントと呼びます。主軸座標系では角運動量は:
$$ \bm{L} = (I_1\omega_1, \, I_2\omega_2, \, I_3\omega_3) $$
オイラーの運動方程式の導出
空間に固定された座標系（慣性系）における角運動量の運動方程式は:
$$ \left(\frac{d\bm{L}}{dt}\right)_{\text{space}} = \bm{\tau} $$
剛体に固定された回転座標系（物体系）との関係は:
$$ \left(\frac{d\bm{L}}{dt}\right)_{\text{space}} = \left(\frac{d\bm{L}}{dt}\right)_{\text{body}} + \bm{\omega} \times \bm{L} $$
ここで右辺第2項はコリオリ項に相当します。物体系では慣性テンソルが時間変化しないため:
$$ \left(\frac{d\bm{L}}{dt}\right)_{\text{body}} = \bm{I}\dot{\bm{\omega}} $$
よって:
$$ \bm{I}\dot{\bm{\omega}} + \bm{\omega} \times \bm{I}\bm{\omega} = \bm{\tau} $$
主軸座標系で各成分を書き下すと:
$$ \begin{align} I_1 \dot{\omega}_1 + (I_3 – I_2)\omega_2\omega_3 &= \tau_1 \\ I_2 \dot{\omega}_2 + (I_1 – I_3)\omega_3\omega_1 &= \tau_2 \\ I_3 \dot{\omega}_3 + (I_2 – I_1)\omega_1\omega_2 &= \tau_3 \end{align} $$
$$ \boxed{I_i \dot{\omega}_i + \epsilon_{ijk}(I_k – I_j)\omega_j\omega_k = \tau_i} $$
これがオイラーの運動方程式です。非線形連立微分方程式であり、$\omega_i$ 同士の積の項が存在するため、一般には解析的に解くことができません。
トルクフリー運動（自由回転）
外部トルクがない場合（$\bm{\tau} = \bm{0}$）:
$$ \begin{align} I_1 \dot{\omega}_1 &= (I_2 – I_3)\omega_2\omega_3 \\ I_2 \dot{\omega}_2 &= (I_3 – I_1)\omega_3\omega_1 \\ I_3 \dot{\omega}_3 &= (I_1 – I_2)\omega_1\omega_2 \end{align} $$
軸対称体の場合（$I_1 = I_2 \neq I_3$）
$I_1 = I_2 = I_\perp$ とすると、第3式より:
$$ I_3 \dot{\omega}_3 = (I_\perp – I_\perp)\omega_1\omega_2 = 0 \implies \omega_3 = \text{const} $$
残りの2式は:
$$ \dot{\omega}_1 = -\Omega\omega_2, \quad \dot{\omega}_2 = \Omega\omega_1 $$
ここで $\Omega = \dfrac{(I_3 – I_\perp)}{I_\perp}\omega_3$ です。
解は:
$$ \omega_1(t) = \omega_\perp \cos(\Omega t + \phi_0), \quad \omega_2(t) = \omega_\perp \sin(\Omega t + \phi_0) $$
角速度ベクトル $\bm{\omega}$ は対称軸の周りを角速度 $\Omega$ で回転します。これが自由歳差運動（フリー・プレセッション）です。
安定性
主軸のうち最大または最小の慣性モーメントの軸まわりの回転は安定ですが、中間の軸まわりの回転は不安定です。これをテニスラケットの定理（ジャニベコフ効果）と呼びます。
地球ゴマ
地球ゴマは3つのジンバルリングの中で高速回転する円盤で、ジャイロ効果を実験的に示す道具です。
高速回転する剛体の角運動量 $\bm{L}$ は大きく、外部トルクが加わっても:
$$ \Delta\bm{L} = \bm{\tau}\Delta t $$
$\bm{L}$ の向きの変化は小さくなります。これがジャイロ効果であり、回転体が姿勢を安定に保つ原理です。
Pythonでの実装
import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
# --- オイラーの運動方程式（トルクフリー） ---
def euler_equations(t, omega, I1, I2, I3):
w1, w2, w3 = omega
dw1 = (I2 - I3) / I1 * w2 * w3
dw2 = (I3 - I1) / I2 * w3 * w1
dw3 = (I1 - I2) / I3 * w1 * w2
return [dw1, dw2, dw3]
fig = plt.figure(figsize=(16, 5))
# --- (1) 軸対称体のトルクフリー運動 ---
ax1 = fig.add_subplot(131, projection='3d')
I1, I2, I3 = 2.0, 2.0, 1.0 # 軸対称
omega0 = [1.0, 0.0, 5.0]
t_span = (0, 20)
t_eval = np.linspace(*t_span, 2000)
sol = solve_ivp(euler_equations, t_span, omega0, t_eval=t_eval,
args=(I1, I2, I3), rtol=1e-10, atol=1e-12)
ax1.plot(sol.y[0], sol.y[1], sol.y[2], 'b-', linewidth=1)
ax1.set_xlabel('$\\omega_1$')
ax1.set_ylabel('$\\omega_2$')
ax1.set_zlabel('$\\omega_3$')
ax1.set_title('軸対称体 ($I_1=I_2=2, I_3=1$)')
# --- (2) 非対称体のトルクフリー運動（ジャニベコフ効果） ---
ax2 = fig.add_subplot(132, projection='3d')
I1, I2, I3 = 1.0, 2.0, 3.0 # 非対称
omega0_unstable = [0.01, 5.0, 0.01] # 中間軸まわり（不安定）
t_span2 = (0, 50)
t_eval2 = np.linspace(*t_span2, 5000)
sol2 = solve_ivp(euler_equations, t_span2, omega0_unstable, t_eval=t_eval2,
args=(I1, I2, I3), rtol=1e-10, atol=1e-12)
ax2.plot(sol2.y[0], sol2.y[1], sol2.y[2], 'r-', linewidth=0.5)
ax2.set_xlabel('$\\omega_1$')
ax2.set_ylabel('$\\omega_2$')
ax2.set_zlabel('$\\omega_3$')
ax2.set_title('非対称体・中間軸回転（不安定）')
# --- (3) 保存量の確認 ---
ax3 = fig.add_subplot(133)
I1_c, I2_c, I3_c = 1.0, 2.0, 3.0
# エネルギーと角運動量の大きさ
T = 0.5*(I1_c*sol2.y[0]**2 + I2_c*sol2.y[1]**2 + I3_c*sol2.y[2]**2)
L2 = (I1_c*sol2.y[0])**2 + (I2_c*sol2.y[1])**2 + (I3_c*sol2.y[2])**2
# 規格化して表示
ax3.plot(sol2.t, T/T[0], 'b-', linewidth=2, label='$T/T_0$（運動エネルギー）')
ax3.plot(sol2.t, L2/L2[0], 'r--', linewidth=2, label='$L^2/L_0^2$（角運動量）')
ax3.set_xlabel('時間')
ax3.set_ylabel('規格化された値')
ax3.set_title('保存量の確認')
ax3.legend()
ax3.grid(True, alpha=0.3)
ax3.set_ylim(0.999, 1.001)
plt.tight_layout()
plt.savefig('euler_rigid_body.png', dpi=150, bbox_inches='tight')
plt.show()
print(f"運動エネルギー: 初期={T[0]:.6f}, 最終={T[-1]:.6f}, 相対変化={abs(T[-1]-T[0])/T[0]:.2e}")
print(f"角運動量の2乗: 初期={L2[0]:.6f}, 最終={L2[-1]:.6f}, 相対変化={abs(L2[-1]-L2[0])/L2[0]:.2e}")
まとめ
本記事では、オイラーの運動方程式と剛体の回転運動について解説しました。
- 剛体の角運動量: $\bm{L} = \bm{I}\bm{\omega}$（慣性テンソル）
- 主軸座標で慣性テンソルが対角化: $I_1, I_2, I_3$
- オイラーの運動方程式: $I_i\dot{\omega}_i + (I_k – I_j)\omega_j\omega_k = \tau_i$
- 軸対称体のトルクフリー運動は自由歳差運動
- 中間軸まわりの回転は不安定（ジャニベコフ効果）
- 運動エネルギーと角運動量の大きさは保存量
次のステップとして、以下の記事も参考にしてください。