黏弹性固体中地下爆炸辐射地震波能量的演化

卢强; 丁洋; 刘赟哲; 唐仕英; 郭志昀; 王占江

doi:10.11883/bzycj-2021-0058

黏弹性固体中地下爆炸辐射地震波能量的演化

doi: 10.11883/bzycj-2021-0058

西北核技术研究所，陕西西安 710024

基金项目: 国家自然科学基金（12072290）

详细信息

作者简介:
卢　强（1984－　），男，博士，副研究员，luqiang@nint.ac.cn

通讯作者:
王占江（1961－　），男，博士，研究员，wangzhanjiang@nint.ac.cn

中图分类号: O382.2；O347.4
计量
- 文章访问数: 670
- HTML全文浏览量: 251
- PDF下载量: 75
- 被引次数: 5
出版历程
- 收稿日期: 2021-02-07
- 修回日期: 2021-04-21
- 网络出版日期: 2021-08-24
- 刊出日期: 2021-09-14

Evolution of the radiated seismic wave energy of underground explosion in visco-elastic solids

Northwest Institute of Nuclear Technology, Xi’an 710024, Shaanxi, China

摘要

摘要: 地下爆炸与介质的能量耦合和介质中的波传播机制是理解地下爆炸源物理的重要基础。为研究地下爆炸辐射地震波能量的传播衰减规律，分析了黏弹性介质中地下爆炸地震波能量的组成。基于无限介质中黏弹性球面波理论，给出了速度、位移、应力、应变等物理量Laplace域的理论解。利用Laplace数值逆求解方法，建立了黏弹性介质中地下爆炸辐射地震波场的计算方法。以干黄土作为典型黏弹性材料，计算给出了地震波能量的传播特征，分析了地下爆炸辐射能量的传播衰减规律。结果表明：（1）在黏弹性介质中，某球面处流入的能量随半径增加而逐渐降低。在理想弹性介质中，某球面处流入的能量在几倍弹性半径外即可稳定到某一定值；（2）在某一固定的有限观测区域内，当观测时间足够长时，势能和耗散能均趋于某一定值，辐射动能趋于零；（3）当有限的观测区域能容纳一个完整波长的地震波时，地震波辐射动能的稳态值随波传播距离的增大而减小，总体上可以用指数函数和幂函数进行分段拟合。
- 地下爆炸 /
- 球面波 /
- 黏弹性 /
- 地震波能量 /
- 辐射动能
Abstract: The energy coupling between the underground explosion and the medium and the wave propagation mechanism in the medium are important bases for understanding the physics of the underground explosion source. In order to study the law of the propagation and attenuation for the seismic wave energy of underground explosion, the composition of the radiated energy of underground explosion in viscoelastic medium was analyzed, and the formulas for calculating inflow energy, outflow energy, radiated kinetic energy, potential energy and dissipation energy in a limited observation region were given. Based on the theory of viscoelastic spherical wave in infinite medium, the theoretical solutions of velocity, displacement, stress and strain in the Laplace domain were given by using the exponential attenuation pressure model and the generalized Maxwell viscoelastic model. The numerical solutions of velocity, displacement, stress and strain were given by using the Laplace numerical inverse method, and the inflow energy, outflow energy, radiated kinetic energy, potential energy and dissipation energy were calculated by these numerical results. The numerical results of different components of seismic wave energy are consistent with the theoretical results, and the correctness of this method is proved. Using the dry loess as typical viscoelastic material, the radial stress and particle velocity at different radii were calculated, and the relationship between the inflow energy at different radii and the wave propagation distance was obtained. The spatial distributions of the radiated kinetic energy of seismic wave were calculated by using the spatial distributions of the radial particle velocity at different times, and the propagation law of the radiated kinetic energy was obtained. The changes of the inflow energy and the radiated kinetic energy with the propagation distance in the limited observation area were analyzed, and the results show as follows: (1) In a viscoelastic medium, the energy flowing into a sphere surface decreases gradually with the increase of radius. In an ideal elastic medium, the energy flowing into a sphere surface at the elastic radius of about several times can be stabilized to a constant value. (2) The potential energy and the dissipative energy tend to constant values when the observation time is long enough in a fixed limited observation region, and the radiated kinetic energy tends to zero. (3) When a limited observation area can hold the seismic waves with complete wave length, the steady-state value of the radiated kinetic energy of seismic waves decreases with the increase of the wave propagation distance. In general, exponential function and power function can be used for piecewise fitting of the attenuation law for the steady-state value of the radiated kinetic energy of the seismic wave.
- underground explosion /
- spherical wave /
- viscoelasticity /
- seismic wave energy /
- radiated energy

HTML全文

研究地下封闭爆炸时通常把爆炸源简化为球形爆源，这样地下爆炸问题可简化为球面波问题。地下爆炸能量与源区介质发生强烈的耦合，最终以波的形式向周围介质辐射地震波能量。为表征地下爆炸或地震激发的震源强度，学术界提出了多种模型^[1-3]。

根据Mueller等^[4-5]提出的二阶震源模型，地下爆炸震源可用作用在弹性半径上的稳态压力叠加一个指数衰减的动态压力组成。在地下爆炸地震效应研究中，一般用折合位移势 $\psi (\tau )$ 或折合速度势 $\gamma (\tau )$ 的振幅谱表征震源强度^[6]，折合位移势稳态值 ${\psi _\infty }$ 和爆后介质弹性半径处的稳态压力（也称准静态压力）呈线性关系。同时，折合位移势稳态值 ${\psi _\infty }$ 还可以和地震学中反映源区不可恢复非弹性形变的地震矩 ${M_0}$ 联系起来。结合Murphy^[7]给出 ${M_0} \propto {\psi _\infty }$ 的线性经验公式，则爆后弹性半径处的稳态压力和地震矩 ${M_0}$ 也呈线性关系。在研究波传播路径上的震动效应时，地下爆炸空腔中的准静态压力和动态压力对震动效应都有贡献，仅用和准静态压力相关的地震矩或折合位移势稳态值表征地下爆炸震源强度是不完备的，不能全面反映爆炸空腔内动态压力在地介质内引起的震动效应。

Choy等^[2]和Boatwright等^[3]利用大量的地震事件，基于能量的概念建立了地震波能量和能量震级之间的关系，以此表征震源强度。袁乃荣等^[8]和李赞等^[9]利用国家地震台网和全球地震台网的宽频地震波数据对地震波的能量震级进行研究，指出地震波对传播路径区域造成的破坏直接与地震释放的地震波能量相关，可见地下爆炸向介质中耦合的地震波能量是表征地下爆炸震源强度的另外一个重要参数。

为研究地下爆炸与介质的能量耦合效应，Sanchidriár等^[10]利用多次野外岩石爆破，结合地震仪记录图、高速摄像机拍摄的岩石爆破面的运动速度以及爆破后岩石碎块尺寸的分布等实验结果，从地震波能量、动能、断裂能等方面研究了爆炸能量的组成。田振农等^[11]以块体离散元计算方法为基础，研究了岩体中爆腔内的压力脉动及爆炸能量分布特征，指出爆源近区岩体结构性质对爆腔内压力脉动规律影响不大，爆腔内压力峰值与炸药量无关，给出了岩石断裂破坏消耗能量、变形能、用于抛掷的能量的比例。郭家豪等^[12]基于反映高程放大的质点振动公式建立了侵彻钻地弹爆炸地震震源能量的计算模型，并结合经验公式推算爆炸总能量，从而评估钻地弹威力。吴亮等^[13]对耦合装药条件下柱状炸药爆炸的能量分布进行了研究，分析了爆炸能量和岩石介质的耦合机制，计算给出了爆炸能量与介质耦合过程中不同力学过程的能量占比。肖卫国等^[14]基于地下爆炸试验实测的地表速度，以折合速度势低频稳态值（即折合位移稳态值）表示的地震耦合强度，研究了黄土、砂砾石和花岗岩场地地下爆炸的地震耦合效应，指出对于相同当量、相当比例埋深下的地下爆炸，地震耦合强度和源区介质特性有着强烈的依赖关系。

本文中，利用广义Maxwell黏弹性体中球面波Laplace域的解，反演计算黄土中地下爆炸的速度场、位移场、应力场和应变场，以任意两个不同半径处球面之间区域作为观测区域，建立地下爆炸自由场辐射地震波能量的计算方法，分析地震波能量的传播演化规律，相关结果对于理解地下爆炸地震波能量在黏弹性介质中的衰减规律以及分析地下爆炸中实测的自由场数据具有一定参考作用。

1. 地下爆炸黏弹性区辐射能量的组成

球面波条件下，某个半径 $r$ 的球面上，其径向应力为 ${\sigma _{{r}}}(r,t)$ ，则径向面力 ${F_{{r}}}(r,t)$ 可以表示为 $4\text{π} {r^2}{\sigma _{{r}}}(r,t)$ 。球面波下材料仅有径向位移 ${u_{{r}}}$ ，则面力在径向位移上的做功 $W$ 可以写为^[10]：

$W(r,t) = \int\nolimits_0^{{u_{{r}}}} {{F_{{r}}}(r,t)d{u_{{r}}}(r,t)} = 4\text{π} \int\nolimits_0^\infty {{r^2}{\sigma _{{r}}}(r,t){v_{{r}}}(r,t)} {\rm{d}}t$

(1)

式中： ${v_r}$ 为径向粒子速度， $r$ 为爆心距。

力 ${F_r}(r,t)$ 在半径 $r$ 的球面上做的功 $W(r,t)$ 即是流入此球面的总能量。根据能量守恒，能量流经 ${r_0} \text{≤} r \text{≤} {r_1}$ 的区域时，流入能量 $W({r_0},t)$ 转化为此区域内所有质点的动能 ${E_{\rm{k}}}({r_0},{r_1},t)$ （即辐射动能）、介质应变能（也称内能或势能） ${E_{\rm{p}}}({r_0},{r_1},t)$ 、材料耗散能 ${E_{\rm{d}}}({r_0},{r_1},t)$ 以及从半径 ${r_1}$ 的球面上流出的能量 $W({r_1},t)$ ，如图1所示。则有：

图 1 黏弹性固体中爆炸地震波能量的组成

Figure 1. The composition of the energy of explosion seismic wave in visco-elastic solids

下载: 全尺寸图片幻灯片

$W({r_0},t) = W({r_1},t) + {E_{\rm{k}}}({r_0},{r_1},t) + {E_{\rm{p}}}({r_0},{r_1},t) + {E_{\rm{d}}}({r_0},{r_1},t)$

(2)

式中： ${E_{\rm{k}}}({r_0},{r_1},t)$ 和 ${E_{\rm{p}}}({r_0},{r_1},t)$ 分别为 $t$ 时刻整个区域 ${r_0} \text{≤} r \text{≤} {r_1}$ 内所有质点动能和应变能密度函数的空间积分，则有：

${E_{\rm{k}}}({r_0},{r_1},t) = 2\rho \text{π} \int\nolimits_{{r_0}}^{{r_1}} {\left( {{r^2}v_{{r}}^2(r,t)} \right){\rm{d}}r}$

(3)

${E_{\rm{p}}}({r_0},{r_1},t) = 4\text{π} \int\nolimits_{{r_0}}^{{r_1}} {{r^2}\left( {\frac{1}{2}{\sigma _{{r}}}(r,t){\varepsilon _{{r}}}(r,t) + {\sigma _{\rm{\theta }}}(r,t){\varepsilon _{\rm{\theta }}}(r,t)} \right){\rm{d}}r}$

(4)

式中： $\ \rho$ 为材料密度， ${\varepsilon _{{r}}}$ 为径向应变， ${\varepsilon _{\rm{\theta }}}$ 为切向应变， ${\sigma _{\rm{\theta }}}$ 为切向应力。

式（2）是基于能量守恒获得的，适用于任意固体中球面波能量传播演化的分析。若要求解球面波能量流经 ${r_0} \text{≤} r \text{≤} {r_1}$ 的区域时材料的耗散能 ${E_{\rm{d}}}({r_0},{r_1},t)$ ，则有：

${E_{\rm{d}}}({r_0},{r_1},t) = W({r_0},t) - W({r_1},t) - {E_{\rm{k}}}({r_0},{r_1},t) - {E_{\rm{p}}}({r_0},{r_1},t)$

(5)

考虑 $t \to \infty$ 的特殊情况，此时波完全离开 ${r_0} \text{≤} r \text{≤} {r_1}$ 区域，则此区域内质点动能 ${E_{\rm{k}}}({r_0},{r_1},\infty ) = 0$ 。按照黏弹性球面波的解，式（4）可以写为^[15-16]：

${E_{\rm{p}}}({r_0},{r_1},\infty ) = \frac{3}{2}\text{π} \int\nolimits_{{r_0}}^{{r_1}} {\left( {\frac{{\sigma _{{\rm{0s}}}^2r_0^6}}{{{G_{\rm{e}}}{r^4}}}} \right){\rm{d}}r} = \frac{{\text{π} \sigma _{{\rm{0s}}}^2r_0^6}}{{2{G_{\rm{e}}}}}\left( {\frac{1}{{r_0^3}} - \frac{1}{{r_1^3}}} \right) = \frac{{\text{π} \sigma _{{\rm{0s}}}^2r_0^3}}{{2{G_{\rm{e}}}}}\left( {1 - \frac{{r_0^3}}{{r_1^3}}} \right)$

(6)

式中： ${\sigma _{{\rm{0s}}}}$ 为稳态空腔压力， ${G_{\rm{e}}}$ 为黏弹性材料的静态剪切模量。从式（6）可以看出，稳态空腔压力 ${\sigma _{{\rm{0s}}}}$ 形成的势能主要分布在几倍空腔半径的范围内，比如当 ${r_1} = 10{r_0}$ 时，此观测区域储存了99.9%的势能。

由式（3）～（6）可给出球面波能量流经 ${r_0} \text{≤} r \text{≤} {r_1}$ 的区域时材料的耗散能 ${E_{\rm{d}}}({r_0},{r_1},\infty )$ ，有：

$\begin{aligned} & {E_{\rm{d}}}({r_0},{r_1},\infty ) = W({r_0},\infty ) - W({r_1},\infty ) - \frac{{\text{π} \sigma _{{\rm{0s}}}^2r_0^6}}{{2{G_{\rm{e}}}}}\left( {\frac{1}{{r_0^3}} - \frac{1}{{r_1^3}}} \right) =\\ &4\text{π} \left[ {\int\nolimits_0^\infty {r_0^2{\sigma _{{r}}}({r_0},t){v_{{r}}}({r_0},t)} {\rm{d}}t - \int\nolimits_0^\infty {r_1^2{\sigma _{{r}}}({r_1},t){v_{{r}}}({r_1},t)} {\rm{d}}t} \right] - \frac{{\text{π} \sigma _{{\rm{0s}}}^2r_0^6}}{{2{G_{\rm{e}}}}}\left( {\frac{1}{{r_0^3}} - \frac{1}{{r_1^3}}} \right) \end{aligned}$

(7)

2. 黏弹性固体中地下爆炸辐射地震波场的计算方法

根据无限介质中黏弹性球面波理论，在Laplace域，以波传播系数表征的黏弹性球面波粒子速度 ${\overline v _{{r}}}(r,s)$ 、粒子位移 ${\overline u _{{r}}}(r,s)$ 、径向应力 ${\overline \sigma _{{r}}}(r,s)$ 、切向应力 ${\overline \sigma _{\rm{\theta }}}(r,s)$ 、径向应变 ${\overline \varepsilon _{{r}}}(r,s)$ 、切向应变 ${\overline \varepsilon _{\rm{\theta }}}(r,s)$ 可写为^[15-16]：

${\overline v _{{r}}}(r,s) = - {\overline \sigma _{{r}}}({r_0},s)\dfrac{{s{\beta ^2}(s)\left[\dfrac{1}{r}\beta (s) + \dfrac{1}{{{r^2}}}\right]{{\rm{e}}^{{\rm{ - }}\beta {\rm{(}}s{\rm{)(}}r{\rm{ - }}{r_0}{\rm{)}}}}}}{{\dfrac{1}{{{r_0}}}\rho {s^2}\left[ {{\beta ^2}(s) + 2{k_{\rm{\mu }}}\left( {\dfrac{1}{{{r_0}}}\beta (s) + \dfrac{1}{{r_0^2}}} \right)} \right]}}$

(8)

${\overline u _{{r}}}(r,s) = - {\overline \sigma _{{r}}}({r_0},s)\dfrac{{{\beta ^2}(s)\left[\dfrac{1}{r}\beta (s) + \dfrac{1}{{{r^2}}}\right]{{\rm{e}}^{{\rm{ - }}\beta {\rm{(}}s{\rm{)(}}r{\rm{ - }}{r_0}{\rm{)}}}}}}{{\dfrac{1}{{{r_0}}}\rho {s^2}\left[ {{\beta ^2}(s) + 2{k_{\rm{\mu }}}\left( {\dfrac{1}{{{r_0}}}\beta (s) + \dfrac{1}{{r_0^2}}} \right)} \right]}}$

(9)

${\overline \sigma _{{r}}}(r,s) = {\overline \sigma _{{r}}}({r_0},s)\dfrac{{\dfrac{1}{r}\left[ {{\beta ^2}(s) + 2{k_{\rm{\mu }}}\left( {\dfrac{1}{r}\beta (s) + \dfrac{1}{{{r^2}}}} \right)} \right]{{\rm{e}}^{{\rm{ - }}\beta {\rm{(}}s{\rm{)(}}r{\rm{ - }}{r_0}{\rm{)}}}}}}{{\dfrac{1}{{{r_0}}}\left[ {{\beta ^2}(s) + 2{k_{\rm{\mu }}}\left( {\dfrac{1}{{{r_0}}}\beta (s) + \dfrac{1}{{r_0^2}}} \right)} \right]}}$

(10)

${\overline \sigma _{\rm{\theta }}}(r,s) = - {\overline \sigma _{{r}}}({r_0},s)\dfrac{{\dfrac{1}{r}\left[ {(1 - {k_{\rm{\mu }}}){\beta ^2}(s) - {k_{\rm{\mu }}}\left( {\dfrac{1}{r}\beta (s) + \dfrac{1}{{{r^2}}}} \right)} \right]{{\rm{e}}^{{\rm{ - }}\beta {\rm{(}}s{\rm{)(}}r{\rm{ - }}{r_0}{\rm{)}}}}}}{{\dfrac{1}{{{r_0}}}\left[ {{\beta ^2}(s) + 2{k_{\rm{\mu }}}\left( {\dfrac{1}{{{r_0}}}\beta (s) + \dfrac{1}{{r_0^2}}} \right)} \right]}}$

(11)

${\overline \varepsilon _{{r}}}(r,s) = {\overline \sigma _{{r}}}({r_0},s)\dfrac{{\dfrac{1}{r}{\beta ^2}(s)\left[ {{\beta ^2}(s) + \dfrac{1}{r}\beta (s) + \dfrac{1}{{{r^2}}}} \right]{{\rm{e}}^{{\rm{ - }}\beta {\rm{(}}s{\rm{)(}}r{\rm{ - }}{r_0}{\rm{)}}}}}}{{\dfrac{1}{{{r_0}}}\rho {s^2}\left[ {{\beta ^2}(s) + 2{k_{\rm{\mu }}}\left( {\dfrac{1}{{{r_0}}}\beta (s) + \dfrac{1}{{r_0^2}}} \right)} \right]}}$

(12)

${\overline \varepsilon _{\rm{\theta }}}(r,s) = - {\overline \sigma _{{r}}}({r_0},s)\dfrac{{\dfrac{1}{r}{\beta ^2}(s)\left[ {\dfrac{1}{r}\beta (s) + \dfrac{1}{{{r^2}}}} \right]{{\rm{e}}^{{\rm{ - }}\beta {\rm{(}}s{\rm{)(}}r{\rm{ - }}{r_0}{\rm{)}}}}}}{{\dfrac{1}{{{r_0}}}\rho {s^2}\left[ {{\beta ^2}(s) + 2{k_{\rm{\mu }}}\left( {\dfrac{1}{{{r_0}}}\beta (s) + \dfrac{1}{{r_0^2}}} \right)} \right]}}$

(13)

${k_{\rm{\mu }}} = \frac{{1 - 2\mu }}{{1 - \mu }}$

(14)

式中： $\ \mu$ 为材料的泊松比， $s$ 为Laplace变量， $\ \beta (s)$ 和 ${\overline \sigma _{{r}}}({r_0},s)$ 分别为波传播系数和弹性半径为 ${r_0}$ 的球腔上压力边界条件的Laplace表示。

由式（8）～（13）可知，若 $\ \rho$ 、 ${r_0}$ 、 $\ \mu$ 、 $\ \beta (s)$ 和 ${\overline \sigma _{{r}}}({r_0},s)$ 已知，则方程组完备。下面确定波传播系数 $\ \beta (s)$ 和弹性空腔半径上的压力边界条件 ${\overline \sigma _{{r}}}({r_0},s)$ 。广义Maxwell黏弹性模型可以用多个Maxwell体并联描述，如图2所示，其Laplace域表示可写为^[15]：

图 2 广义Maxwell体模型

Figure 2. The generalized Maxwell element model

下载: 全尺寸图片幻灯片

$\overline E (s) = \sum\limits_{i = 0}^N {\frac{{{E_i}{\theta _i}s}}{{{\theta _i}s + 1}}}$

(15)

式中： ${\theta _i} = {\eta _i}/{E_i}$ 、 ${\eta _i}$ 和 ${E_i}$ （ $i = 0,1,2, \cdots ,N$ ）分别表示第i个Maxwell体松弛时间、黏度系数和弹性模量， $\overline E (s)$ 为Laplace域的弹性模量，弹性元件的弹性模量 ${E_0}$ 、松弛时间 ${\theta _0} = \infty$ 。

波传播系数 $\ \beta (s)$ 和 $\overline E (s)$ 之间满足^[17]：

$\beta (s) = \sqrt {\frac{{(1 - 2\mu )(1 + \mu )\rho {s^2}}}{{(1 - \mu )\overline E(s)}}}$

(16)

Muller等^[4-5]通过对地下爆炸自由场测量结果的分析，认为球腔上压力边界条件的 ${\overline \sigma _{{r}}}({r_0},s)$ 可写为：

$\begin{aligned}&{\overline \sigma _{{r}}}({r_0},s) = L[{\sigma _{{r}}}({r_0},t) = ({\sigma _{\rm{d}}}{{\rm e}^{ - t/{t_0}}} + {\sigma _{{\rm{0s}}}})H(t)] = \\ &\frac{{{\sigma _{\rm{d}}}}}{{s + 1/{t_0}}} + \frac{{{\sigma _{{\rm{0s}}}}}}{s}\end{aligned}$

(17)

式中： ${\sigma _{\rm{d}}}$ 为动态压力峰值， ${\sigma _{{\rm{0s}}}}$ 为稳态空腔压力， ${\sigma _{\rm{d}}} + {\sigma _{{\rm{0s}}}}$ 为空腔压力峰值， ${t_0}$ 为压力衰减的时间常数， $H(t)$ 为亥维赛阶跃函数（Heaviside function）， $L$ 表示Laplace算符。

至此，把式（16）～（17）代入式（8）～（13），即可给出黏弹性球面波各参数的Laplace域的形式；通过对式（8）～（13）进行Laplace逆变换，即可求解出时域的粒子速度 ${v_{{r}}}(r,t)$ 、粒子位移、 ${u_{{r}}}(r,t)$ 径向应力 ${\sigma _{{r}}}(r,t)$ 、切向应力 ${\sigma _{\rm{\theta }}}(r,t)$ 、径向应变 ${\varepsilon _{{r}}}(r,t)$ 、切向应变 ${\varepsilon _{\rm{\theta }}}(r,t)$ ；把这些时域的解代入式（1）～（5），即可给出某个观测区域 ${r_0} \text{≤} r \text{≤} {r_1}$ 地下爆炸的辐射能量。

3. 地下爆炸辐射地震波能量的计算实例

以黄土中空腔爆炸为例说明地下爆炸辐射地震波能量的特征。参照文献[16]，取空腔半径为0.025 m，空腔边界压力为 ${\sigma _{{{{r}}_{{0}}}}} = - ({{\rm e}^{ - t/{t_0}}} + 1)$ MPa，其中 ${t_0} = {\theta _1}$ 。地震波的观测区域分别为 ${r_0}$ =0.025 m， ${r_1}$ =0.096 5，0.146 5，0.200 0 m，黄土黏弹性参数见表1。

表 1 黄土黏弹性参数^[18]

Table 1. The visco-elastic parameters of loess^[18]

ρ/（kg·m⁻³）	E₀/GPa	θ₀/µs	E₁/GPa	θ₁/µs	µ
1 800	1.60	∞	0.33	21.0	0.25

下载: 导出CSV

| 显示表格

图3给出黄土中地下爆炸辐射地震波能量的传播演化。对于某个特定的 ${r_0}$ 和 ${r_1}$ ，在 ${r_0} \text{≤} r \text{≤} {r_1}$ 的观测区域内， ${r_0}$ 边界上流入的能量最终趋于稳定值 $W({r_0},\infty )$ ；在 ${r_1}$ 边界上，当波到达此处，能量开始流出，最终趋于稳定值 $W({r_1},\infty )$ ；自 ${r_0}$ 边界上开始流入能量后，势能 ${E_{\rm{p}}}({r_0},{r_1},t)$ 和动能 ${E_{\rm{k}}}({r_0},{r_1},t)$ 即逐渐增加，当波开始流出 ${r_1}$ 边界时， ${E_{\rm{p}}}({r_0},{r_1},t)$ 和 ${E_{\rm{k}}}({r_0},{r_1},t)$ 突然减少，最终 ${E_{\rm{p}}}({r_0},{r_1},t)$ 趋于稳定值（如式（6））， ${E_{\rm{k}}}({r_0},{r_1},t)$ 趋于0， ${E_{\rm{k}}}({r_0},{r_1},t)$ 为0时说明波已完整地从 ${r_0} \text{≤} r \text{≤} {r_1}$ 观测区域离开；自 ${r_0}$ 边界上开始流入能量后， ${E_{\rm{d}}}({r_0},{r_1},t)$ 即逐渐增加，最终趋于稳定值 ${E_{\rm{d}}}({r_0},{r_1},\infty )$ 。 ${E_{\rm{p}}}({r_0},{r_1},\infty )$ 和 ${E_{\rm{d}}}({r_0},{r_1},\infty )$ 分别表示波完整通过 ${r_0} \text{≤} r \text{≤} {r_1}$ 区域后留在该区域介质的势能和材料吸收的耗散能。

图 3 黄土中地下爆炸辐射地震波能量随观测区域的变化

Figure 3. The variation of the energy of the radiated seismic wave from underground explosion in loess with the observation region

下载: 全尺寸图片幻灯片

当观测区域变大后，比如 ${r_0}$ 不变，随着 ${r_1}$ 由0.0965 m增大至0.2 m， ${r_1}$ 边界上流出的能量 $W({r_1},\infty )$ 越来越低；势能 ${E_{\rm{p}}}({r_0},{r_1},t)$ 、动能 ${E_{\rm{k}}}({r_0},{r_1},t)$ 和耗散能 ${E_{\rm{d}}}({r_0},{r_1},t)$ 在观测区域重合的部分其变化曲线也是重合的，观测区域越大则获得的能量传播信息越多； ${E_{\rm{d}}}({r_0},{r_1},t)$ 随着观测区域的变大， ${E_{\rm{d}}}({r_0},{r_1},\infty )$ 随之增大。

在地下爆炸实验研究中，一般是通过测量不同爆心距处粒子速度、应力等参量表征应力波的传播特性^{[7, 19-22]}，而地下爆炸辐射的地震波能量不能直接进行测量，但可以通过测得的粒子速度、应力等物理量计算。比如，通过同时测量某球面处径向应力和径向粒子速度波形，计算得到该球面处流入能量 $W(r,t)$ ；通过测量 ${r_0} \text{≤} r \text{≤} {r_1}$ 区域内不同位置的有限个粒子速度波形，利用一定的物理和数学处理方法^[23]，获得该区域的粒子速度场，从而可由式（3）计算得到 ${E_{\rm{k}}}({r_0},{r_1},t)$ ；通过测量 ${r_0} \text{≤} r \text{≤} {r_1}$ 区域内不同位置有限个径向应力、切向应力和径向粒子速度波形，结合一定的物理和数学处理方法构造径向应力 ${\sigma _{{r}}}(r,t)$ 、切向应力 ${\sigma _{\rm{\theta }}}(r,t)$ 、径向应变 ${\varepsilon _{{r}}}(r,t)$ 和切向应变 ${\varepsilon _{\rm{\theta }}}(r,t)$ 的场，从而利用式（4）计算得到 ${E_{\rm{p}}}({r_0},{r_1},t)$ ；在获得 $W(r,t)$ 、 ${E_{\rm{k}}}({r_0},{r_1},t)$ 、 ${E_{\rm{p}}}({r_0},{r_1},t)$ 之后，即可求得耗散能 ${E_{\rm{d}}}({r_0},{r_1},t)$ 。

4. 地下爆炸辐射地震波能量的传播特征

4.1 不同球面处流入能量的传播特征

为研究黏弹性固体中地下爆炸在不同球面处流入能量的传播特征，把观测区域取为 $r$ =0.025～10 m。

图4给出黄土中不同半径球面处最终流入能量 $W(r,\infty )$ 的变化。可看出， $W(r,\infty )$ 随半径 $r$ 增加而减小，并呈现出不同的变化规律。结合数据的变化特点，人为把0.025～10 m的空间范围划分为4个区域，并分段对 $W(r,\infty )$ 进行拟合。可以发现，在图中R1、R3、R4区域可近似用幂函数进行拟合，而在R2区域用指数函数拟合更符合 $W(r,\infty )$ 的变化规律。若用R1～R4这4个区域数据的拟合公式分别反推空腔壁上的能量参数，则预估的能量分别为真实输入能量的0.93、0.60、189、2.71倍，可见若用某局部观测区域的实测数据预估源区或远区波动参数时仍需慎重。

图 4 黄土不同半径球面处最终流入能量

$W(r,\infty )$ 的变化

Figure 4. The changes in the inflow energy

$W(r,\infty )$ at different radii of the sphere in loess

下载: 全尺寸图片幻灯片

在无黏性的理想弹性固体中，其耗散能为零，则可通过式（7）得到理想弹性材料中地下爆炸不同球面处流入能量：

$W(r,\infty ) = W({r_0},\infty ) - \frac{{\text{π} \sigma _{{\rm{0s}}}^2r_0^3}}{{2{G_{\rm{e}}}}}\left( {1 - \frac{{r_0^3}}{{{r^3}}}} \right)$

(18)

图4同时给出了不考虑介质黏性时 $W(r,\infty )$ 的传播图像，可看出，理想弹性固体中不同球面处流入的能量随着波传播距离的增加而逐渐降低，并在传播至几倍空腔半径处即基本保持不变。可见，黏弹性材料中地下爆炸在不同球面处流入能量的传播特征要比理性弹性材料复杂的多。

4.2 地震波辐射动能的传播特征

由公式（3）可知，为求得 ${r_0} \text{≤} r \text{≤} {r_1}$ 区域内的辐射动能 ${E_{\rm{k}}}({r_0},{r_1},t)$ ，需先求得此区域内 ${v_{{r}}}(r,t)$ 的空间分布。按照第2节的方法，图5给出了黄土中不同时刻粒子速度 ${v_{{r}}}(r,t)$ 的空间分布。可以看出，当 $t$ =100 μs时，在空腔壁上激发的应力波还未完全进入观测区域，随着时间的增加，波阵面一直向远离爆心方向传播，同时波长逐渐增加、幅度逐渐降低；当 $t$ =300 μs时，在观测区域内可以观察到完整的应力波。

图 5 黄土中不同时刻粒子速度的空间分布

Figure 5. The spatial distribution of the particle velocity at different times in loess

下载: 全尺寸图片幻灯片

利用式（3）计算给出不同时刻整个空间中总动能的空间分布，如图6所示，图中平台对应的即是地震波辐射动能稳态值 ${E_{\rm{k}}}$ 。当在观测区域内能观测到完整应力波时（取 $t$ =300～9000 μs），把波阵面对应位置 ${r_{{\rm{WF}}}}$ 设为横坐标、地震波辐射动能稳态值 ${E_{\rm{k}}}$ 设为纵坐标，则可以给出地震波辐射动能稳态值 ${E_{\rm{k}}}$ 随波前位置 ${r_{{\rm{WF}}}}$ 的变化曲线，如图7所示。与 $W(r,\infty )$ 类似，在能观测到完整应力波空间分布的区域内，地震波辐射动能稳态值 ${E_{\rm{k}}}$ 随着波传播距离 ${r_{{\rm{WF}}}}$ 的增加而减少。结合数据变化特征，在离爆心较近的区域，地震波辐射动能稳态值随波前位置的变化曲线可近似用指数函数拟合，而在较远区域，可用分段幂函数进行拟合。

图 6 黄土中不同时刻辐射动能的空间分布

Figure 6. The spatial distribution of the radiated kinetic energy at different times in loess

下载: 全尺寸图片幻灯片

图 7 黄土中地震波辐射动能随波前位置的变化

Figure 7. The variation of the radiated kinetic energy of seismic wave with the position of wave front in loess

下载: 全尺寸图片幻灯片

5. 结　论

由上述分析得到以下几点结论：

（1）利用建立的黏弹性固体中地下爆炸辐射地震波场的计算方法，给出了地下爆炸有限黏弹性观测区域中流入（出）能量、动能、应变能、耗散能的结果，数值计算结果和理论分析结果相符，说明本文计算方法的正确性。

（2）在黏弹性介质中，某球面处流入的能量 $W(r,\infty )$ 随半径增加而逐渐降低，同时呈现出复杂的衰减规律，大部分区域 $W(r,\infty )$ 的衰减规律可用幂函数拟合，而部分区域 $W(r,\infty )$ 的衰减更符合指数函数衰减特征；在理想弹性介质中， $W(r,\infty )$ 在几倍弹性半径外即可稳定到某一定值。

（3）在某一固定的有限观测区域内，当观测时间足够长时，势能和耗散能均趋于某一定值，辐射动能趋于零。

（4）当有限的观测区域能完整地容纳地震波时，利用粒子速度的空间分布可以计算给出此时地震波辐射动能的空间分布，地震波辐射动能的稳态值随波传播距离的增加而减少，且衰减规律比较复杂，总体上可以用指数函数和幂函数进行分段拟合。

图 1 黏弹性固体中爆炸地震波能量的组成

Figure 1. The composition of the energy of explosion seismic wave in visco-elastic solids

下载: 全尺寸图片幻灯片

图 2 广义Maxwell体模型

Figure 2. The generalized Maxwell element model

下载: 全尺寸图片幻灯片

图 3 黄土中地下爆炸辐射地震波能量随观测区域的变化

Figure 3. The variation of the energy of the radiated seismic wave from underground explosion in loess with the observation region

下载: 全尺寸图片幻灯片

图 4 黄土不同半径球面处最终流入能量 $W(r,\infty )$ 的变化

Figure 4. The changes in the inflow energy $W(r,\infty )$ at different radii of the sphere in loess

下载: 全尺寸图片幻灯片

图 5 黄土中不同时刻粒子速度的空间分布

Figure 5. The spatial distribution of the particle velocity at different times in loess

下载: 全尺寸图片幻灯片

图 6 黄土中不同时刻辐射动能的空间分布

Figure 6. The spatial distribution of the radiated kinetic energy at different times in loess

下载: 全尺寸图片幻灯片

图 7 黄土中地震波辐射动能随波前位置的变化

Figure 7. The variation of the radiated kinetic energy of seismic wave with the position of wave front in loess

下载: 全尺寸图片幻灯片

表 1 黄土黏弹性参数^[18]

Table 1. The visco-elastic parameters of loess^[18]

ρ/（kg·m⁻³） E₀/GPa θ₀/µs E₁/GPa θ₁/µs µ

1 800 1.60 ∞ 0.33 21.0 0.25

下载: 导出CSV

参考文献(23)

[1]	靳平, 徐果明, 楼为涛. 受低频动态正压力加载的椭球腔的地震矩张量表示及其在无限介质中辐射的地震波 [J]. 地震学报, 1997, 19(5): 447–456. JIN P, XU G M, LOU W T. Seismic moment tensor representation of ellipsoidal cavity loaded with low frequency dynamic positive pressure and seismic wave radiated in infinite medium [J]. Aata Seismologica Sinica, 1997, 19(5): 447–456.
[2]	CHOY G L, BOATWRIGHT J. Global patterns of radiated seismic energy and apparent stress [J]. Journal of Geophysical Research, 1995, 100(89): 18205–18228. DOI: 10.1029/95JB01969.
[3]	BOATWRIGHT J, CHOY G L. Teleseismic estimates of the energy radiated by shallow earthquakes [J]. Journal of Geophysical Research, 1986, 91(B2): 2095–2112. DOI: 10.1029/JB091iB02p02095.
[4]	MUELLER R A, MURPHY J R. Seismic characteristics of underground nuclear detonations: Part Ⅰ: seismic spectrum scaling [J]. Bulletin of the Seismological Society of America, 1971(61): 1675–1692.
[5]	MUELLER R A, MURPHY J R. Seismic characteristics of underground nuclear detonations: Part Ⅱ. Elastic energy and magnitude determinations [J]. Bulletin of the Seismological Society of America, 1971(61): 1693–1704.
[6]	周钟, 王肖钧, 肖卫国, 等. 花岗岩介质中地下爆炸震源函数研究 [J]. 爆炸与冲击, 2007, 27(1): 18–25. DOI: 10.11883/1001-1455(2007)01-0018-08. ZHOU Z, WANG X J, XIAO W G, et al. Study on the main characteristics of underground explosion seismic source function in granite [J]. Explosion and Shock Waves, 2007, 27(1): 18–25. DOI: 10.11883/1001-1455(2007)01-0018-08.
[7]	MURPHY J R. Free-field seismic observations from underground nuclear explosions [M]. Explosion Source Phenomenology: American Geophysical Union Monograph, 1991. DOI: 10.1029/GM065p0025.
[8]	袁乃荣, 刘瑞丰, 李赞, 等. 能量震级及其测定 [J]. 地震地磁观测与研究, 2018, 39(5): 1–7. DOI: 10.3969/j.issn.1003-3246.2018.05.001. YUAN N R, LIU R F, LI Z, et al. Energy magnitude and its determination [J]. Seismological and Geomagnetic Observation and Research, 2018, 39(5): 1–7. DOI: 10.3969/j.issn.1003-3246.2018.05.001.
[9]	李赞, 刘瑞丰, 孔韩东, 等. 中强地震能量震级测定 [J]. 地震学报, 2019, 41(3): 289–301. DOI: 10.11939/jass.20180139. LI Z, LIU R F, KONG H D, et al. Energy magnitude determination of moderate-strong earthquakes [J]. Acta Seismologica Sinica, 2019, 41(3): 289–301. DOI: 10.11939/jass.20180139.
[10]	SANCHIDRIÁR J A, SEGARRA P, LÓPEM L M. Energy components in rock blasting [J]. International Journal of Rock Mechanics & Mining Sciences, 2007, 44: 130–147.
[11]	田振农, 张乐文, 李世海. 岩体中爆腔内压力脉动特征和爆炸能量分布的数值模拟 [J]. 岩土工程学报, 2010, 32(8): 1247–1252. DOI: http://dspace.imech.ac.cn/handle/311007/43766. TIAN Z N, ZHANG L W, LI S H. Numerical simulation of pulsation features of pressure in explosion cavity and distribution of explosive energy in rock blasting [J]. Chinese Journal of Geotechnical Engineering, 2010, 32(8): 1247–1252. DOI: http://dspace.imech.ac.cn/handle/311007/43766.
[12]	郭家豪, 范锦彪. 基于振动的钻地弹爆炸能量计算方法研究 [J]. 振动与冲击, 2020, 39(10): 180–184. GUO J H, FAN J B. Numerical calculation method for earth penetrating weapon explosion energy based on vibration [J]. Journal of Vibration and Shock, 2020, 39(10): 180–184.
[13]	吴亮, 卢文波, 宗琦. 岩石中柱状装药爆炸能量分布 [J]. 岩土力学, 2006, 27(5): 734–738. WU L, LU W B, ZONG Q. Distribution of explosive energy consumed by column charge in rock [J]. Rock and Soil Mechanics, 2006, 27(5): 734–738.
[14]	肖卫国, 王肖钧, 朱号锋, 等. 不同介质地下爆炸的地震耦合效应 [J]. 爆炸与冲击, 2012, 32(3): 267–272. DOI: 10.11883/1001-1455(2012)03-0267-06. XIAO W G, WANG X J, ZHU H F, et al. Experimental study on seismic coupling effects of underground explosions in different materials [J]. Explosion and Shock waves, 2012, 32(3): 267–272. DOI: 10.11883/1001-1455(2012)03-0267-06.
[15]	LU Q, WANG Z J. Studies of the propagation of viscoelastic spherical divergent stress waves based on the generalized Maxwell model [J]. Journal of Sound and Vibration, 2016, 371: 183–195. DOI: 10.1016/j.jsv.2016.02.034.
[16]	卢强, 王占江. 标准线性固体材料中球面应力波传播特征研究 [J]. 物理学报, 2015, 64(10): 108301. DOI: 10.7498/aps.64.108301. LU Q, WANG Z J. Characteristics of spherical stress wave propagation in the standard linear solid material [J]. Acta Physica Sinica, 2015, 64(10): 108301. DOI: 10.7498/aps.64.108301.
[17]	LU Q, WANG Z J, Ding Y. Inversion for the complex elastic modulus of material from spherical wave propagation data in free field [J]. Journal of Sound and Vibration, 2019, 459: 1–18. DOI: 10.1016/j.jsv.2019.114851.
[18]	卢强, 王占江, 李进, 等. 球面波加载下黄土线黏弹性本构关系 [J]. 岩土力学, 2012, 33(11): 3292–3298. LU Q, WANG Z J, LI J, et al. Linear viscoelastic constitutive relation of loess under spherical stress wave [J]. Rock and Soil Mechanics, 2012, 33(11): 3292–3298.
[19]	JAMES K G, JOHN Q E, CARGILE J D. Cavity expansion experiments with spherical explosive charges in concrete [R]. US Army Corps of Engineers, 2009.
[20]	WANG L L, LAI H W, WANG Z J, et al. Studies on nonlinear visco-elastic spherical waves by characteristics analyses and its application [J]. International Journal of Impact Engineering, 2013, 55: 1–10. DOI: 10.1016/j.ijimpeng.2012.12.001.
[21]	MILLER S A, FLORENCE A L. Laboratory particle velocity experiments on rock from a USSR underground nuclear tests site: AD-A223108 [R]. 1990.
[22]	王占江, 李孝兰, 张若棋, 等. 固体介质中球形发散波的实验装置 [J]. 爆炸与冲击, 2000, 20(2): 103–109. WANG Z J, LI X L, ZHANG R Q, et al. An experimental apparatus for spherical wave propagation in solid [J]. Explosion and Shock Waves, 2000, 20(2): 103–109.
[23]	卢强, 王占江, 朱玉荣, 等. 基于波传播系数构建填实爆炸下花岗岩中运动及变形场 [J]. 爆炸与冲击, 2019, 39(8): 1–10. DOI: 10.11883/bzycj-2019-0140. LU Q, WANG Z J, ZHU Y R, et al. Construction of motion and deformation field in granite under tamped explosion using wave propagation coefficient [J]. Explosion and Shock Waves, 2019, 39(8): 1–10. DOI: 10.11883/bzycj-2019-0140.

施引文献

期刊类型引用(2)

1.	周岸峰，李道奎，周仕明，周旋，崔达. 爆炸地冲击载荷计算方法综述. 工程爆破. 2023(05): 38-46+56 . 百度学术
2.	曾援，李剑，马明星，贺斌，刘瑞. 基于自适应极化滤波的地下浅层P/S波解耦方法. 电子测量技术. 2023(23): 181-187 . 百度学术