基于源项的未解LBM–DEM方法用于流体-颗粒流模拟

《Particuology》：A source-term-based unresolved LBM–DEM method for fluid–particle flow simulations

　　赵超|李翔|张家宝|周杜|刘子涵|杨健|周连勇|赵永志•提出了一种新的未解析LBM-DEM方法，该方法能够精确恢复体积平均纳维-斯托克斯方程。•所提方法可本质上消除虚假速度，并实现完全显式计算。•能够可靠模拟固体体积分数在时间和空间上均发生变化的高动态流体-颗粒流。该方法的有效性已通过基准测试案例和流化床实验得到验证。

引言
颗粒流体-固体耦合系统在众多工程应用中十分常见，包括流化床（Alobaid等人，2022；Deivendran等人，2025；Zhao等人，2024）、催化反应器（Chen等人，2020；Wakefield等人，2022；Wang等人，2022）、沉积物输送（Sun与Xiao，2016；Tazaki等人，2023；Zhao，2022）等。准确捕捉分散的固体颗粒与周围流体之间的复杂相互作用，对于理解和优化这些系统至关重要。在各种建模方法中，欧拉-拉格朗日框架，也就是通常所说的CFD–DEM（计算流体动力学-离散元方法）耦合技术（Ma等人，2022；Wang与Shen，2025；Wang等人，2023），因其能够提供详细的颗粒尺度信息而得到了广泛应用。这种方法能够捕捉多种颗粒行为，包括稀薄和密集的颗粒流。根据处理流体-颗粒相互作用的方式，CFD-DEM耦合方法大致可分为两类：解析型（Gallier等人，2014）和未解析型（Zhou等人，2010a）方法。在解析型方法中，流体网格足够细，以至于每个颗粒的界面都能被明确解析，从而可以准确模拟颗粒周围的局部流场，但代价是计算需求显著增加（Washino等人，2023）。相比之下，未解析型方法则利用体积平均纳维-斯托克斯方程（VANSE）（GIDASPOW，1994）以及经验关联式（Difelice，1994；Hoelzer与Sommerfeld，2008；Lu与Gidaspow，2003）来近似处理流体-颗粒相互作用，因此在实际应用中，对于大规模或高密度颗粒系统而言，具有更好的可扩展性。

在传统的CFD-DEM耦合方法中，流体相通常采用有限体积法等隐式方案进行求解（Kuruneru等人，2018；Ma等人，2019；Xue等人，2024；Zhou、Li等人，2026）。尽管这些方法已经十分成熟且被广泛使用，但它们往往涉及复杂的数值实现方式，且计算量很大，尤其是在处理具有复杂几何结构的瞬态流动时。相比之下，格子玻尔兹曼方法（LBM）近年来作为一种颇具吸引力的流体模拟替代方案出现（Sharma等人，2020）。LBM具有两大优势：（1）它是一种天然可并行化的显式算法，能够在瞬态流动模拟中大幅提升计算效率（Chen等人，2014）；（2）它能够更简便地处理复杂边界条件，这在含颗粒或多孔流动系统中尤为重要（Han与Cundall，2013）。许多研究人员通过引入浸没移动边界（IMB）（Noble与Torczynski，1998）方法和浸没边界法（IBM）（Feng与Michaelides，2004）等边界处理技术，开发出了解析型LBM-DEM耦合方法。这些方法已被广泛应用于研究颗粒沉降（Amiri Delouei等人，2022；Qian等人，2023；Zhang等人，2016）、颗粒输送（Ahmadian与Zheng，2024；Sajjadi等人，2017；Zhang等人，2024）、流化现象（Derksen，2022；Ngoma等人，2018；Rosemann等人，2021）等领域。然而，与完全解析型的LBM–DEM相比，未解析型LBM–DEM在方程级一致性、源项处理以及在孔隙率变化剧烈的情况下的稳定性方面仍不够成熟。

两者的关键区别在于颗粒边界的表示方式。在完全解析型模拟中，颗粒表面会被计算网格明确捕获，因此可以利用成熟的边界处理方法来设定复杂的无滑移或移动边界条件。这样一来，就可以结合一些边界处理方法，使用标准LBM来求解流体流动。而未解析型模拟则无法解析单个颗粒表面，流体相必须通过VANSE来建模。不过，标准LBM公式并不能直接得到未解析型欧拉-拉格朗日公式中所需要的体积平均纳维-斯托克斯方程，因此在使用于未解析型LBM–DEM之前，需要额外的建模处理。这一差异在试图将LBM-DEM耦合扩展到未解析型场景时带来了巨大挑战。

一些研究人员为将LBM用于恢复VANSE做出了重要贡献。郭和赵提出了一种专为多孔介质中的流动设计的模型（Guo与Zhao，2002），后来又将其扩展为适用于高压梯度系统的不可压缩模型（Guo等人，2009）。值得注意的是，他们的模型主要针对孔隙率恒定的场景。Song等人（2013）开发了一种LBM–TFM（双流体模型）方法，通过专门设计的碰撞算子避免了隐式计算，同时能够恢复VANSE。Zhang等人（2014）通过将孔隙率直接纳入平衡分布函数，将LBM扩展到了VANSE，并成功模拟了多孔介质中的库埃特流。尽管他们的方法在理论上可行，但后续的一些研究（Blais等人，2015；Bukreev等人，2023；Fu、Hao与Wang，2024；H?cker等人，2018）发现，在存在剧烈空间孔隙率梯度的情况下，会出现虚假速度的问题。为了解决这一问题，Blais等人（2015）修改了碰撞算子，Hocker等人（2018）调整了流动步骤，Bukreeva等人（2023）则修正了零阶矩，以抑制非物理性的速度异常。最近，Fu、Hao和Wang（2024）在郭的不可压缩多孔介质模型的基础上进行了改进，引入了源项以确保其与VANSE的不可压缩形式保持一致。他们通过多个理想固定床基准测试验证了该模型，并证明在模拟高压系统时，其精度优于传统的基于密度的LBM。总体而言，这些宝贵的研究成果已在理想化的流固系统中得到成功应用和验证，为未解析型LBM-DEM耦合提供了坚实的理论基础。

至于将未解析型LBM-DEM应用于实际系统，Xiong等人（2014）将Zhang（2014）的模型与DEM相结合，用于模拟三维鼓泡流化床。他们的结果在床层高度方向上的固体体积分数分布方面，与FVM–DEM模拟结果具有良好的一致性。Rettinger和Rüde（2018）则采用了郭的多孔介质LBM模型（Guo与Zhao，2002），并结合大涡模拟（LES），研究了不同固体体积分数下颗粒组合的沉降行为。尤其值得一提的是，Wang的研究团队在未解析型LBM–DEM的发展方面做出了重要贡献（Fu、Hao与Wang，2024；Fu、Su等人，2024；Fu与Wang，2023；Wang等人，2013）。他们的方法包括先计算每个流体单元内的平均颗粒速度和固体体积分数，然后再利用IMB型方法来评估相互作用力（Wang等人，2013）。随后，这种作用力被用来修正基于经验关联式计算出的阻力，以确保符合牛顿第二定律。在他们最新的研究中，他们进一步改进了这一方法，消除了基于IMB的作用力对松弛时间的依赖，并确保了由IMB计算出的力与基于经验定律计算出的阻力相等（Fu与Wang，2023）。他们还进一步扩展了模型，使其能够用于处理具有复杂内部流场结构的LBM-DEM模拟（Fu、Su等人，2024）。他们的方法已被成功应用于各种气固流动系统，包括鼓泡流化床和循环流化床。

尽管这些研究已经证明了未解析型LBM–DEM的实际价值，但仍有若干问题有待解决。在一些公式中，所恢复的宏观方程并未得到明确阐述。此外，当细胞尺度上的相互作用力是基于平均颗粒信息构建时，其精度可能会受到局部颗粒分布和流动状态的影响（Rettinger与Rüde，2018）。另外，对于那些存在强烈孔隙率梯度以及激烈界面动量交换的实际气固模拟，可能还需要强大的数值处理手段来保证气相LBM求解器的稳定性。这些情况表明，仍需进一步发展一种透明的未解析型LBM–DEM公式，既能恢复目标体积平均方程，又能保持稳健且显式的数值实现方式。

鉴于上述发展，本研究提出了一种基于源项的新的未解析型LBM–DEM耦合方法，用于流体-颗粒流模拟。该方法的主要创新之处在于在格子玻尔兹曼方程中引入独立的质量和动量源项，从而能够在LBE框架内恢复未解析型公式所需的体积平均控制方程。与一些现有的未解析型LBM-DEM和体积平均LBM公式相比，所提方法保留了标准的平衡分布函数，同时也保持了显式的宏观更新流程。这些特点有助于在孔隙率变化的场合抑制虚假速度，同时避免在宏观变量重建过程中采用额外的隐式处理方法。

本文的其余部分结构如下：第2节介绍了所提出的未解析型LBM-DEM耦合方法的发展过程。其中详细推导了用于恢复体积平均纳维-斯托克斯方程的修正格子玻尔兹曼方程（LBE），同时还全面阐述了源项的公式设计及实现方式。该节还介绍了LBM和DEM的基本背景，以及所采用的耦合算法。第3节介绍了所提方法的验证与应用情况，包括在空间孔隙率变化条件下的无虚假速度测试、固定床流动案例，以及一个准二维流化床实验。关于主要用于验证在均匀孔隙率条件下体积平均动量方程恢复情况的孔隙率库埃特流基准测试，则在补充数据中给出。第4节进行了针对特定基准测试的性能研究，包括分析单个时间步内各模块的耗时情况，以及将该方法与针对相同流化床问题的未解析型CFD–DEM基准方法进行对比。第5节总结了本研究内容，并提出了未来工作的展望。