Richtmyer-Meshkov instability in elastic-plastic media

中国科学技术大学硕士学位论文一般状态方程多流体界面数值方法研究姓名：郑建国申请学位级别：硕士专业：流体力学指导教师：孙德军;尹协远20050501中文摘要摘要①本文发展了一类一般状态方程可压缩多流体界面的数值模拟方法，并具体应用到三种不同的非理想气体状态方程，包括ｓｔｉ＆ｎｆ刚性）气体状态方程，ｖａｒｌｄｅｒＷａａｌｓ状态方程以及工程上广泛适用的更一般的Ｍｉｅ—Ｇｒｉｉｎｅｉｓｅｎ状态方程。

此方法主要的特点是：ｆ１）．采用体积分数多流体数学模型，这是在假设多流体交界面两侧压力和速度平衡的基础上根据二相流理论建立的，并引入计算混合流体压力的“状态方程”使系统封闭。

（２）．将高精度、高分辨率的ＰｉｅｃｅｗｉｓｅＰａｒａｂｏｌｉｃＭｅｔｈｏｄ（ＰＰＭ）数值方法推广到多流体问题中，用膨胀激波代替稀琉波，采用双波近似的方法求解多流体Ｒｉｅｍａｎｎ问题。

（３）．使用Ｌａｇｒａｎｇｉａｎ－Ｒｅｍａｐｐｉｎｇ两步法求解模型方程组。

与以往的多流体方法相比，本文的方法具有一些优点。

首先，体积分数多流体数学模型所采用的交界面两侧压力和速度平衡的假设与真实的物理情况比较接近，它消除了交界面上压力的振荡；特别是其模型简单，并且不因为具体的状态方程而改变，便于应用到复杂状态方程的多流体流动问题。

其次，文中推广的多流体ＰＰＭ方法处理交界面问题的效果非常好，它继承了原始ＰＰＭ的高分辨率和能有效抑制间断上压力振荡的优点。

最后，Ｌａｇｒａｎｇｉａｎ—Ｒｅｍａｐｐｉｎｇ形式的ＰＰＭ方法具有Ｌａｇｒａｎｇｅ类方法的特点，它可以有效地处理多流体界面，为了验证方法是否合理有效，进行了大量的数值实验。

一维和二维算例表明本文的方法可以有效地处理一般状态方程的接触间断、激波、激波和接触间断的相互作用以及多维滑移线等物理问题。

从数值结果中可以很明显地看出交界面附近压力无振荡，并能够比其它一般多流体数值方法更糟细地模拟多流体交界面。

本文还研究了柱坐标下内聚激波诱导的Ｐｄｃｈｔｍｙｅｒ—ＭｅｓｈｋｏｖＩｎｓｔａｂｉｌｉｔｙ（ＲＭｔ）ｅ从模拟的结果来看，演化过程中出现的钉状（ｓｐｉｋｅ）和泡状（ｂｕｂｂｌｅ）结构以及后期的蘑菇状交界面都很清晰。

With y y x x R J R J 1,1== and ()y x y x x yk k k J k J Q −+=22, the coefficients i a in the equation (9’) are:()()()()y x y x x y y x y x y x y x y x k k k J k J k k J J J J J J Q k k g a +−−+−++−=2222209991418()()()()(()()()()()())y y x x y x y x x y y x y x y x y x x y y x y x y x y x x y y x y x x y y x y x x y y x y y x x y x k J k J k k k J k J k k J J k k k J k J J J k k J J k J k J J J k J k J J J k J k J J J k J k J J J Q g a ++++++−+−+−+−−+++++=593576286332323223322222/31()()()(()()()())y x y x x y y x y x x y y x y x y x y x y x x y y x y x y x k k k k k J k J k J k J k k J J J J k k k J k J J J k k J J Q g a +++++−−−+−+++=4422222222222330141813()()()()()x y y x y x y x y x x y y x k J k J J J k k k J k J J J Q g a 332/13233+++−+=14=a1. Lord Rayleigh, “Investigations of the character of the equilibrium of an incompressible heavy fluid of variable density” Proc.Lond.Math.Soc. 14, 170 (1883)2. R.M.Davies, G.I.Taylor, “The mechanics of large bubbles rising through extended liquids and through liquids in tubes” Proc.Roy.Soc. A200, 375 (1950)3. D.H.Sharp, “An overview of Rayleigh-Taylor instability” Physica D 12, 3 19844. S.Chandrasekhar, Hydrodynamic and hydromagnetic stability, 3d ed. (Dover Publications Inc., New York, 1981), pp.428-477; K.O.Mikaelian, “Effect of viscosity on Rayleigh-Taylor and Richtmeyer-Meshkov instabilities” Phys.Rev.E 47, 375 (1993)5. K.I.Read, “Experimental investigation of turbulent mixing by Rayleigh-Taylor instability” Physica D 12, 45 (1984); Yu.A.Kucherenko at al, “Experimental study into the asymptotic stage of separation of the turbulized mixtures in gravitationally stable mode” Laser Part. Beams 15, 17 (1997); Izv. Akad. Nauk SSSR, Mekh. Zhidk. Gaza 6, 157 (1978)6. G.Dimonte M.Schneider, “Density ratio dependence of Rayleigh-Taylor mixing for sustained and impulsive acceleration histories” Phys.Fluids 12, 304 (2000);M.Schneider, G.Dimonte, B.Remington, “Large and small scale structures in Rayleigh-Taylor mixing” Phys.Rev.Lett. 80, 3507 (1998)7. F.H.Harlow, J.E.Welch, “Numerical study of large-amplitude free-surface motion” Phys.Fluids 9,842 (1966); D.L.Youngs, “Numerical simulations of turbulent mixing by Rayleigh-Taylor instability” Physica D 12, 32, (1984); R.Menikoff, C.Zemach,“Rayleigh-Taylor instability and the use of conformal maps for ideal fluid flow”p.Phys. 51, 28 (1983); G.R.Baker, D.I.Meiron, S.A.Orszag, “Vortexsimulations of the Rayleigh-Taylor instability” Phys.Fluids 23, 1485 (1980); G.Gardner, J.Glimm, O.McBryan, R.Menikoff, D.H.Sharp, Q.Zhang, “The dynamics of bubble growth for Rayleigh-Taylor instability” Phys.Fluids 31, 447 (1988)8. X.L.Li, “Study of three-dimensional Rayleigh-Taylor instability” Phys.FluidsA5, 1904 (1993); X.L.Li, “A numerical study of three-dimensional bubble merger in the Rayleigh-Taylor instability” Phys.Fluids 8, 336 (1996)9. D.L.Youngs, “Three-dimensional numerical simulations of turbulent mixing by Rayleigh-Taylor instability” Phys.Fluids A 3, 1312 (1991)10. U.Alon, J.Hecht, D.Offer, D.Sharts, “Power-laws and similarity of Rayleigh-Taylor and Richtmyer-Meshkov mixing fronts at all density ratios” Phys.Rev.Lett.74, 534 (1995)11. S.I.Abarzhi, “Stable steady flows in Rayleigh-Taylor instability”Phys.Rev.Lett.89, 1332 (1998); “Non-linear three-dimensional Rayleigh-Taylor instability” Phys.Rev.E, 59, 1729 (1999)12. A.V.Shubnikov and V.A.Koptsik, Symmetry in science and art, Plenum New York, 197413. yzer, “On the instability of superposed fluids in a gravitational field”Astrophys.Jour. 122, 1 (1955); P.R.Garabedian, “On steady-state bubbles generated by Taylor instability” Proc.Roy.Soc. A241, 423 (1957)14. S.I.Abarzhi, “The stationary spatially periodic flows in Rayleigh-Taylor instability: solutions multitude and its dimension” Physica Scripta T66, 238 (1996);“Stationary solutions of the Rayleigh-Taylor instability for spatially periodic flows:questions of uniqueness, dimensionality and universality” Zh.Eksp.Teor.Fiz. 110,1841 (1996) (Sov.Phys.JETP 83, 1012)15. S.I.Abarzhi and N.A.Inogamov, “Stationary solutions in the Rayleigh-Taylor instability for spatially periodic flow” Zh.Eksp.Teor.Fiz. 107,245 (1995) (Sov.Phys.JETP 80,132 (1995)); N.A.Inogamov, “Stationary solutions in the theory of the hydrodynamic Rayleigh-Taylor instability” Pis’ma Zh.Eksp.Teor.Fiz. 55, 505 (1992) (Sov.Phys.JETP Lett. 55, 521 (1992))16. A.Oparin, “Numerical simulations of Rayleigh-Taylor instability using quasi-monotonic grid-characteristic approach” private communication, 1998; M.F.Ivanov,A.M.Oparin, V.G.Sultanov, V.E.Fortov, “Certain features of development of the Rayleigh-Taylor instability in 3D geometry” Doklady Physics 44, 491 (1999) (Doklady Akademii Nauk 367, 464 (1999))17. R.J.Taylor, J.P.Dahlburg, A.Iwase, J.J.Gardner, D.E.Fyfe O.Willi,“Measurement and Simulation of Laser Imprinting and Consequent Rayleigh-Taylor growth” Phys.Rev.Lett.76, 1643 (1996); M.M.Marinak, S.G.Gledinning, R.J.Wallace, B.A.Remington, et al, “Non-linear evolution of a three-dimensional multimode perturbation” Phy.Rev.Lett. 80, 4426 (1998);18. Dalziel SB, Linden PF, Youngs DL, “Self-similarity and internal structure of turbulence induced by Rayleigh-Taylor instability” Jour.Fluid Mechanics, 399, 1 (1999)19. V.E.Zakharov, “Stability of periodic waves of finite amplitude on deep fluid surface” Prikl.Mat.Teor.Fiz 2, 86 (1968) (Sov.Math. Appl. Math. and Theor. Phys. 2, 86 (1968)); E.Fermi and J.von Newman, “Taylor instability of an incompressible liquid”(1951) in book: E.Fermi, Collected papers, The University of Chicago Press, Chicago 1962, v.2, 81620. J.W.Jacobs and I.Catton, “Three-dimensional Rayleigh-Taylor instability. Part I. Weakly nonlinear theory” Jour. Fluid Mech. 187, 329 (1988)21. X.Y.He, R.Y.Zhang, S.Y.Chen, G.D.Doolen, “On three-dimensional Rayleigh-Taylor instability” Phys.Fluids 11, 1143 (1999)22. J.Hecht, U.Alon, D.Shvarts, "Potential flow models of Rayleigh-Taylor and Richtmeyer-Meshkov bubble fronts" Phys.Fluids 6, 4019 (1994)23. V.S.Arpachi, A.Asmaeeli, “On interface dynamics” Phys.Fluids 12, 1244 (2000)Figure captionsFIG.1. Bubbles generated by the Rayleigh-Taylor instability. Flow symmetry group is p2mm , λx,y are spatial periods in the plane (x, y ), rotation axis 2 and mirror planes m are perpendicular to the plane of translations. Black circle marks the positions of bubble top.FIG.2.(a). Plane of parameters R x,y of steady solutions family, finite values of aspect ratio x y k k . Physical region is bounded by the flattened bubbles (“solitary spikes”) ∞→y x R ,and by the edge solutions cr y x R ,. Dashing line relates to bubbles with circular contour ()()2x y y x k k R R=. Bubbles are y -elongated at ()()2x y y x k R R < and x -elongated at ()()2x y y x k k R R >.(b). Low-symmetric steady solutions under the dimensional crossover 0→x y k k : the solutions with ()2x y y x k k R R ≥ are restricted to region of flattened bubbles (“solitary spikes”), while solutions with ()2x y y x k k R R < become solutions of the 2D flow.FIG.3. The real parts of the Lyapunov exponents for low-symmetric steady solutions.(a). Aspect ratio 2/1=x y k k and the radius of curvature ∞≡y R .(b). Limiting case of “square flow” p4mm , 1=x y k k and x y R R =.FIG.4. Stability region for the low-symmetric bubbles in first approximation (dashing area). Finite values of aspect ratio x y k k , dashing line relates to the bubbles with circular contour, ()()2x y y x k k R R =.FIG.5. Stability of steady solutions under dimensional crossover 0→x y k k .Dashing lines are real parts of the Lyapunov exponents of the 2D flow, unstable mode D D w 23− is plotted at 1.0=x y k k and 001.0=δ. cr c R k R ≈=2. Black points are the Hopf bifurcations, which bound the region of stability of the 2D flow at N =2.FIG.6. Secondary instabilities of the low-symmetric bubbles under spatial modulations: merging and splitting.。

封闭空间中不同压力下火焰过孔板加速机理研究赵健福;周磊;钟力嘉;卫海桥【摘要】爆震或超级爆震发生时总会伴随着湍流火焰-冲击波相互作用,对其开展研究是揭示爆震或超级爆震机理的关键,研究火焰加速产生压力波的过程是火焰-压力波相互作用研究的基础性前提.基于自主设计的定容燃烧弹和Converge三维数值模拟方法,对封闭空间中火焰过孔板加速机理及影响因素开展了研究,讨论了初始压力对火焰过孔板加速的影响.依据火焰传播形态与速度,将火焰过孔板加速过程分为3个阶段:层流火焰阶段、射流火焰阶段和湍流火焰阶段.通过分析火焰过孔板过程中的流场情况,发现在火焰未到达孔板前,孔板附近存在强射流,火焰受强射流的驱动而急剧加速;但当火焰穿过孔板之后,火焰锋面前的流场速度沿着远离火焰的方向而逐渐下降,说明开始由火焰驱动未燃气体运动.比较不同压力下的火焰过孔板过程,发现湍流火焰传播速度和缸压振荡均随着初始压力的提高而升高.【期刊名称】《实验流体力学》【年(卷),期】2019(033)004【总页数】10页(P11-20)【关键词】孔板;火焰加速;射流;火焰传播;压力振荡【作者】赵健福;周磊;钟力嘉;卫海桥【作者单位】天津大学内燃机燃烧学国家重点实验室,天津300072;天津大学内燃机燃烧学国家重点实验室,天津300072;天津大学内燃机燃烧学国家重点实验室,天津300072;天津大学内燃机燃烧学国家重点实验室,天津300072【正文语种】中文【中图分类】TK4170 引言发动机小型强化技术被认为是最有前途的点燃式(Spark Ignition，SI)发动机节能减排技术措施之一[1]。

然而，小型强化SI发动机热负荷的增加，导致在其燃烧过程中更容易发生爆震(Knock)[2]、超级爆震(Super- Knock)[3- 4]等不正常燃烧现象，限制了小型强化SI发动机热效率进一步提升。

目前爆震产生的机理尚不明确，末端气体自燃理论由于很好地解释了光学实验的结果而得到广泛的认可。

雷诺数：对于不同的流场，雷诺数可以有很多表达方式。

这些表达方式一般都包括流体性质（密度、黏度）再加上流体速度和一个特征长度或者特征尺寸。

这个尺寸一般是根据习惯定义的。

比如说半径和直径对于球型和圆形并没有本质不同，但是习惯上只用其中一个。

对于管内流动和在流场中的球体，通常使用直径作为特征尺寸。

对于表面流动，通常使用长度。

管内流场对于在管内的流动,雷诺数定义为:式中:•是平均流速(国际单位: m/s)•管直径(一般为特征长度) (m)•流体动力黏度 (Pa·s或N·s/m²)•运动黏度 (ρ) (m²/s)•流体密度(kg/m³)•体积流量 (m³/s)•横截面积(m²)假如雷诺数的体积流率固定，则雷诺数与密度（ρ）、速度的开方（）成正比；与管径（Ｄ）和黏度（ｕ）成反比假如雷诺数的质量流率（即是可以稳定流动）固定，则雷诺数与管径（Ｄ）、黏度（ｕ）成反比；与√速度（）成正比；与密度（ρ）无关平板流对于在两个宽板(板宽远大于两板之间距离)之间的流动,特征长度为两倍的两板之间距离。

流体中的物体对于流体中的物体的雷诺数,经常用Re p表示。

用雷诺数可以研究物体周围的流动情况，是否有漩涡分离，还可以研究沉降速度。

流体中的球对于在流体中的球，特征长度就是这个球的直径，特征速度是这个球相对于远处流体的速度，密度和黏度都是流体的性质。

在这种情况下，层流只存在于Re=0.1或者以下。

在小雷诺数情况下，力和运动速度的关系遵从斯托克斯定律。

搅拌槽对于一个圆柱形的搅拌槽，中间有一个旋转的桨或者涡轮，特征长度是这个旋转物体的直径。

速度是ND,N是转速(周/秒)。

雷诺数表达为:当Re>10,000时，这个系统为完全湍流状态。

[1]过渡流雷诺数对于流过平板的边界层，实验可以确认，当流过一定长度后,层流变得不稳定形成湍流。

对于不同的尺度和不同的流体，这种不稳定性都会发生。

Richtmyer-Meshkov不稳定性强化混合参变机理王兵;卢梦【摘要】在不同参数条件下，计算分析了H2O和N2等混合物界面上激波诱导Richtmyer-Meshkov（R—M）不稳定性过程．采用有限差分方法数值求解了二维可压缩Navier-Stokes方程，对流项以5阶特征紧致．WENO混合格式离散，输运项以6阶对称紧致格式离散，时间方向以3阶显式Runge—Kutta方法推进．研究表明，界面振幅和激波强度增大，均可增强界面附近涡量场，强化混合．【期刊名称】《气体物理》【年(卷),期】2016(001)006【总页数】17页(P5-21)【关键词】激波加速界面;紧致格式;激波强度;涡合;拟涡能【作者】王兵;卢梦【作者单位】清华大学航天航空学院,北京100084【正文语种】中文【中图分类】V211.3存在扰动的两种密度不同的流体交界面,在受到外部作用时,扰动将发展,界面将失稳,导致两种流体发生混合,即界面不稳定性.激波与不同流体界面的相互作用广泛存在于自然界和工程应用领域[1],导致的流动失稳被称为Richtmyer-Meshkov(R-M)不稳定性.它对流体的混合具有促进作用,从而影响流动发展进程.对于激波加速界面导致的流体不稳定性的研究可追溯到20世纪50年代[2].1960年,Richtmyer对激波加速下的界面流动不稳定性给出了严格的理论推导和分析预测; 1969年,Meshkov通过激波管实验验证了Richtmyer的分析预测[3].球形界面和平面激波是研究激波与界面作用的良好模型.Haas等[4]通过在水平激波管中观测弱激波对He和R22的球形界面的作用,研究了激波反射折射及绕射等复杂现象,详细分析了作用过程中的激波速度和界面速度变化规律.球形界面实验研究的不足是需要支杆支撑,支杆会对其附近的界面发展产生影响.柱形界面也是研究R-M不稳性的基本模型之一.Tomkins等[5]运用实验测量了激波与无膜技术形成的气柱作用过程中气体混合的速率,揭示了气体混合的机制,并指出了二次不稳定性在界面发展后期的重要作用.单模态界面在R-M不稳定性研究中也被大量使用[6].例如,Brouillette等[7]采用厚度为0.5 μm的硝化纤维膜进行了相关研究,Sadot等[8]采用厚度为0.1 μm的薄膜形成的单模界面开展了与平面激波作用的实验.通常,分开两种流体的薄膜破裂的碎片会影响不稳定性的发展,且这种影响很难预测.扰动振幅增长率的激波管实验结果与理论预测之间的偏差大小取决于其强度.为了避免薄膜对实验结果的影响,可采用薄金属板分开流体,然后在启动激波管前通过机械装置快速抽出隔板,将抽出隔板时产生的扰动作为界面上初始的扰动.因为分子扩散,这样得到的是一个密度梯度为有限值的连续界面,与密度间断界面相比不稳定性的发展会减慢.然而,采用这种方法形成的初始扰动不具有可重复性.一种改进的生成无膜界面的方法是,让两种介质分别从竖直激波管的两端进入观察段中形成界面,并在实验段位置开设狭缝排出多余的气体,通过振荡实验仪器来得到正弦的界面[9-10]. 虽然实验研究比较直观,但难度大,重复性差,且大多数结果都是定性的流动显示,定量的实验数据少.随着计算机和数值计算技术的发展,越来越多的学者运用数值模拟的方法对R-M不稳定性开展研究工作.Holmes等[11]用前缘追踪数值方法模拟了R-M不稳定性,准确捕捉了物质界面和激波.由于要求追踪所有的波系和界面的发展,该算法的计算量很大.Anuchina等[12]用MAH-3程序得到了R-M不稳定性线性和非线性发展阶段的结果,对比了二维和三维的界面增长率.Latini等[13]用WENO重构方法对R-M不稳定性进行数值模拟,发现低阶方法在不稳定性发展的后期保留了大尺度结构和对称性,而高阶方法复现了小尺度结构,导致对称结构被破坏,从而增强了流体间的混合.Hu等[14]用6阶中心-迎风WENO格式模拟了R-M不稳定性,并与相关实验结果做了对比.近年,Olson等[15]用大涡模拟方法模拟了R-M不稳定性,讨论了数值方法网格精度Reynolds数的影响.关于早期的R-M不稳定性研究,可参考Holmes等[16]的综述,文献系统而详尽地比较了实验数值模拟线性理论非线性理论,以及脉冲理论的结果.在不稳定性发展早期,界面以大尺度结构为主,线性理论通常能较好地与实验和数值模拟结果吻合.在不稳定性发展后期,界面两侧的流体出现速度差,产生二次不稳定性,出现一些小尺度结构,精度越高的数值算法,越能识别出小尺度结构,而实验方法则难以测量.本文针对激波与H2O和N2等混合物界面相互作用而发生的R-M不稳定性流动现象数值模拟,主要分析界面振幅和激波强度的影响特征及规律,为工程应用提供理论指导.因研究的流体由多种组分构成,故基本控制方程由二维可压缩Navier-Stokes方程和组分质量守恒方程组成.不考虑外部体积力做功和辐射换热,其守恒型为第k种组分的黏性系数采用基于Chapman-Enskog理论的公式,第k种组分的导热系数采用改进的Eucken公式,两种组分间的相互扩散系数计算采用Chapman 和Enskog基于Boltzmann方程提出的公式,分别为混合流体的黏性系数计算公式如下:解出守恒变量后,依据守恒变量可以直接得到密度动能总能质量分数.压力和静焓都是温度的函数,可使用Newton迭代法得到温度.改写混合流体状态方程和总能表达式,得到:考虑有限差分法离散,对流项采用5阶特征紧致-WENO混合格式[19],输运项采用6阶对称紧致格式,时间推进采用3阶显式Runge-Kutta方法.为此,改写基本控制方程如下:特征紧致-WENO混合格式是紧致格式和WENO格式通过加权而形成的,用以解决光滑流场中有间断的问题.紧致格式可用较少的模板点高精度逼近导数,但在处理间断时会出现严重的数值振荡; WENO格式可抑制这种振荡,但数值耗散大.于是,通过加权因子使两者结合起来,在光滑流场使用紧致格式,在参数梯度比较大的间断部分使用WENO格式,同时满足了对高精度和捕捉激波的算法要求.省略下标“c”,记基本控制方程x方向的对流通量F的Jacobi矩阵为A=∂F/∂U,则矩阵A有特征值λ(i),对应特征向量L(i),这些特征向量构成特征空间.将网格点xj 的对流通量投影到特征空间，得到整点特征向量wj.在xj+1/2处计算半点守恒变量.对于输运项的离散也以x方向为例.记,u″j分别为(∂u/∂x)j,[∂(μ∂u/∂x)/∂x]j的近似差分,则可由6阶对称紧致格式离散为对于时间推进,一种3步3阶显式Runge-Kutta方法为选取一维Sod激波管中物质界面变化问题.初始条件为选取一维激波与气泡相互作用问题,一道激波沿y方向穿过一个静止气泡.初始条件为5阶WENO格式是经典的被广泛用来捕捉可压缩流动中激波等复杂间断的高精度格式.以上两个典型的算例均对比显示出本文算法的计算精度与其一致,表明所引入的5阶特征紧致-WENO混合格式对于激波等复杂间断的捕捉也具有很强的适用性和有效性.二维物理模型如图3所示.计算区域中,界面上层流体A为H2O和N2的混合物,Ab 表示激波后的流体,Af表示激波前的流体; 界面下层流体B为H2,O2,N2的混合物.平面正激波从流体A中向流体B中运动,在距离激波一定距离处有受单模态扰动的呈一定几何分布规律的界面,a为界面扰动的振幅,λa为界面扰动的波长.计算区域为矩形区域,宽XL为0.0314 m,高YL为0.15 m.正弦型流体间断界面为[20]上下边界采用流动物理量外插边界条件,左右边界采用周期性边界条件.网格为均匀分布正交网格,网格数为100 500.流体A中物质的质量分数分别为:wH2O=0.135,wN2=0.865,密度ρA=0.092kg/m3; 流体B中物质的质量分数分别为:wH2= 0.015,wO2= 0.12,wN2= 0.865,密度ρB= 0.476 kg/m3.上下层流体等压,均为0.5 atm.为便于参变对比分析,设定基本工况并加以计算分析.取激波Mach数Ms=1.3,界面扰动的振幅a=0.005 m.图4为瞬时密度等值云图时间序列，图5为相应的瞬时纹影时间序列图,计算公式为由图可以清晰地观察到激波与界面相互作用的过程.激波后流体的密度变大,激波与界面相遇后继续入射.t=0.02 ms后,激波离开界面,平面入射波发生微弱变形,而界面仍然是比较规则的正弦形.t=0.075 ms 后,激波基本恢复成平面波,界面逐渐发生变形.在t=0.125 ms后可以看到,界面已经完全失去正弦形状,轻流体以“气泡”形式进入重流体,重流体以“漏斗”形式进入轻流体,且两端发生卷起,形成“蘑菇”结构,“蘑菇杆”上出现一些小尺度结构.t=0.25 ms后,“蘑菇”结构非常明显,“蘑菇杆”上的小尺度结构会逐渐消失,流体混合剧烈.图6为瞬时流向速度等值云图时间序列.在t=0.02 ms 以后,界面两侧的流体存在流向速度梯度,所以沿界面切线方向两侧流体也存在流向速度梯度,引起Kelvin-Helmholtz(K-H)不稳定性[21],从而产生涡,使界面发生变形.在t=0.125 ms的界面发展后期,K-H不稳定性的作用更为明显,界面“漏斗”结构两侧在涡的作用下发生卷起,形成“蘑菇”结构.图7为瞬时涡量等值云图时间序列.t=0 ms时,激波与界面没有接触,流场中基本没有涡.激波与界面相遇后,流体界面附近产生涡,涡量有正有负,对称分布.随时间推移,涡量的大小和方向发生变化.“蘑菇”结构的顶部涡量最小,两侧最大.涡量表征流体旋转运动的强度,在涡的作用下,不同流体不断被卷吸掺混.在R-M不稳定性发展的前期,激波对界面的发展起主要作用; 在中后期,界面附近涡的产生和分布对界面变形和演化起关键作用.现分析涡动力学方程:这里,界面两侧流体密度的不同产生密度梯度,激波的作用产生压力梯度.涡线的拉伸和扭曲可以忽略,即第2项.因此,涡量的变化主要取决于斜压项散度项黏性项.图8为t = 0.02,0.075,0.25 ms不同时刻斜压项与散度项的空间分布对比,对应激波与界面相遇阶段界面线性发展阶段界面非线性发展阶段.对比发现,在激波与界面作用的前中后3个时刻,斜压项均比散度项大1个量级,表明在激波与界面作用的过程中,斜压项始终对涡的产生起主导作用.定义拟涡能ε为图9给出了激波与界面作用过程中拟涡能和界面振幅随时间变化的曲线.初始阶段,激波作用于界面,斜压作用剧烈,流场的拟涡能迅速增大; 激波离开界面后,拟涡能在一定时间内减小,然后又开始增大; 界面进入非线性发展阶段后,拟涡能发生振荡,但总体上继续增大.对于界面振幅,先减小后增大,且约在0.12 ms以前增长速率较大,随后增长速率变小.激波向着界面运动,二者没有接触时,界面不会发生明显的变化.激波与界面相遇时,在激波的强压缩作用下,界面振幅有所减小,图中曲线的第1个极值点是激波离开界面的时刻.激波离开界面后,斜压作用在界面附近产生的涡主导界面的变形,使界面振幅增长,增长速率较大且呈线性增长.界面发展到后期,伴随K-H不稳定性引起的界面“漏斗”结构向“蘑菇”结构的变化,界面振幅增长速率有所减小,不再呈线性变化,而是进入非线性增长阶段.基于基本工况a=0.005 m,再增加界面扰动的振幅a=0.0025,0.01 m算例,计算分析界面振幅对R-M不稳定性的影响规律.相对基本工况,增大界面振幅.图10为a=0.01 m 算例瞬时密度等值云图时间序列.激波与界面作用后继续入射,激波后流体密度增大,t=0.02 ms 前界面被压缩,激波离开界面后,界面上出现小尺度结构,界面振幅线性增长,随时间发展小尺度结构消失;t=0.075 ms时出现“蘑菇”结构,界面振幅继续增大,进入非线性发展阶段,“蘑菇”结构变明显,“蘑菇杆”结构变细,流体不断发生混合; t=0.25 ms后混合剧烈,“蘑菇”结构明显被破坏.同图4的基本工况相比,两种情况在界面附近涡的分布发生了变化.图11为a=0.01 m算例瞬时涡量等值云图时间序列.流场中的涡始终聚集在界面附近,t=0.075 ms 前流场中最大涡量较大,分布较集中; t=0.25 ms 后流场最大涡量变小,由于流体混合,分布较分散.同图7的基本工况相比,a=0.01 m时激波与界面作用使界面附近产生的涡量更大,且发展过程前后期最大涡量的变化也更大.这主要是由于a=0.01 m时混合进行得更剧烈,使界面间密度梯度减小,减弱了斜压作用.图12分别给出了a=0.01 m算例激波与界面作用过程中拟涡能和界面振幅随时间变化的曲线.初始时激波作用于界面,斜压作用剧烈,流场的拟涡能迅速增大; 激波离开界面后,拟涡能在一定时间内减小,然后继续增大,发展后期拟涡能发生振荡,但总体上仍保持增大趋势.界面振幅先减小后增大.在增大过程中,约0.08 ms以前振幅增长速率较大,为线性增长阶段; 0.08 ms之后振幅增长速率开始变小,进入非线性发展阶段.同图9的基本工况相比,拟涡能变化规律相似.相对基本工况,减小界面振幅.图13为a=0.0025 m 算例瞬时密度等值云图时间序列.界面首先被激波压缩,激波发生变形; 激波离开界面后,界面较长一段时间仍然保持正弦型,之后界面振幅开始增长; t=0.125 ms时出现“蘑菇”结构,t=0.25 ms 后混合较明显,但与另2种较大界面振幅情况相比,混合较弱.图14为a=0.0025 m算例瞬时涡量等值云图时间序列.流场中涡量集中在界面附近,混合发生得并不剧烈; 随时间发展,涡量最大值增大,发展后期最大涡量明显集中在“蘑菇”结构两侧.图15分别给出了a=0.0025 m算例激波与界面作用过程中拟涡能和界面振幅随时间变化的曲线.拟涡能发展趋势与另2种较大界面振幅情况相似.界面振幅先减小后增大,t=0.14 ms左右增长速率开始有所减小,进入非线性发展阶段.综合对比分析3种界面振幅情况可知,当激波Mach数相同时,a=0.01 m情况的界面最先进入非线性增长阶段,出现细小的“蘑菇”结构; 其次是a=0.005 m 情况; 最后是a=0.0025 m情况.就“蘑菇”结构诱导的混合强度而言,也呈逐渐减弱趋势.也就是说,界面振幅越大,界面线性发展阶段时间越短,进入非线性增长阶段越快,出现特征结构越快,流体发生混合也越剧烈.另一方面,激波与界面相遇时,不同界面振幅对应的拟涡能变化趋势及早期增长速率相差不多,但界面振幅越大,激波与界面相互作用的时间越长; 界面振幅越大,对应的拟涡能越大,拟涡能再增长速率越大.拟涡能的增强,强化了混合.在3种初始界面振幅情况下,界面振幅开始都会有一定程度的减小; 达到最小值后,在涡的作用下逐渐增大.在界面发展前期,界面振幅越大,界面的变化速率越大; 在界面发展中后期,界面振幅越大,界面振幅增长率越早开始减小,也越早进入非线性发展阶段,界面变形发展得越快.基于基本工况Ms=1.3,再增加激波Mach数Ms=1.15,1.5算例,计算分析激波强度对R-M不稳定性的影响规律.相对基本工况,增大激波强度.图16为Ms=1.5 算例瞬时密度等值云图时间序列.激波强度变大,激波前后流体密度比增大,界面运动速度也变大.当激波与界面相遇时,在激波压缩作用下,界面振幅变小,激波发生微弱变形,继续运动,经过一段时间恢复为平面波; 界面在线性增长阶段,维持正弦型特征,振幅增大; 然后进入非线性增长阶段,流场出现“蘑菇”特征结构.“蘑菇”结构两侧的混合最剧烈,混合导致混合区域流场中密度梯度变小,界面上的结构被破坏.图17为Ms=1.5算例瞬时涡量等值云图时间序列,也显示界面附近的涡很明显,最大涡量集中在“蘑菇”结构两侧.图18分别给出了Ms=1.5算例激波与界面作用过程中拟涡能和界面振幅随时间变化的曲线.拟涡能首先在激波作用下迅速增大,激波离开界面后,拟涡能开始减小,然后又开始增大,后期出现振荡,但总体保持增大趋势.界面振幅先减小,约在0.018 ms开始增大,在0.1 ms增长速率变小.相对基本工况,减小激波强度.图19为Ms=1.15 算例瞬时密度等值云图时间序列.由图可以看出,激波后的流体密度变大,但是激波前后流体密度比和界面运动速度相比于另2种较大激波强度情况明显减小,且界面变化减慢,接近t=0.25 ms 时界面才开始出现“蘑菇”结构,混合剧烈程度也小得多.图20为Ms=1.15算例瞬时涡量等值云图时间序列.与另2种较大激波强度情况相比,流场中的涡量明显较小,在界面发展后期,流场中最大涡量同样集中在“蘑菇”结构两侧.图21分别给出了Ms=1.15算例激波与界面作用过程中拟涡能和界面振幅随时间变化的曲线,它们随时间变化的规律与另2种较大激波强度情况基本一致.综合对比分析3种激波强度情况可知,当初始界面振幅相同时,Ms=1.5情况的界面最先从线性增长阶段进入非线性增长阶段,出现细小的“蘑菇”结构,两种流体在“蘑菇”两侧发生混合,界面结构被破坏.其次是Ms=1.3情况,最后是Ms=1.15 情况.激波在与界面作用过程中,若初始界面振幅相同,则激波强度越强,界面线性发展阶段时间越短,进入非线性增长阶段越快,出现特征结构越快,两种流体混合得越剧烈.类似地,激波强度越大,流场中拟涡能越大,且拟涡能的增长速率也越大; 界面振幅均先在激波的压缩作用下一定程度地变小,达到最小值后,在涡的作用下开始增长.在界面发展前期,激波强度越大,界面振幅增长的速率越大; 界面发展中后期,激波强度越大,界面振幅增长率越先开始减小,进入非线性发展阶段后,出现特征结构,此刻振幅增长速率会变小.总体上,激波强度越大,界面变形越快且越早进入非线性发展阶段,拟涡能伴随增大可导致流体混合发生得越剧烈.最后,以基本工况为例,将本文数值模拟得到的界面振幅随时间变化曲线与经典的脉冲模型[23]线性模型[24]非线性模型[25]等理论进行对比(见图22).从图中曲线的对比可得出:各曲线整体趋势一致,数值模拟结果与非线性理论结果定量上吻合较好,差异较小.因此,本文数值模拟得到的界面发展规律是合理的.采用数值模拟方法,详细研究了平面正激波与单模态正弦波界面相互作用导致的R-M不稳定性现象与强化混合流动物理机制.改变界面初始振幅和运动激波强度,比较了不同条件下R-M不稳定性的规律.主要得到以下结论:(1)当激波从轻流体运动到重流体中时,激波与界面作用,激波发生变形.由于激波的压缩作用,界面振幅小幅度减小.激波离开界面后,界面发展进入线性增长阶段,界面振幅增大,发展到后期,轻流体以“气泡”形式进入重流体,重流体以“漏斗”形式进入轻流体,出现“蘑菇”特征结构,界面发展呈现非线性发展规律,界面增速变小.由于斜压项和散度项作用,界面附近产生涡,诱导界面发生变形,斜压作用占主导地位.由涡诱导出界面的速度差,流动出现K-H不稳定性,使流场中产生更多的涡.(2)激波Mach数不变,比较界面初始振幅为a=0.0025,0.005,0.01 m时R-M不稳定性发展规律表明,增大界面初始振幅,能诱导界面附近产生更强的涡,界面变形发展更快,加快R-M不稳定性发展,更有利于流体混合.(3)界面初始振幅不变,比较激波Mach数为Ms=1.15,1.3,1.5时R-M不稳定性发展规律表明,类似地,增大激波强度,能诱导界面附近产生更强的涡,界面变形发展更快,加快R-M不稳定性发展,同样更有利于流体混合.。