INFLUENCE OF LIQUID VISCOSITY ON THE INTERFACIAL INSTABILITY OF CYLINDRICAL DROPLET INDUCED BY A CAVITATION BUBBLE
-
摘要: Rayleigh-Taylor不稳定性存在于爆炸、液滴形成和液体喷雾等工程应用过程中, 是流体力学关注的经典问题之一. 内空泡振荡诱导液滴界面演化问题是其研究中基本模型之一, 空泡振荡作用下液滴界面发生扰动并发展, 其特征形态主要表现为破碎、通气和稳定. 液体黏性是影响界面不稳定性发展的重要因素, 文章通过建立高精度的数值模拟方法, 开展液体黏性对内空泡诱导柱状液滴界面不稳定性的影响研究. 在数值模拟中, 基于开源OpenFOAM框架的多相可压缩求解器直接求解Navier-Stokes方程, 采用isoAdvector的几何流体体积法捕捉界面演化特征. 结果表明, 液体黏性的增加会减缓空泡的收缩, 进而减缓液滴界面扰动的发展, 该影响下通气工况液滴通气发生时间增加, 而稳定工况最大扰动幅值减小. 最大扰动幅值的减小直接影响了液滴的特征形态, 基于一系列数值模拟结果归纳得到液滴不稳定性相图. 在文章讨论的参数范围内, 随着黏性增加, 小液滴(Rd0 < 2 mm)的形态从破碎转变为通气进而变成稳定; 中液滴(2 mm < Rd0 < 3 mm)的形态从通气转变为稳定, 不出现破碎形态; 而大液滴(Rd0 > 3 mm)的形态不随液体黏性改变, 液滴形态会随着液体黏性增加趋于稳定.
-
关键词:
- Rayleigh-Taylor不稳定性 /
- 液体黏性 /
- 空泡振荡 /
- 柱状液滴界面不稳定
Abstract: Rayleigh-Taylor instability is one of the classical problems in fluid mechanics, which is widely used in underwater explosions, droplet formation, liquid sprays, and other engineering applications. The interfacial instability of cylindrical droplet induced by a cavitation bubble is one of the classical models in its study. As the perturbation of the interfacial instability develops, three distinct characteristics of the droplet deformation are obtained: (i) a splashing state; (ii) a ventilating state; and (iii) a stable state. Considering that liquid viscosity is an important factor affecting the development of interfacial instability, the influence of liquid viscosity on the interfacial instability of a cylindrical droplet induced by a cavitation bubble is investigated through a high-precision numerical simulation method in this paper. A direct numerical simulation is set up based on the compressible multiphase solver in the OpenFOAM framework. The interfaces are captured using a geometric fluid volume method named isoAdvector. The results show that within the parameters discussed in the paper, an increase in liquid viscosity slows down the oscillation of the bubble, thus reducing the perturbation growth rate during the bubble collapse. In these cases, the onset time of ventilation increases for a ventilating state, and the maximum perturbation reduces for a stable state. This reduction directly affects the distinct characteristics of the droplets. A phase diagram of the interfacial instability of droplets with different initial radius and liquid viscosity is obtained based on the results of numerical simulations. With the increase of liquid viscosity, the small droplets (Rd0 < 2 mm) change from a splashing state to a ventilating state and then a stable state; the medium droplets (2 mm < Rd0 < 3 mm) change from a ventilating state to a stable state; whereas, the large droplets (Rd0 > 3 mm) does not change and stay only in a stable state. The droplet tends to stabilize as the liquid viscosity increases. -
引 言
Rayleigh-Taylor (RT)不稳定性是一种液体或气体界面上出现的不稳定现象. 当两个不同密度的流体在接触面上形成密度梯度时, 通过重力或其他外部力的作用, 界面发生扰动和不规则变形, 并产生波浪、褶皱或涡旋等[1-2]. RT不稳定性现象在自然界和工程应用中广泛存在[3-5]. 例如, 在海洋中不同密度的海水层之间的相互作用可能导致海浪的形成和传播[6-7]. 在火山喷发过程中, 熔岩和气体之间的界面会出现RT不稳定性, 导致喷发的物质和形状的变化[8-9]. 在惯性约束聚变中, 两种不同密度的激光束界面可能引RT不稳定性, 这可能导致激光能量的损失和散射, 从而影响聚变反应的效率[10-11]. 在航空航天领域, 飞行器在大气中飞行时, 机翼表面的气流分离可能导致RT不稳定性, 从而影响飞行器的性能和稳定性[12]. 因此, 理解与定向控制RT不稳定性对于工程应用来说是至关重要的.
Taylor等[13-14]提出平面几何界面线性不稳定性分析理论, 结合实验总结了RT不稳定性的发展的几个阶段: (1)线性阶段, 扰动振幅远小于扰动波长, 扰动随时间指数增长, 可以用线性理论描述; (2)变形阶段, 扰动振幅发展, 界面开始表现出不对称特征[15]; (3)规则非线性阶段, 低密度流体形成空泡进入高密度流体, 高密度流体形成尖钉进入低密度流体中. Birkhoff[16]在其基础上总结出了规则非线性后的不稳定性发展阶段; (4)不规则非线性发展阶段, 界面发展受到速度差异等多种作用耦合影响; (5)湍流混合阶段, 界面两侧流体会发生掺混.
同时, Bell[17]和Plesset[18]分别给出了柱面和球面几何下的RT不稳定性演化方程, 建立了运动几何效应影响下的不稳定性理论. 这种界面运动轨迹影响界面扰动发展的效应被称为Bell-Plesset (BP)效应. BP效应下界面不稳定性在工程中同样广泛存在, 比如预混燃烧[19]、爆炸[20-21]和液体喷雾[22], 可以简化成内部体积振荡诱导液滴界面变形这一模型. Obreschkow等[23]通过高速摄像机观察了微重力条件下球形液滴中空泡振荡诱导界面扰动、变形等现象. 吕明等[24]基于数值模拟研究了液滴内空泡生长的规律, 提出了空泡溃灭的3个历程, 包括快速溃灭期、缓慢溃灭期以及稳定期. Zeng等[25]基于实验和数值模拟研究了液滴内空泡溃灭产生的射流, 他们提出内部空泡振荡产生的径向加速度诱导了界面RT不稳定性的形成. 为了在更大参数空间评估界面稳定性, 推导了包含液体黏性的球面RT不稳定性扰动增长解析模型, 结果发现发现空泡溃灭后产生的表面射流会随着液体黏性增加而减少. Wang等[26]鉴于实验观测手段难以捕捉三维不稳定性的发展, 将模型简化至柱坐标系, 结合实验、数值模拟和理论全面研究了板间液滴下空泡溃灭诱导的RT不稳定性现象, 他们总结了对于不同液滴初始半径和空泡初始内压工况下液滴出现的破碎、通气和稳定3种形态特征(如图1所示), 同时建立了考虑RT不稳定性和空泡振荡的分析模型, 用于分析扰动的增长和液滴形态的评估. Zhang等[27]对不同空泡和液滴半径比的飞溅情况进行分析, 结果表明随着半径比的增加, 射流飞溅效果会更加显著, 在这种情况下空泡的溃灭时间会变小. Rosselló等[28]对自由落体液滴内空泡溃灭的界面现象展开研究, 发现空泡溃灭诱导液滴表面出现RT不稳定性和Rayleigh-Plateau不稳定性, 并提出了不稳定性产生的机制: 空泡快速膨胀或收缩加速液层导致液滴表面波纹产生, 波纹在界面进一步发展从而观察到不稳定性现象.
在工程应用中, 液体黏性是影响不稳定性产生和发展的关键因素, 对应用的效率和可靠性有着重要的影响. Zeng等[25]提出了考虑黏性的球形液滴内空泡诱导RT不稳定性形成的机理, 但是在黏性对液滴特征形态和不稳定性评估的影响规律尚未研究清楚. 本文采用直接数值模拟和高精度界面捕捉方法, 研究了液体黏性影响下柱状液滴内空泡诱导不稳定性的产生与发展, 以期为惯性约束核聚变等工程应用提供参考.
1. 数值模拟方法
1.1 计算模型
本文基于开源框架OpenFOAM建立直接数值模拟方法, 开展液体黏性对内空泡振荡诱导柱状液滴界面不稳定性的影响规律研究.
计算模型和边界条件如图2所示, 计算域内液滴和空泡分别用蓝色和白色表示. 为了模拟空泡的振荡, 初始空泡初始设置成高温高压状态. 模拟过程忽略轴向的变化而重点关注径向的不稳定性现象. 计算域为三维矩形结构, 选取1/4对称结构开展计算, 其大小为12 mm × 12 mm × 1 mm, 右面和下面为对称面, 上面和左面为出口条件, 前面和后面为空边界. 全区域采用结构性网格均匀划分, 最小网格尺寸Δx = 20 μm, 其中初始空泡内包含20个网格.
1.2 控制方程
数值模拟中, 通过多相可压缩求解器compressibleInterIsoFoam直接求解可压缩的纳维-斯托克斯方程, 假设气相和液相均为可压缩的互不相溶牛顿流体[25], 考虑两相间的热量传输[26,29], 其控制方程为ρ
连续性方程
$$ \frac{\partial \rho }{\partial t} + \nabla \cdot (\rho {{\boldsymbol{U}}}) = 0 $$ (1) 动量方程
$$ \frac{{\partial \rho {{\boldsymbol{U}}}}}{{\partial t}} + \nabla \cdot (\rho {{\boldsymbol{UU}}}) = - \nabla p + \nabla \cdot {{\boldsymbol{\tau}} }{ } + \sigma \kappa {\delta _s}{{\boldsymbol{n}}} $$ (2) 能量方程
$$ \frac{{\partial \rho (e + K)}}{{\partial t}} + \nabla \cdot [\rho (e + K){{\boldsymbol{U}}}] - \nabla \cdot (\beta \nabla e) + \nabla \cdot (\rho {{\boldsymbol{U}}}) = 0 $$ (3) 式中, t为时间, ρ为混合物密度, U为速度矢量, p为压力, τ为应力张量, σ为表面张力系数为0.07 N/m, $\kappa $为界面曲率, δs为特征函数, 在界面处网格为1, 非界面为0, n为界面单位法向量, e为比内能, K为比机械能, β为热扩散系数, 液相热扩散系数为1.4 × 10−7 m2/s, 气相热扩散系数为1.8 × 10−5 m2/s.
考虑到气相和液相之间会发生热交换, 其状态方程采用理想气体和理想流体进行建模求解[26,30], 表达式分别为
$$\qquad\qquad\qquad {\rho _g} = \frac{1}{{RT}}p $$ (4) $$\qquad\qquad\qquad {\rho _l} = \frac{1}{{RT}}p + {p_{l0}} $$ (5) 式中, ρg为气相密度为1.29 kg/m3, ρl为液相密度为1000 kg/m3,R为常数为8.31 J/(mol·K), T为温度, pl0为空泡初始内压.
两相界面通过volume of fluid (VoF) 方法捕捉, 其核心是引入流体体积分数α, 定义为对应相体积与网格体积的比值(气相为0, 液相为1), 通过求解相分数方程求解实现两相界面捕捉[31-32]. 相分数方程定义为
$$ \begin{split} & \frac{{\partial \alpha }}{{\partial t}} + {{\boldsymbol{U}}} \cdot \nabla \alpha + \nabla \cdot \left[\alpha (1 - \alpha ){{{\boldsymbol{U}}}_{r}}\right]= \\ &\qquad - \alpha (1 - \alpha )\left(\frac{1}{{{\rho _l}}}\frac{{{\mathrm{D}}{\rho _l}}}{{{\mathrm{D}}t}} - \frac{1}{{{\rho _g}}}\frac{{{\mathrm{D}}{\rho _g}}}{{{\mathrm{D}}t}}\right) \end{split} $$ (6) 式中, α为体积分数, Ur为界面处两相之间的相对速度, 作用于人工压缩确保界面尖锐. 基于流体体积分数, 两相系统中的密度和黏度定义为
$$ \rho = \alpha {\rho _l} + (1 - \alpha ){\rho _g} $$ (7) $$ \mu = \alpha {\mu } + (1 - \alpha ){\mu _g} $$ (8) 式中, μ为液相的黏度, 具体取值在后文讨论, μg为气相的黏度为1.0 × 10−5 Pa·s.
考虑到不稳定性发展对于精细界面捕捉的要求, 本文采用了isoAdvector的几何VoF方法. 该方法采用了有效等值面的思想, 对等值面处网格点切割并重新进行跨面连接, 形成单元内部新的等值面. 在迭代过程中, 该方法假设界面在子时间步长的间隔期内稳定移动, 来分析计算通过网格表面的体积分数通量, 过程中没有对单元形状作出假设, 因此可以适用于任意形状的网格[33].
求解过程采用PIMPLE算法[34]求解瞬态的可压缩纳维-斯托克斯方程, 采用自适应时间步长基于库朗数进行调整, 最大库朗数设定为0.5. 每个时间步分为外循环和内循环, 外循环首先求解相分数方程, 然后顺序求解连续性方程、动量方程和能量方程, 内循环是压力修正过程, 利用温度通过状态方程更新各相密度和混合物密度实现压力修正, 达到求阶精度后结束循环, 并步进到下一时间步, 直至达到预定求解时间步停止求解. 求解过程时间采用Euler格式进行离散, 相分数方程对流项采用Gauss-vanLeer格式进行离散, 动量方程对流项采用Gauss线性格式进行离散, 其余项均为Gauss迎风格式.
1.3 数值模拟验证
为了验证数值模拟的准确性, 本小节讨论了数值模拟结果的网格无关性, 并与参考文献[26]中的实验结果进行对比. 根据实验的空泡最大半径和溃灭时间, 确定数值模拟空泡的初始条件为Rb0 = 0.2 mm, Pb0 = 45 MPa, Tb0 = 1500 K; 液滴的初始条件与实验一致, 初始液滴半径Rd0 = 6 mm, 初始温度Td0 = 300 K. 网格无关性选取了粗、中、细3套网格, 最小尺寸分别为40, 20和10 μm, 具体参数见表1. 表1给出了3套网格中数值模拟空泡最大半径Rmax的结果, 并与实验进行了对比, 同时给出不同网格对扰动波数n和最大扰动幅值ηmax的计算结果, 最大扰动幅值为扰动最大位置(Rηmax(t))与最小值位置(Rηmin(t))的差: ηmax = 0.5(Rηmax(t) + Rηmin(t)), 结果显示空泡最大半径、扰动波数和最大扰动幅值误差逐级减小. 图3(a)展示了液滴半径(Rd)、空泡半径(Rb)及扰动幅值(η)变化曲线对比结果, 液滴半径通过当前时刻扰动最大位置(Rηmax(t))和最小值位置(Rηmin(t))的平均值计算得到Rd(t) = 0.5(Rηmax(t) + Rηmin(t)). 其中圆形和方形标记分别代表实验的空泡半径和液滴半径, 线条为数值模拟结果, 蓝色虚线、橙色实线和绿色点线分别对应了粗网格、中网格和细网格. 结果表明3套网格计算结果逐步收敛, 粗、中和细网格空泡最大半径与实验结果的误差分别为0.56%, 0.43%和0.37%. 因此采用中网格进行计算是合适并且准确的.
表 1 粗中细网格参数细节Table 1. Details for coarse, medium and fine gridsExp. Coarse Medium Fine cells/103 — 90 360 1440 Δx/μm — 40 20 10 Rbmax/mm 5.378 5.346 5.355 5.358 error Rbmax — 0.56% 0.43% 0.37% wave number n — 13 15 15 ηmax/mm — 2.511 1.849 1.729 图3(b)展示了4个典型时刻液滴和空泡形态的实验与数值模拟对比结果, 其中左半为数值模拟结果, 右半为实验结果, 气液界面为α = 0.5等值线. 过程中, t = 0.1 ~ 0.45 ms为空泡的膨胀阶段, t = 0.45 ~ 0.85 ms为收缩阶段. 从t = 0.1 ms开始, 高温高压的空泡向外膨胀, 并带动液滴膨胀, 此时液滴界面仍然保持光滑. 随后在t = 0.3 ms液滴界面观察到扰动的产生(如黄色方框局部放大图所示), 此时不稳定性开始出现, 随后空泡继续膨胀在t = 0.45 ms达到最大半径, 液滴界面逐渐明显. 收缩阶段, 随着空泡体积迅速收缩, 液滴界面扰动快速增长, 在t = 0.7 ms时能观察到扰动明显变大, 该阶段内不稳定性迅速发展. 整个发展过程实验与数值模拟对比良好, 本文建立的数值模拟方法可以准确描述空泡振荡和液滴界面扰动演化过程. 下文基于该空泡初始条件, 改变液体黏性(μ)和液滴的初始半径(Rd0)讨论液体黏性对柱状液滴界面不稳定性的影响.
2. 结果与讨论
柱状液滴内空泡溃灭诱导不稳定性产生, 随着不稳定性的发展液滴会展现出稳定、通气和破碎3种不同的特征形态, 在液体黏性影响下空泡和液滴的行为都会发生改变.
2.1 液滴破碎形态
图4展示了液滴初始半径Rd0 = 1.5 mm不同液体黏性条件下液滴与空泡的行为, 其中图4(a) ~ 图4(c)分别对应液体黏性为μ = 0.001, 0.05, 0.1 Pa·s, 该工况下液滴表现为破碎特征形态. 图4(a)观察可得, 随着空泡膨胀, 液滴界面开始向外运动, 在膨胀初期t = 0.04 ms观察到液滴表面初始扰动在产生. 随后扰动继续增长, 在空泡和液滴膨胀到最大半径前与外界大气形成通路. 此时处于膨胀阶段, 扰动穿透后空泡和液滴由于惯性继续向外运动, 液滴在扰动增长和空泡膨胀的共同作用下发生破碎(t = 0.15 ms), 不再保持完整形态. 这种在膨胀阶段扰动穿透空泡界面的情况被称为液滴的破碎形态.
图4(b)展示了液体黏性变大后的液滴和空泡行为. 当液体黏性μ = 0.05 Pa·s时, 扰动在t = 0.05 ms时产生, 且扰动幅值变小. 扰动穿透空泡界面仍处于膨胀阶段, 随后液滴发生破碎.
2.2 液滴通气形态与通气发生时间
图5展示了液滴初始半径Rd0 = 2.5 mm不同液体黏性条件下液滴与空泡的行为, 其中图5(a) ~ 图5(c)分别对应液体黏性为μ = 0.001, 0.05, 0.1 Pa·s, 该工况下液滴初始半径相较于图4变大, 液滴表现为通气特征形态. 图5(a)观察可得, 膨胀阶段空泡带动液滴向外膨胀, 在t = 0.2 ms扰动开始产生, 其相比于破碎形态(图4)更晚, 并且此时空泡界面与液滴界面距离更大. 空泡在t = 0.25 ms膨胀到最大半径, 随后进入收缩阶段. 在收缩阶段扰动幅值继续增大, 扰动穿透空泡界面并在空泡和外部大气形成通路(t = 0.35 ms), 定义该时刻为通气发生时间(tv), 此时空泡内部压力小于外部大气压, 气体从外部进入内部发生通气. 这种在收缩阶段扰动穿透空泡界面的情况被称为液滴通气形态.
图5(b)和图5(c)展示了液体黏性变大后的液滴和空泡行为. 液体黏性增大后扰动产生时间延后, μ = 0.05, 0.1 Pa·s分别在t = 0.22, 0.23 ms观察到扰动的产生. 当空泡膨胀到最大半径, 扰动幅值明显小于μ = 0.001 Pa·s黏性工况. 随着扰动发展其通气发生时间也延后, 分别是t = 0.38, 0.4 ms. 随着液体黏性增大, 抑制了扰动的发展, 使得需要更长时间才能穿透空泡界面, 形成通气形态.
在液滴通气形态下, 提取了空泡半径和扰动幅值作为关键特征量开展分析. 图6(a)展示了液滴通气特征形态下(Rd0 = 2.5 mm, μ = 0.001 Pa·s), 液体黏性对空泡半径的影响. 观察发现, 随着液体黏性增大空泡在膨胀阶段行为几乎相同, 这是由于膨胀阶段空泡膨胀会挤压液滴, 在液滴内部形成高压一起向外膨胀, 此时液体黏性并不起主要作用, 液体黏性的影响主要体现在收缩阶段. 在收缩阶段, Rd0 = 2.5 mm, μ = 0.001 Pa·s条件下空泡的最小半径为0.69 mm, 当液体黏性增大至μ = 0.05, 0.1 Pa·s时最小半径分别增大了179%和238%. 收缩阶段, 空泡内部低压会在外部大气压的作用下收缩, 相比与高温高压的初始空泡此时大气压的压差并不大, 液体黏性作用逐渐显著, 形成一个与压差方向相反的力, 因此空泡在液体黏性的作用下会更难收缩. 液体黏性增加使空泡收缩减慢, 最小半径变大.
图6(b)中展示了液滴通气特征形态下(Rd0 = 2.5 mm, μ = 0.001 Pa·s), 液体黏性对扰动幅值发展历程的影响, 定义扰动幅值为η(t), 其通过当前时刻液滴界面扰动的最大位置和最小位置计算得到
$$ \eta (t) = 0.5({R_{\eta \max }}(t) - {R_{\eta \min }}(t)) $$ (9) 观察发现, 扰动幅值随着发展过程逐渐增大, 但在膨胀阶段和收缩阶段行为有所差异. 在膨胀阶段(t < 0.25 ms)扰动幅值较小并且随着空泡膨胀缓慢增长, 明显受到液体黏性抑制, 液体黏性较小时(μ = 0.001 Pa·s)更早观察到扰动. 收缩阶段(t > 0.25 ms)开始后, 扰动随着空泡收缩快速增长, 这是由于空泡收缩带动液滴内界面向内运动造成的. 当扰动增长至穿透空泡界面, 液滴发展到通气形态. 随着液体黏性增加扰动的整体幅值均减小, 但整体趋势仍保持一致. 因此可得, 液体黏性增大扰动幅值会逐渐变小, 但并不改变扰动的发展趋势.
图7展示了液体黏性对不同液滴初始半径通气发生时间(tv)的影响. 由图7可知, 液滴初始半径相同时, 通气发生时间随着液体黏性增大逐渐增加, 对于Rd0 = 2.5 mm的工况从μ = 0.001 ~ 0.1 Pa·s增加了15.9%. 这是由于液体黏性增大扰动幅值变小, 其需要更长时间的增大以形成通气形态. 相比之下, 液滴初始半径的增加也会增加通气开始的时间, 并且其效果比液体黏性的更加明显. 对于μ = 0.1 Pa·s和Rd0 = 1.8 mm的工况, 发生通气仅需要0.22 ms而Rd0 = 2.5 mm则需要0.39 ms, 增加了77.27%. 液滴初始半径增大达到通气形态需要的扰动幅值更大, 因此扰动需要更长的发展时间, 对于继续加大液滴初始半径则液滴形态会变成稳定, 此时扰动在最大值时无法发生通气.
2.3 液滴稳定形态与扰动幅值
图8展示了液滴初始半径Rd0 = 6 mm不同液体黏性条件下液滴与空泡的行为, 其中图8(a) ~ 图8(c)分别对应液体黏性为μ = 0.001, 0.05, 0.1 Pa·s, 此时液滴初始半径更大, 该工况下液滴表现为稳定特征形态. 图8(a)观察可得, 随着空泡膨胀, 液滴界面开始向外运动, 在t = 0.3 ms观察到液滴表面初始扰动在产生, 其扰动产生时间相比于通气形态(图5)更晚, 并且此时空泡界面与液滴界面距离更大. 随后在t = 0.45 ms时膨胀到最大半径, 同时扰动逐渐变大. 在空泡收缩过程中液滴界面扰动继续增长, 在t = 0.85 ms空泡收缩到最小半径, 此时扰动并未穿透空泡界面, 空泡内部与外界气体并未形成通路, 该液滴形态成为稳定特征形态.
图8(b)和图8(c)展示了液体黏性变大后的液滴和空泡行为. 随着液体黏性增加, 出现扰动时间同样延后, μ = 0.05, 0.1 Pa·s分别在t = 0.35, 0.4 ms观察到扰动的产生. 空泡在t = 0.45 ms膨胀到最大半径随后开始收缩, 液体黏性增大的条件下, 收缩阶段扰动幅值的增加明显变缓, 在t = 0.85 ms观察到扰动的幅值随着液体黏性增加逐渐变小, 空泡最小半径则变大. 液滴黏性增加, 扰动幅值变小, 液滴保持稳定形态.
在液滴稳定形态下, 同样对关键特征量开展分析. 图9(a)中展示了液滴稳定特征形态下(Rd0 = 6 mm, μ = 0.001 Pa·s), 液体黏性对空泡半径变化的影响. 与通气形态下(图6(a))空泡行为相似, 液体黏性的影响主要发生在收缩阶段. 液体黏性较小时(μ = 0.001 Pa·s)空泡收缩会更快并且半径会更小, 在收缩到最小半径时为0.56 mm. 当液体黏性增大到μ = 0.05, 0.1 Pa·s时, 其最小半径变为了1.1和1.9 mm, 分别增大了96%和239%. 相比于液滴通气形态, 该工况下空泡周期更长, 这是由于液滴变大, 空泡膨胀需要驱动水的质量增加导致膨胀和收缩都变得更慢.
图9(b)中展示了液滴稳定特征形态下(Rd0 = 6 mm, μ = 0.001 Pa·s), 液体黏性对扰动幅值发展历程的影响. 在膨胀阶段(t < 0.45 ms)液体黏性较小时(μ = 0.001 Pa·s)更早观察到扰动, 此后扰动幅值随着空泡膨胀缓慢增加. 收缩阶段(t > 0.45 ms)开始后, 扰动随着空泡收缩快速增长, 并在t = 0.85 ms到达最大值. 与通气形态下(图6(b))液滴和空泡行为不同, 该条件下液滴初始半径更大不再发生通气, 扰动幅值会随着空泡第2次膨胀而减小, 从而出现先增加后减小的趋势. 此时扰动存在最大幅值, 在液体黏性小时(μ = 0.001 Pa·s)可达到1.85 mm.
图10展示了液滴稳定特征形态下, 液体黏性对不同液滴初始半径最大扰动幅值(ηmax)的影响. 由图可知, 当液滴初始半径相同时, 扰动最大值随着液体黏性增大呈现出逐渐减小的趋势. 对于Rd0 = 4 mm的工况, μ = 0.001 Pa·s时ηmax = 5.36 mm而在μ = 0.1 Pa·s时仅为ηmax = 1.85 mm, 减小了66.2%. 这也是由于液体黏性增大后空泡收缩变慢最小半径变大, 更缓慢的变化导致了不稳定性扰动增长变慢, 最大值变小. 对于不同液滴初始半径其基本趋势保持一致, 但随着液滴初始半径增大, 其扰动幅值也减小, 对于μ = 0.001 Pa·s的工况, Rd0 = 7 mm时仅为ηmax = 2.23 mm, 减小了58.4%. 这是由于液滴初始半径增大后空泡膨胀和收缩变得更加缓慢, 因此扰动的增长会减慢. 相比之下, 从Rd0 = 4 mm增大到Rd0 = 5 mm最大扰动幅值平均减少了29.3%, 而从Rd0 = 6 mm增大到Rd0 = 7 mm最大扰动幅值平均仅减少了15.0%, 证明该效果在液滴初始半径较小时更明显.
2.4 液体黏性影响下的液滴界面不稳定性相图
基于液体黏性和液滴初始半径对液滴内空泡诱导不稳定性形态的影响, 本小节在μ∈[0.001, 100] Pa·s, Rd0∈[1.5, 7] mm两个参数区间上开展了一系列的数值模拟, 并归纳得到液滴不稳定性相图, 如图11所示, 其中蓝色为破碎形态, 橙色为通气形态, 绿色为稳定形态, 图11中蓝色圆形、橙色方形和绿色三角形为数值模拟得到的结果. 由图11可得, 在本文的空泡参数下, 小液滴(Rd0 < 2 mm)条件下, 液体黏性较小时(μ = 0.001 Pa·s)通气发生在膨胀阶段, 此时液滴为破碎形态. 随着液体黏性增加, 扰动幅值减小, 通气发生时间延后. 在μ > 0.1 Pa·s时, 空泡与外界大气连通延后至收缩阶段, 液滴形态从破碎变为通气(如μ = 0.1, 1 Pa·s). 随着液体黏性持续增大(μ = 10 Pa·s), 扰动减小至不能穿透空泡界面, 液滴变为稳定形态. Rd0 = 1.5 mm的典型液体黏性的液滴形态在下方黑色方块中展示, 随着液体黏性增加, 液滴形态从破碎转变为通气再转变为稳定. 中液滴(2 mm < Rd0 < 3 mm)的特征形态受液体黏性影响, 在液体黏性较小时(μ = 0.001, 0.1 Pa·s)扰动随着发展穿透空泡界面, 液滴为通气形态, 随着液体黏性增大最大扰动幅值减小至小于空泡界面和液滴界面的距离, 在μ > 2 Pa·s时不再通气, 液滴从通气转变为稳定形态(μ = 1, 10 Pa·s). Rd0 = 2.5 mm的典型液体黏性的液滴形态在下方紫色方块中展示, 液体黏性增加液滴形态从通气转变为稳定. 随着液滴初始半径变大液滴变得更稳定, 大液滴(Rd0 > 3 mm)的特征形态不再随液体黏性发生变化, 液体黏性增加最大扰动幅值减小, 液滴维持稳定形态. 综上, 液滴初始半径增大和液体黏性增大均能使液滴更加稳定, 液滴会从破碎变为通气, 进而变成稳定形态.
3. 结 论
在内空泡振荡下液滴会发生破碎、通气及保持稳定3种特征形态, 本文基于开源OpenFOAM框架, 采用直接数值模拟和isoAdvector几何VoF方法对液滴界面不稳定性开展研究, 讨论了不同液滴初始半径下液体黏性对液滴界面扰动、通气时间等关键特征量的影响规律, 总结了液体黏性与液滴初始半径参数下的液滴不稳定性相图.
(1)液体黏性同时影响空泡振荡和液滴界面扰动演化. 对于空泡发展, 液体黏性主要影响空泡的收缩运动, 液体黏性的增加会使空泡收缩变得更慢, 且最小半径会更大; 对于扰动发展, 液体黏性会抑制扰动前期增长, 随着液体黏性增加扰动幅值减少, 但扰动的发展趋势不发生改变.
(2)液体黏性和液滴初始半径共同影响了该过程的关键特征量趋势. 对于通气形态, 随空泡发展扰动幅值增大至穿透空泡界面发生通气, 液体黏性或液滴初始半径增加, 通气发生时间增加; 对于稳定形态, 扰动幅值呈现先增加后减小趋势, 随着液体黏性或液滴初始半径增加, 最大扰动幅值减小.
(3)液体黏性和液滴初始半径显著改变柱状液滴的特征形态. 在文章讨论的参数范围内, 随着液滴初始半径增加, 液滴形态从破碎先转变为通气再转变为稳定. 随着黏性增加, 小液滴(Rd0 < 2 mm)的形态从破碎转变为通气进而变成稳定; 中液滴(2 mm < Rd0 < 3 mm)的形态从通气转变为稳定, 不出现破碎形态; 大液滴(Rd0 > 3 mm)的形态不随液体黏性改变, 保持稳定形态. 液滴形态会随着液体黏性或初始半径增加趋于稳定.
本文重点分析了液滴初始半径和液体黏性对液滴界面不稳定性发展的影响, 在今后研究中, 计划基于稳定性分析方法从理论上开展系统参数研究, 推导得出液体黏性作用下的扰动模数及其增长特性.
-
表 1 粗中细网格参数细节
Table 1 Details for coarse, medium and fine grids
Exp. Coarse Medium Fine cells/103 — 90 360 1440 Δx/μm — 40 20 10 Rbmax/mm 5.378 5.346 5.355 5.358 error Rbmax — 0.56% 0.43% 0.37% wave number n — 13 15 15 ηmax/mm — 2.511 1.849 1.729 -
[1] 丘润荻, 王静竹, 黄仁芳等. 改进的物理融合神经网络在瑞利-泰勒不稳定性问题中的应用. 力学学报, 2022, 54(8): 2224-2234 (Qiu Rundi, Wang Jingzhu, Huang Renfang, Du Tezhuan, Wang Yiwei, Huang Chenguang. The application of modified physicsinformed neural networks in Rayleigh-Taylor instability. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(8): 2224-2234 (in Chinese) Qiu Rundi, Wang Jingzhu, Huang Renfang, Du Tezhuan, Wang Yiwei, Huang Chenguang. The application of modified physicsinformed neural networks in Rayleigh-Taylor instability. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(8): 2224-2234 (in Chinese)
[2] 陈璐, 赖惠林, 林传栋等. 多模态瑞利-泰勒不稳定性的离散玻尔兹曼数值研究. 空气动力学学报, 2022, 40(3): 140-150 (Chen Lu, Lau Huilin, Lin Chuandong, et al. Numerical study of multimode Rayleigh-Taylor instability by using the discrete Boltzmann method. Acta Aerodynamica Sinica, 2022, 40(3): 140-150 (in Chinese) Chen Lu, Lau Huilin, Lin Chuandong, et al . Numerical study of multimode Rayleigh-Taylor instability by using the discrete Boltzmann method. Acta Aerodynamica Sinica,2022 ,40 (3 ):140 -150 (in Chinese)[3] Zhang AM, Li SM, Cui P, et al. A unified theory for bubble dynamics. Physics of Fluids, 2023, 35(3): 033323 doi: 10.1063/5.0145415
[4] Jia BQ, Fu QF, Xu X, et al. Spray characteristics of Al-nanoparticle-containing nanofluid fuel in a self-excited oscillation injector. Fuel, 2021, 290: 120057J doi: 10.1016/j.fuel.2020.120057
[5] Wang J, Wang G, Zeng Q, et al. Recent progress on the jetting of single deformed cavitation bubbles near boundaries. Journal of Hydrodynamics, 2023, 35(5): 832-857 doi: 10.1007/s42241-023-0071-6
[6] Fang L, Ma K, Liu C, et al. Ocean mixing induced by three-dimensional structure of equatorial mode tropical instability waves in the Pacific Ocean. Frontiers in Marine Science, 2023, 10: 1228897 doi: 10.3389/fmars.2023.1228897
[7] Augustin D, César MB. Gravity waves’ modulational instability under the effect of drag coefficient in the ocean. Physica Scripta, 2023, 98(12): 125014 doi: 10.1088/1402-4896/ad05a9
[8] Tupper A, Textor C, Herzog M, et al. Tall clouds from small eruptions: the sensitivity of eruption height and fine ash content to tropospheric instability. Natural Hazards, 2009, 51(2): 375-401 doi: 10.1007/s11069-009-9433-9
[9] Fowler SJ, Spera FJ. Phase equilibria trigger for explosive volcanic eruptions. Geophysical Research Letters, 2008, 35(8): 2008GL033665 doi: 10.1029/2008GL033665
[10] 王立锋, 叶文华, 陈竹等. 激光聚变内爆流体不稳定性基础问题研究进展. 强激光与粒子束, 2021, 33(01): 5-64 (Wang Lifeng, Ye Wenhua, Chen Zhu, et al. Review of hydrodynamic instabilities in inertial confinement fusion implosions. High Power Laser and Particle Beam, 2021, 33(01): 5-64 (in Chinese) doi: 10.11884/HPLPB202133.200173 Wang Lifeng, Ye Wenhua, Chen Zhu, et al . Review of hydrodynamic instabilities in inertial confinement fusion implosions. High Power Laser and Particle Beam,2021 ,33 (01 ):5 -64 (in Chinese) doi: 10.11884/HPLPB202133.200173[11] 林祖德, 戴羽, 徐梦飞等. 基于高精度3D打印工艺的ICF调制靶. 强激光与粒子束, 2023, 35(10): 64-70 (Lin Zude, Dai Yu, Xu Mengfei, et al. ICF modulation targets based on high-precision 3D printing technology. High Power Laser and Particle Beam, 2023, 35(10): 64-70 (in Chinese) Lin Zude, Dai Yu, Xu Mengfei, et al . ICF modulation targets based on high-precision 3D printing technology. High Power Laser and Particle Beam,2023 ,35 (10 ):64 -70 (in Chinese)[12] Ligrani PM, Mcnabb ES, Collopy H, et al. Recent investigations of shock wave effects and interactions. Advances in Aerodynamics, 2020, 2(1): 4 doi: 10.1186/s42774-020-0028-1
[13] Taylor G. The Instability of Liquid Surfaces when Accelerated in a Direction Perpendicular to their Planes. I. Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, 1950, 201(1065): 192-196
[14] Lewis DJ. The instability of liquid surfaces when accelerated in a direction perpendicular to their planes. II. Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences, 1950, 202(1068): 81-96
[15] Guo HY, Cheng T, Li YJ. Weakly nonlinear multi-mode Bell–Plesset growth in cylindrical geometry. Chinese Physics B, 2020, 29(11): 115202 doi: 10.1088/1674-1056/ab9c14
[16] Birkhoff G. Stability of spherical bubbles. Quarterly of Applied Mathematics, 1956, 13(4): 451-453 doi: 10.1090/qam/79932
[17] Bell GI. Taylor instability on cylinders and spheres in the small amplitude approximation. Los Alamos National Laboratory, 1951, Rep. LA-1321
[18] Plesset MS. On the stability of fluid flows with spherical symmetry. Journal of Applied Physics, 1954, 25(1): 96-98 doi: 10.1063/1.1721529
[19] 南江浪, 张政, 姚卫等. 分区参数对超声速湍流燃烧动态分区火焰面模拟的影响. 力学学报, 2023, 56(3): 1-11 (Nan Jianglang, Zhang Zheng, Yao Wei, et al. Effect of zoning parameters on dynamic zone flamelet modeling of supersonic turbulent combustion. Chinese Journal of Theoretical and Applied Mechanics, 2024, 56(3): 1-11 (in Chinese) Nan Jianglang, Zhang Zheng, Yao Wei, et al. Effect of zoning parameters on dynamic zone flamelet modeling of supersonic turbulent combustion. Chinese Journal of Theoretical and Applied Mechanics, 2024, 56(3): 1-11 (in Chinese)
[20] 吕可, 邹佳俊, 陈颖等. 近刚性壁面异相双空泡耦合及射流增强效应研究. 力学学报, 2023, 55(8): 1605-1617 (Lyu Ke, Zou Jiajun, Chen Ying, Yan Shuai, Li Shuai, Zhang Aman. Study on the interaction and jet enhancement effect of two out-of phase bubbles near a rigid boundary. Chinese Journal of Theoretical and Applied Mechanics, 2023, 55(8): 1605-1617 (in Chinese) Lyu Ke, Zou Jiajun, Chen Ying, Yan Shuai, Li Shuai, Zhang Aman . Study on the interaction and jet enhancement effect of two out-of phase bubbles near a rigid boundary. Chinese Journal of Theoretical and Applied Mechanics,2023 ,55 (8 ):1605 -1617 (in Chinese)[21] 王平平, 张阿漫, 彭玉祥等. 近场水下爆炸瞬态强非线性流固耦合无网格数值模拟研究. 力学学报, 2022, 54(8): 2194-2209 (Wang Pingping, Zhang Aman, Peng Yuxiang, Meng Zifei. Numerical simulation of transient strongly-nonlinear fluid-structure interaction in near-field underwater explosion based on meshless method. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(8): 2194-2209 (in Chinese) Wang Pingping, Zhang Aman, Peng Yuxiang, Meng Zifei . Numerical simulation of transient strongly-nonlinear fluid-structure interaction in near-field underwater explosion based on meshless method. Chinese Journal of Theoretical and Applied Mechanics,2022 ,54 (8 ):2194 -2209 (in Chinese)[22] Livescu S, Roy RV, Schwartz LW. Leveling of thixotropic liquids. Journal of Non-Newtonian Fluid Mechanics, 2011, 166(7-8): 395-403 doi: 10.1016/j.jnnfm.2011.01.010
[23] Obreschkow D, Kobel P, Dorsaz N, et al. Cavitation bubble dynamics inside liquid drops in microgravity. Physical Review Letters, 2006, 97(9): 094502 doi: 10.1103/PhysRevLett.97.094502
[24] 吕明, 宁智, 孙春华. 单液滴内空化空泡的生长及溃灭研究. 力学学报, 2016, 48(4): 857-866 (Lyu Ming, Ning Zhi, Sun Chunhua. Study on the growth and collapse of cavitation bubble within a droplet. Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(4): 857-866 (in Chinese) doi: 10.6052/0459-1879-15-434 Lyu Ming, Ning Zhi, Sun Chunhua . Study on the growth and collapse of cavitation bubble within a droplet. Chinese Journal of Theoretical and Applied Mechanics,2016 ,48 (4 ):857 -866 (in Chinese) doi: 10.6052/0459-1879-15-434[25] Zeng Q, Gonzalez-avila SR, Voorde ST, et al. Jetting of viscous droplets from cavitation-induced Rayleigh–Taylor instability. Journal of Fluid Mechanics, 2018, 846: 916-943 doi: 10.1017/jfm.2018.284
[26] Wang J, Li H, Guo W, et al. Rayleigh–Taylor instability of cylindrical water droplet induced by laser-produced cavitation bubble. Journal of Fluid Mechanics, 2021, 919: A42 doi: 10.1017/jfm.2021.401
[27] Zhang Y, Zhang X, Zhang S, et al. A review of the dynamics progress of bubble collapse within droplet and droplet splash. Applied Sciences, 2023, 13: 7822 doi: 10.3390/app13137822
[28] Rosselló JM, Reese H, Raman KA, et al. Bubble nucleation and jetting inside a millimetric droplet. Journal of Fluid Mechanics, 2023, 968: A19 doi: 10.1017/jfm.2023.542
[29] Suponitsky V, Plant D, Avital EJ, et al. Pressure wave in liquid generated by pneumatic pistons and its interaction with a free surface. International Journal of Applied Mechanics, 2017, 9(03): 1750037 doi: 10.1142/S1758825117500375
[30] Huang J, Wang J, Huang J, et al. Effects of wall wettability on vortex flows induced by collapses of cavitation bubbles: A numerical study. Physics of Fluids, 2023, 35(8): 087122 doi: 10.1063/5.0164694
[31] Miller ST, Jasak H, Boger DA, et al. A pressure-based, compressible, two-phase flow finite volume method for underwater explosions. Computers & Fluids, 2013, 87: 132-143
[32] Koch M, Lechner C, Reuter F, et al. Numerical modeling of laser generated cavitation bubbles with the finite volume and volume of fluid method, using OpenFOAM. Computers & Fluids, 2016, 126: 71-90
[33] Roenby J, Bredmose H, Jasak H. A computational method for sharp interface advection. Royal Society Open Science, 2016, 3(11): 160405 doi: 10.1098/rsos.160405
[34] Ye B, Wang Y, Huang C, et al. Numerical study of the pressure wave-induced shedding mechanism in the cavitating flow around an axisymmetric projectile via a compressible multiphase solver. Ocean Engineering, 2019, 187: 106179 doi: 10.1016/j.oceaneng.2019.106179