境界層理論(プラントル)を基礎から徹底解説

飛行機の翼はなぜ揚力を生むのでしょうか。自動車のボディはなぜ空気抵抗を受けるのでしょうか。これらの問いに答えるには、物体表面のごく薄い領域で起きている現象を理解しなければなりません。その薄い領域こそが境界層(boundary layer)です。

1904年、ルートヴィヒ・プラントル(Ludwig Prandtl)は「高レイノルズ数の流れでは、粘性の影響は物体表面近傍の薄い層に集中する」という画期的な着想を発表しました。この発見により、それまで理論(非粘性の完全流体)と実験(粘性を含む実際の流れ)の間にあった深刻な乖離 — いわゆるダランベールのパラドックス(非粘性理論では物体に抗力が生じない)— が解消され、流体力学は飛躍的に発展しました。

境界層理論は現在でも、航空機の翼型設計、自動車の空力最適化、タービン翼の冷却設計、気象学における大気境界層の解析など、工学と科学の幅広い分野で基礎となっています。

本記事の内容

  • 境界層の概念とプラントルの着想
  • 境界層方程式の導出(NS方程式からのオーダー解析)
  • ブラジウスの相似解(平板上の層流境界層)
  • 境界層厚さの3つの定義($\delta$, $\delta^*$, $\theta$)
  • 運動量積分法(カルマンの積分方程式)
  • 層流から乱流への遷移
  • 乱流境界層の1/7乗則
  • 摩擦抗力の計算
  • Pythonによる数値計算と可視化

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

境界層の概念

プラントルの着想

一様流の中に置かれた物体を考えましょう。レイノルズ数 $Re = \rho U L / \mu$ が非常に大きい場合、流れの大部分では粘性力は慣性力に比べて無視できるほど小さくなります。しかし物体表面では滑りなし条件(no-slip condition)が成り立ち、流速はゼロです。一方、表面から十分離れた場所では流速はほぼ主流速度 $U$ に等しくなります。

つまり、物体表面近傍には流速が $0$ から $U$ へと急激に変化する薄い領域が存在します。この領域が境界層です。プラントルの核心的な洞察は次の点にあります。

  • 境界層の外側: 粘性を無視できる。ベルヌーイの定理やポテンシャル流れの理論が適用できる
  • 境界層の内側: 粘性力と慣性力が同程度のオーダーを持つ。完全なナビエ・ストークス方程式を解く必要があるが、層の薄さを利用した簡略化(境界層近似)が可能

この「内と外に分けて考える」手法は、流体力学にとどまらず、特異摂動法やマッチドアシンプトティクス展開として数学的にも体系化されています。

境界層の厚さのスケール

境界層がどの程度薄いかを大まかに見積もりましょう。粘性の影響が及ぶ範囲は、粘性拡散の時間スケールから推定できます。流体粒子が物体の先端から距離 $x$ を移動する時間は $t \sim x/U$ です。この時間で粘性が拡散する距離は:

$$ \delta \sim \sqrt{\nu t} \sim \sqrt{\frac{\nu x}{U}} $$

ここで $\nu = \mu / \rho$ は動粘度です。これを $x$ で無次元化すると:

$$ \frac{\delta}{x} \sim \frac{1}{\sqrt{Re_x}} $$

$Re_x = Ux/\nu$ は位置 $x$ に基づくレイノルズ数です。$Re_x = 10^6$ なら $\delta/x \sim 10^{-3}$、つまり境界層の厚さは代表長さの1000分の1程度です。この「薄さ」が境界層近似の正当性を保証します。

ここまでで境界層の物理的描像が掴めました。次に、ナビエ・ストークス方程式から出発して、境界層方程式を数学的に導出しましょう。

境界層方程式の導出

出発点: 2次元非圧縮NS方程式

2次元の定常非圧縮性流れを考えます。$x$ 方向を主流方向、$y$ 方向を壁面法線方向とします。連続の式とNS方程式の $x$ 成分、$y$ 成分は:

連続の式:

$$ \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} = 0 $$

$x$ 方向の運動量方程式:

$$ u\frac{\partial u}{\partial x} + v\frac{\partial u}{\partial y} = -\frac{1}{\rho}\frac{\partial p}{\partial x} + \nu\left(\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2}\right) $$

$y$ 方向の運動量方程式:

$$ u\frac{\partial v}{\partial x} + v\frac{\partial v}{\partial y} = -\frac{1}{\rho}\frac{\partial p}{\partial y} + \nu\left(\frac{\partial^2 v}{\partial x^2} + \frac{\partial^2 v}{\partial y^2}\right) $$

オーダー解析(無次元化)

境界層方程式を導くために、各量のオーダーを見積もります。特性量として主流速度 $U$、代表長さ $L$(物体長さ)、境界層厚さ $\delta$ を使います。

各変数のオーダーを次のように置きます:

$$ u \sim U, \quad x \sim L, \quad y \sim \delta $$

連続の式 $\partial u / \partial x + \partial v / \partial y = 0$ から、各項は同じオーダーでなければなりません。$\partial u / \partial x \sim U/L$ なので:

$$ \frac{\partial v}{\partial y} \sim \frac{U}{L} \quad \Longrightarrow \quad v \sim \frac{U\delta}{L} $$

つまり、法線方向の速度 $v$ は主流速度の $\delta/L$ 倍程度と非常に小さいことがわかります。

次に $x$ 方向の運動量方程式の各項のオーダーを評価します:

オーダー
$u \dfrac{\partial u}{\partial x}$ $\dfrac{U^2}{L}$
$v \dfrac{\partial u}{\partial y}$ $\dfrac{U\delta}{L} \cdot \dfrac{U}{\delta} = \dfrac{U^2}{L}$
$\nu \dfrac{\partial^2 u}{\partial x^2}$ $\dfrac{\nu U}{L^2}$
$\nu \dfrac{\partial^2 u}{\partial y^2}$ $\dfrac{\nu U}{\delta^2}$

左辺の慣性項は $O(U^2/L)$ です。粘性項のうち $\partial^2 u / \partial x^2$ は $O(\nu U / L^2)$ であり、慣性項に比べて $O(1/Re_L)$ だけ小さいので無視できます。一方、$\partial^2 u / \partial y^2$ は $O(\nu U / \delta^2)$ です。

境界層内では粘性力と慣性力が釣り合わなければなりません(これこそが境界層の存在意義です)。したがって:

$$ \frac{U^2}{L} \sim \frac{\nu U}{\delta^2} \quad \Longrightarrow \quad \delta \sim \sqrt{\frac{\nu L}{U}} = \frac{L}{\sqrt{Re_L}} $$

これは先ほどのスケール推定と一致します。

$y$ 方向の運動量方程式の簡略化

$y$ 方向の方程式についても同様のオーダー解析を行うと、全ての項が $x$ 方向の方程式に比べて $O(\delta/L)$ だけ小さいことがわかります。結果として:

$$ \frac{\partial p}{\partial y} \approx 0 $$

つまり、境界層内の圧力は $y$ 方向に一定であり、境界層の外縁で決まる圧力がそのまま壁面まで伝わります。圧力は $x$ のみの関数 $p = p_e(x)$ で、これは外部流れ(非粘性理論)から求めることができます。

外部流れでは粘性を無視できるので、ベルヌーイの定理から:

$$ -\frac{1}{\rho}\frac{dp_e}{dx} = U_e\frac{dU_e}{dx} $$

ここで $U_e(x)$ は境界層外縁の速度です。

プラントルの境界層方程式

以上のオーダー解析により、NS方程式から不要な項を落として得られる方程式系がプラントルの境界層方程式です:

$$ \boxed{u\frac{\partial u}{\partial x} + v\frac{\partial u}{\partial y} = U_e\frac{dU_e}{dx} + \nu\frac{\partial^2 u}{\partial y^2}} $$

$$ \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} = 0 $$

境界条件:

  • 壁面($y = 0$): $u = 0, \; v = 0$(滑りなし条件)
  • 境界層外縁($y \to \infty$): $u \to U_e(x)$

元のNS方程式と比較して、次の簡略化が行われています:

  1. $x$ 方向の粘性拡散項 $\nu \partial^2 u / \partial x^2$ を無視
  2. $y$ 方向の運動量方程式を $\partial p / \partial y = 0$ に置き換え
  3. 圧力項を外部流れの速度 $U_e(x)$ で表現

この方程式は元のNS方程式と異なり、$x$ 方向に放物型であるため、上流から下流へ順番に解いていくことができます。これは数値計算上も大きな利点です。

境界層方程式を導出できました。まずは最も基本的なケースとして、圧力勾配のない平板上の流れ($U_e = U = $ 一定)を考え、ブラジウスの相似解を求めましょう。

ブラジウスの相似解

平板上の流れ

無限に続く平板の前縁に一様流速 $U$ の流れが当たる場合を考えます。外部流速が一定($U_e = U$, $dU_e/dx = 0$)なので、境界層方程式は:

$$ u\frac{\partial u}{\partial x} + v\frac{\partial u}{\partial y} = \nu\frac{\partial^2 u}{\partial y^2} $$

$$ \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} = 0 $$

ブラジウス(Blasius, 1908)は、この問題に相似解が存在することを見出しました。相似解とは、$x$ のどの位置でも速度分布の「形」が同じで、$y$ のスケールだけが変わるという解です。

相似変数の導入

無次元の相似変数 $\eta$ と流れ関数 $\psi$ を次のように定義します:

$$ \eta = y\sqrt{\frac{U}{\nu x}} $$

$$ \psi = \sqrt{\nu x U} \, f(\eta) $$

ここで $f(\eta)$ は未知関数です。流れ関数 $\psi$ は $u = \partial\psi/\partial y$, $v = -\partial\psi/\partial x$ と定義されるので、連続の式は自動的に満たされます。

速度成分を計算します。まず $u$ について:

$$ u = \frac{\partial \psi}{\partial y} = \sqrt{\nu x U}\, f'(\eta) \cdot \frac{\partial \eta}{\partial y} = \sqrt{\nu x U}\, f'(\eta) \cdot \sqrt{\frac{U}{\nu x}} = Uf'(\eta) $$

$v$ についても計算すると:

$$ v = -\frac{\partial \psi}{\partial x} = \frac{1}{2}\sqrt{\frac{\nu U}{x}}\left(\eta f'(\eta) – f(\eta)\right) $$

ブラジウス方程式

$u, v$ を境界層方程式に代入すると、偏微分方程式が $\eta$ のみの常微分方程式に変換されます。代入の計算は煩雑ですが、結果は驚くほど簡潔です:

$$ \boxed{2f”’ + f f” = 0} $$

これがブラジウス方程式です。3階の非線形常微分方程式です。

境界条件は:

  • $\eta = 0$: $f(0) = 0$(壁面で $v = 0$)
  • $\eta = 0$: $f'(0) = 0$(壁面で $u = 0$)
  • $\eta \to \infty$: $f'(\infty) = 1$(外部流速 $u \to U$)

数値解法

ブラジウス方程式は解析的な閉じた形の解を持ちません。そこで射撃法(shooting method)を使って数値的に解きます。

射撃法のアイデアは次の通りです。初期値問題として $f(0) = 0$, $f'(0) = 0$, $f”(0) = \alpha$ を仮定し、$\eta$ を増加させて積分します。$\eta \to \infty$ で $f'(\eta) \to 1$ となるように $\alpha$ を調整します。

数値計算の結果、$f”(0) = 0.33206$ が得られます。この値は境界層理論で繰り返し現れる重要な定数です。

ブラジウスの相似解を数値的に求めるPythonコードを示しましょう。

Python: ブラジウス速度分布の計算

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

# ブラジウス方程式: 2f''' + f*f'' = 0
# y = [f, f', f''] とすると
# y' = [f', f'', -f*f''/2]
def blasius_ode(eta, y):
    f, fp, fpp = y
    return [fp, fpp, -0.5 * f * fpp]

# 射撃法: f''(0) = alpha を仮定し、f'(inf) = 1 を満たすalphaを求める
def shooting(alpha):
    sol = solve_ivp(blasius_ode, [0, 10], [0, 0, alpha],
                    t_eval=np.linspace(0, 10, 1000), rtol=1e-10, atol=1e-12)
    return sol.y[1][-1] - 1.0  # f'(inf) - 1 = 0 を満たすalphaを探索

# Brent法でalphaを求める
alpha_opt = brentq(shooting, 0.1, 1.0)
print(f"f''(0) = {alpha_opt:.5f}")

# 最適alphaで改めてブラジウス方程式を解く
sol = solve_ivp(blasius_ode, [0, 8], [0, 0, alpha_opt],
                t_eval=np.linspace(0, 8, 500), rtol=1e-10, atol=1e-12)

eta = sol.t
f_vals = sol.y[0]
fp_vals = sol.y[1]   # f' = u/U
fpp_vals = sol.y[2]  # f''

# 可視化
fig, axes = plt.subplots(1, 3, figsize=(16, 5))

# (1) 速度分布 u/U = f'(eta)
axes[0].plot(fp_vals, eta, 'b-', lw=2.5)
axes[0].axvline(0.99, color='r', ls='--', alpha=0.6, label='$u/U = 0.99$')
idx_99 = np.argmin(np.abs(fp_vals - 0.99))
axes[0].plot(0.99, eta[idx_99], 'ro', ms=8)
axes[0].annotate(f'$\\eta = {eta[idx_99]:.2f}$',
                 xy=(0.99, eta[idx_99]), xytext=(0.7, eta[idx_99]+0.5),
                 fontsize=11, arrowprops=dict(arrowstyle='->', color='red'))
axes[0].set_xlabel("$u/U = f'(\\eta)$", fontsize=13)
axes[0].set_ylabel('$\\eta = y\\sqrt{U/(\\nu x)}$', fontsize=13)
axes[0].set_title('Blasius velocity profile', fontsize=14)
axes[0].legend(fontsize=11); axes[0].grid(True, alpha=0.3)
axes[0].set_xlim(-0.05, 1.1); axes[0].set_ylim(0, 8)

# (2) f(eta), f'(eta), f''(eta)
axes[1].plot(eta, f_vals, 'b-', lw=2, label="$f(\\eta)$")
axes[1].plot(eta, fp_vals, 'r-', lw=2, label="$f'(\\eta) = u/U$")
axes[1].plot(eta, fpp_vals, 'g-', lw=2, label="$f''(\\eta)$")
axes[1].set_xlabel('$\\eta$', fontsize=13)
axes[1].set_ylabel('Function values', fontsize=13)
axes[1].set_title('Blasius solution components', fontsize=14)
axes[1].legend(fontsize=11); axes[1].grid(True, alpha=0.3)

# (3) 壁面せん断応力分布 (f''(0)を用いた局所摩擦係数)
x_vals = np.linspace(0.01, 2.0, 200)
U_inf = 10.0  # 主流速度 [m/s]
nu = 1.5e-5   # 空気の動粘度 [m^2/s]
Rex = U_inf * x_vals / nu
Cf_x = 0.664 / np.sqrt(Rex)

axes[2].plot(x_vals, Cf_x * 1e3, 'b-', lw=2.5)
axes[2].set_xlabel('$x$ [m]', fontsize=13)
axes[2].set_ylabel('$C_f \\times 10^3$', fontsize=13)
axes[2].set_title('Local skin friction coefficient (laminar)', fontsize=14)
axes[2].grid(True, alpha=0.3)

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

上のコードは射撃法によるブラジウス方程式の数値解を計算し、3つのグラフを生成します。

左図の速度分布 $u/U = f'(\eta)$ を見ると、壁面($\eta = 0$)でゼロから滑らかに立ち上がり、$\eta \approx 5.0$ でほぼ $u/U = 0.99$ に達しています。この $\eta = 5.0$ が、後述する境界層厚さ $\delta$ の定義に使われます。

中央図では $f$, $f’$, $f”$ の3つの関数を示しています。$f”(\eta)$ は壁面での値 $f”(0) = 0.33206$ が最大で、$\eta$ とともに単調減少します。$f”(0)$ は壁面せん断応力に直結する重要な値です。

右図の局所摩擦係数 $C_f$ は $x$ の増加とともに減少します。これは下流に行くほど境界層が厚くなり、壁面での速度勾配が緩やかになるためです。

ブラジウスの解が得られたので、次に境界層の厚さをより精密に定義しましょう。工学的に重要な3つの厚さの概念があります。

境界層厚さ

境界層の「厚さ」には、目的に応じた3つの定義があります。それぞれの物理的意味を理解することが、工学的な応用で不可欠です。

境界層厚さ $\delta$(99%厚さ)

最も直感的な定義です。流速が主流速度の99%に達する位置を境界層の端とします:

$$ u(x, \delta) = 0.99\, U_e $$

ブラジウス解から $f'(\eta) = 0.99$ となる $\eta$ を求めると $\eta \approx 5.0$ です。$\eta = y\sqrt{U/(\nu x)}$ の定義から:

$$ \boxed{\delta = \frac{5.0 x}{\sqrt{Re_x}}} $$

$\delta$ は $x$ の平方根に比例して成長します。つまり、前縁から離れるほど境界層は厚くなりますが、その成長率は次第に鈍化します。

排除厚さ $\delta^*$(displacement thickness)

境界層の存在により、実際の流量は一様流の場合よりも減少します。排除厚さは、その流量の減少分を等価的に壁面を $\delta^*$ だけ外側にずらすことで表現したものです。

イメージとしては、「境界層があることで主流が $\delta^*$ だけ外側に押し出される」ということです。数学的には:

$$ \boxed{\delta^* = \int_0^\infty \left(1 – \frac{u}{U_e}\right) dy} $$

被積分関数 $1 – u/U_e$ は「速度欠損の割合」を表しています。壁面近くでは $u$ が小さいので欠損が大きく、外部に行くほど $u \to U_e$ で欠損はゼロに近づきます。

ブラジウス解の場合、相似変数を用いて積分すると:

$$ \delta^* = \frac{1.7208\, x}{\sqrt{Re_x}} $$

運動量厚さ $\theta$(momentum thickness)

運動量厚さは、境界層による運動量の欠損を表す量です。

$$ \boxed{\theta = \int_0^\infty \frac{u}{U_e}\left(1 – \frac{u}{U_e}\right) dy} $$

被積分関数に $u/U_e$ が掛かっているのは、運動量欠損が流速 $u$ で輸送されることを反映しています。排除厚さが「質量の欠損」を測るのに対して、運動量厚さは「運動量の欠損」を測っています。

ブラジウス解の場合:

$$ \theta = \frac{0.664\, x}{\sqrt{Re_x}} $$

形状係数

排除厚さと運動量厚さの比を形状係数 $H$ と呼びます:

$$ H = \frac{\delta^*}{\theta} $$

ブラジウスの層流境界層では $H \approx 2.59$ です。形状係数は境界層の状態を特徴づける重要なパラメータで、層流では $H \approx 2.3 \sim 3.5$、乱流では $H \approx 1.3 \sim 1.5$ です。$H$ が大きいほど速度分布が壁面近くに偏っていることを意味し、剥離に近い状態では $H$ が増大します。

Python: 境界層厚さの成長

import numpy as np
import matplotlib.pyplot as plt

# パラメータ
U = 30.0        # 主流速度 [m/s]
nu = 1.5e-5     # 空気の動粘度 [m^2/s]
x = np.linspace(0.001, 1.5, 500)  # 平板上の位置 [m]

Re_x = U * x / nu

# 層流境界層の各厚さ(ブラジウス解)
delta = 5.0 * x / np.sqrt(Re_x)
delta_star = 1.7208 * x / np.sqrt(Re_x)
theta = 0.664 * x / np.sqrt(Re_x)

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

# (1) 3つの境界層厚さの比較
axes[0].plot(x * 100, delta * 1000, 'b-', lw=2.5, label='$\\delta$ (99% thickness)')
axes[0].plot(x * 100, delta_star * 1000, 'r--', lw=2.5,
             label='$\\delta^*$ (displacement)')
axes[0].plot(x * 100, theta * 1000, 'g-.', lw=2.5,
             label='$\\theta$ (momentum)')
axes[0].set_xlabel('$x$ [cm]', fontsize=13)
axes[0].set_ylabel('Thickness [mm]', fontsize=13)
axes[0].set_title('Boundary layer thickness growth (laminar)', fontsize=14)
axes[0].legend(fontsize=11)
axes[0].grid(True, alpha=0.3)

# (2) 形状係数と各厚さの比
H = delta_star / theta
axes[1].axhline(2.59, color='k', ls='--', alpha=0.5, label='Blasius: $H = 2.59$')
axes[1].plot(x * 100, H, 'b-', lw=2.5, label='$H = \\delta^* / \\theta$')
axes[1].plot(x * 100, delta / theta, 'r-', lw=2, label='$\\delta / \\theta$')
axes[1].set_xlabel('$x$ [cm]', fontsize=13)
axes[1].set_ylabel('Ratio', fontsize=13)
axes[1].set_title('Shape factor and thickness ratios', fontsize=14)
axes[1].legend(fontsize=11)
axes[1].grid(True, alpha=0.3)
axes[1].set_ylim(0, 10)

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

print(f"x = 1.0 m での各厚さ:")
print(f"  Re_x = {U * 1.0 / nu:.0f}")
print(f"  delta = {5.0 * 1.0 / np.sqrt(U * 1.0 / nu) * 1000:.2f} mm")
print(f"  delta* = {1.7208 * 1.0 / np.sqrt(U * 1.0 / nu) * 1000:.2f} mm")
print(f"  theta  = {0.664 * 1.0 / np.sqrt(U * 1.0 / nu) * 1000:.2f} mm")

左図では、3つの境界層厚さがいずれも $x$ の平方根に比例して成長する様子が見えます。$\delta$ が最も大きく、次いで $\delta^*$、$\theta$ の順で、常に $\delta > \delta^* > \theta$ という関係が成立します。$x = 1.0$ m($Re_x \approx 2 \times 10^6$)では $\delta \approx 3.5$ mm と、主流方向の長さに比べて極めて薄いことが確認できます。

右図では形状係数 $H$ が $x$ によらず一定値 2.59 をとることがわかります。これは相似解の本質的な性質で、「どの $x$ でも速度分布の形が同じ」ということの帰結です。

3つの厚さの定義を理解したところで、次に境界層方程式を直接解くのではなく、積分形式で扱う強力な手法を紹介します。

運動量積分法(カルマンの積分方程式)

積分方程式の導出

境界層方程式を $y = 0$ から $y = \delta$ まで積分する方法は、フォン・カルマン(von Karman, 1921)によって提案されました。この方法は、速度分布の詳細を知らなくても、適切な速度分布の近似形を仮定するだけで境界層の主要な量を求めることができる実用的な手法です。

境界層方程式の $x$ 方向成分を $y = 0$ から $y = \infty$ まで積分します。導出の詳細は技術的ですが、連続の式を使って $v$ を消去し、ライプニッツの積分則を適用すると、次のカルマンの運動量積分方程式が得られます:

$$ \boxed{\frac{d\theta}{dx} + (2 + H)\frac{\theta}{U_e}\frac{dU_e}{dx} = \frac{C_f}{2}} $$

ここで:

  • $\theta$: 運動量厚さ
  • $H = \delta^* / \theta$: 形状係数
  • $C_f = \tau_w / (\frac{1}{2}\rho U_e^2)$: 局所摩擦係数
  • $\tau_w = \mu(\partial u / \partial y)_{y=0}$: 壁面せん断応力

平板(圧力勾配なし)の場合

$U_e = U = $ 一定の場合、$dU_e/dx = 0$ なので運動量積分方程式は:

$$ \frac{d\theta}{dx} = \frac{C_f}{2} $$

この式は物理的に明快です。境界層の運動量厚さの成長率が、壁面での摩擦(運動量の損失率)に等しいということです。

ポールハウゼン法

運動量積分法を使うには、速度分布の近似形を仮定する必要があります。代表的なものとして、4次多項式近似(ポールハウゼン法)があります。

$\eta_\delta = y/\delta$ として:

$$ \frac{u}{U_e} = 2\eta_\delta – 2\eta_\delta^3 + \eta_\delta^4 $$

この近似分布は以下の境界条件を満たすように構成されています:

  • $y = 0$ で $u = 0$(滑りなし条件)
  • $y = \delta$ で $u = U_e$(外部流速との連続性)
  • $y = \delta$ で $\partial u / \partial y = 0$(外部流速との滑らかな接続)
  • $y = \delta$ で $\partial^2 u / \partial y^2 = 0$

この速度分布を運動量厚さと壁面せん断応力の積分に代入すると:

$$ \theta = \frac{37}{315}\delta, \quad \tau_w = \mu\frac{2U}{\delta} $$

$d\theta/dx = C_f/2$ に代入して $\delta$ について解くと:

$$ \delta = \frac{5.84\, x}{\sqrt{Re_x}} $$

ブラジウスの厳密な数値解 $\delta = 5.0\,x/\sqrt{Re_x}$ と比較すると、係数が5.84対5.0で約17%の過大評価ですが、複雑な形状や圧力勾配がある場合にも適用できるという実用上の利点があります。

運動量積分法は管内流れの入口領域の解析などにも広く応用されています。ここまでは層流境界層を扱ってきましたが、実際のエンジニアリングでは境界層が層流から乱流に遷移するケースが非常に多く現れます。次にこの遷移現象を見ていきましょう。

層流から乱流への遷移

遷移の機構

境界層は、上流では層流として発達しますが、ある位置から下流で乱流に遷移します。層流と乱流で解説したように、遷移のメカニズムは次のように段階的に進行します:

  1. 不安定性の発生: 層流境界層に微小な擾乱が加わると、特定の波数を持つ擾乱が増幅される。この不安定波はトルミーン・シュリヒティング波(Tollmien-Schlichting wave, TS波)と呼ばれる
  2. 非線形成長: 増幅されたTS波が十分に大きくなると、非線形効果が支配的になり、3次元的な構造($\Lambda$渦)が発生する
  3. 乱流スポット: 3次元構造が崩壊し、局所的に乱流状態(turbulent spot)が発生する
  4. 完全乱流: 乱流スポットが合体し、境界層全体が乱流状態となる

臨界レイノルズ数

遷移が起こる位置は、局所的なレイノルズ数 $Re_x = Ux/\nu$ で特徴づけられます:

$$ Re_{x,\text{crit}} \approx 3.2 \times 10^5 \sim 5 \times 10^5 $$

この値は一定ではなく、主流の乱れ度(turbulence intensity)、表面粗さ、圧力勾配、壁面の振動などの環境条件に大きく依存します。

条件 $Re_{x,\text{crit}}$ の目安
高乱れ度の風洞 $\sim 10^5$
自由飛行(低乱れ度) $\sim 3 \times 10^6$
典型的な工学計算 $5 \times 10^5$
順圧力勾配(加速流れ) 遷移を遅らせる
逆圧力勾配(減速流れ) 遷移を促進する

遷移領域の取り扱い

実際の工学計算では、遷移領域をある $x$ 位置で瞬間的に起こるものとして扱うことが多いです。つまり、$x < x_{\text{tr}}$ では層流の公式を、$x > x_{\text{tr}}$ では乱流の公式を適用します。遷移位置 $x_{\text{tr}}$ は $Re_{x,\text{tr}} = 5 \times 10^5$ から:

$$ x_{\text{tr}} = \frac{5 \times 10^5 \cdot \nu}{U} $$

遷移後の乱流境界層では、速度分布の形が層流とは大きく異なります。次にその代表的な近似モデルを見ましょう。

乱流境界層の1/7乗則

乱流境界層の特徴

乱流境界層は層流境界層と比較して、以下の特徴があります:

  • 厚い: 同じ $Re_x$ でも乱流境界層の方が厚い。乱流による混合が壁面からの運動量輸送を促進するため
  • 壁面せん断応力が大きい: 壁面近くの速度勾配が急峻。摩擦抗力が増大する
  • 速度分布がフラット: 中央部で速度がほぼ一定に近い。壁面近くでのみ急変する
  • 剥離しにくい: 壁面近くの運動量が大きいため、逆圧力勾配に対する耐性が高い

1/7乗速度分布

乱流境界層の速度分布の近似式として、べき乗則が広く使われます。最も基本的なものが1/7乗則です:

$$ \boxed{\frac{u}{U_e} = \left(\frac{y}{\delta}\right)^{1/7}} $$

この分布は経験的なものですが、管内流れで得られるべき乗速度分布と対応しており、$Re$ が適度な範囲($Re_x \sim 10^5 \sim 10^7$)で良好な近似を与えます。

1/7乗則の速度分布を用いて排除厚さと運動量厚さを計算すると:

$$ \delta^* = \int_0^\delta \left(1 – \left(\frac{y}{\delta}\right)^{1/7}\right) dy = \frac{\delta}{8} $$

$$ \theta = \int_0^\delta \left(\frac{y}{\delta}\right)^{1/7}\left(1 – \left(\frac{y}{\delta}\right)^{1/7}\right) dy = \frac{7}{72}\delta $$

形状係数は:

$$ H = \frac{\delta^*}{\theta} = \frac{1/8}{7/72} = \frac{9}{7} \approx 1.286 $$

層流の $H \approx 2.59$ と比較して著しく小さく、乱流境界層の速度分布が壁面近くに集中していること(フラットな分布)を反映しています。

乱流境界層厚さ

1/7乗則とカルマンの運動量積分方程式を組み合わせ、壁面摩擦係数にブラジウスの管摩擦公式の類推を用いると、乱流境界層の厚さは:

$$ \boxed{\delta = \frac{0.37\, x}{Re_x^{1/5}}} $$

層流の場合($\delta \propto x / Re_x^{1/2}$)と比較すると、乱流では $Re_x$ の指数が $-1/2$ から $-1/5$ に変わっています。これは同じ $Re_x$ でも乱流の方が厚い境界層を持つことを意味し、乱流の混合効果によるものです。

層流と乱流の境界層厚さの違いを理解したところで、工学的に最も重要な量の一つである摩擦抗力について考えましょう。

摩擦抗力の計算

局所摩擦係数

壁面せん断応力 $\tau_w$ を無次元化した局所摩擦係数 $C_f$ は、流れの状態に応じて以下の公式で与えられます。

層流境界層(ブラジウス解):

$$ \boxed{C_f = \frac{0.664}{\sqrt{Re_x}}} $$

これは $\tau_w = \mu(\partial u / \partial y)_{y=0}$ から直接導けます。ブラジウス解で $f”(0) = 0.33206$ を用いると:

$$ \tau_w = \mu U \sqrt{\frac{U}{\nu x}} f”(0) = 0.33206 \, \mu U \sqrt{\frac{U}{\nu x}} $$

これを $C_f = \tau_w / (\frac{1}{2}\rho U^2)$ で無次元化すると上式が得られます。

乱流境界層

局所摩擦係数にはいくつかの経験式があります。代表的なものを挙げます:

$$ C_f = \frac{0.0592}{Re_x^{1/5}} \quad (Re_x < 10^7) $$

$$ C_f = \frac{0.370}{(\log_{10} Re_x)^{2.584}} \quad (\text{Schlichting}) $$

平均摩擦係数(摩擦抗力係数)

平板全体の摩擦抗力を求めるには、局所摩擦係数を板の長さ $L$ にわたって積分します:

$$ C_D = \frac{1}{L}\int_0^L C_f \, dx $$

全層流の場合:

$$ C_D = \frac{1}{L}\int_0^L \frac{0.664}{\sqrt{Ux/\nu}} dx = \frac{1.328}{\sqrt{Re_L}} $$

全乱流の場合:

$$ C_D = \frac{0.074}{Re_L^{1/5}} $$

混合境界層(前縁から $x_{\text{tr}}$ まで層流、それ以降乱流)の場合、シュリヒティングの公式:

$$ \boxed{C_D = \frac{0.074}{Re_L^{1/5}} – \frac{A}{Re_L}} $$

ここで $A$ は遷移レイノルズ数に依存する定数です。$Re_{x,\text{tr}} = 5 \times 10^5$ の場合、$A = 1741$ がよく使われます。

摩擦抗力の物理的意味

摩擦抗力は、流体の性質である粘性に起因して物体表面に作用する力です。片面の摩擦抗力は:

$$ D_f = \frac{1}{2}\rho U^2 \cdot L \cdot b \cdot C_D $$

ここで $b$ は平板のスパン(幅)です。重要なのは、乱流境界層の摩擦係数は層流の数倍以上大きいということです。たとえば $Re_L = 10^6$ の場合:

  • 全層流: $C_D = 1.328 / \sqrt{10^6} = 0.00133$
  • 全乱流: $C_D = 0.074 / (10^6)^{0.2} = 0.00467$

乱流の方が約3.5倍大きい値です。航空機の翼型設計で自然層流翼(Natural Laminar Flow, NLF)が追求される理由がここにあります。

摩擦抗力の理論を実際に数値で確認してみましょう。

Python: 層流 vs 乱流の摩擦係数

import numpy as np
import matplotlib.pyplot as plt

Re_L = np.logspace(4, 9, 500)

# 平均摩擦係数(平板全長にわたる積分値)
CD_laminar = 1.328 / np.sqrt(Re_L)
CD_turbulent = 0.074 / Re_L**0.2
CD_mixed = 0.074 / Re_L**0.2 - 1741 / Re_L  # 遷移 Re = 5e5

# 局所摩擦係数
Re_x = np.logspace(4, 9, 500)
Cf_laminar = 0.664 / np.sqrt(Re_x)
Cf_turbulent = 0.0592 / Re_x**0.2

fig, axes = plt.subplots(1, 2, figsize=(14, 6))

# (1) 平均摩擦係数
axes[0].loglog(Re_L, CD_laminar, 'b-', lw=2.5, label='All laminar')
axes[0].loglog(Re_L, CD_turbulent, 'r-', lw=2.5, label='All turbulent')
axes[0].loglog(Re_L[CD_mixed > 0], CD_mixed[CD_mixed > 0], 'g--', lw=2.5,
               label='Mixed ($Re_{tr}=5 \\times 10^5$)')
axes[0].axvline(5e5, color='gray', ls=':', alpha=0.5)
axes[0].set_xlabel('$Re_L$', fontsize=13)
axes[0].set_ylabel('$C_D$ (average friction coefficient)', fontsize=13)
axes[0].set_title('Average friction coefficient on flat plate', fontsize=14)
axes[0].legend(fontsize=11)
axes[0].grid(True, which='both', alpha=0.3)
axes[0].set_xlim(1e4, 1e9)
axes[0].set_ylim(1e-4, 1e-1)

# (2) 局所摩擦係数
axes[1].loglog(Re_x, Cf_laminar, 'b-', lw=2.5, label='Laminar: $C_f = 0.664/Re_x^{0.5}$')
axes[1].loglog(Re_x, Cf_turbulent, 'r-', lw=2.5,
               label='Turbulent: $C_f = 0.0592/Re_x^{0.2}$')
axes[1].axvline(5e5, color='gray', ls=':', alpha=0.5, label='$Re_{x,crit}$')
axes[1].set_xlabel('$Re_x$', fontsize=13)
axes[1].set_ylabel('$C_f$ (local friction coefficient)', fontsize=13)
axes[1].set_title('Local skin friction coefficient', fontsize=14)
axes[1].legend(fontsize=11)
axes[1].grid(True, which='both', alpha=0.3)
axes[1].set_xlim(1e4, 1e9)
axes[1].set_ylim(1e-4, 1e-1)

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

左図の平均摩擦係数では、全乱流の $C_D$ が全層流の $C_D$ を全域にわたって大きく上回っていることがわかります。混合境界層の曲線は両者の間にあり、$Re_L$ が大きくなるほど全乱流の値に近づきます。これは、平板が長くなるほど乱流部分の占める割合が増えるためです。

右図の局所摩擦係数を見ると、層流の $C_f$ は $Re_x^{-1/2}$ で減少するのに対し、乱流の $C_f$ は $Re_x^{-1/5}$ とより緩やかに減少します。遷移点 $Re_x = 5 \times 10^5$ 付近で両者の差が特に顕著で、遷移により $C_f$ が跳ね上がる様子が読み取れます。

以上で摩擦抗力の基本理論を押さえました。最後に、これまでの結果を統合して、平板上の境界層を包括的に可視化するPythonコードを作成しましょう。

平板上の境界層の可視化

Python: 境界層の成長と遷移

ここでは、平板上の境界層が前縁から発達し、層流から乱流へ遷移する全体像を1つのコードで描画します。

import numpy as np
import matplotlib.pyplot as plt

# 物理パラメータ
U = 20.0        # 主流速度 [m/s]
nu = 1.5e-5     # 空気の動粘度 [m^2/s]
L = 2.0         # 平板の長さ [m]

# 遷移位置
Re_tr = 5e5
x_tr = Re_tr * nu / U
print(f"遷移位置: x_tr = {x_tr:.3f} m")

# x座標
x_lam = np.linspace(0.001, x_tr, 300)
x_turb = np.linspace(x_tr, L, 300)

# 層流領域の境界層厚さ(ブラジウス解)
Re_x_lam = U * x_lam / nu
delta_lam = 5.0 * x_lam / np.sqrt(Re_x_lam)

# 乱流領域の境界層厚さ
# 遷移点での仮想原点を考慮: 等価乱流長さ x_eff
delta_tr = 5.0 * x_tr / np.sqrt(Re_tr)  # 遷移点での層流境界層厚さ
# 乱流境界層: delta = 0.37 * x_eff / Re_x_eff^(1/5)
# 遷移点で連続になるよう仮想原点 x0 を求める
# delta_tr = 0.37 * x0 / (U*x0/nu)^(1/5) => x0 を数値的に求める
from scipy.optimize import brentq

def find_x0(x0):
    return 0.37 * x0 / (U * x0 / nu)**0.2 - delta_tr

x0 = brentq(find_x0, 0.001, 2.0)
x_eff = x_turb - x_tr + x0
Re_x_eff = U * x_eff / nu
delta_turb = 0.37 * x_eff / Re_x_eff**0.2

fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# (1) 境界層の全体像
ax = axes[0, 0]
ax.fill_between(x_lam * 100, 0, delta_lam * 1000, alpha=0.3, color='blue',
                label='Laminar BL')
ax.fill_between(x_turb * 100, 0, delta_turb * 1000, alpha=0.3, color='red',
                label='Turbulent BL')
ax.plot(x_lam * 100, delta_lam * 1000, 'b-', lw=2)
ax.plot(x_turb * 100, delta_turb * 1000, 'r-', lw=2)
ax.plot(x_lam * 100, -delta_lam * 1000, 'b-', lw=2)
ax.plot(x_turb * 100, -delta_turb * 1000, 'r-', lw=2)
ax.fill_between(x_lam * 100, 0, -delta_lam * 1000, alpha=0.3, color='blue')
ax.fill_between(x_turb * 100, 0, -delta_turb * 1000, alpha=0.3, color='red')
ax.axvline(x_tr * 100, color='green', ls='--', lw=2, label='Transition')
ax.axhline(0, color='k', lw=2)  # 平板
ax.set_xlabel('$x$ [cm]', fontsize=12)
ax.set_ylabel('$\\delta$ [mm]', fontsize=12)
ax.set_title('Boundary layer development on flat plate', fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)

# (2) 層流 vs 乱流の速度分布比較
ax = axes[0, 1]
y_norm = np.linspace(0, 1, 200)

# ブラジウス風の近似(3次多項式)
u_lam = 2 * y_norm - 2 * y_norm**3 + y_norm**4
# 1/7乗則
u_turb_profile = y_norm**(1.0/7.0)

ax.plot(u_lam, y_norm, 'b-', lw=2.5, label='Laminar (polynomial approx.)')
ax.plot(u_turb_profile, y_norm, 'r-', lw=2.5, label='Turbulent (1/7 power law)')
ax.fill_betweenx(y_norm, u_lam, u_turb_profile, alpha=0.15, color='gray')
ax.set_xlabel('$u / U_e$', fontsize=12)
ax.set_ylabel('$y / \\delta$', fontsize=12)
ax.set_title('Velocity profiles: laminar vs turbulent', fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3)

# (3) 各境界層厚さの比較(層流vs乱流)
ax = axes[1, 0]
x_all = np.linspace(0.01, L, 500)
Re_x_all = U * x_all / nu

# 層流
d_lam_all = 5.0 * x_all / np.sqrt(Re_x_all)
ds_lam_all = 1.7208 * x_all / np.sqrt(Re_x_all)
th_lam_all = 0.664 * x_all / np.sqrt(Re_x_all)

# 乱流
d_turb_all = 0.37 * x_all / Re_x_all**0.2
ds_turb_all = d_turb_all / 8
th_turb_all = 7 * d_turb_all / 72

ax.semilogy(x_all * 100, d_lam_all * 1000, 'b-', lw=2, label='$\\delta$ (laminar)')
ax.semilogy(x_all * 100, d_turb_all * 1000, 'r-', lw=2, label='$\\delta$ (turbulent)')
ax.semilogy(x_all * 100, th_lam_all * 1000, 'b--', lw=1.5, label='$\\theta$ (laminar)')
ax.semilogy(x_all * 100, th_turb_all * 1000, 'r--', lw=1.5, label='$\\theta$ (turbulent)')
ax.set_xlabel('$x$ [cm]', fontsize=12)
ax.set_ylabel('Thickness [mm]', fontsize=12)
ax.set_title('Boundary layer thicknesses: laminar vs turbulent', fontsize=13)
ax.legend(fontsize=9, ncol=2)
ax.grid(True, which='both', alpha=0.3)

# (4) 局所摩擦係数
ax = axes[1, 1]
Cf_lam_all = 0.664 / np.sqrt(Re_x_all)
Cf_turb_all = 0.0592 / Re_x_all**0.2

# 混合: 層流部 + 乱流部
Cf_mixed = np.where(x_all < x_tr, Cf_lam_all, Cf_turb_all)

ax.semilogy(x_all * 100, Cf_lam_all * 1000, 'b-', lw=2, label='All laminar')
ax.semilogy(x_all * 100, Cf_turb_all * 1000, 'r-', lw=2, label='All turbulent')
ax.semilogy(x_all * 100, Cf_mixed * 1000, 'g--', lw=2.5, label='Mixed')
ax.axvline(x_tr * 100, color='green', ls=':', alpha=0.7)
ax.set_xlabel('$x$ [cm]', fontsize=12)
ax.set_ylabel('$C_f \\times 10^3$', fontsize=12)
ax.set_title('Local skin friction coefficient', fontsize=13)
ax.legend(fontsize=10)
ax.grid(True, which='both', alpha=0.3)

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

この4パネルの図から、境界層理論の主要な要素を包括的に理解できます。

左上の境界層全体像では、層流部(青)から乱流部(赤)への遷移で境界層の成長率が急激に変化する様子が明瞭です。層流部では $\delta \propto x^{1/2}$ のゆっくりした成長ですが、乱流部では $\delta \propto x^{4/5}$ とより急速に成長します。

右上の速度分布比較は、層流と乱流の本質的な違いを示しています。乱流の1/7乗則は壁面近くでの速度勾配が層流よりも急峻です。これが乱流の方が壁面せん断応力(摩擦抗力)が大きい理由です。一方、外部領域では乱流の方が速度が大きく、境界層の「壁」が薄い。このため乱流境界層は逆圧力勾配に対して強く、剥離しにくいのです。

左下の厚さ比較では、同じ $x$ 位置において乱流の $\delta$ が層流の数倍に達していることがわかります。$\theta$ についても同様の傾向が見られます。

右下の摩擦係数では、遷移点での $C_f$ の跳ね上がりが特に印象的です。層流から乱流への遷移が摩擦抗力に与える影響の大きさを如実に示しています。

Python: ブラジウス速度分布の平板上での発達

最後に、平板上の各位置での速度分布を空間的に描画し、境界層内の流速場を可視化します。

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

# ブラジウス方程式を解く
def blasius_ode(eta, y):
    f, fp, fpp = y
    return [fp, fpp, -0.5 * f * fpp]

def shooting(alpha):
    sol = solve_ivp(blasius_ode, [0, 10], [0, 0, alpha],
                    t_eval=np.linspace(0, 10, 500), rtol=1e-10, atol=1e-12)
    return sol.y[1][-1] - 1.0

alpha_opt = brentq(shooting, 0.1, 1.0)
sol = solve_ivp(blasius_ode, [0, 8], [0, 0, alpha_opt],
                t_eval=np.linspace(0, 8, 500), rtol=1e-10, atol=1e-12)
eta_bl = sol.t
fp_bl = sol.y[1]  # u/U = f'(eta)

# パラメータ
U = 15.0       # 主流速度 [m/s]
nu = 1.5e-5    # 動粘度 [m^2/s]

fig, ax = plt.subplots(figsize=(14, 6))

# 平板を描画
ax.plot([0, 150], [0, 0], 'k-', lw=3)
ax.annotate('Flat plate', xy=(75, -1.5), fontsize=12, ha='center')

# 各x位置での速度分布を描画
x_positions = [0.05, 0.15, 0.30, 0.50, 0.80, 1.20]
colors = plt.cm.viridis(np.linspace(0.2, 0.9, len(x_positions)))

for i, x_pos in enumerate(x_positions):
    Re_x = U * x_pos / nu
    delta = 5.0 * x_pos / np.sqrt(Re_x)

    # y座標(物理座標)
    y_phys = eta_bl * np.sqrt(nu * x_pos / U)

    # 速度分布 u(y) = U * f'(eta)
    u_profile = U * fp_bl

    # 描画用にスケーリング(速度を空間的に表示)
    scale = 15.0  # 速度の表示スケール [cm/(m/s)]
    x_display = x_pos * 100 + u_profile * scale / U
    y_display = y_phys * 1000

    # 境界層内の速度分布のみ描画
    mask = y_display < delta * 1000 * 1.5
    ax.plot(x_display[mask], y_display[mask], '-', color=colors[i], lw=2)
    ax.fill_betweenx(y_display[mask], x_pos * 100, x_display[mask],
                     alpha=0.15, color=colors[i])

    # 境界層の厚さを点で表示
    ax.plot(x_pos * 100, delta * 1000, 'o', color=colors[i], ms=6)
    ax.annotate(f'$x={x_pos*100:.0f}$ cm\n$\\delta={delta*1000:.2f}$ mm',
                xy=(x_pos * 100, delta * 1000),
                xytext=(x_pos * 100 - 5, delta * 1000 + 2),
                fontsize=8, color=colors[i])

# 境界層外縁のエンベロープ
x_env = np.linspace(0.005, 1.3, 300)
Re_env = U * x_env / nu
delta_env = 5.0 * x_env / np.sqrt(Re_env)
ax.plot(x_env * 100, delta_env * 1000, 'r--', lw=1.5, alpha=0.7,
        label='$\\delta(x)$ envelope')

# 主流の矢印
for x_arrow in [10, 40, 70, 100]:
    ax.annotate('', xy=(x_arrow + 8, 18), xytext=(x_arrow, 18),
                arrowprops=dict(arrowstyle='->', color='blue', lw=1.5))
ax.text(55, 20, '$U_\\infty$', fontsize=14, color='blue', ha='center')

ax.set_xlabel('$x$ [cm]', fontsize=13)
ax.set_ylabel('$y$ [mm]', fontsize=13)
ax.set_title('Blasius boundary layer: velocity profiles along flat plate', fontsize=14)
ax.set_xlim(-5, 140)
ax.set_ylim(-3, 25)
ax.legend(fontsize=11, loc='upper left')
ax.grid(True, alpha=0.3)
ax.set_aspect('auto')

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

この図は、平板の前縁から下流に向かって境界層が成長する様子を空間的に描画したものです。各位置での速度分布を横向きの棒グラフのように表示しています。

前縁近く($x = 5$ cm)では境界層が薄く壁面近くの速度勾配が急ですが、下流に行くほど境界層は厚くなり、壁面近くの速度勾配はなだらかになっていきます。これは局所摩擦係数 $C_f$ が下流で減少する事実と対応しています。赤い破線は境界層外縁($\delta(x)$)のエンベロープで、$\sqrt{x}$ に比例する成長が見て取れます。

また、どの $x$ 位置でも速度分布の「形」が同じ(相似変数 $\eta$ に対して同一曲線になる)ことが視覚的にも確認できます。これがブラジウスの相似解の本質です。

工学的な応用と発展

圧力勾配のある境界層

ここまでは圧力勾配のない平板を主に扱ってきましたが、実際の翼型や機体では圧力勾配が存在します。

  • 順圧力勾配($dp/dx < 0$、流れが加速): 境界層が薄くなり、安定化する。層流を維持しやすい
  • 逆圧力勾配($dp/dx > 0$、流れが減速): 境界層が厚くなり、不安定化する。やがて境界層剥離(flow separation)が起きる

境界層剥離は、壁面近くの流体がもはや逆圧力勾配に逆らって前進できなくなり、逆流が始まる現象です。剥離が起こる条件は、壁面でのせん断応力がゼロになること:

$$ \tau_w = \mu\left(\frac{\partial u}{\partial y}\right)_{y=0} = 0 $$

剥離は航空機の失速、ディフューザの効率低下、建造物まわりの渦度構造の形成など、多くの工学的問題に関わります。

熱境界層

速度境界層と同様に、温度場にも熱境界層が存在します。壁面温度と主流温度が異なる場合、温度が壁面値から主流値に変化する薄い層が形成されます。熱境界層の厚さと速度境界層の厚さの比はプラントル数 $Pr = \nu / \alpha$($\alpha$ は温度拡散率)で決まります。

  • $Pr > 1$(オイル等): 熱境界層 < 速度境界層
  • $Pr < 1$(液体金属等): 熱境界層 > 速度境界層
  • $Pr \approx 1$(空気): 両者がほぼ同じ厚さ

境界層理論のまとめと限界

プラントルの境界層理論は、高レイノルズ数の流れにおいて非常に強力な解析手法ですが、以下のような限界もあります:

  • 低レイノルズ数($Re \lesssim 100$)では境界層の「薄さ」の仮定が崩れる
  • 大規模な剥離が起こると、境界層近似(放物型の仮定)が成り立たなくなる
  • 3次元境界層、非定常問題、圧縮性流れではさらなる拡張が必要

現代の計算流体力学(CFD)では、直接数値シミュレーション(DNS)やラージエディシミュレーション(LES)により境界層内部の詳細な乱流構造まで解像できますが、境界層理論の知見は物理的理解、初期設計、結果の検証において今なお不可欠です。

まとめ

本記事では、プラントルの境界層理論を基礎から体系的に解説しました。

  • 境界層の概念: 高 $Re$ の流れでは、粘性の影響は物体表面近傍の薄い層($\delta/x \sim 1/\sqrt{Re_x}$)に集中する。プラントルの着想により、流れ場を境界層内(粘性重要)と外(非粘性)に分けて解析できる

  • 境界層方程式: ナビエ・ストークス方程式にオーダー解析を適用し、$x$ 方向の粘性拡散と $y$ 方向の圧力変化を無視して得られる。放物型であり上流から下流へ順次解ける

  • ブラジウス解: 平板上の層流境界層の相似解。射撃法で数値的に解く。$f”(0) = 0.33206$ が壁面せん断応力を決定する

  • 境界層厚さ: $\delta$(99%厚さ)、$\delta^*$(排除厚さ: 質量欠損)、$\theta$(運動量厚さ: 運動量欠損)の3種類。形状係数 $H = \delta^*/\theta$ は層流で $\approx 2.59$、乱流で $\approx 1.29$

  • カルマンの運動量積分方程式: 境界層方程式の積分形。速度分布の近似形を仮定すれば、$\delta$, $C_f$ を求められる実用的手法

  • 層流→乱流遷移: $Re_{x,\text{crit}} \approx 5 \times 10^5$(環境条件で変動)。TS波の増幅、非線形成長、乱流スポットの発生を経て遷移

  • 乱流境界層: 1/7乗則で速度分布を近似。層流より厚く、摩擦抗力が大きいが、剥離に強い

  • 摩擦抗力: 層流 $C_f = 0.664/\sqrt{Re_x}$、乱流 $C_f = 0.0592/Re_x^{0.2}$。乱流は層流の数倍の抗力を生む

境界層理論は流体力学の中核であり、翼型の空力設計、熱伝達、環境流体など多くの分野の基盤です。

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