1绪论
1.1 课题研究目的
1.1.1 课题背景
冲击波效应是多种弹箭的主要破坏力之一。持续时间大于0.25ms压力值大于3.91KPa时就会对人员听觉器官造成损伤[1]。当冲击波的动压高于70KPa时对于性能最好的C类地面野战通信设备结构和大部分零件破坏, 失去修复价值[2]。超压大于100KPa 时各类飞机完全破坏, 当超压为50KPa-110KPa时, 能使地雷引爆, 无线电雷达站和各种轻武器受损或破坏[3]。所以冲击波场超压测试是评价武器系统、工程爆破的有效手段。爆炸产生的冲击波特性研究包括炸药爆轰过程以及其相应的物理效应、爆炸的特点、能量释放以及冲击波荷载形成的机理。武器系统在研制、改进、定型、验收等环境都需要进行冲击波超压的测试, 以获得威力评估参数。冲击波特点为冲击波上升沿陡峭, 信号带宽, 爆炸为瞬间单次过程, 爆炸时伴随高温高压强光强电磁干扰, 测试系统工作环境恶劣, 对测试系统的要求高[4]。
电测法是一种日渐成熟的压力测试方法，随着现代数字存储测试技术的发展，它已广泛用于工农业生产，特别是武器工业。传统的电测试方法一般采用“压力传感器——放大器——瞬态信号记录仪——计算机”的形式构成测试系统。这种测试系统存在几个不足之处：第一是信号传输电缆会产生“电缆效应”，有时严重干扰测试信号；第二是各级组成体分离不便于系统的标定和校准；第三是系统的测试设备较多，价格高，体积大，机动性差。还有受测试环境影响较大，操作麻烦等不足的地方。鉴于以上不足的地方，现在已经有了使用存储测试技术的微型测试系统。存储测试系统一般是把传感器、适配放大器、A/D转换器、存储器、控制电路、接口电路以及电源集合在一起的微型化测试系统，它具有微小体积，能耐高的冲击加速度，较高的环境温度及环境压力，不需引线，在被测体工作过程中把信号记录下来，然后回收装置，用计算机读出和处理数据。体积微小和测量时不受外界的电磁辐射是他们的最突出特点，使得过去不便于测量的场合能够进行直接测量[5]。
1997年，华北工学院测试技术实验室基于动态测试理论，提出一种冲击波超压的存储式测试方法。通过试验获取了相关数据，并对实验曲线进行了分析和处理，用系统辨识的方法给出系统的传输特性曲线，对于测试系统动态特性的补偿问题作了探讨。1998年3月，华北工学院测试技术实验室在冲击波测试方法上进一步进行了探索，根据测试
系统存在超调而无法给出冲击波超压峰值误差的难题，设计了兼顾被测信号的有效宽度和测试系统的特性的滤波器，能够得到被测信号误差很小的估计值。1998年8月，北京理工大学提出了基于燃料空气炸药（FAE）爆炸压力场的基本特征及爆轰波、冲击波压力测量理论，系统地研究了FAE的爆炸压力测试技术。在综合分析研究多种FAE爆炸压力测试方法的基础上，提出了以内装电压放大器的测压系统代替电荷放大器的测压系统及数字存储式的压力测试系统为主要的两种测试方法。1999年，中国工程物理研究院也对存储测试方面提供了一些相关的看法。2002年，北京理工大学针对爆炸冲击波存储测试问题，提供了一种新型数字存储式测压系统设计思想，建立了存储测试系统信息传输模型及测试通道信息有效性准则，引入了信息的真和率失真函数概念进行测试通道设计[6]。
随着专用集成电路ASIC和用户现场可编程门列阵（FPGA）的飞速发展，针对冲击波测试的存储测试系统也向着微型化、大容量、模块化和智能化的方向发展。华北工学院测试技术实验室近年以来一直从事冲击波测试技术的研究，目前在冲击波的存储测试方面已经取得了一定的发展,采用了直接测试和间接测试一体化的设计，系统功能上也达到了智能化和自适应的特点。无线通信在存储测试技术的应用，进一步简化测试操作，提高测试的可靠性。
1.1.2课题的国内外发展状况
国内在爆炸问题方面的研究起步较晚，在此方面目前还比较匮乏。在理论研究方面，我国的学者李翼祺对爆炸问题进行了详细的分析与探讨，在其编著的《爆炸力学》一书中，作者陈列了各个不同学者对于爆炸问题的研究成果，并且将其与中国的国防规范进行了对比。阐述了不同地域对于爆炸问题的不同确定方法[7]。李夕兵，凌同华，张义平等人对于爆炸的信号采集及其处理进进行了详细的研究[8]，并对比了傅里叶变换，小波变换，以及目前很流行的一种新的信号处理方法EMD-HHT法，即对于爆炸信号进行自适应的EMD分解，并对其进行Hilbert-Huang变换，得到其能量谱，书中还分别指出了三种方法对于爆炸问题分析的优缺点以及其适用范围。
由于数值模拟可以很好的模拟爆炸问题的全过程，在国内，更多的学者选择使用数值模拟的方法进行科学研究。武汉大学潘超研究了并使用了多物质欧拉算法[9]，证明了该算法适合用于复杂流体的运动分析。采用多物质欧拉算法，对近水面爆炸进行模拟，可以很好的计算爆轰产物、水、空气的相互作用过程，证明了该方法可用来分析复杂流
体的运动可以得到较好的结果。武汉理工大学的申祖武对爆炸冲击波的环流效应进行了数值模拟[10]，利用任意拉格朗日欧拉(ALE)算法和炸药爆轰产物的JWL状态方程，针对有障碍物阻挡时空气冲击波的环流现象进行了数值模拟，得到了在爆源周围有障碍物时爆炸场初始发展、环流情况及整个传播过程,并对空气冲击波的环流规律进行了分析。研究结果表明，该种方法可以很好地描述空气冲击波环流流场参量分布与变化规律，数值模拟结果也符合客观物理规律。中北大学张广福对爆炸冲击波在无限空气领域传播的数值模拟进行了研究[11]，利用LS-DYNA程序，研究了爆炸发生后爆炸冲击波在无地面障碍物阻挡时的传播规律,并将数值模拟结果和经验公式的计算结果进行了对比。模拟结果表明,炸药近地爆炸产生的峰值超压与经验公式的计算结果基本一致。
国外学者在爆炸问题的理论研究方面开展较早，并在一些方面已经获得了较系统的研究成果：在爆炸理论的研究方面，Johansson详细的论述了炸药爆轰的整个过程及其中产生的相应的物理效应。Fickett对爆轰的特点进行了严密的数学推导，并对其中的物理现象和化学现象两者构建了严格的数学模型。BrodeH.J.[12]根据理论推到，计算了实际空气中的点爆炸的数值解。Brode对球形冲击波的传播进行了计算，假设在初始条件取为点源爆炸和等温球的情况下，对空气冲击波的传播规律进行了研究。W.E.Baker 在《空中爆炸》一书中对空中爆炸的理论和现象进行了详细的阐述，并在书中对爆炸实验的操作与设计方法进行了详细的讲解，同时给出了大比例尺的比例空中爆炸波参考图。Henrych综述了炸药爆炸现象及爆炸对各种介质和结构系统的破坏效应，并且通过大量的试验给出了不同比例距离下TNT炸药的入射超压，反射超压以及正压作用时间等参数的经验公式。Edward J.Conrath等人对于爆炸荷载作用下的各种结构形式上的荷载及力学行为进行了研究，并同时对结构构件的安全问题做了详细的分析与研究，如对安全门或者安全窗等一些构件进行了关于爆炸方面的研究[13]。
国外一些大的软件公司都在各自的软件中集成了可以进行爆炸问题分析的动力学模块，如ANSYS，ABAQUS等大型商用软件都拥有各自的显式分析功能，除此之外，美国军事部门还基于大量的试验所测试数据，编制了CONWEP和BLASTX程序以用于计算爆炸冲击载荷，其中BLASTX程序不仅能考虑入射角度对壁面反射压力的影响，而且Mehdi SotudehChafi等人[14]对比了TNT与C4在开敞空间的爆炸的数值模拟，并对RHA金属板的爆炸冲击进行了模拟。证明了利用数值模拟方法可以准确的计算爆炸对于结构形成的冲击波荷载作用，以及冲击波在介质中的传播规律。
1.2 研究内容、步骤
1.2.1 研究内容
根据空气中TNT爆炸峰压计算的经验公式，得出实验的理论值。利用ANSYS软件中的Fluent模块，建立压力测试装置外流场分布的二维有限元模型，模拟稳态下平面流场中，不同粘性系数。不同流动模型，不同流速等初始条件下两种测试装置表面及其附近位置的流场分布。重点在于测试装置表面压力分布，压力传感器安装位置的压力变化规律，绕流典型剖面流速值、压力值和尾流流场结构等。将两种测试装置所得模拟结果与实验理论值比较，分析它们的测量精确度。
1.2.2 研究步骤
（1）分析爆炸压力场特点，主要是爆炸冲击波超压的求解及其特点，查阅相关文献，了解国内外压力测试装置结构特点、绕流流场数值模拟和实验的研究成果。大致了解本课题的研究内容、措施和目的。
（2）复习流体力学有关圆柱绕流的相关基础，学习计算流体动力学（CFD）的基础知识；学习Fluent软件，包括前处理（建模，划分网格，设置湍流模型、材料参数及边界条件）、求解（设置残差监测曲线、监测点），及后处理（应用相关后处理软件，比如tecplot、origin处理模拟得到的图像和数据）。
（3）利用Fluent软件，建立两种压力测试装置外流场分布的二维有限元模型，进行稳态数值仿真。模拟不同流速等初始条件下，测试装置表面压力分布、压力传感器安装位置的压力变化规律、绕流典型剖面流速值、压力值和尾流流场结构等。
（4）将两种压力测试装置的数值模拟结果分别与理论值比较，得出结论。尝试指出两者之间的区别联系或区别的原因；结合现有的理论知识，说明分析结果；得出具有创新意义的线索或理论。
（5）总结实验数据的处理和分析，认真分析试验中出现的问题，概括实验中的心得体会，给出由建设性和预见性的意见。完成论文初稿的撰写、修改以及论文的最后定稿和提交。
2 爆炸冲击波理论和Fluent数值模拟基础
2.1 爆炸冲击波基本理论
2.1.1 爆炸冲击波的产生和传播规律
当炸药在空气中或水中爆炸时，其周围物质直接受到具有高温、高速、的爆炸产物作用，它强烈地压缩着相邻的介质，使其压力、密度、温度等突跃式的升高，形成初始冲击波。在空气中，爆炸发生的初始时刻，炸药产生的能量，在瞬时转化为高温高压状态的气态爆轰产物，此气体急剧膨胀，并迫使周围空气离开它原来占据的位置，于是在此气体的前端形成了一层压力梯度变化巨大的压缩空气层，此时外层的压缩空气层就是我们常说的爆炸波。事实上，爆炸产生的全部能量几乎都转化为爆炸波的能量。随着冲击波不断地膨胀，气态爆轰产物的压力逐渐减小直到等于大气压，此时开始爆炸波不再由气态爆轰产物支持，与其脱离开并相对独立的向前传播。气态爆轰产物的质点由于惯性作用继续向前运动，当其压力降低到周围大气压以下时，跟随在爆炸压力波后面传播，形成爆炸稀疏波。由于周围空气的压力值较高，气爆轰产物将逐渐停止前进并开始向后运动；由于惯性作用，它们的压力又逐渐增加，直至略微超过大气压，并重新形成气态爆轰产物膨胀的条件。由以上所述可以看出，整个过程就是“气态爆轰产物-空气”系统的一种自由振[15]。图2-1给出了理想情况下炸药在无限空气中爆炸时的压力随时间变
图2-1空气中炸药爆炸理想P-t曲线
化曲线。冲击波以球形向外扩张，随距离起爆点距离的增大，波阵面的面积也迅速增大，
因此，即使没有能量损耗，单位面积上的能量携带量也会不断减少，直至恢复到空气的初始大气压；并且，冲击波的传播并不是等熵的过程，在波阵面上熵是不断增加的，因此传播过程中始终存在着空气受到冲击绝热压缩而产生不可逆的能量损失，冲击波越强，这种能量损耗越大。综上所述，冲击波在传播过程中波阵面压力是迅速衰减的，并且初始阶段衰减较快，后期衰减逐渐缓慢。
在超声速流场中，激波是最普遍的现象，在管道内或在飞行器绕流环境之中几乎无处不在。飞行器以超声速飞行，对周围环境的扰动以压力波的形式传播，而压力波的传播速度为当地声速，可想而知以超声速运动的物体速度超过了压力波的传播速度，从而引起多重压力波的追赶与汇集，在飞行器前端形成了超强的波前，称之为激波。爆炸冲击波也是超声速流动，在爆炸波传播过程中遇到障碍物时也会产生激波。激波的外形取决于绕过物体的外形和来流马赫数，因此描述激波形状的方程式如果已确定，那么根据相应的一阶导数，则可确定激波角度和对应的激波强度。最简单的两种激波为正激波与斜激波。正激波与来流方向垂直，波后马赫数小于1，为亚音速流；斜激波与来流方向不垂直，倾斜一定角度（β角），称为激波角，通常称为超声速流，马赫速大于1，但在激波角接近900时，波后也可能出现马赫数小于1的情况。
高速运动物体（相对）表面气体的黏性作用不可忽略，但是在流场的雷诺数相当高时，黏性影响层很薄，在层内由于黏性而引起的剪切力和热交换不可忽略，定义为边界层。在平板上可压缩边界层可视为二维流动，是最简单的边界层流动。在亚声速、超声速、高超声速流中，平板边界层流动特性存在着明显的差别。对于亚/超声速来流，在平板上的边界层厚度，沿流向增长缓慢，对于平板压力梯度影响较小，因此不计这种影响（本文应用测压装置测试爆炸超压正是考虑了这种假设）。
2.1.2 爆炸冲击波参数
人们对高爆炸药爆炸产生的入射冲击波的传播规律主要用冲击波压力、超压峰值、冲量、持续时间等冲击波参数来描述。目前，比较简单，但却是公认的、常用的方法，是用比例距离表达冲击波的各种参数。比例距离的定义为：
Z=
式中：R为测点与爆心之间的距离（m）；W等效TNT药量（kg）。
等效TNT 药量计算公式：i T i T
Q W W Q =? 式中：i W 为所用炸药的药量（kg ）；i Q 为所用炸药的爆热（kCal kg ）；T Q 为TNT 爆热
（ kCal kg ）；T W 为i W 折算的TNT 当量[16]。
不同的研究者对空气中入射冲击波的参数进行了研究，并提出了相应的经验预测方
法。下面简要叙述一下常见的几个经验公式[17][18][19]。
2.1.2.1 冲击波峰值超压
Brode （1955年）建议高爆炸药爆炸冲击波峰值超压（MPa ）的表达式为：
3230.670.110.09750.14550.5850.00190.011
SO SO SO P Z P P Z Z Z ?+>??=??++-≤≤?? Henrych （1979年）建议空气中冲击波的峰值超压（MPa ）表达式为：
234
23231.407170.553970.035720.006250.050.30.619380.032620.213240.310.06620.4050.3288110
SO Z Z Z Z Z P Z Z
Z Z Z Z Z Z ?+-+≤≤???=-+≤≤???++≤≤?? Mills （1987年）介绍了高爆炸药冲击波峰值超压（MPa ）的一种表达式：
230.1080.114 1.772SO P Z Z Z
=-+ Crawford 和Karagoziam （1995年）给出了求解峰值超压SO P （MPa ）的表达式：
212
22240.4810(1434)(9.77)(10.55)SO O P R P R R R +=??+-??
式中O P 为大气压力（MPa ）。
文献[16]认为TNT 球形炸药在无限空气中爆炸时的冲击波峰值超压（MPa ）还可
按下式计算：
230.0840.270.7SO P Z Z Z
=++ 萨多夫斯基根据模型相似律理论建立公式，由试验确定系数，得到高爆炸药冲击波
超压峰值（MPa ）的表达式：
3231.070.110.0760.2550.65
115
SO Z Z P Z Z Z Z ?-≤??=??++<≤?? 由于爆炸问题的复杂性，冲击波超压峰值的各种经验公式预测方法难免会带来一些
偏差。小比利距离范围内，在直角坐标中不同预测方法的冲击波超压值与比例距离的关系，见图2-2所示。
图2-2经验公式的峰值超压对比
由图2-2可以看出，当比例距离Z 大于113m kg 时，各个公式预测的结果比较接近，
其中Mills 和Wu C ＆Hao H 预测的结果较其他四个公式稍微偏高。随着比例距离的减小，各个公式给出的结果的偏差逐渐增大，Mills 公式的值偏高，而Henrych 公式的值偏低，且在Z 小于0.813m kg 的情况下其结果较其他公式偏差较大。在小比例范围内，Brode 公式和文献[16]推荐的公式比较吻合。
通过对比分析可知，将文献[16]和萨多夫斯基推荐的公式综合起来能比较好的描述
冲击波超压与比例距离的关系，其表达式为：
23230.0840.270.7+10.0760.2550.65
115
SO Z Z Z Z P Z Z Z Z ?+≤??=??++<≤?? (a) 根据公式(a)的本文用于数值模拟的10组理论值：
2.1.2.2 爆炸冲击波速度D 与峰值压力SO P 的关系[16]：
由于空气冲击波在传播过程中逐渐均匀化的结果，不同装药形状的药包的作用与同
重量的球形装药相近。对空中球形TNT 爆炸而言，冲击波速度D 和峰值压力SO P 存在如
下关系式：
2
1002
07(1)6So D P P P P C =-=?-? 式中，0P 为未扰动的空气初始压力；0C 为未扰动的空气声速，对于不同温度有：
020.1C =（m s ），而0T 为未扰动的空气的初始温度（K ）。
注：同等重量的TNT 炸药在地面爆炸相当于2倍等效质量的炸药在空气中爆炸。
2.2 CFD 数值模拟理论基础
计算流体动力学(Computational fluid dynamics ，简称CFD)是通过计算机数值计
算和图像显示，对包含有流体流动和传热等相关物理现象的系统所做的分析[20]。CFD 的基本思想可以归结为：把原来在时间和空间上连续的物理量的场，如速度场和压力场，用一系列有限个离散点上的变量值的集合来代替，通过一定的原则和方式建立起关于这些离散点上场变量之间关系的代数方程组，然后求解代数方程组获得场变量的近似值。
2.2.1 流体力学控制方程
流体力学研究流体的宏观运动，它是在远远大于分子运动尺度的范围里考察流体运
动，而不考虑分子的个别运动，因此我们可以把流体视为连续介质[21]，它具有如下性质：
①流体是连续分布的介质，它可以无限分割为具有均布质量的宏观流体；②在不发生化学反应和离散等非平衡热力学过程的运动流体中，微元体内流体状态服从热力学关系；
③除了特殊面外，流体的力学和热力学状态参数在时空中连续分布的，并且通常认为是无限可微的。连续介质是一种力学模型，它适用于所考察的流体运动尺度L远远大于流体运动平均自由程l的情况，即：
L l?
1
流体流动要受物理守恒的支配，基本的守恒定律包括：质量守恒定律、动量守恒定律、能量守恒定律。如果在流体流动过程中出现了湍流情况，此时系统还应该将湍流输运方程考虑进去。控制方程是这些守恒定律的数学描述。下面就这些基本的守恒定律所对应的控制方程及控制方程的通用形式做如下介绍：
（1）质量守恒方程
在流体力学的分析中，要标识和跟踪固定的研究对象，通常是不方便的，所以采用场论的观点，在空间中标识一个明确的固定区域或体积，称为控制体。流场中任取一封闭曲面所围的控制体如图2-3所示，其体积为 ，表面面积为A，n为微元面积矢量dA 外法线方向上的单位矢量。
图2-3流场中的控制体
根据雷诺输运定理：s CV CS
dN d V dA dt t ηρτηρ?=+???????；体系内特征量N 关于时间t 的总变化率等于控制体内该特征量关于t 的变化率CV
d t ηρτ?????与单位时间内该变量通过控制面的净流出量CS
V dA ηρ???之和。这样，流体体系内某力学量随流体流动的变化情况可
以用该力学量对控制体的积分形式表出。借此，我们导出连续性方程对控制体的积分形式：
0CV CS
d V d A t ρτρ?+?=??????u r u r 即单位时间内，控制体中流体质量的减少量恰等于流体通过控制面所流出的总质量，满足质量守恒定律。
（2）动量守恒方程
动量守恒定律是任何流动系统都必须满足的基本定律。根据理想流体的Euler 运动
方程，微元体中流体动量对时间的变化率等于外界作用在该微元体上的各种力之和。
CV CS CV CS
V d V V d A Fd pnd A t ρτρρτ?+?=-???????????u r u r u r u r u r u r 相对应的微分形式的动量守恒方程为
1()grad V V V p F t ρ
?+??=-+?u r uu r u r u r 直角坐标系下的展开形式为：
111x y z u u u u p u v w F t x y z x v v v v p u v w F t x y z y w w w w p u v w F t
x y z z ρρρ??????+++=-+?????????????+++=-+?????????????+++=-+??????? （3）能量守恒方程
能量守恒定律是包含有热交换的流动系统必须满足的基本定律。该定律可表述为：
微元体中能量的增加率等于进入微元体的净热流量加上体力与面力对微元体所做的功。该定律实际是热力学第一定律。
()()E DE f V P V k T S Dt
ρρ=?+???+???+u r u r u r u r 或
()()()()E E EV f V P V k T S t
ρρρ?+??=?+???+???+?u r u r u r u r u r 本构方程：
'2()()23ij i j j
i ij j i P pI V I S S S g g u u S x x μμμ?=-+-??+???=?????=+????
u r r u r r u r u r u r u r 'μ为第二粘性系数，也称为膨胀粘性系数。
表征流体所处的热力学状态的宏观物理量叫做流体的热力状态参数，如p 、ρ、T 。
状态方程是描述流体宏观热力学参数之间关系的方程，其通式可以表示为：
(,)p p T ρ=
对于完全气体有，
m R p RT T M
ρρ== 其中，R 为特定气体常数，P V R C C =-，P V C γ=，m R MR =为通用气体常数，
m R =8.314()J K mol ?。
p
h e ρ=+；V de C dT =；P dh C dT =
2.2.2 Fluent 使用介绍[22][23]
Fluent 在当前应用广泛，是一个用于模拟和分析复杂几何区域内流体流动与热交换
的专用计算流体动力学软件。Fluent 提供了灵活的网格特性，用户可方便的使用结构网格和非结构网格对各种复杂区域进行网格划分。Fluent 还允许用户根据求解规模、精度和效率等因素，对网格进行整体和局部的细化或粗化。对于具有较大梯度的流动区域，Fluent 提供的网格自适应特性可让用户在很高的精度下得到流畅的解。Fluent 软件还推出了多种优化的物理模型，如定常和非定常流动、层流（包括各种非牛顿流模型）、湍流（包括最先进的湍流模型）、不可压缩和可压缩流动、传热、化学反应等等。针对
每种物理问题的流动特点，Fluent都有合适的数值解法，用户可对显示或隐式进行选择，以期在计算速度、稳定性和精度方面达到最佳。
Fluent使用C语言开发完成，支持基于MPI的并行环境。它通过交互的菜单界面与用户进行交互，用户可通过多窗口的方式随时观察计算过程和计算结果。计算结果可用云图、等值线图、矢量图、XY散点图等多种方式显示、存储和打印，也可传送给其他CFD或FEM软件。Fluent提供了用户编程的接口，让用户定制或控制相关的计算或输入输出。利用Fluent软件进行流体流动计算的流程：首先利用GAMBIT构建流动区域的几何形状，生成计算网格并设置边界类型，输出用于Fluent求解器的格式；然后用于Fluent求解器对流动区域进行求解计算，并进行计算结果的后处理。
Fluent可以计算二维和三维流动问题，主要使用范围包括：可压缩与不可压缩流动，瞬态和稳态流动，无粘流、层流及湍流，牛顿流及非牛顿流，对热交换（包括自然对流和混合对流），导热与对流换热耦合，辐射换热，惯性坐标系和非惯性坐标系下的流动问题模拟，用Lagrangian轨道模型模拟稀疏相（颗粒，水滴，气泡等），一维风扇、热交换器性能计算，两相流，复杂表面形状下的自由面流动。对于模拟复杂流场结构的不可压缩和可压缩流动，Fluent是理想的软件。
应用Fluent进行数值模拟，在确定所要解决问题的特征之后，需按以下几个基本步骤求解问题：①创建网格，输入并检查网格；②运行合适的解算器（2D或3D）；③选择需要解的基本方程：层流或湍流（或无粘流），化学组分或化学反应；④确定其他需要的模型，如传热模型等；⑤确定流体的材料物理性质；⑥确定边界类型及其边界条件；
⑦调节解的控制参数；⑧初始流场；⑨检查结果，保存结果；⑩必要的话，细化网格，改变数值和物理模型，重复⑧⑨步骤。Fluent是非结构解法器，采用非结构网格能缩短生成网格所需时间，简化几何外形和网格产生过程，具有使网格适应流场的特点。Fluent 可以在2D流动中处理三角网格和四边形网格，在3D流动中处理四面体网格、六面体网格、金字塔网格以及楔形网格。Fluent可以接受单块和多块网格，以及二维混合网格和三维混合网格，另外还接受由悬挂节点的网格。
2.2.3 湍流模型
可压缩流计算通常用于气流流动计算，液体流动在绝大多数情况下被认为是不可压的。衡量流体可压缩性的指标是马赫数。在马赫数小于0.3时，气体通常被认为是不可压流体；在马赫数大于0.3时，气体的压缩性影响逐渐增强，此时计算中必须考虑压缩
性影响。可压缩计算中各项参数的设置过程与其他流动计算大同小异，主要区别是计算中引入了粘性加热效应，同时需要在物质属性设置理想气体。可压缩流计算中需要设定的边界条件：
?流场入口：
压强入口条件：给定入口总温、总压，在超音速入口时，还需设定静压；
质量流入口：给定入口质量流和总温。
?流场出口
压强出口条件：给定出口静压。
因为湍流现象是高度复杂的，所以至今没有一种方法能够全面、准确地对所有流动问题中的湍流现象进行模拟。在涉及湍流的计算中，都要对湍流模型模型的模拟能力以及计算所需系统资源进行综合考虑后，再选择合适的湍流模型进行模拟。Fluent中采用的湍流模型包括Spalart-Allmaras模型、standard（标准）kε
-模型、RNG（重整化群）
-模型、RSM（Reynolds Stress Model，kε-模型、Realizable（现实）kε-模型、2v f
雷诺应力模型）模型和LES（Large Eddy Simulation，大涡模拟）方法。
Spalart-Allmaras模型最早被用于有壁面限制情况的流动计算中，特别在存在逆压梯度的流动区域内，对边界层的计算效果较好，因此经常被用于流动分离区附近的计算，后来在涡轮机械的计算中也得到广泛应用。最早的Spalart-Allmaras模型是用于低雷诺数流计算的，特别是在需要准确计算边界层粘性影响的问题中效果较好。Fluent对Spalart-Allmaras进行了改进，主要改进是可以在网格精度不高时使用壁面函数。在湍流对流场影响不大，同时网格较粗糙时，可以选用这个模型。
标准kε
-模型由Launder和Spalding提出，模型本身具有的稳定性、经济性和比较高的计算精度使之成为湍流模型中应用范围最广的一个模型。标准kε
-模型通过求解湍流动能（k）方程和湍流耗散率（ε）方程，得到k和ε的解，然后再用k和ε的值计算湍流粘度，最终通过Boussinesq假设得到雷诺应力的解。因标准kε
-模型假定湍流为各向同性的均匀湍流，所以在旋流（swirl flow）等非均匀湍流问题的计算中存在较大误差，因此后来又发展出很多kε
-
-模型的改进模型，其中包括RNG（重整化群）kε模型和Realizable（现实）kε
-模型在形式上类似于标准
-模型等衍生模型。RNG kε
kε-模型，但是在计算功能上强于标准kε-模型，其改进措施主要有：①在ε方程中
增加了一个附加项，使得在计算速度梯度较大的流场时精度更高；②模型中考虑了旋转效应，因此对强选择流动计算精度也得到提高；③模型中包含了计算湍流Prandtl数的解析公式，而不像标准kε
-模型仅用用户定义的常数；④标准kε
-是一个高雷诺数模型，而RNG kε
-模型在对近壁区进行适当处理后可以计算低雷诺数效应。Realizable kε-模型与标准kε-模型的主要区别是：①Realizable kε-模型中采用了新的湍流粘度公式；②ε方程是从涡量扰动量均方根的精确输运方程推导出来的。Realizable kε
-模型满足对雷诺应力的约束条件，可以在雷诺应力上保持与真实湍流的一致，这一点是标准kε
-模型都无法做到的，因此其能够更精确地模拟平面和圆形-模型和RNG kε
射流的扩散速度，同时在旋转流计算、带方向压强梯度的边界层计算和分离流计算等问题中，计算结果更符合真实情况。
-模型中包含了低雷诺数影响、可压缩性影响kω
-模型也是二方程模型，标准kω
和剪切流扩散，适用于尾迹流动计算、混合层计算、射流计算，以及受到壁面限制的流动计算和自由剪切流计算。剪切应力输运kω
-
-模型，综合了kω
-模型，简称SST kω
模型在近壁区计算的优点和kε
-模型
-模型在远场计算的优点，将kω
-模型和标准kε
都乘以一个混合函数后再相加就得到这个模型。在近壁区，混合函数的值等于1，因此在近壁区等价于kω
-模型；在远离壁面的区域混合函数的值则等于0，因此自动转换为标准kε
-模型中增加了横向耗散导数项，同时-模型。与标准kω
-模型相比，SST kω
在湍流粘度定义中考虑了湍流剪切应力的输运过程，模型中使用的湍流常数也有所不同。这些特点使得SST kω
-模型的适用范围更广，比如可以用于带逆压梯度的流动计算、翼型计算、跨音速激波计算等等。
2
v f
-模型中考虑到了壁面附近湍流的各向-模型与kε-模型比较类似，但是2v f
异性问题和非局部的压强与应变的关系。2v f
-模型属于低雷诺数湍流模型，其适用范围从自由流区一直延伸到壁面，并且无需适用壁面函数，2v f
-模型主要用于边界层计算和分离流计算。雷诺应力模型（RSM）中没有采用涡粘度的各向同性假设，从理论上说比湍流模式理论要精确得多；雷诺应力模型不采用Boussinesq假设，而是直接求解雷诺平均NS方程中的雷诺应力项，同时求解耗散率方程，因此二维问题中需要求解5个附加方程，在三维问题中则需要求解7个附加方程。从理论上说，雷诺应力模型应该比一方程模型和二方程模型的计算精度更高，但实际上雷诺应力模型的精度受限于模型
的封闭形式，因此雷诺应力模型在实际应用中并没有在所有的流动问题中都体现出其优势。只有在雷诺应力明显具有各向异性的特点时才必须使用雷诺应力模型（RSM），比如龙卷风、燃烧室内流动等带强烈旋转的流动问题。
湍流模型计算速度的快慢与计算量成反比，即计算量大则计算速度慢，需要的时间也长。湍流模型计算中的工作量主要取决于方程的数量和方程中函数项的多少。如果不考虑大涡模拟方法，湍流模型计算从总体上说，一方程模型（Spalart-Allmaras模型）计算最快，二方程模型（kε
-模型）次之，雷诺应力模型（RSM）
-模型、2v f
-模型、kω
最慢。
2.2.4 边界条件
边界条件就是流场变量在计算边界上应该满足的数学物理条件。边界条件与初始条件一起并称为定解条件，只有在边界条件和初始条件确定后，流场的解才存在，并且是唯一的。Fluent的初始条件是在初始化过程中完成的，边界条件则需要单独进行设定。边界条件大致可以分为下列几类：
?流场进出口条件：包括压强入口、速度入口、质量入口、吸气风扇、入口通风、压强出口、出口流动、出口通风和排气风扇等条件；
?壁面条件：包括固壁条件、对称轴（面）条件和周期性边界条件；
?内部单元分区：包括流体分区和固体分区；
?内面边界条件：包括风扇、散热器、多孔介质阶跃和其他内部壁面边界条件。内面边界条件在单元边界面上设定，因而这些面没有厚度，只是对风扇、多孔介质膜等内部边界上流场变量发生阶跃的模型化处理。
3 数值模拟与数据处理
3.1 测压装置测压原理
冲击波参数测试法（电测法）:它主要通过传感器将冲击波信号转换成电信号，通
过适配放大，再通过记录装置记录下来。通过测试结果来分析冲击波压力的大小和分布情况，而且能完整记录冲击波的传播过程。这种方法一般分为直接测试和间接测试两种。
① 压力-时间曲线直接测试
采用压力传感器将超压的压力信号转换成电信号，再通过信号放大、滤波等调理方
式，然后通过记录仪将信号记录下来，根据系统传输特性从获取信号中提取被测信号。这就是超压直接测试原理。
② 压力峰值间接测试
根据美国陆军实验规程，测量与爆心在一条直线上相距若干距离的两点的冲击波到
达时间，可获取对应的冲击波传播的平均速度。有了冲击波通过测试点的速度，可以通过兰基涅-胡果尼方程得到测试点的峰值压力。
20211r V K p p r C ??±????=?-??? ? ?+????????
式中：p 为峰值超压；r 为波前和波后的空气比热比；V 为冲击波的波前速度；K 为波前垂直的风矢量分量；C 为测试环境温度下的未扰动的空气声速；0p 为测试环境的未扰
动大气压力。
本文中测压装置的实际测压正是依据电测法思想，如图3-1是用于安装传感器的圆
盘形和楔形测压装置及其尺寸标注：
图3-1测压装置外形和尺寸
3.1.1 二维圆盘和楔形测压装置建模
GAMBIT软件是Fluent公司提供的前处理器软件，它包含功能较强的几何建模能力和强大的网格划分工具，可以划分出包含边界层等CFD特殊要求的高质量网格。GAMBIT 可以生成FLUENT6、FLUENT5.5、FIDAP、PLOYFLOW等求解器所需要的网格。用户可以直接使用Gambit软件建立复杂的实体模型，也可以从主流的CAD/CAE系统中直接读入数据。Gambit软件高度自动化，可以生成包括结构化和非结构化的网格，也可以生成多种类型组成的混合网格。Gambit生成网格的步骤一般为：①创建几何结构图：创建区域的节点，由节点连成直线，由边线创建面；②创建网格：启用边线网格设置操作分别对测压器模型边界及其外部流场边线划分成等距和非等距网格，由线网格生成面网格；③设置边界类型并输出网格文件。
依照Gambit生成网格的步骤，分别对圆盘形和楔形测压装置的流场区域按如下步骤进行：
?圆盘形测压装置：第1步，建立测压装置和外流场的几何形状，测压装置和外流场坐标：A（0,0），B（16,0），C（9.6，-1.9），D（6.4，-1.9），C1（9.8，-1.4），C2（9.6，-1.4），C3（9.8，-1.2），D1（6.2，-1.4），D2（6.4，-1.4），D3（6.2，-1.2），E（-80,0），F（-80,90），G（-80，-90），H（6.4，-90），I（9.6，-90），J（176，-90），K（176,0），L（176,90），M（16,90），N（0,90）。第2步，创建网格，测压装置划分等距网格（Successive Ratio为1，Interval count为节点数）：AB为36，CD为12，AD3、BC3为20，CC1、DD1为8，C1C2、C2C3、D1D2、D2D3为4，与测压装置对应的MN 划分等距36个节点，HI等距12个节点；外流场划分非等距网格（First Length项和Interval count项）：EF、EG、AN、BM、KL、KJ为0.1和92，AE、NF为0.15和80，CI、DH为0.1和80，HG为0.15和104，ML、BK为0.15和150，IJ为0.15和174，由线网格生成面网格。第3步，设置边界类型，inlet：EF、EG、NF、MN、ML、HG、HI、IJ，类型为pressure-far-field；outlet：KJ、KL，类型为pressure-far-field；wall1、wall2、wall3、wall4、wall5、wall6分别为AB、DC、AD3、BC3、[D3D2、D2D1、D1D]、[CC1、C1C2、C2C3]，类型为wall。最后得流场区域大小：X:-80—176（cm），Y：-90—90（cm），并将整个流场区域划分成6个不同大小的子流场区域，其中总节点数（Nodes）71976个，总单元数（Cells）71360个，总网格数（Faces）143336个。划分网格图3-2如下示：
图3-2圆盘测压装置网格图
(2)楔形测压装置：第1步，创建测压装置和外流场的几何形状，，测压装置和外流场坐标：A（0，0），B(5.459，0.974)，C（40.7，0.974），D（40.7，-0.974），E（5.459，-0.974），F（-160，0），G（-160，200），H（-160，-200），I（5.459，-200），J（40.7，-200），K（400，-200），L（400，-0.974），M（400，0.974），N（400，200），O（40.7，200），P（5.459，200），将各点连成线，右边线创建面。第2步，创建网格，测压装置划分等距网格（Successive Ratio为1，Interval count为节点数）：AB、AE为10，BC、ED、PO、IJ为60，CD、LM为6；外流场划分非等距网格（First Length项和Interval count项）：FG、FH、BP、CO、EI、DJ、MN、LK为0.1和110，AF为0.15和90，PG、IH 为0.15和100，ON、CM、DL、JK为0.15和160，由线网格生成面网格。第3步，设置边界类型，inlet：FG、FH、PG、PO、ON、IH、IJ、JK，类型为pressure-far-field；outlet：MN、ML、LK，类型为pressure-far-field；wall1、wall2、wall3、wall4、wall5分别为AB、BC、AE、ED、CD，类型为wall。最后得流场区域大小：X：-160—400，Y：-200—200，将整个流场区域划分成7个不同大小的子流场区域，其中总节点数（Nodes）71976个，总单元数（Cells）71360个，总网格数（Faces）143336个。划分网格数如图3-3所示：
图3-3 楔形测压装置网格图
3.1.2 两测压装置二维流场边界条件设置
参照文献[22]、[23] 中关于可压缩外部绕流的实例，应用Fluent 求解可压缩外部
绕流并发生边界层分离和激波现象问题的一般步骤：①对可压缩流动建模（密度使用理想气体定律）；②对外部绕流设置成无穷远边界条件；③选用合适的湍流模型；④使用耦合隐式求解器进行求解计算；⑤使用力和表面点监测器检查解的收敛性；⑥通过y+分布曲线检查网格；⑦外部绕流流动情况的后处理。
根据求解问题的一般步骤，对数值模拟圆盘测压装置和楔形测压装置外部绕流问题
进行具体设置。流场区域大小应根据被放置其中的障碍物尺寸设置，要求流域尾部边界应不小于10倍障碍物的横向尺寸，前流域大小约为5倍障碍物横向尺寸，由此创建流域大小如前面网格图所示。读入网格文件、网格检查、确定长度单位，因用于实验模拟的测压装置的数量级均为cm ，故具体实验设置时网格文件的长度单位也采用cm 。选择耦合、隐式求解器（Solve 中选择Coupled ，Formulation 中选择Implicit,其他保持默认设置）；选择热传导能量方程求解（即激活Energy Equation ）；在湍流模型设置对话框中选择标准k ε-模型。考虑到压缩性以及热物理特性随温度的变化，需对缺省设置进行修改：密度项（Density ）设置为理想气体（ideal-gas ），流体粘性项（Viscosity ）设置成适用于高速可压缩流动的气体黏度的Sutherland 定律，因为当流体密度和黏度均与温度有关时，Cp 和热传导率就不是常数了，而对于高速可压缩流动，一般来说都要考虑流体物理属性对热的依赖性，又由于温度梯度小，设Cp 和热传导率为常数的计算精度可以接受。
流场的边界可分为物理边界和人工边界两种，边界条件的设定十分重要，不仅会影
响模拟的精度，而且关系到求解过程的稳定和计算的收敛。鉴于理想气体状态方程中的压强为绝对压强，本论文中数值模拟的工作压强设置中操作压强（Operating Pressure ）设为0Pa ，不考虑重力影响。除测压装置模型外的其他流场边界均采用压强远场边界（pressure-far-field ）类型。因压强远场边界条件是在局部一维假设下通过引入黎曼不变量来定义的一类无反射边界条件，对于亚音速流动，存在两个黎曼不变量，分别对应于进入和离开流场的波：
21i i ni c R v γ=+-；21
n c R v γ∞∞∞=--