排水系统是城市建设的重要组成部分，它的稳定运行对城市环境改善具有重要意义。根据现有研究，污水管道内沉积状况在我国普遍存在，西安市城区有80%的污水管道存在沉积，广州老城区近75%的排水管道存在不同程度的沉积现象^[1,2]。沉积物淤积会导致排水管道堵塞、水力能力损失、污染物质储存等问题，极大程度地影响排水系统的运行^[3]。

以往对排水系统内沉积物的研究多关注于管网整体，通常是根据监测数据，确定整个管网内沉积物容易聚集的位置，主要依据统计经验与物理规律定性分析。近年来有学者尝试结合数值模型，研究不同流动状态下，沉积物在管网中的输送运动，以支撑排水管涵的设计^[4]。为了更好地了解沉积物聚集规律，通过计算流体力学（Computational Fluid Dynamics,CFD）方法构建管内沉积物输送模型，进而数值模拟沉积物运动过程，通过结果分析，为排水系统的设计与维护管理提供参考。计算流体动力学模型构建流程如图1所示。

CFD是通过计算机数值模拟和可视化处理，对流体运动及相关物理现象进行计算机数值分析和研究的一门学科^[5]。通过CFD技术，可以利用计算机分析显示流场现象，并在较短时间内预测流场。

CFD克服了试验研究和理论分析的一些缺点，它成本低、耗时短、数据易于收集，可在不同流动参数（如马赫数、雷诺数、气体性质和模型尺度等）取值范围内研究流体力学的问题，给出流场参数的定量结果，这些常常是试验研究和理论分析难以做到的。因此CFD仿真与理论研究、试验分析方法相互补充，已经普遍应用到各个领域。

以SH市某段典型箱涵断面为例，利用CFD建模软件建立排水箱涵CFD模型，如图2所示。箱涵宽2.7m，高2.1m，长度取12.0m；检查井井筒高2.0m，直径0.7m，位于箱涵模拟长度的后1/3处。令断面左下角为原点，x轴与断面底边平行，y轴与断面垂直边平行，z轴与箱涵长边平行。假设箱涵充满度为0.5，即气液两相入口处各占1.05m，模型坡度为0（计算中坡度通过设置z轴重力分量体现）。

使用多块网格生成技术，将排水箱涵模型划分成多个矩形块，并在检查井井筒处使用O型切分法，将井筒处划分成类似于O型的矩形块。网格划分中将在边界层以及流动变化剧烈的位置加密网格，在变化不剧烈的区域尽量粗化网格，以减少网格规模，提高计算速度，故对箱涵壁面、水气交界面等处加密网格，最终获得网格数量为194万。划分结束后检查网格质量，发现所有网格质量在[0.65,1.0]区间内，达到CFD仿真计算网格质量要求。

使用流体体积函数（VOF）模型模拟多相流，不考虑气、液、固的温度差异；启用稳态计算，即不考虑气水两相流随时间的变化；开启离散相（DPM）模型模拟离散相运动特性。双向耦合中先采用稳态计算，得到颗粒的迹线图；再使用瞬态计算，时间步长取0.005s，计算颗粒群在不同时间点的位置。启用冲蚀/沉积子模型，其他模型参数采用默认值，求解方法采用SIMPLE压力-速度耦合算法，离散相源松弛因子设为0.5。

排水管涵中沉积物的来源主要有2方面：一是旱季积累在城市不同汇水面的地表固体颗粒物质，雨季随雨水径流冲刷进入排水管涵；二是污水管涵中悬浮颗粒物的沉降。排水系统内沉积物运动可以用DPM模型模拟。

在使用DPM模型计算之前，需要确定DPM颗粒的物理属性，包括初始位置、速度、粒径分布、颗粒密度和质量流量等。采用罗森-拉姆勒分布表示颗粒的分布方式，假设颗粒粒径d与颗粒粒径大于d的质量分数Y_d存在式（1）中指数关系：

式中———颗粒平均粒径，μm;

n———颗粒分布宽度常数，n值越小表示分布范围越宽。

采用文献[7]提供的颗粒物数据，颗粒物密度假设为1 034kg/m³。利用整理的d～Y_d对应关系，采用最小二乘法适线，得到粒径分布参数，n=4，珚d=810μm，计算结果如图3所示。

当颗粒与模型物理边界（例如壁面或入口）接触时，使用离散相边界条件确定颗粒轨迹在壁面应满足的条件。可以在每个边界面单独定义模型中的离散相边界条件，可选的边界条件如下：反射边界条件：颗粒在此处反弹而发生动量变化，颗粒反射前后法向及切向的速度变化与用户指定的动量恢复系数相关；捕获边界条件：颗粒在边界处捕获，颗粒的轨迹计算终止；逃逸边界条件：颗粒在此处穿过边界至计算域外，颗粒的轨迹计算终止；内部边界：颗粒在此处穿过内部流动区域。

为简化VOF与DPM耦合模型的计算，在VOF模型流场边界条件的基础上，从液相入口加入面DPM颗粒入射源，入射速度为1.5m/s，质量流量为50g/s，颗粒密度为1 034kg/m³；粒径采用罗森-拉姆勒分布，令最小粒径d_min=350μm，平均粒径珔d=810μm，最大粒径d_max=1 400μm，颗粒分布宽度n=4，粒径采样点为11。考虑到连续相紊流对颗粒的运动效应，启用入射源紊流颗粒弥散子模型。

液相入口、气相入口采用反射边界条件，即颗粒接触这些边界时反射回计算域；液相出口、气相出口采用逃逸边界，即颗粒从这些边界逃逸出计算域；箱涵模型底部壁面采用捕获模型，即颗粒在该边界沉积；箱涵模型侧壁、顶板及井筒采用反射边界条件，即颗粒接触这些壁面即刻反射回计算域，动量恢复系数取默认值，即完全弹性碰撞。

VOF与DPM耦合稳态计算至结果收敛后，使用CFD后处理软件读取计算结果，生成颗粒轨迹线，并以色阶及圆球大小区分颗粒粒径，如图4所示。

由图4可知，面颗粒入射源以颗粒均匀分布于整个入射面的形式，颗粒在液相紊流的作用下，受曳力、重力及浮力的影响，作上下起伏运动。由于颗粒密度与水密度大致相当，而颗粒轨迹线近乎平行于渠底。靠近渠底的颗粒群，在上下起伏的运动过程中，与壁面碰撞，因受壁面水流黏滞力的影响，颗粒动量大大消耗，因此沉积在箱涵底板。本次计算中，DPM颗粒计算结果为：“共追踪22 968个颗粒包，其中逃逸17 531个颗粒包，捕获5 437个颗粒包”，计算结果中的“逃逸”即指从箱涵模型液相出口运动出去的颗粒数目；“捕获”是指沉积在箱涵底板的颗粒数目。在一次CFD仿真计算中，从面源释放的颗粒有20%左右沉积在箱涵底板，80%左右随水流运动至下游。

将沉积在箱涵底板以及从液相出口逃逸的颗粒作颗粒粒径与频率统计，得图5。由图5可知，箱涵底板的沉积颗粒和液相出口的颗粒粒径分布趋势恰好相反。在箱涵底板上，颗粒粒径越大，沉积数量越多。因为颗粒越大，受到的流体阻力越大，动量损失则越大，因此越容易沉积在箱涵底板。而液相出口中，小颗粒物质占比较多，可能因为小颗粒物质受到的流体阻力较小，随流体作悬浮运动。

将沉积在箱涵底板的颗粒位置作图，并以色阶区分颗粒粒径，所得箱涵底板颗粒位置分布见图6。

由图6可知，沉积在箱涵底板的颗粒中，大颗粒占比较多，颗粒沉降的位置主要在箱涵底板两侧，箱涵底板中心位置沉积颗粒数量较少。这是因为在箱涵底板两侧，颗粒受箱涵侧壁的影响，水流黏滞力较大，颗粒动量损失越大，因此沉积在底板两侧的颗粒较多。而在箱涵底板中心位置，水流流速较大，挟带颗粒作迁移运动，因此沉积在底板中心位置的颗粒较少。

根据图6可以预测，由于排水管涵中存在颗粒物质，排水管涵几乎总是存在不同程度的污染物沉积现象。沉积物在排水管涵中受重力、浮力及曳力等力的作用下，存在冲刷及累积等运动过程。在排水箱涵中，沉积物厚度随位置不同而略有不同，其中靠近侧壁的沉积物较厚，位于箱涵底板中心位置的沉积物较薄。由于箱涵底板存在着沉积物，一定程度上保护了底板免受水流冲刷及污水中酸性物质对混凝土的腐蚀，因此，箱涵底板腐蚀程度较侧壁腐蚀程度低。

使用VOF稳态与DPM瞬态耦合计算至计算结果收敛，其中瞬态步长取0.005s，面颗粒入射源入射起始时间为第0s，终止时间为第1s，即入射总质量为50g的颗粒进入计算域。取时间点0s、1.5s、3.0s、4.5s、6.0s和7.5s的颗粒位置作图，如图7所示。

由图7可知，颗粒进入流场中，随水流挟带向前运动，由于壁面边界层的存在，靠近壁面处水流速度较低，靠近中心处速度较高。颗粒受到水流作用的曳力及浮力，不同颗粒的动量耗散并不一致。运动一定时间以后，存在一定的速度差，因此颗粒群产生了明显的尾迹现象。其中，靠近箱涵中心的颗粒速度快、位置靠前；靠近壁面的颗粒速度慢、位置靠后。这些速度较慢的颗粒在运动一定时间后有可能沉积在箱涵底板，从而形成沉积层，根据颗粒群尾迹的数量分布可知，靠近侧壁处沉积层较厚。

利用CFD数值模拟技术，研究了排水管涵内沉积物输送规律。根据模拟结果，颗粒在液相紊流的作用下，受曳力、重力及浮力的影响，作上下起伏的运动。由于颗粒密度与水密度相当，故而颗粒的轨迹线几乎平行于渠底。靠近渠底的颗粒群运动过程中，与壁面碰撞，因受壁面水流黏滞力的影响，颗粒动量损失较大，从而沉积在箱涵底板。沉积在箱涵底板的颗粒中，大颗粒占比较多，颗粒主要沉降在箱涵底板两侧，箱涵底板中心颗粒沉积较少。因此在排水箱涵中，沉积物厚度随位置不同而略有不同，其中靠近侧壁的沉积物较厚，底板中心的沉积物较薄。方形排水箱涵底边两侧沉积物积累速度要超过中心位置，在此作用下排水系统的水力性质也会改变。通过模型的构建验证了CFD在沉积物运动模拟方面的可行性，后续可以根据坡度、管径、断面形状以及不同的粒径性质，进行更深入的研究，以服务于排水管涵设计运营维护。