水滴为何是圆的，杯沿的水为何会爬升 — 表面张力与接触角

1881年，Agnes Pockels在厨房水槽里，用一颗纽扣、一根线和一个白铁盆测量了水的表面张力。她无法上大学，洗碗时观察水面上漂浮的油膜，并用自制装置测出了面积与张力的关系。本文沿着这个厨房实验所涉及的物理 — 表面张力与接触角 — 一直讲到CFD如何在网格上再现这种力（CSF模型）。最后，我们用真实数值验证毛细上升。

把接触角滑块在 30 度与 120 度之间移动，同样体积的液滴会贴在底面上，或立成球。

把 $\theta$ 降到30度以下，液滴会摊平；推到120度以上，它几乎像球一样隆起。接触线上黄色弧线所指的角，正是Young平衡。

表面张力是分子在表面上损失的能量

液体内部的分子四面都被邻居包围。表面的分子有一侧是空的。它正好缺了这一部分的结合能。

自然为了减少这种损失而最小化表面积。所以失重中的水滴会成为同体积下表面积最小的形状，也就是球。

表面张力 $\sigma$（单位长度上的力，N/m）表示"撑大表面所需的能量"。常温下水–空气界面约为 $\sigma \approx 0.072$ N/m。

由表面张力主导的世界用无量纲数来区分。Weber数（惯性/表面张力）$We = \rho U^2 L / \sigma$ 小，则液滴不会破碎。Bond数（重力/表面张力）$Bo = \rho g L^2 / \sigma$ 小，则液滴不被重力压扁而保持圆形。在毛细管中，这两个数都变小。

Young–Laplace — 弯曲的面向内挤压#

面一旦弯曲，表面张力就产生向内的合力。结果是压力跨越界面发生跳变。

$\Delta p = \sigma\left(\frac{1}{R_1} + \frac{1}{R_2}\right)$

这里 $\Delta p$ 是内侧更高的压差，$R_1,R_2$ 是两个主曲率半径。球形液滴时 $R_1=R_2=R$，所以 $\Delta p = 2\sigma/R$。

关键在于：半径越小，内部压力越高。一个5 µm的雾滴内部比外界约高29 kPa。这就是小气泡被大气泡吞并的原因。

接触角 — 三种张力的拔河

液滴落在固体上时，三个界面沿一条线相遇：固–气（$\sigma_{sg}$）、固–液（$\sigma_{sl}$）、液–气（$\sigma_{lg}$）。它们在水平方向达到力平衡的角度就是接触角 $\theta$。

$\sigma_{sg} = \sigma_{sl} + \sigma_{lg}\cos\theta$

这是Young方程。整理后 $\cos\theta = (\sigma_{sg} - \sigma_{sl})/\sigma_{lg}$。

若 $\theta < 90^\circ$，液体喜欢固体而铺展（润湿，亲水）。若 $\theta > 90^\circ$，液滴收缩成珠（不润湿，疏水）。荷叶上的水珠 $\theta$ 超过150度。

CSF — 把面上的力变成网格上的体积力#

问题出在CFD。表面张力作用在零厚度的面上。但网格单元有体积。怎样把面上的力放进单元里？

Brackbill的CSF（Continuum Surface Force）模型把面力转换成在界面附近铺开的体积力。

$\mathbf{f}_{sv} = \sigma\,\kappa\,\hat{\mathbf{n}}\,\delta_s$

$\kappa$ 是界面曲率，$\hat{\mathbf{n}}$ 是界面法线，$\delta_s$ 是只在界面附近接近1的表面delta函数。法线和曲率由相指示函数 $c$（例如VOF的体积分数）得到。

$\hat{\mathbf{n}} = \frac{\nabla c}{|\nabla c|},\qquad \kappa = -\nabla\cdot\hat{\mathbf{n}}$

曲率的精度决定这个模型的成败。直接用 $\nabla c$ 会给曲率引入噪声，产生寄生电流。所以要先用距离加权插值或最小二乘平面拟合把法线磨光滑，再算曲率。

用壁面附着强制接触角

接触壁面的界面要多受一个条件：壁面处的接触角必须等于我们指定的 $\theta$。

CSF通过强制倾斜近壁单元的法线来实现这一点。用壁面法线 $\hat{\mathbf{n}}_w$ 和切线 $\hat{\mathbf{t}}_w$ 重写界面法线。

$\hat{\mathbf{n}} = \hat{\mathbf{n}}_w\cos\theta + \hat{\mathbf{t}}_w\sin\theta$

倾斜法线会改变曲率，而曲率又生成表面张力，把液滴拽向指定的接触角。这就是wall adhesion（壁面附着）模型。要处理三相相遇线的运动，还需要额外的模型。

Python — 用数值验证Jurin定律#

把细管浸入水中，水会自己往上爬。这就是毛细上升。平衡高度由Young–Laplace与静水压的平衡得出。

$h = \frac{2\sigma\cos\theta}{\rho g r}$

$r$ 是管半径，$\theta$ 是接触角。管越细（$r$ 越小），爬得越高。

import numpy as np
 
SIGMA = 0.072    # N/m, 水-空气 (常温)
RHO   = 1000.0   # kg/m^3
G     = 9.81     # m/s^2
 
def young_laplace_dp(sigma, R1, R2):
    """跨越弯曲界面的压力跳变 dp = sigma*(1/R1 + 1/R2)."""
    return sigma * (1.0 / R1 + 1.0 / R2)
 
def jurin_capillary_rise(r, theta_deg, sigma=SIGMA, rho=RHO, g=G):
    """半径 r 毛细管中的平衡上升高度 (m)."""
    theta = np.radians(theta_deg)
    return 2.0 * sigma * np.cos(theta) / (rho * g * r)
 
# 球形水滴 (R1=R2=R): 半径越小内部压力越高
for R_um in (5, 50, 500):
    R = R_um * 1e-6
    print(f"R={R_um:4d} um -> dp = {young_laplace_dp(SIGMA, R, R):8.1f} Pa")
 
# Jurin定律: 管越细爬得越高
print()
for r_mm in (0.1, 0.5, 1.0):
    h = jurin_capillary_rise(r_mm * 1e-3, theta_deg=30)
    print(f"r={r_mm:.1f} mm -> h = {h*1000:6.1f} mm")

输出:

R=   5 um -> dp =  28800.0 Pa
R=  50 um -> dp =   2880.0 Pa
R= 500 um -> dp =    288.0 Pa
 
r=0.1 mm -> h =  127.1 mm
r=0.5 mm -> h =   25.4 mm
r=1.0 mm -> h =   12.7 mm

改变管半径和接触角，水柱会按公式走过输出里的 127.1 mm、25.4 mm 和 12.7 mm。

把接触角推到90度以上，$\cos\theta$ 变为负值，液面反而下降。这正是水银在玻璃管中被压低（凸液面）的同一行为。

半径缩小到十分之一，压力跳变和上升高度都增大十倍。这正是植物用毛细与蒸腾把水提升到几十米高的机理的一部分。

输出的三行，就是把代码里的 0.072 N/m 代进 $\Delta p = 2\sigma/R$。 5 µm 是 28800 Pa，50 µm 是 2880 Pa，500 µm 是 288 Pa。 半径变成 10 倍、100 倍，压力跳变就变成 1/10、1/100。 28800 除以 288 等于 100，500 µm 除以 5 µm 也是 100。 前面说的约 29 kPa，是把 28800 Pa 四舍五入之后的说法。

Jurin 高度用的是同一组常数。30 度、0.1 mm 给出 127.1 mm，0.5 mm 给出 25.4 mm，1.0 mm 给出 12.7 mm。 127.1 除以 5 等于 25.42，所以 25.4 是保留一位小数的结果。 127.1 除以 10 等于 12.71，与 12.7 一致。 荷叶上的 150 度只把余弦的符号反过来，同一根 0.1 mm 的管子，高度是向下的 127.1 mm。 把管半径代进 $Bo = \rho g L^2 / \sigma$，0.1 mm 约为 0.00136，1 mm 约为 0.136。 两个都小于 1。这三根管子的高度由表面张力决定。

只问平衡高度时不要打开 CSF#

有管半径和接触角，Jurin 就够了。0.1 mm，30 度，127.1 mm。 这几个数用不到网格上的曲率。只问平衡高度时，CSF 只是在加噪声。

液滴在壁上铺展、合并之后，计算就换了。壁面附着只是把法向扳到指定角度，曲率一粗糙，静止液滴也会流。 下一步的限制是表面张力的时间步，不是接触角滑块。 用隐式打开的时间步窗口在那个位置。