3D打印夹芯复合材料模拟冰型设计与分析

彭锦峰，吴东润，崔为运，蔡登安，周光明

南京航空航天大学 机械结构力学及控制国家重点实验室，南京 210016

飞机在湿润的环境中巡航时，空气中的水滴则撞击并吸附在飞机迎风部件上。由于液体在固化的过程中会释放本身热量，当大量水滴聚集在飞机表面上时，水滴温度会急剧降低至降至结冰水平，使部分水滴在飞机蒙皮表面凝结成冰[1-2]。飞机不同部位结冰对飞行安全有着不同的影响，例如机翼前缘结冰厚度达到1.3厘米时，机翼升力下降50%，阻力增加30%，飞机的失速速度也会大幅度提高。结冰厚度继续增加，使得飞机整体重心前移，使飞机形成一个低头力矩，加大对飞机操作的难度，严重影响飞行安全[3]。故在美国联邦航空法规23部分1419号“防冰”条款中指出“飞机必须通过连续最大结冰型件和最大间断结冰状态的试验来验证飞机能够安全运行”[4]。此条款说明在适航认证的中，获得飞机在适航条例规定环境下的结冰形状，从而对飞机表面结冰的情况进行安全评估的环节必不可少[5]。

根据中国民用航空条例25部规定，一款取得适航认证的飞机必须进行自然结冰的试验[6]，但是自然结冰下的飞行测试成本较高、难度较大，同时自然结冰形状多样化，单次飞行试验难以达到实际试验的要求，试验结果不具备普遍性。因此通过安装模拟冰型试验不仅减少自然结冰飞行试验的次数降低了试验成本，提高了适航测试效率，还可以提前发现飞机性能的缺陷，从而对飞机进行再改进设计。在飞机机翼的关键位置安装相应的临界模拟冰型(临界冰型指导致飞机飞行性能下降最严重情形下的积冰形状)，在干空气条件下进行飞行试验，得到飞机结冰状态下的气动特性。

在模拟冰型的研制的进程中，Raab等[7]基于机型Phenom 300进行研究，通过测试软件对实际飞行过程中形状不同的模拟冰型对机翼的气动性能影响进行了分析与研究，验证安装模拟冰型的飞机在风洞中进行实验的结果和气动仿真计算的结果。Broeren等[8]在NACA 23012机翼前缘上安装全尺寸模拟冰铸件，通过冰风洞试验对整个机翼进行气动性能的研究，并总结出短脊冰对机翼流场较高脊冰对机翼流场影响较大。对比国外的研究进展，中国在模拟冰型方面相对应的研究较少。ARJ21-700作为中国独立自主研制的第一款支线客机，其模拟冰型采用以红松芯材为主体，树脂基玻璃纤维增强材料作为基体的复合材料夹层结构[9]，通过有限元计算预测了不同载荷状态下的破坏载荷，并通过地面破坏试验验证了有限元模型的准确性，通过对比技术要求以及气动载荷，验证了这种模拟冰型的可行性。

相比于传统机械加工如钻孔、切削等材料去除技术，3D打印技术通过液体、粉末或薄板材料的结合，逐层形成实体，是一种自下而上的材料增加技术[10]。北京航空航天大学的王华明教授以及他的团队于2013年利用3D打印技术，在全世界范围内首次制造出飞机大型主承力构件[11]，从而验证了3D打印技术的可靠性，说明3D打印技术可以适用于高精端的航空航天领域并且能够满足一定的结构性能要求。根据工作原理的差异，3D打印一般有以下几种形式：光固化成型、选择性激光烧结、分层实体制造、数字化光处理。由于采用光固化成型技术的3D打印方式能够提高材料的利用率，完成结构复杂、尺寸精度要求高的模型制造[12]。

以某民用飞机机翼自然结冰过程为研究对象，确定飞机的临界冰型，基于3D打印夹芯材料，设计模拟冰型结构，通过CFD仿真得到带冰模型在极端工况下所承受的气动载荷，利用有限元分析得到模拟冰型的失效载荷及破坏模式。进行了3D打印夹芯复合材料模拟冰型的制备与试验，验证了分析方法的有效性及设计方案的可行性。

1 模拟冰型结构设计

研究对象为某民用飞机的结冰机翼，其机翼弦长3.314 m。在结冰情况下，该翼型可能形成如图1所示的临界模拟冰型。

图1 临界冰型截面形状Fig.1 Critical ice cross-sectional shape

在ARJ21-700飞机模拟冰型复合材料夹层结构的设计基础上[13]，对整个模型结构进行改进处理，通过粘贴铝胶带代替原结构内层玻璃钢，采用3D打印技术制作空心光敏树脂芯材。模型构成如图2所示。

图2 改进设计后模拟冰型夹层结构Fig.2 Simulated ice sandwich structure after improved design

在安装模拟冰型之前，铺设一层铝胶带对机翼表面蒙皮进行保护，然后采用环氧树脂将模拟冰型与蒙皮表面的铝胶带进行粘接，等待环氧树脂固化之后在模拟冰型外表面铺设环氧树脂基玻璃纤维增强面板(俗称玻璃钢)，用于对模拟冰型进行保护，同时设置玻璃钢的方法可以为整体模型结构提供一定的强度与刚度。

采用光固化成型的3D打印模拟冰型夹芯结构。由于光敏树脂密度大于明冰密度，使得等体积模拟冰型的质量大于真实结冰情况下冰块的质量，在质量上不具备相似性。因此将模拟冰型设计为空腔结构，这样设计虽然有效地降低了模拟冰型质量，使得模拟冰型在质量上与真冰更具相似性。但是这样会使得冰型整体刚度和强度的降低，可能导致模拟冰型在带冰飞行测试过程中出现损伤或者脱落情况发生。由于在真实结冰状况下，冰块的各个部位形状变化多端，个别位置的扭曲程度较大，故设计将沿展向横截面变化不大的实际冰块划分成展长250 mm的冰型，在进行等截面处理后进行分析。因此，如图3所示设计壁厚为1 mm的模拟冰型进行分析，验证该模拟冰型能否满足适航飞行测试的要求。

图3 沿中心界面剖开的壁厚1 mm的模拟冰型Fig.3 Simulated ice with a wall thickness of 1 mm cut along central interface

2 结冰翼型气动载荷分析

气动分析根据连续最大结冰型件，对一般民用飞机处于跨音速飞行的极端飞行状态进行气动仿真计算，进而研究粘接于机翼前缘的冰型所承受气动载荷的强度。

最大结冰型件指出飞机巡航区域位于处于距离海平面6700 m，环境满足最大云层垂直厚度1980 m，云层覆盖的水平标准距离32.22 km时，飞机处于高风险结冰状态[14]。表1给出了仿真工况用于计算处于高风险结冰区域飞机巡航中结冰机翼的环境参数。

表1 高风险结冰区域环境参数Table 1 Environmental parameters of high-risk icing areas

利用Fluent对带冰翼型在高风险结冰区域环境内，以不同迎角进行飞行的状态进行仿真，得到不同飞行迎角下带冰翼型的气动参数，进而研究粘接于机翼前缘的冰型所承受的气动载荷。由于机翼结冰使飞机整体的气动性能大幅下降，导致飞机的失速迎角减小，严重威胁飞机的飞行安全[15]。设置的迎角变化范围为-4°至6°。图4所示为结冰翼型处于2°至6°迎角下的压强系数分布和流线分布。

图4 不同迎角下结冰翼型表面压强分布Fig.4 Pressure distribution of icing airfoil surface at different angles of attack

如图4所示，在高速飞行下，结冰对翼型绕流影响较大，积冰的角状突起处压力系数变化剧烈，改变了翼型前缘的压强分布。随着迎角的增大，积冰迎风面处压力激增，积冰背风面与前缘处压强变化更加复杂，翼型上下表面的逆压梯度增大，使得气体流动分离提前，翼型的气动特性受到了严重的影响。

图5所示为高风险结冰区域内跨音速飞行时6°迎角结冰翼型表面的流线分布。冰型的迎风面出现大面积高压区，从而产生严重轴向载荷的根本原因。其次冰型对绕流影响比较严重，冰型上下背风面均出现驻涡，上下驻涡的强度变化复杂，同时冰型与翼型前缘粘接处的上表面形成一片负压区，从而产生严重法向载荷的根本原因。

图5 6°迎角下结冰翼型表面流线分布Fig.5 Distribution of streamlines on iced airfoil surface at 6° angle of attack

如图6所示，将气动载荷分解成沿翼型弦线方向的轴向载荷与垂直于弦线方向的法向载荷进行研究，记轴向载荷为F1，法向载荷为F2。

图6 气动载荷分解Fig.6 Aerodynamic load decomposition

通过Fluent后处理得到高风险结冰区域内以大迎角跨音速飞行下冰型所承受的轴向载荷与法向载荷情况见表2。随着迎角的变化，冰型所受轴向载荷的变化较小，法向载荷变化较大，且法向载荷随着迎角的增大而大幅增加。

表2 气动载荷与迎角变化关系

考虑模拟冰型的安装方式，法向载荷是影响冰型与机翼粘接强度的主要因素；而轴向载荷的作用实际上是将冰型牢牢压在翼型前缘上，对冰型与机翼胶接面的影响程度小于法向载荷。因此确认6°迎角时结冰翼型承受最大气动载荷。

由于气动仿真模型经过归一化处理，其真实弦长为3.314 m，因此该模型实际展向线载荷为3 503.46 N/m，最大法向线载荷为1 755.92 N/m。同时模型沿展向延伸250 mm，所以实际等效的最大法向载荷F2 max=438.98 N。

3 模拟冰型失效分析

如图7(a)所示，采用ABAQUS软件建立有限元分析模型，网格单元总数为260 780个。其中机翼蒙皮、玻璃钢以及模拟冰型均采用C3D8R单元(8节点六面体线性减缩积分单元)，单元数量为227 780个，胶层采用COH3D8单元(8节点三维内聚力单元)，单元数量为330 00个，离散模型如图7(b)所示。

图7 有限元模型Fig.7 Finite element model

有限元模型中，光敏树脂密度为1 100 kg/m3，弹性模量为3.64 GPa，泊松比为0.29，铝合金弹性模量为73 GPa，泊松比为0.33。表3给出了环氧树脂基玻璃纤维增强面板的9个工程弹性常数[16]，其中E1、E2、E3分别为复合材料3个方向的拉伸模量；G12、G23、G13分别为3个面内剪切模量；μ12、μ23、μ13分别为3个面内主泊松比。界面元应用于对环氧树脂粘接层损伤的模拟[17-19]，设定其厚度为0.2 mm。采用了非耦合本构关系，二次应力判据作为损伤起始判据，以BK准则为损伤扩展判据[20]，内聚力单元的刚度系数和强度[9]由表4给出，表中：Knn、Kss和Ktt分别为界面法向与两个切向的内聚力单元刚度系数；tn、ts和tt分别为平拉试验测定的界面拉伸强度与双剪切试验测定的剪切强度。选取的损伤扩展判据BK准则涉及到的Ⅰ型裂纹能量释放率GⅠC为0.218 N/mm、Ⅱ型裂纹能量释放率GⅡC为0.56 N/mm、Ⅲ型裂纹能量释放率GⅢC为0.56 N/mm[21]。

表3 玻璃纤维增强复合材料参数Table 3 Material parameters of glass fiber reinforced composites

表4 界面元刚度系数与强度

对铝合金蒙皮内表面全部固支约束，同时利用等效计算给外层玻璃钢表面施加位移量，使得模型产生变形，从而达到破坏胶层的目的。通过ABAQUS仿真计算，在约束条件下得到胶层部件上的位移-载荷趋势图，从而获得加载情况下模拟冰型的胶接强度。由图8可知，胶层结构破坏的极限载荷为2 539.23 N，对比仿真计算得到的最大法向气动载荷，说明胶层结构破坏载荷大于最大法向气动载荷，从理论上验证在飞行安全允许的前提下即使处于最恶劣的飞行环境中，模拟冰型能够胜任带冰飞行测试的任务。

图8 胶层部件的位移-载荷趋势图Fig.8 Displacement-load trend diagram of adhesive layer components

根据图8可知，在初始阶段的载荷位移曲线基本呈线性增加的趋势。但对于1 mm壁厚的模型而言，在位移U=0.007 mm，载荷F=506.889 N时，载荷加载方向的界面元与外层玻璃钢交界处开始出现损伤，损伤位于交界处的两端边缘，但是损伤程度较小，损伤区域也不明显，如图9所示。

图9 模型在U=0.007 mm时的界面损伤云图Fig.9 Interface damage cloud of model at U=0.007 mm

随着位移进一步增加，当位移U=0.031 mm时，此时载荷F=1 843.43 N，曲线斜率开始下降，中间界面元发生大面积受到损伤，同时界面元单元只在与外层玻璃钢交界处边缘位置发生明显破坏，在其余位置损伤程度并不明显，因此整体结构仍能继续承载，但是其承载能力开始下降，U=0.031 mm时的界面损伤云图如图10所示。

图10 模型在U=0.031 mm时的界面损伤云图Fig.10 Interface damage cloud of model at U=0.031 mm

当位移U=0.061 mm时，载荷达到峰值2539.23 N，此时在位移加载方向的界面元与外层玻璃钢交界处的单元几乎全部发生破坏，随着位移进一步增大，载荷反而开始下降，说明此时结构已经无法继续承载，该模拟冰型结构在升力作用下的破坏载荷Fmax=2539.23 N，其损伤云图如图11所示。

图11 模型在U=0.061 mm时的界面损伤云图Fig.11 Interface damage cloud of model at U=0.061 mm

在破坏载荷作用下，需对模拟冰型其他部位强度进行校核以保证模拟冰型不会发生破坏，尤其是光敏树脂冰型是空心结构。因此需对模拟冰型的芯材部分进行强度校核。

图12为极限载荷作用下中间冰型的应力分布云图。由于壁厚较薄的因素，冰型在玻璃钢与胶层结合处区域内应力水平较高，产生明显的应力集中，并且向上逐渐减小，产生明显的阶梯分布。冰型上表面的应力分布较为复杂，平均应力高于其余部位的应力，这是由于空心结构传力性能较差，因此除冰型上表面外，其余部位应力相对处于低水平状态。作为冰型材料，各向同性材料的3D打印光敏树脂拉伸强度为45.76 MPa[22]。根据仿真计算结果，冰型最大应力17.62 MPa，说明即使粘接模拟冰型和机翼蒙皮表面的胶层在因外力作用而失效的情况下，光敏树脂模拟冰型也不会发生破坏。

图12 模拟冰型应力分布云图Fig.12 Cloud map of simulated ice stress distribution

4 静力试验

如图13所示，加工外侧曲面与机翼外形尺寸一致的铝模，用作试验基座来模拟机翼前缘外形，同时采用3D打印技术制作壁厚1 mm的模拟冰型。

图13 试验部件Fig.13 Test components

将铝胶带铺设在铝模表面，并用塑料刮子将其刮平。将模拟冰型与粘有铝胶带的铝模配合，确定冰型边线并且在上下翼面的冰型边线处分别向外延伸50 mm作为玻璃钢边线玻璃钢由玻璃纤维、环氧树脂胶液配置，其中环氧树脂胶液由环氧树脂618、苯二甲胺固化剂和邻苯二甲酸二丁酯按照质量比为10:2:1的比例配置而成，同时将胶液充分搅拌并涂均匀涂抹在芯材内表面，之后安装在相应位置的铝胶带上。待胶液初步固化，冰型不在重力作用下滑动则表示安装结束。

试验如图14所示，采用帆布袋杠杆系统，帆布带设计为展向长度230 mm，弦向长度100 mm。通过手动液压车进行加载，采用拉压载荷传感器和载荷记录仪采集载荷数据。使用DH3816型静态电阻应变仪采集应变数据，试验过程中每300 N采集一次应变数据以监控试验过程中结构表面的应变水平，连续加载至模拟冰型破坏或模拟冰型从模具表面脱落为止。

图14 试验设备Fig.14 Test equipment

根据有限元计算结果以及帆布带分布，在模拟冰型上表面以及迎风表面区域内布置双向应变片用于采集加载过程中模拟冰型表面的应变数据，双向应变片采用0°/90°粘贴方向与玻璃钢中纤维方向相同，应变花位置如图15所示，其中#1粘贴于玻璃钢迎风面的中心处，用于采集迎风面处的应变数据；参考图12冰型的应力云图，发现位于上角点区域内的应力水平较高，故#2与#3分别对称粘贴于模拟冰型上角点截面处，用于采集冰型受载荷作用时的应变。

图15 应变片分布位置Fig.15 Distribution of strain gauge

如图16所示，当模拟冰型结构试验件被破坏时，呈现出环氧树脂胶层与铝胶带之间分离，同时粘接于机翼表面的玻璃钢发生部分脱落，但是整体结构并未发生完全破坏。与有限元仿真结果类似，模拟冰型外形完整，模型结构未受到严重破坏，验证采用光敏树脂的空腔模拟冰型结构能够应用于飞机带冰适航测试。同时机翼铺设铝胶带表面凹凸不平，并且部分区域内的环氧树脂出现脱胶的现象发生。记录试验件破坏载荷以及模拟冰型有效粘接面积具体结果见表5。

图16 模拟冰型结构试验件破坏结果Fig.16 Failure results of simulated ice structure test pieces

由表5可知，有限元计算结果与试验结果相对误差为4.72%，验证上述的有限元仿真的有效性，同时采用的有限元模型也能较为准确地计算实际模拟冰型的破坏载荷以及预测模拟冰型的破坏形式。根据模拟冰型试验件的破坏形式可以看到，界面破坏顺序首先为靠近冰型端面处的环氧树脂出现损伤，然后损伤逐渐由端面向内部进行扩散。

表5 模拟冰型静力破坏试验结果Table 5 Static failure test results of simulated ice

理论上2号与3号应变片测量的应力值相同，但考虑实际粘贴过程中应变片定位存在一定偏差，所以根据试验实际测量的应变数据，将二者应变取平均值，记为2号应变片的数据进行分析。如图17所示，1号应变片用于采集迎风面玻璃钢的应变数据，在法向载荷作用下应变较小，且主要受力形式为受拉，最大应变发生在0°方向，约为440 με。2号应变片用于采集模拟冰型芯材端面靠近加载面区域的应变数据，此处为应力集中区域，且主要受力形式为受拉，同时最大应变出现在0°方向，约为1 250 με，符合有限元仿真的结果。

图17 应变片应变载荷曲线Fig.17 Strain-load curve of strain gauge

根据气动载荷计算结果，该段展长为250 mm的带冰机翼在5倍安全系数下的极限法向载荷为2.19 kN。通过静力拉伸试验得到模拟冰型试验件的破坏载荷2.42 kN。因此该模拟冰型能够在5倍升力载荷作用下不发生脱落，说明设计的采用光敏树脂材料设计壁厚1 mm的模拟冰型是可以满足实际适航测试的需要。

5 结 论

1) 对结冰翼型进行了气动仿真分析，阐明了翼型积冰对翼型的升阻力特性带来了不利影响。由于明冰出现在翼型前缘处，故翼型绕流在冰型处易产生紊乱流动，影响整个流场，使流动提前分离，导致翼型过快地进入失速状态，对飞行安全造成很大威胁。

2) 对模拟冰型夹层结构进行有限元仿真，分析了模拟冰型与机翼前缘表面之间胶层失效情况，得到了胶层破坏的极限载荷。通过对比极端气动载荷与胶层破坏载荷，在理论上验证民机即使处于极限飞行环境下，测试用的模拟冰型不会出现破坏或者脱落的情况发生。进一步说明采用3D打印技术制备的空腔模拟冰型能够满足测试所需要的强度要求。

3) 对比有限元仿真与地面静力试验结果，在可接受的误差范围之内，验证了仿真计算结果的有效性，确定了适航试飞用3D打印夹芯复合材料模拟冰型结构的可行性。