《Water Resources Research》：A Semi-Resolved VOF-DEM Coupling Method for Bedload Transport of Non-Spherical Particles
编辑推荐：
粗颗粒泥沙的推移质输移涉及复杂的颗粒-流体相互作用，这些相互作用控制着泥沙的运动和河床表面的演变。在本研究中，开发了一种半解析CFD-DEM方法来模拟推移质泥沙动力学，其中使用高斯核函数计算每个颗粒的背景流体信息。明确考虑了颗粒形状效应和自由表面动力学。采用粘
粗颗粒泥沙的推移质输移涉及复杂的颗粒-流体相互作用，这些相互作用控制着泥沙的运动和河床表面的演变。在本研究中，开发了一种半解析CFD-DEM方法来模拟推移质泥沙动力学，其中使用高斯核函数计算每个颗粒的背景流体信息。明确考虑了颗粒形状效应和自由表面动力学。采用粘结球模型构建圆盘状非球形颗粒，并执行了一系列验证算例以评估所提方法的可靠性。在不同水流条件下对球形和非球形颗粒进行了比较研究，包括单颗粒沉降、沉降过程中的双颗粒相互作用以及推移质输移模拟。此外，所提方法能够模拟床面自组织结构的演变，包括颗粒聚群形态和叠瓦状结构。结果表明，考虑颗粒尺度的水流扰动对于准确复现床面自组织结构至关重要，而这是传统的非解析CFD-DEM方法无法捕捉的。
推移质输移是山区河流中泥沙运动的主导模式。在砾石河床河流中，推移质的输移和沉积不仅对河床演变和河流地貌形态产生显著影响，还在维持河流生态环境中发挥关键作用。因此，准确模拟和预测山区河流中砾石的推移质输移对河流系统的保护和管理至关重要。目前广泛应用于泥沙输移的计算流体力学(CFD)方法是Euler-Euler模型，该模型无法解析流场内部的详细局部相互作用或捕捉颗粒-颗粒之间的碰撞，而这些对于准确描述推移质泥沙动力学至关重要。基于Euler-Lagrange框架的CFD-DEM模型被越来越多地用于高保真模拟复杂流体-颗粒系统。根据流场解析度的不同，耦合策略通常分为全解析、半解析和非解析方法。全解析方法计算成本极高，限制了其在大尺度计算域的应用；非解析方法无法捕捉单个颗粒与周围流体之间的近场流体动力学相互作用；而半解析方法在准确性和效率之间提供了折中方案。然而，半解析CFD-DEM方法尚未经过系统性验证或应用于推移质泥沙输移模拟。此外，传统CFD-DEM模拟通常将泥沙颗粒表示为球形，这限制了其复现真实推移质动力学的能力。自然砾石的形状效应增强了颗粒间摩擦，可能主导整体泥沙输移率，并影响自组织床面结构的形成。因此，在真实的推移质输移模拟中纳入颗粒形状效应至关重要。
针对上述问题，研究人员开发了一种半解析CFD-DEM耦合方法，用于模拟粗颗粒的推移质输移过程，并明确考虑了自由表面动力学和颗粒形状效应的影响。研究得出结论，考虑颗粒尺度的局部流场结构对于准确复现泥沙沉积过程和床面自组织结构至关重要，该成果在计算效率与解析精度之间取得了平衡，为数值研究山区河流推移质输移提供了一种新颖且实用的方法，对河床演变和河流地貌形态过程的研究具有重要意义。该论文发表在《Water Resources Research》上。
在关键技术方法方面，研究人员采用基于流体体积(VOF)方法的两相流模型来处理流体相计算，并捕捉气液自由表面动力学。在流体求解中，采用Navier-Stokes(NS)方程并引入气液体积分数和固液体积分数以考虑相间作用。颗粒相的运动基于离散元方法(DEM)求解，为表征非球形颗粒，采用了粘结球模型，每个颗粒由四个呈正方形排列的球形子单元构建而成，颗粒相的碰撞计算采用Hertz-Mindlin接触模型，并引入虚拟质量力修正。核心的半解析技术在于使用基于高斯核函数的加权平均方法计算颗粒的背景流场信息：以颗粒中心为原点建立半径为3倍颗粒直径的扩展域，利用核函数为流场网格赋予0到1的权重，通过加权平均提取背景流体参量，从而实现对大于网格尺寸颗粒的局部流场解析。流体-颗粒相互作用力主要由曳力、浮力和升力组成，其中通过Gidaspow模型计算曳力，并结合Saffman模型计算升力。模型基于开源C++框架CFDEM实现开源计算流体力学软件OpenFOAM与离散元软件LIGGGHTS的耦合求解，流体相方程采用PIMPLE算法求解。样本队列来源方面，高浓度输移验证采用Kidanemariam等设计的明渠几何参数进行测试；低浓度聚类及叠瓦结构验证则参考了Zhang等和Strom等的实验工况及水槽观测数据。
研究结果如下：
验证与单颗粒沉降：通过单颗粒在静水中的沉降模拟验证了球形与非球形颗粒的曳力系数。对于球形颗粒，模拟得到的曳力系数与实验值和理论相关性高度一致。对于非球形颗粒，基于形状因子的模拟结果与实验数据具有可接受的一致性。相同体积的球形颗粒比非球形颗粒能达到更高的最终沉降速度，验证了方法捕捉颗粒穿越气液界面时曳力变化的准确性。
验证与双颗粒沉降：通过沉降过程模拟重现了经典的牵引-接触-翻滚(DKT)现象。球形颗粒清晰地呈现了DKT过程，而非球形颗粒在尾流影响下表现出显著的偏转和快速加速。传统的非解析方法无法解析详细的尾流结构，因此无法捕捉DKT现象。该半解析方法在未显著增加计算成本的前提下成功捕捉了DKT现象。
高浓度条件下的泥沙输移率验证：在明渠高浓度条件下计算了球形和非球形颗粒的输移率。结果表明，非球形颗粒床面具有更大的粗糙度，会诱发更强烈的近床湍流和更高的摩擦速度。模拟的整体输移率与经典的Meyer-Peter实验结果吻合较好。研究指出，在相同水流条件下球形颗粒表现出最高的机动性，若在模拟中仅将泥沙颗粒视为球体，可能会高估推移质输移率。
低浓度条件下的颗粒聚群：在低浓度条件下模拟了少量颗粒的输移和床面演变。当Shields数(Sh)在1.25至2*Sh_cr之间时，床面形成明显的颗粒聚群；当Sh超过2.25*Sh_cr时，聚群开始解体。非解析方法无法在床面发展出类似实验中观测到的条带状聚群结构，颗粒主要聚集在通道侧壁附近。半解析方法成功捕捉了由颗粒尺度尾流屏蔽效应引起的减速和沉积过程，从而自发形成顺下游倾斜的条带状自组织聚群。对于非球形颗粒，模拟成功复现了扁平颗粒相互叠靠、最大投影面向上游倾斜的显著叠瓦状结构，其倾角的概率密度分布(PDF)峰值约为25°，与以往水槽实验和现场观测一致。
总结讨论部分，研究结论如下：本研究采用半解析CFD-DEM耦合方法对非球形推移质颗粒的输移进行了数值模拟，明确考虑了自由表面的影响。与非解析模型相比，该方法有效捕捉了单个颗粒与局部流场之间的详细相互作用，从而更准确地描述了粗颗粒推移质动力学。通过包括单颗粒沉降曳力系数计算、DKT现象模拟、高浓度泥沙输移验证以及床面自组织结构演变在内的一系列验证算例，证明了该方法具有稳健的计算可靠性。模拟结果进一步表明，颗粒尺度的水流扰动在准确复现床面颗粒聚群过程中发挥着不可忽视的作用，这是非解析方法无法表征的。该研究的核心贡献在于能够在不显著增加计算成本的情况下，考虑到影响泥沙输移和沉积的关键颗粒尺度流场特征。当前方法在处理颗粒取向变化引起的受力调整时进行了简化，未来工作可纳入多种形状的颗粒以更贴近自然泥沙形态，该方法可为研究和管理河流生态系统提供有价值的数值支持。
