单液滴撞击弯曲壁面的相场动力学模拟

尹宗军，苏 蓉，王清清，陈倩楠，洪 雨

(安徽信息工程学院 机械工程学院，安徽 芜湖 241100)

在过去的几十年中，对湿润现象的动力学分析以及数值模拟一直是流体力学、工程热物理学科前沿的研究主题[1-3]，如喷雾冷却、喷墨打印、液体雾化，以及等离子喷涂。对喷墨打印，精细的摊铺工艺有利于提高打印质量。在内燃机内部，液滴撞击汽缸壁时，油滴如能反弹则不会导致汽缸壁过热。对于喷雾冷却，高速的撞击靶表面或液膜，则有利于快速散热。对于农药喷洒，叶面表面绒毛和腊质与农药的液滴相互作用将直接影响农药能否滞留及其杀虫效能。

在介观、宏观领域，两种互不相容的流体(一般为液气两相)在固体表面形成移动的三相接触线称为移动接触线(MCL,Moving Contact Line)。Huh[4]在无滑移边界条件下，分析了在连续介质力学框架下(宏观足够小、微观足够大)三相接触线处流场的应力和压强具有r-1奇异性，能量耗散率具有对数奇异性。后来，这种现象被称为Huh-Scriven佯谬。三相接触线作为气、液、固界面汇合的极限区域，跨越“原子-纳米-微米-宏观”等尺度，是浸润问题的核心，它不仅影响着流体系统宏观行为规律，同时也是纳米尺度分子间近、长程力和跨尺度关联研究的重要载体。为了消除该奇异性，需要借助触线处的微观信息，引入微观特征尺度用以消除接触线处的奇异性的一些机制被相继提出，例如滑移边界、分子动理论、前驱膜、以及扩散界面等[5-6]。Ren[7]基于唯象理论认为滑移边界模型应该满足三个的本构关系：(1)在体相区，应力张量P和变形速度张量S满足牛顿黏性定律，即P=-pI+2μS，I为单位张量，μ为粘度；(2)在远离移动接触线的固液界面处，壁面剪切应力τs和壁面相对滑移速度us满足Navier滑移关系，即τs=-βus，β为滑移系数；(3)在移动接触线处，Young应力τY与接触线速度uCL存在线性本构关系，即τY=-βCLuCL，βCL为三相摩擦系数。最后，分析得到滑移边界模型体系的能量耗散是有限且非增的。Sibley[8]根据各物理场的阶次分析得到滑移边界模型的压强场仍存在O(log(r))对数奇异性。事实上，当MCL附近液膜厚度达到纳微观尺度时，分子间的相互作用，即分离压力对其流动特性的影响不能被忽略。液体分子在固体表面上替代气体分子发生在前驱膜前缘。由于前驱膜的存在，使得“移动接触线”相对于固体基体产生了相对滑移，前驱膜理论是解决Huh-Scriven佯谬的最重要的途径和机制之一。另一方面，分子动理论认为三相接触线的速度uCL是由Young应力τY导致接触线附近液体分子振动频率κ的改变。Blake[9-10]认为液体分子在固体表面铺展湿润时，液体分子必须克服壁面势阱H后才能取代临近的气体分子。若液滴受到Young应力驱动时，则液滴向前移动速度uCL可以表示为uCL=2λκsinh(λ2τY/2kBT)，其中，kB波尔兹曼常数，λ为分子驱替尺度，T为热力学温度，κ为平衡态振动频率。MCL在连续介质力学框架下的应力奇异性是由于大尺度水动力学和局部小尺度分子动力学的交互作用。基于此，界面扩散模型认为当物质内存在化学位梯度时，由于热运动而导致质点的定向迁移，如果此过程发生在彼此分隔的两相界面上，则称为界面扩散。相比较尖锐界面模型而言，有限扩散界面模型(DIM,Diffusion Interface Model或者称为相场理论)采用了各物理量在有限厚度垂直界面的方向上具有连续变化梯度的思想，从而消除了接触线处的奇异性[11-12]。

结合水平管降膜蒸发器中液滴对弯曲换热管表面产生质热交换的实际应用场所，当前的研究是以毫米尺寸水滴撞击弯曲固体壁面作为研究载体，利用移动接触线奇异问题和相场动力学理论模拟了液滴的铺展动力学，着重分析了四种撞击速度下水滴撞击后的运动形态演化，从而为工业设备中的质热交换强化提供了模拟应用依据。

1 相场动力学的控制方程和能量律

1.1 相场动力学的控制方程

在流体力学中，界面模型应该需要满足三个控制条件：(1)能够捕捉相界面并能够计算表面张力；(2)能够消除MCL附近流场奇异性问题；(3)能够处理好固体表面润湿性质[13]。1958年，Cahn和Hilliard在研究热力学中两种物质之间相互扩散现象时提出一类重要的四阶非线性Cahn-Hilliard扩散方程[14]。若考虑两种互不相容的流体(假设其中一相为气相)在固体表面湿润过程，扩散界面模型提供了一种自然的方式追踪界面拓扑结构的变化，它引入了相场函数φ(r,t)，其在液相内取值为1，在气相内取值为-1，在有限扩散界面层内迅速地由-1变化到1，相界面由φ(r,t)=0确定，如图1所示。两相界面处的相场函数的Cahn-Hilliard扩散方程表示如下[14-15]：

(1)

μ=F′(φ)-ε2▽2φ

(2)

式中，ε为相场界面厚度参量，M(φ)=(1-φ2)表示该方程为具变迁移率Cahn-Hilliard扩散方程，若M(φ)=γλ/ε则表示该方程为具常迁移率Cahn-Hilliard扩散方程,μ表示化学势，γ表示迁移系数，λ表示混合能量密度，F(φ)=0.25(φ2-1)2是单一均质单位体积的自由能，亦被称为双阱位势。此时，系统体相流动的控制方程为Navier-Stokes方程组[15-16]：

▽·(ρu)=0

(3)

η(▽u+(▽u)T)+F+ρg

(4)

其中，ρ=ρ1(1+φ)/2+ρ2(1-φ)/2，η=η1(1+φ)/2+η2(1-φ)/2，在液相中ρ=ρ1，η=η1，而在气相中ρ=ρ2，η=η2，F为表面张力项，它是由于表面张力引起的Laplace附加压力。另外，在尖锐模型中F=σκδSn，式中δS为Dirac算子，为了避免对曲率κ的直接计算，引入毛细应力张量T=-σ▽·(I-n⊗n)δS，(I-n⊗n)为投影张量，并且有等式关系如下：

F=σκδSn= -▽·T

(5)

因此可以将其类推到有限扩散界面模型中，则F=σκεn=σε▽·(I-n⊗n)，另外，n=▽φ/∣▽φ∣。考虑二维平衡状态下，则液-液间的Ginzburg-Landou自由能应该满足最小能量化，即等于系统的表面张力σ：

(6)

根据变分法：

(7)

求得：

(8)

因此平衡态时，相场函数的分布φ(x)为：

(9)

故，表面张力、混合能量密度以及相场界面厚度参量有如下关系：

(10)

图1 液滴的相场动力学模型

1.2 相场动力学的能量律

从能量的观点来看，润湿体系的总的Ginzburg-Landou自由能泛函为[17-20]：

E(φ,▽φ)=

(11)

其中：

(12)

静态平衡时，Fw(1)=σw1，Fw(-1)=σw2，σ为液气界面单位自由能，σw1为固液界面单位自由能，σw2为固气界面单位自由能，θ为动态接触角，θE为静态接触角，并且满足Young关系σw2-σw1=σcosθE。现在考虑液滴铺展时的能量耗散率：

(1)流体混合自由能的耗散率EL

(13)

(2)壁面-流体相互作用自由能的耗散率Ew

=uCL(σw1-σw2+σcosθ)=σuCL(cosθ-cosθE)

(14)

因此总的能量耗散为：

(15)

在这里σ(cosθ-osθE)被称为Young应力，总的能量耗散递减且有界。在壁面处，接触角可以由下式计算：

(16)

2 单液滴撞击弯曲壁面的相场动力学

2.1 建模与参数设置

单液滴撞击弯曲壁面的几何模型如图2(a)所示，该物理模型具有对称特性，因此我们在COMSOL软件中建立了一个轴对称计算域，计算域的大小能够使得液滴自由的铺展。液滴的初始直径D0=6mm，弯曲壁面的初始直径Ds=24mm，初始时液滴的质量中心距离弯曲壁面上顶点的距离h0=6mm。h0过大会造成计算浪费，也不能过小，过小则导致初始时界面网格扭曲严重。除左侧面设置计算条件为轴对称和弯曲壁面为固壁以外，其余边界均为压力出口条件。液滴周围充满气体，以一定碰撞速度V0垂直冲击弯曲壁面，壁面静的静态接触角为θE=90°，不考虑动态接触角的影响。在常温20°下，我们设置水的表面张力为γ=71×10-3N·m-1，粘度μ=2.98×10-3pa·s，密度为ρ=1000kg·m-3。在COMSOL中，在相场模型里面将界面厚度控制参数设置为εpf=3.2×10-4(一般都为最大网格尺寸的五分之一)，迁移调整参数为χ=10，这两个参数与表面张力σ计算相关。计算域采用非结构化网格，选用较细化，由COMSOL自带的物理场控制网格自动生成计算网格，如图2(b)所示，从局部放大图可以看出网格在界面处较细化，没有大的扭曲。为保证计算稳定性，计算时间步长选择1.0×10-3ms，残差由COMSOL自动缩放。

图2 液滴的相场动力学模型：(a)计算模型；(b)网格图

图3 四种不同撞击速度下液滴的运动形态演化：a)0.5m/s;b)1.0m/s;c)1.5m/s;d)2.0m/s

2.2 不同速度下液滴的铺展动力学

在单液滴撞击弯曲壁面铺展过程中，液滴的运动形态主要因撞击速度、液滴的属性以及壁面特性的变化而发生改变，而撞击速度则是引起液滴形态变化、液相传热的最主要因素。本节将就不同速度下液滴的铺展动力学进行研究。

如图3(a)所示，当V0=0.5m/s时，初始时刻两液滴以相同速度相互靠近，我们设置液滴刚开始接触时时间节点为t=0.0ms。随后液滴开始铺展，液滴的初始动能开始转化为粘性耗散和表面能。在t=30.0ms时，在液滴铺展的前端出现两个突起物，由于表面张力的作用，该突起物并没有在重力影响下继续下滑；但到了t=40.0ms时，由于后续的流体继续流入到这个突起物内，使得重力的影响大于表面张力，突起物由于夹断而脱离，在重力作用下继续向下铺展。在t=70ms时，可以看见两个突起物在弯曲壁面的底部汇合而滴落[21-22]。

如图3(b)所示，当V0=1.0m/s时，液滴也在底部汇聚后滴落。其变形过程基本与图3(a)相似。不过，当速度增大后，液滴形变过程加快，克服粘性力和表面张力的共同束缚的能力增大，液滴在惯性作用下铺展的也较远。由于中间薄膜流体流入边缘液体突起物，突起物夹断时刻的润湿位置较远，时间也较早。

如图3(c)所示，当速度增大到V0=1.5m/s时，由于速度较大时，液滴所具有的动能较大。在t=10.0ms时，液滴的润湿位置以及液滴的形态，明显较图3(a)和3(b)有较大改变。在t=30.0ms时，在弯曲壁面底部不远处出现边缘液体突起物，注意这两个突起物并不是重力造成的滑落，而是由于惯性力较大液滴直接铺展到此。此时液滴的初始动能足以克服粘性力和表面张力的共同束缚，两个突起物直接撕裂与液膜的联系向弯曲壁面底部滑落；在t=36.9ms时，我们还可以看见两个更小的液滴突起物，这是由于剩余液膜内的液体在惯性力和重力作用，被表面张力约束所形成。在t=30.0ms时，两个较大液滴突起物在没有汇聚之前，由于惯性力的作用，液滴摆脱了壁面的黏附而发生了飞溅。

如图3(d)所示，当速度增大到V0=2.0m/s时，液滴接触弯曲壁面后立即铺展，随着时间的推移，液滴与弯曲壁面的形成薄液层，液滴沿径向的快速铺展使得被液层拉伸的很薄。由于表面张力的约束，薄液层的两端形成边缘液滴。在t=30.0ms时，薄液层出现了表面不稳定现象，边缘液滴在拉伸过程中逐渐断裂，形成分离液滴[23-25]。这是由于具有更大的动能，克服了表面张力和粘性力的作用而发生撕裂分离。在t=39.9ms时，分离液滴再次脱离壁面而自由运动。由于脱体后的两个分离液滴相向运动，在t=40.0ms时，它们再次相遇融合。

从上述现象可以看出，从能量的角度来分析，液滴的撞击动能越大，液滴具有克服液体粘性做功的能力越大，液滴形变生成的表面越大，铺展越远。当速度较小时，液滴易在重力作用下发生底部滴落，而当速度较大时，液滴容易出现脱离飞溅。当然壁面的直径也会有影响，例如液滴有可能在较低速度时，液滴撞击小直径壁面，动能还没有耗散完就已经达到壁面底部而滴落或飞溅。

3 结论

针对液滴撞击弯曲表面的实际物理现象，当前的研究以相场动力学为研究手段，探讨了相场理论在处理移动接触线奇异型的基本原理，然后对单水液滴撞击弯曲壁面的运动形态演化进行了模拟仿真，讨论了四种不同撞击速度下液滴的铺展特性。仿真结果精确捕捉到了液滴的形态演化，得到以下结论:随着撞击速度的增大，液滴在球面铺展越快，在本文撞击速度0.5m/s-2m/s内,液滴撞击润湿球面后出现了铺展、滴落、分裂、飞溅等现象。