液滴密度和粘度对液滴撞击水平壁面振荡特性的影响

尹宗军，苏 蓉，彭闪闪，陈姗姗，黄俊玲

（安徽信息工程学院 机械工程学院，安徽 芜湖 241100）

液滴撞击过程有着重要的基础研究意义和工程应用场景，例如喷漆沉积、防冰、冷喷涂和燃油喷射等［1-2］.液滴撞击的结果在很大程度上取决于液滴特性（例如液滴密度、表面张力、液滴粘度）、撞击表面的类型（例如固体表面、薄液和深液池）和其他影响条件（例如碰撞速度、液滴大小、碰撞角度、表面粗糙度和壁温）.描述液滴碰撞的主要无量纲数包括：韦伯数We、雷诺数Re以及奥内佐格数Oh［3］.

近年来，液滴碰撞现象已经吸引了许多研究者的兴趣.李晓宇等在平整聚乙烯表面上构建的具有矩形、半圆形和三角形宏观结构的超疏水性纳米线结构，均可造成高速撞击液滴分裂，并且分裂回弹接触时间最低可达到4.8 ms［4］.雷基林等研究了液滴撞击低温壁面的动态行为［5］.解芳芳利用分子动力学模拟方法，研究了超疏水表面上液滴聚并自弹跳及撞击行为［6］.焦丽丽设计了楔形化学非均匀表面用于精确操控撞击液滴定向反弹研究，研究了不同表面温度和撞击条件下液滴撞击三角脊冷表面的动力学行为［7］.王宇航等模拟了不同的液体物性参数下液滴聚并自弹跳行为，认为黏度变化与真实弹跳速度关联性更高，密度变化与真实接触时间关联性更高［8］.刘晓华等对液滴撞击超疏水管状壁面的动态特性展开了研究，讨论了粘性系数和表面张力系数对撞击过程的影响［9］.熊伟分析了基底润湿性、颗粒润湿性、颗粒尺寸分布、环境空气和基底粗糙度等因素对液滴碰撞动力学行为的影响［10］.陈明镜等在绝热条件下开展了液滴在不同惯性力作用下撞击不同润湿性以及不同微观结构表面的动态润湿研究［11］.张莹等制备了亲水区域与疏水区域共存的复合表面，采用高速相机记录了液滴撞击水平与倾斜表面的动力过程［12］.朱阳研究了液滴撞击球面的动力学行为以及液滴与圆球的流固耦合问题，他们通过改变液滴的We数和圆球与液滴的相对大小发现了4种模态，即液滴回弹模态、圆球穿透液滴模态、液滴包裹圆球模态以及液滴脱离球面模态［13］.李玉杰等发现撞击韦伯数、密度比及动力黏性比、圆柱半径等对液滴撞击后沿圆柱内表面的铺展均有一定影响［14］.邹磊等对液滴撞击亲疏水间隔条纹壁面的过程进行了模拟，发现在撞击后容易分裂出卫星液滴［15］.马学虎等综述了核化过程、液滴多尺度润湿、液滴弹跳和液滴多向迁移的固液界面作用机理［16］.李栋等对液滴撞击冷表面的动态行为进行了实验，对比了单液滴撞击普通冷表面以及超疏水冷表面的动力学特性［17］.梁刚涛等探讨了液滴冲击加热壁面时中间射流及宝塔状气泡的形成机理［18］.毕菲菲观测了蒸馏水、无水乙醇以及甘油液滴撞击水平、倾斜以及弯曲壁面的流动规律［19］.

本文中，VOF相追踪模拟被用于研究液滴在光滑表面上前进/后退的振荡特性.动态自适应网格也被用于自适应跟踪液滴铺展时的微小流动结构，这有利于提高气体、下落液滴和固体之间相互作用的计算精度.模拟结果可靠性通过与前人实验报告进行比对验证.通过改变液滴密度ρ和粘度μ来研究这2因素对液滴振荡特性（用铺展因子β=D/D0随时间的变化来表征，D为当前湿润长度，D0为液滴直径）的影响，并通过最大铺展因子βmax、平均振幅比χ以及平均振荡周期τ来分析液滴振荡特性的一般规律.

1 控制方程与参数设置

1.1 流体流动的控制方程

由于液滴撞击速度较低，气相（ρV和mV分别表示空气的密度和粘度）和液相（ρL和mL表示液滴的密度和粘度）均假定为不可压缩、等温、不互溶地流动.在本研究中，采用VOF相界面追踪法描述液滴撞击壁面过程的气液界面.通过引进相体积分数这一变量c，不同的流体组分共用一个相体积分数方程，实现对每一个计算单元相界面的追踪；在气相中，c=0，在液相中，c=1，在气液界面接口处，c的值介于0和1之间.在每一个计算单元中，密度ρ可以表示为ρ=cρL+（1-c）ρV，粘度为μ=cμL+（1-c）μV.因此，液滴撞击壁面的流动控制方程包括连续性方程、带表面张力的动量方程以及基于VOF相界面追踪的相分数方程［20］：

式中，u≡（u，v）为速度矢量，p为压强，D定义为Dij≡（∂iuj+∂jui）/2，γLV为表面张力系数，g为重力加速度，δ为迪拉克算子，κ表示界面的平均曲率，n表示从液相流出的界面的单位法线.根据连续表面张力（CSF）近似模型，δn≈∇c，κ≈∇（=∇c/‖∇c‖）.采用具有二阶精度的交错时间离散法对上式（1）至（3）中体积分数/密度和压力进行数值离散，得到：

式中，ρi、ui、μi、Di以及ci为对应时间节点i（i可取n、n+1/2 和n+1）的值.

若指定壁面接触角为θE时，则壁面单元处界面单位向量的法向量可以表示为：

其中和分别表示与壁面垂直和相切的单位向量.本文通过GERRIS软件（一个用于求解描述流体流动的偏微分方程的开源软件）对上述离散方程进行了数值求解.

1.2 建模与参数设置

图1显示了液滴撞击壁面过程中的初始液滴和自适应网格图.图1（a）中计算域的右侧、左侧以及上侧边界为压力出口边界，底部为滑移边界条件（满足u=u0+λ∂u/∂y，λ取0.1 mm）.计算域是一个二维模型Lx×Ly=2L×L，其中L为8 mm，计算域的大小足以保证流体可以自由流动.在模拟开始时，设定两相中的初始体积分数：空气体积分数cV=0，液体体积分数cL=1，设定流体域的初始压力p为0，设定初始时间对液滴施加沿 -y方向撞击速度为V0.在计算中，水和空气的密度分别设定为ρL=9.98×102kg·m-3、ρV=1.2 kg·m-3，水和空气的粘度分别设定为μL=1.003×10-3N·s·m-3、μV=1.8×10-5N·s·m-3.气 液 两 相 的 表 面 张 力 设 置 为γLV=7.5×10-2N·m-1.

为精确捕捉铺展时的微小流动，自适应网格细化被应用去提高气液界面的网格分辨率.在液滴界面和液膜界面处局部细化网格，其余地方可采用合适的粗化网格：在液滴界面和液膜界面处局部细化网格，液滴界面局部区域内加密等级选择为10；除液滴界面外的计算域加密等级选择为6.计算单元相应尺寸为L/210和L/26（即 7.812 5×10-3mm 和 1.25×10-1mm）.GERRIS 软件利用动态网格自适应手段在相界面处自动加密可以减少不必要的计算.图1（b）显示了液滴撞击水平壁面后铺展的自适应笛卡尔网格图.

图1 液滴的铺展动力学模型

1.3 验证

如图2所示，Gunjal等介绍了水滴（直径D0=2.5 mm，冲击速度V0=0.3 m·s-1）撞击聚四氟乙烯表面（静态接触角θE=110°）时，水滴高度随时间的变化，其被用来与我们的模拟进行定量比较［21］.需要注意的是，Gunjal等文献中的实验数据（即绿线，h*=h′/D0，h′表示液滴瞬时形状的最大高度）与他们的仿真模拟结果（即红线，h*=h/D0，h表示液滴中心处的高度）之间的差异可能是由于测量方法的差异［21］.在不包括第一个振荡周期的100.0 ms内，在他的实验中液滴铺展和回缩引起的振荡（超过465.0 ms）得到的平均振幅比χ为1.032，平均振荡周期τ为18.12 ms.而本文计算中平均振幅比和平均振荡周期分别为1.051和18.33 ms，这与Gunjal等的工作具有一致性.此外，VOF仿真模拟中第一个振荡循环期间的最大无量纲高度为1.48，这也与Gunjal等的实验值（即1.48）一致.

图2 Gunjal的实验结果与当前VOF仿真模拟结果验证比较

2 结果与讨论

2.1 液滴密度ρL对液滴振荡的影响

液滴密度与液滴惯性力密切相关.在进行VOF仿真模拟实验时，选择直径D0=1.6 mm的类水液滴来进行撞击水平光滑壁面（θE=110°），选择只改变液滴密度ρL分别为 700、800、900、998、1 100 kg·m−3（这些密度值的分别对应汽油、煤油、矿物油、水和乙二醇的密度值），其余物理属性与水一样，来模拟液滴密度ρL对液滴振荡的影响.

图3（a）显示了液滴密度ρL=700 kg·m−3对液滴振荡的影响.时间序列t=0.0 ms设置为液滴刚与固体表面接触时刻；随后，液滴发生铺展，当液滴高度连续降低时，液滴的下部形成一个盘状结构（t=2.0 ms），并随时间变化而变大；在t=7.4 ms，液滴达到最大润湿半径，液滴外围边缘开始向中心收缩.在t=24.2 ms时，液滴上升到某一高度；惯性流动不足以将液滴从固体中升起，液滴开始回落.在那之后，滴水会经历一系列的周期振荡，最后静止在固体表面上（t=400.0 ms）.图3（b）~（e）显示了其余4种不同液滴密度的液滴运动形态演化，它们的形态变化几乎相似，但是所对应的时间有所不同.

图3 5种不同密度ρL 的液滴运动形态演化

图4展示了不同液滴密度值的类水滴铺展因子β时间变化趋势.如图4（a）所示，在第1个铺展回缩循环下，所有液滴密度值的扩展因子β变化趋势相似，但峰值有所不同.对于较大的液滴密度（例如，ρL=1 100 kg·m−3），最 大 扩展 因 子βmax显 著 较 高.图4（b）为 100 ms 内 的 类 水滴铺展因子β时间变化趋势图.其振荡参数如平均振幅比χ分别为1.052、1.053、1.055、1.055、1.061，平均振荡周期τ分别为14.1、15.3、16.3、17.3、18.3 ms.可以发现，密度值从 700 kg·m−3变化至 1 100 kg·m−3时，平均振幅比χ和平均振荡周期τ随着液滴密度增加而增加，这是因为由于惯性能作用，较高液滴密度（例如ρL=1 100 kg·m−3）会使液滴有较高初始动能.同时可发现，在其他条件相同情况下，较小液滴密度在早期阶段有助于加速达到平衡；但反过来，较大液滴密度会使液滴运动受到额外阻尼，从而加速振荡结束时的振荡衰减.

2.2 液滴粘度μL对液滴振荡的影响

选择只改变液滴的粘度值μL分别为 0.6、0.8、1.003、2.5、4.0 N·s·m−3，其余物理属性与水一样，来探究液滴密度μL对液滴振荡的影响.

图5（a）显示了液滴粘度μL为 0.6 N·s·m−3对液滴振荡的影响.此时，当液滴撞击水平壁面时，液滴发生扩散并创造新的表面积，液滴动能的一部分也消散在克服粘性力上（t=2.0 ms），直到液滴动能耗散完（t=9.6 ms），作用于界面的表面张力阻止了液滴下降并开始回缩液滴.由于惯性流动，液滴的动能和高度在回缩过程中逐渐增加，直到所有动能转化为势能（t=21.2 ms）为止.如果惯性发展过程中惯性流动大到足以将液滴从固体中升起，液滴会反弹；否则，液滴升到最大高度后会再次下落.这样的循环铺展和回缩将持续相当长的一段时间（t=400.0 ms）.由于表面张力和壁面润湿性的不同，液滴振荡行为可能表现出不同的变化.图5（b）~（e）显示了类似的形态变化，区别在于对应的时间有所不同.

图5 5种不同粘度μL的液滴运动形态演化

如图6所示，5种不同液滴粘度μL下液滴的铺展因子变化趋势图相似.图6（a）显示最大铺展因子βmax随液滴粘度的降低而增加.这是由于粘度值高（如μL=4.0 N·s·m−3）会 降低液滴变形速率，阻滞液滴运 动.图6（b）表 明，粘 度对液滴振荡有显著影响.平均振幅比χ分别为1.066、1.060、1.057、1.090、1.094，平均振荡周期τ分别为 18.1、18.0、17.3、15.3、15.1 ms.当液滴粘度值从 0.6 N·s·m−3变化至 4.0 N·s·m−3，可以观察到液滴粘度越小（如μL=0.6 N·s·m−3），振幅越大，相应的平均振荡周期越长.粘度值为 4.0 N·s·m−3，振幅小，平均振荡周期小，振荡速度明显加快.总之，图6表明液滴具有相对较高的粘度（例如μL=4.0 N·s·m−3）往往会抑制液滴在壁面上振荡，尤其会导致振幅衰减.

图6 不同粘度μL的液滴撞击过程中铺展因子的变化曲线图

3 结论

本文对液滴垂直撞击在水平壁面上的振荡行为进行了二维VOF数值模拟，采用动态自适应网格细化技术捕捉液滴铺展的微小流动.借助于最大铺展因子βmax，平均振幅比χ以及平均振荡周期τ，定量描述了液滴密度和粘度对液滴振荡特性的影响，研究表明当液滴密度值从 700 kg·m−3变化至 1 100 kg·m−3时，最大扩展因子βmax逐渐增大，平均振幅比χ和平均振荡周期τ也逐渐增大；液滴密度越大，惯性能越大，越会促进液滴的扩散.另一方面，较大的液滴粘度会阻碍液滴变形，将导致平均振荡周期变小和振幅衰减加快.