冰内球形颗粒周围预融薄膜动力学模拟的数值方法

《International Journal of Multiphase Flow》：A numerical approach to simulate the dynamics of premelted film around a spherical particle in ice

　　

研究人员提出了一种策略来解决热再凝固（regelation）问题模拟中的困难，在一个统一的框架内求解了流动动力学、传热和相变，并将结果与球形颗粒迁移的润滑理论进行了验证，该方法可直接应用于复杂颗粒和颗粒团簇。

论文解读：《冰内球形颗粒周围预融薄膜动力学模拟的数值方法》

研究背景与意义

多孔介质中的水冻结会引发一系列微观与宏观耦合的动力学现象，其中最具代表性的便是土壤膨胀、冰分凝以及由此产生的地下水迁移（即冻吸作用），统称为“冻胀（frost heave）”。冻胀不仅会破坏土木基础设施的稳定性，还参与了岩石破裂和图案地面形成等地质过程。研究表明，冻胀的主要驱动力并非水结冰时的体积膨胀，而是源于预融（premelting）动力学。即使在远低于体相熔点的温度下，固体颗粒与冰基体界面处也会形成一层极薄的液态水膜（预融薄膜）。在温度梯度下，较暖侧增强融化而较冷侧重新冻结，驱动颗粒相对于周围冰基质发生迁移，这一过程被称为“热再凝固（thermal regelation）”。然而，预融薄膜的厚度通常仅为几十纳米，而典型土壤颗粒的尺寸则从几十微米到毫米级不等。这种跨尺度的巨大差异，使得在传统数值模拟中同时解析薄膜内的流动以及相间的流固耦合（Fluid-Structure Interaction, FSI）变得极具挑战。为了弥合这一尺度鸿沟，研究人员开发了一种新型计算模型，旨在单一模拟中复现所有热再凝固现象，为理解冻胀的微观机制提供工具。

主要关键技术方法

为实现上述目标，研究人员采用了以下几个关键技术方法：

1. 1.
   参数缩放策略：研究人员并未直接模拟真实的纳米级薄膜，而是通过调整材料属性中的自由参数 λ，使预融薄膜厚度达到颗粒直径的0.05至0.20倍。这种理想化的流体属性使得在保留微观物理现象的前提下，能够显式解析整个系统及薄膜内的流动。
2. 2.
   相场法（Phase-Field Method）：采用相场方程来捕捉融化-冻结界面的位置，并结合不可压缩Navier-Stokes方程和能量方程，在一个统一框架内处理了涉及移动边界、液相流动和多相传热的强耦合问题。
3. 3.
   浸没边界法（Immersed Boundary Method, IBM）：利用基于移动最小二乘（moving-least-squares）实现的IBM来处理球形颗粒的存在。该方法在动量方程中引入体积力，确保了颗粒与液体界面的速度连续性，并强制执行颗粒内部的刚体速度分布，从而高效处理了流固耦合。
4. 4.
   无量纲化分析：研究以颗粒直径 d、特征温度差 d?T、特征速度 gαd2?T和密度 ρ作为尺度基准，将控制方程无量纲化，重点关注雷诺数（Reynolds number, Re）、佩克莱特数（Peclet number, Pe）和斯特藩数（Stefan number, Ste）等关键无量纲参数。

研究结果

3.1 流动动力学与极限行为

研究发现，在预融薄膜厚度趋近于零的极限下，薄膜内的流动处于润滑（lubrication）状态。此时，粘性合力（Fμ）相较于流体动力压力合力（Fp）可以忽略不计，颗粒迁移速度主要由分离压（disjoining pressure）与流体动力压力的平衡决定。此外，研究验证了重力效应在地球条件下完全可以忽略，因为1 K/m的温度梯度产生的热分子驱动力相当于约700 m/s2的重力加速度。在有限薄膜厚度下，数值解与解析润滑解表现出良好的一致性，但由于有限厚度带来的粘性应力贡献，两者的偏差也在预期之内。

3.2 长期特征与能量学

随着颗粒向更温暖的区域迁移，赤道处的预融薄膜厚度会增加，分离压减小，导致颗粒的平移速度随时间缓慢增加，直至到达底部边界。通过对全局能量平衡的分析发现，穿过边界的热交换主要通过潜热（latent heat）进行冰的融化，并伴随界面能的变化，其余的能量贡献则通过粘性耗散掉。这表明该系统是一个由热梯度驱动的、伴随着相变的动态非平衡过程。

3.3 变材料属性的影响

当考虑不同相和材料的真实热力学属性（如热导率、比热容和密度的空间变化）时，系统的动力学受到显著影响。以粘土矿物和高石英颗粒为例，较高的热导率和较低的容积热容会削弱颗粒内部的局部温度梯度。根据理论公式，这种梯度的降低直接导致颗粒平移速度的下降。这一发现为解释自然界中不同土质（如导热性好的石英砂与导热性差的粘土）抗冻胀性能的差异提供了微观物理解释。

讨论与结论

总结讨论部分

研究人员指出，尽管当前的模型尚未包含冰透镜生长和大尺度水分传输等冻胀的关键过程，但其揭示了热梯度下预融薄膜介导的水分迁移微观物理机制。颗粒平移速度的降低与冻胀过程中驱动水分向冻结锋面迁移的低温吸力减弱直接相关。此外，不同材料的热物理性质差异是导致石英砂抗冻胀能力优于粘土的重要原因之一。未来的工作重点将扩展至颗粒团簇的研究，因为在相邻颗粒形成的孔隙空间中，将涌现出孔隙尺寸、水分传输与预融层厚度相互作用的新动力学，从而建立与冻胀现象更直接的联系。

研究结论翻译

本文提出并验证了一种旨在模拟热再凝固现象的新型计算模型。通过设定理想的流体属性以获得有限的预融薄膜厚度 δ，同时保留与无穷小厚度相关的分离压，成功规避了固体颗粒尺寸与预融液层厚度之间的跨尺度难题。结果表明，该模型的行为与理论一致，即预融层内的流动根据润滑极限演化，并通过平衡分离压和流体动力压力达到平衡。研究人员证实，在 δ趋近零且固体颗粒尺寸较小的极限下，重力和流体热膨胀的影响相较于分离压均可忽略不计。相反，各相和材料的不同热学性质会显著改变系统动力学；特别是，具有较大热导率和较小容积热容量的固体夹杂物会导致其平移速度相较于均匀属性参考情况有所降低。该模型通过浸没边界法处理流固耦合，使其能够像模拟球形颗粒一样轻松地模拟复杂几何形状，为后续研究真实土颗粒乃至颗粒团簇的再凝固现象奠定了基础。