質点の運動方程式と解法

運動方程式 $\bm{F} = m\bm{a}$ は、力が与えられたときに物体の運動を予測するための基本方程式です。実際の問題では、力を具体的に同定し、座標系を設定し、微分方程式として解くスキルが必要です。

本記事の内容

  • 運動方程式の成分表示
  • 典型的な問題の解法(自由落下、斜面、抵抗力)
  • 微分方程式としての運動方程式
  • Pythonでの数値積分

前提知識

運動方程式の成分表示

3次元デカルト座標では:

$$ m\ddot{x} = F_x, \quad m\ddot{y} = F_y, \quad m\ddot{z} = F_z $$

極座標 $(r, \theta)$ では:

$$ m(\ddot{r} – r\dot{\theta}^2) = F_r, \quad m(r\ddot{\theta} + 2\dot{r}\dot{\theta}) = F_\theta $$

解法の手順

  1. 物体を特定し、着目する質点を決める
  2. 座標系を設定(問題の対称性に合わせる)
  3. 全ての力を同定して自由体図を描く
  4. 運動方程式を各成分について書き下す
  5. 初期条件を設定して微分方程式を解く

具体例

速度に比例する抵抗力

$$ m\dot{v} = mg – bv $$

変数分離で解くと:

$$ v(t) = \frac{mg}{b}\left(1 – e^{-bt/m}\right), \quad v_{\infty} = \frac{mg}{b} $$

2次元斜方投射

初速 $v_0$、仰角 $\alpha$ での投射:

$$ \ddot{x} = 0, \quad \ddot{y} = -g $$

$$ x(t) = v_0\cos\alpha \cdot t, \quad y(t) = v_0\sin\alpha \cdot t – \frac{g}{2}t^2 $$

飛距離 $R = v_0^2\sin(2\alpha)/g$ は $\alpha = 45°$ で最大。

Pythonでの実装

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

g = 9.8

# 斜方投射(空気抵抗なし)
fig, axes = plt.subplots(1, 2, figsize=(13, 5))

v0 = 30  # m/s
for alpha_deg in [15, 30, 45, 60, 75]:
    alpha = np.radians(alpha_deg)
    t_flight = 2*v0*np.sin(alpha)/g
    t = np.linspace(0, t_flight, 200)
    x = v0*np.cos(alpha)*t
    y = v0*np.sin(alpha)*t - 0.5*g*t**2
    axes[0].plot(x, y, label=f'$\\alpha={alpha_deg}°$')

axes[0].set_title('斜方投射(空気抵抗なし)')
axes[0].set_xlabel('x [m]')
axes[0].set_ylabel('y [m]')
axes[0].legend()
axes[0].grid(True, alpha=0.3)
axes[0].set_ylim(0, None)

# 速度に比例する抵抗力付き落下
m, b = 1.0, 0.5
def drag_ode(t, state):
    x, y, vx, vy = state
    v = np.sqrt(vx**2 + vy**2 + 1e-10)
    ax = -(b/m)*vx
    ay = -g - (b/m)*vy
    return [vx, vy, ax, ay]

v0 = 30
alpha = np.radians(45)
sol = solve_ivp(drag_ode, [0, 8], [0, 0, v0*np.cos(alpha), v0*np.sin(alpha)],
               t_eval=np.linspace(0, 8, 500), events=lambda t,y: y[1])

# 空気抵抗なし比較
t_nd = np.linspace(0, 2*v0*np.sin(alpha)/g, 200)
x_nd = v0*np.cos(alpha)*t_nd
y_nd = v0*np.sin(alpha)*t_nd - 0.5*g*t_nd**2

axes[1].plot(x_nd, y_nd, 'b--', label='抵抗なし')
mask = sol.y[1] >= 0
axes[1].plot(sol.y[0][mask], sol.y[1][mask], 'r-', linewidth=2, label='抵抗あり')
axes[1].set_title('空気抵抗の影響($\\alpha=45°$)')
axes[1].set_xlabel('x [m]')
axes[1].set_ylabel('y [m]')
axes[1].legend()
axes[1].grid(True, alpha=0.3)
axes[1].set_ylim(0, None)

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

まとめ

  • 運動方程式を成分に分解し、各方向独立に解く
  • 抵抗力がある場合は1階ODEとして変数分離で解ける
  • 斜方投射は $\alpha=45°$ で飛距離最大(空気抵抗なし)
  • 空気抵抗により飛距離は減少し、軌道は非対称になる

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