基于MDPD方法的液滴动态润湿研究

孟欣,刘汉涛,李海桥

(中北大学能源动力工程学院，山西 太原 030051)

表面微观结构对液滴动态润湿性能的影响，具有重要的科学意义和工程应用价值，液滴润湿固体壁面这一过程它不仅存在于自然界，还广泛存在于日常生活中[1]，如荷叶表面具有疏水性和自清洁性[2]，主要是其表面的微米结构凸起与水之间相互作用决定的[3]。在工业领域，表面微结构的这种特殊属性能够引领仿生材料的开发与应用[4]。深入揭示固壁表面微结构的浸润特性对推动新型材料的研发，以及特殊功能器件的开发具有非常重要的指导意义。因此，表面微结构上的浸润问题一直是当前科学研究的热点与难点之一。

针对表面微结构液滴浸润，很多学者进行了研究。Zhao等[5]研究了浸没微柱表面上液滴结合的自发润湿转变，从Wenzel态转变到Cassie-Baxter态，探讨了一些由低能态向高能态转变的方法；Zhao等[6]还进行了基于多体耗散粒子动力学(many-body dissipative particle dynamics,MDPD)的液滴浸润研究，对MDPD方法的理论及其应用做了详细介绍，如在化学非均匀的微结构表面上，对亲水条纹和疏水条纹相间情况下液滴的碰撞和扩散动力学行为进行研究，用幂次定律验证模拟结果，分析了液滴的撞击和扩散全过程；Kwon等[7]在实验和模拟的基础上，通过设置不同表面微结构验证Cassie-Baxter方程[8]，证明Cassie-Baxter方程可以说明接触线对表观接触角的影响。

研究表明，表面微观结构对液滴动态润湿性能能够产生显著影响。目前，在宏观和微观尺度下的研究较多，介观尺度下的研究较少。MDPD方法作为典型的介观尺度方法被应用于液滴动态浸润研究[9]，研究内容涉及液滴运动、接触角滞后、液滴扩散，以及液滴在不同表面的碰撞、液滴润湿状态的转变等。本文利用MDPD方法构建不同的微结构粗糙表面和液滴，与文献[10]的实验结果和文献[11]的模拟结果进行对比验证，在此基础上统计液滴平衡接触角和铺展变化，通过固液作用力对微表面结构上的浸润问题进行深入研究，为控制液滴亲、疏水性提供介观尺度的理论基础。

1 MDPD方法

MDPD方法是一种基于随机粒子的方法，该方法中单个粒子i代表实际的大量原子，这一粗粒化过程能够有效获得所需物理学的主要特征，在提高计算效率的同时降低了计算成本[12]。这里的每个MDPD粒子i的运动都要受牛顿第二运动定律的约束，运动方程为：

(1)

式中：ri和vi分别为单个粒子i的位置矢量和速度矢量；m为MDPD粒子的质量，所有MDPD粒子的质量相等；Fi为施加在粒子i上的合力。合力Fi包括粒子间内部的相互作用力和外部对粒子的作用力，可以表示为：

(2)

式中：Fij为第i个粒子和第j个粒子之间相互作用力的总和；Fi,ext为外部对第i个粒子施加力的总和。

(3)

(4)

(5)

(6)

(7)

(8)

式中：γ和σ分别为耗散力和随机力强度的系数；ωD(rij)和ωR(rij)为耗散力和随机力的权重函数；vij=vi-vj，为颗粒i对于颗粒j的相对速度；θij是一个具有高斯分布和单位方差的随机变量。为了维持恒定的比动能，MDPD与DPD系统有着相似的涨落耗散定理。

(9)

式中：kBT为系统的能量，通常取值为单位值1。其中kB为玻尔兹曼常量，T为系统热力学温度。

另外，在本文中固-液粒子的相互作用尤为重要，适当的作用力可以有效地防止液滴粒子穿透固体壁面，也能明显改变液滴的浸润形态，结合SPH钟形权函数的吸引力和排斥力，构造一种短程排斥远程吸引的固液作用函数，类似于式(4)的形式。

2 模型验证

为了验证本文模型的准确性和可行性，与文献[10]的实验结果和文献[11]的模拟结果进行对比。文献[10]的实验中液体为蒸馏水，初始液滴直径为50.0 μm，稳定后的接触角为90.0°。将实验数据和数值计算结果进行简单的无量纲化处理，包括铺展因子(d/D)、无量纲高度(h/D)与无量纲时间(tV/D)，d和h分别为液滴的瞬时铺展直径和高度，D为液滴初始直径，V为反映液滴振荡的特征速度，其定义为：

(10)

式中：ρ，Γ分别为液滴的密度和表面张力。

文献[11]利用MDPD计算液滴动态浸润问题时，液滴的直径为384.0 nm，奥内佐格数(Ohnesorge number)为0.19。本文所采用的MDPD方法与文献[11]构建的方法之间的最大差异体现在，本文的粒子间采用近程排斥、远程吸引的相互作用力，壁面设置为光滑，而文献[11]采用调整保守力系数来实现固壁间相互作用。针对液滴在固体壁面上润湿的模型，本文与文献[10]、[11]结果间的对比如图1所示，从铺展因子和无量纲高度随时间的变化可以看出，本文所采用的方法更容易达到稳定状态。为了进一步证明本文所采用方法的有效性，利用最小二乘法对两组计算结果与实验结果之间的误差进行了评估，接触角为90.0°时的铺展因子和无量纲高度的误差均很小，证明本文所使用的方法精确度更高。

图1 模型验证

3 结果与讨论

3.1 设置初始化模型

本文以液滴动态浸润过程为研究对象，研究壁面微结构参数对固体表面润湿性的影响。本文使用的液滴和固体壁面模型是简化模型，液滴对应真实水滴的流体参数见表1，固体壁面数密度设置为液滴密度的2倍，固体壁面上的微柱阵列是微结构粗糙表面的简化模型。

表1 模拟液滴的参数对照

液滴半径为4(单位对照见表1)，且各模型中作用力保持不变，通过接触角来观测微柱高度对浸润过程的影响。第一组设置微柱高度从0～7，其具体的物理模型如图2所示；第二组微柱高1，微柱的俯视图为正方形，每个小单元的边长为10，本文的面积分数是每个微柱俯视图的面积占小单元面积的比值，微柱面积分数分别设置为0.04,0.16,0.36,0.64；第三组微柱高2，微柱面积分数设置同第二组。

图2 液滴浸润初始化模型

对于液滴浸润，用无量纲数奥内佐格数Oh表示黏性力与表面张力之间相互作用的比值：

(11)

式中：μ为液滴的黏度。本文的液滴半径为4，奥内佐格数为0.263。以MDPD为单位，并基于本文模拟的流体物理特性，用奥内佐格数将模拟的MDPD系统与以水为流体的真实系统相匹配，计算得到水滴的直径为200.952 nm。本文模拟系统的单位1相当于25.12 nm。在构建的固液模型中，计算区域为30×30×30(单位对照见表1)。

3.2 计算结果分析

第一组的固液作用力系数，吸引项为-40，排斥项为20，在液滴铺展过程中，不同微柱高度上的接触角随时间的变化曲线如图3(a)所示，由图可知，液滴开始铺展后，接触角一直在变化，在表面张力作用下，液滴的外轮廓变得光滑。在0.1～0.5时，液滴铺展过程逐步减缓，表面张力趋于稳定，在耗散力、随机力和保守力的共同作用下，液滴形状出现明显的波动，对应的接触角在此时也出现了明显波动；0.5之后，在耗散力的作用下，液滴内部趋于稳定，液滴的形态和接触角也到达平衡状态。在光滑表面上，液滴接触角稳定在27.48°，固壁表现为亲水性。当存在微柱阵列时，液滴铺展受限，液滴接触角在100°～140°，固壁表现为疏水性。不同微柱高度下铺展因子变化曲线如图3(b)所示，由图可知，随着面积分数的增大，固体的表面积增大，固体表面与液体之间的相互作用增强，阻碍了液滴的铺展，当微柱高度增加到一定量时，此处的阈值为2，固体表面与液体间的接触面积基本不变，液滴的铺展因子基本保持不变。

图3 不同微柱高度下的计算结果

第二组的接触角和铺展因子变化曲线如图4所示，液滴平衡后的浸润状态如图5所示。在0.763之后，液滴内部趋于稳定，固壁都表现为亲水性。面积分数为0.04和0.16的条件下，液滴形成Wenzel态，即液滴底面与固壁完全接触；面积分数为0.36和0.64条件下，液滴部分粒子进入微柱之间，形成不完全润湿。面积分数为0.04时接触角最小，稳定在40.20°；面积分数增加后，对浸润状态的影响较小，液滴接触角在80°～86°，面积分数为0.64时铺展更为充分，铺展因子稳定在1.64。研究表明，微柱高度较低时，表面分数(宽度)对液滴浸润状态影响有限，固壁表现为亲水性。

图4 微柱高1的计算结果

图5 微柱高1不同面积分数下稳定后的液滴浸润图

第三组的接触角和铺展因子变化曲线如图6所示，平衡后的浸润状态如图7所示。在0.763之后，液滴内部趋于稳定，宽度对浸润产生了较大影响。面积分数为0.04和0.16的条件下，接触角稳定在137.40°、149.00°，固壁表现为疏水性，液滴外观呈球形，液滴粒子进入微柱之间，发生明显的“钉扎”效应。面积分数为0.36和0.64的条件下，接触角稳定在96.04°、67.94°，且液滴外观呈半球形，液滴粒子均处于微柱之上。研究表明，当微柱高度增加到一定量时(此处阈值为2)，不同的表面分数会对固壁表面润湿性产生显著影响，即随着面积分数的增加，固壁表面从Wenzel态的疏水界面逐渐转变为Cassie-Baxter态的亲水界面。

图6 微柱高2的计算结果

图7 微柱高2不同面积分数下稳定后的液滴浸润图

4 结论

本文基于MDPD方法研究表面微结构参数对液滴动态浸润的影响，结果表明：

1)微柱高度对于液滴浸润状态的影响并不是一直存在，光滑壁面条件下液滴平衡后的接触角为27.48°，表现为亲水性。当微柱高度大于2时，固液间的相互作用不会随着微柱高度的增加而发生改变，固壁表现为疏水性，液滴不完全铺展。

2)微柱高度为1的情况下，微柱的面积分数对浸润状态影响不大，液滴粒子会进入微柱间形成Wenzel态，固壁表现为亲水性，液滴普遍呈半球形；微柱高度为2的情况下，微柱的面积分数可以对浸润状态发生显著影响，随着面积分数的增大，固壁表面由疏水变为亲水，液滴粒子不再进入微柱间，由Wenzel态转变为Cassie-Baxter态。