【www.zhangdahai.com--其他范文】
张 健， 陈秀华， 陈 勇， 范 寅
(上海交通大学 a. 机械与动力工程学院；
b. 航空航天学院，上海 200240)
复合材料由于强度高、刚度大、质量轻和抗疲劳等优点，在航空航天、能源和交通等领域应用广泛[1].随着热塑性预浸料纤维浸渍技术及低成本制造技术的发展和成熟，热塑性复合材料由于高韧性、高损伤容限、生产效率高、易修复和可回收等众多优点，在商用航空、汽车和风力发电叶片等领域逐渐得到应用[2-4].为了设计出一个能满足要求的热塑性复合材料结构件，在设计分析阶段需要一个合适的本构方程来描述热塑性复合材料在不同载荷条件下的力学行为.对于热塑性复合材料，已有试验研究表明AS4/PEEK、C/PPS和AS4/PEKK等热塑性复合材料在剪切方向会由于塑性变形而表现出显著的非线性行为[5-7].另外，由于不同组分的力学性能存在较大差异，复合材料内部的应力分布极不均匀，其损伤模式极为复杂，包括纤维断裂、基体断裂和纤维/基体脱黏等，这些损伤模式可以单独或者同时发生[8].在加载的过程中，复合材料将会出现上述损伤并逐渐演化直至失去承载能力.因此，为了准确预测热塑性复合材料的非线性力学行为和极限强度，在本构方程中应同时考虑塑性和损伤的影响.
连续介质损伤力学为复合材料的有限元分析提供了判断损伤起始和损伤演化的框架，国内外学者在此基础上建立了众多的弹性损伤模型[9-17].弹性损伤模型假设复合材料在破坏前表现为线弹性，该模型适用于描述在加载过程中应力-应变曲线表现为线弹性的复合材料.但对于在加载过程中塑性行为较为明显的热塑性复合材料，弹性损伤模型的模拟效果并不理想.例如Maa等[18]使用弹性损伤模型分析AS4/PEEK开孔层压板的极限拉伸载荷，模拟结果与试验结果误差在10%以上，最大的误差达到了25.83%，误差较大的原因为模型所假设的纤维方向为脆性断裂并不能准确模拟实际的损伤演化.上述研究表明，对于在加载过程中产生了塑性变形的热塑性复合材料，弹性损伤模型将产生较大的误差；
Xiao[19]的研究表明，未考虑塑性效应会低估复合材料结构的能量吸收能力，从而导致较大的误差.
为了同时考虑复合材料的塑性行为和损伤累积，以往研究在弹性损伤模型中加入塑性项，即在塑性力学的框架内建立弹塑性损伤模型.Ladeveze等[20]提出一种复合材料弹塑性损伤模型，该模型不仅计入了塑性变形，同时也考虑了由基体开裂和纤维/基体脱黏而造成的刚度退化，但模型不考虑沿纤维方向的损伤，与试验结果对比显示，该模型能够准确模拟复合材料的塑性力学行为和破坏行为；
陈静芬[8]提出一种基于单参数塑性模型[21]和Hashin准则[22]的二维弹塑性损伤本构模型.为解决有限元隐式计算中由于应变软化而导致的计算困难问题，作者对损伤变量进行了黏性正则化.结果表明该模型能够准确预测存在显著塑性力学行为的复合材料的破坏载荷，但该模型忽略了层间应力的影响，因而对于预测包含分层损伤的复合材料强度问题显得不够准确.
本文对不同角度的单向AS4/PEEK层压板进行了偏轴拉伸试验和相应的数值模拟.通过对断裂面的观察，分析了断裂面角度与铺层角度存在差异的原因.对不同角度单向层压板的应力-应变曲线的分析表明，在破坏前纤维方向仅有线弹性行为，将Sun等[21]提出的三维塑性势函数进行了简化，塑性势函数中的参数通过信任域反射算法计算得到.结合LaRC05准则[23]和裂纹带理论[24]建立了三维弹塑性损伤本构模型，通过编写基于Abaqus的用户自定义材料子程序VUMAT并将其应用于不同铺层角度的偏轴拉伸试验的数值模拟中，与试验获得的结果进行对比，验证了该模型的正确性.本文提出的三维弹塑性损伤模型考虑了除纤维方向外各方向的塑性变形并以及断裂面角度对复合材料损伤的影响，与文献[8,18]相比，该模型能更加准确地预测试验件的拉伸强度.
1.1 试验设置
试验采用单向AS4/PEEK层压板，共有7种铺层，分别为[0]8、[15]20、[30]20、[45]20、[60]20、[80]20和[90]16.其中，[0]8铺层的试验件尺寸为250 mm×15 mm×1.12 mm，[90]16铺层的试验件尺寸为250 mm×25 mm×2.24 mm，其余试验件的尺寸为250 mm×25 mm×2.8 mm.试验件两端粘有加强片，[0]8铺层加强片的尺寸为25 mm×15 mm×1.5 mm，其余试验件加强片的尺寸均为 25 mm×25 mm×1.5 mm，斜削角均为90°.每种铺层有5个试验件.
图1 偏轴拉伸试验的试验设备Fig.1 Test equipment for off-axis tensile test
试验使用UTM5105X微机控制电子万能力学试验机进行加载，如图1所示.该试验机最大载荷为100 kN，准确度精度为0.5级；
使用应变片测量试验件的应变，应变片型号为BX120-3AA，精度等级为A.试验根据ASTM D 3039/D 3039M-17[25]标准进行，采用位移加载，加载速度为5 mm/min.
图2 偏轴拉伸应力-应变曲线Fig.2 Stress-strain curves of off-axis tensile
1.2 试验结果及分析
采用ASTM D3039/D3039M-17[25]提出的方法对试验数据进行处理，拉伸应力的计算公式为
(1)
式中：σi为第i个点的拉伸应力；
Pi为第i个点的力；
A为试验件的横截面积.
偏轴拉伸试验的应力-应变(σ-ε)曲线如图2所示.从应力-应变曲线可知，当偏轴角度为α=0°时，试验件在破坏前的力学行为都表现为线弹性；
随着偏轴角度的逐渐增加，试验件的力学行为首先表现为线弹性，在应力超过一定范围后表现出明显的非线性行为；
当偏轴角度在80°～90°之间时，试验件的力学行为只有一小段的线弹性行为，随后发生断裂，原因是这些偏轴角度纤维已经不能承受足够的载荷，在表现出塑性变形前基体已经发生了断裂破坏.
不同铺层的破坏模式、破坏部位和断裂面角度总结如表1所示.对比不同试验件的断裂面角度和偏轴角度，断裂面角度与偏轴角度的差异在 -1°～2°之间.破坏后的试验件图片如图3(a)所示，断裂面的显微照片如图3(b)所示.断裂面内观察到的主要是基体的断裂，在基体的断裂处可以观察到裸露的纤维，纤维只有少数发生了断裂，未发生断裂的纤维没有观察到明显的偏转.对断裂面的观察和分析表明，断裂面角度与偏轴角度存在差异的原因为在试验件的制造过程中， 层压板的铺层角度与设定的角度存在一定的误差.
表1 不同铺层的破坏模式、破坏部位及断裂面角度[25]
图3 试验结束后的试验件及断裂面Fig.3 Test specimens and fracture surface after test
2.1 三维弹塑性损伤模型
由试验结果可知，单向AS4/PEEK层压板在偏轴拉伸载荷作用下存在显著的由塑性变形引起的非线性力学行为及破坏行为.为了准确模拟上述行为，需要建立一个包含弹性模型、塑性模型和损伤起始及演化模型的三维弹塑性损伤模型.
2.1.1弹性模型 在线弹性阶段，本构关系采用三维模型.对于正交各向异性层压板，其典型的应力-应变关系为
(2)
式中：
(3)
(4)
(5)
(6)
(7)
(8)
C44=G12
(9)
C55=G23
(10)
C66=G13
(11)
(12)
式中：E1、E2和E3为弹性模量；
G12、G23和G13为剪切模量；
ν12、ν13和ν23为泊松比；
σij(i,j=1, 2, 3)为材料坐标系下的应力；
εij(i,j=1, 2, 3)为材料坐标系下的应变.
2.1.2塑性模型 本模型采用的塑性屈服函数为
(13)
Fp(σ)=
(14)
式中：a66为复合材料各向异性程度的参数，假定其值与材料有关.
等效应力定义为
(15)
采用关联流动法则确定塑性应变增量为
(16)
式中：dλ为比例因子.
为保证单位体积塑性功率(dWp)相等，有
(17)
(18)
强化函数采用Sun等[21]提出的幂函数形式的各向同性强化函数，为
(19)
式中：β和n为材料参数，可通过拟合等效应力-等效塑性应变曲线得到.
塑性势函数中的参数a66可以通过偏轴拉伸试验得到[21, 27].如图4所示，x-y坐标系为整体坐标系，1-2坐标系为材料坐标系，x轴和1轴之间有一个大小为θ的夹角，在平行于x轴的拉伸载荷的作用下，复合材料层压板的受力状态视为平面应力状态.将层压板在整体坐标系的应力转换到材料坐标系，转换公式为
图4 偏轴拉伸试验示意图Fig.4 Schematic diagram of off-axis tensile test
σ11=cos2θσx
(20)
σ22=sin2θσx
(21)
σ12=-sinθcosθσx
(22)
式中：σx为x方向的应力.将式(20)～式(22)代入式(15)可得等效应力为
(23)
式中:
(24)
由塑性应变增量的坐标转换可得
(25)
(26)
假设材料具有显著的塑性效应，且满足小变形假设，则总应变增量可以分解为弹性应变增量和塑性应变增量两部分[8]，即
(27)
故x方向的塑性应变增量可以表示为
(28)
(29)
通过上述变换，可将偏轴拉伸试验获得的应力-应变曲线转换为等效应力-等效塑性应变曲线.对于不同角度的偏轴拉伸试验，经过转换后的等效应力-等效塑性应变曲线应该是一致的，故可以通过选取适当的a66的值使得不同角度的等效应力-等效塑性应变曲线重合.对于a66的选取，目前有试错法[21]和信任域反射算法[28]等计算方法.本文采用Jang等[28]提出的信任域反射算法进行计算.具体计算流程为：首先设定a66的初值，此处设定为0.5，然后采用式(20)～式(29)将试验获得的不同角度的应力和应变数据变换为等效应力和等效塑性应变数据，将等效应力视为等效塑性应变的函数，通过最小化目标函数来获得a66.目标函数的形式为
(30)
(θ1,θ2=15°, 30°, 45°, 60°, 80°, 90°)
2.1.3损伤起始及演化模型 对于纤维拉伸损伤，本模型采用最大应力准则来预测，即
(31)
式中：Xt为纤维方向的拉伸强度.当fft≥1时，发生纤维拉伸损伤.
为了预测基体损伤，本模型采用LaRC05准则[23].该准则认为基体损伤的发生只与潜在断裂面上的正应力(σn)和两个剪应力(σnT和σnL)有关，表达式为
(32)
(33)
(34)
式中：ϑ为潜在断裂面角度，范围为0°～180°；
ST为横向剪切强度；
SL为纵向剪切强度；
μT为横向摩擦系数；
μL为纵向摩擦系数；
YT为基体抗拉强度；
Yc为基体抗压强度；
φ0为纯横向压缩破坏的断裂面的角度，对于碳纤维增强复合材料，试验[29]测定其值范围为φ0=53°±2°.当基体最大应力fmt≥1时，发生基体拉伸损伤.
在使用失效准则对复合材料的损伤模式进行预测后，需要采用适当的损伤变量来模拟复合材料发生损伤后的力学行为.本模型使用的损伤变量包括纤维拉伸损伤变量(dft)和基体拉伸损伤变量(dmt)，损伤变量从损伤起始时的0逐渐演化为完全退化时的1[30].损伤变量的更新采用Donadon等[31]提出的形式，即
(35)
(36)
式中：εf为最终失效时的应变.
对于纤维的损伤，当前增量步的等效应变对应于1方向上的应变分量，损伤开始时的应变和最终失效时的应变的计算公式为
(37)
(38)
式中：Gft为纤维拉伸破坏断裂能；
L为特征长度.
对于基体的损伤，由于在断裂面上的拉应力会促进基体的损伤，损伤开始时的应变和最终失效时的应变的计算相对更加复杂.当满足式(33)的失效准则时，基体发生损伤，此时的等效应力和等效应变为
(39)
(40)
(41)
对于混合失效情形，基体失效的断裂能为
(42)
(43)
当层压板发生基体损伤时，通过式(32)将应力旋转到断裂面上，对断裂面上的应力进行软化，为了考虑反向载荷对裂纹闭合的影响，断裂面上的正应力只在其为拉伸应力时进行软化，断裂面上的应力软化式为
(44)
Matzenmiller等[32]认为纤维的断裂、屈曲和压缩破坏会对基体的损伤演化产生影响，故纤维损伤导致的应力软化计算公式为
(45)
2.2 数值模拟设置
依据三维弹塑性损伤模型的数值算法编写基于Abaqus的用户自定义材料子程序VUMAT，其计算流程如图5所示.
计算过程采用位移增量逐级加载的方式，在第n+1步的具体计算步骤如下.
(2) 假定当前增量步的塑性应变张量和等效塑性应变与上一增量步一致，计算试应力为
(46)
(3) 检查屈服条件：
(47)
(4) 采用Newton-Raphson方法迭代求解非线性方程：
(48)
表2 AS4/PEEK性能参数Tab.2 Properties of AS4/PEEK
表3 AS4/PEEK的强度及断裂能Tab.3 Strength and fracture energy of AS4/PEEK
首先初始化参量：
(49)
其次计算第k+1次迭代的比例因子增量为
(50)
并更新第k+1次迭代的塑性应变张量和等效塑性应变为
(51)
再次计算试应力为
(52)
Ftol(=10-6)
(5) 判断是否发生损伤起始.
首先计算fft，如果fft≥1，则计算纤维损伤变量dft；
否则dft=0.
然后计算fmt，如果fmt≥1，则计算基体损伤变量dmt；
否则dmt=0.
(6) 更新应力，进入下一增量步.
AS4/PEEK材料的性能参数、强度及断裂能[33-35]如表2和表3所示.其中，Gfc为纤维压缩压缩破坏断裂能，Xc为纤维压缩破坏断裂能.为了减少计算量，有限元模型不考虑加强片部分，模型长度为200 mm，其余尺寸与试验件尺寸一致，模型的铺层包含试验件的7种铺层.有限元模型的单元类型为C3D8R，厚度方向的单元尺寸为0.14 mm，长度和宽度方向的单元尺寸为1 mm.在加载过程中，如图1(b)所示试验件的上边界在x、y和z的3个方向的位移保持一致，对试验前后的试验件的上边界横截面尺寸进行测量，上边界横截面尺寸变形微小，可以将试验件的上边界视作刚性平面，故有限元模型的边界条件为一端约束x、y和z方向的位移Ux、Uy和Uz，另一端与参考点耦合， 在参考点上施加x方向位移并约束y和z方向的位移.将一组节点耦合至参考点并在参考点施加位移载荷已在文献[10,34,36]中使用广泛.包含边界条件的有限元模型如图6所示.
图5 VUMAT计算流程图Fig.5 Flow chart of subroutine VUMAT
图6 包含边界条件的有限元模型Fig.6 Finite element model with boundary conditions
3.1 塑性参数的确定
利用MATLAB软件编写程序进行塑性参数的求解.模型假设纤维方向是线弹性的，故[0]8铺层的数据未考虑，对其余铺层获得的应力-应变曲线进行计算，最优化目标函数后获得的塑性参数为a66=1.50，拟合等效应力-等效塑性应变曲线获得的硬化函数的参数值为β=292.67，n=0.134 6.选取上述值绘制的等效应力-等效塑性应变曲线如图7所示，其中ξ为铺层角度.拟合得到的等效应力-等效塑性应变曲线与试验值存在误差的原因有两个：一是试验测量获得的应力和应变数据存在一定的测量误差；
二是在进行数值优化过程中存在的舍入误差.
图7 等效应力-等效塑性应变曲线Fig.7 Equivalent stress-equivalent plastic strain curve
3.2 数值模拟结果
数值模拟获得的不同角度的应力-应变曲线与试验获得的应力-应变曲线的对比如图8所示，可以看出数值模拟获得的应力-应变曲线与试验获得的应力-应变曲线吻合较好，表明本模型能够较为准确地模拟单向AS4/PEEK层压板的塑性效应.
试验和数值模拟预测的拉伸强度见表4，表中实验值为5件试验件的拉伸强度的平均值.从表4可知，在试验件厚度相同的情况下，随着偏轴角度的增加，拉伸强度逐渐下降.数值模拟预测的拉伸强度与实验值的最大误差为19.05%，不超过20%，表明该模型适用于预测单向AS4/PEEK层压板偏轴拉伸的拉伸强度.
表4 拉伸强度及误差的对比Tab.4 Comparison of tensile strength and error
[45]20铺层数值模拟和试验获得的应力-总应变(εtotal)曲线对比如图9所示.结合图10的塑性应变(εp)-总应变曲线可知，当总应变在 1.5×10-3以下时，各个方向都不存在塑性变形，此时层压板表现为线弹性；
随着总应变的增大，除纤维方向外其余方向的塑性变形逐渐增大，使得层压板表现出明显的非线性力学行为；
当总应变达到 22.7×10-3时，应力出现急剧下降，表明此时层压板已经完全破坏.从图10可以看出，塑性变形主要发生在1-2平面内和3方向，且1-2平面内的剪切塑性应变约为2方向和3方向的塑性应变的3倍，表明剪切变形在塑性变形中处于主导地位.
图8 数值模拟与试验的应力-应变曲线对比Fig.8 Comparison of stress-strain curves between numerical simulation and test
图9 [45]20铺层试验与数值模拟应力-总应变曲线对比Fig.9 Comparison of stress-total strain curves between experiment and numerical simulation of [45]20 specimens
图10 塑性应变-总应变曲线Fig.10 Plastic strain-total strain curve
图11 [45]20铺层损伤变量Fig.11 Damage variables of [45]20 specimens
最大载荷和分析结束时的损伤变量(δ)云图如图11所示，当损伤变量值为0时表示未发生破坏，损伤变量值为1时表示已经完全破坏.模拟结果表明，当载荷达到7.34 kN时，层压板首先在端部的边缘发生基体破坏；
随着载荷的增加，在基体的破坏部位附近发生纤维破坏，然后基体破坏和纤维破坏开始向未破坏的部位传播；
当载荷达到最大值9.66 kN时，基体破坏已经扩展至层压板端部的另一边缘；
计算结束时，基体破坏和纤维破坏已经扩展到试验件的中部，此时的层压板已完全失去承载能力.
通过对[0]8、[15]20、[30]20、[45]20、[60]20、[80]20和[90]16共7种铺层进行偏轴拉伸试验，获得了相应的应力-应变曲线及拉伸强度.对不同铺层角度的断裂面角度进行测量，测得的断裂面角度与铺层角度的差异在 -1°～2°之间，对断裂面的观察和分析表明导致此差异的原因为AS4/PEEK层压板在制造过程中铺层角度和设计角度存在误差.
对于AS4/PEEK层压板在加载过程中表现出显著的非线性力学行为和破坏行为，提出一个同时考虑塑性效应和损伤累计的三维弹塑性损伤本构模型，并开发了基于Abaqus的用户材料子程序VUMAT.塑性模型中的参数使用信任域反射算法计算得到，其值为a66=1.50，β=292.67和n=0.134 6.
将三维弹塑性损伤本构模型应用于单向AS4/PEEK层压板的偏轴拉伸数值模拟，预测得到的应力-应变曲线与试验获得的应力-应变曲线吻合良好，表明该模型能够较为准确地模拟单向AS4/PEEK层压板的塑性效应.对加载过程中的各方向的塑性应变的分析表明，1-2平面内的剪切塑性应变在塑性变形中占主导地位.对比数值模拟与试验测量的拉伸强度，在相同的厚度下，随着偏轴角度的增大，拉伸强度逐渐下降，预测的拉伸强度与实际的拉伸强度的最大误差不超过20%，表明该模型适用于预测单向AS4/PEEK层合板的偏轴拉伸强度.
猜你喜欢层压板铺层塑性基于应变梯度的微尺度金属塑性行为研究九江学院学报(自然科学版)(2022年2期)2022-07-02硬脆材料的塑性域加工装备制造技术(2020年2期)2020-12-14含表面裂纹的复合材料层压板剩余强度的工程算法民用飞机设计与研究(2020年1期)2020-05-21铍材料塑性域加工可行性研究制造技术与机床(2019年12期)2020-01-06复合材料风扇叶片铺层设计方法研究航空发动机(2018年3期)2018-06-23不同工艺参数对自动铺带碳纤维层压板的性能影响纤维复合材料(2018年3期)2018-04-25典型双曲率变厚度层压板结构自动铺带工艺研究纤维复合材料(2018年3期)2018-04-25准静态压痕力作用下复合材料层压板损伤分析方法纤维复合材料(2018年3期)2018-04-25卫星天线复合材料框架的铺层优化设计北京航空航天大学学报(2017年5期)2017-11-23复合材料轴结构力学性能预测及铺层方案设计材料科学与工程学报(2016年2期)2017-01-15本文来源：http://www.zhangdahai.com/shiyongfanwen/qitafanwen/2023/0917/655704.html
