蓄意恐怖袭击事件和意外爆炸事故引发的建筑倒塌等灾害可导致巨大的人员伤亡和财产损失，因此，开展工程结构在爆炸荷载作用下的动态响应和损伤破坏研究尤为重要。其中，爆炸试验是重要的研究手段之一，目前，工程结构抗爆试验主要包括野外炸药爆炸试验^[1]和实验室爆炸波模拟试验^[2-3]。野外爆炸试验成本较高，存在一定的危险性，且试验数据受环境和测量因素的影响较大，具有一定的离散性。相比较而言，爆炸波模拟器可对爆炸冲击波荷载峰值和持时进行控制和调节，安全性较高，主要包括利用作动器加载的爆炸波模拟器^[2]和利用平面压力波加载的激波管^[3]。

20世纪90年代以来，为满足防护结构及武器装备足尺和大比尺的试验需求，国内外军工科研单位相继建造了大型爆炸波模拟器装置，如：法国Centre d'Etudes de Gramat试验场^[4-5]、美国白沙导弹靶场^[4-6]、英国Foulness试验场^[7-8]、德国Reiter Alpe Range试验场^[5]以及我国原总参谋部工程兵科研三所^[5]等。其中，英国Foulness试验场的大型激波管加载段内径达10.2 m，可进行框架等结构的缩尺模型抗爆试验。近年来，已有较多科研机构的激波管设备用于研究构件和复合材料在爆炸荷载下的响应，如渥太华大学^{[3, 9]}、卧龙岗大学^[10]、美国陆军工程研究与发展中心^[11]、米兰理工大学^[12]和挪威科技大学^[13]等。

激波管的驱动方式有高压气体驱动、炸药或氢氧爆轰驱动等，本研究主要关注高压气体驱动方式。高压气体驱动的激波管主要由驱动段、膨胀段和膜片系统组成。膜片系统将驱动段和膨胀段隔开。试验时，在驱动段内充入高压气体，膜片在两侧的压力差达到临界值时破裂，高压气体瞬间从驱动段流入膨胀段，可在膨胀段末端模拟爆炸时产生的超压荷载（峰值约为数百千帕，持续时间数十毫秒），并作用于构件或结构。

Lloyd^[9]基于渥太华大学的激波管开展了超压荷载与构件响应试验，研究了钢筋混凝土（reinforced concrete，RC）激波管驱动段的长度和超压与膨胀段末端反射超压峰值和持时的关系，并进一步对有/无抗震设计的RC柱开展了激波管试验。Ismail等^[14]采用商用有限元分析软件ABAQUS开展了数值模拟，分析了激波管驱动段的长度和超压对膨胀段末端超压荷载特性的影响，并通过激波管试验得到了驱动段不同长度和超压下反射超压峰值及反射冲量对应的半球形TNT装药的质量及爆炸距离。Jacques^[3]采用渥太华大学的激波管设备开展了未加固和纤维增强复合材料（FRP）布加固RC柱和板在爆炸荷载下的动力响应试验。

已有试验研究表明，激波管能够较好地复现均布的高峰值、短持时荷载，可用于工程构件和结构的抗爆试验，然而已有数值模拟和激波管设计的相关工作较少。本研究基于有限元分析软件ANSYS/LS-DYNA^[15]，对Jacques^[3]开展的RC板激波管试验进行数值模拟，通过对比膨胀段末端反射超压时程（作用于RC板的荷载）和RC板的挠度时程，对超压荷载的均布性以及采用的激波管建模和数值模拟方法的适用性进行讨论。进一步通过开展激波管不同设计参数（驱动段长度、直径、超压和膨胀段角度）对末端反射超压峰值和持时的影响分析，提出基于膨胀段末端反射超压峰值和持时的激波管设计方法，并开展相关的验证工作。

1.   激波管试验数值模拟

对Jacques^[3]开展的激波管试验进行精细化数值模拟，通过对比试验数据与预测结果（反射超压和RC板挠度时程），对激波管建模、参数取值和数值模拟方法的适用性进行验证，为激波管设计提供可靠的分析方法。

1.1   试验简介

渥太华大学激波管（图1）驱动段长度可根据试验荷载需求在305～5 180 mm范围内进行调节，可以产生反射超压为8～104 kPa、反射冲量为217～2 690 Pa·s的爆炸荷载。膨胀段长6 096 mm，其截面由直径597 mm的圆截面扩展为边长2 032 mm的方形截面，末端的刚性框架用于安装构件。

图  1  激波管示意图

Figure  1.  Schematic diagram of shock tube

下载: 全尺寸图片 幻灯片

Jacques^[3]基于上述激波管开展了RC板（尺寸为2 440 mm×2 440 mm×75 mm）的爆炸波加载试验。RC板内，双层双向布置11根直径为6.3 mm的钢筋，混凝土保护层厚度为6 mm，钢筋屈服强度为580 MPa，弹性模量为206 GPa。试验中两种RC板边界支撑条件分别为单向简支和双向固支，如图2所示，对应的混凝土单轴压缩强度分别为60.0和49.5 MPa。

图  2  边界支撑条件

Figure  2.  Boundary conditions

下载: 全尺寸图片 幻灯片

1.2   有限元模型

针对上述试验，考虑到对称性，建立如图3所示的1/4有限元模型。模型由激波管壳体、驱动段、膨胀段、激波管外围空气以及RC板组成。其中空气采用欧拉网格，RC板和激波管壳体采用拉格朗日网格，二者之间通过任意拉格朗日-欧拉流固耦合算法相互作用，由关键字*CONSTRAINED_ LAGRANGE_IN_SOLID实现。在数值模拟中，忽略膜片对气体流出的影响，通过合并驱动段和膨胀段连接界面网格节点对膜片进行简化，使气体理想地从驱动段向膨胀段流出。为避免在此过程中出现气体泄漏现象，激波管内部的空气网格外边界与激波管壳体网格内边界重合。在变截面膨胀段部分，激波管壳体和空气划分为渐变网格。RC板中混凝土为实体单元，钢筋为梁单元，混凝土与钢筋采用关键字*CONSTRAINED_BEAM_IN_SOLID相互耦合。支座与RC板之间的接触分别由关键字*CONTACT_AUTOMATIC_NODES_TO_SURFACE和*CONTACT_ AUTOMATIC_SURFACE_TO_SURFACE实现。

图  3  激波管试验有限元模型

Figure  3.  Finite element model of shock tube test

下载: 全尺寸图片 幻灯片

激波管外部有大量的横向和纵向的钢肋，刚度较大，可将激波管壳体以及支撑简化为刚体，其材料模型为*MAT_RIGID，并约束其空间平动和转动。激波管外围和膨胀段气体（空气）及驱动段高压气体均采用*MAT_NULL材料模型，二者均可视作理想气体，其状态方程为

式中：p为理想气体压强；${c_0}$～${c_6}$为与气体性质有关的常数；${e_0}$为初始单位体积气体内能；$\mu = \rho /{\rho _0} - 1$表示相对体积，$\, \rho$和${\,\rho _0}$分别为理想气体的密度和初始密度。对于理想气体，有${c_0} = {c_1} = {c_2} = {c_3} = {c_6} = 0$，${c_4} = {c_5} = 0.4$^[15]。对于空气和高压气体，${\,\rho _0}$和${e_0}$的取值见表1，其中n为高压气体压强与空气压强的比值。

表  1  理想气体材料参数

Table  1.  Parameters of gas

Gas type	p/MPa	${\rho {_0}}$/(kg·m−3)	${e{_0}}$/(MJ·m−3)
Air	0.1	1.29	0.25
High-pressure gas	0.1n	1.29n	0.25n

下载: 导出CSV 

| 显示表格

混凝土采用K&C模型（*MAT_CONCRETE_DAMAGE_REL3），该模型包括初始屈服面、极限强度面和残余强度面，考虑了强度硬化、软化效应和应变率效应^[14]，已广泛应用于爆炸荷载下混凝土结构动态响应数值模拟中。钢筋采用*MAT_PLASTIC_KINEMATIC模型，应变率效应采用Cower-Symonds模型描述，计算公式为

式中：f_y,d和f_y分别为钢筋的动态和静态屈服强度；$\dot \varepsilon$为应变率；C₀和P为应变率参数^[16]，分别取40.4 s⁻¹和5。钢筋密度为7 800 kg/m³，弹性模量为206 GPa，泊松比取0.3。

1.3   网格敏感性分析

针对试验中驱动段长度和超压为1 830 mm和80.7 kPa的工况，分别将高压气体离散为10、20和40 mm的实体单元网格，对应的膨胀段末端空气网格尺寸分别为35、70和140 mm，模型网格数量分别为586 000、80 000和10 000。试验和数值模拟得到的作用于RC板的反射超压时程如图4所示。可以看出，3种网格尺寸下，数值模拟结果与试验数据均吻合较好，荷载峰值误差分别为14.4%、18.7%和12.1%。为避免膨胀段末端空气网格与RC板网格尺寸差别过大，驱动段内高压气体网格尺寸取10 mm，此时膨胀段空气网格为10～35 mm的渐变网格。RC板中混凝土和钢筋的单元尺寸基于网格敏感性分析选取为12.5 mm，单元总数为660 000，如图3所示。

图  4  不同网格尺寸下膨胀段末端反射超压时程

Figure  4.  Reflected overpressure-time histories with different mesh sizes

下载: 全尺寸图片 幻灯片

图5为数值模拟得到的不同时刻激波管内部的压力传播云图，可以看出，高压气体由驱动段向膨胀段传播时可保持较好的平面度。图6进一步给出了膨胀段末端不同测点（图3）的反射超压时程，可以看出，不同测点的超压荷载时程基本重合，验证了膨胀段末端作用于构件的超压荷载为平面波。

图  5  不同时刻的压力云图

Figure  5.  Instantaneous pressure contours

下载: 全尺寸图片 幻灯片

图  6  作用于RC板的反射超压时程

Figure  6.  Reflected overpressure-time histories acted on RC slab

下载: 全尺寸图片 幻灯片

1.4   结果对比

1.4.1   超压荷载

为进一步验证激波管有限元模型预测膨胀段末端反射超压荷载的准确性，选取Jacques^[3]试验中4种典型工况（表2）进行数值模拟。相应的数值模拟结果与试验数据的对比如图7所示，反射超压峰值的试验值与预测值的相对误差均小于20%，表明建立的激波管有限元模型、参数取值和数值模拟方法可以准确预测激波管产生的超压荷载，其中由于试验中的采样误差导致图7(d)中出现幅值约为150 kPa的压力尖峰，文献[3]指出该工况的压力幅值为87.5 kPa，本研究中的数值模拟结果为91.38 kPa。此外，在S-4试验中，膨胀段末端采用封板封闭，而本数值模拟末端采用刚性封板。在数值模拟过程中为理想封闭条件，气体不会发生泄露，因而会出现冲击波的反射叠加效应，导致末端在40 ms后出现比较明显的二次峰值；而试验过程中封板出现了变形，气体外泄，无明显的二次峰值，如图7(d)所示。从试验结果与数值模拟结果的对比可以进一步看出，膨胀段末端的边界条件对超压的第一个峰值影响不大。

表  2  试验工况

Table  2.  Test cases

Case	Overpressure of driver section/kPa	Length of driver section/mm
S-1	61.4	305
S-2	80.7	1 830
S-3	697.1	2 745
S-4	637.1	4 880

下载: 导出CSV 

| 显示表格

图  7  数值模拟与试验得到的反射超压荷载曲线对比

Figure  7.  Comparisons of reflected overpressure-time histories obtained by simulation and experiment

下载: 全尺寸图片 幻灯片

1.4.2   RC板挠度响应

为进一步验证建立的有限元模型预测RC板位移响应的适用性，选取Jacques^[3]试验中3种不同RC板加载工况（表3）开展数值模拟。图8给出了试验和数值模拟得到的RC板中心点的挠度时程曲线对比，可以看出，3种工况下最大挠度误差分别为9.5%、2.3%和9.3%，进一步表明本研究建立的有限元模型可较好地预测激波管试验中RC板的动态响应。

表  3  RC板试验工况

Table  3.  Test cases of RC slabs

Case	Boundary conditions	Thickness of RC slabs/mm	Concrete strength/MPa	Overpressure of driver section/kPa	Length of driver section/mm
S-5	Simply-supported	75	60.0	80	1 830
S-6	Simply-supported	75	60.0	133	1 830
S-7	Fixed	75	49.5	125	1 830

下载: 导出CSV 

| 显示表格

图  8  数值模拟与试验得到的RC板挠度时程曲线对比

Figure  8.  Comparisons of deflection-time histories of RC slab obtained by simulation and experiment

下载: 全尺寸图片 幻灯片

2.   激波管参数对超压荷载的影响分析

如图9所示，激波管驱动段超压p_d、长度l、直径d和膨胀段角度θ均会影响激波管膨胀段末端的反射超压荷载。考虑到我国建筑层高约为3 m，梁、板和柱等构件跨度约3 m，为满足开展足尺试验的需求，并保证荷载的均布性，设计膨胀段末端为3 m×3 m的对称结构。基于第1节验证的激波管有限元模型、参数取值和数值分析方法，建立激波管有限元模型，并针对驱动段超压p_d和激波管几何参数（l、d和θ）对末端反射超压荷载峰值Δp和持时Δt的影响展开参数分析。上述参数的变化范围分别为：0.1 MPa≤p_d≤2.0 MPa，1 m≤l≤5 m，0.4 m≤d≤0.8 m，8°≤θ≤15°。

图  9  激波管参数

Figure  9.  Parameters of shock tube

下载: 全尺寸图片 幻灯片

2.1   驱动段超压

当激波管驱动段直径为0.6 m，长度可在1～5 m范围内变化，膨胀段角度为11°时，改变驱动段超压p_d进行数值模拟。图10给出了不同驱动段超压下末端反射超压时程曲线，结果表明，随着驱动段超压的增大，高压气体的传播速度升高，即高压气体到达末端时间提前。图11进一步给出了不同驱动段超压p_d下膨胀段末端反射超压峰值Δp和持时Δt对比。可以看出：反射超压峰值随驱动段超压增加而不断增加，而持时在p_d≥0.5 MPa后不再发生明显变化。需要指出的是：当驱动段长度l为1、2、3和4 m时，膨胀段末端的持时Δt为第一个正向压力持续时间（即超压起始时刻至超压降至零的时间段），而当驱动段长度为5 m时，无负压区，该持时为超压起始时刻到超压降至最小值的时间段。由于l=5 m时超压持时取值标准发生变化，所取持时较l=4 m时减少约2.7 ms。

图  10  不同驱动段超压下末端反射超压时程曲线

Figure  10.  Reflected overpressure-time histories at the end of driven section under different overpressures of driver section

下载: 全尺寸图片 幻灯片

图  11  不同驱动段超压下末端反射超压峰值和持时

Figure  11.  Peak reflected overpressure and positive time duration at the end of driven section under different overpressures of driver section

下载: 全尺寸图片 幻灯片

2.2   驱动段长度

当激波管驱动段超压在0.1～2.0 MPa范围内变化，直径为0.6 m，膨胀段角度为11°时，改变驱动段长度l进行数值模拟。图12给出了不同驱动段长度下膨胀段末端的反射超压时程曲线。可以看出：当驱动段超压为0.1 MPa时，随着驱动段长度的增加，末端反射超压的峰值不变，持时不断增加；在本研究模拟的其他工况下，即驱动段超压为0.5～2.0 MPa，随着驱动段长度的增加，膨胀段末端的超压反射峰值和持时有所提高。

图  12  不同驱动段长度下末端反射超压时程

Figure  12.  Reflected overpressure-time histories at the end of driven section under different driver section lengths

下载: 全尺寸图片 幻灯片

图13给出了不同驱动段长度下膨胀段末端反射超压峰值Δp和持时Δt。可以看出：随着驱动段超压的增加，其长度的增加对末端反射超压峰值的影响程度趋于显著，而对持时的影响逐渐减弱。结合图10可以看出，当驱动段超压增加时，驱动段和膨胀段压差增大，会导致：(1) 驱动段内高压气体流入膨胀段的速度增大，到达膨胀段末端形成超压的时刻不断提前；(2) 末端超压到达峰值后超压衰减较快，使超压结束时刻提前。当驱动段较短，如l=1 m时，前者的影响较大，超压起始时刻提前量更大，持时随着驱动段超压的增加而增加；当驱动段较长，如l=5 m时，后者的影响较大，超压结束时刻提前量更大，持时随着驱动段超压的增加而降低。

图  13  不同下驱动段长度下末端反射超压峰值和持时

Figure  13.  Peak reflected overpressure and positive time duration at the end of driven section under different driver section lengths

下载: 全尺寸图片 幻灯片

2.3   驱动段直径

由于激波管末端尺寸固定为3 m×3 m，膨胀段角度随着驱动段直径的增大而减小。当激波管驱动段超压为0.1～2.0 MPa，驱动段长度为3 m时，改变驱动段直径为0.4、0.6和0.8 m，其对应的膨胀段角度分别为11.9°、11.0°和10.1°。图14为相应的末端反射超压峰值Δp和持时Δt的计算结果，可以看出：在本研究模拟的工况下，当驱动段直径变大时，末端反射超压峰值不断增大，而相应的持时变化并不明显。

图  14  不同驱动段直径下末端超压峰值和持时

Figure  14.  Peak reflected overpressure and positive time duration at the end of driven section under different driver section diameters

下载: 全尺寸图片 幻灯片

2.4   膨胀段角度

进一步在驱动段超压为0.1～2.0 MPa，长度为3 m，直径为0.6 m条件下，改变膨胀段角度（15°、13°、11°、9°和8°）进行数值模拟。图15为不同膨胀段角度对应的末端反射超压峰值Δp和持时Δt。可以看出，在本研究模拟的工况下，随着膨胀段角度增大，膨胀段长度相应减小，导致末端反射超压峰值变大和持时降低。

图  15  不同膨胀段角度下末端反射超压峰值和持时

Figure  15.  Peak reflected overpressure and positive time duration at the end of driven section under different driven section angles

下载: 全尺寸图片 幻灯片

根据2.1节～2.4节共55个工况的计算结果，对末端反射超压峰值Δp和持时Δt与驱动段气体超压p_d、长度l、直径d和膨胀段角度θ的关系进行二次多项式拟合，表达式为

式中：A～O、A'～O'为拟合参数，如表4所示。需要指出的是，各变量（p_d、l、d和θ）的适用范围与参数影响分析中一致。

表  4  拟合系数

Table  4.  Fitted coefficients

A	B	C	D	E	F	G	H
–8.09	–16.88	405.54	–2.90	41.19	–93.76	9.75	542.06
I	J	K	L	M	N	O
–2.09	11.73	–11.78	–231.52	–326.53	31.10	–49.57
A′	B′	C′	D′	E′	F′	G′	H′
–0.61	1.22	–26.59	–0.20	–2.10	3.79	0.41	–11.95
I′	J′	K′	L′	M′	N′	O′
–0.34	2.57	1.62	13.90	5.97	0.99	1.24

下载: 导出CSV 

| 显示表格

3.   激波管设计算例

通过设计反射超压峰值和持时分别为130 kPa、17 ms和350 kPa、15 ms的两种目标工况对式(3)的适用性进行验证。基于式(3)选取的激波管参数和式(3)的计算结果如表5所示。进一步对设计的激波管进行数值模拟，得到的超压荷载峰值和持时同样在表5中给出。结果表明：算例D-1中，数值模拟与式(3)计算的反射超压峰值和持时的相对误差分别为6.3%和1.6%；算例D-2中，相应的相对误差分别为1.4%和7.6%。因此，式(3)可以用于激波管设计。

表  5  设计工况

Table  5.  Designed cases

Designed case	Δp/kPa				Δt/ms			l/m	d/m	p/MPa	θ/(°)
	Designed	Calc.	Sim.		Designed	Calc.	Sim.
D-1	130	136.57	127.97		17	17.99	17.70	3	0.7	0.5	10.55
D-2	350	352.25	357.10		15	16.08	14.85	3	0.7	1.5	10.55

下载: 导出CSV 

| 显示表格

4.   结　论

在验证激波管模型、参数取值和有限元分析方法适用性的基础上，进一步设计了膨胀段末端尺寸为3 m×3 m的激波管，开展了激波管几何参数和驱动段超压对膨胀段末端反射超压峰值和持时影响的参数分析。得到如下结论。

(1) 随着驱动段超压、长度、直径和膨胀段角度的增大，膨胀段末端超压峰值增加。

(2) 驱动段长度较长时，膨胀段末端超压持时随驱动段超压增加而缩短，超压持时受驱动段超压和长度的耦合影响；膨胀段末端超压持时随驱动段直径、膨胀段角度变化不大。

(3) 建立的激波管末端反射超压峰值和持时与驱动段超压、直径、长度和膨胀段角度的关系式可为激波管设计提供初步参考。

本研究只讨论了激波管几何尺寸和驱动段内超压对膨胀段末端超压荷载的影响，对于材料选型、结构强/刚度计算分析、连接件的强度校核等详细设计应进一步分析讨论。