NUMERICAL STUDY OF THE INTERACTION BETWEEN SHOCK WAVES AND STAGGERED LIQUID COLUMNS
-
摘要: 采用时空守恒元与解元方法对平面激波冲击交错排列的3液柱形态演化过程进行研究. 展示了平面激波冲击3液柱后的波系结构及界面形态演化过程, 发现液柱间的通道增强了激波间的相互作用, 使下游液柱迎风面上压力分布不均, 导致其迎风面在较短时间内从弧面演化为平坦表面. 测量了下游液柱的迎风面极点位移和宽度, 观察到上游两液柱间通道的变窄趋势使得气流对下游液柱迎风面极点的冲击减弱, 降低了其迎风面极点位移速度. 通过测量3液柱整体的无量纲宽度、高度和质心位移, 发现第二通道宽度在气流与界面的相互作用中发挥了关键作用. 随着第二通道宽度减小, 气流对界面的作用能力增强, 上游液柱后部凹陷的演化速率加快, 第一通道侧的尖端逐渐消失, 3液柱整体的无量纲宽度由减小到增大的转折点提前, 而后期的无量纲质心位移速度减慢. 在第二通道宽度较大的算例中, 上游液柱仍保持了较为对称的结构特征. 而在第二通道宽度较小的算例中, 上游液柱在波系演化阶段便显现出显著的非对称性, 不再展现出对称的结构特征.
-
关键词:
- 激波 /
- 交错排列液柱 /
- 波系演化 /
- 界面变形 /
- 时空守恒元与解元方法
Abstract: The morphology of three liquid columns impacted by a plane shock is investigated by using the space-time conservation element and solution element method. The evolution of the shock system and interface deformation of the three liquid columns are demonstrated. The flow inside the channels between the liquid columns enhance the interaction between the shock waves. It leads to uneven pressure distributions on the windward surface of the downstream liquid column. And it also causes its windward surface to evolve from a curved surface to a flat one within a short period of time. The windward pole displacement and width of the downstream liquid column are measured. The narrowing trend of the channel between the two upstream liquid columns weakens the impact of the air flow on the windward pole of the downstream liquid column. It also reduces the velocity of the windward pole displacement. The dimensionless widths, heights, and center-of-mass displacements of the entire three liquid columns are measured. It is found that the width of the second channel plays a key role in the interaction between the airflow and the liquid interface. As the width of the second channel decreases, the ability of the airflow to interact with the interface strengthens and the deformation rate of the dent at the back of the upstream liquid column accelerates. The tip on the first channel side disappears. Then, the turning point of the dimensionless width of the three liquid columns from decreasing to increasing occurs earlier. Additionally, the dimensionless centroid displacement slows down at a later stage. In the case of large second channel width, the upstream liquid column maintains a symmetric shape. However, in the case of small second channel width, the upstream liquid column shows significant asymmetry at the wave system evolution stage and no longer exhibits symmetric features. -
引 言
燃料射流在超声速气流中的破碎行为受到广泛关注. 射流的破碎时间、破碎后燃料液滴的尺寸分布等参数对发动机的工作效率具有重要影响[1]. 在爆炸防护实验中, 多液柱形成的水幕帘能有效阻挡爆炸冲击波, 通过激波的反射和透射作用显著降低其能量强度, 从而有效减轻爆炸造成的伤亡和损失[2]. 激波与液柱相互作用过程中, Kelvin-Helmholtz (KH)不稳定性[3]和Rayleigh-Taylor (RT)不稳定性[4]起着关键作用, 液柱的变形是这两种不稳定性共同作用的结果. 深入研究液滴、液柱破碎现象有助于加深对超声速气流中液体变形行为的理解.
关于液滴破碎的研究已取得一系列研究成果. 早期的研究重点是对单个液滴破碎形态和破碎机制进行归纳总结. Hinze[5]最早提出We数是表征液滴破碎的主要参数并提出了临界We数的概念. Pilch等[6]则根据We数的不同对液滴的破碎模式进行分类: 振荡破碎[7](vibrational breakup, We≤12)、袋状破碎[8](bag breakup, 12≤We≤50)、袋蕊破碎[9](bag-stamen breakup, 50≤We≤100)、剪切破碎[10](sheer breakup, 100≤We≤350)和毁灭性破碎[11-12](catastrophic breakup, We≥350). Theofanous等[11]利用激光诱导荧光技术(laser induced fluorescence, LIF)获得了清晰的图像, 并揭示了液滴破碎的两个基本机制: We数较小时, 由于RT不稳定性的发展, 使得液滴发生变形破碎, 其对应机制为RT穿刺破碎(Rayleigh-Taylor piercing, RTP); We数较高时液滴表面的液雾多为边缘剪切生成, 对应机制是剪切诱导的剥离破碎(shear induced entrainment, SIE). Thoufanous等[11-13]对毁灭性破碎的真实性提出质疑, 认为该模式并非是实际的物理现象. 他们认为之前的实验忽略了液滴破碎过程中气体出现的高温现象且纹影图中液体突然消失时观察到的波的波长比相应加速度下增长最快的RT波波长高一个数量级. Nicholls等[14]和Liu等[15]针对液滴破碎过程中其表面产生的液雾分别提出了“剪切剥离”(shear stripping)机制和薄层细化机制(sheet-thinning). Guildenbecher等[16]讨论了大量关于液滴二次破碎在实验方法、破碎形态、破碎时间、破碎液滴尺寸和速度分布的相关研究成果, 通过对比分析, 阐明了液滴破碎的物理机制.
传统的实验方法在研究液滴或液柱破碎过程时, 需要配备高精度的高速摄影设备, 以确保能够清晰捕捉液滴破碎的细节. 随着计算流体力学的不断发展, 越来越多的学者采用数值模拟的方式研究此类问题. 相较于传统实验, 数值模拟能够提供更全面的数据信息且不受实验设备或环境限制, 方便进行参数研究, 进而更全面地探究不同因素对液滴破碎的影响. Meng等[17]采用高保真的数值方法研究了水柱在激波冲击后流动初期的破碎现象, 研究结果证实了先前实验中观察到的剪切剥离破碎模式. 宋家喜等[18]对高马赫数和极高马赫数下液滴变形破碎的动力学过程进行全三维数值模拟研究, 揭示了液滴破碎过程中的循环破碎机制. Xu等[19]采用理论射线分析和高分辨率数值模拟的方法, 研究了不同弯曲激波作用下二维水柱内的波汇聚现象, 建立了描述激波在水柱内部波系结构演化特征的解析方法. Ling等[20]通过数值模拟研究了中等We数下的液滴破碎现象, 讨论了液滴变形为前向袋状过程中的不同阶段和主导机制.
目前的研究已不再局限于波后气流与单个液滴或液柱的相互作用, 而是结合实际情况进行数值模拟以期得到全面深入的理解. Xu等[21]研究了含颗粒液滴在气流中的破碎行为, 发现颗粒的存在会加速液膜的破碎, 并促进局部快速穿透, 从而影响液滴的破碎模式. 而在实际应用中, 液滴和液柱的分布也更为复杂. Igra等[22-23]强调了对二维液柱研究的重要性, 研究了激波冲击下双液柱的形态演化过程并指出上游的水柱会影响下游的流场, 从而影响下游水柱的变形. Chen等[24]通过数值模拟研究了平面激波与两个水柱相互作用的复杂过程, 分别讨论了水平排列和垂直排列两种配置的影响, 发现垂直排列的水柱间形成的高速区域对漩涡的产生和演化起着关键作用. 吴润龙等[25]利用守恒清晰界面方法深入研究了激波作用下并排液滴变形的规律, 发现液滴间距与液滴间通道内的压力峰值有关. 随着并排液滴间距离的缩短, 通道侧界面波动的幅度增大, 这表明液滴间距对并排液滴尤其是在通道侧界面形态的演变过程中的影响至关重要. 研究多液滴或液柱在激波冲击下的波系及流动界面演化, 对于深入理解其背后的物理机制具有重要意义.
本文采用数值方法研究了平面激波冲击交错排列3液柱的物理过程, 分析了该过程中的波系结构演化和液柱的形态演化, 研究了3液柱之间的通道对下游液柱界面变形形态的影响, 探讨了上游液柱与下游液柱之间的通道宽度对3液柱整体变形形态的影响及作用机理.
1. 计算域设置与数值方法
1.1 计算域设置
本文计算域设置如图1所示.
激波马赫数为1.935, 波前空气处于静止状态, 密度为1.174 kg/m3. 3个液柱A, B和C处于静止状态, 其初始半径R均为
0.00175 m, 密度为1000kg/m3, 液柱位置如图1所示. 本文研究的模型是平面对称模型, 为节省计算资源, 计算域设置在对称线之上, 其大小为16R × 8R. 参照Guan等[26]的研究, 文中液滴半径上的网格数为150. 计算域左侧边界为入口边界, 右侧边界为出口边界, 顶部边界设为开放边界, 底部边界则设为对称边界.1.2 数值方法
入射激波波后气流流速高, 气动力在液柱破碎过程中占主导地位, 因此忽略了液柱的表面张力及黏性. 本文用于研究激波与液柱相互作用的控制方程为Euler方程, 方程的无量纲形式为
$$ \frac{{\partial \left( {{y^n}{\boldsymbol{U}}} \right)}}{{\partial t}} + \frac{{\partial \left( {{y^n}{\boldsymbol{F}}} \right)}}{{\partial x}} + \frac{{\partial \left( {{y^n}{\boldsymbol{G}}} \right)}}{{\partial y}} = {\boldsymbol{0}} $$ (1) 对于二维问题n = 0. 控制方程中, U表示守恒量, F和G表示相对应的对流通量矢量为
$$ {\boldsymbol{U}} = \left( \begin{gathered} \rho \\ \rho u \\ \rho v \\ E \\ \end{gathered} \right),\quad {\boldsymbol{F}} = \left( \begin{gathered} \rho \\ \rho {u^2} + p \\ puv \\ \left( {E + p} \right)u \\ \end{gathered} \right),\quad {\boldsymbol{G}} = \left( \begin{gathered} \rho \\ \rho uv \\ p{v^2} + p \\ \left( {E + p} \right)v \\ \end{gathered} \right) $$ (2) 其中, ρ表示流体密度, u和v分别为x和y方向上流体的速度, p表示流体压强, E表示流体总能量.
使用基于体积分数的5方程模型来描述气液界面的演化, 表达式为
$$ \frac{{\partial \alpha }}{{\partial t}} + {\boldsymbol{V}} \cdot \nabla \alpha = 0 $$ (3) $$ \frac{{\partial {\rho _s}{\alpha _s}}}{{\partial t}} + \nabla \cdot \left( {{\rho _s}{\alpha _s}{\boldsymbol{V}}} \right) = 0,\quad s = 1,2 $$ (4) $$ \frac{{\partial \rho {\boldsymbol{V}}}}{{\partial t}} + \nabla \cdot \left( {\rho {\boldsymbol{V}} \otimes {\boldsymbol{V}} + p} \right) = {\boldsymbol{0}} $$ (5) $$ \frac{{\partial E}}{{\partial t}} + \nabla \cdot \left[ {{\boldsymbol{V}}\left( {E + p} \right)} \right] = 0 $$ (6) 式中, α表示流体的体积分数, ρs表示流体$s$的密度, ρ表示混合物的密度, V表示速度矢量.
采用加强气体状态方程来封闭系统
$$ p = \left( {\gamma - 1} \right)\left( {E - \frac{1}{2}\rho {\boldsymbol{V}} \cdot {\boldsymbol{V}}} \right) - \gamma \textit{π} $$ (7) 其中
$$ \frac{1}{{\gamma - 1}} = \sum {\frac{{{\alpha _i}}}{{{\gamma _i} - 1}}} ,\quad \frac{{\gamma \textit{π} }}{{\gamma - 1}} = \sum {\frac{{{\alpha _i}{\gamma _i}{\textit{π} _i}}}{{{\gamma _i} - 1}}} $$ (8) 混合物的总密度和声速分别为
$$ \rho = \sum {{\alpha _i}{\rho _i},\quad c = \sqrt {\gamma \left( {p + \textit{π} } \right)/\rho } } $$ (9) 对于空气, 比热比γ = 1.4, π = 0 Pa使加强气体状态方程降为理想气体方程. 对于水, 采用经验参数γ = 7.15和π = 3.31 × 108 Pa.
本文采用二维结构网格对计算域进行划分, 使用满足极大值准则的迎风时空守恒元与解元(conservation element and solution element, CESE)方法求解方程. CESE方法基于双曲微分方程的时空积分形式, 将空间和时间统一对待. Shen 等[27]将其改进为迎风格式. 该格式保留了中心格式低耗散的优势, 同时确保了局部和全局的守恒性, 有效地抑制了数值震荡, 并且其耗散程度不受 CFL数的影响. 对网格划分以及CESE方法的更多描述见参考文献[26-30]. 数值程序的正确性也已在前期研究中得到验证, 如图2所示, Guan等[26]使用该数值程序模拟了马赫数为1.39的激波冲击液滴直径为3.03mm液滴的过程, 将数值结果与Yi等[31]的实验结果进行对比后发现, 该程序对液滴形态演化和雾化(AT)区域的刻画与实验结果基本吻合.
![]()
2. 数值结果分析
2.1 初期波系结构演化
以入射激波IS与液柱A和B接触时刻为0时刻. 在1.74 ${\text{μs}}$时, 液柱表面产生了向上游传播的反射激波RS1, 液柱内部产生向下游传播的透射激波TS1, 如图3(a)所示. 由于气-液两相声阻抗不同, 透射激波TS1将先与液柱下游壁面作用产生内部反射膨胀波, 膨胀波汇聚处产生负压区[32]. IS与RS1在液柱边界连接, 当入射激波与界面夹角由小于脱离角变为大于脱离角时, 规则反射RR转变为马赫反射MR, 产生马赫杆MS1, 且IS, RS1和MS1交于3波点TP1, 如图3(b)所示.
![]()
图 3 初期波系结构演化纹影图(上)和压力云图(下)Figure 3. Schlieren (top) and pressure (bottom) contours of the wave system evolution at the early stage当IS经过液柱A和B赤道时, 液柱A和B之间形成的通道(后文称为第一通道)内部RS1相互作用发生反射使上游的入射激波IS弯曲并与MS1融合为弯曲激波BS1, RS1在通道内部与液柱经过多次作用, 形成向上游传播的“鱼鳞状”波系结构(fish scale wave), 并产生马赫杆MS2. 如图3(d)所示, 弯曲激波BS1与液柱C相互作用, 产生反射激波RS2, 液柱内部产生透射激波TS2. 由于液柱A和B的存在, 气流经第一通道沿第二通道(液柱A与C和液柱B与C之间形成的通道)向下游传播, RS2波阵面逐渐趋于平缓. BS1与MS1相互作用产生反射激波RS3和RS4, RS2与RS3作用产生两道反射激波, 其中一道与RS4融合形成弯曲激波BS2, 另一道则与液柱C发生碰撞, 流场中出现马赫杆MS3, 如图3(e)所示. 绕射激波DS1相交, 产生马赫杆MS4.
2.2 液柱形态演化
在激波冲击交错布置的3液柱的波系演化过程中, 液柱周围流场在不断变化. 如图3(e)和图3(f)所示, 液柱A和B整体变形幅度较小, 液柱C阻碍了从第一通道流出的气流, 导致第一通道出口形成激波. 液柱A和B界面上均出现不同程度的波动, 其中赤道面非通道侧较通道侧的不稳定波发展迅速, 表现为非通道侧迎风面赤道附近凹坑更显著. 液柱A和B后部通道侧分别受到来自第二通道的气流作用剥离出大量液雾并沿通道向非通道侧移动. 气流夹带液雾向非通道侧移动使得液柱A和B背风面附近的流场演化更为复杂, 背风面在气流作用下发生变形, 变形后的界面又反过来影响流场. 在气流与界面相互作用的过程中, 液柱A和B界面形态的不对称演化趋势显现.
图3(e)和图3(f), 红色框线所示区域中的气流由于经过了一道强激波, 导致该区域压强较高. 因此液柱C迎风面在内外压差的作用下, 呈现扁平化趋势. 对于液柱C, 除了迎风面的扁平化趋势更显著以外, 其赤道侧及背风面也发生了微小形变. 图3(f)和图4(a)中可以观察到, 液柱C赤道侧流场变化迅速, 在气流作用下, 该区域的界面也出现略微的扁平化趋势. 而液柱C背风面的流场则较为稳定, 在两个漩涡的作用下, 背风面极点处的液体向两侧延伸, 使之逐渐变为平坦表面.
![]()
图 4 液柱破碎中期压力图(上)和涡量云图(下)Figure 4. Pressure (top) and vorticity (bottom) contours at the middle stage of liquid column breakup在t = 29.52 ${\text{μs}}$时, 液柱A和B迎风表面极点附近仍表现得较为光滑, 由极点向赤道附近延伸, 出现了小幅度褶皱, 非通道侧形成较为明显的凹坑P1(pit1), 背风面突起T1(tip1)在气流持续作用下向下游方向延伸. 从液柱A和B剥离产生的液雾在气流的裹挟下沿第二通道向下游传播, 液柱C迎风面弧度减小变为平坦表面, 在平坦表面与弧面交界处即液柱C外部高压区域的边缘, 产生了凹坑P2. 由于液柱A和B通道侧和非通道侧的流场不同, 导致液柱两侧涡量分布不均匀. 液柱受到的气流的剪切作用不一致, 使得液柱赤道两侧界面的不对称演化加剧. 从涡量图中可以看出, 以红色点划线(与水平方向的夹角α = 28.90°)为分界线, 液柱B两侧形态相似, 两侧涡量分布也相似(方向相反). 因此, 液柱B的形态演化在气流作用下呈“顺时针扭转”的趋势. 相对应地, 液柱A的形态演化则呈现“逆时针扭转”的趋势.
在t = 61.26 μs时, 液柱A和B的不对称演化趋势更为显著, 液柱迎风面极点处也出现了细碎褶皱, 凹坑P1受气流的剪切作用变得平缓狭长. 液柱非通道侧液柱主体迎风面受到气流的直接冲击, 背风面受到回流的冲击, 在两者的共同作用下液柱出现了平台状结构Pf1 (platform1). 选取液体体积分数αl为0.5的物质线作为液柱主体的轮廓线, 发现该平台结构Pf1两端出现突起. 在液柱内部流动与外部流场的共同作用下, 背风面非通道侧出现大幅度隆起, 出现平台结构Pf2和Pf3. 液柱C两侧出现了尖端T2, T3和T4, 这与单液滴受激波冲击发生变形破碎过程中出现的典型现象相符[33].
在t = 87.26 μs时, 液柱A和B通道侧Z1 (zone1)区域出现凹槽结构, 位于通道侧的液柱外层液体受剪切作用形成不同程度的突起, 位于非通道侧的平台结构Pf1变得细长. 液柱C迎风面在气流作用下由平面演变为凹面. 此时, 向内凹陷的迎风面的导流作用增强, 气流沿第二通道传播, 使得通道两侧界面的KH不稳定波逐渐发展, 通道出口处的尖端T2拉长并沿着第二通道向上游移动. 如图4(d)所示, 在t = 118.71 μs时, 液柱A和B背风面产生的大幅度隆起变得平缓. 液柱C迎风面被拉长, 极点处出现凹坑P3, 其余部分上的不稳定波的波数增加, 背风面逐渐在y轴方向上延伸形成平坦表面. 此时, 液柱A和B主体部分仍保持较好的完整性. 从涡量图中可以看出, 产生的涡量大多集中于第二通道中以及液柱C的赤道侧, 因此以上位置的界面演化更为迅速.
2.3 液柱C的特征尺度
从Igra等[23]的激波冲击双液柱实验中可知, 在上游液柱的影响下, 下游液柱的变形速率较为缓慢. 为了探究激波冲击3液柱的过程中, 液柱之间通道内部流动对下游液柱变形形态的影响, 本节主要对液柱C的迎风面极点的位移s以及宽度w进行定量分析. 根据单液柱的演化图像绘制了迎风面极点的位移s以及宽度w示意图, 见图5.
![]()
图 5 液柱迎风面极点位移与宽度示意图Figure 5. Definition of windward pole point drift and width of the liquid column对液柱C暴露在波后气流的时间进行无量纲化
$$ \tau = t\frac{{{u_g}}}{D}\sqrt {\frac{{{\rho _g}}}{{{\rho _l}}}} $$ 其中, t为弯曲激波BS1接触到液柱C的时间, D为液柱C的初始直径, ug, ρg和ρl分别为入射激波IS的波后气流速度、气流密度和液柱密度.
单液柱与液柱C迎风面极点位移和液柱宽度随时间变化的对比分别如图6(a)和图6(b)所示. 单液柱的破碎过程如图7所示. 根据迎风面极点位移图像将液柱的变形过程分为3个阶段.
![]()
图 6 液柱迎风面极点位移和宽度随时间变化Figure 6. Time variation of windward pole point drift and width of the liquid column将τ < 0.058定义为第一阶段(Ⅰ), 该阶段波系结构尚未发展完全, 液柱C界面变形可以忽略. 液柱C与单液柱的迎风面极点位移s和宽度w在该阶段均未发生显著变化.
通过比较数值结果发现τ = 1.009时, 二者的极点位移相等. 因此, 将0.058 < τ < 1.009定义为第二阶段(Ⅱ), 该阶段液柱C迎风面的演化形态与单液柱情况差异显著. 液柱C迎风面前部存在的高压区域使得迎风面由弧面更快地演化为平坦表面. 在气流的持续作用下, 液柱C的迎风平面逐渐演化为凹面, 并呈现出向来流方向倾斜的趋势, 极点附近受冲击产生凹坑, 如图4(d)所示. 气流从第二通道流出, 增强了对通道两侧界面即液柱A和B的背风面以及液柱C迎风面上的KH不稳定性的发展. 气流在第二通道中剥离迎风面上因KH不稳定性所产生的突起, 还使得液柱C的迎风面被不断拉伸, 导致其迎风直径迅速增大. 在液柱C背风区形成的两个强度较大的漩涡作用下, 见图4(a), 其背风面由弧面变为平面的演化速率较单液柱更为迅速, 液柱C背风面极点向迎风面极点靠近速度更快. 因此该阶段液柱C迎风面极点位移s较单液柱略大, 而宽度w的减小速率明显快于单液柱.
将τ > 1.009定义为第三阶段(Ⅲ), 该阶段液柱C前部的凹陷持续扩大, 第二通道两侧界面逐渐平行. 但第一通道有逐渐变窄的趋势, 使得气流对液柱C迎风面极点的冲击减小, 见图4(d). 因此, 该阶段液柱C迎风面极点位移s的增长速率有所减缓. 由于平台结构Pf1逐渐变得细长, 且液柱C迎风面上的尖端T2被不断拉伸, 第二通道出口处的气流与非通道侧来流相互作用, 使得液柱C背风区的涡量减少, 界面受到的冲击减弱, 从而使液柱宽度w的变化速率减缓.
对单个液柱而言, 波后气流直接与液柱作用, 使液柱发生变形, 但这种变形对气流的影响较本文研究的3液柱工况下小. 3液柱工况下, 气流与通道两侧界面的相互作用更为显著. 随着通道两侧界面演化, 气流在通道内部的流动受到的阻碍减小. 在液柱C迎风面的尖端T2被拉伸至超过赤道侧最高点之前, 如图4(b)所示, 液柱C背风区回流作用均得到加强, 显著改变了位于下游的液柱C的界面变形速率以及界面演化形态.
3. 液柱位置对变形的影响
本节通过调整液柱C距离液柱A和B圆心连线的距离Xc, 在其他参数保持不变的情况下, 对液柱破碎过程进行数值仿真, 探究液柱位置对液柱变形过程的影响. 3个算例的设置见表1. 对数值结果进行初步分析发现改变液柱位置对早期波系演化过程影响较小, 因此本节将主要分析A2和A3工况下波系演化过程中出现的不同于A1工况的现象.
表 1 数值模拟算例的参数设置Table 1. Parameters setup for the numerical simulation casesCase Ms Xc/m A1 1.915 0.004546 A2 1.915 0.003500 A3 1.915 0.002778 液柱C向液柱A和B方向的靠近直接减小了第二通道的宽度, 在A2和A3工况下, 入射激波在第二通道内多次反射, 在第一通道处形成了一道向上游传播的强激波, 流场中形成了面积更大、持续时间更长的高压区域, 如图8(a)和图8(b)所示. 同时, 液柱C向来流方向的靠近使得其内部的高压区域不仅局限于迎风面极点附近, 还使得高压区域在液柱外部向其赤道两侧延伸. 通常, 当激波冲击液柱时, 液柱内部会产生透射激波. 透射激波与液柱下游的壁面作用, 产生向上游传播的膨胀波, 称为第一反射波[19]. 由于液柱下游界面的凹面形式, 膨胀波在向上游传播的过程中会聚焦于一点. 第一反射波在液柱内部聚焦后继续传播, 在迎风面发生反射, 形成第二反射波. 在3个液柱共存的情况下, 激波在通道内部发生多次反射. 由于反射激波与液柱接触点、接触时刻以及反射激波强度的不同, 液柱内部产生的多重透射激波必然会与生成的多道第一反射波和第二反射波相互作用, 使液柱内部的压力波动更为剧烈. 液柱内部的高压区域还会影响其内部流动, 进而影响到液柱外部的凸起产生[26]. 随着初始通道宽度的减小, 液柱在外部持续高压与内部剧烈压力变化的共同作用下, 在波系演化阶段就已经发生了较大的形变. 如图8(b)所示, 在黑色框线所示区域中, 液柱A和B的通道侧已经出现了明显的凹坑, 界面在液柱内外的压力差作用下变平.
![]()
图 8 算例A2 (左)和A3 (右)的纹影图(上)和压力云图(下)Figure 8. Schlieren (top) and pressure (bottom) contours of case A2 (left) and case A3 (right)液柱的位置分布对于其变形破碎过程有着显著影响. 尽管在液柱变形的波系演化阶段, 3个算例中都出现了相似的液滴形态特征, 但在整体上的变形却存在较大差异. 随着气流的持续作用, 这些差异逐渐累积和放大, 导致气液界面发生较大的形变. 而这种界面形变又反过来影响气流, 在气液两相的相互作用下, 3个液柱的变形趋势更加复杂.
气流经第一通道向下游传播时,由于第二通道宽度的减小, 气流对该通道两侧界面的作用显著增强. 如图8(c)所示, 在A2算例中, 由于气流从第二通道流出时受到的阻碍增强, 气流对第二通道入口处界面的冲击也相应增强, 这导致液柱A和B背风面靠近第二通道入口的部分由平面变为一个内部较为光滑的凹陷. 液柱C极点附近的迎风面在更强的气流冲击下变为平面, 第二通道出口处两侧均拉伸出更长的液膜; 如图8(g)所示, 在算例A3中, 气流受到的阻碍作用较算例A2更为显著, 因此, 液柱A和B背风面的凹陷进一步扩大, 且凹陷内部被卷起形成褶皱. 液柱C极点附近的迎风平面演化更为迅速, 平面中部已经形成微小凹陷. 靠近第二通道出口处的迎风面上也产生了较大的凹陷, 两侧界面拉伸出的液膜相互靠近, 导致该通道出口被封闭. 在算例A2和A3中, 液柱C后部的凹陷均进一步加深, 导致第二通道的形态演化趋势转变为向下游延伸.
为了研究通道对3液柱整体变形形态的影响, 通过测量3液柱整体的无量纲宽度W/W0、无量纲高度H/H0和无量纲质心位移(S − S0)/R对3液柱整体的变形过程进行定量研究, 如图9所示. 该整体的质心位置定义为
![]()
图 9 3液柱整体无量纲宽度、无量纲高度和无量纲质心位移随时间变化Figure 9. Time variation of dimensionless width and dimensionless height and dimensionless centroid displacement of the three liquid columns as a whole$$ \bar x = \frac{{\displaystyle\int {{\alpha _l}} {\rho _l}x{\text{ }}{\mathrm{d}}x{\mathrm{d}}y}}{{\displaystyle\int {{\alpha _l}} {\rho _l}{\text{ }}{\mathrm{d}}x{\mathrm{d}}y}}{\text{ }}\quad ({\alpha _l} \geqslant 0.5) $$ (10) 其中, x和y分别为计算域中参与计算的网格横向坐标和纵向坐标.
随着通道宽度减小, 气流流出通道所受到的阻碍增强, 气流在通道内部被压缩, 使得通道内部形成的高压区域作用时间更长, 作用于液柱的表面积更广. 从图9可以观察到, 对于A1, A2和A3算例, 其3液柱整体无量纲宽度变化趋势均为先减小后增大. 由于液柱变形破碎过程中, 液柱之间的通道被显著压缩, 液柱也被压扁, 其整体的无量纲宽度减小. 但是3个算例在发展后期, 其整体宽度均有显著增长. 在算例A1中, 整体无量纲宽度增长的原因是液柱C在气流作用下前部凹坑P3持续扩大, 背风面由向内凹陷演化为向外突出的弧面, 如图9中液滴轮廓所示. 但外凸弧面形成较晚, 算例A1中整体无量纲宽度由减小到增大的转折点则相应出现较晚.
而在算例A2和A3中, 其整体宽度的增长趋势是因为第二通道出口附近的下游液柱C的赤道侧在气流作用下形成尖端继而拉伸出较长的液膜, 如图8(g)和图8(h)中红色框线区域所示. 第二通道宽度的减小增强了气流与液柱界面的相互作用, 因此在第二通道宽度最小的算例A3中, 液柱C赤道侧形成尖端向下游移动的速度最快, 整体无量纲宽度变化趋势最早出现增长的转折.
3个算例中的整体无量纲高度随时间变化的趋势差异较小, 这是因为气流在通道内部与界面的作用对液柱A和B后部以及液柱C的形态演化影响较为显著, 而对液柱A和B赤道侧尖端的影响较小. 第二通道宽度的减小, 使得液柱A和B后部在气流冲刷和通道内部高压的共同作用下形成凹陷, 液柱A和B主体部分逐渐向各自非通道侧移动, 如图8(g)和图8(h)所示. 3液柱整体宽度和高度在增长过程中的骤然下降则是由于液滴表面“脱涡”现象的出现所致.
初期阶段各算例的整体无量纲质心位移差异不大. 随着时间推进, 算例A1的整体无量纲质心位移逐渐超越算例A2和A3. 在算例A1和A2中, 液柱A和B的第一通道侧均形成了向下游移动的尖端(见图4(d)和图8(g)的黑色框线所示区域), A2和A3算例的液柱A和B后部在气流作用下形成较大尺度的凹陷(见图8(g)和图8(h)红色弧线所示区域), 该凹陷尺寸逐渐扩大, 凹陷附近的回流现象导致液柱A和B后部通道侧的部分向下游移动的速率减小, 而非通道侧的部分则呈现出在y轴方向移动的趋势. 因此, 在液柱演化后期, 算例A3的整体无量纲质心位移相较于A2要小. 相比之下, 算例A1中液柱后部并未形成较大尺度凹陷, 液柱A和B后部y轴方向的移动趋势不明显, 这使得算例A1的无量纲质心位移最大. 从能量传递的角度来看, 随着第二通道宽度的减小, 气流与液柱界面之间的相互作用增强, 从气相传递到液相的能量显著提升, 液柱之间形成的通道类似于3液柱作为一个整体的内部空腔结构, 在这种结构中, 气体在空腔内被进一步压缩, 增加了其对液柱界面的作用能力, 因而界面变形速率随之加快.
4. 结 论
本文采用时空守恒元与解元方法, 对平面激波冲击交错排列3液柱的物理过程进行研究. 通过分析波系结构的演化过程, 发现液柱间的通道显著增强了激波间的相互作用, 展示了平面激波冲击3液柱后, 液柱的变形过程, 探究了液柱位置对液柱变形过程的影响. 研究揭示了下游液柱(液柱C)在激波冲击下的动态行为. 数值结果表明, 液柱C迎风面上压力分布不均匀导致其迎风面从弧面演化为平坦表面, 迎风面向来流方向倾斜, 极点附近成为凹坑. 液柱C宽度减小速率始终快于单液柱, 但τ > 1.009阶段, 第二通道出口处的气流与非通道侧来流相互作用, 减少了液柱C背风区的涡量, 减弱了对背风面的冲刷作用, 从而使液柱宽度w的变化速率减缓. 研究发现第二通道宽度在气流与界面的相互作用中发挥了关键作用. 随着第二通道宽度减小, 气流对界面的作用能力增强, 液柱A和B后部凹陷的演化速率加快, 3液柱整体的无量纲宽度由减小到增大的转折点提前, 后期的无量纲质心位移速度减小. 在第二通道宽度较大的算例中, 尽管液柱A和B在其第一通道侧与非通道侧存在气流流动的不对称性, 但两液柱仍保持了较为对称的结构特征. 而在第二通道宽度较小的算例中, 由于第二通道变窄导致气流的冲刷作用增强, 液柱A和B在波系演化阶段显现出显著的非对称性. 本文对工程应用中出现的多液柱共存情况进行了初步探索, 研究对象为规则的交错排列3液柱, 在丰富性和完善性方面尚显不足, 后续将对激波冲击非规则的交错排列多液柱过程开展研究.
-
![]()
![]()
图 3 初期波系结构演化纹影图(上)和压力云图(下)
Figure 3. Schlieren (top) and pressure (bottom) contours of the wave system evolution at the early stage
![]()
图 4 液柱破碎中期压力图(上)和涡量云图(下)
Figure 4. Pressure (top) and vorticity (bottom) contours at the middle stage of liquid column breakup
![]()
图 5 液柱迎风面极点位移与宽度示意图
Figure 5. Definition of windward pole point drift and width of the liquid column
![]()
图 6 液柱迎风面极点位移和宽度随时间变化
Figure 6. Time variation of windward pole point drift and width of the liquid column
![]()
图 8 算例A2 (左)和A3 (右)的纹影图(上)和压力云图(下)
Figure 8. Schlieren (top) and pressure (bottom) contours of case A2 (left) and case A3 (right)
![]()
图 9 3液柱整体无量纲宽度、无量纲高度和无量纲质心位移随时间变化
Figure 9. Time variation of dimensionless width and dimensionless height and dimensionless centroid displacement of the three liquid columns as a whole
表 1 数值模拟算例的参数设置
Table 1 Parameters setup for the numerical simulation cases
Case Ms Xc/m A1 1.915 0.004546 A2 1.915 0.003500 A3 1.915 0.002778 
下载: 导出CSV
-
[1] Reitz RD, Bracco FV. Mechanism of atomization of a liquid jet. Physics of Fluids, 1982, 25(10): 1730-1742 doi: 10.1063/1.863650
[2] 刘贵兵, 侯海量, 朱锡等. 液滴对爆炸冲击波的衰减作用. 爆炸与冲击, 2017, 37(5): 844-852 (Liu Guibing, Hou Hailiang, Zhu Xi, et al. Attenuation of shock wave passing through liquid droplets. Explosion and Shock Waves, 2017, 37(5): 844-852 (in Chinese) doi: 10.11883/1001-1455(2017)05-0844-09 Liu Guibing, Hou Hailiang, Zhu Xi, et al. Attenuation of shock wave passing through liquid droplets. Explosion and Shock Waves, 2017, 37(5): 844-852 (in Chinese) doi: 10.11883/1001-1455(2017)05-0844-09
[3] Marmottant P, Villermaux E. On spray formation. Journal of Fluid Mechanics, 2004, 498: 73-111 doi: 10.1017/S0022112003006529
[4] Joseph DD, Beavers GS, Funada T. Rayleigh-Taylor instability of viscoelastic drops at high Webber numbers. Journal of Fluid Mechanics, 2002, 453: 109-132 doi: 10.1017/S0022112001006802
[5] Hinze JO. Critical speeds and sizes of liquid globules. Flow, Turbulence and Combustion, 1949, 1(1): 273-288 doi: 10.1007/BF02120335
[6] Pilch M, Erdman CA. Use of breakup time data and velocity history data to predict the maximum size of stable fragments for acceleration-induced breakup of a liquid drop. International Journal of Multiphase Flow, 1987, 13(6): 741-757 doi: 10.1016/0301-9322(87)90063-2
[7] Wierzba A. Deformation and breakup of liquid drops in a gas stream at nearly critical Weber numbers. Experiments in Fluids, 1990, 9(1-2): 59-64 doi: 10.1007/BF00575336
[8] Lane WR. Shatter of drops in streams of air. Industrial and Engineering Chemistry, 1951, 43(6): 1312-1317 doi: 10.1021/ie50498a022
[9] Hanson AR, Domich EG, Adams HS. Shock tube investigation of the breakup of drops by air blasts. Physics of Fluids, 1963, 6(8): 1070-1080 doi: 10.1063/1.1706864
[10] Chou WH, Hsiang LP, Faeth GM. Temporal properties of drop breakup in the shear breakup regime. International Journal of Multiphase Flow, 1997, 23(4): 651-669 doi: 10.1016/S0301-9322(97)00006-2
[11] Theofanous TG, Li GJ, Dinh TN, et al. Aerobreakup in disturbed subsonic and supersonic flow fields. Journal of Fluid Mechanics, 2007, 593: 131-170 doi: 10.1017/S0022112007008853
[12] Joseph DD, Belanger J, Beavers GS. Breakup of a liquid drop suddenly exposed to a high-speed airstream. International Journal of Multiphase Flow, 1999, 25(6-7): 1263-1303 doi: 10.1016/S0301-9322(99)00043-9
[13] Theofanous TG, Li GJ. On the physics of aerobreakup. Physics of Fluids, 2008, 20(5): 052103
[14] Nicholls JA, Ranger AA. Aerodynamic shattering of liquid drops. AIAA Journal, 1969, 7(2): 285-290 doi: 10.2514/3.5087
[15] Liu Z, Reitz RD. An analysis of the distortion and breakup mechanisms of high speed liquid drops. International Journal of Multiphase Flow, 1997, 23(4): 631-650 doi: 10.1016/S0301-9322(96)00086-9
[16] Guildenbecher DR, Lopez-Rivera C, Sojka PE. Secondary atomization. Experiments in Fluids, 2009, 46(3): 371 doi: 10.1007/s00348-008-0593-2
[17] Meng JC, Colonius T. Numerical simulations of the early stages of high-speed droplet breakup. Shock waves, 2015, 25: 399-414 doi: 10.1007/s00193-014-0546-z
[18] 宋家喜, 潘书诚. 高马赫数下激波液滴相互作用的数值模拟研究. 力学学报, 2022, 54(9): 2419-2434 (Song Jiaxi, Pan Shucheng. Numerical investigation of shock-droplet interaction with high-Mach numbers. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(9): 2419-2434 (in Chinese) Song Jiaxi, Pan Shucheng. Numerical investigation of shock-droplet interaction with high-Mach numbers. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(9): 2419-2434 (in Chinese)
[19] Xu S, Fan WQ, Wu WX, et al. Analysis of wave converging phenomena inside the shocked two-dimensional cylindrical water column. Journal of Fluid Mechanics, 2023, 964: A12 doi: 10.1017/jfm.2023.239
[20] Ling Y, Mahmood T. Detailed numerical investigation of the drop aerobreakup in the bag breakup regime. Journal of Fluid Mechanics, 2023, 972: A28 doi: 10.1017/jfm.2023.708
[21] Xu ZK, Wang TY, Che ZZ. Breakup of particle-laden droplets in airflow. Journal of Fluid Mechanics, 2023, 974: A42 doi: 10.1017/jfm.2023.847
[22] Igra D, Takayama K. Numerical simulation of shock wave interaction with a water column. Shock Waves, 2001, 11(3): 219-228 doi: 10.1007/PL00004077
[23] Igra D, Takayama K. Experimental investigation of two cylindrical water columns subjected to planar shock wave loading. Journal of Fluids Engineering, 2003, 125(2): 325-331 doi: 10.1115/1.1538628
[24] Chen H, Liang SM. Flow visualization of shock/water column interactions. Shock Waves, 2008, 17(5): 309-321 doi: 10.1007/s00193-007-0115-9
[25] 吴润龙, 李祝军, 丁航. 平面激波冲击并排液滴的三维数值模拟研究. 力学学报, 2022, 54(11): 2958-2969 (Wu Runlong, Li Zhujun, Ding Hang. Impact of a planar shock onto side-by-side droplets: A 3D numerical study. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(11): 2958-2969 (in Chinese) Wu Runlong, Li Zhujun, Ding Hang. Impact of a planar shock onto side-by-side droplets: A 3D numerical study. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(11): 2958-2969 (in Chinese)
[26] Guan B, Liu Y, Wen CY, et al. Numerical study on liquid droplet internal flow under shock impact. AIAA Journal, 2018, 56(9): 3382-3387 doi: 10.2514/1.J057134
[27] Shen H, Wen CY, Parsani M, et al. Maximum-principle-satisfying space-time conservation element and solution element scheme applied to compressible multifluids. Journal of Computational Physics, 2017, 330: 668-692 doi: 10.1016/j.jcp.2016.10.036
[28] Shen H, Wen CY, Liu K, et al. Robust high-order space-time conservative schemes for solving conservation laws on hybrid meshes. Journal of Computational Physics, 2015, 281: 375-402 doi: 10.1016/j.jcp.2014.10.023
[29] Shen H, Wen CY, Zhang DL. A characteristic space-time conservation element and solution element method for conservation laws. Journal of Computational Physics, 2015, 288: 101-118 doi: 10.1016/j.jcp.2015.02.018
[30] Guan B, Yang HS, Yang HW, et al. On the irregular jet formation of shock-accelerated spherical heavy gas bubbles. Physics of Fluids, 2022, 34(12): 126111
[31] Yi XY, Zhu YJ, Yang JM. On the early-stage deformation of liquid drop in shock-induced flow//30th International Symposium on Shock Waves, Israel, 2015-7-19-24. Switzerland: Springer, 2017: 1269-1273
[32] Sembian S, Liverts M, Tillmark N, et al. Plane shock wave interaction with a cylindrical water column. Physics of Fluids, 2016, 28(5): 056102
[33] Liu N, Wang ZG, Sun MB, et al. Numerical simulation of liquid droplet breakup in supersonic flows. Acta Astronautica, 2018, 145: 116-130 doi: 10.1016/j.actaastro.2018.01.010
计量
- 文章访问数: 128
- HTML全文浏览量: 23
- PDF下载量: 30
