三角形蜂窝的共面冲击动力学

孙德强， 宫 凯， 李国志， 李剑玲， 邢月卿， 巩桂芬， 赵郁聪

(陕西科技大学 设计与艺术学院， 陕西 西安 710021)

0 引言

作为一种二维缓冲多孔材料，三角形蜂窝在航空、包装、军事、建筑、交通等领域内广泛应用.结构如图1所示，其周期性特征单元包含两个长为l的斜边.壁厚为t，孔深为b，底边与斜边的夹角θ称为扩展角.t/l称为壁厚边长比.假设ρ*和ρs分别为三角形蜂窝及其壁材的密度，则三角形蜂窝的相对密度为[1]：

(1)

关于二维多孔材料静动态力学行为已有大量的理论、试验和有限元方面的研究.Hohe 等人研究了相关弹性力学性能[1, 2].Grenestedt给出二维多孔材料的弹性模型[3].文献[4]和[5]系统总结了二维多孔材料静力学性能.Zhu和Mills对正六边形蜂窝静压性能进行了理论和实验研究[6].为了克服实验样品尺寸的局限性，有限元法应用越来越普遍.Papka和Kyriades利用试验和有限元法对比研究了六边形和圆形蜂窝在共面单双轴压缩载荷下的响应[7-10].Hönig和Stronge探讨了六边形蜂窝共面坍塌带的产生和冲击波的扩展规律[11, 12].Ruan等人采用有限元法研究了壁厚和冲击速度对正六边形蜂窝共面变形模式和峰应力的影响[13].后来Zheng等人[14]、Li等人[15]、Ali等人[16]和Sun等人[17, 18]借助有限元法研究了各种二维多孔材料的共面变形模式、峰应力和能量吸收.Hu和Yu推导了正六边形蜂窝共面峰应力的公式，并用有限元计算结果进行验证[19].上述研究证明，微单元结构参数和冲击速度决定着二维多孔材料的共面宏观力学性能.然而这些研究将重点放在了六边形蜂窝材料上.后来Liu和Zhang利用有限元法模拟了三角形和四边形蜂窝的共面冲击响应，并指出特征单元排列方式对其能量吸收的影响[20].然而，他们使用的蜂窝单元是正三角形的，而且关于结构参数和冲击速度对其共面动态力学性能的影响并没有深入研究.

图1 常用三角形蜂窝的三维结构示意图

本论文将借助于软件ANSYS/LS-DYNA来模拟研究三角形蜂窝在3～250 m/s的冲击速度下的共面(如图1的x1方向)冲击性能.不同结构参数下，三角形蜂窝的变形模式和动态峰应力将以图和经验公式的形式来表示，以此揭露结构参数和冲击速度对其共面冲击力学性能的影响.

1 共面冲击分析的有限元方法

1.1 有限元模型

图2 三角形蜂窝共面冲击分析有限元模型

图3 典型的三角形蜂窝共面冲击响应曲线 (l=3 mm， t=0.25 mm, θ=60°， v=70 m/s)

类似于Ruan[13]、Zheng[14]、Li[15]、Ali[16]和Liu与Zhang[20]等人的模型，三角形蜂窝共面冲击分析有限元模型示于图2.样品夹在两个钢性板间，异面方向(图1的x3方向)上的位移自由度限制为0以模拟共面应变状态.支撑板P2固定，冲击板P1以恒速v向下匀速冲击至样品密实，v处于3～250 m/s.采用ANSYS/LS-DYNA软件建模.使用5个积分点的Belytschko-Tsay壳单元Shell163进行网格划分.整个模型定义无摩擦的self-contact型接触，样品与两板间定义摩擦系数为0.02的surface-to-surface型接触.如Papka和Kyriakides[7]所用材料，壁材为双线性硬化材料[21, 22].如Hönig与Stronge[11, 12]、Ruan[13]、Zheng[14]、Li[15]、Ali[16]、Liu与Zhang[20]和Tan[23]等人的处理，壁材视为应变率不敏感的.对于所有样品，l=3 mm.样品中x1和x2上单元数量应在8×10以上，孔深b≥8 mm.

1.2 分析方法

采用后处理软件LSPREPOSTD对计算结果做进一步处理，可以得到相关响应曲线，如图3所示.u为冲击板位移，F为样品与冲击板间的冲击力，BF为样品与支撑板间的所用力，K、U、E分别为样品动能、内能和总能量吸收，W是冲击力F的功.则E=K+U，W表示为

(2)

计算结果表明，三角形蜂窝共面冲击响应遵循二维多孔材料的基本响应规则.如图3(a)，线弹性阶段是一瞬态过程，冲击力快速增加至初始峰冲击力F0，相应小的初始位移为u0.因屈服使冲击力突降，紧随其后是一长的平台变形，冲击力绕一水平上下波动.相应地在W-u曲线上对应一条近乎直的斜线，而后曲线突然上拐密实化开始.当大部分相邻单元壁相互接触时，密实化开始.对应W-u曲线上关键点的位移为uD，称为“密实化位移”.u0和uD对应的功分别为W0和WD.定义(WD-W0)/(uD-u0)为动态峰冲击力FP.设样品沿冲击方向上横截面的面积为A0，定义三角形蜂窝的共面动态峰应力σp*为

(3)

2 结果分析

对正三角形蜂窝样品(l=4.5 mm，t=0.5 mm，θ=60°)，在冲击速度分别为10、35和70 m/s下采用上述有限元模型进行计算.我们也得到了类似于Liu 和Zhang[20]给出的变形模式，同时得到的相应动态峰应力与此模拟结果相吻合.证明了本有限元模型的可靠性.

2.1 共面变形模式

2.1.1 冲击速度的影响

三角形蜂窝变形模式不单单有Liu 和Zhang[20]所述的，完整三角形蜂窝的变形模式可以分3类：均匀模式、过渡模式和动态模式，如图4～6所示.相应F和BF曲线如图7所示.

低速冲击时的变形为均匀模式(图4).均匀变形过后，一个随机倾斜的局部坍塌带出现在冲击板或支撑板端；随着初始坍塌带的逐渐密实，更多临近的单元逐步加入该坍塌带，直至整个蜂窝密实.随着冲击速度的增加，发生靠近冲击板和支撑板端先后出现“一”字形坍塌的过渡模式(图5)，Liu 和Zhang[20]所述变形模式就属此种.高速冲击条件下出现动态模式(图6)，靠近支撑板的坍塌带消失，“一”字形坍塌带仅存于冲击板端.

图4 三角形蜂窝共面冲击均匀变形模式(l=3 mm， t=0.15 mm, θ=70°， v=3 m/s)

图5 三角形蜂窝共面冲击过渡变形模式(l=3 mm， t=0.15 mm, θ=70°， v=50 m/s)

图6 三角形蜂窝共面冲击动态变形模式(l=3 mm， t=0.15 mm, θ=70°， v=250 m/s)

与Liu 和Zhang[20]的有限元模型相比，本模型中单元层数少，这仅导致了应力曲线上峰值点个数的减少，变形模式的本质是一样的.

由图7可见，支撑板上的应力应变曲线独立于冲击速度，而冲击板上的应力应变曲线是速度敏感的.均匀模式下，冲击板和支撑板上的冲击应力应变曲线几乎是完全重合的(图7(a))，表明此时冲击过程中的力平衡状态.随着冲击速度的增加，冲击板上的应力在平台区段的幅值和振幅越来越大(图7(b)和(c)).同时，支撑板上应力开始点推迟越来越大.之所以出现这样的规律，是因为支撑板上的应力随着冲击波的反射和相互作用而不断增强.低速冲击时冲击波有足够的时间在冲击板和支撑板之间反复反射而达到平衡；然而高速冲击时样品靠近冲击板端在相对较短的时间内很快发生密实，甚至冲击波还未到达支撑板并反射时就已发生.故随着冲击速度的增加，支撑板上由冲击波所引起的应力开始时间越来越滞后.

图7 三角形蜂窝共面冲击应力应变曲线(l=3 mm， t=0.07 mm, θ=60 ° (a)均匀模式;(b)过渡模式；(c)动态模式

2.1.2 结构参数的影响

选6组具有不同结构参数的样品，前三组用来研究扩展角的影响，壁厚分别固定为0.07、0.15和0.3 mm，每组8个样品的扩展角分别为10 °、20 °、30 °、40 °、50 °、60 °、70 °和80 °.后三组用来研究壁厚边长比的影响，扩展角分别固定为20 °、40 °和60 °，每组8个样品的壁厚分别为0.03、0.07、0.1、0.15、0.2、0.25、0.3和0.4 mm.

从不同结构参数和冲击速度下的变形模式来看，三种变形模式间的转换速度是由结构参数决定的，由均匀模式向过渡模式转化的冲击速度相对于由过渡模式向动态模式转换的冲击速度来说，关于结构参数不够敏感，一般处在15～35 m/s之间.而后者的转换速度随着扩展角和壁厚边长比的增加而明显增加.

2.2 共面动态峰应力

2.2.1 冲击速度的影响

Reid和Peng[24]给出了木材动态峰应力公式

σp*=σp0*+Av2

(4)

该公式被证明同样适用于规则蜂窝(Harrigan[25]、Ruan[13]、Ali[16]和Sun[17])、不规则蜂窝(Zheng[14]和Li[15])和泡沫(Tan[23]).该式同样也适用于过渡和动态模式下的三角形蜂窝.σp0*为相应共面静态峰应力，A为关于冲击速度的关系系数，定义为

A=ρ*/εD0

(5)

εD0为静态载荷下的密实化应变，A可以表示为

A=ρ*(1-ρ*/ρs)

(6)

从模拟结果很容易看出，在3 m/s的冲击速度下，所有样品密实化之前动能在总能量吸收中的比重都小于2%，所以此时的动态峰应力称为“类静态峰应力”，可视为σp0*[26].根据6组样品共面动态峰应力的计算结果拟合可得系数A的值，A关于t/l或θ的关系曲线见图8.由图8可见，随着t/l的增加，有限元结果拟合的A值相比式(6)的理论解变得越来越小，偏离理论值越来越大.当扩展角变小或单元壁厚变大时，相邻单元行间的距离变小，相对密度增加，此时单元行间较小的位移就会使单元发生密实，此时在固定冲击速度下密实化应变会变得越来越大.但是，在动态冲击载荷条件下密实化应变受冲击速度的影响更大.根据有限元计算结果发现，在v≤50 m/s时三角形蜂窝的密实化应变在0.7左右；在50 m/s

图8 三角形蜂窝共面动态峰应力关系系数A关于结构参数的拟合和理论关系曲线

样品编号θ(°)t/mmσp*(v=3 m/s)(MPa)σp*(v=20 m/s)(MPa)σp*(v=50 m/s)(MPa)σp*(v=70 m/s)(MPa)σp*(v=100 m/s)(MPa)σp*(v=150 m/s)(MPa)σp*(v=200 m/s)(MPa)σp*(v=250 m/s)(MPa)A(102 kg/m3)1100.0750.320.624.087.7918.6044.3659.8081.2313.442200.1441.882.564.417.1512.1323.3541.9761.829.643300.2004.835.488.069.6914.0427.0143.3461.479.154400.2417.659.5613.6715.3718.5028.5248.9366.259.235500.25910.6213.0916.9719.1825.4436.5747.8392.9511.836600.25011.7513.4618.7423.2626.7038.4166.7477.6910.807700.20712.2116.1421.9923.9230.0143.3558.90105.6013.468800.1264.979.5416.3822.0528.4541.1858.6477.2610.92

2.2.2 壁厚边长比的影响

由式(1)可见，三角形蜂窝相对密度随壁厚边长比的增加而增加.尽管A-t/l曲线不再遵循公式(6)，但由图8(a)、(c)和(e)可以看出A仍有随t/l增加而增加的趋势.A的增加意味着在给定冲击速度下有更高的共面动态峰应力.因此共面动态峰应力仍随t/l的增加而增加.Ruan等人[13]发现正六边形蜂窝的共面动态峰应力在给定冲击速度和其它结构参数保持不变的情况下与壁厚边长比成幂指数关系.根据4～6组样品的计算结果，即使每一组样品都有不同的扩展角，对于三角形蜂窝来说也出现同样的力学规律，关系式如下

σp*=σysB(t/l)k

(7)

B和k为关系系数，取决于三角形蜂窝的扩展角和冲击速度.相应4～6组样品，三角形蜂窝共面动态峰应力与t/l的幂指数关系曲线如图9所示.

对第四组样品来说，v分别为3、20、50、70、100、150、200和250 m/s时，由最小二乘法得到的关系式分别为

711(t/l)2、347(t/l)1.67、236(t/l)1.32、255(t/l)1.21、371(t/l)1.14、631(t/l)1.09、1 040(t/l)1.07和

1 709(t/l)1.09；对第5组样品来说分别为

1 977(t/l)2.16、1 503(t/l)1.99、842(t/l)1.67、708(t/l)1.52、580(t/l)1.35、656(t/l)1.21、951(t/l)1.18和1 174(t/l)1.13；

对第6组样品来说分别为

3 399(t/l)2.25、2 803(t/l)2.09、1 789(t/l)1.81、

1 359(t/l)1.64、1 146(t/l)1.5、1 007(t/l)1.31、

1 532(t/l)1.29和1 665(t/l)1.22.随着冲击速度的增加关系曲线逐步变成近似直线关系.

2.2.3 扩展角的影响

如图8(b)、(d)和(f)所示，不考虑t/l的影响，A与θ的关系仍旧满足式(6).由于扩展角趋于0和90°时，三角形蜂窝的相对密度急剧增加，所以在给定冲击速度下此时动态峰应力会迅速增加.由1～3组样品的计算结果可以看出，只要t/l固定，共面动态峰应力与扩展角都满足一致的关系规律.那就是，均匀变形模式下三角形蜂窝共面动态峰应力与扩展角的关系可用二次曲线来拟合

σp*=σys(Aθ(θ)2+Bθ(θ)+Cθ)

(8)

其中，Aθ、Bθ和Cθ为关系系数，取决于t/l和v.当三角形蜂窝处于过渡和动态变形模式时，共面动态峰应力与扩展角的关系可用如下关系式拟合

σp*=σysk1/(cosθ(k2+sinθ))

(9)

其中，k1和k2为关系系数，取决于t/l和v.

按照该关系规律对计算结果利用最小二乘法进行拟合，得到1～3组样品在不同冲击速度下共面动态峰应力与扩展角的关系曲线，如图10所示.同时还发现，当三角形蜂窝处于过渡或动态变形模式时，每一条曲线都有一个最小值点，相应最小值点的扩展角处于20 °～50 °之间.当v≤50 m/s时，相应扩展角约为20°左右；当50 m/s

图9 三角形蜂窝共面动态峰应力在不同冲击速度下与壁厚边长比的关系曲线图10 三角形蜂窝共面动态峰应力在不同冲击速度下与扩展角的关系曲线

3 结论

本论文采用有限元软件ANSYS/LS-DYNA，研究了t/l、θ和v对三角形蜂窝共面变形模式和动态峰应力的影响规律.

三角形蜂窝在共面冲击载荷作用下，会表现出3种变形模式：均匀、过渡和动态.由均匀模式转向过渡模式的转换速度相对于由过渡模式转向动态模式的转换速度来说，对结构参数不够敏感.而由过渡模式转向动态模式的转换速度随着扩展角和壁厚边长比的增加而增加.

计算结果显示，当三角形蜂窝单元结构参数固定时，共面动态峰应力与冲击速度的平方呈线性关系.在给定冲击速度下，共面动态峰应力与壁厚边长比成幂指数关系，与扩展角也可以用较复杂的关系曲线进行拟合.随着冲击速度的增加动态峰应力与壁厚边长比的关系曲线会由幂指数关系变成近似直线关系.当三角形蜂窝处于过渡和动态变形模式时，动态峰应力与扩展角关系曲线上会出现极小值点，该极小值点相对应的扩展角会随着冲击速度的增加而增加.

[1] Hohe J,Becker W.Effective elastic properties of triangular grid structures[J].Comp. Struct.,1999,45(2):131-145.

[2] Hohe J,Beschorner C,Becker W.Effective elastic properties of hexagonal and quadrilateral grid structures[J].Comp. Struct.,1999,46(1):73-89.

[3] Grenestedt JL.Effective elastic behavior of some models for ‘perfect’ cellular solids[J].Int. J. Solids Struct., 1999,36(10):1 471-1 501.

[4] Gibson LJ,Ashby MF.Cellular solids: structures and properties, 2nd ed.[M].Cambridge: Cambridge University Press,1997.

[5] Wang AJ,McDowell DL.In-plane stiffness and yield strength of periodic metal honeycombs[J].J. Eng. Mater. Technol.,2004,126(2):137-156.

[6] Zhu HX, Mills NJ.The in-plane non-linear compression of regular honeycombs[J].Int. J. Solids Struct.,2000,37(13):1 931-1 949.

[7] Papka SD,Kyriakides S.In-plane compressive response and crushing of honeycombs[J]. J. Mech. Phys.Solid,1994, 42(10):1 499-1 532.

[8] Papka SD,Kyriakides S.Experiments and full-scale numerical simulations of in-plane crushing of a honeycombs[J].Acta. Mater., 1998,46(8):2 765-2 776.

[9] Papka SD,Kyriakides S.Biaxial crushing of honeycombs. Part I: Experiments[J].Int. J. Solids Struct.,1998,36(29):4 367-4 396.

[10] Papka SD,Kyriakides S.In-plane biaxial crushing of honeycombs. Part II: Analysis[J].Int. J. Solids Struct,1998,36(29):4 397 -4 423.

[11] Hönig A,Stronge W J.In-plane dynamic crushing of honeycomb. Part I: crush band initiation and wave trapping[J].Int.J.Mech.Sci.,2002,44(8):1 665-1 696.

[12] Hönig A,Stronge WJ.In-plane dynamic crushing of honeycomb. Part II: application to impact[J].Int. J. Mech. Sci.,2002,44(8):1 697-1 714.

[13] Ruan D,Lua G,Wang B,et al.In-plane dynamic crushing of honeycombs-a finite element study[J].Int. J. Impact Eng.,2003,28(2):161-182.

[14] Zheng ZJ,Yu JL,Li JR.Dynamic crushing of 2D cellular structures-A finite element study[J].Int. J. Impact Eng.,2005,32(1-4):650-664.

[15] Li K,Gao X-L,Wang J.Dynamic crushing behavior of honeycomb structures with irregular cell shapes and non-uniform cell wall thickness[J].Int. J. Solids Struct.,2007,44(14-15):5 003-5 026.

[16] Ali M,Qamhiyah A,Flugrad D,et al.Theoretical and finite element study of a compact energy absorber[J].Adv. Eng. Soft.,2008,39(2):95-106.

[17] Sun D,Zhang W.Mean in-plane plateau stresses of hexagonal honeycomb cores under impact loadings[J].Comp. Struct.,2009,91(2):168-185.

[18] Sun D,Zhang W,Wei Y.Mean out-of-plane dynamic plateau stresses of hexagonal honeycomb cores under impact loadings[J].Comp Struct,2010,92(11):2 609-2 621.

[19] Hu LL, Yu TX.Dynamic crushing strength of hexagonal honeycombs[J].Int. J. Impact Eng.,2010,37(5):467-474.

[20] Liu Y, Zhang XC.The influence of cell micro-topology on the in-plane dynamic crushing of honeycombs[J].Int. J. Impact Eng.,2009,36(1):98-109.

[21] Gao X-L.An exact elasto-plastic solution for a thick-walled spherical shell of elastic linear-hardening material with finite deformations[J].Int. J. Pres. Ves. and Piping.,1994,57(1):45-56.

[22] Gao X-L.Strain gradient plasticity solution for an internally pressurized thick-walled cylinder of an elastic linear-hardening material[J].Z. angew Math. Phys.,2007,58(1):161-173.

[23] Tan PJ,Reid SR,Harrigan JJ,et al.Dynamic compressive strength properties of aluminum foams. Part II-‘shock’ theory and comparison with experimental data and numerical models[J].J. Mech. Phys. Solids,2005,53(10):2 206-2 230.

[24] Reid SR,Peng C.Dynamic uniaxial crushing of wood[J].Int. J. Impact Eng.,1997,9(5-6):531-570.

[25] Harrigan JJ,Reid SR, Peng C.Inertia effects in impact energy absorbing materials and structures[J].Int. J. Impact Eng.,1999,22(9):955-979.

[26] 王青春,范子杰.利用Ls-dyna计算结构准静态压溃的改进方法[J].力学与实践,2003,25(3):20-23.