多孔材料内部孔隙占比可达30%，孔隙结构对其物理性能有着重要的影响^[1]。尤其是在力学过程中，多孔材料变形过程伴随着孔隙坍塌行为，发生显著的体积变化，因而具有优异的能量吸收特性^[2]。在冲击防护与吸能领域应用广泛^[3]。同时强动载冲击作用下多孔材料孔隙闭合过程中更容易形成局部热点，与传统密实材料相比更容易形成高温响应区，是高能量密度物理领域关注的焦点^[4]问题之一。实现对多孔材料冲击压缩特性参数的可靠预测对于指导该类材料设计与制备具有重要意义。

早期研究中，研究人员对多孔材料的孔隙坍塌行为主要采用“雪崩”模型^[5]进行简化，假设在零压时多孔材料已经压实，进而采取Grüneisen物态方程沿等容路径建立疏松材料的物态方程。Oh等^[6]沿等压路径构建了多孔材料物态方程。Wu等 ^[7-8]假设材料在恒压吸热过程中定压比热容保持常数，获得了与Grüneisen型物态方程形式相近的热力学物态方程。Boshoff等 ^[9]通过对比Wu-Jing型物态方程与其余几种疏松材料的物态方程的计算结果后，认为Wu-Jing方程的计算结果最好，且在此基础上，进一步提出以shchetinin方程描述冷压与冷能，使用Slater方程描述格林爱森系数，改进了Wu-Jing方程的求解流程；Geng等^[10-11]讨论了非谐振效应对材料比热容的影响，考虑等压路径下热电子效应以及计算冲击温升的方法，并给出了Wu-Jing参量的统计力学形式，阐明了其微观物理图像；Zhang等^[12]结合冷能叠加原理建立了混合物冷能曲线，建立了多孔混合物的冲击压缩特性计算方法，扩展了Wu-Jing型物态方程的计算对象；Zhou等^[13]沿等压线引入了材料在高温与室温下的比内能差值，建立了求解具有初始温度的材料冲击压缩特性方法；Wu等^[14]通过考虑孔隙内气体对比内能的贡献，建立了求解孔隙内含气体的疏松材料冲击压缩特性计算方法，计算结果与分子动力学模拟结果相符合。

多孔材料冲击试验中，研究者很早就观测到具有三波结构的波剖面，如1969年Boade^[15]在研究多孔钨的冲击波剖面时发现材料内存在具有稳定传播速度的两个前驱波及一个后继冲击波。类似的波剖面结构同样可以在其他多孔金属^[16-19]、聚合物^[20]和多孔陶瓷^[21-22]等多种不同多孔材料的冲击试验中观察到。但是，对多孔材料冲击压缩初始阶段存在的弹性前驱波、塑性波及冲击波“三波结构”^{[15-16, 18-19]}的理论研究较少，这可能导致多孔材料在低压区压缩特性计算结果误差较大。亟待对多孔材料孔隙坍塌行为与冲击压缩特性的联系展开理论分析并建立相应计算方法。

综上所述，本文中针对孔隙坍塌行为对多孔材料冲击物态特性的影响开展了研究，分析多孔材料中孔隙坍塌行为与冲击波结构间的内在联系，基于冲击波分离和冲击波追赶模型构建多孔材料冲击波结构的低压单波模式、双冲击波模式、高压单波模式理论模型。围绕3类冲击波模式建立与其相容的多孔材料冲击压缩特性计算方法，讨论了孔隙坍塌行为对多孔物质冲击压缩特性的影响规律。

1.   多孔材料冲击压缩作用过程分析

图1所示为Wang等^[18]在研究多孔铝的动态压实与断裂行为时使用VISAR记录到的具有“三波结构”的波剖面。压缩过程中，波剖面包含两个强度较低的前驱波和一个强度较高的后继冲击波。试验结果表明^{[15-16, 19]}，冲击过程中第1道前驱冲击波的波速与通过超声波测试测得的弹性纵波波速相近，其应视为弹性波在多孔材料内的传播结果，但第2道前驱冲击波以及后继冲击波的产生机制并不明确。

图  1  典型的具有三波结构的波剖面^[18]

Figure  1.  Typical wave profile with three-wave structure

下载: 全尺寸图片 幻灯片

Carroll等^[23]指出多孔材料在压缩过程中材料内部应力并非均匀分布，而是孔隙内表面附近的材料首先进入屈服。孔隙坍塌过程分为3个阶段：弹性阶段、弹塑性阶段、塑性阶段，仅在进入塑性阶段后材料孔隙度才会发生剧烈的变化。

结合上述试验现象以及孔隙演化相关理论，图2所示为具有三波结构的多孔材料冲击波波剖面示意图，第1道前驱冲击波对应着孔隙坍塌行为中材料未进入屈服的弹性阶段，而第2道前驱冲击波及第3道后继冲击波分别对应着孔隙坍塌行为中材料开始部分进入屈服的弹塑性阶段和材料孔隙度开始快速下降的塑性阶段。下文将具体分析孔隙坍塌的不同阶段对于冲击波演化及三波结构的联系。

图  2  具有三波结构的多孔材料冲击波波剖面示意图

Figure  2.  Schematic diagram of shock wave profile of porous material with three-wave structure

下载: 全尺寸图片 幻灯片

1.1   多孔材料冲击波演化过程分析

在多孔材料的冲击压缩过程中，基体材料比容变化与孔隙度变化均可引起多孔材料的比容变化，材料压缩曲线的性质受控于孔隙的坍塌行为。因此，孔隙坍塌行为将直接影响冲击波的演化过程。下面结合图3所示活塞模型分析材料内部冲击波演化过程。

图  3  不同压缩曲线的材料在压缩过程中的波系特征示意图

Figure  3.  The schematic diagram of wave system characteristics of materials with different compression curves in the process of compression

下载: 全尺寸图片 幻灯片

已有研究证明^[24]，材料压缩曲线的二阶导数决定了材料声速的变化趋势，若二阶导数为正，则材料声速随着压缩度的提高而增加，若二阶导数为负，则材料声速随着压缩度的提高而减小。如图3(a)所示，常规材料的压缩曲线表现为上凹型，二阶导数恒大于零，在活塞加速运动条件下，压缩波将追上先产生压缩波，其特征线将发生汇合，由此在波阵面上将形成间断的力学参量，即冲击波。而对于具有上凸型压缩曲线的物质，其在受到冲击时波系特征如图3(b)所示，其特征线是发散的，在冲击加载时只产生连续的压缩波而无法产生间断的冲击波。

如图4所示，在冲击压缩过程中多孔材料的压缩曲线在弹塑性阶段发展的过程中B点以及弹塑性阶段向塑性阶段转化时C点处具有上凸结构，并因此产生了两道分离的前驱冲击波以及后继分离的冲击波。

图  4  多孔材料冲击绝热线示意图

Figure  4.  Schematic diagram of Hugoniot curve of porous material

下载: 全尺寸图片 幻灯片

B点上凸结构的产生过程如下，在弹性前驱波OE的波后，对材料的进一步压缩将使得多孔材料的孔隙坍塌行为由弹性阶段进入弹塑性阶段，孔隙内侧部分材料进入塑性，材料的纵波波速趋近于塑性材料纵波波速，此时纵波波速随材料密度提升而降低，多孔材料的压缩曲线表现为上凸型的特征，压缩过程中压缩波特征线发散，无法形成冲击波，此时的波系结构为弹性前驱波及后方特征不交汇的简单压缩波。随着加载压力进一步提高，对已进入塑性状态的材料，其体积模量随压缩度的上升而上升，因此随着体胞中塑性区占比的进一步提升，存在一个时刻，即B点，纵波波速随塑性部分材料密度的提升而提升，材料压缩曲线重新表现为上凹型特征，压缩波特征线交汇，物质在压缩过程中能够以B点参量为参考形成冲击波，此冲击波即为试验中观测到的第2道前驱冲击波。

C点上凸结构的产生过程如下：当驱动压力位于C点上方时，驱动压力超过了临界压力，使得孔隙的坍塌行为进入了塑性区。此时，多孔材料仍然为上凹型，但是由于孔隙此时随着压力迅速坍塌，材料的声速在C点迅速减小，此时C点可以看作是一段极小的上凸型曲线，在C点两侧发生了冲击波分离的现象，形成了第三道后继冲击波。

综上，由于孔隙坍塌过程中先后经历了弹性阶段、弹塑性阶段以及塑性阶段，多孔材料在压缩过程中形成了由第1道弹性前驱波、第2道前驱冲击波、第3道后继冲击波构成的“三波结构”。

1.2   多孔材料冲击波追赶问题及冲击波模式划分

在图4所示的压缩曲线中，随着终态压力进一步提高，达到H点时，此时瑞利线CH的斜率大于瑞利线BC，后继冲击波波速高于前驱冲击波波速。因此在终态压力足够高时，多孔材料中需要考虑冲击波的追赶问题。

为了便于讨论后继冲击波对前驱冲击波时的追赶问题，且考虑到弹性前驱波强度较低^{[15, 16]}，进一步对上一小节内的模型进行简化，即仅在描述孔隙的坍塌行为时考虑材料的弹塑性行为，对材料的变形行为采用流体模型处理，不再考虑材料的弹性前驱波，获得了图5所示的“双波结构”简化模型。如图5所示，材料压缩曲线在弹塑性阶段表现为上凹型，在压力逐渐提高并达到临界压力时，由于孔隙度的快速下降，发生冲击波分离现象，冲击波表现为双波模式，由一个由临界压力决定的前驱冲击波和一个后继冲击波组成。当后继冲击波的波速大于前驱波波速时，即终点压力高于M点时，材料中将发生冲击波追赶现象。

图  5  多孔材料简化模型冲击绝热线示意图

Figure  5.  Schematic diagram of Hugoniot curve of porous material in Simplified mode

下载: 全尺寸图片 幻灯片

如图6所示，在前驱冲击波的作用下，物质从未被扰动的O区进入E区，再经过后继冲击波的作用下进入H'区。由于后继冲击波波速较快，其将在某一时刻追上前驱冲击波并在两者接触后，形成一个前向冲击波，其波前是未被扰动的O区介质，波后是H区，同时物质内还会形成后传冲击波或者后传稀疏波，以使H'区介质进入H''区，H'区与H''区之间通过接触间断面隔开，接触间断面两侧比容、温度等参量间断，而介质压力、质点速度一致，因此在物质坐标上看，接触间断面不向前传播。

图  6  冲击波追赶问题示意图

Figure  6.  Characteristic line diagram of shock wave overtaking problem

下载: 全尺寸图片 幻灯片

假设冲击加载的速度足够快，总是在形成前驱波之后在很短的时间内形成一个强度稳定的后继冲击波，即图6中的${x_0}$足够小。前驱冲击波在${x_{\text{1}}}$处完成追赶，形成前向冲击波与驻定的接触间断面，在后续冲击波传播中大部分材料中处于H区，此时可以认为多孔材料内仅有一个单一冲击波在传播。因此，基于上述分析，可以将多孔材料的冲击响应如图7所示划分为低压单波模式、双冲击波模式、高压单波模式。

图  7  多孔材料中的三种冲击波模式

Figure  7.  Three shock wave modes in porous materials

下载: 全尺寸图片 幻灯片

2.   多孔材料冲击压缩特性计算模型

2.1   考虑冲击波结构的冲击压缩特性理论计算模型

为了充分考虑材料的孔隙坍塌行为，采取可沿等压路径求解的Wu-Jing方程建立多孔态材料和密实态材料冲击绝热线之间的关系，其形式如式(1)、式(2)所示。

式中：R为Wu-Jing参数，其是压力p的函数。${V_{\text{H}}}$、${V_{\text{c}}}$、${H_{\text{H}}}$、${H_{\text{c}}}$和$V_{\text{H}}'$、${V_{\text{c}}}$、$H_{\text{H}}'$、${H_{\text{c}}}$分别代表着密实材料和多孔材料的冲击比容、零温比容、冲击比焓、零温比焓。

比焓的定义为H=E+pV，根据冲击波结构的不同，并假设初始比内能与压力可以忽略，多孔物质的比内能及比焓可以表达为如下形式，其中式(3)对应双冲击波模式后继冲击波波后状态、式(4)对应单冲击波模式冲击波波后状态。

式中：下标为E的物理量为前驱冲击波波后的物理量，${V_{00}}$为多孔材料的初始比容，$E_{\text{c}}'$为多孔材料的冷能。

假设压力相等时密实态物质和多孔态物质中的Wu-Jing参量相等，密实材料与多孔材料的冷能相等，即有$E_{\text{c}}' = {E_{\text{c}}}$，通过联立上述关系式，获得冲击比容表达式，形式如下，其中式(5)对应双冲击波模式的后继波波后状态，式(6)对应单冲击波模式的波后状态：

式(5)与Wu^[7]等提出的冲击比容计算式是一致的，区别在于本文模型中明确了前驱波是由孔隙坍塌引起的第2前驱冲击波。式(6)与Geng等^[10]基于Wu-Jing方程给出的简化冲击比容计算式相似，其形式如式(7)所示：

式(6)与简化冲击比容计算式的主要差别在于方程的第三项，若无该项，简化冲击比容计算式的起点将为密实态的初始比容点，其物理假设与“雪崩”模型相似，忽略了多孔材料的孔隙坍塌行为。同时，式(6)使得求解式(5)中的${V_{\text{E}}}$时通过对该式带入${p_{\text{E}}}$即可求解得到${V_{\text{E}}}$，而在原有的求解模型中，需要引入${V_{\text{E}}}{\text{ = }}{V_{{\text{00}}}}$的近似条件，因此在本文中模型在分析前驱冲击波的波速时具有更明确的物理意义和计算精度。

处于双波模式时，后继冲击波的波前是已受前驱冲击波扰动的物质，其Hugoniot关系的计算式要发生变化。根据冲击波基本关系式，Hugoniot关系的计算式如下：

式中：$D$为冲击波速度，$U$为质点速度，${P_{\text{H}}}$、${V_{\text{H}}}$分别为冲击绝热线上的压力与比容，下标为E的物理量为前驱冲击波波后的物理量。其中式(8)对应单冲击波结构，式(9)对应双冲击波结构，式(10)对应前驱冲击波的参量。

2.2   孔隙坍塌$p{\text{-}}\alpha$模型

为了求解式(5)～(6)，需要计算多孔材料在孔隙坍塌过程中的冷容。针对孔隙坍塌时材料比容的变化，相关学者^{[23, 25, 26]}已经开发了许多物理模型，其中以p-α类模型以压力p与孔隙度$\alpha$为变量，物理意义明确且计算简便。在选用p-α类模型时，多孔材料冷容$V_{\text{c}}'$与孔隙度的关系为$V_{\text{c}}' = \alpha {V_{\text{c}}}$。

p-α模型可以较好的反映材料孔隙坍塌现象，随着加载压力的提高，孔隙度在弹性阶段、弹塑性阶段逐渐减小，并在塑性阶段快速减小。Carroll等^[23]在采用静态球面模型分析时指出孔隙度在弹性区、弹塑性区的变化较小，在容许一定误差的情况下，可以忽略材料弹性阶段以及弹塑性阶段引起的孔隙度变化，认为孔隙度可用下式描述：

式中：${\sigma _{\text{y}}}$为基体材料屈服强度、$\alpha$为材料当前状态孔隙度、${\alpha _0}$为初始孔隙度，${p_{\text{E}}}$为孔隙坍塌临界压力，其计算式为：

Johnson等^[26]指出在多孔材料的受力单元中，孔隙的存在降低了材料的受力面积，实际作用于孔隙坍塌模型上的压力与受力分析时的单元压力之间应有如式(13)指出的关系：

式中：$\bar p$为作用于静态球面模型上的压力。

基于静态球面模型，忽略多孔材料的孔隙拓扑结构或孔隙形貌的影响，在坍塌过程采用以下孔隙坍塌模型描述材料的p-α关系：

式中：$G$为基体材料剪切模量、${\alpha _{\text{E}}}$为弹塑性状态分界点孔隙度。该修正孔隙坍塌模型中，孔隙度线性地减小到多孔材料弹塑性阶段与塑性阶段的分界，并考虑了孔隙坍塌模型中的压力与受力单元上的压力之间的关系，因此相比Wu^[7]等模型中采取的简化p-α关系式(11)，式(14)～(15)充分考虑了多孔材料的弹塑性阶段的孔隙演化行为，有利于正确反映孔隙坍塌行为。

3.   计算结果与分析

冲击压缩特性计算中涉及了铜、钨、铁等3种材料，相关材料参数如表1所示，表中：Q、q分别为描述材料冷压、冷能的玻恩-迈耶势的所需的参数，使用最小二乘法对相应材料的冲击压缩数据拟合获得，Wu-Jing参数的计算方法参考Geng的工作^[10]。

表  1  计算密实材料Hugoniot数据所需的材料常数

Table  1.  The material constants required for dense materials to calculate Hugoniot data

材料	ρ0/(g∙cm−3)	Q/GPa	q	σy/GPa	G/GPa
铜	8.924	52.90	10.0696	0.6[27]	47.7[27]
钨	19.2	185.27	7.1086	2.2[27]	160[27]
铁	7.856	34.94	11.9378	0.68[28]	89[28]

下载: 导出CSV 

| 显示表格

本节在模型验证、规律分析时引用了前人的试验工作^[10]，试验及相关条件如下：多孔钨、多孔铜的试验数据来源于Boade，多孔材料的孔隙尺寸分别为小于1 μm与小于10 μm，试验中采取气枪发射与炸药透镜的方式实现不同的加载压力，最大加载压力可达105 Gpa；多孔铁的相关试验数据来源于Butcher，多孔材料的孔隙尺寸为20 μm，试验中采取气枪发射实现加载，最大加载压力5 Gpa左右。

3.1   模型验证

计算材料Hugoniot数据时考虑3种计算模型：(1)考虑孔隙坍塌行为而不考虑双波模式的Viljoen模型；(2) 既不考虑孔隙坍塌行为也不考虑双波模式的简化Wu-Jing模型,即式(7)；(3)根据本文中方法，判断冲击波模式并选择对应的方程进行计算。

Viljoen的模型的关系式为：

采用Boade等^[15]对初始密度为12.64 g/cm³多孔钨与Butcher等^[29]对初始密度为5.76 g/cm³的多孔铁开展的冲击压缩试验进行对比，冲击压缩特性计算结果与试验结果的对比如图8所示。在计算中3种方法参考的密实态材料冲击绝热线、冷压线与冷能线是一致的。

图  8  三种计算方式下Hugoniot数据计算结果与试验结果对比

Figure  8.  Comparison of the experimental results and computation of Hugoniot data under three calculation methods

下载: 全尺寸图片 幻灯片

从图中可以看出，采取简化的Wu-Jing模型与Viljoen模型进行计算时，计算的冲击波速度与试验数据偏差较大，而本文方法的冲击波波速的计算结果与试验数据符合较好，能够反映在多孔材料冲击压缩数据中低压段的反常的低冲击波速度。这说明考虑双波模式及孔隙坍塌行为对计算多孔材料冲击压缩特性是必要的。

根据式(13)与式(24)可以看到，多孔材料的双波模式的准确计算依赖于对前驱冲击波参数的准确计算。根据简化孔隙坍塌模型(式(11))与本文孔隙坍塌模型(式(14))计算了双波模式下的前驱冲击波波速并与试验值^{[15, 16]}进行了对比，计算结果如图9所示。计算结果显示，简化孔隙坍塌模型计算得到的结果与前驱冲击波波速的试验结果偏差过大。而本文得到的修正孔隙坍塌关系计算的前驱冲击波波速与试验值较为接近。

图  9  不同的孔隙坍塌模型下前驱冲击波波速与试验结果^[15-16]的对比

Figure  9.  Comparison of precursor shock wave velocity with experimental results under different pore collapse models

下载: 全尺寸图片 幻灯片

综上所述，本文考虑材料孔隙坍塌行为的计算模型能够较好地反映多孔材料的冲击波结构，可以准确计算得到相关冲击波参数。

3.2   孔隙坍塌行为对冲击波参数的影响分析

多孔材料的基体强度、初始疏松度、剪切模量均可能影响材料的前驱冲击波波速。下面以多孔钨为例，分析不同屈服强度多孔材料的剪切模量、疏松度对前驱冲击波波速的影响。结果如图10所示，在不同屈服强度下的多孔材料，其前驱冲击波波速均随着剪切模量的提升而提升，随孔隙度的提升而下降，且屈服强度的提高有利于前驱冲击波波速的提高。

图  10  多孔钨的不同屈服强度下前驱冲击波波速与剪切模量及初始孔隙度的关系

Figure  10.  The relationship between precursor shock wave velocity and shear modulus and initial porosity under different yield strength

下载: 全尺寸图片 幻灯片

本文中所选取的孔隙坍塌模型中，孔隙度的变化主要发生在基体材料进入塑性区后，由式(14)可知，材料的孔隙坍塌行为主要由材料基体强度与初始孔隙度控制。对不同基体屈服强度下初始密度为6.052 g/cm³多孔铜的冲击压缩试验特性进行了计算并与试验结果^[16]作了对比，计算结果如图11所示。计算结果显示在冲击压力较低时，基体屈服强度对多孔材料Hugoniot数据的计算结果有较明显的影响，随着基体屈服强度的提高，冲击波速度起点增大，反映阻抗匹配特性的冲击压力-质点速度曲线提高，当冲击压力较高时，不同基体屈服强度的Hugoniot数据趋于一致。

图  11  不同屈服强度下多孔铜Hugoniot数据计算结果与试验结果对比

Figure  11.  Comparison of the experimental results and computation of Hugoniot data of porous copper under different yield strength

下载: 全尺寸图片 幻灯片

不同初始孔隙度多孔铁的Hugoniot数据计算结果如图12所示。计算结果显示初始孔隙度对Hugoniot数据有较明显的影响。多孔材料的冲击波速度-质点速度曲线表现出明显的非线性特点，并且随着初始孔隙度的提高，冲击波波速的起点降低，冲击压力-质点速度曲线变得平缓，反映出较高孔隙度的多孔材料在冲击时具有更低的界面压力，具有更好的吸能特性。

图  12  不同初始孔隙度下多孔铁Hugoniot数据计算结果与试验结果对比

Figure  12.  Comparison of the experimental results and computation of Hugoniot data of porous iron under different initial porosity

下载: 全尺寸图片 幻灯片

4.   结　论

(1) 孔隙坍塌行为导致了多孔材料具有反常压缩曲线，导致多孔材料的冲击波形成过程中将出现1个前驱冲击波和1个后继冲击波，且在压力进一步提高时发生后继冲击波对前驱冲击波的追赶。因此，多孔材料中存在3种冲击波模式，分为低压单波模式、双冲击波模式、高压单波模式。

(2) 基于Carroll的工作建立了修正孔隙坍塌关系，结合Wu-Jing物态方程与冲击波基本关系，建立了与多孔材料中的不同冲击波模式相容的Hugoniot数据计算模型。其中与单波模式相容的计算模型解决了低压时Wu-Jing方程不能计算得到临界比容的问题。本文的计算模型与已有试验结果吻合较好，反映了多孔材料的前驱冲击波波速和低压区的反常冲击波波速。

(3)基于本文模型分析了孔隙坍塌关系中基体屈服强度、剪切模量、初始孔隙度对冲击波参量的影响规律。前驱冲击波波速随基体屈服强度、剪切模量增大而提高，随初始孔隙度增大而降低；初始孔隙度、基体屈服强度均可影响材料的后继冲击波Hugoniot数据，冲击波速度与冲击压力的计算结果与基体强度正相关、与初始孔隙度负相关；尽管终态压力提高，不同基体屈服强度的Hugoniot数据趋于一致，但是对多孔材料的低压区Hugoniot数据的计算不应忽略基体屈服强度的影响。