混凝土是最普遍应用的结构工程材料，由水泥石、标准砂和集料构成，水泥石是从微观到细观结构复杂的复合物，集料尺寸则可达厘米以上，它们构成的水泥基体是一个多组元多相多孔、含有大量微孔洞微裂纹的复合体，具有跨尺度非均质多组元多相复杂结构。其动力学特性表现出强烈的应变率效应、尺寸效应、加载路径效应、静水压效应等特性^[1-4]，应用分离式霍普金森压杆（split Hopkinson pressure bar, SHPB）进行混凝土动态性能研究绕不开这些效应的影响。相关研究较多，但仍存在较大争议。

冲击下材料内部波动的传播使得材料内部处于非平衡状态。作为SHPB实验的基本假定之一的均匀性假定^[5]，也允许试件两端存在一定的载荷差异，例如5%。而在试样受到冲击的过程中，试样内部沿冲击方向存在纵向惯性力，垂直冲击方向存在横向惯性力或横向运动的动能，试样与杆的接触端面存在摩擦力，并且试样本身可能存在一定的黏性。纵向惯性效应和黏性效应要求试样的厚度不能太厚；横向惯性效应和接触端面的摩擦作用要求试样的直径不能太粗^[5-7]。但是，对于混凝土试样而言，由于水泥石、集料等大颗粒结构成分的存在，为使得试样具有较好的代表性，试样尺寸需要足够大^[1]，例如，直径为50 mm^[8-10]、74 mm^[1]、80 mm^[11]的圆柱形试样，甚至更大尺寸的试样都进行了尝试。试样尺寸越大，带来的问题越多，突出的问题是：此时如何正确认识测试得到的强度特性。

另外，混凝土试样SHPB冲击实验结果表明：应变率高于10～10² s⁻¹时，动态强度会迅速提高^[10,12]，这一般归结于黏性机制控制的应变率效应。Li等^[13]基于有限元法和Drucker-Prager (DP)模型研究静水压力对混凝土材料动态强度的影响。他们采用的是一种应变速率不敏感的材料模型，其计算结果表明：材料宏观动态强度的增强是由横向惯性约束引起的，而不是材料的应变率敏感性。而如果将这种动态强度的增强归因于应变率效应，会造成极大的计算误差。他们后续的工作从实验和数值分析方面进行了更深入的研究^[14-15]，通过引入圆环试样表明引入动摩擦模型的必要性^[16]。杨茨等^[17] 基于金属圆环试样进行数值分析，结果发现：在界面摩擦因数的变化过程中，惯性效应起主要作用，应变率效应起次要作用；存在着一个表征圆环试样属性的临界摩擦因数。但由于圆环状混凝土试样的加工难度较大，进行相关的尝试较少。这些数值模拟和实验的尝试已经可以说明，横向惯性效应含有应变率的作用。

此时，问题又回到原点，即：如何消除惯性效应的影响，尤其是大直径试样的横向惯性效应的影响。进一步，如何评估大直径试样动态强度中应变率效应的影响所占的比重。Forrestal等^[18]推导了大直径试样冲击过程中横向应力的应变率表达式，Xu等^[10]给出了受约束试样横向应力的应变率表达式，结合实际冲击实验可初步估计试样的横向约束效应。在实验技术方面，Grote等^[11]、Forquin等^[19]引入了一种金属环限制的高应变率下的混凝土实验装置。金属环的限制可以显著减小试样的径向变形，以减小横向运动的动能，因此，在一定程度上达到了减小横向惯性效应的目的。进一步，Grote等^[11]估计出应变率在高达10⁴ s⁻¹的条件下，横向惯性约束引起的强度升高约为52%。这比Xu等^[10]基于真三轴静载的混凝土冲击实验得到的应变率50 s⁻¹下的6%和应变率100 s⁻¹下的14.5%高得多，呈现显著的应变率效应。由于Grote等^[11]的高应变率实验中试样处于一维应变状态，因此，静水压效应带来的横向约束也是一个不可忽视的影响因素。

应变率效应、尺寸效应和静水压效应在混凝土动态强度的表征上都起到了一定的作用，并且各个效应之间具有一定的耦合作用，如何对它们进行评估存在较大的困难。本文中，将基于真三轴应力状态（σ_xs≠σ_ys≠σ_zs）作用下混凝土冲击实验结果^[10]，分别结合考虑应变率效应的HJC模型^[20-21,23]和考虑静水压效应的Drucker-Prager (DP)模型^[22]进行数值分析，对这3种效应进行探讨。

1.   混凝土冲击实验与数值分析模型

1.1   真三轴静载下混凝土冲击实验装置

实验设备如图1所示，包括2部分，即：（1）真三轴静载施加系统，分别由3个方向的液压缸和对应的反力支架组成，可对立方体试件施加3向不等的压应力；（2）撞击杆发射和信号测试系统，主要由冲击方向（x方向）的高压气炮、入射方杆、支撑方杆、垂直于冲击方向的水平y方向左支撑方杆、y方向右支撑方杆以及垂直于冲击方向的z方向下支撑方杆和z方向上支撑方杆等组成。图中，实验控制台控制液压系统和发射系统，液压站提供伺服控制过程3个液压缸中的液压油。图1(b)中，ε_i、ε_r和ε_t分别为x方向杆上入射、反射和透射应变信号；ε_y-left和ε_y-right分别为y方向2根杆上的应变信号。

图  1  真三轴静载混凝土冲击实验设备

Figure  1.  The experimental device for concrete specimens under true tri-axial confinement

下载: 全尺寸图片 幻灯片

首先，通过液压系统和各方向方杆对立方体试件3个方向分别施加预定的静载σ_xs、σ_ys和σ_zs，各个方向施加的静载值可直接通过控制界面读取。然后，x方向的高压气体驱动撞击杆，撞击入射方杆，杆中产生应力波作用于试件。记录撞击杆速度和6根杆上的波动信号，由此通过与传统Hopkinson杆相同的实验思路可得到冲击过程中试件6个面的动态响应，计算3个方向的扰动载荷σ_xd、σ_yd和σ_zd，得到试件3个方向的应力-应变关系。最后，由3向应力σ_x=σ_xs+σ_xd，σ_y=σ_ys+σ_yd和σ_z=σ_zs+σ_zd，估计混凝土材料的动态强度。关于混凝土试件的详细实验结果见文献[8-10]。

1.2   数值分析模型

真三轴静载冲击实验系统的数值模拟模型如图2所示：6根方杆的截面尺寸均为50 mm×50 mm，x、y、z方向各有2根长度相同的钢杆夹持试样，杆长分别为2.5、2.0、1.5 m；试件为50 mm×50 mm×50 mm的正方体。模型尺寸与实验尺寸保持一致。

图  2  有限元计算模型

Figure  2.  The finite element calculation model

下载: 全尺寸图片 幻灯片

在3个方向，杆的一端分别设置80 mm×10 mm的挡板为该方向施加静载时提供反力约束。杆与混凝土之间采用自动面面接触，界面摩擦因数取为0.1。

钢杆和挡板均采用线弹性模型：密度为7 850 kg/m³，弹性模量E=210 GPa，泊松比ν_steel=0.3。为对比研究应变率效应和静水压效应的影响，分别运用HJC模型和DP模型来描述混凝土材料。

HJC本构模型包含强度方程、损伤演化方程和状态方程。

（1）强度方程为：

式中：σ^*=σ/ f′，σ为材料真实应力，f′为材料单轴抗压强度；A为特征化黏性强度因数，D为损伤因子，B为特征化压力硬化指数；p^*=p/ f′，p为真实压力；N为压力硬化指数，C为应变率影响参数；${\dot \varepsilon ^*} = {{\dot \varepsilon } \mathord{\left/ {\vphantom {{\dot \varepsilon } {{{\dot \varepsilon }_0}}}} \right. } {{{\dot \varepsilon }_0}}}$，$\dot \varepsilon$为真实应变率，${\dot \varepsilon _0}$为参考应变率； S_max为无量纲化最大等效强度。

（2）损伤演化方程为：

式中：$\Delta {\varepsilon _{\text{p}}}和\Delta {\mu _{\text{p}}}$分别为一个计算步长内的等效塑性应变的改变量和塑性体积应变的改变量；$\varepsilon _{\text{p}}^{\text{f}} + \mu _{\text{p}}^{\text{f}}$为混凝土材料在压力p作用下的断裂塑性应变，$\varepsilon _{\text{p}}^{\text{f}} + \mu _{\text{p}}^{\text{f}} = {D_1}{({p^*} + {T^*})^{{D_2}}}$，D₁、D₂为损伤常数，T^*=T/ f′，T为材料抗拉强度。

（3）状态方程为：

式中：K_e=p_c/μ_c，p_c、μ_c为弹性极限下的静水压力和体积应变；μ为体积应变；K_l=(p_l−p_c)/(μ_l−μ_c)，p_l、μ_l分别为压垮时的静水压力和体积应变；K₁、K₂、K₃为压力常数，$\bar \mu = {{(\mu - {\mu _{\text{l}}})} \mathord{\left/ {\vphantom {{(\mu - {\mu _{\text{l}}})} {(1 + {\mu _{\text{l}}})}}} \right. } {(1 + {\mu _{\text{l}}})}}$。

混凝土HJC模型的具体参数^[20-23]见表1，表中ρ₀为密度，G为剪切模量。运用*MAT_ADD_ EROSION关键字中的剪切应变破坏作为混凝土单元破坏的判据。

表  1  混凝土HJC本构模型参数^[20-23]

Table  1.  Parameters of the HJC model for cement mortar^[20-23]

ρ0/(kg·m−3)	G/GPa	A	B	C	N	$f'$/MPa
2 440	14.875	0.69	1.50	0.007	0.61	48
T/MPa	${\dot \varepsilon _0}/{{\text{s}}^{ - 1}}$	εf,min	Smax	pc/MPa	μc	pl/MPa
4	1	0.01	7	594.3	0.03	800
μl	D1	D2	K1	K2	K3
0.1	0.04	1	85	−171	208

下载: 导出CSV 

| 显示表格

ABAQUS对经典的Drucker-Prager模型进行了扩展，包括线性Drucker-Prager模型、双曲线Drucker-Prager模型和指数Drucker-Prager模型，本文中采用线性Drucker-Prager模型作为混凝土材料的本构。线性Drucker-Prager模型包括屈服面方程和塑性势面方程。

（1）屈服面方程为：

式中：q为Mise应力，k′为三轴拉伸与压缩屈服应力之比，r为偏应力不变量，p为等效围压应力，ϕ为摩擦角，d为材料黏聚力。

（2）塑性势面方程为：

式中：ψ为膨胀角。

混凝土DP模型的具体参数见表2，表中ρ₀为密度，E为弹性模量，ν_con为泊松比。

表  2  混凝土DP模型参数

Table  2.  Parameters of the DP model for cement mortar

ρ0/(kg·m−3)	E/GPa	νcon	ϕ/(o)	k′	ψ/(o)
2440	35.7	0.3	30	0.78	30

下载: 导出CSV 

| 显示表格

2.   应变率效应和静水压效应

2.1   有限元模型

分别采用表1中的HJC模型参数和表2中的DP模型参数对实验波形进行模拟，结果如图3所示。实验对应的冲击速度为13.2 m/s，三向静载[σ_xs, σ_ys, σ_zs]为[15 MPa, 6 MPa, 10 MPa]。由于混凝土材料的非均匀性，从图3(c)～(d)可以看出，在试样同方向的2个面上的应力时程曲线有一定差异。采用描述混凝土骨料分布的细观有限元方法可以反映这种差异，但由于会增加更多的不确定参数，本文中将采用宏观等效的连续介质力学的有限元进行分析。模拟结果表明：采用上述2种模型参数，入射波形态相似，计算得到的波形与各个杆上的测试波形基本一致。图4为无侧限、双向侧限和三向侧限等3种情况下试样的破坏形态，其中的数值模拟均采用HJC模型。在无侧限冲击下，试样发生严重碎裂，如图4(a)所示；图4(a)中的有限元计算结果有相似的碎裂趋势。x轴和y轴双向侧限、x方向冲击下， z方向处于自由状态，试件破坏呈现出类似于成层剥离的特点，碎块剥离层面与z轴垂直，如图4(b)所示；图4(b)中的计算结果表现出相似的剥离趋势。三向侧限、x方向冲击下，试件不产生明显的宏观破坏，其失效为材料内产生局部细微观破裂，如图4(c)中的白色条带所示；图4(c)中的计算结果表现出相似的破坏形态。

图  3  数值模拟波形

Figure  3.  Simulated wave profiles in three directions

下载: 全尺寸图片 幻灯片

图  4  不同侧限条件下混凝土试样的破坏形态

Figure  4.  Failure patterns of concrete samples under different confinement conditions

下载: 全尺寸图片 幻灯片

2.2   应变率效应和静水压效应

图5为2种强度模型模拟的应力-应变关系曲线。图5(a)对应3种冲击速度下的模拟，图5(b)对应冲击速度为15 m/s和3种摩擦角下的模拟，施加的静载均为[10 MPa, 10 MPa, 10 MPa]。图5(a)为HJC模型计算的结果，随着应变率从40 s⁻¹升高到80 s⁻¹，混凝土试件的动态强度迅速提高，破坏应变也迅速增大。图5(b)为DP模型计算结果，随着内摩擦角ϕ从20°增大到40°，混凝土试件的动态强度迅速提高，但由于静水压随之提高，其破坏应变迅速减小。这说明，提高应变率和静水压均可使混凝土试件的动态强度升高。

图  5  2种模型模拟得到的应力-应变关系

Figure  5.  Stress-strain relations simulated by two strength models

下载: 全尺寸图片 幻灯片

混凝土类材料常用的强度准则有Mohr-Coulomb (MC)准则和Drucker-Prager (DP)准则，二者有一定的联系。MC准则和DP准则的表达式分别为：

式中：τ_n、σ_n分别为剪切面上的切向应力和法向应力，c为内聚力，ϕ为内摩擦角；I₁为第一应力不变量，J₂为第二应力偏量不变量，α和k为材料的强度参数。且：

根据式(7)，I₁和$\sqrt {{J_2}}$呈线性关系，对应用HJC模型计算得到的数据点进行线性拟合，如图6所示。由此得到不同应变率所对应的α和k，该参数随应变率的变化如图7(a)所示。随着应变率的升高，α逐渐变小、k逐渐变大，并且k近似与应变率呈线性关系。根据式(9)和式(10)，计算得到不同应变率对应的强度参数c和ϕ，如图7(b)所示。参数c和ϕ的发展规律与α和k的变化相对应，即内聚力c随着应变率的升高而升高，摩擦角ϕ随着应变率的升高而降低。这与文献[10]的实验结果趋势一致，但数值模拟得到的α和摩擦角ϕ偏小、k和内聚力c偏大。这表明，可能实际实验混凝土试样内部的非均匀程度非常高。因此，在数值模拟中应用均匀化的等效连续介质模型对混凝土进行量化分析的方法需要进一步提高。

图  6  不同应变速率下$\sqrt {{J_2}} \text{-} {I_1}$平面的强度统计

Figure  6.  Statistics of the strength in the $\sqrt {{J_2}} \text{-} {I_1}$ plane at different strain rates

下载: 全尺寸图片 幻灯片

图  7  计算强度参数的应变率效应

Figure  7.  Strain rate effect of simulated strength parameters

下载: 全尺寸图片 幻灯片

不过，从图6～7的计算结果来看，随着应变率的升高，α逐渐变小、k逐渐变大；由式(7)和图6可得出I₁和$\sqrt {{J_2}}$也应该相应增大；而等效应力${\sigma _{\text{e}}} = \sqrt {3{J_2}}$。因此，随着应变率的升高，混凝土的强度会提高，而且，这种强度的提高，也有一部分原因是I₁的增大所导致的。综上可以说明，混凝土试件的应变率效应含有静水压作用。反过来，类似Li等^[13]的工作，利用DP模型进行计算，可以得到混凝土试件的静水压效应也含有应变率的作用，也即应变率效应和静水压效应具有较强的耦合作用。

3.   惯性效应与应变率效应和静水压效应的关系

3.1   横向惯性效应

轴向冲击时，由于横向惯性效应会产生横向应力的分布。Forrestal等^[18]推导了大直径试样冲击过程中横向应力σ_r的应变加速度表达式为：

式中：r为以中心为起点的横向坐标，a₀为试样半径，ν为泊松比，ρ₀为初始密度，${\ddot \varepsilon _x}$为冲击方向的应变加速度。对于有侧限的实验条件，初始条件和边界条件为：

式中：p₀为作用于试件外表面的静水压，u¹为横向位移。可以得到受约束试样横向位移和应力的应变加速度表达式分别为^[10]：

式(14)～(15)表明，试件内部沿垂直于冲击方向的横向位移和横向应力的分布与冲击方向的应变加速度相关。Forquin等^[19]使用金属环对试件进行限制，以减少横向位移；式(14)表明，这种限制可以降低混凝土试样的横向惯性效应。Li等^[16]使用管状试件来减小横向惯性。对于管状结构而言，r和a₀较小，可以降低试样的横向惯性效应，这可以由式(15)进行解释。基于此，本文中结合实际实验中测试得到的冲击方向的应变时程曲线，可利用式(15)进行横向应力分布的计算，计算结果如图8所示。其中，在应变率50和200 s⁻¹下，在三轴围压值分别为0、7.5、45.0 MPa下进行了计算；在应变率100 s⁻¹下，在三轴围压值分别为0、7.5、30.0 MPa下进行了计算。由图8(a)可见，横向应力随着围压和应变率的升高有升高的趋势，其分布为中心大两边小。为明确冲击过程的响应特性，将不考虑三轴静载时横向应力的结果列入图8(b)。由图8(b)可见，随着应变率和围压的升高，横向应力分布的幅值均有所增大，反映出一定的应变率和围压耦合特性。

图  8  横向应力分布

Figure  8.  Distribution of the transverse stress

下载: 全尺寸图片 幻灯片

3.2   尺寸效应

波在较大尺寸试件中传播时会产生几何弥散效应。Rayleigh推导了考虑几何弥散的谐波近似解^[5]，得到与圆频率ω对应的谐波的相速度C为：

式中：r_g为截面对冲击轴的旋转半径，λ为与圆频率ω对应的波长。

不同频率的谐波由于其相速度不同在传播的过程中会分散开来，产生几何弥散现象；在大直径杆中传播时会产生较大的局部波形振荡^[24]，从而造成测试强度偏高。同时，较大尺寸对试样宏观力学性能也有一定的影响。控制应变率和围压状态不变，进行边长分别为37、50、74 mm等3种尺寸的正方体试样的数值模拟。图9(a)～(b)为三向静载[0, 0, 0]、应力脉冲幅值为100 MPa、应变率为40 s⁻¹的情况下，试样纵向应力σ_x和横向应力σ_y沿试样横向的分布情况。图9(c)为三向静载[5 MPa, 5 MPa, 5 MPa]、应力脉冲幅值为140 MPa、应变率为50 s⁻¹的情况下，试样纵向应力σ_x沿试样横向的分布情况。由此可见，沿试样横向，σ_x和σ_y表现出明显的非均匀分布特征：中心应力最高，试样边界处应力最低；其幅值随试样尺寸的增大而减小。由于三向静载的存在，图9(c)中的σ_x分布与图9(a)中的 σ_x分布存在一定差异；与图8中基于实验测试的理论计算结果相比，σ_x的数值计算结果要高一些，实际混凝土试样的非均匀性和数值模拟的材料均匀化处理可能是导致这一现象的原因。

图  9  σ_x和σ_y的分布

Figure  9.  Distributions of σ_x and σ_y

下载: 全尺寸图片 幻灯片

3.3   应力三轴度

为了探讨横向惯性带来的强度提升效果，采用下式定义参数ξ：

以对x轴冲击方向最大应力σ_x和等效应力${\sigma _{\text{e}}}{\text{ = }}\sqrt {3{J_2}}$进行比较分析。

（1）试样初始静载为三向等压，即σ_xs=σ_ys=σ_zs，此时x方向冲击产生的动态扰动σ_yd=σ_zd，σ_y=σ_z，由式(17)则有：

特别地，初始静载为零时，$\xi = {{{\sigma _{y{\text{d}}}}} \mathord{\left/ {\vphantom {{{\sigma _{y{\text{d}}}}} {{\sigma _{x{\text{d}}}}}}} \right. } {{\sigma _{x{\text{d}}}}}}$，为y方向与x方向动态扰动载荷之比。在平板撞击实验对应的一维应变条件下，${\sigma _{y{\text{d}}}} = {\sigma _{z{\text{d}}}} = \dfrac{\nu }{{1 - \nu }}{\sigma _{x{\text{d}}}}$。当泊松比v从0.15变化到0.25时，σ_yd/σ_xd由17.6%变化到33.3%，与应变率无关。在本实验中为测试此比值，将侧向杆与试件紧密接触，但不施加静载，可得到应变率50 s⁻¹下ξ=6%和应变率100 s⁻¹下ξ=14.5%^[10]，如图10(a)所示。这比基于Grote等^[11]在应变率10⁴ s⁻¹下的实验结果估计的ξ=52%低很多。其原因在于：相对于平板撞击一维应变条件，本文实验中应变率低很多，而且试样侧向也具有较好的波传播性能。

图  10  ξ的应变率效应

Figure  10.  Strain rate effect of ξ

下载: 全尺寸图片 幻灯片

（2）试样初始静载三向不相等，即σ_xs≠σ_ys≠σ_zs，此时x方向冲击产生的动态扰动σ_yd≠σ_zd，令：

由式(18)则有：

模型计算的ξ如图10(a)所示，表现出强烈的应变率效应：在低应变率下，ξ的计算结果与实验结果一致；但是应变率较高时，实际试样的非均匀性使得实验结果比数值计算结果低很多，需要进一步改进数值分析模型。图10(b)为0、10、20 MPa等3种静围压作用下的计算结果，有明显的静水压和应变率效应。同时，式(15)、式(17)和图10表明：此比值ξ在某种程度上可以反映x方向冲击产生的横向惯性效应。另外，对边长分别为20、37、50、74、100 mm等5种尺寸的正方体试样进行计算，其结果如图10(c)所示，计算得到的ξ有明显的尺寸效应。

进一步，在如图11所示的主应力空间引入应力三轴度η的定义^[25]，即：

图  11  主应力空间

Figure  11.  The space of principal stresses

下载: 全尺寸图片 幻灯片

式中：$p = - {\sigma _{\text{m}}} = - {{{I_{\text{1}}}}/ 3}$，$q = \bar \sigma = \sqrt {3{J_2}}$。

根据文献[10]中表2中的实验数据和本文的数值计算结果，整理得到ξ与应力三轴度η的关系如图12所示。由此可见：不同应变率情况下的数据点都集中在一条曲线上，具有较好的一致性，ξ与η的关系表现出应变率不敏感特性。此关系为混凝土材料应变率效应和尺寸效应的评估提供了一种新的思路。值得注意的是，图12中圈出的部分双向侧限的数据点，有些偏离数据点集中的曲线。原因在于，其初始剪切应力较大，虽然没有达到混凝土材料的剪切强度，但已经足以产生局部剪切破坏。

图  12  ξ与η的关系

Figure  12.  Relationship between ξ and stress triaxiality η

下载: 全尺寸图片 幻灯片

4.   结　论

基于真三轴应力状态作用下混凝土材料的冲击性能的实验结果，进行了Holmquist-Johnson-Cook (HJC)模型和Drucker-Prager (DP)模型数值模拟等方法的对比研究，得到的主要结论如下。

（1）基于考虑应变率效应的HJC模型，数值计算得到的混凝土试件动力学性能具有明显的静水压效应；基于考虑静水压效应的DP模型，数值计算得到的混凝土试件动力学性能具有明显的应变率效应。因此可得出，混凝土材料的应变率效应和静水压效应具有较强的耦合作用。

（2）基于受约束试样横向位移和应力的应变率表达式，结合基于HJC模型的数值计算方法，探讨了混凝土试样的横向惯性效应和尺寸效应，由此提出了一种对x冲击方向最大应力σ_x和等效应力σ_e进行比较分析的参数ξ。参数ξ有明显的静水压效应、应变率效应和尺寸效应。值得注意的是，ξ与应力三轴度η的关系曲线表现出应变率不敏感特性，此现象可为混凝土材料的应变率效应和尺寸效应的研究提供一种新的方法。