0800中国力学学会.DOC

水下爆炸气泡膨胀及冲击波和气泡群相互作用的数值模拟 曹乐 张亚军 张竹鹤 徐胜利 （中国科学技术大学力学和机械工程系， 合肥 230026） 摘 要 采用2阶迎风TVD求解欧拉方程得到流场解，采用5阶WENO求解Level set方程追踪多流体界面，采用GFM方法[1]（ghost fluid method）处理流体内界面。本文研究水下爆炸气体产物膨胀（气泡群）以及冲击波和气泡群相互作用问题进行了数值模拟，给出了不同时刻压力和密度云图分布和指定点的压力时间曲线。结果表明：本文方法和计算软件可很好地捕捉多个水下爆炸气泡及其演化图象。 关键词：Level set 方法，迎风TVD格式，MGFM 方法 引 言 流体运动界面模拟在环境灾害、水下爆炸等领域中有广泛应用。高压比、高密度比界面模拟及其与冲击波相互作用问题的难度更大。其中，多爆源的水下爆炸问题就是一类典型情况。该问题包含：（1）流场计算和运动激波捕捉。已有成熟的高精度、高分辨率激波捕捉算法，如TVD[1]、WENO[2] 。（2）界面模拟。精确捕捉界面方法对界面附近流场和物理量分布影响很大[3]。结合不需要捕捉界面位置GFM（ghost fluid method）思想[4]，Fedkiew提出level set方法。为适合现有的计算格式，在界面两侧设定2-5个节点组成的窄带。伪流体（ghost fluid）和真实流体在窄带中并存。合理设定界面边界条件，即指定伪流体参数，可保证界面附近不产生伪物理振荡和波系。先后出现基于等压[5]和等熵[4]固定GFM方法。为适合水下爆炸场计算，LIU T G提出了MGFM（modified ghost fluid method）方法[6] [7]。本文利用MGFM方法，对水下多爆源的爆炸气体产物膨胀（气泡群）以及冲击波和气泡群的相互作用进行数值研究。 1.控制方程和数值方法 1.1控制方程 水下爆炸中，水为可压缩介质，采用无粘Euler方程描述其流场 其中，矩阵列向量含义同文献【8】。爆炸气体产物采用理想气体状态方程，但比热比取3。 水的状态方程采用Tait方程。 其中，，A1=105Pa，B1=3.31×108Pa， 1.2界面方程 采用Level－set方程组描述多爆源爆炸气体产物界面，该方法已成功用于气体冲击波和多界面相互作用问题的数值模拟。 其中，，t为时间，(u，v，w)为空间坐标x，y，z对应的流场速度分量。为第n个界面的距离函数（定义为界面到流场中任一点距离）。等于0表示界面位置。界面边界条件采用MGFM方法，见文献[7]。 1.3 数值方法 本文流场计算采用TVD格式[1] 本文采用5阶WENO格式求解界面方程[8] 2计算结果分析及讨论： 2.1 水下运动激波与串列气泡群相互作用 计算条件说明如下：计算域为，Levelset函数定义为 激波波前气体参数： ρ＝1000.0，u=0.0，v=0.0，p=1 .0 激波波后气体参数： ρ=1187.487，u=，v=0.0，p=8000.0 气泡参数： 　　　 ρ=1.0， u=0.0，v=0.0，p=1.0，γ=1.4。 边界条件：给定指定Mach数的运动激波波后参数，右边界取外推边界条件，上、下边界取外推边界条件。计算网格数取112x82。 针对水下左行运动激波和串列空气泡群相互作用，图1给出不同时刻密度与压力云图分布。图1（a）表明冲击波和气泡群尚未发生相互作用。图1（b）表明冲击波和气泡群的左边第一个气泡发生了相互作用，气泡形状发生畸变。类似与激波和圆柱发生相互作用，冲击波阵面形状也发生了畸变。在气泡区域，产生了反射激波。图1（d）表明冲击波已扫过右边的第二个气泡，但由于气泡群区域的冲击波强度下降，因此，第二个气泡形状尚未发生变化。图1（d）表明冲击波已传出计算域右边界，此时第二个气泡才发生变形，第一个气泡收到压缩，体积已近似为0，压力出现汇聚。严格地说，气泡消失不满足气泡质量守恒，只是气泡体积小于网格尺度。图1（e）表明第二个气泡受压后又产生反射冲击波。图1（f）表明冲击波第二个气泡的反射冲击波继续向左传播，该气泡受到进一步的压缩。 密度 压力 （a）nt=0 密 度 压 力 （b）nt=50 密度 压力 （c）nt=100 密度 压力 （d）nt=150 密度 压力 （e）nt=200 密