竖直振动激励下半球形液滴界面失稳特性分析1)

黎一锴 朱 铭 席 儒 王东方† 杨子明 吴 坤,2)

*(北京理工大学机械与车辆学院,北京 100081)

† (北京理工大学重庆创新中心,重庆 401120)

** (中国科学院力学研究所高温气体动力学国家重点实验室,北京 100190)

†† (昆明理工大学云南省内燃机重点实验室,昆明 650500)

引言

竖直振动板上的附着液滴动力学广泛出现在工业和科学领域,如蒸发结晶、超声雾化、风力涡轮机叶片除冰和内燃机二次雾化等[1-6].当振动幅度超过一定阈值时,首先在液滴表面观察到对称的规则图案[7].随着振幅的进一步增加,表面波变形加剧,这会产生尖锐的液线,并从液体射流的表面产生二次液滴[8].

研究具有不同物理性质的液体、不同接触壁面条件和不同振动模式下的液滴的表面变形是液体动力学的一个焦点.James 等[9]观察到液滴出现轴对称驻波(纬向波)、静止然后缓慢旋转的方位角波(经向波)、随时间变化的弹坑和液线,以及在逐渐增加的激励振幅下的微小二次液滴的快速喷溅过程.Chang 等[10]将水滴放置在水平壁面,并通过施加一系列频率的外部激励来观察水滴的共振行为.液滴的表面形态分为3 种模式: 纬向模式、经向模式和镶嵌模式.通过引入球谐函数,根据液滴表面的几何形态图案对应识别每一个模态.Singla 等[11]将水银放置在球形玻璃上,并采用不同的振幅和频率对其施加垂直扰动,同时用高速摄像机记录水银液滴的振动.实验观察到,水银液滴受到垂直扰动时液滴的振动频率会达到扰动频率的一半,不同的振动液滴的频率对应于不同的表面模态.王凯宇等[12]研究发现超疏水表面上液滴的体积对共振振幅、模式区间、共振频率等振动特性影响较大.Noblin 等[13]研究了有重力条件下聚苯乙烯振动基板上液滴的振动特性.通过调整频率和振幅,观察到液滴的两种轴对称振荡的状态,一种为接触线保持固定,主要出现在低振幅条件下,另一种接触线处于振荡状态,主要出现在临界振幅之上.Shao 等[14]通过实验研究了圆柱形容器中机械振动激发下的空气水界面表面波,研究表明界面与容器接触角对界面表面波波形有着很大的影响.

在垂直振动的液层表面上触发周期性表面波的机制通常归因于Faraday 不稳定性,这是Faraday[15]在研究垂直振动作用下液层表面的响应时提出的,他观察到,表面驻波以驱动频率一半的频率振荡.Kumar 等[16]将液体黏性纳入了考虑范畴,修正了Mathieu 方程并给出了黏性液体表面不稳定性参考表,并且与实验结果[17]吻合良好.Bronfort 等[18]研究了系统在垂直摇晃下气泡-液体界面的Faraday波,研究表明存在一个明确的加速度阈值,超过阈值气泡-液体界面出现次谐波.Chen[19]采用二维数值模拟研究了不同模式下Faraday 波的非线性动力学特性.Li 等[20]通过三维模拟研究了Faraday 不稳定下液层表面液线形成及其随后断裂相关的基本非线性动力学.Dinesh 等[21]指出,在机械和静电引起的Faraday 不稳定性情况下,利用无黏理论都可以预测黏性流体共振时的模态响应.在接触线迟滞现象对Faraday 不稳定性的研究中,Yuan 等[22]发现迟滞会使固有频率增大,并通过计算得到了接触角范围与接触线位置之间的线性关系.Bestehorn 等[23]用相场方法对Faraday 不稳定性进行了研究,进一步分析得到不同流体黏度下的非线性解.Dong 等[24]研究了扁平液滴的上表面在垂直激励下由表面参数不稳定性引起的星形振荡,并推导了一个结合平面Faraday波模态和液滴方位振荡模态的三维振荡模型.

在径向周期正弦外加激励作用下,球形液滴表面也会出现Faraday 不稳定性.Adou 等[17]在受径向激励的球形Faraday 不稳定性分析中,通过将Floquet方法[25]应用于球坐标系,解决了线性稳定性问题.姚慕伟等[26]基于线性小扰动理论研究了受径向振荡体积力的黏弹性液滴表面波的不稳定性.康宁[27]基于二维轴对称数值仿真模型,并结合水平集与流体体积耦合算法(CLSVOF 方法)[28],对液滴在 Faraday不稳定性主导下的二次雾化过程及雾化机理做了深入研究.

现有的研究主要集中在Faraday 不稳定性引起的径向力作用下液滴的界面不稳定性[9].当液滴受到垂直加速度时,表面也会出现次谐波.很少有研究探讨竖直力产生的次谐波与Faraday 不稳定性之间的关系.以前的许多研究工作都致力于研究液体层而不是液滴的不稳定性.关于Faraday 不稳定性下球形液滴的表面波与竖直激励下的表面波之间的联系的文献很少.此外,现有的模拟研究主要通过二维轴对称数值模型研究液滴失稳产生的表面波的物理特性.然而,当发生非轴对称经向表面波时,二维轴对称数值模型不适用.因此,有必要通过三维数值模型进一步研究不同工况下液滴表面波的物理特性和机理.

本工作采用液滴振荡实验结合三维数值仿真,研究了振动平板上液滴的振荡特性,主要讨论了表面经向波和纬向波的演变过程.并将竖直作用力下的表面失稳特性与径向作用力下由Faraday 不稳定性引起的表面失稳进行了比较,阐明了液滴在竖直加速度下发生经向失稳机理.丰富了液滴表面失稳形式的理论研究.

1 研究手段

1.1 实验设置

实验装置如图1 所示,通过该装置可以独立控制振动平板的激励振幅和频率.液滴发生器是一个尾端连接步进电机的注射器,注射器出口靠近平板.开启步进电机可产生一定体积的液体(本文选用去离子水)使其聚集在平板中心位置,液滴直径为7.5 ±0.02 mm.相关液滴物理特性如表1 所示,实验中测量得到的液滴在铝基板上的静态接触角均值为87°.

表1 水滴在室温、标准大气压下的物理特性Table 1 Physical properties of water at atmospheric temperature and pressure

图1 液滴振动实验系统Fig.1 Schematic diagram of the experimental setup

UTG9005C 型函数信号发生器串联SA-PA010型功率放大器,用于产生一定振幅和频率的正弦电压信号.铝平板通过M5 螺纹连接到JZ-002 型激振器的顶端,并以与激励电压信号相同的频率在竖直方向上振动.Ds1054 示波器清晰地显示出频率和电压.1200 W 的镝灯为实验光源,均光片放置在观测液滴与光源中间使光线均匀分布.实验中使用PhantomV7.3 高速相机,使用焦距为180 mm 的Tamron 微距镜头.拍摄帧率为11 000 fps,曝光时间为88 μs,拍摄视窗为512×512 像素.

图2 显示了在不同频率下,平板振幅和激励电压之间的线性相关性.因此,可以通过分别控制驱动频率和平板振幅来观察液滴在不同激励参数下的形态变化和表面波特性.实验条件根据表2 和表3设置.

表2 实验工况: 激励振幅的影响 (300 Hz)Table 2 Experimental condition: influence of forcing amplitude(300 Hz)

表3 实验工况: 激励频率的影响Table 3 Experimental condition: influence of forcing frequency

图2 不同频率下的激励振幅与激励电压的关系Fig.2 The relationship between forcing amplitude and excitation voltage of the plate under different frequencies

1.2 数值模型

在竖直振动激励下液滴动力学特性的数值模拟中,将做出以下假设: 流体视为不可压缩,气体和液体彼此不混溶.表面张力和惯性力被视为体积力源项[29],并添加到Navier-Stokes 方程中.假设液滴的初始状态是半球形的.因此,具体的控制方程为

式中,u=(u,v,w) 为流体速度矢量(其中,速度分量v为平板振动方向),p为压力,ρ 为流体密度,µ 为流体动力黏度,混合物中的密度和黏度作为各相的体积分数的函数来计算

式中,下标l 和g 分别代表液相和气相,c为基于流体体积(VOF)法的流体体积函数c(x,t),描述网格单元中某一种流体所占的体积

Fs为表面张力项,在Gerris 中采用改进优化后的CSF 近似方法计算表面张力,该方法最早由Duan等[30]提出,可以表示为

式中,σ 为表面张力系数,κ 为曲率,δs为界面的狄拉克函数,n为界面的法向量,为平滑处理后的体积率函数.A为惯性加速度,其表达式为

式中,Δ 为位移振幅,Ω 为角频率,即 Ω=2πf,f为驱动频率,j为竖直方向(y轴正方向)的单位矢量.

特征参量的定义如下: 特征长度L等于液滴直径D,特征时间t=1/f,特征速度u=Df,特征密度ρ=ρL,特征动力黏度 µ=µL,特征压力p=ρLD2f2,对控制方程式(1)和式(2)进行无量纲化,无量钢化后的控制方程如下

图3(a)为竖直作用力下的物理模型示意图.在一个正方体的计算域中,液滴位于底部中心位置处,正方体计算域边长L=2D.坐标系原点与液滴中心位置重合.

图3 垂直惯性力和径向惯性力作用物理模型Fig.3 Physical model with vertical inertia force and radial inertia force

底部的边界条件设置为滑移边界,实验中测量工况内所得动态接触角范围约81°～128° (300 Hz,振幅小于100 μm),为了便于计算,模型中设置接触角恒定为90°.因为在比较高频率的激励下,波长较小,接触线处滑移距离很小,对于液体内部影响可以忽略不计.且在仿真工作中只选取了施加激励的前有限个周期的时间段,液体半径变化不大.

针对控制方程式(10)与式(11),采用时间交错法可获得时间上的二阶精度离散方程组,通过使用时间分割投影方法和多重网格泊松方程求解器求解方程[31].数值模拟中使用了可变时间步长,受Courant-Friedrichs-Lewy (CFL)数0.3 的约束,以确保计算的稳定性.仿真工况的参数设置见表4,工况S-A 和S-B对应于振动频率300 Hz,振幅分别为40 μm 和80 μm.由于液滴的惯性加速度远大于重力加速度,因此模拟中不包括重力的影响.

表4 仿真工况设置Table 4 Simulation operating condition parameters

另外,本文还探究了液滴径向加速度下失稳机制,图3(b)是径向加速度下液滴计算域示意图,相比于前文中的数值模型,此模型仅仅将竖直加速度改变为径向加速度,其他计算条件均保持一致,并且工况仍为工况S-B,此时加速度形式为

式中,i,j,k分别为直角坐标系下x,y,z方向的单位矢量,α,θ 分别为径向加速度与i,j的夹角.

1.3 数值验证

本文采用Gerris 仿真软件[32]进行计算,用于求解具有两相界面流的不可压缩Navier-Stokes 方程.Gerris 在三维模拟中使用八叉树形式的自适应规则将域离散为不同级别的计算网格.本文数值模拟的关键是捕捉液滴气液两相界面在不同工况下的变形,因此需要细化两相界面的网格.在本文中,(Minlevel,Midlevel,Maxlevel)的组合用于指定数值模拟中的网格细化级别,即每个部分细分为2Minlevel,2Midlevel或 2Maxlevel个单元.流体介质分为液体、气体和气-液界面3 部分,其中气-液两相界面要求最高的网格分辨率,设置为Maxlevel;液体和气体介质分别用Midlevel 和Minlevel 表示.本文使用自定义函数来指定基于流体体积函数c及其梯度的细化级别.这确保了以可承受的计算成本准确地解决两相界面.采用3 种网格离散标准(4,5,6),(4,6,7)和(4,6,8)进行网格无关性研究.图4 为实验中液滴顶点从平衡位置到最高位置的位移与对应工况S-A 下采用不同网格离散标准的仿真结果对比图.考虑到实验结果一致性以及仿真计算时间成本,因此采用网格离散标准(4,6,7)作为本文的网格自适应标准.

图4 液滴顶端位移的仿真和实验结果对比Fig.4 Comparison between computational results of the droplet vertex with the experiment measurement

图5 显示了与实验结果相比,液滴界面表面波形态的时间演变.在液滴宏观表面波形态特征方面,模拟结果与实验观测结果基本一致.这表明所建立的三维模拟模型和网格自适应准则在再现外部振动下液滴界面的演变方面是准确的.

图5 实验(工况E-A)和仿真(工况S-A)表面波形态的比较Fig.5 Comparison of surface wave morphology of experiment (case E-A) and simulation (case S-A)

2 结果和分析

2.1 实验结果

2.1.1 不同激励振幅下液滴界面的失稳特性

为研究不同激励振幅下液滴在振动平板上的振荡特性,实验中选取3 组典型的工况进行对比.实验中驱动频率保持为300 Hz,激励振幅按照表2 中选取.图6(a)～图6(c)分别为E-A～E-C 实验工况下的液滴的表面波形态.

图6 不同激励振幅下表面波形态的时间演化Fig.6 Time evolution of surface wave morphology under different forcing amplitudes

图6(a)为激励振幅为40 μm 下液滴的表面形态.定义液滴的顶端在其最高点时为初始时刻,记录表面波形态在一个驱动周期T(T=3.33 ms)中发生的变化.液滴的表面波主要表现为轴对称的纬向表面波.t=0 时,除了顶端表面波外,液滴还存在3 条轴向对称的纬向表面波.之后,液滴顶端开始下降,表面波的波峰和波谷在0.25T处返回到平衡位置.0.5T时,液滴的顶端处于其最低位置,表面形成4 个纬向表面波.然后,液滴在0.75T时返回平衡界面.在1T时,液滴顶端将再次返回最高点,形成顶端表面波加上3 个轴对称纬向表面波的形式.因此,可以将表面波变化周期确定为1T.

图6(b)中,激励振幅增加到160 μm.出现了显著的经向表面波,这意味着轴对称的表面波被经向不稳定性打破.以经向表面波最明显时,即处于波峰位置时作为起始时刻,很容易识别出经向表面波的产生源于和板的接触线.0 ＜t＜0.5T时,经向表面波从波峰逐渐回落到平衡位置,而在纬向表面波的作用下,经向表面波将逐渐向上移动直至消失.t=1T时,在液滴的底部再次产生清晰的经向表面波.然而,此时的经向表面波的位置是在t=0 的波谷位置.比较t=2T时液滴经向表面波的位置,可以发现液滴经向波的周期为 2T.t=0.5T和t=1.5T处的液滴顶端位置表明,在经向波出现后,纬向波仍然存在,并保持1T的周期.

进一步增加振幅到200 μm,从图6(c)可以看出,在更大的惯性力作用下,液滴经向和纬向的表面波振幅显著增大.两种表面波之间也存在相互干涉,使得液滴表面由经向和纬向表面波组成胞格状的驻波.此时,液滴的纬向和经向表面波的形状演变将趋于混乱,但基于0T和2T处液滴形状的相似性,推测表面波演变周期仍然是2T.

2.1.2 不同驱动频率下液滴界面的失稳特性

在本节,改变驱动频率,并适当调节激励电压使得振动平板的位移振幅保持相对一致,以研究不同驱动频率下液滴所产生的表面波特性.实验工况的参数设置见表3.

图7(a)～图7(e)对应不同驱动频率下液滴的表面波形态图.其中图7(a)～图7(d)以液滴顶端位于最高点作为初始参考对象,分别给出5 个不同位置,即最高位置、平衡位置、最低位置、平衡位置和最高位置作为参考位置.从图中可以看出,不同频率下液滴的纬向表面波周期与对应的驱动周期一致.液滴此时的表面波始终为简谐波.另外,不同频率下液滴产生的表面波也有区别.高驱动频率下液滴表面波模态数明显增多,这与Pang 等[33]的实验规律一致.图8 为不同驱动频率下对应的液滴纬向表面波的波长,波长定义为沿球面经向所有波长的平均值.从图中可以看出,随着频率的增加,纬向表面波的波长首先快速降低,当频率大于400 Hz 时,波长降低的速率减缓.当频率较大时,如图7(e)所示,此时液滴表面产生2 倍于驱动周期的经向表面波,在Pang等[33]关于圆形和星形振荡之间相互转换的临界振幅研究中也观察到了类似现象.这可以解释为高频率下为了保持振幅一致而增加的电压促使液滴发生了经向失稳,液滴此时的经向波为亚简谐波.

图7 不同驱动频率下表面波形态的时间演化Fig.7 Time evolution of surface wave morphology under different forcing frequencies

图8 液滴纬向表面波平均波长随频率的变化Fig.8 Variation of the average wavelength of droplet latitudinal surface waves with frequency

根据以上分析,通过实验手段在宏观层面上观测到在激振作用下产生的经向波和纬向波,并探索了不同激振频率和激励振幅对两种波形的影响.具体来讲,不同激励振幅下液滴存在着不同的界面失稳特性.在低振幅下,液滴界面主要发生纬向失稳从而产生具有轴对称性质的纬向表面波,而纬向表面波的周期等于驱动周期1T.在高振幅下,液滴界面不仅存在纬向表面波,还会发生经向失稳产生经向表面波,经向表面波周期为驱动周期的两倍,即为2T.而当激励振幅足够大时,表面波的形状演变将趋向于混沌,但演变周期仍然是2T.频率的变化也会引起振动模式的转变,频率增加,表面波模态数增加,相应波长减小.同振幅下,频率继续增加,波形会发生从只存在纬向波模式向纬向波叠加经向波模式的转变.然而仅仅通过实验手段不能说明表面波演化的微观内在机理.

2.2 仿真结果

为了探明表面波演化的内在机理,通过仿真的手段对液体内部速度流场和压力波动进行研究.

2.2.1 竖直惯性力作用下的纬向失稳特性

图9 是液滴在 工况S-A 下纬向表面波的生成及发展过程.在初始时刻,液滴处于静止状态.然后,在惯性力的作用下液滴整体向正向(y轴正方向)移动,因此在表面张力的作用下液滴底部向内收缩.当t=0.25T时惯性力反向,但是液滴在惯性的作用下仍然有正向移动的趋势,因此液滴底部继续向内收缩.当t=0.55T时,液滴底部停止收缩,底部处于波谷位置,底部附近形成完整表面波的一半.随后液滴在惯性力的作用下开始加速向y轴负向移动,底层液滴受到上方液层的挤压将会抵抗表面张力逐渐向外扩张.t=1T时,液滴底部扩张到最大,形成波峰.此时,底部波峰和向上传播的波谷形成一组完整的表面波.随后液滴又将开始向正向移动并重复这样的过程.当表面波接近液滴顶端时(t=2T),液滴顶端到液滴底部已经形成了两组表面波.液滴底部首先产生表面波,且波的传播方向是朝向液滴顶端的.

图9 液滴纬向表面波演化仿真结果(工况S-A)Fig.9 Simulation results of droplet latitudinal surface wave generation process (case S-A)

由这个过程可以得知,液滴在竖直惯性加速度下立即产生纬向方向的失稳.本文在其他工况下皆出现此现象,Vukasinovic 等[34]在灵敏度更高的实验系统中也观察到了同样的现象,他们将此现象归结于单自由度系统在小强迫振幅下的行为,这也验证了本文仿真模型的正确性.

图10 是工况S-A 下液滴纬向表面波在一个稳定周期内随时间的演化过程.本文定义与顶端表面波相邻的表面波为次级表面波.以液滴顶端表面波处于波峰位置时作为周期的起始点,此时次级表面波处于波谷状态.根据图10 可知,液滴纬向表面波演化周期与所施加惯性力的周期一致,这意味着液滴纬向表面波的频率等于驱动频率,即此时液滴纬向表面波是简谐波.除此之外,发现液滴顶端位移与惯性力相位相差1/4 个周期.

图10 稳定循环周期液滴顶端位移和惯性力随时间的演变(工况S-A)Fig.10 Evolution of droplet vertex displacement and inertial force with time during the stabilization cycle (case S-A)

为了进一步研究液滴纬向表面波变化规律的内在机理,图11 给出了液滴纬向表面波各关键时刻下液滴xoy截面的压强场及速度场图,其中黑色粗曲线为液滴的表面轮廓线,黑色小箭头为速度矢量,压强为相对压强.液滴顶端处于波峰状态的时刻被选择为循环的开始.右侧视图给出了三维状态下液滴表面的压强场分布.液滴的波峰和波谷下方的液体区域分别是高压区域和低压区域.从整体来看,液滴表面各位置波峰和波谷是周期性交替变化的,而对应的高压区和低压区也会相应地周期性地交替变化.图12 为t=0 时刻液滴表面从顶端到底端的压强分布示意图.液滴顶端高压区的压力最高,大约为124 Pa.而顶端之外的高压区的压强大小相似,都低于顶端高压区的压强.

图11 工况S-A 在xoy 截面上不同时刻的液滴压力场和速度矢量场Fig.11 Droplet pressure field and velocity vector field at different moments for case S-A in the xoy cross section

图12 液滴在x 方向压力分布(工况S-A,t=0)Fig.12 Pressure distribution from top to bottom surface of droplet(case S-A,t=0)

0 ＜t＜0.25T时,液滴顶端表面波由波峰位置向平衡位置移动.而此时液滴顶点处的表面张力也是向下的,所以表面张力对液滴顶端表面波做正功.而根据图10 可知,此时惯性力方向也是向下的,因此惯性力对顶端表面波也做正功,这将进一步增大顶端表面波下方液体区域的动能,加速顶端表面波向平衡位置偏移.顶端液体区域向下流动时会逐渐流向次级表面波的低压区,原本低压区处的压力逐渐升高,而表面波也从波谷位置向平衡位置处移动.其余位置处的表面波也将发生类似的流动及能量转化.由图11(b)可知,当t=0.25T时,液滴各表面波处于平衡位置.图11(b)图片右侧是左侧顶端部分的放大图,红色箭头为液体流动方向,绿色箭头为惯性力的方向.由图可知,此时液滴顶端下方流体仍具备较高的动能往次级表面波位置处流动.

当 0.25T＜t＜0.5T时,惯性力向下并且逐渐减小,因此惯性力仍然对顶端表面波做正功,表面张力对顶端表面波做负功.由于刚开始表面张力小于惯性力,因此一开始合力仍然对顶端表面波做正功,这将使得其压力势能转化为动能.顶端下方流体向周围流动,使得次级表面波由平衡位置向波峰位置移动,此时表面张力和惯性力都对其做负功,这将增加对应区域的压力势能,逐渐形成高压区.由于顶端表面波向波谷位置偏移时,液滴界面曲率将逐渐增大,表面张力也会相应增大,当表面张力大于惯性力时,液滴顶端表面波下行的速度将会降低.当t=0.5T时,在表面张力和惯性力的共同作用下,液滴顶端表面波行至波谷位置,由于下方流体区域的流动,其能量逐渐传递给次级表面波区域,所以形成顶端低压区,而次级表面波接受能量,压力势能增加,逐渐形成高压区.此高压区的压强大约为55 Pa.图13 为液滴顶端及次级表面波下压强随时间的变化规律.其中,液滴顶端及次级表面波压强取点位置相距半个表面波长.图11(c)右侧视图为三维状态下液滴表面的压强场分布,这与t=0 时正好相反.进一步说明液滴纬向表面波处于稳定周期性变化时,其压强场也有类似于表面波的驻波特性.

图13 顶端表面波和次级表面波作用下压力随时间的变化(工况S-A)Fig.13 The variation of pressure with time under the top surface wave and secondary surface wave (case S-A)

当 0.5T＜t＜0.75T时,液滴顶端表面波由波谷位置向平衡位置恢复,此时惯性力方向由负向转为正向,因此此时表面张力做正功,惯性力也做正功,液滴顶端表面波动能逐渐增加.当t=0.75T时,图11(d)右侧图片为液滴顶端部分的放大图.高压点出现在液滴的顶部.对速度矢量场的分析表明,高压点的出现是由于流体在圆周方向上的碰撞以及动能转化为压力势能.

在整个周期的最后一部分,即 0.75T＜t＜1T时,液滴顶端表面波由平衡位置向波峰方向移动,此时惯性力方向为正,做正功,表面张力做负功.随着液滴顶端的曲率逐渐增大,表面张力也在迅速增大,而惯性力在逐渐减小.当表面张力超过惯性力时,液滴顶端的移动将被减速.并且由于表面张力始终做负功,导致顶端表面波动能逐渐转化为压力势能,因此液滴顶端的高压区逐渐形成.当t=1T时,液滴整体形态、压强场和速度矢量场逐渐回到周期起始的状态.表明纬向波周期是1T,这与实验相对应.

分析可知液滴形成轴对称的周期性纬向表面波机理: 液滴表面张力和惯性力的共同作用下,通过液滴表面波完成不断低能量转化和传递,逐渐形成频率与驱动频率一致的周期性表面波.

2.2.2 竖直惯性力作用下的经向失稳特性

根据以上分析,液滴在惯性力和表面张力的共同作用下发生纬向失稳并形成稳定的纬向表面波.当激励增大致使惯性力和表面张力无法再通过纬向波进行能量的转化、传递而达到平衡状态时,就会产生经向表面波.本节仿真工况为S-B.

在图14 中,顶视图用于更好地呈现液滴的经向表面波.t=0 时液滴处于静止状态,随着惯性力的作用,液滴发生纬向失稳而产生纬向表面波,纬向波发展到 1 5.25T时基本保持稳定且只有纬向波存在,波形保持驻波特性.而当t=19.05T时,液滴周向出现了经向表面波,如图红色箭头所示,此时经向表面波的振幅较小.经向表面波的出现意味着液滴发生了经向失稳,但同时注意到经向表面波在液滴底部位置较明显,靠近液滴顶部位置仍然为较为平滑的纬向表面波.当t=21.35T时,液滴的经向表面波振幅增加,此时纬向表面波和经向表面波共存.对比工况S-A 和工况S-B,说明 ΔD数的增大会增加液滴经向失稳的趋势.此外,施加激励后的较长时间后才会发生经向不稳定性.

图14 液滴的表面波演化(工况S-B)Fig.14 Surface wave evolution of liquid droplets (case S-B)

当t=23.35T时,液滴经向波振幅最大,并且存在波峰和波谷.在下一时刻,处于波峰的表面波和处于波谷位置的表面波都将向平衡位置移动.因此,当t=23.85T时,液滴经向表面波处于平衡位置.而当t=24.35T时,原先处于波峰位置的表面波变成波谷,而原先波谷位置的表面波变成了波峰.随后当t=24.85T时液滴又回到经向平衡位置,而当t=25.35T时完成完整的液滴经向表面波的演化周期.液滴经向表面波演化周期大约为2T,即液滴经向表面波的频率为驱动频率的一半,这与实验结论也相对应.经向表面波的次谐波相应表明,液滴在经向方向上的界面不稳定性是一种类似于Faraday 波的参数驱动不稳定性.同时观察液滴顶部纬向表面波可以发现,此时液滴纬向表面波的周期仍然接近 1T,即纬向表面波仍然为简谐波.

图15 是t=23.35T时刻的液滴形态及其表面处压强分布图.可以看出,处于波峰位置处的经向表面波处产生了压力集中形成经向高压区,且同一条经波压强从底部到顶部逐渐减小.由经向波形态和压强分布推断,促使液滴经向失稳形成经向表面波的能量集中在液滴底部,并且越往顶部逐渐降低.

图15 经向波波形和表面压力分布(工况S-B,t=23.35T)Fig.15 Longitudinal wave form and surface pressure distribution(case S-B,t=23.35T)

由上文分析可得液滴底部的状态在液滴界面演化过程中具有重要的意义.图16 是经向波半个周期内不同时刻下液滴底面的压强场和速度矢量场图.压强变化周期为1T.当t=23.35T时,液滴底部中心压强最高,并向四周逐渐降低.从速度矢量图中可以得知,液滴靠近底部的界面发生径向运动,接触线向内收缩,这将进一步增大底部中心压强.

图16 液滴底部的压力场和速度矢量场的轮廓(工况S-B,y=0)Fig.16 Contour of the pressure field and velocity vector field at the bottom of the droplet (case S-B,y=0)

当 2 3.35T＜t＜23.75T时,液滴底部中心高压强将推动液体反向运动.波谷处因曲率更大而拥有更高的压力势能,因此当接触线向外移动时压力势能转化成动能产生更高的速度.当t=23.75T时,底部液体流向远端导致远端压强高于内部,且之前的波谷处压强高于两侧.当 2 3.75T＜t＜24.35T时,之前处于波峰的表面波变成了波谷,波谷变成了波峰,而液滴底面中心处压强完成了周期性的变化,而在下个压强周期完成时,经向表面波完成完整的周期.

综上分析,竖直激励振幅较大时,首先立即出现纬向不稳定性波,施加激励后的较长时间后才会发生经向不稳定性.并且经向表面波的频率为驱动频率的一半,即次谐波效应.因此猜测竖直激励下液滴在经向方向上的界面不稳定性是一种类似于Faraday波的参数驱动不稳定性.

2.2.3 径向惯性力作用下的经向失稳特性

由以前的研究可知[15],Faraday 波具有次谐波效应.为了探明竖直惯性力下的经向波和由Faraday不稳定产生的表面波之间的联系,对液滴施加垂直于表面,也就是径向方向的激励.这一工作通过实验手段较难实现,所以用建立的径向加速度模型进行仿真研究,工况为S-B.

图17 是工况S-B 径向加速度下由Faraday 不稳定性主导的液滴界面失稳时表面波演化状态,此时液滴直接发生经向失稳产生经向表面波而不会产生纬向表面波,经向表面波具有次谐波的特性,周期为2T.此时液滴界面失稳而产生的表面波与竖直加速度下有着密切的联系.由于加速度垂直于液滴界面,因此液滴底部界面将完全处于径向运动状态,在惯性加速度作用下发生失稳形成经向表面波.

图17 径向加速度下Faraday 不稳定性引起的液滴界面不稳定性(工况S-B)Fig.17 Droplet interface instability dominated by Faraday instability under radial acceleration (case S-B)

在仿真结果中,波长被定义为液滴上所有表面波的平均波长.根据线性理论,Faraday 表面波的波长可以近似计算为 Λ=2πr0/l,式中,l为模态数,r0为液滴半径.从图17 中可以看出,液滴模态数为8,可得波长直径比约为0.418.表5 是工况S-B 不同加速度形式下液滴经向表面波物理特征值.

表5 不同加速度形式下经向波的物理特征值Table 5 Physical characteristic values of longitudinal surface waves in different acceleration forms

因此,在同一工况下,竖直加速度下液滴经向表面波的波长直径比与经向加速度下及线性理论结果下的波长直径比几乎一致,这进一步说明液滴在竖直加速度下的经向失稳是由Faraday 不稳定性主导的.而其振幅直径比的差异可以归因于在同一 ΔD下,一部分竖直惯性力用于维持液滴发生纬向失稳.

结合章节2.2.2 及2.2.3 的分析.竖直激励下液滴界面发生经向失稳的机制可以归结为: 竖直惯性力通过液滴的几何特征使液滴底部径向运动,这种径向运动会在接触线处产生水平方向的径向力,其方向垂直于接触线处的液面,从而产生Faraday 不稳定性.当惯性力较大超过临界条件时,一部分径向力用于维持纬向失稳产生纬向表面波,另一部分将促进形成新的失稳特征,即液滴的经向失稳,并且经向表面波的频率为驱动频率的一半.

3 结论

本文基于平板上液滴振荡实验系统结合三维数值仿真模型研究了竖直振动平板上液滴的表面波的演化特性,得到主要结论如下.

(1) 基于平板上液滴振荡实验系统,分别研究了液滴界面在不同激励电压和不同激励频率下液滴失稳特性.发现在激励电压较小,即惯性力较小时,液滴表面发生纬向失稳并且形成轴对称的纬向表面波,并且纬向表面波的频率等于驱动频率.增大激励电压液滴纬向表面波振幅增加.继续增加激励电压到一定程度时,液滴表面会发生经向失稳并产生经向表面波,经向表面波频率为驱动频率的一半.当激励电压足够高,即惯性力足够大时,液滴纬向表面波和经向表面波会相互干涉,此时液滴表面形状演化趋于混沌.频率的变化也会引起振动模式的转变,频率增加,表面波模态数增加,相应波长减小.同振幅下,频率继续增加,波形会发生从只存在纬向波模式向纬向波叠加经向波模式的转变.

(2) 建立了滴液受迫振动的三维数值仿真模型,对低 ΔD下液滴仿真结果中的物理场及液滴顶点位移与惯性力相位关系的研究,阐明了液滴形成轴对称的周期性纬向表面波机理: 液滴表面张力和惯性力的共同作用下,通过液滴表面波完成不断低能量转化和传递,逐渐形成频率与驱动频率一致的周期性表面波.并且表面张力趋向于稳定液滴,惯性力趋向于促进液滴失稳.

(3) 在高 ΔD下,通过对液滴经向表面波形态及液滴内部物理场进行分析,并且结合Faraday 理论与同工况径向加速度下由Faraday 不稳地性主导的液滴界面失稳产生的表面波物理特征进行对比,阐明液滴在竖直加速度下发生经向失稳机理: 由于液滴的几何特征使得液滴接触线处产生径向力,径向力法向地作用于液滴底部界面.当增加惯性力达到临界条件时,随之增加的周期性径向力使得液滴底部界面发生变形,并由Faraday 不稳定性主导,最终形成Faraday 波即经向表面波,并且其频率为驱动频率的一半.