ジェット旅客機が巡航高度1万メートルを時速900 kmで飛んでいるとき、機体の先端付近では空気がぎゅっと押し縮められています。一方、水道管の中を流れる水や、そよ風が吹く穏やかな午後の空気は、ほとんど密度が変化しません。同じ「流体の流れ」でも、速度が大きくなると流体は圧縮され、密度が場所ごとに変わるという全く別の物理が顔を出します。この「密度が変わる流れ」を扱うのが圧縮性流体力学（compressible fluid dynamics）であり、その中心的なパラメータがマッハ数（Mach number）です。
日常の低速流れではベルヌーイの定理や連続の式をそのまま使えばよいのですが、高速気流ではそれだけでは不十分です。音速を超えると衝撃波が発生し、温度・圧力・密度が不連続に跳ぶなど、低速流れでは想像もつかない現象が起こります。圧縮性流体力学を理解することで、以下のような幅広い応用分野に対処できるようになります。
- 航空宇宙工学: 遷音速・超音速域での翼型設計、ラバルノズルによるロケットエンジンの推力最適化
- ターボ機械: ガスタービンの圧縮機・タービン翼列における衝撃波の制御
- 風洞実験: 超音速風洞のノズル設計とテストセクション条件の計算
- 配管・プラント工学: 高圧ガスの配管内流れにおける臨界流量の評価
- 天体物理学: 超新星爆発の衝撃波伝播、太陽風と地球磁気圏の相互作用
本記事では、「流体の圧縮性とは何か」「それをどう判定するのか」という基本的な問いから出発して、音速の導出、マッハ数の定義と分類、等エントロピー流れの基本関係式、よどみ点の性質、臨界条件、面積-速度関係まで、圧縮性流体力学の基礎を一貫して解説します。最後にPythonで等エントロピー関係式やノズル流れを可視化し、理論の理解を数値で裏付けます。
本記事の内容
- 圧縮性の判定基準（マッハ数0.3が分水嶺となる理由）
- 理想気体中の音速の導出と物理的意味
- マッハ数の定義と流れの分類（亜音速・遷音速・超音速・極超音速）
- 等エントロピー流れの基本関係式（温度比・圧力比・密度比）
- よどみ点（全温・全圧・全密度）の意味と使い方
- 臨界条件（M = 1 でのスロート条件）と臨界比
- 面積-速度関係と収束-発散ノズルの物理
- Pythonによる4種類の可視化と計算ツール
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
- ベルヌーイの定理 — 非圧縮性流体のエネルギー保存則。圧縮性流れではこれを拡張する必要があります
- 連続の式 — 質量保存の定式化。圧縮性流れでは密度変化を陽に含めます
- 熱力学第1法則 — 内部エネルギー・仕事・熱の関係。流体のエネルギー方程式の出発点です
- エントロピー — 等エントロピー流れの前提条件を理解するために必要です
- 理想気体の状態方程式 — $pV = nRT$ と理想気体の基本性質
- 比熱と比熱比 — 比熱 $C_v$, $C_p$ と比熱比 $\gamma$ の定義と関係
圧縮性とは何か
非圧縮性流体と圧縮性流体の違い
流体の性質で学んだように、流体には「圧縮性」という性質があります。ここでいう圧縮性とは、圧力の変化に応じて流体の密度が変化する性質のことです。
日常的な感覚で考えてみましょう。注射器の先を指で塞いで、ピストンを押してみてください。中に空気が入っていれば、力を入れるとピストンは少し動きます。空気は圧縮されて体積が減り、密度が上がったのです。一方、中に水が入っていたらどうでしょうか。どんなに力を入れてもピストンはほとんど動きません。水はほとんど圧縮されないのです。
しかし、流体力学における「圧縮性」と「非圧縮性」の区別は、流体そのものの性質だけでは決まりません。流れの速度が本質的に重要です。空気は圧縮性のある流体ですが、そよ風程度の低速流れであれば密度変化はごく小さく、非圧縮性として扱っても実用上問題ありません。逆に、空気であっても高速で流れると密度変化が無視できなくなり、圧縮性を考慮しなければならなくなります。
つまり、圧縮性が重要かどうかは「流体の種類」と「流れの速度」の両方で決まるのです。この判定基準を与えるのがマッハ数であり、その閾値がおよそ $M = 0.3$ です。
圧縮性の定量的な判定
それでは、なぜ $M = 0.3$ が分水嶺になるのでしょうか。これを定量的に確認してみましょう。
非圧縮性のベルヌーイの定理を思い出してください。定常・非粘性・非圧縮の流れでは、流線に沿って次の関係が成り立ちます。
$$ \begin{equation} p + \frac{1}{2}\rho V^2 = p_0 \end{equation} $$
ここで $p$ は静圧、$\rho$ は密度、$V$ は流速、$p_0$ はよどみ点圧力（全圧）です。この式を変形すると、圧力変化の相対的な大きさは次のように書けます。
$$ \begin{equation} \frac{p_0 – p}{p} = \frac{\rho V^2}{2p} \end{equation} $$
理想気体の状態方程式 $p = \rho R T$ と音速 $a = \sqrt{\gamma R T}$（後述）を使うと、$\rho / p = 1/(RT) = \gamma / a^2$ となるため、次が得られます。
$$ \begin{equation} \frac{p_0 – p}{p} = \frac{\gamma}{2} M^2 \end{equation} $$
ここで $M = V/a$ がマッハ数です。一方、圧縮性を正しく考慮した等エントロピー関係式（後のセクションで導出します）を使うと、全圧と静圧の比は次のようになります。
$$ \begin{equation} \frac{p_0}{p} = \left(1 + \frac{\gamma – 1}{2} M^2\right)^{\gamma/(\gamma – 1)} \end{equation} $$
空気（$\gamma = 1.4$）で $M = 0.3$ のとき、この式に数値を代入すると $p_0/p \approx 1.0644$ です。一方、非圧縮のベルヌーイの式からは $p_0/p = 1 + \gamma M^2/2 = 1.063$ となります。両者の差はわずか0.13%程度です。つまり $M \leq 0.3$ の範囲では、密度変化を無視しても圧力の計算誤差は1%未満に収まるのです。$M = 0.3$ を超えると誤差は急速に大きくなり、$M = 0.5$ では約2.5%、$M = 0.8$ では約8%にまで成長します。
工学的には、$M < 0.3$ の流れは非圧縮性として扱ってよいというのが広く受け入れられた目安です。$M > 0.3$ の流れでは圧縮性効果を無視できず、本記事で扱う圧縮性流体力学の枠組みが必要になります。
では、この判定基準の核心にある「音速」とは何なのか、そしてなぜ音速がこれほど重要な役割を果たすのか、次のセクションで掘り下げていきましょう。
音速の導出
音速の物理的意味
音速とは文字通り「音が伝わる速さ」ですが、流体力学においてはもっと本質的な意味を持っています。音速とは、微小な圧力の擾乱（じょうらん）が流体中を伝播する速度です。
手を叩くと空気中に圧力の波（音波）が生じます。この波が伝わる速さが音速です。重要なのは、「流体の状態が変化したという情報」が伝わる速度が音速であるという点です。たとえば、パイプの中を流れる気体の下流端で弁を急に閉じると、「弁が閉まった」という情報は音速で上流に伝わります。もし流れの速度が音速を超えていたら、この情報は上流に伝わることができません。この事実は、超音速流れの振る舞いを理解する上で決定的に重要です。
微小擾乱の伝播速度としての導出
音速を熱力学的な量で表す式を導出しましょう。非常に細い管の中を音波が伝わる状況を考えます。管の中の静止した気体（圧力 $p$、密度 $\rho$）に微小な圧力擾乱が速度 $a$ で右に伝播していくとします。
この波と一緒に動く座標系に移ると、気体が速度 $a$ で左から右に流れてきて、波面を通過した後、速度が $a + dV$、圧力が $p + dp$、密度が $\rho + d\rho$ に変化するという定常問題になります。
まず連続の式（質量保存則）を適用します。管の断面積を $A$ とすると、波面の前後で質量流量が等しいので次が成立します。
$$ \begin{equation} \rho a A = (\rho + d\rho)(a + dV) A \end{equation} $$
右辺を展開し、微小量の二次の項 $d\rho \cdot dV$ を無視すると次のようになります。
$$ \begin{equation} \rho \, dV + a \, d\rho = 0 \end{equation} $$
ここから $dV = -a \, d\rho / \rho$ が得られます。
次に、波面に運動量の法則を適用します。波面に作用する正味の圧力力は $(p – (p+dp))A = -dp \cdot A$ であり、これが質量流量 $\rho a A$ に速度変化 $dV$ を掛けたものに等しいので次が成立します。
$$ \begin{equation} -dp \cdot A = \rho a A \cdot dV \end{equation} $$
$A$ で割り、先ほど得た $dV = -a \, d\rho / \rho$ を代入すると次のようになります。
$$ -dp = \rho a \cdot \left(-\frac{a \, d\rho}{\rho}\right) = -a^2 d\rho $$
したがって、音速 $a$ は次の式で定まります。
$$ \begin{equation} a^2 = \frac{dp}{d\rho} \end{equation} $$
この式は一般的に成り立つ結果であり、流体の種類を問いません。音速は「圧力が密度変化にどれだけ敏感に応答するか」を表すものだとわかります。密度変化に対して圧力が大きく変化する（硬い）流体ほど音速が大きいのです。水が空気よりも音速が大きい（水中: 約1500 m/s、空気中: 約340 m/s）のは、水の方がはるかに「硬い」からです。
理想気体中の音速
音波の伝播は非常に速い現象であり、微小振幅の音波では各流体要素と周囲との間に熱のやりとりがほとんどありません。したがって、音波の伝播過程は等エントロピー過程（断熱可逆過程）として扱えます。
エントロピーが一定の理想気体では、次の関係が成り立ちます。
$$ \begin{equation} p \rho^{-\gamma} = \text{const} \end{equation} $$
ここで $\gamma = C_p / C_v$ は比熱比です。この関係式の両辺を $\rho$ で微分すると次が得られます。
$$ \frac{dp}{d\rho} = \gamma \frac{p}{\rho} $$
理想気体の状態方程式 $p = \rho R T$（$R$ は気体のガス定数: 比気体定数 $R = R_u / M_w$）を代入すると、理想気体中の音速は次のようになります。
$$ \begin{equation} a = \sqrt{\gamma R T} \end{equation} $$
この結果は非常にすっきりしています。理想気体中の音速は温度だけで決まり、圧力や密度には直接依存しません。高度10,000 mの上空では気温が約-50°C（223 K）なので、空気（$\gamma = 1.4$, $R = 287$ J/(kg K)）の音速は $a = \sqrt{1.4 \times 287 \times 223} \approx 299$ m/s です。海面上の標準大気（15°C = 288 K）では $a \approx 340$ m/s となります。
ニュートンは音速を等温過程で計算し $a = \sqrt{RT}$（$\gamma$ なし）という式を導きましたが、これは実験値より約15%低い値を与えます。ラプラスが等エントロピー過程を正しく適用して $\gamma$ を含む式を導き、実験値との一致を得ました。この歴史的エピソードは、音波の伝播が等温ではなく等エントロピーであることの重要な証拠です。
音速の表式が得られたところで、いよいよ圧縮性流体力学の中心的なパラメータであるマッハ数を正式に定義し、流れの分類を見ていきましょう。
マッハ数の定義と流れの分類
マッハ数の定義
マッハ数 $M$ は、流速 $V$ とその地点の局所音速 $a$ の比として定義されます。
$$ \begin{equation} M = \frac{V}{a} \end{equation} $$
注意すべきは、音速 $a = \sqrt{\gamma R T}$ は温度に依存するため、場所によって異なるという点です。マッハ数は流速だけでなく局所温度も含んだ無次元数であり、流れの圧縮性の度合いを総合的に表現しています。
マッハ数が物理的に意味するところを直感的に捉えてみましょう。$M < 1$ の場合、流体は音速（情報伝達速度）よりも遅く流れています。このとき、下流で起きた変化の情報は上流に伝わることができます。ちょうど、川の流れよりも速く泳げる魚が上流に向かって情報を伝達できるようなものです。一方、$M > 1$ では流体が音速を超えているため、下流の情報は上流に届きません。川の流れが速すぎて、魚がどんなに頑張っても上流に泳ぎ戻れないのと同じです。
この「情報が上流に伝わるかどうか」の違いが、亜音速流れと超音速流れの決定的な違いを生み出します。
流れの分類
マッハ数に基づいて、流れは以下のように分類されます。
| 分類 | マッハ数の範囲 | 特徴 |
|---|---|---|
| 非圧縮性流れ | $M < 0.3$ | 密度変化は無視可能。ベルヌーイの定理がそのまま適用可能 |
| 亜音速流れ | $0.3 \leq M < 0.8$ | 圧縮性効果が顕在化。密度変化を考慮する必要あり |
| 遷音速流れ | $0.8 \leq M < 1.2$ | 流れの一部が音速を超え、局所的に衝撃波が発生し得る |
| 超音速流れ | $1.2 \leq M < 5$ | 流れ全体が音速を超過。衝撃波・膨張波が支配的 |
| 極超音速流れ | $M \geq 5$ | 衝撃波背後の温度上昇で気体の解離・電離が生じる |
遷音速域（$0.8 \leq M < 1.2$）は特に厄介な領域です。主流のマッハ数が1未満であっても、物体（たとえば翼の上面）の周りでは流れが局所的に加速されてマッハ数が1を超え、衝撃波が発生することがあります。ジェット旅客機の巡航マッハ数が0.78~0.85程度に設定されているのは、この遷音速域の衝撃波による抗力増大（波動抗力）を避けるためです。
極超音速域（$M \geq 5$）では、衝撃波背後の温度が数千度に達します。このような高温では、空気の分子が解離したりイオン化したりするため、比熱比 $\gamma$ が一定という理想気体の仮定そのものが崩れます。大気圏再突入カプセルやスクラムジェットエンジンの設計では、こうした高温気体力学（real gas effects）の考慮が不可欠です。
マッハ数による流れの分類を押さえたところで、次に圧縮性流体力学の中核となる等エントロピー流れの関係式を導出していきます。この関係式は、流れの中の温度・圧力・密度をマッハ数だけで結びつけてくれる非常に強力なツールです。
等エントロピー流れの基本関係式
等エントロピー流れの条件
等エントロピー流れ（isentropic flow）とは、エントロピーが一定に保たれる流れのことです。これは以下の2つの条件が同時に満たされるとき実現します。
- 断熱: 流体と外部との間に熱のやりとりがない
- 可逆: 粘性による摩擦散逸や衝撃波による不可逆性がない
現実の流れは完全に等エントロピーではありませんが、衝撃波のない滑らかな流れの領域では非常に良い近似になります。特に、ノズルやディフューザの内部流れの解析で等エントロピー関係式は中心的な役割を果たします。
エネルギー方程式とよどみ点
等エントロピー関係式を導出する出発点は、定常流のエネルギー方程式です。熱力学第1法則を定常流れに適用すると、流線に沿って次の関係が成り立ちます。
$$ \begin{equation} h + \frac{V^2}{2} = h_0 = \text{const} \end{equation} $$
ここで $h$ は比エンタルピー（単位質量あたりのエンタルピー）、$V$ は流速、$h_0$ は全エンタルピー（よどみエンタルピー）です。この式は、流体の運動エネルギーとエンタルピーの合計が流線に沿って一定であることを表しています。
大雑把に言うと、これは圧縮性流体版のベルヌーイの定理です。非圧縮性のベルヌーイの定理では $p + \frac{1}{2}\rho V^2 = \text{const}$ ですが、圧縮性流れでは圧力の代わりにエンタルピー（圧力と内部エネルギーの両方を含む量）が登場します。
よどみ点（全温・全圧）の物理的意味
よどみ点（stagnation point）とは、流れが断熱的に速度ゼロまで減速されたときの仮想的な状態です。よどみ点での温度 $T_0$、圧力 $p_0$、密度 $\rho_0$ をそれぞれ全温（total temperature）、全圧（total pressure）、全密度（total density）と呼びます。
イメージとしては、流れに正面から突き刺さるように置いた管（ピトー管のようなもの）の先端で流体が完全に停止する点がよどみ点です。流体が持っていた運動エネルギーの全てがエンタルピー（温度）の上昇に変換されます。
理想気体では $h = C_p T$ なので、エネルギー方程式は次のように書けます。
$$ \begin{equation} C_p T + \frac{V^2}{2} = C_p T_0 \end{equation} $$
両辺を $C_p T$ で割ると次が得られます。
$$ 1 + \frac{V^2}{2 C_p T} = \frac{T_0}{T} $$
ここで $C_p = \gamma R / (\gamma – 1)$ と $a^2 = \gamma R T$ を使うと、$V^2/(2 C_p T) = V^2 / (2 \cdot \frac{\gamma R}{\gamma-1} \cdot T) = \frac{(\gamma-1) V^2}{2 \gamma R T} = \frac{(\gamma-1) V^2}{2 a^2} = \frac{\gamma – 1}{2} M^2$ と変形できます。
したがって、全温比（温度比）は次のようになります。
$$ \begin{equation} \frac{T_0}{T} = 1 + \frac{\gamma – 1}{2} M^2 \end{equation} $$
これは等エントロピー関係式の中で最も基本的な式です。物理的な意味は明快です。マッハ数が大きいほど、流体が停止したときの温度上昇が大きくなります。超音速旅客機コンコルド（$M \approx 2$）の機体表面温度は巡航中に100°C以上に達しましたが、これはまさにこの関係式が示す現象です。$\gamma = 1.4$ の空気で $M = 2$ とすると、$T_0/T = 1 + 0.2 \times 4 = 1.8$ となり、外気温が $-57$°C（216 K）でも、よどみ点温度は $T_0 = 389$ K（116°C）にも達します。
圧力比と密度比の導出
等エントロピー過程では、エントロピー一定の条件から温度と圧力、温度と密度の間に次の関係が成り立ちます。
$$ \frac{p_0}{p} = \left(\frac{T_0}{T}\right)^{\gamma/(\gamma-1)}, \qquad \frac{\rho_0}{\rho} = \left(\frac{T_0}{T}\right)^{1/(\gamma-1)} $$
先ほど導いた全温比の式を代入すると、全圧比と全密度比が得られます。
$$ \begin{equation} \frac{p_0}{p} = \left(1 + \frac{\gamma – 1}{2} M^2\right)^{\gamma/(\gamma – 1)} \end{equation} $$
$$ \begin{equation} \frac{\rho_0}{\rho} = \left(1 + \frac{\gamma – 1}{2} M^2\right)^{1/(\gamma – 1)} \end{equation} $$
これら3つの式（全温比・全圧比・全密度比）が、等エントロピー流れの基本関係式です。マッハ数と比熱比 $\gamma$ さえわかれば、流れの中の任意の点での温度・圧力・密度を、よどみ点の値を基準にして求めることができます。
$M = 0$ を代入すると、すべての比が1になり、静止状態ではよどみ点の値と一致することが確認できます。$M$ が大きくなると、温度比・圧力比・密度比はいずれも急速に大きくなりますが、その増大速度は $\gamma$ の指数の違いにより異なります。圧力比が最も速く増大し、次に密度比、最も穏やかなのが温度比です。
この等エントロピー関係式がどのような振る舞いを見せるか、ここでPythonを使って可視化してみましょう。
Pythonで等エントロピー関係式を可視化する
等エントロピー関係式（温度比・圧力比・密度比）がマッハ数に対してどのように変化するかを視覚的に把握するために、グラフを描画してみます。
import numpy as np
import matplotlib.pyplot as plt
# パラメータ設定
gamma = 1.4 # 空気の比熱比
M = np.linspace(0, 5, 500) # マッハ数 0〜5
# 等エントロピー関係式
# 温度比 T/T0、圧力比 p/p0、密度比 rho/rho0
T_ratio = (1 + (gamma - 1) / 2 * M**2) ** (-1)
p_ratio = (1 + (gamma - 1) / 2 * M**2) ** (-gamma / (gamma - 1))
rho_ratio = (1 + (gamma - 1) / 2 * M**2) ** (-1 / (gamma - 1))
# グラフ描画
fig, ax = plt.subplots(figsize=(10, 6))
ax.plot(M, T_ratio, label=r"$T/T_0$（温度比）", linewidth=2, color="#2196F3")
ax.plot(M, p_ratio, label=r"$p/p_0$（圧力比）", linewidth=2, color="#F44336")
ax.plot(M, rho_ratio, label=r"$\rho/\rho_0$（密度比）", linewidth=2, color="#4CAF50")
# M=1 の目印
ax.axvline(x=1.0, color="gray", linestyle="--", linewidth=0.8, alpha=0.7)
ax.text(1.02, 0.92, "M = 1", fontsize=10, color="gray")
# M=0.3 の目印（圧縮性判定境界）
ax.axvline(x=0.3, color="orange", linestyle=":", linewidth=0.8, alpha=0.7)
ax.text(0.32, 0.92, "M = 0.3", fontsize=9, color="orange")
ax.set_xlabel("Mach number $M$", fontsize=13)
ax.set_ylabel("Ratio to stagnation value", fontsize=13)
ax.set_title("Isentropic Flow Relations ($\\gamma = 1.4$)", fontsize=14)
ax.legend(fontsize=12, loc="center right")
ax.set_xlim(0, 5)
ax.set_ylim(0, 1.05)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
このグラフから、以下の重要な特徴が読み取れます。
- $M = 0.3$ 以下では全ての比がほぼ1に近い — $T/T_0 \approx 0.982$, $p/p_0 \approx 0.939$, $\rho/\rho_0 \approx 0.956$ であり、よどみ点の値からほとんど変化していません。これが「$M < 0.3$ では圧縮性効果が小さい」という判定基準の数値的な根拠です。
- 圧力比が最も急激に減少する — 赤い曲線（圧力比）が最も速く0に近づいています。$M = 3$ では $p/p_0 \approx 0.027$ とよどみ点圧力のわずか2.7%しかありません。これは圧力比の指数が $\gamma/(\gamma-1) = 3.5$ と最も大きいことに対応しています。
- 温度比の減少は最も穏やか — 青い曲線（温度比）は他の2つよりもなだらかに減少します。$M = 5$ でも $T/T_0 \approx 0.167$ と、まだ17%程度の温度を保っています。
これらの関係式は、ノズル流れの設計や風洞試験の条件設定など、圧縮性流体力学のあらゆる場面で繰り返し使われる基本ツールです。
次に、等エントロピー関係式の中で特別な意味を持つ $M = 1$（臨界条件）について詳しく見ていきましょう。
臨界条件（M = 1）
臨界状態とは
等エントロピー流れにおいて $M = 1$ となる状態を臨界状態と呼び、そのときの温度・圧力・密度を臨界温度 $T^*$、臨界圧力 $p^*$、臨界密度 $\rho^*$ と表します（上付きのアスタリスク $*$ が臨界状態を示す記号です）。
臨界状態は、ラバルノズルのスロート（最小断面積部）で実現される状態であり、流量の上限を決定するという点で工学的に非常に重要です。なぜ臨界状態が流量の上限を決めるのかは、面積-速度関係のセクションで詳しく見ていきます。
臨界比の計算
等エントロピー関係式に $M = 1$ を代入すると、臨界比が得られます。
$$ \begin{equation} \frac{T^*}{T_0} = \frac{2}{\gamma + 1} \end{equation} $$
$$ \begin{equation} \frac{p^*}{p_0} = \left(\frac{2}{\gamma + 1}\right)^{\gamma/(\gamma – 1)} \end{equation} $$
$$ \begin{equation} \frac{\rho^*}{\rho_0} = \left(\frac{2}{\gamma + 1}\right)^{1/(\gamma – 1)} \end{equation} $$
空気（$\gamma = 1.4$）の場合に数値を計算してみましょう。
$$ \frac{T^*}{T_0} = \frac{2}{2.4} = 0.8333 $$
$$ \frac{p^*}{p_0} = \left(\frac{2}{2.4}\right)^{3.5} = 0.5283 $$
$$ \frac{\rho^*}{\rho_0} = \left(\frac{2}{2.4}\right)^{2.5} = 0.6340 $$
臨界圧力比 $p^*/p_0 \approx 0.528$ は特に重要な数値です。これは、ノズル出口圧力と全圧の比がこの値以下になると、スロートでマッハ数が1に達し、流れがチョーク（choking）するということを意味します。チョーク状態では、ノズルの上流側から見た質量流量が最大に達し、下流の圧力をさらに下げても流量は増えなくなります。
この $p^*/p_0 \approx 0.528$ という値は、圧力容器やボイラーの安全弁の設計、ガスタービンのタービンノズルの設計など、多くの工学的場面で使われます。安全弁が開いたときにスロート部がチョークするかどうかの判定に、まさにこの臨界圧力比が用いられるのです。
臨界速度
$M = 1$ での流速（臨界速度）$V^* = a^*$ も重要な量です。臨界温度 $T^* = 2T_0/(\gamma+1)$ を音速の式に代入すると次が得られます。
$$ \begin{equation} a^* = \sqrt{\gamma R T^*} = \sqrt{\frac{2\gamma R T_0}{\gamma + 1}} \end{equation} $$
全温 $T_0$ が与えられれば臨界速度が一意に定まります。臨界速度は、等エントロピー流れの中で「流速がこの値を上回れば超音速、下回れば亜音速」という基準線を与えます。
臨界条件はノズルのスロートで実現されるものですが、なぜスロートでなければ $M = 1$ が実現できないのでしょうか。これを理解するためには、管路の断面積と流速の関係を調べる必要があります。
面積-速度関係
準一次元流れの仮定
ノズルやディフューザのように断面積が軸方向に変化する管路を考えます。流れの方向を $x$ 軸に取り、断面積を $A(x)$ とします。各断面における流れの状態量（速度、圧力、密度など）は断面内で一様と仮定します。これを準一次元流れの仮定と呼びます。実際のノズル内の流れは二次元・三次元の複雑な構造を持ちますが、管路の断面積変化がなだらかな場合、準一次元近似は良い精度を与えます。
面積-速度関係の導出
定常・等エントロピーの準一次元流れに対して、連続の式は次のように書けます。
$$ \begin{equation} \rho V A = \text{const} = \dot{m} \end{equation} $$
ここで $\dot{m}$ は質量流量です。この式の両辺の対数を取って微分すると次のようになります。
$$ \begin{equation} \frac{d\rho}{\rho} + \frac{dV}{V} + \frac{dA}{A} = 0 \end{equation} $$
一方、粘性と体積力を無視した一次元運動量方程式（オイラーの方程式）は次の通りです。
$$ \begin{equation} dp + \rho V \, dV = 0 \end{equation} $$
等エントロピー過程では $dp = a^2 d\rho$ が成り立つので、これを運動量方程式に代入すると次が得られます。
$$ a^2 d\rho + \rho V \, dV = 0 $$
$d\rho$ について解くと $d\rho/\rho = -V \, dV / a^2 = -M^2 \, dV / V$ です。これを連続の式の微分形に代入すると次の関係が得られます。
$$ -M^2 \frac{dV}{V} + \frac{dV}{V} + \frac{dA}{A} = 0 $$
整理すると、面積-速度関係（area-velocity relation）が得られます。
$$ \begin{equation} \frac{dA}{A} = (M^2 – 1) \frac{dV}{V} \end{equation} $$
この式は非常にコンパクトですが、圧縮性流体力学で最も重要な結論の一つを含んでいます。
面積-速度関係の物理的意味
面積-速度関係の含意を場合分けして整理しましょう。
亜音速（$M < 1$）の場合: $M^2 – 1 < 0$ なので、$dA$ と $dV$ は逆符号になります。つまり、断面積が減少すると速度が増加し、断面積が増加すると速度が減少します。これは非圧縮性流体と同じ直感に合う振る舞いです。ホースの先を絞ると水が勢いよく出る、あの現象です。
超音速（$M > 1$）の場合: $M^2 – 1 > 0$ なので、$dA$ と $dV$ は同符号になります。つまり、断面積が増加すると速度が増加し、断面積が減少すると速度が減少します。これは非圧縮性の直感に反する、圧縮性流体に特有の振る舞いです。
なぜ超音速では断面積を広げると加速するのでしょうか。密度変化を考えると理解できます。超音速流れでは流れが加速すると密度が急激に低下します。連続の式 $\rho V A = \text{const}$ において、密度の低下が速度の増加を上回るため、質量流量を保つには断面積を増やす必要があるのです。
音速（$M = 1$）の場合: $M^2 – 1 = 0$ なので、$dV \neq 0$ であれば $dA = 0$ でなければなりません。つまり、$M = 1$ は断面積が極値を取る点（$dA/dx = 0$）でしか実現しないのです。断面積の極値は通常、最小値（スロート）か最大値のどちらかですが、流れが加速して音速に達する状況を考えると、スロートで $M = 1$ が実現することがわかります。
この結果が、収束-発散ノズル（ラバルノズル）の動作原理の核心です。亜音速流れを超音速まで加速するには、まず断面積を縮小して亜音速流れを加速し、スロートで $M = 1$ を実現してから、断面積を拡大してさらに超音速まで加速するという二段階の過程が必要なのです。単に断面積を縮小するだけの収束ノズルでは、スロートで $M = 1$ に達した時点で加速が止まってしまいます。
面積比とマッハ数の関係
準一次元等エントロピー流れでは、任意の断面積 $A$ とスロート面積 $A^*$（$M = 1$ の断面積）の比がマッハ数の関数として表せます。連続の式 $\rho V A = \rho^* a^* A^*$ から次を得ます。
$$ \frac{A}{A^*} = \frac{\rho^* a^*}{\rho V} = \frac{\rho^*}{\rho} \cdot \frac{a^*}{V} $$
等エントロピー関係式を用いて $\rho^*/\rho$ と $a^*/V$ を $M$ で表すと、次の面積-マッハ数関係が得られます。
$$ \begin{equation} \frac{A}{A^*} = \frac{1}{M} \left[\frac{2}{\gamma + 1}\left(1 + \frac{\gamma – 1}{2} M^2\right)\right]^{(\gamma + 1) / [2(\gamma – 1)]} \end{equation} $$
この関係は $M = 1$ で $A/A^* = 1$（最小値）を取り、$M < 1$ でも $M > 1$ でも $A/A^* > 1$ となります。つまり、一つの面積比 $A/A^*$ に対して、亜音速解と超音速解の2つのマッハ数が対応します。ノズル流れの設計では、収束部は亜音速解、発散部は超音速解を使います。
面積-マッハ数関係をPythonで可視化して、この二重値の構造を確認してみましょう。
Pythonでマッハ数と面積比の関係を可視化する
面積-マッハ数関係がどのような曲線を描くかを確認し、亜音速解と超音速解が存在することを視覚的に把握します。
import numpy as np
import matplotlib.pyplot as plt
gamma = 1.4
# 亜音速域と超音速域を分けて計算
M_sub = np.linspace(0.01, 1.0, 300) # 亜音速
M_sup = np.linspace(1.0, 5.0, 500) # 超音速
def area_ratio(M, gamma=1.4):
"""面積比 A/A* を計算"""
exponent = (gamma + 1) / (2 * (gamma - 1))
term = (2 / (gamma + 1)) * (1 + (gamma - 1) / 2 * M**2)
return (1 / M) * term**exponent
A_sub = area_ratio(M_sub)
A_sup = area_ratio(M_sup)
fig, ax = plt.subplots(figsize=(10, 6))
ax.plot(M_sub, A_sub, color="#2196F3", linewidth=2.5, label="Subsonic branch")
ax.plot(M_sup, A_sup, color="#F44336", linewidth=2.5, label="Supersonic branch")
ax.plot(1.0, 1.0, "ko", markersize=8, zorder=5, label="Throat (M = 1)")
# 面積比 = 2 の水平線で二重値を表示
A_target = 2.0
ax.axhline(y=A_target, color="gray", linestyle="--", linewidth=0.8, alpha=0.6)
ax.text(0.05, A_target + 0.1, f"$A/A^* = {A_target}$", fontsize=10, color="gray")
ax.set_xlabel("Mach number $M$", fontsize=13)
ax.set_ylabel("Area ratio $A / A^*$", fontsize=13)
ax.set_title("Area-Mach Number Relation ($\\gamma = 1.4$)", fontsize=14)
ax.legend(fontsize=12)
ax.set_xlim(0, 5)
ax.set_ylim(0, 6)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
このグラフから、以下の重要な特徴が読み取れます。
- $M = 1$ で面積比が最小値1を取る — これはスロートに対応します。面積比はスロートで最小（$A/A^* = 1$）であり、亜音速側にも超音速側にも面積比は単調に増加します。
- 一つの $A/A^*$ に対して2つのマッハ数が存在する — 灰色の破線 $A/A^* = 2$ と曲線の交点が2つあり、亜音速側（$M \approx 0.31$）と超音速側（$M \approx 1.69$）のそれぞれに解が存在します。ノズル流れでどちらの解が実現するかは、上流・下流の境界条件（圧力比）で決まります。
- 超音速側では面積比が急速に増大する — $M = 5$ では $A/A^* \approx 25.0$ となり、スロート面積の25倍もの出口面積が必要です。これが超音速ノズル（ロケットエンジンのベルノズルなど）が末広がりの大きな形状をしている理由です。
ここまでで等エントロピー流れの基本関係式と面積-マッハ数関係を見てきました。次に、等エントロピー関係式を使った計算を手軽に行えるPythonの計算ツールを作成し、さまざまなマッハ数での状態量を数値的に確認してみましょう。
Pythonで等エントロピー関係式の計算ツールを作成する
等エントロピー関係式を使って、マッハ数を入力すると温度比・圧力比・密度比・面積比をまとめて計算するツールを作成します。代表的なマッハ数での数値を表形式で確認します。
import numpy as np
def isentropic_relations(M, gamma=1.4):
"""
等エントロピー関係式の各比を計算する
Parameters:
M : float or array - マッハ数
gamma : float - 比熱比
Returns:
dict : T/T0, p/p0, rho/rho0, A/A* の各比
"""
factor = 1 + (gamma - 1) / 2 * M**2
T_ratio = factor ** (-1)
p_ratio = factor ** (-gamma / (gamma - 1))
rho_ratio = factor ** (-1 / (gamma - 1))
# 面積比 A/A*
exponent = (gamma + 1) / (2 * (gamma - 1))
A_ratio = (1 / M) * ((2 / (gamma + 1)) * factor) ** exponent
return {
"T/T0": T_ratio,
"p/p0": p_ratio,
"rho/rho0": rho_ratio,
"A/A*": A_ratio,
}
# 代表的なマッハ数で計算
mach_numbers = [0.0, 0.3, 0.5, 0.8, 1.0, 1.5, 2.0, 3.0, 5.0]
print(f"{'M':>6s} {'T/T0':>10s} {'p/p0':>10s} {'rho/rho0':>10s} {'A/A*':>10s}")
print("-" * 52)
for M_val in mach_numbers:
if M_val == 0.0:
# M=0 では A/A* が発散するため別扱い
print(f"{M_val:6.1f} {'1.0000':>10s} {'1.0000':>10s} {'1.0000':>10s} {'inf':>10s}")
else:
result = isentropic_relations(M_val)
print(f"{M_val:6.1f} {result['T/T0']:10.4f} {result['p/p0']:10.4f} "
f"{result['rho/rho0']:10.4f} {result['A/A*']:10.4f}")
実行すると、以下のような出力が得られます。
M T/T0 p/p0 rho/rho0 A/A*
----------------------------------------------------
0.0 1.0000 1.0000 1.0000 inf
0.3 0.9823 0.9395 0.9564 2.0351
0.5 0.9524 0.8430 0.8852 1.3398
0.8 0.8865 0.6560 0.7400 1.0382
1.0 0.8333 0.5283 0.6340 1.0000
1.5 0.6897 0.2724 0.3950 1.1762
2.0 0.5556 0.1278 0.2300 1.6875
3.0 0.3571 0.0272 0.0762 4.2346
5.0 0.1667 0.0019 0.0113 25.0000
この数値表から、いくつかの重要な観察ができます。
- $M = 0.3$ では密度変化はわずか4.4% — $\rho/\rho_0 = 0.956$ なので、密度はよどみ点から約4.4%しか減少していません。これが $M < 0.3$ を非圧縮性として扱える根拠です。
- $M = 1$（臨界状態）で $p/p_0 = 0.528$ — 先ほど解析的に計算した臨界圧力比と一致しています。スロートでの圧力はよどみ点圧力の約53%まで低下します。
- $M = 5$ では圧力がよどみ点のわずか0.19% — 超音速流れでは圧力が劇的に低下します。また、面積比が25に達し、スロートの25倍の断面積が必要です。
このツールは、ノズルの設計計算や風洞実験の条件設定など、圧縮性流体力学の実務で日常的に利用されます。
数値による確認ができたところで、次に視点を変えて、非圧縮性の取り扱いと圧縮性の正しい取り扱いの間にどの程度の差が生じるかを、定量的に比較してみましょう。
圧縮性 vs 非圧縮性の比較
理論的な差異
非圧縮性のベルヌーイの定理では、密度 $\rho$ を一定として全圧を次のように求めます。
$$ \begin{equation} p_0^{\text{incomp}} = p + \frac{1}{2}\rho V^2 \end{equation} $$
この式を静圧 $p$ で割り、$\rho V^2 / (2p)$ を $\gamma M^2 / 2$（先ほど導いた関係）で置き換えると、非圧縮の全圧比は次のようになります。
$$ \begin{equation} \frac{p_0^{\text{incomp}}}{p} = 1 + \frac{\gamma}{2} M^2 \end{equation} $$
一方、等エントロピーの正しい全圧比は次の通りです。
$$ \frac{p_0}{p} = \left(1 + \frac{\gamma – 1}{2} M^2\right)^{\gamma/(\gamma – 1)} $$
両者の差を「圧縮性補正」と呼びます。この差がどの程度あるかをPythonで可視化してみましょう。
Pythonで圧縮性と非圧縮性の比較を可視化する
非圧縮性のベルヌーイの定理による全圧比と、等エントロピー関係式による正しい全圧比を比較し、マッハ数が上がるにつれて両者の乖離がどれだけ大きくなるかを確認します。
import numpy as np
import matplotlib.pyplot as plt
gamma = 1.4
M = np.linspace(0.01, 2.0, 500)
# 等エントロピー（圧縮性あり）
p0_comp = (1 + (gamma - 1) / 2 * M**2) ** (gamma / (gamma - 1))
# ベルヌーイ（非圧縮性）
p0_incomp = 1 + gamma / 2 * M**2
# 相対誤差 (%)
relative_error = (p0_incomp - p0_comp) / p0_comp * 100
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))
# 左: 全圧比の比較
ax1.plot(M, p0_comp, label="Compressible (isentropic)", linewidth=2, color="#F44336")
ax1.plot(M, p0_incomp, label="Incompressible (Bernoulli)", linewidth=2,
color="#2196F3", linestyle="--")
ax1.axvline(x=0.3, color="orange", linestyle=":", linewidth=0.8, alpha=0.7)
ax1.text(0.32, 1.1, "M = 0.3", fontsize=9, color="orange")
ax1.set_xlabel("Mach number $M$", fontsize=12)
ax1.set_ylabel("$p_0 / p$", fontsize=13)
ax1.set_title("Total Pressure Ratio: Compressible vs Incompressible", fontsize=12)
ax1.legend(fontsize=11)
ax1.set_xlim(0, 2)
ax1.grid(True, alpha=0.3)
# 右: 相対誤差
ax2.plot(M, relative_error, linewidth=2, color="#9C27B0")
ax2.axhline(y=1, color="gray", linestyle="--", linewidth=0.8, alpha=0.5)
ax2.axhline(y=5, color="gray", linestyle="--", linewidth=0.8, alpha=0.5)
ax2.axvline(x=0.3, color="orange", linestyle=":", linewidth=0.8, alpha=0.7)
ax2.text(0.32, 6, "M = 0.3", fontsize=9, color="orange")
ax2.text(0.02, 1.5, "1% error", fontsize=9, color="gray")
ax2.text(0.02, 5.5, "5% error", fontsize=9, color="gray")
ax2.set_xlabel("Mach number $M$", fontsize=12)
ax2.set_ylabel("Relative error (%)", fontsize=13)
ax2.set_title("Error of Incompressible Approximation", fontsize=12)
ax2.set_xlim(0, 2)
ax2.set_ylim(0, 30)
ax2.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
このグラフから、以下の重要な結論が読み取れます。
- $M = 0.3$ での相対誤差は約0.13% — 左のグラフで2本の曲線はほぼ重なっています。右のグラフでも誤差は1%のラインをはるかに下回っています。これにより、$M < 0.3$ では非圧縮ベルヌーイの式で十分な精度が得られることが定量的に確認できました。
- $M \approx 0.6$ で誤差が1%を超える — 工学的に1%以上の誤差が問題になる場合は、$M = 0.6$ あたりから圧縮性補正が必要です。
- $M = 1$ での誤差は約5.3% — 音速近傍では非圧縮近似は完全に破綻し始めています。$M = 2$ では誤差が約27%にまで達し、非圧縮のベルヌーイの式は全く使えないことがわかります。
- 非圧縮近似は常に全圧を過小評価する — 非圧縮の式（破線）は常に圧縮性の正しい値（実線）よりも小さい値を与えます。これは、圧縮性効果により流体が圧縮されてさらにエネルギーが蓄えられる効果を、非圧縮近似が見落としているためです。
ここまでの議論で、等エントロピー関係式を軸にした圧縮性流体の基本的な振る舞いを理論とPythonで確認してきました。次に、面積-速度関係をさらに掘り下げ、収束-発散ノズル内のマッハ数分布を計算・可視化してみましょう。
収束-発散ノズルの流れ
ノズル形状とマッハ数分布
ラバルノズルは、収束部（断面積が減少する部分）と発散部（断面積が増加する部分）を持つノズルです。面積-速度関係で見たように、亜音速流れを超音速まで加速するにはこの収束-発散形状が不可欠です。
ノズル内のマッハ数分布は、面積比 $A/A^*$ から逆算できます。面積-マッハ数関係式は陽に $M$ について解くことができないため、各位置の $A/A^*$ に対して数値的に $M$ を求めます。
ここで重要なのは、ノズルの出口圧力（背圧）によって流れのパターンが変わるという点です。代表的なパターンは以下の通りです。
- 完全に亜音速の流れ: 背圧が十分高い場合、流れはノズル全体で亜音速のままです。収束部で加速し、スロートで最大速度（$M < 1$）に達し、発散部で減速します。
- スロートでちょうどチョーク: 背圧がちょうど臨界値の場合、スロートで $M = 1$ に達しますが、発散部では再び亜音速に戻ります（亜音速解）。
- 超音速流れ（設計条件）: 背圧が十分低い場合、スロートで $M = 1$ に達した後、発散部で超音速に加速し続けます（超音速解）。
- 衝撃波を含む流れ: 背圧が2と3の間の場合、発散部の途中に衝撃波が発生し、流れが不連続に亜音速に戻ります。
以下のPythonコードでは、設計条件（完全超音速流れ）と亜音速解のマッハ数分布を、与えられたノズル形状に対して数値的に計算し、可視化します。
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import brentq
gamma = 1.4
def area_mach_relation(M, gamma=1.4):
"""面積比 A/A* をマッハ数から計算"""
exp = (gamma + 1) / (2 * (gamma - 1))
return (1 / M) * ((2 / (gamma + 1)) * (1 + (gamma - 1) / 2 * M**2)) ** exp
def mach_from_area_ratio(A_ratio, gamma=1.4, subsonic=True):
"""面積比 A/A* から対応するマッハ数を数値的に求める"""
def equation(M):
return area_mach_relation(M, gamma) - A_ratio
if subsonic:
return brentq(equation, 1e-6, 1.0)
else:
return brentq(equation, 1.0, 30.0)
# ノズル形状の定義（放物線型）
# x: 0 → 1 (入口 → スロート → 出口), スロートは x = 0.4
x = np.linspace(0, 1, 500)
x_throat = 0.4
# 面積分布: スロートで A/A* = 1、入口・出口で拡大
# 入口面積比 = 3, 出口面積比 = 2.5
A_inlet = 3.0
A_exit = 2.5
A_throat = 1.0
A_over_Astar = np.where(
x <= x_throat,
A_throat + (A_inlet - A_throat) * ((x_throat - x) / x_throat) ** 2,
A_throat + (A_exit - A_throat) * ((x - x_throat) / (1 - x_throat)) ** 2,
)
# 各位置でのマッハ数を計算（亜音速解と超音速解の両方）
M_subsonic = np.zeros_like(x)
M_supersonic = np.zeros_like(x)
for i, A_r in enumerate(A_over_Astar):
if A_r < 1.0001: # スロート付近
M_subsonic[i] = 1.0
M_supersonic[i] = 1.0
else:
M_subsonic[i] = mach_from_area_ratio(A_r, subsonic=True)
M_supersonic[i] = mach_from_area_ratio(A_r, subsonic=False)
# 設計条件: 収束部は亜音速、スロート以降は超音速
M_design = np.where(x <= x_throat, M_subsonic, M_supersonic)
# 亜音速のみ（発散部も亜音速）
M_allsub = M_subsonic.copy()
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(10, 8), sharex=True)
# 上段: ノズル形状
ax1.fill_between(x, -A_over_Astar / 2, A_over_Astar / 2,
alpha=0.15, color="#607D8B")
ax1.plot(x, A_over_Astar / 2, color="#607D8B", linewidth=1.5)
ax1.plot(x, -A_over_Astar / 2, color="#607D8B", linewidth=1.5)
ax1.axvline(x=x_throat, color="gray", linestyle="--", linewidth=0.8, alpha=0.5)
ax1.text(x_throat + 0.01, 1.3, "Throat", fontsize=10, color="gray")
ax1.set_ylabel("Nozzle shape $A/A^*$", fontsize=12)
ax1.set_title("Converging-Diverging Nozzle: Shape and Mach Distribution", fontsize=13)
ax1.set_ylim(-2, 2)
ax1.grid(True, alpha=0.3)
# 下段: マッハ数分布
ax2.plot(x, M_design, label="Design condition (supersonic exit)",
linewidth=2.5, color="#F44336")
ax2.plot(x, M_allsub, label="Fully subsonic",
linewidth=2.5, color="#2196F3", linestyle="--")
ax2.axhline(y=1.0, color="gray", linestyle=":", linewidth=0.8, alpha=0.5)
ax2.text(0.02, 1.05, "M = 1", fontsize=10, color="gray")
ax2.axvline(x=x_throat, color="gray", linestyle="--", linewidth=0.8, alpha=0.5)
ax2.set_xlabel("Axial position $x / L$", fontsize=12)
ax2.set_ylabel("Mach number $M$", fontsize=13)
ax2.legend(fontsize=11)
ax2.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
このグラフから、収束-発散ノズルの2つの典型的な流れパターンを視覚的に理解できます。
- 設計条件（赤の実線）: 収束部ではマッハ数が亜音速で増加し、スロートで $M = 1$ に達した後、発散部で超音速に加速し続けています。出口マッハ数は約2.4に達しています。これがロケットエンジンなどで実現される理想的な超音速ノズル流れです。
- 完全亜音速（青の破線）: 収束部ではマッハ数が増加し、スロートで最大値（$M = 1$）に達しますが、発散部では減速して再び亜音速に戻っています。これは背圧が高い場合に実現される流れパターンであり、発散部はディフューザ（減速・昇圧装置）として機能しています。
- スロートの特殊性: 両方のパターンでスロート位置（$x/L = 0.4$）で $M = 1$ が実現しています。面積-速度関係で導いた「$M = 1$ は $dA = 0$ の位置でしか実現しない」という理論的結論がグラフ上で確認できます。
質量流量とチョーク現象
質量流量の式
ノズルの設計で最も重要な量の一つが質量流量 $\dot{m}$ です。スロートでの質量流量を等エントロピー関係式で表してみましょう。
スロートの断面積を $A^*$、よどみ点での温度と圧力を $T_0$, $p_0$ とすると、スロートで $M = 1$ のときの質量流量は次のようになります。
$$ \dot{m} = \rho^* a^* A^* $$
等エントロピー関係式と理想気体の状態方程式を使って $\rho^*$ と $a^*$ を $p_0$, $T_0$ で表すと、次の式が得られます。
$$ \begin{equation} \dot{m}_{\max} = A^* p_0 \sqrt{\frac{\gamma}{R T_0}} \left(\frac{2}{\gamma + 1}\right)^{(\gamma + 1) / [2(\gamma – 1)]} \end{equation} $$
この式から読み取れる物理を整理しましょう。
質量流量は $p_0$ に比例し、$\sqrt{T_0}$ に反比例します。よどみ点圧力を2倍にすると質量流量は2倍になりますが、よどみ点温度を2倍（絶対温度で）にすると質量流量は $1/\sqrt{2} \approx 0.71$ 倍に減少します。温度が高い気体は密度が低い（同じ圧力でも体積あたりの質量が少ない）ため、質量流量が減るのです。
チョーク現象
スロートでマッハ数が1に達すると、これ以上背圧を下げても質量流量は増加しません。この状態をチョーク（choking）と呼びます。チョーク状態ではスロートがボトルネックとなり、流量はスロート面積 $A^*$ とよどみ点条件 $p_0$, $T_0$ だけで決まります。
日常的な例で言えば、高圧ガスボンベのバルブを開けるとき、バルブの開度を大きくしても（下流の圧力を下げても）、ある時点から出てくるガスの量が増えなくなる現象がチョークです。バルブの最小断面積（スロート）で流れが音速に達しているため、「もっとガスを出せ」という信号（圧力波）が上流のボンベに伝わらないのです。
チョーク現象は安全工学の観点からも重要です。圧力容器の安全弁はチョーク流量を基準に設計されており、上述の質量流量式はそのための基本ツールです。
ここまでで、等エントロピー流れの基本的な理論体系を一通り見てきました。最後に、マッハ数を変化させたときの状態量（温度・圧力・密度）の変化を三次元的に可視化して、等エントロピー流れの全体像を俯瞰してみましょう。
Pythonで比熱比の影響を可視化する
等エントロピー関係式は比熱比 $\gamma$ に依存します。空気（$\gamma = 1.4$）、単原子分子ガス（ヘリウムなど、$\gamma = 5/3$）、ロケット燃焼ガスの典型値（$\gamma = 1.2$）について、圧力比のマッハ数依存性を比較してみましょう。異なる気体で等エントロピー関係式がどの程度変わるかを確認します。
import numpy as np
import matplotlib.pyplot as plt
M = np.linspace(0.01, 4.0, 500)
gammas = [1.2, 1.4, 5/3]
labels = [r"$\gamma = 1.2$ (rocket exhaust)",
r"$\gamma = 1.4$ (air)",
r"$\gamma = 5/3$ (monatomic, He)"]
colors = ["#FF9800", "#F44336", "#2196F3"]
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))
for gamma, label, color in zip(gammas, labels, colors):
p_ratio = (1 + (gamma - 1) / 2 * M**2) ** (-gamma / (gamma - 1))
T_ratio = (1 + (gamma - 1) / 2 * M**2) ** (-1)
ax1.plot(M, p_ratio, label=label, linewidth=2, color=color)
ax2.plot(M, T_ratio, label=label, linewidth=2, color=color)
ax1.set_xlabel("Mach number $M$", fontsize=12)
ax1.set_ylabel("$p / p_0$", fontsize=13)
ax1.set_title("Pressure Ratio for Different $\\gamma$", fontsize=13)
ax1.legend(fontsize=10)
ax1.set_xlim(0, 4)
ax1.set_ylim(0, 1.05)
ax1.grid(True, alpha=0.3)
ax2.set_xlabel("Mach number $M$", fontsize=12)
ax2.set_ylabel("$T / T_0$", fontsize=13)
ax2.set_title("Temperature Ratio for Different $\\gamma$", fontsize=13)
ax2.legend(fontsize=10)
ax2.set_xlim(0, 4)
ax2.set_ylim(0, 1.05)
ax2.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
このグラフから、比熱比 $\gamma$ の違いが等エントロピー関係式に与える影響を読み取ることができます。
- 圧力比（左図）では $\gamma$ の影響が顕著 — $M = 3$ のとき、$\gamma = 1.2$ では $p/p_0 \approx 0.052$、$\gamma = 5/3$ では $p/p_0 \approx 0.013$ と、約4倍の差があります。$\gamma$ が大きいほど圧力の低下が急激です。これは指数 $\gamma/(\gamma-1)$ が $\gamma$ の増大とともに大きくなるためです。
- 温度比（右図）では $\gamma$ の影響は相対的に小さい — 指数が $-1$ で固定されているため、$\gamma$ による違いは分母の $(\gamma-1)/2$ の項の違いのみです。$\gamma = 5/3$（単原子分子）は $(\gamma-1)/2 = 1/3$ と大きいため、温度低下がやや急になります。$\gamma = 1.2$ の場合は $(\gamma-1)/2 = 0.1$ と小さく、温度低下が緩やかです。
- ロケット推進工学における $\gamma$ の重要性 — ロケットエンジンの燃焼ガスは $\gamma \approx 1.2$ であり、空気（$\gamma = 1.4$）と比べて同じマッハ数での圧力低下が穏やかです。これはノズル出口の圧力が相対的に高いことを意味し、推力や比推力の計算に大きな影響を与えます。ラバルノズルの設計において $\gamma$ の正確な値が重要とされる理由がここにあります。
エネルギー方程式の詳細
全温の保存
等エントロピー流れだけでなく、より一般的に断熱流れ（粘性散逸や衝撃波を含む場合でも）では、全エンタルピー $h_0 = h + V^2/2 = C_p T_0$ が保存されます。つまり、全温 $T_0$ は衝撃波を通過しても変化しません。
これは直感的には次のように理解できます。衝撃波は流体を減速させ（運動エネルギーが減少）、圧縮して温度を上げます（内部エネルギーが増加）。衝撃波の前後で熱のやりとりはないので、運動エネルギーの減少分がそのまま内部エネルギー（温度）の増加に使われ、合計（全エンタルピー）は変わらないのです。
ただし、全圧は衝撃波を通過すると減少します。衝撃波はエントロピーを増大させる不可逆過程であり、全圧の減少はこのエントロピー増大に対応しています。全圧の損失は衝撃波の強さ（マッハ数）が大きいほど大きく、これが超音速飛行における波動抗力の源泉の一つです。
速度と温度のトレードオフ
エネルギー方程式 $C_p T + V^2/2 = C_p T_0$ は、流体が持つ「温度の形のエネルギー」と「運動の形のエネルギー」がトレードオフの関係にあることを示しています。
流体が加速する（$V$ が増える）と温度 $T$ が下がり、減速すると温度が上がります。超音速風洞では、高温高圧のガスをノズルで膨張・加速させることで、テストセクションに高マッハ数の低温流れを作り出します。たとえば $T_0 = 500$ K のガスを $M = 3$ まで加速すると、$T = T_0/(1 + 0.2 \times 9) = 500/2.8 = 178.6$ K（$-94.5$°C）にまで冷却されます。
この温度低下は凝縮の問題を引き起こすことがあります。空気中の水蒸気がテストセクションの低温で凝縮し、流れの特性を変えてしまうのです。これを防ぐため、超音速風洞では除湿が重要な課題になります。
このように、圧縮性流体力学では温度・圧力・密度・速度が互いに密接に結びついており、一つの量の変化が他の全ての量に波及します。ナビエ・ストークス方程式の圧縮性バージョンではこれらの全てを連立して解く必要があり、非圧縮性の場合よりもはるかに複雑な問題になります。
マッハ円錐とマッハ角
マッハ円錐の形成
超音速で飛行する物体は、移動するにつれて次々と球面状の圧力擾乱を放出しますが、物体の速度が音速を超えているため、これらの擾乱は物体の後方にしか伝わりません。こうしてできる擾乱の包絡面がマッハ円錐（Mach cone）です。二次元ではこれをマッハ波と呼びます。
マッハ円錐の半頂角 $\mu$（マッハ角）は、幾何学的に次のように求まります。物体が時間 $\Delta t$ の間に $V \Delta t$ だけ進む一方、その間に発せられた音波は半径 $a \Delta t$ の球面まで広がります。物体が音速より速く移動しているため、音波の球面は物体の後方に取り残されます。マッハ円錐の半頂角は次の通りです。
$$ \begin{equation} \sin \mu = \frac{a \Delta t}{V \Delta t} = \frac{a}{V} = \frac{1}{M} \end{equation} $$
すなわち、マッハ角はマッハ数の逆正弦です。$M = 1$ ではマッハ角が90°（垂直な平面波）、$M = 2$ では30°、$M$ が大きくなるほどマッハ円錐は細くなります。
超音速機が発生するソニックブーム（衝撃音）は、このマッハ円錐が地上に到達する現象です。マッハ数が大きいほどマッハ円錐が鋭くなり、ソニックブームの到達範囲（地上に影響する帯の幅）が変化します。
マッハ角の意味
マッハ角 $\mu$ は超音速流れの解析で基本的な役割を果たします。超音速流れ中に置かれた微小な擾乱源（たとえば壁面の小さな突起）は、マッハ角の方向にマッハ波を放射します。このマッハ波は微弱な圧縮波であり、流れの方向変化や圧力変化の情報を伝えます。
有限の大きさの物体（翼やくさび）がある場合は、マッハ波が集積して衝撃波（斜め衝撃波）を形成します。斜め衝撃波の角度はマッハ角よりも大きく、流れの転向角とマッハ数から決定されます。衝撃波と膨張波の詳細な理論については、発展記事で扱います。
圧縮性流体力学の体系的な位置づけ
ここまでの内容を振り返り、圧縮性流体力学の体系の中での各概念の位置づけを整理しましょう。
基本パラメータ: マッハ数 — 流速と音速の比であるマッハ数が、流れの圧縮性の度合いを支配する最も重要な無次元数です。マッハ数は流れの分類（亜音速・遷音速・超音速・極超音速）の基準であり、等エントロピー関係式を通じて全ての状態量を結びつけます。
基本仮定: 等エントロピー — 断熱かつ可逆という理想化のもとで、マッハ数だけから温度比・圧力比・密度比・面積比が求まります。この等エントロピーの仮定は、衝撃波のない滑らかな流れの領域で成り立ち、ノズル・ディフューザの設計における出発点となります。
キーとなる関係式: – 全温比 $T_0/T = 1 + (\gamma-1)M^2/2$ — エネルギー方程式から直接得られる – 全圧比・全密度比 — 等エントロピー条件を組み合わせて得られる – 面積-マッハ数関係 — 連続の式と運動量方程式から導かれる
衝撃波が登場する場面 — 等エントロピーの仮定が破れるのは、衝撃波が存在する場合です。衝撃波を横切ると圧力・温度・密度が不連続に変化し、エントロピーが増大します。衝撃波の前後関係を記述するのがランキン-ユゴニオの関係式であり、これは等エントロピー関係式と並んで圧縮性流体力学の二大ツールです。
まとめ
本記事では、圧縮性流体力学の基礎を音速とマッハ数を軸に一貫して解説しました。
- 圧縮性の判定基準: マッハ数 $M = 0.3$ を境に、非圧縮性近似の誤差が工学的に無視できなくなります。Pythonによる定量比較で、この閾値の妥当性を確認しました
- 音速の物理: 音速 $a = \sqrt{\gamma RT}$ は微小擾乱の伝播速度であり、理想気体では温度のみに依存します。音速は「情報伝達速度」として、亜音速と超音速の本質的な違いを生み出す根源です
- マッハ数の分類: 亜音速（$M < 0.8$）、遷音速（$0.8 \leq M < 1.2$）、超音速（$1.2 \leq M < 5$）、極超音速（$M \geq 5$）の4つに分類され、それぞれ異なる物理現象が支配的になります
- 等エントロピー関係式: 全温比 $T_0/T$、全圧比 $p_0/p$、全密度比 $\rho_0/\rho$ の3つの基本関係式が、マッハ数と比熱比だけから全ての状態量を結びつけます
- よどみ点の性質: 全温 $T_0$ は断熱流れで保存され（衝撃波を含む場合でも）、全圧 $p_0$ は等エントロピー流れでのみ保存されます
- 臨界条件とチョーク: $M = 1$ はスロートでしか実現できず、チョーク状態では質量流量が最大値に達します。臨界圧力比 $p^*/p_0 \approx 0.528$（空気）は安全弁やノズルの設計に直結する重要な数値です
- 面積-速度関係: 亜音速では断面積の縮小が加速をもたらし、超音速ではその逆になるという、圧縮性流体に特有の非直感的な関係が成り立ちます。これが収束-発散ノズルの動作原理です
これらの概念は、航空宇宙工学・ガスタービン工学・推進工学など、高速気流を扱うあらゆる分野の基盤です。
次のステップとして、以下の記事も参考にしてください。
- 衝撃波と膨張波 — 等エントロピーの仮定が破れる場面。垂直衝撃波と斜め衝撃波の理論
- ラバルノズル — 収束-発散ノズルの設計理論と推力の計算
- ナビエ・ストークス方程式 — 粘性を含む圧縮性流れの支配方程式