颗粒物沉降估算法:从斯托克斯公式到工程实践

颗粒物沉降估算法:从斯托克斯公式到工程实践 简介《空气中颗粒物沉降估算法》课件是一份面向环境科学、化学工程及大气污染控制方向学生与从业者的专业学习教案。内容围绕颗粒物在空气中的沉降过程展开系统梳理了颗粒直径、密度、流体性质等因素对沉降速度的影响并重点讲解层流区的斯托克斯定律、过渡区的阿伦公式及湍流区的牛顿公式既有公式推导也有流型判断与试差法、判据法说明。对于非球形颗粒、不均匀颗粒以及干扰沉降、壁效应、端效应等工程实际问题也给出了校正和处理思路。课件后段还结合降尘室、沉降槽、离心沉降设备等典型装置分析了停留时间、生产能力、设备尺寸计算等关键设计参数能帮助读者建立从理论到工程应用的完整知识链。压缩包内共1个pptx文件仅442KB内含约30页内容结构紧凑、逻辑清晰适合课程复习、教学备课或自学入门使用。目前已有74人学习下载。1. 颗粒物沉降估算法先弄清“多久落地”再谈浓度控制空气中颗粒物到底多久才会落回地面一个反直觉的结论是按矿物尘密度估算10 微米颗粒在 1 米高的静止空气里沉降大约需要 2 分钟而 1 微米颗粒要悬停 3 个小时以上。这种量级上的巨大差异正是“颗粒物沉降估算法”要解决的核心问题。它不研究气流里复杂的浓度分布而是先把重力、阻力和颗粒粒径之间的关系算明白再用沉降速度去估算沉积通量、清除时间以及表面污染风险。洁净室风量设计、机房防尘、传感器进风口选型、生产车间职业暴露评估都会用到这套估算思路。本文从斯托克斯沉降公式出发逐步落到可复现的 Excel 模板、Python 脚本和实测校验手段适合需要把理想公式转成工程结论的从业者阅读。2. 用斯托克斯沉降公式估算颗粒物终端沉降速度的建模过程2.1 颗粒物沉降的受力平衡与终端沉降速度公式由来颗粒在静止空气中下落受力并不复杂重力向下空气浮力和阻力向上。颗粒刚释放时速度为零阻力也是零于是加速速度变大后阻力按线性关系上升很快加速度归零颗粒进入匀速下落状态。这个匀速速度就是终端沉降速度符号记为 v_t。阻力项在低雷诺数条件下满足斯托克斯定律F_d 3 · π · μ · d_p · v其中 μ 是空气动力粘度20 ℃ 时约 1.81×10⁻⁵ Pa·sd_p 是颗粒直径v 是当前下落速度。把重力、浮力和阻力三者平衡整理后得到v_t (d_p² · ρ_p · g) / (18 · μ)这里的 ρ_p 是颗粒物真实密度g 是重力加速度。空气密度只有颗粒密度的千分之一左右浮力项可以直接忽略不会给工程估算带来可感知的误差。以 10 微米、密度 2600 kg/m³ 的石英粉尘为例代入上式v_t (10×10⁻⁶)² · 2600 · 9.8 / (18 × 1.81×10⁻⁵) ≈ 0.0078 m/s换算成工程习惯单位约为 28 米每小时也就是在 1 米高的房间里自由落体大约 2 分钟。1 微米颗粒的直径只有前者的十分之一但由于公式里直径是平方项速度直接缩到不足百分之一变成了 0.28 米每小时。这里有一个容易误用的前提斯托克斯沉降公式只适用于雷诺数 Re 小于 1 的粘性绕流区。计算 Re 的表达式为 Re ρ_air · d_p · v_t / μ当粒径达到 100 微米以上时Re 已经接近或超过 1实际阻力会比线性模型偏大公式给出的沉降速度也会偏快。做工程估算时超过 100 微米的粗颗粒建议查阅阻力系数曲线修正而不是硬套这个公式。2.2 亚微米颗粒的滑移修正Cunningham 因子怎么取斯托克斯定律的前提是“颗粒表面处空气速度为零”也就是把空气当成连续介质。但当粒径小到微米级以下、和空气分子平均自由程可比时颗粒会在气体分子之间发生“滑移”受到的阻力比连续介质模型更小实际沉降速度会高于公式计算值。工程上引入 Cunningham 滑移修正因子 C_c把原公式改写成v_t (d_p² · ρ_p · g · C_c) / (18 · μ)C_c 的经验表达式为C_c 1 Kn · [1.257 0.4 · exp(-1.1 / Kn)]其中 Kn 2λ / d_p 是努森数λ 为空气分子平均自由程标准条件下约 0.065 微米。下表给出不同粒径下的典型取值。粒径 (μm)Knudsen 数 KnCunningham 因子 C_c对沉降速度的影响0.0113约 22不修正会严重低估0.11.3约 2.86速度差约 3 倍10.13约 1.16修正约 16%100.013约 1.017可忽略从表中能直接看到工程上最常用的口径粒径大于 10 微米时可以不乘 C_c误差在 2% 以内1 微米附近建议乘上 1.16而 0.1 微米以下的亚微米颗粒如果不修正估算速度会比实际低 60% 以上。这也是很多新手用公式算细颗粒沉降时“总觉得太慢”的原因。需要特别提醒的是C_c 的数值与温度和压力有关因为空气分子平均自由程会随环境变化。0 ℃ 和 40 ℃ 之间λ 的变化会给 0.1 微米颗粒带来百分之十几的速度偏差。严谨的估算应写明环境温度而不是默认“常温”。2.3 一个可复现的估算模板Excel 公式与 Python 函数手算只适合单个粒径实际工作里往往要同时看 0.1、1、2.5、10 微米几条曲线这时候用表格或脚本更可靠。Excel 的复用方式很直接A1 填粒径微米B1 填颗粒密度kg/m³C1 填空气粘度Pa·sD1 计算努森数E1 计算修正因子F1 输出沉降速度m/s。D1 0.13 / A1 E1 1 D1 * (1.257 0.4 * EXP(-1.1 / D1)) F1 (A1 * 10^-6)^2 * B1 * 9.8 / (18 * C1) * E1如果嫌 Excel 公式太长或者要批量生成对比表Python 更适合。下面是一个最小可运行版本。import math MU_AIR 1.81e-5 # 20°C 空气动力粘度, Pa·s LAMBDA_M 0.065e-6 # 空气分子平均自由程, m G 9.8 def cunningham(dia_m): kn 2 * LAMBDA_M / dia_m # 努森数 return 1 kn * (1.257 0.4 * math.exp(-1.1 / kn)) def settling_velocity(dia_um, density2600.0, slipTrue): d dia_um * 1e-6 # 微米转米 vt d * d * density * G / (18 * MU_AIR) if slip: vt * cunningham(d) return vt for dia in [0.1, 1, 2.5, 10, 50]: print(fd {dia:6.1f} um, v_t {settling_velocity(dia):.6f} m/s)逻辑上settling_velocity 函数先做单位换算再套斯托克斯公式最后根据参数 slip 决定是否做滑移修正。cunningham 函数单独拆出来便于在多个公式间复用。参数上的关键点有三个dia_um 是几何粒径或空气动力学当量直径不是质量中值直径density 应填颗粒真密度不要填堆积密度返回值的单位是 m/s要换成 cm/s 或 m/h 时直接乘系数即可。3. 颗粒物沉降通量与室内清除时间的估算从公式到工程参数3.1 沉降通量公式与监测点位置的设定算出终端沉降速度只完成了一半工程上更关心“单位时间有多少颗粒落到某个表面上”这对应沉降通量 JJ C · v_t其中 C 是颗粒物的体浓度v_t 是沉降速度。通量单位常写成 μg/(m²·s) 或 mg/(m²·day)。举一个直观例子室内 TSP 平均浓度为 100 μg/m³粒径 10 微米的颗粒沉降速度取 0.0078 m/s则单位面积上的沉降通量约为 0.78 μg/(m²·s)一天下来接近 67 mg/m²。如果表面是朝上的桌面或设备顶板这个值会成为可观的积尘量。这个公式常被用来评估监测点位设置是否合理。颗粒物浓度传感器如果紧挨着通风死角或设备散热面传感器内部的进风通道会在几天内积累一层细灰导致采样流量下降。用沉降通量估算可以在部署前判断该位置是否需要加装过滤网或缩短清洁周期而不是等设备报警后再处理。需要注意的是J C·v_t 描述的是“颗粒垂直于表面方向净通量”的重力分量。真实表面还有气流湍流扩散、热泳和静电作用因此工程上更常用“有效沉积速度 v_d”替代 v_t把一切非重力输运效应折算成等效沉降速度。3.2 室内浓度衰减模型与有效沉降速度 v_d室内颗粒物浓度随时间的变化可以写成一阶衰减方程V · dC/dt - (v_d · A Q) · C其中 V 是房间体积A 是有效沉积表面积Q 是通风体积流量。解这个微分方程得到C(t) C0 · exp[-(v_d · A / V k) · t]k 是换气次数单位 1/s工程上习惯用 h⁻¹ 表示后除以 3600 换算。这个模型的价值在于它把“重力沉降”和“通风稀释”放进了同一个指数项里可以直观比较谁占主导。v_d 的取值不能直接用 v_t因为颗粒到达表面的方向不一样气流状态也不一样。下表是常见工程参考量级以 v_t 为基准按表面朝向划分。表面类型v_d 参考范围说明水平朝上地板、桌面、机柜顶0.51.0 倍 v_t重力直接叠加在通量上垂直表面墙面、侧板0.10.3 倍 v_t主要靠湍流输运携入水平朝下吊顶底面0.010.1 倍 v_t颗粒需逆重力扩散到达实际取数时还要考虑房间气流速度。空调风口直吹的区域垂直表面 v_d 可能提升到 0.5 倍 v_t完全无风的环境中水平朝下的沉积甚至可以忽略。建议先用上表做敏感性分析而不是只取一个“标准值”。3.3 洁净室与机房的沉降量和清除时间计算示例用一个典型的小型机房来演示整套计算。房间尺寸 6 m × 6 m × 3 m体积 V 108 m³有效沉积面积只算地面 A 36 m²初始颗粒物浓度 C0 100 μg/m³换气次数 k 0.5 h⁻¹。假设颗粒密度 2600 kg/m³v_d 取 0.5 倍 v_t。计算逻辑可直接用以下 Python 片段import math def clear_time(dia_um, c0100.0, ratio0.1, v108.0, a36.0, ach0.5): vt settling_velocity(dia_um, slipTrue) vd 0.5 * vt k_s ach / 3600.0 # 换气次数换算为 1/s lam vd * a / v k_s # 总衰减常数 t_s math.log(c0 / (c0 * ratio)) / lam return lam, t_s / 60.0 for dia in [1, 2.5, 10]: lam, t_min clear_time(dia) print(f{dia} um: λ{lam:.6f} 1/s, 降到10%约{t_min:.1f} min)逻辑说明lambda 为沉降和换气两部分的叠加衰减到 10% 需要用初始浓度除以目标浓度取自然对数再除以总衰减常数。参数方面v_d 按地面积取的 0.5 倍 v_t实际多表面场景应加权各类表面积ach 传入的是每小时换气次数代码里统一转成秒避免单位错位。三个粒径的计算结果对比如下粒径v_t (m/s)v_d (m/s)衰减常数 λ (1/s)降到 10% 用时1 μm7.8×10⁻⁵3.9×10⁻⁵1.52×10⁻⁴约 4.2 h2.5 μm4.9×10⁻⁴2.4×10⁻⁴2.22×10⁻⁴约 2.9 h10 μm7.8×10⁻³3.9×10⁻³1.44×10⁻³约 27 min从这个表可以看得很清楚10 微米粒径下沉降项贡献是换气项的近十倍关不关新风影响不大但 1 微米颗粒恰恰反过来换气对清除的贡献占据了主要地位。在设计监测或净化方案时如果不先做这一步估算很容易把资源投到错误的方向上。4. 把颗粒物沉降估算法固化成一个可交互的筛选工具4.1 用标准库实现多工况批量估算的完整脚本单点计算适合教学演示工程选型时却要同时比较多种粉尘和多种粒径。可以把上面两个函数整合成一个批量估算脚本输出可直接贴进表格的结果。import math MU_20C 1.81e-5 LAMBDA_M 0.065e-6 G 9.8 def cunningham_factor(dia_m): kn 2 * LAMBDA_M / dia_m return 1 kn * (1.257 0.4 * math.exp(-1.1 / kn)) def settling_velocity(dia_um, density2600.0): d dia_um * 1e-6 vt d * d * density * G / (18 * MU_20C) return vt * cunningham_factor(d) materials { 石英尘: 2600.0, 燃煤飞灰: 2200.0, 水泥尘: 3100.0, 木材尘: 700.0, } diameters [0.1, 0.5, 1.0, 2.5, 5.0, 10.0, 50.0] for mat, rho in materials.items(): row [f{mat:6s}] for d in diameters: row.append(f{settling_velocity(d, rho):.4g}) print( .join(row))这段脚本有三点设计意图。第一材料密度放在一个字典里便于扩展新粉尘密度值尽量使用真密度而不是表观密度。第二粒径序列跨度从 0.1 到 50 微米覆盖了从细颗粒到粗颗粒的完整区间适合直接画 log-log 图。第三输出用格式化字符串对齐运行时复制到表格或 CSV 都方便。密度对沉降速度的影响是线性的因此同样的粒径下水泥尘比木材尘快约 4.4 倍。这种差异在职业危害评估里意义重大同样是 10 微米颗粒不同材质的沉降行为完全不同笼统说“PM10 会落下来”没有工程价值。4.2 粒度、温度、密度三个参数对估算结果的灵敏度沉降速度公式里最敏感的参数是粒径因为它以平方关系出现。粒径偏差 10%v_t 偏差 21%粒径偏差 20%v_t 偏差 44%。颗粒物测量时粒径分布往往有一定宽度用中位粒径代替整个分布相当于人为引入至少 20% 的沉降速度误差。温度的影响通过空气粘度间接体现。粘度随温度升高而增大沉降速度随之下降。常见工程温度范围内可以参考这张表环境温度空气粘度近似值 (Pa·s)相对 20 ℃ 的 v_t 偏差0 ℃1.71×10⁻⁵约 5.8%20 ℃1.81×10⁻⁵基准40 ℃1.91×10⁻⁵约 -5.2%密度是最容易出错的参数。如果颗粒物成分未知很多人直接用水的密度 1000 kg/m³这会让石英粉尘的估算速度偏小 60% 以上。建议先用 X 射线荧光或灰分分析确认成分再做估算临时评估时宁可给密度一个区间也要避免用单一猜测值。三个参数的偏差可以用一阶近似合并dv_t/v_t ≈ 2·δd/d δρ/ρ - δμ/μ。这个公式不是严格误差传递但能快速判断哪项误差源主导。多数场景下粒径测量误差占了大头密度和温度反而在其次。4.3 输出到培训 PPT 的图表规范对数坐标与原数据表把估算结果做成培训用图表时最常见的错误是线性坐标压缩了细颗粒数据。1 微米和 10 微米的速度差了两个数量级线性坐标下小粒径几乎贴在地轴上。规范的呈现方式是双对数坐标横轴为粒径1100 μm纵轴为 v_tm/s不同材料密度画成一组平行直线组斜率固定为 2。import matplotlib.pyplot as plt fig, ax plt.subplots() d_um [0.1 * i for i in range(1, 100)] # 0.1 到 10 μm d_um [10 i for i in range(0, 91)] # 10 到 100 μm for rho in [700, 2200, 2600, 3100]: vt [settling_velocity(d, rho) for d in d_um] ax.loglog(d_um, vt, labelfρ{rho} kg/m³) ax.set_xlabel(粒径 (μm)) ax.set_ylabel(终端沉降速度 (m/s)) ax.legend()画图代码本身不复杂真正要在 PPT 里写清楚的是边界条件空气温度 20 ℃、压力 101.325 kPa、颗粒为球形、忽略壁面效应。尤其要把“球形”三个字写进注释里否则非球形颗粒会按流体动力学直径处理实际沉降速度比球形假设偏慢。这类培训图表通常还要附带完整数据表不要只留曲线而丢了数值。将 4.1 节的运行结果放进附页读者才能对照曲线做插值而不是重新推导公式。5. 用实测沉降数据反向校验估算结果的三种技巧5.1 双高度采样反推实际沉降速度理论估算再严谨也必须接受现场数据校验。最简单的方法是用两个不同高度的采样点同时测浓度一个放在 0.5 m 高度一个放在 2.5 m 高度。如果颗粒物以沉降为主低处浓度应高于高处。把两个测点的浓度差和采样时间记录下来结合房间体积与表面积可以反推出实际有效沉降速度。具体操作是用表面皿或无尘培养皿做被动沉降板。称重干净皿片的初始质量放在水平台面上暴露 24 小时同步记录平均浓度 C_avg。暴露结束后再次称重用增重除以皿口面积和暴露时间得到实测沉降通量 J_actual再除以 C_avg 就是实测 v_d。将这个值和 0.5 倍 v_t 的理论值对比偏差在 30% 以内说明估算口径合理超过 50% 就要检查是不是有再悬浮或局部气流干扰。5.2 浓度衰减曲线拟合提取有效沉降参数比沉降板更省事的方式是利用室内浓度自然衰减数据做线性拟合。关闭新风和净化设备释放一次示踪颗粒物或等待背景浓度稳定然后连续记录浓度 C(t)。对 ln(C/C0) 关于时间 t 作图斜率的绝对值就是总衰减常数 λ。import numpy as np # t_h: 时间序列, 单位小时 # ln_c: ln(C/C0) 序列 t_h [0.0, 0.5, 1.0, 1.5, 2.0] ln_c [0.0, -0.31, -0.57, -0.88, -1.20] slope, _ np.polyfit(t_h, ln_c, 1) # lambda_s 单位 1/s lam_s -slope / 3600.0 # 若换气次数 k_h 已知分离出沉降项 # vd (lam_s - k_h / 3600.0) * V / A逻辑在于一阶衰减模型里 λ v_d · A / V k。只要换气次数已知把拟合得到的 λ 减掉 k剩余部分就能换算成 v_d。前提条件是两个一是测量期间没有新颗粒物源二是浓度足够高测量噪声不会掩盖衰减趋势。数据点少于 6 个、拟合 R² 低于 0.9 时结果不建议采用。5.3 估算偏差超过 30% 时的参数排查顺序实测数据和估算对不上先不要怀疑公式按下面顺序排查第一检查粒径口径。计算用的是几何粒径还是空气动力学当量直径如果是 D50 质量中位径细颗粒占比高会造成整体沉降偏慢。第二检查密度取值。粉尘成分未知时把 2600 kg/m³ 换成 1600 kg/m³ 这类错误操作非常常见。第三检查环境温度冬季和夏季粘度差约 10%足够让估算值移动 5% 左右。第四看气流状态。估算值偏小往往是因为房间湍流增强了表面沉积偏大则可能是气流把颗粒重新卷扬起来。也有两种特殊情况需要考虑高浓度细颗粒物会发生凝并使粒径分布向粗端移动实际沉降速度高于单颗粒估算非球形颗粒如纤维状粉尘阻力系数比球形大实际速度低于球形假设。把这些检查项按顺序过一遍基本能把偏差来源定位到具体参数。校验通过的 v_d 值才能写进下一版培训教案和选型清单里。本文还有配套的精品资源点击获取