爆炸地冲击载荷计算方法综述

周岸峰,李道奎,周仕明,周 旋,崔 达

(1.国防科技大学空天科学学院,长沙 410073;2.空天任务智能规划与仿真湖南省重点实验室,长沙 410073)

爆炸地冲击主要是指地下或地面发生爆炸时周围土石介质经冲击波的压缩和推动产生以加速度、速度、位移为表征的强烈地运动,也称爆炸地震或爆炸冲击振动[1-2]。与天然地震相比,爆炸地震对土石介质的作用以球对称压缩为主,主要成分为纵波,震源时间函数表现为单脉冲形式且能量急速释放的物理过程[3],具有幅值大、衰减快、频率高、持时短等特点[2]。在隧道、核电站和人防工程等建筑设施的设计中,通常需要考虑爆破施工、爆炸事故以及导弹攻击等极端环境[4-5];为了对结构在极端环境下的动态响应与损伤进行预测和评估,必须实现爆炸地冲击载荷的准确计算。

目前,国内外研究爆炸地冲击载荷的主要方法有理论分析、实验研究和数值仿真[6]。理论分析通常采用比较理想化的近似模型,只能处理一些简单问题。爆炸实验具有瞬时性和破坏性,对测量精度要求非常高,一些超常规的爆炸实验因国际公约无法进行。数值仿真则具有很强的灵活性,可以处理许多复杂的问题,甚至模拟实验达不到的条件,且方便进行参数灵敏度分析,从而指导理论分析和实验研究。因此,大部分研究将上述方法相互结合、相互验证,提出了许多关于爆炸地冲击载荷的计算方法,包括爆炸地冲击峰值载荷的预测[2,7],反应谱的设计[8-9]以及时域过程的模拟[10-15]。

由于现有爆炸地冲击载荷的计算方法种类繁多,各方法的适用条件和应用范围不尽相同,为此对现有研究成果进行梳理和汇总,对爆炸地冲击峰值载荷的经验公式进行比较分析,对爆炸地冲击反应谱的设计以及时域过程的模拟进行归纳总结,从而为进一步的相关研究提供参考与建议。

1 爆炸地冲击峰值载荷的预测

爆炸地冲击峰值载荷主要指峰值位移Dm、峰值速度vm和峰值加速度Am,其大小与爆炸当量、爆心距、爆源埋深、爆炸传播的地质条件等参数密切相关。目前常用经验公式和人工智能算法来表示上述参数之间的复杂关系,从而实现对爆炸地冲击峰值载荷的预测。

1.1 经验公式

爆炸地冲击峰值载荷的经验公式一般为特定的基函数形式,以爆炸当量Q和爆心距R为主要变量,将其它参数的影响作用归纳到比例系数k和衰减系数α中去,然后通过实验数据识别出k和α。

在矿产开采、隧道开挖等国民经济建设领域,峰值速度是工程爆破安全评估的重要准则[16],相关经验公式最多,一些比较经典、常用的峰值速度经验公式[2]如表1所示。

表1 峰值速度经验公式

美国陆军编写的《常规武器防护设计原理》[17]给出的地冲击峰值载荷的经验计算公式为

(1)

式中:C为地冲击耦合系数,封闭爆炸时C=1,其余情况取值可参考文献[1];cp为纵波波速;α为衰减系数;g为重力加速度。cp、α的取值如表2所示。

表2 岩土特性参数

《美国空军防护结构设计与分析手册》[10]给出了不同条件下核爆地冲击峰值载荷的估算方法

(2)

式中:kd、kv和ka为与地质条件和爆炸方式相关的常系数(见表3)。

表3 公式(2)相关系数

此外,Murphy和Lahoud[18]研究了99起地下核实验的地冲击振动记录,分析了爆心距为0.25～600 km之间的500多个观测点的数据,采用以下基函数形式来拟合峰值载荷经验公式

(3)

式中:Y为待拟合的物理量(峰值位移、峰值速度或峰值加速度);k、m、n为待定常系数,具体如表4所示。

表4 Murphy-Lahoud经验公式相关系数

综合比较上述经验公式可知,虽然各个公式的构造形式接近,但由于不同学者对地质条件归类、实验数据处理和基函数定义的方法不尽相同,导致各经验公式的系数取值存在差异,实际计算结果的离散性较大。因此,在一些实际工程中,仍然需要根据爆炸现场的实测数据,对传统经验公式的常系数重新进行识别,或增加其它变量进行多元回归分析,从而提出符合特定工程环境的爆炸地冲击峰值载荷预测公式[19-20]。

1.2 人工智能算法

爆炸地冲击峰值载荷相关参数之间具有很强的非线性,因而采用具体的函数形式来描述这种非线性关系是十分困难的。随着计算机技术的发展,人工智能(artificial intelligence,AI)算法由于非常适合处理复杂的非线性关系,也逐渐应用于爆破工程领域[7]。研究表明,AI算法在预测爆炸地冲击峰值速度时比经验公式和多元回归模型的性能都要好,但由于算法本身的局限性、训练样本来源的单一性以及建模的差异性(见表5),故AI算法还只能针对特定的工程爆破问题进行峰值速度的预测。

表5 爆炸地冲击峰值载荷预测AI算法

未来的研究还需要:①理清影响爆炸地冲击峰值载荷的主要因素,如爆炸当量(或炸药质量)、爆心距、地质条件、埋深、炮孔设置等,建立适用于不同工程问题的AI算法预测模型;②搜集不同地质条件的爆炸现场测量数据,建立一个更广泛、更具多样性的数据库,以对AI算法预测模型进行充分的训练,扩大其适用范围。

2 爆炸地冲击反应谱的设计

与峰值载荷这类单一判据相比,反应谱由于考虑了结构自身固有动力特性对地冲击响应的放大作用,以此进行隔震设计一般来说更加安全。目前常用于爆炸地冲击反应谱设计的方法主要是分量谱法、三系数法和特征频率法[8,9,15]。其中,分量谱法是美国《抗核武器结构设计手册》提出的方法,绘制图谱的过程十分繁琐,需要给定参数较多且难以确定,实际使用存在一定困难,而三系数法和特征频率法相对简单方便,下面主要介绍这两种方法。

2.1 三系数法

三系数法是根据爆炸地冲击峰值加速度、峰值速度、峰值位移,并基于大量实测加速度信号的归一化反应谱计算结果,确定加速度放大系数Ka、速度放大系数Kv、位移放大系数Kd,并按下式确定最大加速度响应Sa、最大速度响应Sv和最大位移响应Sd:

Sa=KaAm,Sv=KvVm,Sd=KdDm

(4)

然后根据Sa、Sv、Sd在三对数坐标图上分别作加速度、速度、位移坐标的等值线(见图1中实线),从而构成爆炸震动现场环境的反应谱包络线。

图1 三系数法Fig.1 Three-coefficient method

三系数法形式简洁,构建方便,但由于加速度、速度、位移峰值的经验公式结果在运动学上难以协调,容易导致反应谱形状畸变,出现三角形图谱和沿频率轴偏移等现象(见图1中虚线)。

2.2 特征频率法

特征频率法是根据爆炸地冲击加速度的衰减规律、加速度放大系数及反应谱的特征频率,通过量纲分析及统计得出的反应谱设计方法,无需速度、位移参数即可在三对数坐标系中同时绘制加速度、速度和位移反应谱,其中加速度反应谱(见图2中实线)可以使用分段函数表示为[15]

图2 特征频率法Fig.2 Specific frequency method

式中:f为频率;F1和F2分别定义为第一和第二特征频率。特征频率计算公式和加速度放大系数可参考文献[8]和文献[9],此处不一一列举。

年鑫哲等[15]认为式(5)确定的设计反应谱在高频部分存在较大误差,以此合成爆炸地冲击时程曲线时高频分量将被不合理地放大;当频率较高(刚度较大)时,系统的最大动力响应是趋近激励峰值的,即动力放大系数应趋近于1(见图2中虚线)。因此,对式(5)设计反应谱曲线的水平阶段(即图2中f>F3时)进行如下修正:

(6)

式中:F3和F4为按照F1、F2的方式定义的第三和第四特征频率。

综上可知,频率特征值法由于控制了特征频率F1和F2,且只需要计算峰值加速度,故克服了三系数法容易产生畸变的缺点[8]。但现有的特征频率经验公式仅仅来源于少量实测数据的回归分析,置信度不高;不同学者采用的经验公式基函数不同,导致计算结果差异巨大,通用性不强;而且文献[15]也没有给出计算特征频率值的详细说明。因此,为建立适用范围更广的爆炸地冲击反应谱模型,还需结合爆炸现场的实测数据对反应谱的图谱特征进行合理的统计分析,一方面可以对现有特征频率的经验计算公式进行修正和改进,另一方面可以尝试采用人工智能算法建立合适的特征频率的预测模型。

3 爆炸地冲击时程的模拟

爆炸地冲击的峰值载荷和反应谱,应用于工程设计时十分方便快捷,但不能分析结构在地冲击作用下的动态响应规律,而爆炸实验研究的成本太高、危险过大;因此,准确模拟地冲击载荷的时程曲线对结构的动力响应分析和减振设计十分必要。目前,计算爆炸地冲击时程过程的方法主要有经验曲线、数值仿真技术以及非平稳随机过程方法。

3.1 经验曲线

《常规武器防护设计原理》给出的爆炸地冲击速度时程经验曲线公式为[17]

(7)

式中:t为时间;β为修正系数,一般可取0.4,t0为冲击波到达时间,t0=R/cp;式(7)对时间求导即为加速度时程曲线。

严东晋等[11]对《美国空军防护结构设计与分析手册》提出的两类典型波形进行了曲线拟合,给出了对应的加速度时程经验公式;其中,类型I加速度波形的表达式为

(8)

类型II加速度波形的表达式为

(9)

分析可知,文献[17]的经验公式十分简单,却忽略了地冲击载荷的爬升过程;而严东晋等[11]提出的经验公式虽然能快速确定不同工况条件下的地冲击近似波形,并为隔震设计提供一定依据,但由于实际地质条件的复杂性,上述近似波形尚不能完全反应实际爆炸地冲击时域过程的复杂变化[26],容易丢失一些高频信息。因此,将上述经验公式确定的近似地冲击波作为结构动力分析的激励输入,本身就存在较大误差。

3.2 数值仿真技术

爆炸冲击波在岩土介质中的传播过程是一类典型的非线性瞬态动力学问题。理论分析可以清晰地认识问题的本质,但只有在岩土动态本构模型、载荷和边界条件都比较简单时,才可能得到一些解析解[26-27];随着数学和计算机技术的飞速发展,有限元法、离散元法、物质点法、光滑粒子流体动力学法等数值方法被广泛应用于爆炸冲击问题的研究[6,28-32]。

由于岩土中的爆炸问题比空气和水中的爆炸复杂得多,为提高数值仿真精度,除了开发高精度的算法格式以外,还需要建立科学合理的岩土介质动态本构模型。目前,常用的岩土动态本构模型有泡沫模型、帽盖模型等弹塑性模型,以及HJC模型、TCK模型等损伤断裂模型[12]。此外,北京理工大学的宁建国团队[12]还给出了基于细观力学理论的混凝土材料动态本构模型以及微裂纹演化与汇合的本构模型。然而,上述方法和本构模型都建立在连续介质力学的框架内,而实际地表的岩土介质存在大量断层、裂隙等特征,爆炸后岩土块的旋转以及相互间的滑移、挤压和分离会导致岩土介质位移场的间断[5],因此,在未来的研究中,仍然需要努力解决以下问题:①多介质界面相互作用的高精度计算[6,28,31];②爆炸冲击波强间断现象的高精度计算[29-30];③高速、高压、高温等极端条件下材料动态本构关系以及损伤演化规律的数学表征[32]。

3.3 非平稳随机过程方法

由于爆轰产物、气体、岩土等介质相互作用的复杂性,即使进行多次相同条件的爆炸实验,同一观测点记录的地冲击波一般也不会完全重现,具有一定的随机性。因此,可采用随机振动分析的方法,将爆炸地冲击时程按非平稳随机过程进行模拟[13-15]。目前常使用包络函数、相位差谱或时变功率谱表征随机过程的非平稳特性,分述如下。

1) 基于包络线函数的非平稳随机过程模拟。假定爆炸地冲击加速度时程为调幅非平稳随机过程为

A(t)=L(t)X(t)

(10)

式中:L(t)为幅值包络线,是一确定性函数;X(t)为一均值为0的平稳高斯随机过程。

关于幅值包络线函数L(t),其数学模型较多,目前常用的有单指数型[33]、改进的单指数型[14]、双指数型[34]以及分段函数型[35]。其中,单指数型包络线函数忽略了幅值的增长过程,实际应用较少。改进的单指数型和双指数型包络线函数不仅弥补了单指数型的不足,而且常系数物理意义明确,计算简单方便;但两者的衰减系数与幅值峰值时刻相关,当保持幅值峰值时刻相等时,两种包络线在衰减阶段会出现明显差异,故只适用于特定工况。分段函数型的平稳段持续时间和衰减系数则具有一定的可调性,适用范围比较大,但这些参数还缺乏明确的计算方法,需要根据实验数据进行确定。

对于平稳随机过程X(t)的模拟,目前常用的是Shinozuka和Jan[36]提出的加权振幅谐波叠加法,其表达式为

(11)

Δω=(ωu-ωl)/N,ωk=ωl+(k-1/2)Δω

(12)

式中:Sxx(ω)为目标功率谱密度函数;N为三角级数项数;ωk为各离散化的子频带的中心圆频率;φk为[0～2π]之间服从均匀分布的随机相位,且φk与φj相互独立;ωu、ωl分别为Sxx(ω)在正频率范围内的上、下限,即X(t)的频带宽度为(ωu-ωl)。

显然,基于谐波叠加法进行非平稳随机过程模拟的先决条件是给定合适的目标功率谱。国胜兵等[14]和周鹏等[37]提出了能够反应装药量和爆距影响的爆炸地冲击功率谱模型,但该模型部分系数的物理意义不明确,必须借助实测地冲击信号进行确定且计算过程比较繁琐。杨佑发等[13]和年鑫哲等[15]则根据反应谱转功率谱的方法,实现了爆炸地冲击信号的模拟,但目标反应谱的设计仍然依赖于实测地冲击信号,因而不能普遍应用。

2)基于相位差谱的非平稳随机过程模拟。一些强天然地震(持时较短,能量集中的强震)的研究[38-40]表明,采用包络线函数来反映地震动的非平稳性具有一定的随意性,且不能表征频率非平稳特性;实际地震动的非平稳性是由其相位差谱决定的,而相位差谱与相位谱具有一一对应的关系,即地震动的非平稳性包含在其相位差谱之中。因此,基于相位差谱进行非平稳随机过程模拟时,无需进行强制调幅处理,即

(13)

式中:各参数的物理意义与上文相同,但随机相位之间的差值应符合特定分布规律。随机相位之间的差值定义为

(14)

目前,针对不同地区的天然强震,不同学者对大量实测数据进行统计分析,得到了相关相位差谱模型,如正态分布、对数正态分布、均匀分布、Beta分布等[39-40],并实现了天然强震的非平稳随机过程模拟。但对于爆炸冲击振动,目前尚无相位差谱的统计分析模型,该方法能否准确模拟爆炸冲击振动有待进一步研究和验证。

3)基于时变功率谱的非平稳随机过程模拟。除相位差谱外,还可以采用时变功率谱来模拟地震动的时频非平稳特性,其数学表达式为

(15)

式中:Gxx(t,ωk)为t时刻ωk频率点的时变功率谱值。

由式(15)可知,建立合理的时变功率谱模型是模拟爆炸地冲击时程的关键。为此,宋浩等[41-42]先后以Kanai-Tajimi模型和Kameda模型作为目标时变功率谱对爆炸地冲击加速度时程进行了模拟;结果表明该模拟信号的反应谱与实际地冲击时程的反应谱整体趋势比较一致,但局部存在明显差异。因此,未来的研究应该以大量的爆炸现场的实测信号为基础,深入分析爆炸冲击振动的时频特征,提出更符合其客观物理过程的功率谱或时变功率谱模型。

4 总结与展望

对现有爆炸地冲击载荷计算方法进行梳理发现,当前计算方法适用范围有限、通用性较差,准确性有待进一步提高。未来需开展以下研究工作:

1)充分收集和整理已有的实验数据,建立一个广泛的数据库;并详细研究爆炸当量、爆心距、地质条件、埋深、炮孔设置等因素的影响作用,提高人工智能预测模型的适用性。

2)基于爆炸现场的实测数据,对反应谱的图谱特征进行合理的统计分析,从而对现有的特征频率法进行修正和改进;此外,还可以尝试采用人工智能算法建立合适的特征频率预测模型。

3) 研究更高效、更精确的计算方法和更准确的岩土介质动态本构模型,以实现对爆炸全过程的高精度仿真。同时,加强对实验数据的统计和分析,提出符合爆炸地冲击的客观物理规律的时变功率谱模型或者相位差谱模型,进一步提高非平稳随机过程模拟的准确性。