发泡聚乙烯最大加速度-静应力曲线快速获取方法研究

宋卫生，薛阳

发泡聚乙烯最大加速度-静应力曲线快速获取方法研究

宋卫生，薛阳

（河南牧业经济学院，郑州 450046）

研究快速、准确预测最大加速度-静应力曲线的方法。首先利用落锤冲击试验机获取了5个不同高度下，5种不同厚度的发泡聚乙烯的最大加速度-静应力曲线。在此基础上，分析对比文中3种不同的改进拟合法与已有的动应力与动能量多项式拟合法的区别。研究发现，当不区分高度的情况下，以最大加速度因子为函数值，以跌落高度、衬垫厚度、静应力为变量进行拟合时，其代表预测精度2的平均为0.835，相比动应力与动能量多项式拟合法的0.299 6要高。但曲线右侧的预测精度偏低。引入以静应力为变量的多项式作为修正因子后，2的平均值为0.934。预测精度有所提高，右侧的预测偏差减小，但仍存在。在区分高度的情况下，以带有修正因子的公式进行预测时，2的平均值为0.984，曲线向右侧预测偏差逐渐增大的现象明显改善。区分高度情况下，利用带修正因子的预测公式可以快速且较准确地预测最大加速度-静应力曲线，可以为冲击防护设计及相关软件的开发提供一定的帮助。

最大加速度-静应力曲线；应力能量法；预测精度；发泡聚乙烯；多项式拟合

缓冲包装设计国家标准中[1]，冲击防护设计方法有缓冲系数-最大应力曲线法和最大加速度静应力曲线法。缓冲系数-最大应力曲线设计法结果偏差较大，最大加速度-静应力曲线法设计结果精确度高。最大加速度-静应力曲线是在动态冲击条件下获取的，因此又称为动态缓冲特性曲线。近些年来，国内外学者对缓冲材料的动态缓冲特性持续进行了研究。金强维等[2]基于4种不同密度的聚丙烯材料，研究发现当密度和冲击能对最大接触力、最大位移、最大应变、比吸能的影响规律。美国Hewlett-Packard公司研发中心的Daum[3]提出了使用应力-能量法来快速预测泡沫塑料的缓冲曲线。研究发现，动应力与能量密度之间存在指数函数。Marcondes等[4]运用应力-能量法获得缓冲材料的动态缓冲特性曲线，并对其准确性进行了验证。Kendalyn等[5]通过研究，认为关于应力-能量方法预测加速度值的准确性如何，以及它是否在统计上可与ASTM D1596相当，仍然是一个问题。史岩等[6]的研究指出应力-能量法是以假定缓冲材料在受到冲击过程中没有能量损失为前提，且该法较易受人为因素的影响，故利用应力-能量法的动能量公式计算得到的能量不够精确。王金梅等[7]以EPS材料为研究对象，基于应力-能量法，分别采用多项式函数以及二参数、三参数指数函数进行拟合，发现多项式函数拟合精度最高，但是三者都存在一定的误差。陈满儒等[8]基于常见材料EPS在3种不同密度下的动态压缩实验数据，结合应力-能量法得到了缓冲材料在不同密度下的动态缓冲曲线以及应力-能量方程。李光等[9]利用发泡聚乙烯，采用多项式拟合方法拟合缓冲曲线，但由于预测方法没有考虑实际试验过程中的能量耗散，导致曲线部分区域误差较大。陈磊[10]通过动态缓冲对不同发泡倍数和不同尺寸大小的EPP泡沫材料进行了性能测试及分析研究，发现应力-能量法可用于EPP动态缓冲特性的确定。邢紫玉等[11]基于应力-能量法，通过实验最终获得了不同跌落高度、不同厚度下的3种不同夹芯结构的聚丙烯板缓冲曲线模型，为今后选用缓冲材料提供了理论依据。此外，EPE的动态冲击特性还受到材料密度的影响[12]。李超等[13]基于柱形空气衬垫缓冲材料，通过试验、数据分析和曲线拟合，分别得到了应力-能量方程和缓冲曲线，证实了应力-能量法同样可以在柱形空气衬垫缓冲材料中使用。文中仅研究在一定密度条件下EPE的动态冲击特性。

综上分析，应力-能量法在获取最大加速度-静应力曲线方面是一个快捷的方法，但曲线预测精度有待提高。本文通过不同方法的对比分析，研究提高最大加速度-静应力曲线预测精度的方法，为其在缓冲包装设计中的应用及相关优化设计、软件开发提供参考。

1 实验

1.1 材料与设备

文中所用材料为发泡聚乙烯，密度为25 kg/m3。其应力应变曲线如图1所示。由于包装缓冲材料冲击试验机构造原理的原因，无法精确地在较小的静应力条件下进行实验，进而无法得到较为完整的最大加速度-静应力曲线。Piatkowski等[14]研究并验证了有限元仿真在缓冲材料冲击实验模拟的可靠性与有效性。因此，在ANSYS LS-DYNA中采用有限元仿真的方法开展了落锤冲击实验，获得了最大加速度-静应力曲线。

图1 EPE的应力应变曲线

1.2 实验方法

1）通常取5个点来拟合1条曲线，但当拟合点增加时，将更加精确[15]。因此，文中至少取5种静应力及对应的最大加速度来拟合曲线。鉴于有限元分析的稳定性，每个静应力条件下，不再多次实验求平均值。

2）分别选择厚度为0.02、0.03、0.04、0.05和0.06 m的衬垫进行冲击。

3）分别选择0.31、0.51、0.66、0.81和0.97 m作为冲击高度。

1.3 实验结果

图2a～e分别为0.31、0.51、0.66、0.81、0.97 m高度下，缓冲材料的厚度分别为0.02、0.03、0.04、0.05、0.06 m的最大加速度-静应力曲线。通过实验获取的最大加速度-静应力曲线整体呈开口向上的“U”形抛物线，加速度值先快速降低后缓慢上升。当砝码质量较小时，缓冲材料的变形较小，冲击时间较短，因而加速度较大；随着砝码质量的增加，缓冲材料的变形逐渐增大，冲击时间变长，缓冲材料吸收的冲击能量也随之增多，加速度逐渐减小；当砝码质量增大到一定程度时，缓冲材料会达到密实化，变形空间变小，导致加速度再次增加。

图2 最大加速度-静应力曲线

2 几种曲线预测方法的对比分析

2.1 动应力与动能量多项式拟合法

王金梅等[7]已经通过实验证实Matthew Daum提出的指数法拟合精度低，而多项式函数拟合精度高。动应力与动能量之间的多项式函数最早由张波涛[16]提出，具体如式（1）所示。

式（1）中动能量实则是冲击过程中的缓冲材料单位体积的能量吸收。动应力和动能量的获取方法如下：

1）设置最大与最小动能量值。动能量计算式如式（2）所示。

最小动能量值意味着最小的静应力、最小的冲击高度以及最大的材料厚度；最大动能量值意味着最大的静应力、最大的冲击高度以及最小的材料厚度。

2）等分能量域。设置能量域步长，选择不同的动能量值进行实验。

4）重复步骤3，得到不同动能量值条件下的动应力均值。

在参数、、、确定的条件下，只需将衬垫厚度、跌落高度代入式（3）中，就可以获得不同静应力对应的加速度因子，进而可以获取不同跌落高度、不同缓冲材料厚度下的最大加速度-静应力曲线。

按照以上方法采用17.22 kN/m2为最小动能量值，172.25 kN/m2为最大动能量值，17.23 kN/m2为能量域步长。冲击高度分别选用0.97、0.66、0.31、0.81 m，对应的衬垫厚度分别采用0.07、0.04、0.05、0.03 m，实验后得到每个动能量对应的平均动应力如表1所示。

表1 动能量与平均动应力

Tab.1 Dynamic energy and average dynamic stress

使用式（1）拟合后得到参数、、、的值分别为111.684 8、−0.562 55、0.078 098、−0.000 2。将得到的、、、具体值代入式（3），按照上述方法即可得到预设定冲击高度和衬垫厚度条件下，不同静应力对应的最大加速度因子。

为了对不同的预测方法的精度进行比较，文中使用了实验获取的数据作为基础数据进行公式预测精度的对比分析，并采用判定系数2为评价指标。

以图2中各高度情况下，衬垫厚度分别为0.03、0.04、0.05 m的实验数据作为基础数据，以获取的预测数据作为对比数据，获得冲击高度为0.31、0.51、0.66、0.81、0.97 m时预测得到的2分别为−0.674、0.297、0.514、0.618和0.743。5个高度2的平均值为0.299 6。2越接近1，说明预测精度越高。以上得到各高度的2说明，随着高度的增加，预测精度逐渐提升，但是整体的预测精度都不是很高。在由动能量获取动应力过程中，需要先在不同的静应力、冲击高度和衬垫厚度条件下开展实验，获取最大加速度因子，然后求平均，得到平均最大加速度因子，最后与静应力的乘积即为平均动应力。分析其过程数据发现，不同实验参数组合条件下，获得的最大加速度因子有一定的偏差，这就决定了后期预测会有一定的误差，因此这个预测方法需要进一步改进。

2.2 改进的预测方法

为了减小误差，提高预测精度，数学模型延续多项式形式，而拟合数据则直接采用实验过程中的冲击高度、衬垫厚度、静应力与最大加速度因子。基础的拟合函数形式如式（3）所示，式中以最大加速度因子为函数值，以冲击高度、衬垫厚度、静应力为变量。

1）第1种改进方法。为了减少实验工作量，第1种改进方法的思路是以较薄和较厚的衬垫厚度所得到的数据进行拟合，并预测处于中间厚度的数据。以图2所有高度情况下衬垫厚度分别为0.02 m和0.06 m的数据进行拟合，得到的、、、的值分别为12 786.45、3.991 6、1.51E−05、−1.24E−11。利用式（3）对0.03、0.05和0.05 m衬垫厚度的数据进行预测，获得冲击高度为0.31、0.51、0.66、0.81、0.97 m时预测得到的2值分别为0.658、0.857、0.88、0.911和0.913。5个高度2的平均值为0.835。图3为0.31 m冲击高度情况下的实验与预测曲线对比，其他高度的情况类似，不再一一列出。可以发现该方法实验的工作量比动应力与动能量多项式拟合法的有一定增加，但仍在可接受范围内，不过预测精度有了明显提升。但从图3中预测曲线与实验曲线的对比可以发现，随着静应力的逐渐增加，误差值越来越大。因此式（3）中有必要引入一个以静应力为变量的多项式作为修正因子，用以提高预测精度。

图3 改进方法1下冲击高度为0.31 m时的最大加速度-静应力曲线对比

2）第2种改进方法。根据第1种改进方法的分析结果，以式（3）为基础，添加了一个修正因子，具体如式（4）所示。式中、为待拟合参数，其他符号含义与式（3）相同，不再一一解释。

同方法1，以图2所有高度情况下0.02 m和0.06 m衬垫厚度的数据进行拟合，来预测处于中间厚度的数据。得到的a、b、c、d、e、f的值分别为42 647、13.38、4.03E−05、0.297 1、8.42E−6，利用式（3）对0.03、0.05和0.05 m衬垫厚度的数据进行预测，获得冲击高度为0.31、0.51、0.66、0.81、0.97 m时预测得到的R2分别为0.942、0.991、0.972、0.915和0.873，5个高度R2的平均值为0.934。图4为0.31 m冲击高度情况下的实验与预测曲线对比图。很显然，预测效果比第1种改进方法高。不但平均R2比第1种改进方法高，而且曲线右侧的预测效果也明显好了很多，说明了修正因子的有效性。但是从图4仍然能看到右侧有一定的误差存在，而且较为明显。

分析原因，可以考虑是由于用于拟合的数据范围选择过大造成的。在运输包装的后期测试中，其测试高度是要按照一定的标准选择，比如上文提到的ISTA 3A标准，当然也有其他的标准。高度选择基本都是与包装件的质量有关，一个质量范围选择一个高度。ISTA 3A标准使用的高度值分别为0.31、0.51、0.66、0.81和0.97 m。也就是说要预测的曲线的高度并不是随机，只需要根据测试需要获取某些高度的最大加速度-静应力曲线即可。因此，根据这一思路，就可以缩小拟合所用的数据范围了。

3）第3种改进方法。根据第2种改进方法的分析结果，仍然采用式（4）进行拟合，不同之处是第3种改进方法的拟合数据和预测数据仅是针对某个高度，比如以0.31 m高度为例，选择该高度情况下0.02 m和0.06 m衬垫厚度的数据进行拟合，来预测0.31 m高度情况下0.03、0.04以及0.05 m衬垫厚度的数据。拟合得到、、、、、分别为2 385.25、0.965 7、4.68E−07、−2.4E−13、4.524 2、0.000 119 2。预测得到的2值为0.976。在0.51、0.66、0.81、0.97 m冲击高度时，进行同样的操作，拟合得到、、、、、均有变化。当高度为0.51 m时，、、、、、分别为42 588、14.91、2.22E−05、−6.82E−11、0.271 49、1.04E−05；当高度为0.66 m时，、、、、、分别为11 368.76、3.562 6、7.95E−06、−9.23E−12、1.140 7、3.45E−05；当高度为0.81 m时，、、、、、分别为−11 811、−3.598、−5.39E−6、1.16E−11、−1.094 8、−5.76E−05；当高度为0.97 m时，、、、、、分别为1 147、0.304 6、9.58E−7、−1.76E−12、12.7、0.000 487 8；而得到的2分别为0.989 7、0.979 9、0.982和0.994。5个高度2的平均值为0.984，相比第2种改进方法有了明显的提高。图5为在0.31、0.51、0.66、0.81、0.97 m这5种不同冲击高度时的实验与预测最大加速度-静应力曲线的对比图。从图5中可以看出预测效果均有明显提升，特别是曲线右侧的偏差明显减小，预测曲线与实验曲线的整体贴合度也比较高，说明了第3种改进方法的可行性。

图5 5个不同冲击高度的最大加速度-静应力曲线对比

3 结语

从以上改进过程的对比分析可以得出结论如下：

1）基于动应力和动能量的多项式拟合法虽然相对指数法在预测精度上有了一定的提升，但是由于在单个动能量的多次实验过程中，得到用于平均的动应力有一定的偏差。导致了最终预测精度不够高，代表预测精度的平均2为0.299 6。

2）以最大加速度因子为函数值，以跌落高度、衬垫厚度、静应力为变量的多项式拟合法虽然实验量有一定的增加，但是预测精度比基于动应力和动能量的多项式拟合法要高，平均2为0.835。存在的问题是曲线右侧的预测精度较差，因此需要在预测公式中添加以静应力为变量的多项式作为修正因子，结果表明预测精度有一定的提升，但是右侧预测偏差仍较明显。

3）结合测试标准中高度值的选取方法，缩小数据范围，单独使用某个高度的数据进行拟合，并用于预测该高度条件下，其他厚度的曲线。结果表明，2明显提高，平均值达到0.984，各个高度的预测效果均较好，曲线向右侧预测偏差逐渐增大的现象得到明显改善。

4）由于目前最常用的动态缓冲设计方法及跌落测试均是基于一定跌落高度的前提下开展的。因此，从设计人员角度来看，第3种改进方法更适合实际需要。使用该方法对式（4）进行拟合时，只需要某个高度的实验数据即可，大大减小了获取基础数据的实验测试工作量，而且提高了曲线预测的精度。从包装材料生产企业角度来看，第2种改进方法相对能够更加全面地描述缓冲材料的动态缓冲性能，但是预测精度相对会有所降低。

[1] 国家质量监督检验检疫总局, 中国国家标准化管理委员会. 缓冲包装设计: GB/T 8166—2011[S]. 北京: 中国标准出版社, 2012: 6-8.

General Administration of Quality Supervision, Inspection and Quarantine of the People's Republic of China, Standardization Administration of the People's Republic of China. Package Cushioning Design: GB/T 8166-2011[S]. Beijing: Standards Press of China, 2012: 6-8.

[2] 金强维，孙德强，李国志．密度对聚丙烯泡沫材料动态冲击性能的影响[D]. 西安: 陕西科技大学, 2019: 15-30.

JIN Q W, SUN D Q, LI G Z. Influence of Density upon Dynamic Crushing Behaviour of Polypropylene Materials[D]. Xi'an: Shaanxi University of Science and Technology, 2019: 15-30.

[3] DAUM M. A Simplified Process for Determining Cushion Curves: The Stress-Energy Method[D]. East Lansing: Michigan State University, 1999: 1-8.

[4] MARCONDES P, BATT G, DARBY D, et al. Minimum Samples Needed to Construct Cushion Curves Based on the Stress Energy Method[J]. Journal of Applied Packaging Research, 2008, 2(3): 191-198.

[5] KENDALYN P. Statistical Analysis of the Stress- Energy Methodology Applied to Cushion Curve Determination[D]. South Carolina: Clemson University, 2012, 2(4): 81-89.

[6] 史岩, 李光. 确定缓冲曲线的应力-能量法研究[J]. 包装学报, 2014, 6(3): 35-39.

SHI Y, LI G. Research on Stress-Energy Method for Determining Cushioning Curve[J]. Packaging Journal, 2014, 6(3): 35-39.

[7] 王金梅, 刘乘. 应力-能量法求取泡沫塑料缓冲曲线时函数模型的研究[J]. 包装工程, 2014, 35(5): 79-82.

WANG J M, LIU C. Function Models of the Stress-Energy Method in Determining Cushion Curves of Foam Materials[J]. Packaging Engineering, 2014, 35(5): 79-82.

[8] 陈满儒, 刘兰. 基于应力-能量法的缓冲包装材料性能研究[J]. 包装工程, 2018, 39(5): 44-47.

CHEN M R, LIU L. Properties of Cushion Packaging Materials Based on Stress-Energy Method[J]. Packaging Engineering, 2018, 39(5): 44-47.

[9] 李光, 要志雯. 基于能量密度确定缓冲材料性能曲线的研究[J]. 包装工程, 2020, 41(23): 103-108.

LI G, YAO Z W. Determining Performance Curves of Cushioning Material Based on Energy Density[J]. Packaging Engineering, 2020, 41(23): 103-108.

[10] 陈磊. EPP包装材料静态压缩与动态缓冲性性能研究[J]. 今日印刷, 2021, 20(4): 43-47.

CHEN L. Research on Static Compression and Dynamic Buffering Properties of EPP Packaging Materials[J]. Print Today, 2021, 20(4): 43-47.

[11] 邢紫玉, 宋海燕, 王立军. 聚丙烯夹芯板动态缓冲性能的研究[J]. 包装与食品机械, 2021, 39(6): 26-30.

XING Z Y, SONG H Y, WANG L J. Study on Dynamic Cushioning Characteristics of Polypropylene Sandwich Panel[J]. Packaging and Food Machinery, 2021, 39(6): 26-30.

[12] PHAM T M, CHEN W S, KINGSTON J, et al. Impact Response and Energy Absorption of Single Phase Syntactic Foam[J]. Composites Part B, 2018(150): 226-233.

[13] 李超, 张雄飞. 应力-能量法在柱形空气衬垫中的应用研究[J]. 绿色包装, 2023(10): 26-30.

LI C, ZHANG X F. Research on the Application of Stress-Energy Method in Cylindrical Air Cushion[J]. Green Packaging, 2023(10): 26-30.

[14] PIATKOWSKI T, OSOWSKI P. Modified Method for Dynamic Stress-Strain Curve Determination of Closed-Cell Foams[J]. Packaging Technology and Science, 2016, 29: 337-349.

[15] NAVARRO J P, GARCIA R M M, CLOQUELL B V, et al. Evaluation of Two Simplified Methods for Determining Cushion Curves of Closed Cell Foams[J]. Packaging Technology and Science, 2012, 25(4): 217-231.

[16] 张波涛. 应力-能量法在测定泡沫塑料缓冲曲线中的应用[J]. 包装工程, 2008, 29(1): 59-60.

ZHANG B T. Application of the Stress-Energy Method in Determining Cushion Curves of Foam Materials[J]. Packaging Engineering, 2008, 29(1): 59-60.

Rapid Acquisition Method of Maximum Acceleration-Static Stress Curve for Foamed Polyethylene

SONG Weisheng, XUE Yang

(Henan University of Animal Husbandry and Economy, Zhengzhou 450046, China)

The work aims to study a method that can predict the maximum acceleration-static stress curve quickly and accurately. Firstly, the maximum acceleration-static stress curves of foamed polyethylene of 5 different thickness at 5 different heights were obtained by a drop hammer impact testing machine. On this basis, the differences between the three different improved fitting methods proposed and the existing dynamic stress and dynamic energy polynomial fitting methods were analyzed and compared. It was found that when the maximum acceleration factor was used as a function value and the drop height, pad thickness, and static stress were used as variables for fitting without distinguishing heights, the average2value representing prediction accuracy was 0.835, which was much higher than the value of 0.299 6 got by the polynomial fitting method of dynamic stress and dynamic energy. However, the prediction accuracy on the right side of the curve was still low. After a polynomial with static stress as the variable was used as the correction factor, the average value of2was 0.934, indicating a significant improvement in prediction accuracy. The prediction deviation on the right side was reduced, but it still existed. When a formula with a correction factor was used for prediction while heights were distinguished, the average value of2was 0.984, and the phenomenon of gradually increasing prediction deviation towards the right side of the curve was significantly improved. Under different heights, the use of prediction formulas with correction factors can quickly and accurately predict the maximum acceleration-static stress curve, which can provide certain assistance for impact protection design and related software development.

maximum acceleration-static stress curve; stress energy method; prediction accuracy; foamed polyethylene; polynomial fitting

TB485.1

A

1001-3563(2024)05-0309-06

10.19554/j.cnki.1001-3563.2024.05.037

2023-11-08

河南省科技攻关项目（222102210267）；河南省高校重点科研项目（24B140004）