1.1　复合结构几何构型及有限元模型

复合结构几何构型示意图如图1所示。该结构由2层厚度为的UHPC面板和蜂窝-钢隔板复合芯层组成，结构内部包含20个蜂窝块，分布于2层，每层10个。单个蜂窝块的长×宽×高为150 mm×150 mm×150 mm的立方体，蜂窝压缩强度为，蜂窝的面外方向与y轴方向一致。水平与垂直隔板的厚度均为。UHPC面板主要发挥耐腐蚀、耐磨损作用，面板与芯层之间未设置连接结构。芯层中，蜂窝为主要吸能组件，隔板用于分隔蜂窝块，水平隔板与蜂窝胶接以提高协同变形性能。

图1全屏查看

复合结构几何构型

基于LS-DYNA建立的有限元模型如图2所示，其中，锤头、UHPC面板、蜂窝采用SOLID单元建模，隔板采用SHELL单元进行建模。水平与纵向隔板以共节点方式连接，水平隔板与蜂窝通过TIE_SURFACE_TO_SURFACE模拟胶接，其余各部件间的接触通过定义CONTACT_AUTOMATIC_SURFACE_TO_SURFACE关键字进行模拟。设置动摩擦因数为0.15，静摩擦因数为0.20^[3]，结构底部设置固定约束。采用精细的均匀网格划分各部件，以更好地模拟冲击下的大变形，UHPC面板网格尺寸为4 mm(沿厚度方向设置4层)，蜂窝与格栅网格尺寸为5 mm。荷载设置方法如下：以MASS_PART关键字定义锤头总质量为，并施加竖直向下的速度。锤头以刚性材料模拟，同时考虑冲击过程中重力的影响。

图2全屏查看

复合结构有限元模型

采用CSCM本构模型来描述混凝土的屈服和流动，屈服面本构方程如下：

(1)

其中：为应力张量第一不变量；、分别为偏应力张量第二和第三不变量；为盖帽硬化系数^[16]；为Rubin折减系数；和分别为剪切破坏参数和盖帽硬化参数，屈服面本构方程中各参数的设置方法见文献[16]。UHPC材料的弹性模量、抗压强度、抗拉强度如表1所示，使用直径为75 mm的霍普金森压杆测试UHPC动态应变率，结果如表2所示。钢纤维质量分数为2%，长度为 12 mm，直径为0.25 mm，抗拉强度为2 800 MPa，采用蒸汽养护7 d。

隔板材料为Q235，密度为7 850 kg/m³，弹性模量为210 GPa，屈服强度为235 MPa，泊松比为0.3。通过Cowper-Symonds方程来描述其黏塑性行为：

(2)

其中：为屈服强度；为初始屈服应力；为应变率；为硬化参数；为等效塑性应变；为塑性硬化模量。

1.2　蜂窝实体等效本构模型

1.2.1　蜂窝实体等效屈服面

蜂窝作为一种周期性薄壁结构，胞元边长通常为毫米级^[18]。对于大尺度的场景，采用壳单元精细化建模会造成网格数激增，计算效率很低。LS-DYNA显动力学软件提供了一种可用实体单元等效模拟蜂窝各向异性行为的MAT_126材料模型，能有效减小计算时间。

目前，应用该材料模型集中于单轴压缩工况的研究^[19]，针对船桥碰撞这类涉及蜂窝压剪混合受载的局部冲击问题，材料本构的设置方法与适用性有待进一步研究。为此，首先设置实体蜂窝的密度，蜂窝的等效密度可据下式计算^[18]：

(3)

其中：为蜂窝胞元壁厚；为蜂窝胞元边长；为材料密度。

该模型可采用3个屈服面来表达，其中，第一屈服面需要分别定义所有法向应力和剪切应力，可通过试验测得，这些参数在求解过程中是完全解耦的。

在MAT_126本构中，激活第二屈服面需定义蜂窝等效屈服强度、强轴硬化应力和弱轴硬化应力，该屈服面在横向上呈各向同性，可考虑加载角度的影响。用于定义屈服应力与加载角度之间的关系，和用于表征蜂窝密实后的应力强化效果，单轴极限应力由下式计算得到：

(4)

其中：为加载角度；为结构的体积应变。为使极限应力随加载角度减小，需确保

(5)(6)

压缩强度为1.1 MPa的蜂窝在面外方向单轴加载下的应力-应变曲线如图3所示。采用预压工艺减小了初始峰值力。单个蜂窝块的长×宽×高为150 mm×150 mm×150 mm，蜂窝胞元边长为4 mm，壁厚为0.062 mm。由应力-应变曲线的平台段确定，极限应力随角度的退化关系见文献[18]。与由应力-应变曲线的密实段确定，工程应变与对数应变的转换关系为

(7)

其中：为对数应变；为工程应变。

图3全屏查看

蜂窝力学性能测试结果

需分别定义LCAB、LCBC、LCCA这3条损伤曲线，以表征蜂窝的各向剪切行为。曲线的自变量为应变，因变量为剪切应力系数。在模拟中，应力会在每个时间步长随损伤值逐步更新。未损伤时，剪切应力系数记作1，损伤后，剪切应力系数定义为(0，1)区间的值。可根据图3中剪切曲线确定。

在激活第二屈服面时，为避免结构由于剪切抗性低、在小应力作用下出现严重的剪切畸变等情况，需要引入第三屈服面，在不影响单轴性能的情况下，进一步规定材料的静水强度和剪切强度。在单轴应力作用下，极限应力由式(4)计算得到；在压力状态下，极限应力由式(8)计算得到；在剪切作用下，极限应力由式(9)计算得到。

(8)(9)

1.2.2　蜂窝实体等效方法验证

设计如图4和图5所示的蜂窝压缩试验，用于验证实体等效模型的模拟结果，并对激活和未激活第二、第三屈服面两种情况均进行对比分析。在单元类型的选择上，参考LS-DYNA手册，选取全积分单元和非线性弹簧单元两种。试验基于INSTRON 1342试验机开展，分别用于验证蜂窝的单轴压缩与局部受载下的压剪混合行为。平面锤头是为长×宽150 mm×150 mm的方形，弧面锤头是半径为75 mm的半圆构型。在2组试验中，锤头均以5 mm/min匀速压缩。蜂窝块结构与材料参数与1.2.1节中的一致。

图4全屏查看

平面荷载下变形模式对比

图5全屏查看

弧面荷载下变形模式对比

从图4可见：在单轴压缩下，本构屈服面模型和单元类型均能得到较好的模拟结果；在弧面锤头荷载下，蜂窝不仅受到沿轴向传递的力，而且受到剪切作用和面内的侧向荷载，在这种复杂应力耦合作用下，采用非线性弹簧单元的蜂窝无法在仿真中反映出合理剪切抗力，产生了明显的单元畸变；采用不激活第二和第三屈服面的设置方法，蜂窝在非轴向荷载下出现了整体向外扩张的现象。上述结果均与试验中的变形结果不符。总体而言，全积分单元与激活第二和第三屈服面结合的建模方法具有最优适应性，无论是在单轴荷载下还是在弧面复合荷载下，所得结果均与试验结果高度吻合。

应力-应变曲线试验结果与数值模拟结果如 图6所示。为确保可比性，硬件设备以及所采用的CPU线程数均保持一致。结果表明，采用实体等效模拟方法所得应力曲线不仅与试验的应力曲线吻合度较高，而且可大幅度减小计算时间，相比于壳单元精细建模，计算时间减小了72.75%。在本文后续研究中，对蜂窝均按该方法进行实体建模。

图6全屏查看

不同数值模拟方法应力曲线与计算效率对比

1.3　复合结构落锤冲击试验与验证

基于中南大学高能级落锤试验平台，开展UHPC-蜂窝复合结构冲击试验。试验实物图如图 7(a)所示，锤头呈半圆弧构型，半径为375 mm，质量为1.6 t。试验时，锤头距试样上表面抬升高度为80 mm，对应撞击初速度为3.5 m/s。试样总长×宽×高为750 mm×300 mm×332 mm，试样上下UHPC板厚度为16 mm，蜂窝面外压缩强度为 0.5 MPa，格栅厚度为0.2 mm。试样横向置于均力板上，均力板底部沿试样对称设置2个200 t量程力传感器，通过高速摄影装置采集复合结构冲击变形响应模式。

图7全屏查看

UHPC-蜂窝复合结构落锤冲击试验与有限元验证

基于LS-DYNA开展显式动力学仿真计算，荷载与边界条件定义方法与1.1节中的相同，设置冲击初速度 m/s，锤头质量 t。对于UHPC面板，CSCM本构采用强度包络和损伤累积机制模拟破坏行为；对于钢材，则定义失效应变为0.2。在计算过程中，不采用仿真软件中的质量缩放算法，实际最小计算时间步长为2.9×10^-7 s。试验与仿真变形模式见图7(b)与图7(c)。从7(b)与图7(c)可见：在锤头冲击下，UHPC顶部呈现出弧形的大挠度屈曲，出现微裂纹但未发生完全断裂，试验与仿真所得的面板整体变形模式相吻合。根据图7，位置1与1′处的格栅产生屈曲，并向蜂窝面内方向挤压。相比于试验结果，在仿真中格栅屈曲的波长更小。位于冲击中心位置的蜂窝块其变形以面外方向轴压为主，位于中心两侧的蜂窝逐渐演变为压剪混合受载模式，如位置2所示；而位于最外侧的蜂窝则在格栅挤压以及斜向荷载的耦合作用下，出现屈曲行为，以面内变形为主(位置3)。

试验与仿真所得的冲击力-时间曲线对比见 图8。从图8可见：在冲击初期出现1个初始峰值；由于格栅的屈曲，曲线呈现振荡形状；弧面构型锤头在冲击过程中与样件的接触面积不断增大，使得撞击力随之提升，并逐渐达到最大峰值；最后，由于样件的速度减小，锤头回弹并迅速卸载。结果表明：冲击力试验曲线与仿真曲线高度吻合，其中，撞击初始峰值力相对误差为7.8%，最大峰值力的相对误差为6.4%，撞击响应时间的相对误差为7.0%。

图8全屏查看

UHPC-蜂窝复合结构落锤冲击力-时间曲线对比

2.1　冲击响应与能量分配关系

在15 mm， MPa， m/s，mm， t下，UHPC-蜂窝复合结构在冲击荷载作用下的动态响应过程见图9。其中为冲击响应时间。从图9可见在局部弧面荷载的作用下，UHPC面板出现屈曲，并且挠度随冲击位移提高不断增加；同时，芯层的变形特征更加多样，体现在：在冲击的初期阶段(=0.02 s)，芯层的变形由中心蜂窝的轴向压缩主导，同时，伴随着隔板的初个屈曲折叠形成，此时段的变形局限于一个相对较小的区域；随着冲击荷载进一步作用，塑性变形区域开始向外扩散(=0.04 s)，位于中心两侧的蜂窝结构也开始参与能量吸收。该位置的蜂窝处于复杂应力状态，表现出多种变形模式，在其顶部与UHPC面板的接触区域，蜂窝呈压剪混合变形模式，而在中部，由于受到格栅屈服的挤压作用，蜂窝产生面内局部凹痕；在0.06 s时，中心蜂窝依然表现为轴压主导，而中部两侧的蜂窝因侧向力作用不断加强出现整体屈曲，并进一步向最外侧蜂窝传递面内荷载。这一现象揭示了结构响应在时序和空间上的依赖性，表明结构的不同部分对整体耐撞性的影响各不相同。

图9全屏查看

冲击响应过程

图10(a)所示为整个冲击响应历程的力-位移曲线，可反映响应模式的时空依赖性以及各个部件在碰撞吸能中的贡献。从图10(a)可见：在冲击响应第一阶段，由于锤头与复合结构的接触区域较小，变形集中在中心附近；冲击力随隔板的屈曲而呈现波形变化；随着冲击位移逐步增加，变形区域不断扩大，与之相应的冲击力也逐渐升高；同时，曲线在第二阶段末期出现最大峰值，反映了结构的吸能效率最高。

图10全屏查看

复合结构抗冲击特性及时空特征

不同组件的吸能能力与时序性如图10(b)所示。从图10(b)可见UHPC的吸能主要集中于第一阶段。图10(c)中的挠度-位移曲线和图10(d)中的吸能量-位移图进一步阐释了这一现象，其中，x为锤头位移， 代表UHPC沿跨度方向的位置坐标，当锤头位移为30 mm时，UHPC面板的挠度增量最大，且在该时刻面板的吸能占全过程吸能的61%，表明面板的塑性变形与吸能集中于冲击响应的前期和中期。

蜂窝与隔板在整个冲击历程中的吸能量均显著高于UHPC，且第二阶段的吸能量大于第一阶段，说明塑性区域的扩展强化了芯层的吸能能力。同时，蜂窝作为复合结构主要吸能单元，对总吸能量的贡献率达52.36%。蜂窝与隔板的吸能贡献与空间位置的关系见图10(e)。从图10(e)可见：蜂窝被划分为5列，从左至右用自然数标记，位于蜂窝之间的隔板相应划分为4部分；芯层的吸能量随空间位置呈现正态分布。由于存在最高的压缩量以及面外的单轴压缩主导变形模式，仅中心一列的蜂窝吸能量贡献了全蜂窝46.37%的吸能量；中心两侧的两列蜂窝处于压剪混合与面内凹陷的复杂受载状态，贡献了40.57%的吸能量。类似地，处于冲击中心位置的隔板产生了更多的塑性吸能量，占比为73.68%。表明在局部受载情况下，塑性波由撞击中心向外传播，芯层的吸能效果由中心向两侧递减。

2.2　结构参数-冲击特性关联性研究

2.2.1　Spearman秩相关系数分析

选取UHPC厚度(10~20 mm)、蜂窝壁厚(0.04~0.06 mm)与隔板厚度(0.5~1.0 mm)作为关键结构参数开展研究，蜂窝边长为4 mm。由于采用实体等效方法对蜂窝进行模拟，在仿真模型中直接定义强度，不需要建立精细的正六边形蜂窝胞元。蜂窝壁厚与蜂窝强度之间的关系式如下^[18]：

(10)

其中：为材料屈服强度，对应蜂窝强度区间为0.5~1.0 MPa。基于拉丁超立方采样设计20组工况，如表3所示。

表3全屏查看

拉丁超立方采样结果

序号	/mm	/MPa	mm	序号	mm	MPa	mm
1	12.50	0.94	0.64	11	13.29	0.65	0.62
2	17.40	0.54	0.66	12	12.40	0.51	0.88
3	16.74	0.69	0.54	13	18.49	0.78	0.84
4	10.26	0.6	0.86	14	18.94	0.99	0.91
5	14.56	0.73	0.94	15	14.60	0.58	0.58
6	16.36	0.76	0.81	16	14.01	0.57	0.75
7	10.83	0.64	0.50	17	15.72	0.96	0.57
8	11.58	0.85	0.98	18	19.02	0.70	0.95
9	17.77	0.83	0.78	19	13.73	0.91	0.72
10	19.97	0.81	0.70	20	11.47	0.88	0.69

展开更多

利用Spearman秩相关系数^[20]评估UHPC-蜂窝复合结构设计参数与其冲击性能指标之间的非线性关联性。关键冲击特性指标如下。

1) 初始峰值力。为冲击初期遭遇的瞬时高力冲击力，其在特定冲击场景中需要控制在一个较小的值。例如，在冲击物为柔性体时，过高的初始峰值力可能导致冲击物或支撑结构损伤。

2) 冲击响应过程中的最大峰值力。反映结构的最大吸能效率，但在某些场景下也应予以限制以防产生过大的结构损伤。

3) 撞击深度。用于衡量结构在吸收全部冲击动能时所需的位移，它直接反映了结构在冲击全程的综合吸能效率。

Spearman秩相关系数可通过以下计算式得到：

(11)

其中：通过对变量排序，并计算其中第i个数据点在排名上的差值得到；为数据点的数量。基于Python中的相关函数计算显著度水平(p)。p反映了在零假设下计算的相关系数或更极端情况下的概率，计算结果如表4所示。

表4全屏查看

Spearman秩相关系数分析结果

输入变量	输出变量		p
		0.660	0.001 53
		0.269	0.251 00
		-0.163	0.489 00
		0.624	0.003 27
		0.911	2.330 00×10-8
		-0.924	5.530 00×10-9
		0.592	0.005 91
		0.294	0.207 00
		-0.315	0.174 00

展开更多

从表4可见：、与均对初始峰值力具有较高的正相关影响力度，且p较小表明这种影响较显著。其中，与初始峰值力之间的相关性最强。蜂窝强度与最大峰值力和撞击深度具有最强的相关性，表明在冲击过程中，蜂窝对于承载和吸能起到了决定性作用。其较大的Spearman秩相关系数和较小的p表明，增加蜂窝强度能显著提高结构的吸能效率和抗变形能力。

2.2.2　结构几何参数-冲击特性相关性

为进一步直观展示结构参数对冲击特性的影响，采用高斯过程回归模型作为建模工具，建立多维非线性的代理模型。高斯过程作为一种贝叶斯非参数统计方法，它能够预测连续的输出变量，并且给出预测的不确定性。

基于代理模型进一步分析结构参数-冲击特性相关性，结构参数作为自变量输入，自变量区间与2.2.1节中的一致，为UHPC厚度(10~20 mm)、蜂窝壁厚(0.04~0.06 mm)与隔板厚度(0.5~1.0 mm)。冲击特性关键指标作为输出，预测结果如图11(a)~(c)所示。图11(a)~(c)与表4中的Spearman秩相关系数一致，初始峰值力和最大峰值力与结构参数之间存在正相关关系，而撞击深度则与结构参数呈负相关。随着UHPC厚度增加，在所示区间内，初始峰值力的增值较显著，但撞击深度和最大峰值力的变化较平缓。预测结果进一步佐证了2.1节中的结论，表明UHPC面板在冲击响应的前中期发挥关键作用。分析结果表明：增加UHPC厚度会显著提高初始峰值力，但由于受UHPC本身的吸能量限制，其对撞击深度的调节能力较弱。在整个冲击过程中，芯层具有主导性的吸能作用，其中，蜂窝结构相较于隔板结构提供了更高的吸能量，因此，对最大峰值力和撞深的影响也更显著。相对而言，芯层中隔板结构在影响初始峰值力方面具有更突出的效果。

图11全屏查看

关键结构参数-冲击特性指标关联性分析

3.1　参数优化方法

以UHPC厚度、蜂窝强度与隔板厚度这3个结构参数为优化对象，采用基于非支配排序的遗传算法(NSGA-II)^[21]搜索最优的结构设计参数。

为实现结构的耐撞性与冲击阻抗提升，提出如下优化策略：最小化撞击深度，表示结构具有更强吸能能力；最小化最大峰值力，减少冲击全历程中冲击物与支撑所承受的最大荷载；初始峰值力与最大峰值力的比值接近1，使最小化，以进一步减小冲击损伤。目标函数为：

(12)

(13)

对于每一代种群中的每个个体，计算它在各个目标函数上的排名和拥挤距离，然后，根据这些指标进行选择、交叉和变异。NSGA-II的选择过程可以表示为

(14)

其中：为当前种群的个体集合，为选择后的新种群个体集合，集合包含每个个体的目标函数值和决策变量；为sort函数；为基于非支配排序的个体排名；为拥挤距离指标。交叉操作和变异操作分别定义如下：

(15)

(16)

其中：为函数；为函数；为父代个体的解向量；为变异个体的解向量；为交叉概率，设置为0.7；变异概率设置为0.1，以防止算法过早收敛陷入局部最优解。在优化过程中，种群大小设定为200，进化代数为200代。

3.2　优化结果分析

在多目标优化过程中，可计算得到一系列支配解和非支配解^[21]。Pareto是前沿由中互不支配的解构成的集合。在目标函数空间中，前沿解是在一个目标上改进的同时不会在另一个目标上劣化的解，因此，这些解被认为是最优的。在前沿解集中，采用一种基于最短归一化距离的方法得到最佳的权衡的平衡解：

(17)

其中：为解和解之间的归一化距离；和分别为解和解的目标函数值；为目标总数。对于每个解，计算得到它到其他解的距离之和，作为它的满意度分数：

(18)

其中：为前沿中解的总数，满意度分数为解i到其他所有解的距离总和。得到Pareto前沿以及平衡最优解如图12所示，其中最优解用红色标记。

图12全屏查看

Pareto前沿及最优解

对所设计结构，优化得到的最优结构参数取值如下：UHPC面板厚度为12.40 mm，蜂窝强度为0.84 MPa，格栅壁厚为0.81 mm。计算上述最优解的数值模拟结果，并与优化预测结果进行对比，结果如表5所示。从表5可见：多目标优化结果与数值仿真结果具有较高的吻合度，且采用较小厚度的UHPC与较高强度的蜂窝进行匹配，能得到更优的冲击特性。以雷达图的形式对比优化前后的结构(优化前结构参数与2.1节中的一致)在关键指标上的表现，除、和外，进一步采用最大峰值力与平均压溃力的比值进行评估。为更直观体现各指标的优劣及相对关系，对结果进行归一化处理。在优化前后结构中，将效果更优的指标记为1，另一指标则相应转换为(0，1)区间内的值。如图13所示，对所设计结构，优化前后的与接近，但优化后最大峰值力下降8.32%，撞击深度减小23.51%。证明优化后的结构在不增加对冲击物与支撑结构损伤的前提下，减小了最大峰值力，同时显著提升了吸能能力。

表5全屏查看

冲击阻抗增强设计优化结果

方法	/mm	/MPa	/mm	/kN	/kN	/mm
相对误差/%	0	0	0	-4.89	+2.93	+4.50
优化预测	12.40	0.84	0.81	117.54	225.25	68.84
数值模拟	12.40	0.84	0.81	111.79	231.84	71.94

展开更多

图13全屏查看

优化前后关键指标对比