爆炸冲击波作用下颅脑损伤机理的数值模拟研究

栗志杰; 由小川; 柳占立; 杜智博; 张仡; 杨策; 庄茁

doi:10.11883/bzycj-2018-0348

爆炸冲击波作用下颅脑损伤机理的数值模拟研究

doi: 10.11883/bzycj-2018-0348

栗志杰^1,,
由小川^1, ,,
柳占立¹,
杜智博¹,
张仡¹,
杨策²,
庄茁¹

1.
清华大学航天航空学院，北京 100084
2.
陆军军医大学大坪医院，重庆 400038

详细信息

作者简介:
栗志杰（1987－　），男，硕士，博士研究生，lizhijie_20082006@163.com

通讯作者:
由小川（1977－　），男，博士，副研究员，youxiaochuan@tsinghua.edu.cn

中图分类号: O383.1
计量
- 文章访问数: 5327
- HTML全文浏览量: 2057
- PDF下载量: 169
- 被引次数: 33
出版历程
- 收稿日期: 2018-09-14
- 修回日期: 2018-10-12
- 网络出版日期: 2019-10-25
- 刊出日期: 2020-01-01

Numerical simulation of the mechanism of traumatic brain injury induced by blast shock waves

1.
School of Aerospace Engineering, Tsinghua University, Beijing 100084, China
2.
Daping Hospital, Army Medical University, Chongqing 400038, China

摘要

摘要: 爆炸引起的颅脑损伤已经成为现代战场单兵的主要致伤形式，而相关的致伤机理尚未完全阐明。本文中，针对头部在爆炸冲击波作用下的动态响应及相关致伤机理进行了数值模拟研究。首先，基于颅脑的核磁共振切片建立了人体头部三维数值模型，该模型真实地反映了颅脑的生理特征与细节构造；利用该模型对人体头部碰撞实验进行数值模拟，模拟结果与实验结果吻合良好，验证了头部模型的有效性。在此基础上，基于欧拉-拉格朗日耦合(Euler-Lagrangian coupling method，CEL)方法发展了爆炸冲击波-头部流固耦合模型，对头部受到爆炸冲击波正面冲击工况进行了数值模拟，分别从流场压力分布、脑组织压力、颅骨变形与加速度等方面对头部动态响应过程进行了分析。爆炸冲击波峰值压力在流固耦合作用下增大为入射波的3.5倍，致使受到直接冲击处的颅骨与脑组织发生高频振动，相应的振动频率高达8 kHz，这与碰撞载荷下的脑组织动态响应是完全不同的。同时，该处颅骨的局部弯曲变形会沿着颅骨进行“传播”，影响着整个颅骨的变化构型，从而决定了脑组织压力与损伤的演化过程。
- 三维头部模型 /
- 爆炸冲击波 /
- 创伤性脑损伤 /
- 流固耦合
Abstract: Blast-induced traumatic brain injury (b-TBI) is a signature injury in the current military conflicts. However, the relevant mechanism of injury has not been fully elucidated. In this paper, numerical simulation study is carried out to investigate the dynamic response of brain injury mechanics during the blast loading. Firtstly, the 3D numerical head model is established based on magnetic resonance imaging (MRI) of the human head, whose physiological characteristics and detailed structures are included. The numerical model is adopted to simulate the head collision and the results are in good agreement with the experimental data, demonstrating the validity of the numerical model. Based on the coupled Eulerian-Lagrangian (CEL) theory, a fluid-solid coupling model of explosive shock wave-head is developed. The coupled model is used to simulate the situation of head subjected to frontal impacts by explosion shock wave. The dynamic response of the head is analyzed from the pressure distribution of flow field, brain pressure, skull deformation and acceleration. The peak pressure of explosion shock wave increases 3.5 times as much as that of incident wave under fluid-structure interaction, resulting in high-frequency vibration of skull and brain tissue at the site of direct shock. The corresponding vibration frequency is as high as 8 kHz, which is completely different from the dynamic response of brain tissue under head collision. At the same time, the local bending deformation will “propagate” along the skull, affecting the whole skull configuration, which determines the evolution process of brain tissue pressure and injury.
- 3D head model /
- blast shock waves /
- traumatic brain injury /
- fluid-structure interaction

HTML全文

现代局部战争战伤分析显示，由爆炸导致的单兵伤亡约占总伤亡人数的70%，其中冲击波直接作用引起的伤害约占60%^[1]。同时，国际上爆炸恐怖袭击日益猖獗，严重危害民众的生命安全。而作为人体中枢神经系统的头部，是爆炸冲击波的重要靶器官。据相关报告显示，美军中受到创伤性脑损伤（traumatic brain injury, TBI）影响的军人为总服役人数的20%，这其中由于爆炸冲击波所引起的脑损伤(blast-induced traumatic brain injury, b-TBI)占比高达40%～60%^[1]。爆炸冲击波引起颅脑损伤的途径主要包括冲击波对头部的直接作用、飞起破片对头部的冲击作用、头部与其他物体的撞击作用、冲击波对头部的电磁热等作用^[2-3]。随着爆炸伤的日渐频发，爆炸冲击波作用下的颅脑致伤机理与对应的防护策略成为研究热点。关于颅脑致伤机理，当前的研究主要集中在冲击波对头部直接作用所造成的颅脑损伤。

虽然，经过大量的动物实验、病理生理学分析以及数值模拟研究，建立了颅脑碰撞损伤与头部加速度之间的联系，并确定了相应的损伤准则^[4-6]，而相关理论并不能应用于爆炸冲击波直接作用下的颅脑损伤分析。Cernak等^[7]在生物激波管实验中，将老鼠进行固定，保证头部在冲击波条件下为稳定状态，但是实验结果表明，老鼠的脑组织仍然出现了器质性病变。费周等利用不同尺寸的微型球形炸药，对大鼠进行不同损伤程度的b-TBI实验，并从神经功能、宏观病理学、组织病理学及各种生物标志物表达水平方面，建立了可靠且可重复的b-TBI模型及损伤特征指标^[8]。该实验将大鼠放置在泡沫垫块上，并将球形炸药悬空挂于大鼠头部上面，大鼠头部在爆炸冲击波作用下是相对稳定的，但是仍然出现了脑组织损伤。这些动物实验均表明，爆炸冲击波对颅脑的损伤可以通过直接作用产生。而爆炸冲击波对头部直接作用的时间通常在1～3 ms，是应力波主导的物理现象。因此，基于头部加速度建立的颅脑碰撞损伤理论无法解释上述实验现象。

针对爆炸冲击波直接作用下的颅脑损伤，学者们做了大量的研究工作，并提出了多种致伤机理。Cloots等^[9-11]将微观体积单元（critical volume element, CVE）与常规有限元单元进行结合，建立了脑组织的多尺度模型，从而得到脑组织细胞层次的力学响应，并结合生物细胞的力学损伤实验，建立了相应的损伤准则。Radovitzky等^[12]、Jean等^[13]发展了脑组织的多尺度生物力学模型，研究了爆炸冲击波作用下颅脑在组织、细胞等不同尺度的生物响应特征，并建立了评估b-TBI风险的动物与人之间的尺度律。Goeller等^[14]、Salzar等^[15]提出脑脊液空穴致伤机制，通过激波管实验与数值模拟，对脑脊液的空穴效应进行了研究，研究结果表明，脑脊液在冲击波作用下会发生空化与湮灭，形成局部点高压，从而造成脑组织损伤。Franck^[16]基于类脑软材料，研究了微型空穴效应在b-TBI当中的作用，并利用实时成像技术对不同应变率下神经元细胞的损伤模式进行了探究。Bhattacharjee^[17]提出主动脉冲压致伤理论，冲击波通过压缩主动脉血液，从而使其形成高压液体，并经毛细血管对脑组织形成冲击，造成脑损伤。同时，Courtney等^[18]设计了动物实验，验证了冲击波可以通过作用躯干部位间接对脑组织造成损伤。Moss等^[19]通过数值模拟研究发现，颅骨在爆炸冲击波作用下会发生局部弯曲，从而使脑组织压力出现较大的压力梯度，颅骨的弯曲变形与颅内压力的分布密切相关。但是，颅骨局部弯曲变形与脑组织压力之间的内在联系有待于进一步确立，相应的颅脑损伤机理仍需进一步研究。

基于上述研究现状，本文中针对颅骨局部弯曲变形这个致伤机制，开展进一步的深入研究。首先，建立充分反映颅脑生理结构的三维头部数值模型，并利用头部碰撞试验数据对头部模型的有效性进行验证。然后，基于欧拉-拉格朗日耦合法（Euler-Lagrangian coupling method，CEL）建立爆炸冲击波-头部流固耦合模型，对头部在正面冲击工况下的动态响应过程进行数值模拟，确立颅骨局部弯曲变形与脑组织压力之间的内在联系，初步揭示脑组织的损伤机理。

1. 三维头部数值模型与验证

1.1 三维头部有限元模型

在颅脑结构中，脑组织包括大脑、小脑与脑干3部分，其表面有3层被膜包裹，由内向外依次为软脑膜、蛛网膜和硬脑膜。软脑膜薄而富含血管，紧贴脑表面并深入脑的沟裂中。蛛网膜是半透明结缔组织薄膜，它与软脑膜之间是蛛网膜下腔，内含脑脊液（cerebral spinal fluid, CSF）、蛛网膜小梁和较大血管。硬脑膜呈套状包被脑，在颅腔内分化为大脑镰与小脑幕，两者将颅腔分为左、右与下3个腔体，限制脑组织相对于颅骨的过大滑移。颅骨与硬脑膜外侧紧密结合，其内、外为刚度与密度较大的骨密质，中间为骨松质，为典型的类三明治结构。图1显示了人体头部的解剖结构。

图 1 头部解剖结构示意图^[20]

Figure 1. Schematic diagram for anatomical structures of head ^[20]

下载: 全尺寸图片幻灯片

基于颅脑的核磁共振切片，建立了三维头部有限元模型，该模型包括脑组织（大脑、小脑与脑干），膜结构（软脑膜、蛛网膜、硬脑膜、大脑镰与小脑幕），包围脑组织的脑脊液以及颅骨、皮肤与颈部。该头部模型共有303 588个单元与530 954个节点，单元尺寸为2～3 mm，充分考虑了人体头部的生理结构与相关细节，如图2所示。

图 2 三维头部有限元模型

Figure 2. 3D finite element model of human head

下载: 全尺寸图片幻灯片

1.2 本构模型与材料参数

对于大脑、小脑与脑干等脑组织，采用黏超弹性本构模型进行表征^[21]。利用Mooney-Rivlin模型表征其超弹性行为，对应的应变能函数定义为：

$U = {C_{10}}({J_1} - 3) + {C_{01}}({J_2} - 3) + \frac{1}{{{D_1}}}{({J_{{\rm{el}}}} - 1)^2}$

(1)

式中： $U$ 为应变势能， ${C_{10}}$ 、 ${C_{01}}$ 与 ${D_1}$ 为温度相关的材料参数（ $2/{D_1}$ 为体积模量 $K$ ）， ${J_{{\rm{el}}}}$ 为弹性体积比， ${J_1}$ 和 ${J_2}$ 分别为第一与第二偏应变不变量。

对于脑组织的黏弹部分，第二Kirchhoff应力 ${S_{ij}}$ 通过卷积形式可表示为：

${S_{ij}}{\rm{ = }}\int_0^t {{G_{ijkl}}(t - \tau )} \frac{{\partial {E_{kl}}}}{{\partial \tau }}{\rm d}\tau$

(2)

$G(t) = {G_0}\left(1 - \sum\limits_{i = 1}^N {{g_i}(1 - {{\rm e}^{ - t/{\tau _i}}})}\right )$

(3)

${G_\infty } = {G_0}\left(1 - \sum\limits_{i = 1}^N {{g_i}}\right){\rm{ = }}{G_0}{\rm{ - }}\sum\limits_{i = 1}^N {{G_i}}$

(4)

式中： ${E_{kl}}$ 为Green应变张量分量， $G(t)$ 为t时刻的剪切松弛模量，采用Prony级数进行展开。 ${G_0}$ 、 ${G_\infty }$ 分别为瞬态与稳态剪切模量， ${g_i}$ 和 ${\tau _i}$ 分别为无量纲剪切松弛模量参数与松弛时间， ${G_i}$ 为对应于 ${\tau _i}$ 的剪切松弛模量。对应的Cauchy应力 ${\sigma _{ij}}$ 可表示为：

${\sigma _{ij}} = JF_{ik}^{\rm T} {S_{km}} {F_{mj}}$

(5)

式中：F为变形梯度矢量，J为变换矩阵的雅克比行列式。由于脑组织接近不可压缩，不考虑体积变形相关的黏弹性特性，相应的无量纲体积模量松弛参数 ${k_i}$ 均设置为0。关于上述模型中的具体材料参数，参照文献[21-26]，采用文献[21]中的相关参数，见表1。

表 1 黏超弹性材料参数

Table 1. Material properties for hyper-viscoelastic model

密度/(kg·m⁻³)	C₀₁/Pa	C₁₀/Pa	K/GPa
1 060	3 102.5	3 447.2	2.19

阶数i	g_i	k_i	τ_i/s
1	0.528 262	0	0.008
2	0.301 899	0	0.150

下载: 导出CSV

| 显示表格

脑脊液与水的性质较为接近，在蛛网膜下腔内流动；该腔是由蛛网膜小梁支撑的腔体，能够传递剪切变形，如图1所示。为了便于建模，头部模型中的脑脊液是物理脑脊液、蛛网膜小梁与较大血管的综合体。对于脑脊液CSF的本构模型，文献[22-23]中直接采用剪切模量相对很小的线弹性模型来进行表征。为了真实反映脑脊液与蛛网膜小梁的力学特性，模型中采用Mie-Grüneisen状态方程表征脑脊液的不可压缩特性，同时采用较小的剪切模量（22 kPa）表征蛛网膜小梁的剪切传递特性^[24-26]。Mie-Grüneisen状态方程表示为：

$p = \frac{{{\rho _0}c_0^2(1 - {\rho _0}/\rho )}}{{{{(1 - s(1 - {\rho _0}/\rho ))}^2}}},\quad\quad{u_{\rm s}} = {c_0} + s {u_{\rm p}}$

(6)

式中：p为压力， ${c_0}$ 为声速，取为1 450 m/s； ${\rho _0}$ 和 $\rho$ 分别为初始密度与当前密度，其中 ${\rho _0}$ 取为1 040 kg/m³； $s$ 为无量纲材料参数，建立应力波速度 ${u_{\rm s}}$ 与粒子速度 ${u_{\rm p}}$ 之间的关系，取为1.921^[25]。

头部模型中的软脑膜、蛛网膜、硬脑膜、大脑镰、小脑幕、头骨、颈部与皮肤，均采用线弹性本构模型进行表征，并参考相关文献[23-26]，确定相应的材料参数，见表2。

表 2 线弹性材料参数

Table 2. Material properties for liner elastic model

部位	密度/(kg·m⁻³)	弹性模量/MPa	泊松比	部位	密度/(kg·m⁻³)	弹性模量/MPa	泊松比
大脑镰	1 130	31.5	0.45	骨密质	2 000	15 000	0.22
硬脑膜	1 130	31.5	0.45	骨松质	1 300	1 000	0.24
小脑幕	1 130	31.5	0.45	下颌骨	1 710	5 370	0.19
蛛网膜	1 130	22.0	0.45	脖颈	1 250	354	0.30
软脑膜	1 130	11.5	0.45	皮肤	1 200	16.7	0.42

下载: 导出CSV

| 显示表格

1.3 接触关系与边界条件

基于颅脑的生理结构并参考相关文献[23-28]，三维头部数值模型中的各部分组织之间的接触关系定义为绑定连结。同时，采用状态方程的脑脊液实现了脑组织与颅骨之间的有限滑动，符合颅脑的生理结构。而对于头部模型中颈部的边界条件，至今没有明确统一的结论。Ganpule等^[23]直接采用固定边界条件，对头部碰撞实验^[29]进行数值模拟，两者的结果具有相对较好的吻合性。Chen等^[27]对比分析了自由边界条件与固定边界条件两种工况下的模拟结果，虽然均与实验结果存在一定的差异，但是自由边界对应的模拟结果更能真实反映头部的动态响应过程。躯体对颅脑动态响应的影响是与颈部的肌肉状态密切相关的。当肌肉处于紧绷状态，颈部较为“刚硬”，躯体对颅脑动态响应的影响较大；而当肌肉处于放松状态时，颈部则较为柔软，这种影响作用相对较小。通常情况下，颈部肌肉是处于放松状态的。本文中基于显式有限元程序，对包含躯体的完整模型与自由边界条件下的头部模型在撞击载荷作用下的动态响应过程进行了数值模拟。其中，颈部采用较小的弹性模量（3.54 MPa），表征颈部肌肉的放松状态。图3对比显示了两种数值模型对应的前额处与枕部处脑组织的压力曲线。由图3可知，采用自由边界条件的头部模型，其模拟结果能够较好地吻合完整模型的结果。这表明：当颈部处于放松状态时，躯体对颅脑动态响应过程的影响在碰撞等冲击载荷条件下（通常为几个毫秒）基本上是可以忽略的，此时的颅脑动态响应过程呈现出明显的波动效应。同时，两种数值模型对应的前额处脑组织压力在t=7 ms时开始产生差异，并且呈现出不断扩大的趋势。在这个过程中，颅脑的动态响应会由应力波主导转化为整体结构响应，受到躯体对颈部约束作用的影响。

图 3 前额处与枕部处脑组织的压力曲线

Figure 3. Pressure curues of brain tissue at the forehead and occiput

下载: 全尺寸图片幻灯片

1.4 头部模型验证与数值模拟

Nahum等^[29]进行的头部碰撞实验，已成为头部数值模型有效性验证的经典实验。该实验利用运动块体撞击前额部位，利用传感器记录特征部位（前额、颅顶、枕部与颅后窝处）的颅内压力变化过程。由于该头部碰撞实验没有提供与头部发生碰撞块体的材料力学参数，参考Ganpule等^[23]进行头部数值模型验证时所使用的方法：将实验中测量所得到的碰撞接触力转化为平均压力作用于前额部位，作用面积为1 470 mm²。对应的压力曲线如图4所示。基于上述头部碰撞模型，可获取颅内压力的动态响应过程，与头部碰撞实验数据进行对比分析。为减小采样位置对模拟结果与实验数据对比分析的影响，在前额、颅顶、枕部与颅后窝处各取3个特征点（P1、P2与P3），特征点位置如图5所示。

图 4 碰撞压力曲线

Figure 4. Impact pressure curve

下载: 全尺寸图片幻灯片

图 5 颅内压力测点位置

Figure 5. Locations of measurement points for intracranial pressure

下载: 全尺寸图片幻灯片

本文中提取特征点处颅内压力的模拟结果，与头部碰撞实验数据进行对比，两者具有较好的吻合性，如图6所示。由碰撞实验结果可知，前额与颅顶处的颅内压力在整个碰撞过程中基本上均为正压，表明这两处脑组织处于受压缩状态。同时，前额处的颅内压力要明显高于颅顶处的相应值，这与头部简化模型的理论解吻合^[30]。相比而言，枕部与颅后窝处的颅内压力在碰撞前期阶段为负压，表明脑组织处于拉伸状态，相应的压力幅值较为接近。但是，这两处的颅内压力在后期阶段会发生“反转”，由负压转变为正压，正压幅值相比于负压幅值较小。基于三维头部数值模型所得到的模拟结果，能够较好地预测颅内压力在碰撞条件下的演化过程，与实验数据吻合较好，从而验证了该头部数值模型的有效性。

图 6 颅内压力模拟结果与实验结果的对比

Figure 6. Comparison of intracranial pressure between experiments and calculations

下载: 全尺寸图片幻灯片

为了对比分析头部动态响应在碰撞载荷与爆炸冲击波载荷作用下的差异，提取了头颈部在特征时刻的压力云图，如图7所示。由于颅骨与脑组织的惯性力差异，在碰撞载荷作用下，颅骨加速度明显大于脑组织加速度。因此，前额处（撞击处）颅骨与脑组织之间产生相对径向压缩位移，使得脑组织处于压缩状态，对应的压力为正；而枕部处（对撞处）两者之间产生相对径向拉伸位移，脑组织处于拉伸状态，对应的压力为负。随着撞击载荷的持续作用，颅内正、负压的作用范围持续扩大，发展成为靠近前额的前半部分脑组织受压而靠近枕部的后半部分脑组织受拉的受力状态。随着外部载荷的减弱和头颈部之间相对变形的增大，头部运动会在某一特定时刻发生转向，颅骨与脑组织的惯性差异使得两者之间的相对位移发生转变，致使脑组织压力的状态发生“反转”。但是，“反转”后的脑组织压力，其幅值则相对较小。因此，颅骨与脑组织间的相对径向位移决定脑组织压力，撞击处前半部分脑组织主要处于压缩状态，而对撞处后半部分脑组织则主要处于拉伸状态。

图 7 前额撞击时脑组织压力云图

Figure 7. Nephogram of brain pressure for the forehead collision

下载: 全尺寸图片幻灯片

2. 爆炸冲击波-头部流固耦合模型

为了获得头部在爆炸冲击波下的动态响应规律，需要建立爆炸冲击波-头部流固耦合模型。头部在常规爆炸冲击波作用下的变形相对较小，可采用拉格朗日单元进行网格划分。而空气需要模拟爆炸冲击波的产生、传递以及与头部相互作用的过程，对应区域应当采用欧拉网格，避免单元过度变形。爆炸冲击波与头部在接触边界上的相互作用，需要通过求解流固耦合方程进行确定。通用有限元软件ABAQUS通过非线性瞬态程序与欧拉-拉格朗日耦合法（Euler-Lagrangian coupling method，CEL），对整体系统的三大守恒偏微分方程（动量、质量与能量）同时进行求解。欧拉区域能够对导致严重网格扭曲的高度动态事件进行模拟（如爆炸冲击波），并为拉格朗日区域提供了压力边界条件；拉格朗日区域可完全或部分位于欧拉区域内，拉格朗日表面为欧拉区域提供边界条件，该边界条件不允许在其表面的法线方向上有流动产生。同时，增强浸入边界方法提供欧拉区域和拉格朗日区域之间的耦合作用。因此，固定欧拉网格与基于增强浸入边界法的流固界面模型，能够模拟爆炸冲击波与头部之间耦合作用，获得头部的动态响应规律。

2.1 有限元离散

模拟爆炸冲击波的欧拉区域采用正方体形状，边长为600 mm，这个尺寸的选择充分考虑了冲击波与头部的作用时间，防止冲击波边界反射对头部产生非真实的二次冲击作用。经过验证的三维头部数值模型，在爆炸冲击波-头部流固耦合模型中可以直接使用，并放置在欧拉区域的中心位置。欧拉区域使用缩减积分六面体欧拉单元（EC3D8R）进行网格划分，并采用分块划分单元网格的方法：与头部发生接触的区域采用较小的单元尺寸，利于模拟爆炸冲击波与头部之间的流固耦合作用；没有与头部发生作用的区域则采用相对较大的单元，利于节约计算成本。同时，针对单元尺寸对模拟结果的影响进行参数化研究。研究表明：当中心区域的单元尺寸为3 mm时，相应的模拟结果能够收敛，与文献[23]结论基本一致。因此，将单元尺寸定义为3 mm，欧拉区域共有8 524 800个单元。欧拉区域采用理想气体状态方程来表征空气，空气密度为1.18 kg/m³，气体常数为287.04 J/(kg·K)，比热容为1 000 J/(kg·K)。最终，建立的爆炸冲击波-头部流固耦合模型如图8所示。

图 8 爆炸冲击波-头部的流固耦合模型

Figure 8. Fluid-solid coupling model of explosive blast wave-head

下载: 全尺寸图片幻灯片

2.2 荷载与边界条件

根据我国工程设计规范中的设计准则，爆炸冲击波的超压值（以大气压力为参考零值）会随时间呈现指数型衰减。其理论表达式为：

$\Delta p(t) = \Delta {p_ + }\left(1 - \frac{t}{{{\tau _ + }}}\right)\exp\left ( -\left (\frac{1}{2} + 10\Delta {p_ + }\left(1.1 - (0.13 + 2\Delta {p_ + })\frac{t}{{{\tau _ + }}}\right)\right)\frac{t}{{{\tau _ + }}}\right)$

(7)

式中： $\Delta {p_ + }$ 为爆炸冲击波的峰值压力，单位为MPa； ${\tau _ + }$ 为冲击波的正压持续时间，时间为s。超压峰值采用上述规范的理论公式，而正压时间则采用萨多夫斯基的理论公式。两个公式分别表示如下：

$\Delta {p_ + } = \frac{{0.082}}{{\overline r }} + \frac{{0.265}}{{{{\overline r }^2}}} + \frac{{0.686}}{{{{\overline r }^3}}}\quad\quad\overline r = r/\sqrt[3]{w}$

(8)

${\tau _ + }{\rm{ = }}B{r^{1/2}}{w^{1/6}}\quad\quad B = (1.0 {\text{～}} 1.5) \times {10^{ - 3}}$

(9)

式中： $\bar r$ 为对比距离，单位为m/kg^1/3；r为距离爆心的距离，单位为m；w为炸药当量，单位为kg；B为爆炸环境参数，这里取为1.2×10⁻³，单位为s/(m^1/2· kg^1/6)。

基于上述理论公式，可获得炸药当量为7 kg且距离爆心3.8 m处的压力曲线。将理论结果与实爆实验结果进行对比，两者具有较好的吻合性，如图9所示。

图 9 爆炸冲击波压力曲线的理论结果与实验结果对比

Figure 9. Comparison of pressure curves of explosion blast between theoretical and experimental results

下载: 全尺寸图片幻灯片

本文中，采用这个典型爆炸冲击波压力曲线作为欧拉区域的输入载荷，直接作用在入口处。欧拉场的出口处采用无反射边界条件，冲击波可以直接从这里逸出，不再对欧拉场产生影响。而欧拉场的环向4个边界则采用法向速度约束边界条件，防止冲击波从侧面逸出，模拟平面型爆炸冲击波与头部的相互作用，可真实反映现场环境。

2.3 接触与初始条件

由爆炸冲击波-头部流固耦合模型的理论基础可知，两者之间的相互作用需要借助基于增强浸入边界法的流固界面模型来实现。在程序中，通用接触（general contact）采用增强浸入边界法能够较好地处理拉格朗日单元与欧拉单元之间的接触问题。因此，通用接触是进行流固耦合分析的必备条件。

对于初始条件，系统中有两个方面需要进行合理设置。（1）欧拉区域中空气的初始温度与参考压力。这里设定初始时刻的空气温度为303 K，以标准大气压力作为参考值，即欧拉区域的压力初始值为零。（2）空气材料的初始分布。这里可以使用体积分数工具（volume fraction tool）来设定空气在离散欧拉场中的分布，保证每个位置点处只有一种材料。

3. 结果与分析

基于上述爆炸冲击波-头部流固耦合模型，可以获得头部在爆炸冲击波作用下的动态响应规律，进而揭示颅脑的损伤机理。

3.1 流场分析

爆炸冲击波从正面作用于头部时，流场压力的整个演化过程如图10所示。在t=0.375 ms时，平面冲击波与头部即将发生接触，此时的冲击波压力峰值为170 kPa。随后，冲击波与头部发生相互作用，反射叠加效应使其压力峰值急剧上升；当t=0.450 ms时，前额处的压力峰值达到600 kPa，是入射冲击波峰值的3.5倍。在t=0.720 ms时，冲击波经过前额进入颅顶，曲面型的颅顶结构使其在该处发生边界层分离现象。因此，冲击波对颅顶处的压力峰值相对较小。当t=0.960 ms时，由于头部阻挡作用而分散的冲击波会在后脑部位重新汇集，该处的压力值迅速上升，从而对头部造成二次冲击。相比于面部首次冲击的压力峰值，二次冲击的对应峰值出现了显著下降。

图 10 正面冲击时流场压力分布

Figure 10. Pressure distribution of flow field in blast frontal impact

下载: 全尺寸图片幻灯片

图11显示了t=0.480 ms时眼部周边处的流场压力分布。由图11可见，眼部处的流场压力峰值达到757 kPa，是入射平面波峰值（170 kPa）的4.5倍。眼部的凹陷结构使其容易成为靶器官，是爆炸冲击波防护的重点。

图 11 正面冲击时眼部的流场压力分布

Figure 11. Pressure distribution of flow field around the eyes in blast frontal impact

下载: 全尺寸图片幻灯片

3.2 脑组织压力

对于正面冲击，先提取脑组织特征点处的压力时程曲线，如图12所示。由图12可知，各特征点处的脑组织压力均呈现出明显的波动特性，并且随着时间的推移，这种波动特性会逐步弱化。同时，从前额至颅顶处，特征点1～4的脑组织压力在前期阶段呈现出明显的周期特性，其中前额处的频率高达8 kHz；而特征点5～8则没有明显的周期特性。在爆炸冲击波下，脑组织压力的高频波动特性与碰撞载荷下的脑组织响应是完全不同的，这就决定了针对两种不同工况需要采用不同的颅脑防护策略。结合脑组织压力的动态演化过程与损伤判断准则可知，多处脑组织的压力将超过重度损伤阈值（235 kPa）^[31]，造成脑组织出现多处重度损伤。

图 12 正面冲击时特征点处脑组织的压力曲线

Figure 12. Pressure curves of brain tissue at feature points in blast frontal impact

下载: 全尺寸图片幻灯片

基于模拟结果，可提取颅脑在特征时间点处的压力分布云图和颅骨变形云图，如图13所示。在颅骨变形图中，带有网格的颅骨为初始形态，没有网格的颅骨为当前形态，采用的放大系数为500。由压力云图可以看出，前额处脑组织在冲击波作用下起初处于受压状态（t=0.465 ms），而后转化为受拉状态（t=0.520 ms）；在后续过程中，脑组织压力会在压、拉两种状态之间进行转化。紧邻颅骨内表面处的脑组织，其压力沿着颅骨内表面出现正、负压交替现象，形成较大的压力梯度；而远离颅骨内表面的脑组织，其压力梯度则相对较小。同时，具有典型三明治结构的头部颅骨，利用材料阻抗不匹配特性截断并分散冲击处的能量，从而为脑组织提供防护。

图 13 正面冲击时颅脑压力云图与颅骨变形云图

Figure 13. Nephogram of brain pressure and skull displacement in blast frontal impact

下载: 全尺寸图片幻灯片

从流场压力的演化过程可知，爆炸冲击波主要作用于前额等面部结构，而颅顶等部位由于存在边界层脱离效应，受到冲击波的作用相对较小，在进行颅骨变形分析时可以忽略这个影响。当t=0.465 ms时，流固耦合作用使得爆炸冲击波达到最大值，冲击波直接挤压皮肤与颅骨，前额处颅骨的“平面”变形造成该处脑组织处于挤压状态，体积压力为正压。随着前额颅骨变形的进一步增大，当t=0.520 ms时，该处颅骨的“平面”变形转变为“内凹”弯曲，致使该处与颅骨内表面紧密连结的脑脊液与脑组织处于受拉状态，体积压力为负压。而紧邻前额处的颅骨则相当于发生了“外凸”弯曲，对该处的脑脊液与脑组织进行挤压，体积压力为正压；变形的连续性使得颅顶处颅骨产生较小的“内凹”弯曲，体积压力为负压。当t=0.585 ms时，冲击波作用的减弱与颅骨的弹性边界使得前额处颅骨的“内凹”弯曲在达到最大值后发生回弹，由“内凹”弯曲转变为“外凸”弯曲，并主导着颅骨其他部位弯曲类型的改变。对比t=0.520 ms与t=0.585 ms时的云图可知，各位置处的脑组织的拉压状态发生转变。对比分析后续特征时刻的脑组织压力与颅骨弯曲变形，可以得到相同的结论：理想状态下，颅骨发生“内凹”弯曲，该处脑组织处于受拉状态；颅骨发生“外凸”弯曲，脑组织则为受压状态。而真实条件下，颅骨内表面处脑组织在特定位置处的受力状态同时会受到周边颅骨变形的影响。

本文中进一步提取了前额处颅骨的加速度时程曲线，并将该处特征时刻的加速度标注在该曲线上，如图14所示。其中实心圆点表示该处脑组织压力为正时的加速度，空心圆点则表示压力为负时的加速度。由图14可见，爆炸冲击波致使前额颅骨发生局部高频振动，其频率与该处脑组织压力的波动频率相一致。前额颅骨最值加速度与颅骨弯曲变形模式、脑组织状态均具有较好的相关性：颅骨受冲击方向加速度最大时，颅骨为“外凸”弯曲模式，脑组织处于受压状态；而颅骨受冲击方向加速度最小（冲击反方向最大）时，颅骨为“内凹”弯曲模式，脑组织处于受拉状态。

图 14 前额处颅骨加速度曲线

Figure 14. Acceleration curves of forehead skull

下载: 全尺寸图片幻灯片

对比分析脑组织压力、颅骨局部弯曲变形与加速度的整个演化过程可知，前额处颅骨在爆炸冲击波作用下发生弯曲变形并引起局部振动，该处的局部变形会沿着颅骨进行“传播”，影响着颅骨的变化构型，从而使脑组织在外载作用下由前额处的“单一波源”逐步发展到整个颅骨的“多处波源”，这些波源共同决定脑组织压力的演化过程。颅骨局部弯曲致使脑组织沿着颅骨内表面出现拉、压交替状态，形成较大的压力梯度。同时，颅骨局部弯曲的振动频率较高，使得脑组织压力具有高频波动特性。

4. 结　论

建立了具有详细解剖学结构和较高生物仿真度的三维头部数值模型，并基于Nahum等^[29]的碰撞实验验证了该头部模型的有效性。同时，基于CEL方法建立了爆炸冲击波-头部流固耦合模型，对头部受到爆炸冲击波正面冲击工况进行了数值模拟，分析了脑组织的动态响应规律，揭示了相应的损伤机理。主要结论如下。

（1）在碰撞载荷下，颅骨与脑组织间的相对径向位移决定脑组织压力，撞击处前半部分脑组织主要处于压缩状态，而对撞处后半部分脑组织则主要处于拉伸状态。

（2）流固耦合作用使得爆炸冲击波峰值压力增大为入射波的3.5倍，具有生理凹陷结构的眼部则使该值增大到4.5倍，是冲击波的重要靶器官。

（3）爆炸冲击波在颅顶处会发生边界层分离效应，该处的压力峰值相对较小；冲击波从头部撞击处周边发生绕流，并在撞击的对侧处汇合，对颅脑造成相对较小的二次冲击。

（4）正面作用时前额处颅骨的振动频率高达8 kHz，从而使脑组织的压力具有高频波动特性。这与碰撞载荷下的脑组织动态响应完全不同。

（5）前额处颅骨（正面冲击）的局部弯曲变形，会沿着颅骨进行“传播”，影响着颅骨的变化构型，从而决定了脑组织压力的演化过程。同时，颅骨最值加速度与颅骨弯曲变形模式、脑组织状态均具有较好的相关性。

图 1 头部解剖结构示意图^[20]

Figure 1. Schematic diagram for anatomical structures of head ^[20]

下载: 全尺寸图片幻灯片

图 2 三维头部有限元模型

Figure 2. 3D finite element model of human head

下载: 全尺寸图片幻灯片

图 3 前额处与枕部处脑组织的压力曲线

Figure 3. Pressure curues of brain tissue at the forehead and occiput

下载: 全尺寸图片幻灯片

图 4 碰撞压力曲线

Figure 4. Impact pressure curve

下载: 全尺寸图片幻灯片

图 5 颅内压力测点位置

Figure 5. Locations of measurement points for intracranial pressure

下载: 全尺寸图片幻灯片

图 6 颅内压力模拟结果与实验结果的对比

Figure 6. Comparison of intracranial pressure between experiments and calculations

下载: 全尺寸图片幻灯片

图 7 前额撞击时脑组织压力云图

Figure 7. Nephogram of brain pressure for the forehead collision

下载: 全尺寸图片幻灯片

图 8 爆炸冲击波-头部的流固耦合模型

Figure 8. Fluid-solid coupling model of explosive blast wave-head

下载: 全尺寸图片幻灯片

图 9 爆炸冲击波压力曲线的理论结果与实验结果对比

Figure 9. Comparison of pressure curves of explosion blast between theoretical and experimental results

下载: 全尺寸图片幻灯片

图 10 正面冲击时流场压力分布

Figure 10. Pressure distribution of flow field in blast frontal impact

下载: 全尺寸图片幻灯片

图 11 正面冲击时眼部的流场压力分布

Figure 11. Pressure distribution of flow field around the eyes in blast frontal impact

下载: 全尺寸图片幻灯片

图 12 正面冲击时特征点处脑组织的压力曲线

Figure 12. Pressure curves of brain tissue at feature points in blast frontal impact

下载: 全尺寸图片幻灯片

图 13 正面冲击时颅脑压力云图与颅骨变形云图

Figure 13. Nephogram of brain pressure and skull displacement in blast frontal impact

下载: 全尺寸图片幻灯片

图 14 前额处颅骨加速度曲线

Figure 14. Acceleration curves of forehead skull

下载: 全尺寸图片幻灯片

表 1 黏超弹性材料参数

Table 1. Material properties for hyper-viscoelastic model

密度/(kg·m⁻³)	C₀₁/Pa	C₁₀/Pa	K/GPa
1 060	3 102.5	3 447.2	2.19

阶数i	g_i	k_i	τ_i/s
1	0.528 262	0	0.008
2	0.301 899	0	0.150

下载: 导出CSV

表 2 线弹性材料参数

Table 2. Material properties for liner elastic model

部位	密度/(kg·m⁻³)	弹性模量/MPa	泊松比	部位	密度/(kg·m⁻³)	弹性模量/MPa	泊松比
大脑镰	1 130	31.5	0.45	骨密质	2 000	15 000	0.22
硬脑膜	1 130	31.5	0.45	骨松质	1 300	1 000	0.24
小脑幕	1 130	31.5	0.45	下颌骨	1 710	5 370	0.19
蛛网膜	1 130	22.0	0.45	脖颈	1 250	354	0.30
软脑膜	1 130	11.5	0.45	皮肤	1 200	16.7	0.42

下载: 导出CSV

参考文献(31)

[1]	TANIELIAN T, JAYCOX L H. Invisible wounds of war [R]. Santa Monica, CA: RAND Corporation, 2008. DOI: 10.1037/e527802010-001.
[2]	DEPALMA R G, BURRIS D G, CHAMPION H R, et al. Blast injuries [J]. New England Journal of Medicine, 2005, 352(13): 1335–1342. DOI: 10.1056/NEJMra042083.
[3]	MOORE D F, RADOVITZKY R A, SHUPENKO L, et al. Blast physics and central nervous system injury [J]. Future Neurology, 2008, 3(3): 243–250. DOI: 10.2217/14796708.3.3.243.
[4]	VERSACE J. A review of severity index [C] // Proceedings of the 15th Stapp Car Crash Conference. San Diego: Society of Automotive Engineers, 1971: 771−796. DOI: 10.4271/710881.
[5]	NEWMAN J A. A generalized acceleration model for brain injury threshold (GAMBIT) [C] // Proceedings of International IRCOBI Conference. 1986.
[6]	NEWMAN J A, SHEWCHENKO N. A proposed new biomechanical head injury assessment function: the maximum power index [R]. SAE Technical Paper, 2000.
[7]	CERNAK I, WANG Z G, JIANG J X, et al. Ultrastructural and functional characteristics of blast injury-induced neurotrauma [J]. The Journal of Trauma: Injury, Infection, and Critical Care, 2001, 50(4): 695–706. DOI: 10.1097/00005373-200104000-00017.
[8]	LIU M D, ZHANG C, LIU W B, et al. A novel rat model of blast-induced traumatic brain injury simulating different damage degree: implications for morphological, neurological, and biomarker changes [J]. Frontiers in Cellular Neuroscience, 2015, 9: 168. DOI: 10.3389/fncel.2015.00168.
[9]	CLOOTS R J H, VAN DOMMELEN J A W, KLEIVEN S, et al. Traumatic brain injury at multiple length scales: relating diffuse axonal injury to discrete axonal impairment [C] // IRCOBI conference. 2010.
[10]	CLOOTS R J H, VAN DOMMELEN J A W, GEERS M G D. A tissue-level anisotropic criterion for brain injury based on microstructural axonal deformation [J]. Journal of the Mechanical Behavior of Biomedical Materials, 2012, 5(1): 41–52. DOI: 10.1016/j.jmbbm.2011.09.012.
[11]	CLOOTS R J H, VAN DOMMELEN J A W, KLEIVEN S, et al. Multi-scale mechanics of traumatic brain injury: predicting axonal strains from head loads [J]. Biomechanics and Modeling in Mechanobiology, 2013, 12(1): 137–150. DOI: 10.1007/s10237-012-0387-6.
[12]	RADOVITZKY R, SOCRATE S, TABER K, et al. Investigations of tissue-level mechanisms of primary blast injury through modeling, simulation, neuroimaging and neuropathological studies [R]. Massachusetts Institute of Technology Cambridge, 2012. DOI: 10.21236/ada573887.
[13]	JEAN A, NYEIN M K, ZHENG J Q, et al. An animal-to-human scaling law for blast-induced traumatic brain injury risk assessment [J]. Proceedings of the National Academy of Sciences of the United States of America, 2014, 111(43): 15310–15315. DOI: 10.1073/pnas.1415743111.
[14]	GOELLER J, WARDLAW A, TREICHLER D, et al. Investigation of cavitation as a possible damage mechanism in blast-induced traumatic brain injury [J]. Journal of Neurotrauma, 2012, 29(10): 1970–1981. DOI: 10.1089/neu.2011.2224.
[15]	SALZAR R S, TREICHLER D, WARDLAW A, et al. Experimental investigation of cavitation as a possible damage mechanism in blast-induced traumatic brain injury in post-mortem human subject heads [J]. Journal of Neurotrauma, 2017, 34(8): 1589–1602. DOI: 10.1089/neu.2016.4600.
[16]	FRANCK C. Microcavitation: the key to modeling blast traumatic brain injury? [J]. Concussion, 2017, 2(3): CNC47. DOI: 10.2217/cnc-2017-0011.
[17]	BHATTACHARJEE Y. Shell shock revisited: solving the puzzle of blast trauma [J]. Science, 2008, 319: 406–408. DOI: 10.1126/science.319.5862.406.
[18]	COURTNEY A C, COURTNEY M W. A thoracic mechanism of mild traumatic brain injury due to blast pressure waves [J]. Medical Hypotheses, 2009, 72(1): 76–83. DOI: 10.1016/j.mehy.2008.08.015.
[19]	MOSS W C, KING M J, BLACKMAN E G. Skull flexure from blast waves: a mechanism for brain injury with implications for helmet design [J]. Physical Review Letters, 2009, 103(10): 108702. DOI: 10.1103/PhysRevLett.103.108702.
[20]	FELTEN D L, O’BANION M K, MAIDA M S. Netter’s atlas of neuroscience [M]. Elsevier Health Sciences, 2015: 49.
[21]	CHAFI M S, KARAMI G, ZIEJEWSKI M. Biomechanical assessment of brain dynamic responses due to blast pressure waves [J]. Annals of Biomedical Engineering, 2010, 38(2): 490–504. DOI: 10.1007/s10439-009-9813-z.
[22]	KLEIVEN S. Predictors for traumatic brain injuries evaluated through accident reconstructions [J]. Stapp Car Crash Journal, 2007, 51: 81–114. DOI: 10.12783/dtcse/wcne2017/19838.
[23]	GANPULE S, ALAI A, PLOUGONVEN E, et al. Mechanics of blast loading on the head models in the study of traumatic brain injury using experimental and computational approaches [J]. Biomechanics and Modeling in Mechanobiology, 2013, 12(3): 511–531. DOI: 10.1007/s10237-012-0421-8.
[24]	CHAFI M S, GANPULE S, GU L X, et al. Dynamic response of brain subjected to blast loadings: influence of frequency ranges [J]. International Journal of Applied Mechanics, 2011, 3(4): 803–823. DOI: 10.1142/S175882511100124X.
[25]	WANG C, PAHK J B, BALABAN C D, et al. Biomechanical assessment of the bridging vein rupture of blast-induced traumatic brain injury using the finite element human head model [C] // ASME 2012 International Mechanical Engineering Congress and Exposition. American Society of Mechanical Engineers, 2012: 795−805. DOI: 10.1115/imece2012-88739.
[26]	MOORE D F, JÉRUSALEM A, NYEIN M, et al. Computational biology:modeling of primary blast effects on the central nervous system [J]. Neuroimage, 2009, 47: T10–T20. DOI: 10.1016/j.neuroimage.2009.02.019.
[27]	CHEN Y, OSTOJA-STARZEWSKI M. MRI-based finite element modeling of head trauma: spherically focusing shear waves [J]. Acta Mechanica, 2010, 213(1/2): 155–167. DOI: 10.1007/s00707-009-0274-0.
[28]	ZOGHI-MOGHADAM M, SADEGH A M. Global/local head models to analyse cerebral blood vessel rupture leading to ASDH and SAH [J]. Computer Methods in Biomechanics and Biomedical Engineering, 2009, 12(1): 1–12. DOI: 10.1080/10255840802020420.
[29]	NAHUM A M, SMITH R, WARD C C. Intracranial pressure dynamics during head impact [C] // Proceedings of 21st Stapp Car Crash Conference. Pennsylvania: Society of Automotive Engineers, 1977: 339−366. DOI: 10.4271/770922.
[30]	BENEDICT J V, HARRIS E H, VON ROSENBERG D U. An analytical investigation of the cavitation hypothesis of brain damage [J]. Journal of Basic Engineering, 1970, 92(3): 597–603. DOI: 10.1115/1.3425083.
[31]	WARD C, CHAN M, NAHUM A. Intracranial pressure: a brain injury criterion [R]. SAE Technical Paper, 1980. DOI: 10.4271/801304.

施引文献

期刊类型引用(22)

1.	贾时雨，王成，徐文龙，马东，齐方方. 环形复合内衬头盔冲击波防护性能研究. 兵工学报. 2025(01): 60-69 . 百度学术
2.	黄浩，崔海林，田晓丽，吴浩. 多孔结构对冲击波的衰减影响研究. 机械设计与制造工程. 2024(01): 11-15 . 百度学术
3.	田金，刘少宝，卢天健，徐峰. 持续性高过载下人脑的多孔弹性响应. 应用数学和力学. 2024(06): 691-709 . 百度学术
4.	杨昆，谭向龙，吴艳青，李梦阳，张钊，曾商鉴. 爆炸冲击波作用于生物体损伤的数值仿真研究进展. 兵器装备工程学报. 2024(09): 75-81 . 百度学术
5.	康越，马天，王俊龙，张逸之，张文博，韩笑，栗志杰. 不同海拔高度炮口冲击波动态演化特性数值模拟研究. 爆炸与冲击. 2024(12): 57-73 . 本站查看
6.	Rui Yuan，Yaoke Wen，Weixiao Nie，Dongxu Liu，Zhouyu Shen，Haoran Xu. Dynamic response of armor-piercing bullets to blunt and penetration with protective gelatin. Theoretical & Applied Mechanics Letters. 2024(04): 270-279 . 必应学术
7.	范志强，常瀚林，何天明，郑航，胡敬坤，谭晓丽. 基于PVDF复合压电效应的低强度冲击波柔性测量. 爆炸与冲击. 2023(01): 73-85 . 本站查看
8.	康越，马天，黄献聪，庄茁，柳占立，曾繁，黄超. 颅脑爆炸伤数值模拟研究进展：建模、力学机制及防护. 爆炸与冲击. 2023(06): 3-38 . 本站查看
9.	黄安，曹国鑫. 爆炸冲击波作用下均质颅骨模型有效性研究. 力学学报. 2023(08): 1774-1787 . 百度学术
10.	喻伯牙，高俊宏，王鸿，卢青，范小琳，李亮，李晓. 爆炸冲击波所致的肺损伤与脑损伤. 中国工业医学杂志. 2023(04): 332-335 . 百度学术
11.	鲁菁，屈媛媛，邵玉莹，郭述豪，冯楚文，孙维伯，李彬彬，孙冬玮，杨添淞. 创伤性颅脑损伤动物模型研究概况. 神经损伤与功能重建. 2023(09): 534-538+542 . 百度学术
12.	杜宁，赵梓淇，熊玮，刘闯，张先锋. 壳体厚度对装药爆炸冲击波特性影响研究. 弹道学报. 2023(03): 72-77 . 百度学术
13.	王博，温垚珂，徐诚，刘东旭. 钨合金破片侵彻防弹插板和明胶复合机理研究. 兵器装备工程学报. 2023(11): 38-46+96 . 百度学术
14.	蔡志华，贺葳，汪剑辉，王幸，张磊. 爆炸波致颅脑损伤力学机制与防护综述. 兵工学报. 2022(02): 467-480 . 百度学术
15.	郭建峤，王言冰，田强，任革学，胡海岩. 人体肌骨的多柔体系统动力学研究进展. 力学进展. 2022(02): 253-310 . 百度学术
16.	张文超，王舒，梁增友，覃彬，卢海涛，陈新元，卢文杰. 爆炸冲击波致颅脑冲击伤数值模拟研究. 北京理工大学学报. 2022(09): 881-890 . 百度学术
17.	聂伟晓，温垚珂，董方栋，覃彬，罗小豪，童梁成. 破片侵彻戴防弹头盔头部靶标钝击效应数值模拟. 兵工学报. 2022(09): 2075-2085 . 百度学术
18.	沈周宇，温垚珂，闫文敏，董方栋，张俊斌，李颖. 手枪弹撞击戴防弹头盔人体头颈部靶标的钝击效应. 兵工学报. 2022(09): 2101-2112 . 百度学术
19.	熊漫漫，覃彬，徐诚，安硕，伍杨. 冲击波作用有/无防护颅脑靶标动态响应规律. 兵工学报. 2022(09): 2182-2189 . 百度学术
20.	张文超，王舒，梁增友，覃彬，卢海涛，陈新元，卢文杰. 基于空气流场压力分析的头盔冲击波防护效能研究. 爆炸与冲击. 2022(11): 66-78 . 本站查看
21.	王小峰，陶钢，徐宁，王鹏，李召，闻鹏. 冲击波诱导水中纳米气泡塌陷的分子动力学分析. 物理学报. 2021(13): 283-301 . 百度学术
22.	康越，张仕忠，张远平，柳占立，黄献聪，马天. 基于激波管评价的单兵头面部装备冲击波防护性能研究. 爆炸与冲击. 2021(08): 179-191 . 本站查看