水流崩散的瞬间 — 雷诺临界数与湍流边界层的多层结构

把第一层网格放进对数区，粘性底层就从计算里消失。壁面摩擦只剩截距B一个数。

慢慢拧开水龙头,水流像玻璃棒一样光滑落下。再多拧一点的瞬间,光滑的水柱突然碎成嘈杂散乱的水花。同样的水、同样的喷嘴,改变的只有速度。是什么决定了这一阈值?1883 年,曼彻斯特的雷诺以一缕染料给出了答案,这一答案最终凝结成一个无量纲数。这里要追的是阈值为什么落在 Re ≈ 2300 附近,以及过了之后边界层怎样分成粘性底层、缓冲层和对数律区。Re 滑块停在 1500 时染料丝是直的,推向 8000 时颜色铺满整根管子。

Reynolds number Re = U·d/ν: 1500Laminar (Re_crit ≈ 2300)

越过 2300 的时候,先晃的只有末端。到 4000 之前,涡一块一块地出现又消失。用 Pause 停在那一帧,转捩并不是一幅图就结束的。

染色丝突然紊乱的临界点

雷诺将染料注入玻璃管中央。他改变管径 $d$ 、平均流速 $u$ 与运动粘度 $\nu$ ,只问一个问题:染料丝何时崩散?

无论组合如何变化,紊乱的开始几乎都发生在同一个无量纲数的相同值附近:

\mathrm{Re} = \frac{u\, d}{\nu}

其中 $u$ 是平均流速, $d$ 是特征长度(管径), $\nu = \mu / \rho$ 是运动粘度(动力粘度除以密度)。Re 越过临界值,染料丝便失去平衡。

Re ≈ 2300 究竟说了什么#

Re 不只是实验常数。把动量方程无量纲化,粘性项前正好乘以 $1/\mathrm{Re}$ 。所以 Re 是 惯性力与粘性力之比。

Re 小 → 粘性占主导,小扰动迅速被均匀化 → 层流维持。
Re 大 → 惯性占主导,小扰动经非线性项自我放大 → 转捩为湍流。

在管内流动,经验上的临界 Re 约 2300。但这一数值并不稳定。光洁打磨的管子、无振动的支架、洁净的入口条件下,曾有报告将层流维持到 Re ≈ 50,000。因此「临界 Re ≈ 2300」是工程实用阈值,转捩的本质在于扰动是否拥有足以放大的 能量通道。同一个无量纲数怎样把昆虫的游动和船的尾迹分开,写在雷诺数如何改写流动尺度里。

平板边界层(boundary layer,壁面附近粘性主导的薄层)的尺度则不同。以前缘距离 $x$ 定义的 $\mathrm{Re}_x = U_\infty x / \nu$ 临界值约为 $3 \times 10^5$ ;以边界层厚度 $\delta$ 定义的 $\mathrm{Re}_\delta$ 在 2700 左右失稳。同样的转捩物理,在不同形状下以不同的无量纲数重新现身。

普朗特看见的边界层 — 外层与内层

雷诺给出转捩律的次年(1904),普朗特提出另一个问题。粘性虽小,为何在壁面附近摩擦仍然不可忽视?

他的回答是:把粘性起作用的区域显式地分离出来。他把无滑移条件下从静止加速至自由流速度的薄层命名为 边界层,边界层之外则按势流处理。靠这一分隔,机翼上方 2.5–25 mm 的薄片就足以追踪粘性的影响。

湍流边界层不是一整块板。把壁面距离 $y$ 适当无量纲化,普遍结构便显现:

y^+ = \frac{y\, u_\tau}{\nu}, \qquad u_\tau = \sqrt{\tau_w / \rho}

$u_\tau$ 是 壁面摩擦速度(wall friction velocity,壁面剪切应力 $\tau_w$ 除以密度后开方)。 $\nu / u_\tau$ 是粘性输运动量的距离 — 粘性长度。也就是说, $y^+$ 就是「以粘性长度度量的壁距」。

壁面定律 — 三个区域

平板湍流边界层内部速度分布以 $u^+ = u / u_\tau$ 无量纲化后,呈现三个区域:

区域	范围	速度分布
粘性底层 (viscous sublayer)	$0 < y^+ < 5$	$u^+ = y^+$
缓冲层 (buffer layer)	$5 < y^+ < 30$	两式光滑过渡
对数律区 (log-law / inertial)	$30 < y^+ < 0.2\, \delta^+$	$u^+ = \frac{1}{\kappa} \ln y^+ + B$

其中 $\kappa \approx 0.41$ 是冯·卡门常数, $B \approx 5.0$ 对应光滑壁。

紧贴壁面处,涡块的运动被壁面束缚,分子扩散主导。在剪切应力近似恒定的区域积分 $\tau_w = \mu \, du/dy$ ,正好得到 $u^+ = y^+$ 。给个尺度感:风速 20 m/s、距前缘 2 m 处, $y^+ = 1$ 约对应 0.02 mm — 粘性底层只有头发丝的五分之一。

对数律区源自普朗特的 混合长理论(mixing length,涡块与周围同化前的平均移动距离)。在壁面附近设混合长正比于距离 $y$ ,即 $\ell = \kappa y$ ,再结合恒定剪切应力假设积分,对数式自然浮现。

下面的图中可调整对数律常数 $B$ 。

Log-law constant B: 5.00 smooth buffer blend shade zones

把 $B$ 从 5.0 调到 5.5,对数直线整体平行上移 — 粗糙壁面会降低 $B$ ,有抽吸(suction)时则升高。整面壁的粗糙状态压缩到一条直线的截距里。

只把对数式放在 y+=30,u+ 就是 (1/0.41)×ln(30)+5.0。ln(30) 约为 3.40,除以 0.41 约为 8.30,加上 B 约为 13.3。同一高度的线性式是 u+=30,两者相差近 17。y+=5 处线性式等于 5。把对数式硬套上去,ln(5) 约为 1.61,u+ 约为 8.9。5 和 8.9 之间的空档就是缓冲层。用一条式子接起来,这一段速度会被算两次。

风速 20 m/s、距前缘 2 m 时 y+=1 是 0.02 mm,所以 y+=5 是 0.10 mm。线性律只活在这 0.10 mm 里。壁函数用的 30 到 300,y+ 增大十倍时 ln 增加约 2.30,除以 0.41 之后 u+ 只上升约 5.6。B 从 5.0 移到 5.5,整个对数层平行移动 0.5。这 0.5 大约是一个十倍程上升量 5.6 的 9%。粗糙度让截距偏了 0.5,壁函数给出的摩擦速度会一起偏。

混合长理论 — 把涡块当作分子

层流中负责输运动量的是分子。牛顿粘性律 $\tau = \mu\, du/dy$ 即结果。湍流中,涡块(eddy,各种尺度的组织化流体团)取代了这一角色。

普朗特借用气体分子运动论中平均自由程决定分子摩擦的思路,把涡块的平均移动距离命名为 混合长 $\ell$ 。脉动速度 $u'$ 近似为该距离上的平均速度差,雷诺剪切应力变为

\overline{u' v'} \approx -\ell^2 \left(\frac{d\bar{u}}{dy}\right)^2

由此定义 涡粘性(eddy viscosity) $\nu_t$ :

\nu_t = \ell^2 \left|\frac{d\bar{u}}{dy}\right|

在壁附近, $\nu_t$ 可比分子粘性大 $10^4$ 倍。空气分子运动粘度 $\nu \approx 1.5 \times 10^{-5}$ m²/s,而地面附近大气涡粘性达 0.01–0.1 m²/s,大尺度大气现象中达 1–100 m²/s。

混合长理论是一个假说,但它支撑了零阶 RANS、 $k$ - $\varepsilon$ 、 $k$ - $\omega$ 的发想骨架。把涡块当作分子的简化在 1925 年提出,百年后仍在工业 CFD 求解器中作为基本选项。

NumPy 中的层流与湍流剖面#

管内层流有精确解 — Hagen–Poiseuille 旋转抛物面。湍流以幂律近似,在同一无量纲断面上要平坦得多。

import numpy as np
import matplotlib.pyplot as plt
 
def laminar_poiseuille(r, R, U_max):
    """Hagen–Poiseuille: 旋转抛物面速度分布。"""
    return U_max * (1.0 - (r / R) ** 2)
 
def turbulent_power_law(r, R, U_centerline, n=7):
    """1/n 幂律湍流近似 (n=7 在 Re ~ 1e5 附近吻合较好)。"""
    return U_centerline * (1.0 - np.abs(r) / R) ** (1.0 / n)
 
def law_of_the_wall(y_plus, kappa=0.41, B=5.0):
    """三区合成: 粘性底层 / 光滑缓冲 / 对数律。"""
    y_plus = np.asarray(y_plus, dtype=float)
    u_visc = y_plus
    u_log = (1.0 / kappa) * np.log(np.maximum(y_plus, 1e-12)) + B
    in_buffer = (y_plus > 5) & (y_plus < 30)
    t = np.clip((y_plus - 5) / 25.0, 0.0, 1.0)
    u_blend = u_visc * (1 - t) + u_log * t
    return np.where(y_plus <= 5, u_visc,
           np.where(y_plus >= 30, u_log, u_blend))
 
R = 1.0
r = np.linspace(-R, R, 201)
u_lam = laminar_poiseuille(r, R, U_max=1.0)
u_tur = turbulent_power_law(r, R, U_centerline=1.0, n=7)
 
bulk_lam = np.trapz(u_lam * np.abs(r), r) / np.trapz(np.abs(r), r)
bulk_tur = np.trapz(u_tur * np.abs(r), r) / np.trapz(np.abs(r), r)
print(f"laminar  bulk/centerline = {bulk_lam:.3f}")
print(f"turbulent bulk/centerline = {bulk_tur:.3f}")
# laminar  bulk/centerline = 0.500
# turbulent bulk/centerline = 0.817
 
y_plus = np.logspace(-1, 3, 200)
u_plus = law_of_the_wall(y_plus)
 
fig, axes = plt.subplots(1, 2, figsize=(10, 4))
axes[0].plot(r, u_lam, label="laminar")
axes[0].plot(r, u_tur, label="turbulent (n=7)")
axes[0].set(xlabel="r/R", ylabel="u / U_max", title="pipe profile")
axes[0].legend()
axes[1].semilogx(y_plus, u_plus)
axes[1].set(xlabel="y+", ylabel="u+", title="law of the wall")
plt.tight_layout()

结论压缩为两行。层流的平均流速恰为中心线速度的一半 — 旋转抛物面与外接圆柱体积之比。湍流的平均流速为中心线的 81%。要达到相同平均流速,湍流只需更低的中心线速度,反过来说 湍流将动量更均匀地分布到整个截面。这是涡粘性比分子粘性大一万倍的可视化表达。

转捩在工业上的含义

层流与湍流的分界并非纯理论问题。

设计 自然层流(NLF) 翼型,延长层流区可使摩擦阻力降低近 50% — 这是 NLF 的核心动机。
热交换器有意激发湍流以提高热传递系数。湍流的动量均匀化同时也是热量的均匀化。
管道输油在临界 Re 附近泵功急剧上升。避开转捩区选择运行点直接影响运营成本。
LES 或 RANS 中,首层格子若落在 $y^+ \approx 1$ 即可解析粘性底层;落在 $30 < y^+ < 300$ 则用壁函数补全。一位数的取值决定了网格预算。

第一层网格放进对数层,粘性底层就没了

第一层单元落在 y+ 30 附近时,u+=y+ 的那段不会被算到。壁面剪切就靠对数式的截距 B。光滑壁的 B 是 5.0。粗糙壁仍用这个数,上面那个 0.5 的平移就会进到摩擦速度里。代码里的 0.817 是 n=7、Re 约十万的光滑管近似。中心速度相同,湍流流量是层流 0.500 的 1.63 倍。把这个比贴到截距不对的管子上,流量就会偏。

这个比只在速度剖面已经充分发展时成立。从进口到抛物线,或者到 1/7 次方,管子还要走多远,发展长度与 Hagen–Poiseuille单独在算。漏掉进口长度、只看出口截面,层流 0.500 和湍流 0.817 的比较就会用到一个两种剖面都还没有的位置上。