ニュートンの運動の3法則は、古典力学の根幹をなす最も基本的な原理です。1687年にアイザック・ニュートンが著書「プリンキピア（自然哲学の数学的諸原理）」で発表して以来、物理学と工学のあらゆる分野の基盤として機能しています。

航空宇宙工学、ロボット工学、自動車工学、制御工学、建築工学など、「ものが動く」現象を扱うすべての分野は、ニュートンの運動法則の上に成り立っています。特に宇宙工学では、人工衛星の軌道計算やロケットの推力設計は運動方程式なしには成立しません。大学の理工系学部で最初に学ぶ内容ですが、その重要性はいくら強調しても足りません。

本記事の内容

• ニュートンの運動の3法則の定義と意味
• 運動方程式 $\bm{F} = m\bm{a}$ の微分方程式としての理解
• 作用・反作用の法則と内力・外力の区別
• 運動方程式の解法例（自由落下、空気抵抗付き落下）
• Python による空気抵抗付き自由落下のシミュレーション

第1法則（慣性の法則）

定義

外力が作用しない限り、静止している物体は静止し続け、等速直線運動をしている物体は等速直線運動を続ける。

数式で書くと、外力の合力が $\bm{0}$ のとき、

$$ \bm{F} = \bm{0} \quad \Longrightarrow \quad \bm{v} = \text{const.} $$

が成り立つということです。ここで $\bm{v}$ は速度ベクトルを表します。

慣性系とは

第1法則は単なる「当たり前のこと」を述べているのではなく、慣性系の存在を主張する法則です。

慣性系とは、この法則が成り立つ座標系のことです。例えば、地面に固定された座標系は近似的に慣性系とみなせますが、加速するエレベーターの中の座標系は慣性系ではありません。加速するエレベーターの中では、物体に外力を加えなくても加速して見えるからです。

慣性系に対して等速直線運動をしている座標系もまた慣性系です。これはガリレイ変換のもとで運動法則が不変であることを意味し、ガリレイの相対性原理と呼ばれます。つまり、ある慣性系 $S$ に対して速度 $\bm{V}$ で等速直線運動する座標系 $S’$ においても、第1法則は同じ形で成り立ちます。

第2法則（運動方程式）

定義

物体に力が作用すると、物体はその力の方向に加速度を生じる。加速度の大きさは力に比例し、物体の質量に反比例する。

これを数式で表したものが、有名な運動方程式です。

$$ \bm{F} = m\bm{a} $$

ここで、

• $\bm{F}$：物体に作用する合力 $[\text{N}]$
• $m$：物体の質量 $[\text{kg}]$
• $\bm{a}$：加速度 $[\text{m/s}^2]$

微分方程式としての理解

加速度 $\bm{a}$ は位置 $\bm{r}$ の時間に関する2階微分です。

$$ \bm{a} = \frac{d\bm{v}}{dt} = \frac{d^2\bm{r}}{dt^2} $$

よって、運動方程式は2階常微分方程式として書けます。

$$ m\frac{d^2\bm{r}}{dt^2} = \bm{F} $$

力 $\bm{F}$ が位置 $\bm{r}$、速度 $\bm{v}$、時刻 $t$ の関数として与えられれば、この微分方程式を初期条件のもとで解くことにより、物体の運動（軌道）を決定できます。

初期条件としては、初期位置 $\bm{r}(0) = \bm{r}_0$ と初期速度 $\bm{v}(0) = \bm{v}_0$ の2つが必要です。2階微分方程式なので、2つの積分定数を決めるために2つの初期条件が必要になるのです。

運動量による表現

ニュートンが元々述べた形は、運動量 $\bm{p} = m\bm{v}$ を使った表現です。

$$ \bm{F} = \frac{d\bm{p}}{dt} = \frac{d(m\bm{v})}{dt} $$

質量 $m$ が一定の場合、$\bm{F} = m\frac{d\bm{v}}{dt} = m\bm{a}$ となり、先ほどの形に帰着します。ロケットのように質量が変化する系では、運動量の形で考える必要があります。

第3法則（作用・反作用の法則）

定義

物体Aが物体Bに力を及ぼすとき、物体Bは物体Aに大きさが等しく方向が反対の力を及ぼす。

数式で表すと、

$$ \bm{F}_{AB} = -\bm{F}_{BA} $$

ここで $\bm{F}_{AB}$ は物体Aが物体Bに及ぼす力、$\bm{F}_{BA}$ は物体Bが物体Aに及ぼす力です。

内力と外力の区別

第3法則を正しく適用するには、内力と外力を区別することが重要です。

複数の物体からなる系（質点系）を考えるとき、

• 内力: 系の内部の物体同士が及ぼし合う力
• 外力: 系の外部から系の物体に作用する力

と分類されます。第3法則より、内力は必ず作用・反作用の対として現れるため、系全体の運動量への寄与は打ち消し合います。

2つの質点からなる系を例にとると、質点1と質点2の運動方程式はそれぞれ次のように書けます。

$$ m_1 \bm{a}_1 = \bm{F}_1^{\text{ext}} + \bm{f}_{12} $$

$$ m_2 \bm{a}_2 = \bm{F}_2^{\text{ext}} + \bm{f}_{21} $$

ここで $\bm{F}_i^{\text{ext}}$ は外力、$\bm{f}_{12}$ は質点2から質点1への内力、$\bm{f}_{21}$ は質点1から質点2への内力です。第3法則より $\bm{f}_{12} = -\bm{f}_{21}$ なので、両式を足すと、

$$ m_1 \bm{a}_1 + m_2 \bm{a}_2 = \bm{F}_1^{\text{ext}} + \bm{F}_2^{\text{ext}} $$

内力が消えて、系全体の運動は外力のみで決まることがわかります。外力がゼロであれば、系全体の運動量が保存されるという運動量保存則が導かれます。

$$ \bm{F}^{\text{ext}} = \bm{0} \quad \Longrightarrow \quad \bm{p}_1 + \bm{p}_2 = \text{const.} $$

運動方程式の解法例

例1：自由落下（空気抵抗なし）

質量 $m$ の物体を高さ $h$ から初速度 $0$ で落下させる場合を考えます。鉛直下向きを正とすると、運動方程式は次のようになります。

$$ m\frac{d^2 y}{dt^2} = mg $$

両辺を $m$ で割ると、

$$ \frac{d^2 y}{dt^2} = g $$

これを時刻 $t$ で積分します。

$$ \frac{dy}{dt} = gt + C_1 $$

初期条件 $v(0) = 0$ より $C_1 = 0$ なので、

$$ v(t) = \frac{dy}{dt} = gt $$

もう一度積分すると、

$$ y(t) = \frac{1}{2}gt^2 + C_2 $$

初期条件 $y(0) = 0$（落下開始点を原点とする）より $C_2 = 0$ なので、

$$ y(t) = \frac{1}{2}gt^2 $$

$$ v(t) = gt $$

速度は時間に比例して増加し、変位は時間の2乗に比例して増加します。

例2：空気抵抗付き落下

現実の落下では空気抵抗が存在します。速度に比例する空気抵抗を考えると、運動方程式は次のようになります。

$$ m\frac{dv}{dt} = mg – bv $$

ここで $b$ は空気抵抗係数 $[\text{kg/s}]$、鉛直下向きを正としています。$mg$ が重力、$-bv$ が速度に逆らう空気抵抗です。

この微分方程式を解きましょう。変数を整理すると、

$$ \frac{dv}{dt} = g – \frac{b}{m}v $$

$\gamma = \frac{b}{m}$ と置くと、

$$ \frac{dv}{dt} = g – \gamma v $$

これは1階線形常微分方程式です。$u = g – \gamma v$ と置換すると、$\frac{du}{dt} = -\gamma \frac{dv}{dt} = -\gamma u$ となるので、

$$ u(t) = u(0) e^{-\gamma t} $$

初期条件 $v(0) = 0$ より $u(0) = g$ なので、

$$ g – \gamma v(t) = g e^{-\gamma t} $$

$$ v(t) = \frac{g}{\gamma}\left(1 – e^{-\gamma t}\right) = \frac{mg}{b}\left(1 – e^{-\frac{b}{m}t}\right) $$

$t \to \infty$ のとき、$v \to \frac{mg}{b}$ に収束します。この極限速度 $v_{\infty} = \frac{mg}{b}$ を終端速度（terminal velocity）と呼びます。終端速度では、重力と空気抵抗がつり合い、加速度がゼロになります。

位置 $y(t)$ は速度を積分して得られます。

$$ \begin{align} y(t) &= \int_0^t v(t’)\, dt’ \\ &= \frac{mg}{b}\left[t + \frac{m}{b}\left(e^{-\frac{b}{m}t} – 1\right)\right] \end{align} $$

Python で空気抵抗付き自由落下をシミュレーション

空気抵抗付き落下の運動方程式を scipy.integrate.solve_ivp を使って数値的に解き、解析解と比較します。

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

# パラメータ設定
m = 1.0       # 質量 [kg]
g = 9.8       # 重力加速度 [m/s^2]
b = 0.5       # 空気抵抗係数 [kg/s]
t_end = 15.0  # シミュレーション時間 [s]

# 終端速度
v_terminal = m * g / b
print(f"終端速度: {v_terminal:.2f} m/s")

# --- 解析解 ---
t_anal = np.linspace(0, t_end, 500)
v_anal = (m * g / b) * (1 - np.exp(-b / m * t_anal))
y_anal = (m * g / b) * (t_anal + (m / b) * (np.exp(-b / m * t_anal) - 1))

# --- 数値解（solve_ivp） ---
# 状態変数: [y, v]（位置と速度）
def equations(t, state):
    y, v = state
    dydt = v
    dvdt = g - (b / m) * v
    return [dydt, dvdt]

# 初期条件: y(0)=0, v(0)=0
sol = solve_ivp(equations, [0, t_end], [0, 0],
                method='RK45', t_eval=t_anal, rtol=1e-10)

# --- 空気抵抗なしの場合（比較用） ---
v_free = g * t_anal
y_free = 0.5 * g * t_anal**2

# --- 可視化 ---
fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# 速度のプロット
axes[0].plot(t_anal, v_anal, 'b-', linewidth=2, label='Analytical (with drag)')
axes[0].plot(sol.t, sol.y[1], 'r--', linewidth=1.5, label='Numerical (solve_ivp)')
axes[0].axhline(y=v_terminal, color='gray', linestyle=':', linewidth=1,
                label=f'Terminal velocity = {v_terminal:.1f} m/s')
# 空気抵抗なし（制限して表示）
t_mask = t_anal < 3.0
axes[0].plot(t_anal[t_mask], v_free[t_mask], 'g-.', linewidth=1.5,
             label='No drag (free fall)')
axes[0].set_xlabel('Time [s]')
axes[0].set_ylabel('Velocity [m/s]')
axes[0].set_title('Velocity vs Time')
axes[0].legend()
axes[0].grid(True, alpha=0.3)
axes[0].set_ylim(0, v_terminal * 1.2)

# 位置のプロット
axes[1].plot(t_anal, y_anal, 'b-', linewidth=2, label='Analytical (with drag)')
axes[1].plot(sol.t, sol.y[0], 'r--', linewidth=1.5, label='Numerical (solve_ivp)')
axes[1].plot(t_anal[t_mask], y_free[t_mask], 'g-.', linewidth=1.5,
             label='No drag (free fall)')
axes[1].set_xlabel('Time [s]')
axes[1].set_ylabel('Position [m]')
axes[1].set_title('Position vs Time')
axes[1].legend()
axes[1].grid(True, alpha=0.3)
axes[1].set_ylim(0, max(y_anal) * 1.1)

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

上記のコードを実行すると、次のことが確認できます。

• 速度のグラフ: 空気抵抗なしの場合は速度が直線的に増加し続けますが、空気抵抗ありの場合は終端速度 $v_{\infty} = \frac{mg}{b} = 19.6\,\text{m/s}$ に漸近していく様子がわかります
• 位置のグラフ: 空気抵抗なしの場合は放物線的に変位が増加しますが、空気抵抗ありの場合は終端速度に達した後は直線的（等速）に変位が増加します
• 解析解と数値解の一致: solve_ivp による数値解（赤い破線）が解析解（青い実線）と正確に重なっていることが確認できます

数値解法は、解析的に解けない複雑な力（速度の2乗に比例する空気抵抗など）に対しても適用できるため、工学において非常に重要な手法です。

まとめ

本記事では、ニュートンの運動の3法則について解説しました。

• 第1法則（慣性の法則）: 外力がなければ物体は等速直線運動を続けます。慣性系という基準座標系の存在を主張する法則です
• 第2法則（運動方程式）: $\bm{F} = m\bm{a}$ は2階常微分方程式であり、力と初期条件が与えられれば物体の運動を完全に決定できます
• 第3法則（作用・反作用の法則）: 力は必ず対で現れ、内力は系全体の運動量を変化させません。これにより運動量保存則が導かれます
• 空気抵抗付き落下: 運動方程式を解くと終端速度 $v_{\infty} = \frac{mg}{b}$ の存在が示され、数値計算によってもそれを確認できます

次のステップとして、以下のトピックもぜひ学んでみてください。