晶向对铝合金构件微动疲劳初始起裂的影响研究

王 建, 阮星源

(西南科技大学 土木工程与建筑学院, 四川 绵阳 621010)

微动疲劳是在材料疲劳和接触部位微动磨损协同作用下,构件接触界面附近发生的渐进损伤及开裂过程。微动疲劳现象广泛存在于航空工业领域内各种机械构件中,加速结构构件接触表面及浅表层裂纹的萌生与扩展,进而将造成发动机涡轮叶片榫头断裂、机身蒙皮脱落等事故,是复杂应力环境中长期服役的民用航空飞行器灾难性事故的主要原因[1]。由于微动疲劳构件的相对位移幅度很小(微米级)、造成裂纹萌生及断裂的循环次数很大(1万次以上的高周疲劳)、裂纹发生的位置很隐秘,且裂纹发生后,难以采取有效措施抑制其继续发展,因而微动疲劳又被视作难以处置的“工业癌症”[2]。

超高强铝合金AA7075作为一种航空铝材,广泛用于飞机的机身、机翼等各个部件中,并以铆钉、螺栓、薄板等形式存在[3]。在航空飞行过程中,铝合金板材和构件将面临着微动疲劳导致的结构性能劣化和断裂失事风险[4]。微动疲劳可分为初始起裂和裂纹扩展直至断裂两个阶段;而初始起裂部位的辨识及裂纹萌生寿命的预测,在铝合金微动疲劳研究中,具有重要的现实意义和迫切的实际需求。

微动疲劳影响因素可分为接触状态(接触几何形态、表面状态及受力状态等)、环境因素(温度、湿度及环境介质等)和材料特性(材料成分、结构及宏细观力学特性等)等三个方面[5]。在不改变接触状态和环境因素等条件下,材料特性作为微动疲劳的唯一可控因素,引起了国内外相关研究者的广泛关注。研究者发现,可以通过不同的材料加工方法及表面处理技术改进接触构件表层材料的性状以提高其微动疲劳抗力[6],进而要求确立微动疲劳与材料细观织构之间的关系,以更好地进行抗微动疲劳设计。而在材料加工与设计领域,对于材料细观织构对构件微动疲劳性能的影响研究尚在探索中。Minaii 等[7]利用晶体塑性有限元研究了晶粒尺寸和规则参数等细观结构特征参数对初始起裂寿命及起裂部位的影响。 McCarthy等[8]利用晶体塑性有限元研究了晶向对316LSS不锈钢构件微动疲劳导致的磨损、塑性形变及应力状态的影响。国内相关研究方面,包括田红亮等[9],关海达等[10]和刘丽艳等[11]分别利用室内试验、数值计算对结构部件微动部疲劳、磨损特性及机制所做的探索工作。当前,尚未见晶向对铝合金构件微动疲劳影响的相关报道。

本文拟利用晶体塑性有限元模拟铝合金试件在循环荷载作用下的微动疲劳过程,探讨晶向对铝合金晶体试样初始起裂部位及裂纹萌生寿命的影响情况。文章主体内容分布于1～3章:第1章为微动疲劳模型设计;第2章为晶体塑性有限元模拟方法介绍及材料参数选择,加载方式及荷载选择,以及接触模型验证等;第3章为晶向对晶体试样弹模及屈服强度的影响、晶向对初始起裂部位及裂纹萌生寿命的影响分析。

1 微动疲劳模型

服役期间交变荷载作用下飞机机身铆钉-平板结构单个搭接接头可简化为图1所示的微动疲劳计算模型,这也是常见的单夹具微动疲劳试验的代表。平板接头所受的横向远程振动荷载由同步、同相加载的轴向循环应力σA和轴向反作用应力σR分别施加于试样两侧。

在微动垫模型上部边界施加法向荷载P, 并采用多点约束(MPC),以模拟铆钉的平稳挤压作用。在微动疲劳工况下,微动垫-试样接触保持部分滑动状态[12],即试样顶部所受的切向作用力Q满足条件Q满足条件Q<μP,μ为摩擦因数,Q可表示为σA和σR的函数

(1)

式中,As为试样的横截面积。

为克服单元尺寸差异过大对数值计算所造成的困难,以实施微动疲劳宏-细观耦合多尺度计算,我们采用了全局模型-子模型耦合的方法计算交变荷载作用下试样晶粒内部的弹塑性响应。计算模型尺寸依据Vazquez 等[13]室内微动疲劳试验模型尺寸取值。全局模型中试样长度和宽度分别为 30 mm 和 5 mm。 微动垫接触面圆弧部分半径为 100 mm。除了圆弧部分外,微动垫模型的宽度和高度都被截断为 20 mm,以满足圣维南原理中关于边界效应的消除原则。选择三节点线性平面应变三角形单元 CPE3对微动垫和试样进行网格划分。为了解决试样尺寸、接触半宽尺寸的受限问题,我们在晶体塑性有限元微动疲劳模拟中对试样采用了子模型-全局模型耦合技术,如图1所示。子模型与微动垫毗邻,其尺寸为10 mm × 2.5 mm。子模型由NEPER软件[14]随机生成741颗直径为150 μm的晶粒。晶粒同样由CPE3单元剖分,单元长度设为50 μm,其单元总数为9 968个。以每一时步全局模型-子模型界面处全局模型位移计算结果Ux和Uy作为子模型左右侧和底部动力边界条件;子模型与上方微动垫的接触特性及加载方式均与全局模型一致。在子模型所在区域范围内,试样全局模型单元长度为0.1 mm;在子模型所在区域范围外,试样全局模型单元长度逐渐由0.1 mm增加到左、右边界处的0.2 mm。微动垫模型单元长度由接触面上的0.05 mm逐渐过渡到顶面的0.5 mm。试样模型和微动垫模型分别划分了23 800个和39 682个CPE3单元。微动垫和试样这两个组件的材料均为Al 7075-T651 合金。对于弹性接触而言,试样和微动垫的力学参数均设置为杨氏模量E= 71 GPa 和泊松比ν= 0.33。就接触设置而言,微动垫的圆弧表面定义为接触副的主表面,而试样上表面定义为从表面。 采用拉格朗日乘子摩擦公式描述接触区域的切向行为,根据试验测试结果,其摩擦因数μ=0.72。 此外,对于几何特性而言,平面外表面厚度设置为 7 mm。

2 晶体塑性有限元方法

2.1 晶体塑性公式

晶体塑性有限元法可以模拟循环荷载作用下晶粒内部滑移带错动所引发的微观塑性变形,因而可用于模拟微动疲劳现象。本文选用了Marin 等[15]开发的解释弹性和塑性变形的晶体弹黏塑性本构模型,并利用 Armstrong-Frederick 运动硬化规则对其进行了修改,以适应循环背应力(backstress)的作用。该本构模型假设晶体滑移是塑性变形的唯一机制,而忽略其他机制,例如孪晶、晶界滑动和扩散。 单晶的弹塑性变形梯度F可以乘法分解为弹性分量Fe和塑性分量Fp,弹性变形梯度Fe也可以分解为对称弹性拉伸张量Ve和正交旋转张量Re,即弹塑性变形梯度F可以表示为如下形式

F=FeFp=VeReFp

(2)

塑性速度梯度Lp定义为所有滑移系上的剪应变速率之和[16],即

(3)

第α个滑移系上的剪切应变率用式(4)所示的经验幂律方程流动法则[17-18]计算

(4)

(5)

2.2 材料特性和晶向参数校正

为了与微动疲劳模型一致,本文拟采用二维代表性体积单元(2D RVE)进行参数识别和校正。2D RVE具有两个基本假设:平面应变条件和二维尺度。Vuppala 等[20]通过研究发现,平面应变假设对织构演变几乎没有影响。此外,2D模型的部分缺陷可通过增加晶粒数量克服[21];并且2D模型计算结果也可映射到3D模型[22]。

本文所选用的AA7075晶体本构关系所需参数源于先前工作基础[23]。晶体织构直接基于电子背散射试验(electron back scatter diffraction, EBSD)试验结果,而结晶本构方程的参数通过反分析识别。 弹性材料 C11、C12 和 C44 采用 Hosford试验结果[24]。主要参数如表1 所示。

表1 AA7075材料晶体本构模型参数Tab.1 Calibrated parameters of the crystalline constitutive law for AA7075

除了晶粒形态外,晶体各向异性的其他主要来源是织构或晶向。晶粒初始织构通过电子背散射衍射(EBSD)测试获得,并导入到晶体塑性有限元模型中。如图2所示,本文选取了AA7075材料试验拟合晶向(简称试验晶向),随机分布晶向(简称随机晶向)和随机偏转±10°的(111)晶向(简称(111)晶向)进行对比分析。其中,试验拟合晶向是对Li等[25]电子背散射测试(EBSD)结果的拟合。

(a) 电子背散射测试结果

(b) 试验拟合晶向

(c) 随机晶向

(d) (111)晶向图2 AA7075晶体织构赤平极射投影图Fig.2 Pole figures of AA7075

2.3 模型验证

根据Hertz弹性接触理论[26],微动垫-试样接触面的归一化接触应力如式(6)所示

(6)

式中:p0为接触表面的最大应力;a为接触半宽。

本文所用的二维模型微动垫-试样相互接触与圆柱体与半平面之间的接触等效。根据Nowell等[27]的研究成果,圆柱体与半平面接触面上的归一化剪切牵引力如式(7)所示

(7)

式中:c为黏着区的接触半宽。c由式(8)所示

(8)

e为黏着区中心与接触中心之间的偏移量,可由式(9)所示

(9)

图3为微动垫-试样模型接触界面处归一化接触压力和归一化剪切牵引力的解析解和 FEM 解之间的比较。Nowell等[27]研究发现,当试样的半宽大于或等于接触半宽的10倍时,则数值模型满足Hertz解的假设条件;当试样半宽与接触半宽之比等于3时[28],数值计算结果与解析解之间的偏差高达20%。尽管本文试样模型的宽度不满足完全弹性半无限空间假设,但归一化接触压力和归一化剪切牵引力的解析解和数值解之间吻合得很好,表明本文所采用的模型获得了很好的验证。

(a) 归一化接触压力

(b) 归一化剪切牵引力图3 解析解与有限元数值解的比较Fig.3 Comparison of analytical and FEM solutions

2.4 加载方式

根据时序关系,微动疲劳模拟过程可分为两个不同的阶段:初始加载阶段和微动循环加载阶段,如图4所示。在初始加载阶段(即图4中的时步0),逐渐将法向荷载P施加到设定值,使微动垫和试样之间建立稳定的接触关系,随后在保持P值恒定。在微动循环加载阶段,在试样端部同步施加相位差为0的轴向体积应力σA和轴向反作用应力σR。循环荷载加载过程中,轴向应力比(最小轴向应力与最大轴向应力之比)和切向载荷比(最小Q与最大Q之比)均等于-1。为了将这些加载预设值准确发送到 Abaqus求解器,每个循环划分为四个阶段:加载阶段(如时步1)、卸载阶段(如时步2)、反向加载阶段(如时步3)和反向卸载阶段(如时步4)。为与微动疲劳试验一致,循环加载频率设置为 10 Hz。

图4 轴向应力σA和轴向反力σR相位差为0时的一个完整循环的加载次序Fig.4 Loading sequence of one complete cycle for 0° phase difference between σA and σR

在本文的计算中,我们模拟了300个循环次数的微动疲劳作用,加载参数设定值如表2所示。

表2 微动疲劳模拟加载工况Tab.2 Loading conditions selected for fretting fatigue simulation

3 晶体塑性有限元方法

3.1 晶向对晶体弹模及屈服强度的影响

不同晶向试样单轴拉伸应力应变模拟曲线如图5所示。其中,试验晶向试样,随机晶向试样和(111)晶向试样的赤平极射投影图分别见图2(b)～(d)。单轴拉伸试样采用图1中的子模型尺寸,单轴拉伸速率与微动疲劳模拟中循环荷载轴向加载速率相同,即应变率均为0.001 s-1。

图5 单轴拉伸条件下晶体试样应力应变曲线Fig.5 Stress-strain simulation curves of specimens with different crystallographic orientations under uniaxial tension

在弹性域内,试验晶向单轴拉伸模拟所得的弹性模量为69.4 GPa, 与Boller等[29]的试验结果71 GPa相差仅约2%; (111)晶向试样的弹性模量最大,为77.5 GPa; 随机晶向试样的弹性模量居中,为76.4 GPa。此外,试验晶向试样屈服强度和硬化强度均最小; (111) 晶向试样的屈服强度大于随机晶向试样,而其硬化强度低于随机晶向试样。由此可见,(111) 晶向试样的宏观弹模及屈服强度最大,随机晶向试样次之,而试验晶向试样最小。

3.2 晶向对微动疲劳初始起裂部位的影响

300次微动疲劳循环加载后,试验晶向, 随机晶向和 (111) 晶向这三种晶向试样等效塑性应变云图,如图6所示。三种晶向的接触面初始起裂部位均发生于接触后缘附近的E1714单元或与E1714单元相邻的E5993单元内部。由此可见,晶向分布对接触面初始起裂部位影响不大。

(a) 试验晶向

(b) 随机晶向

(c) (111)晶向图6 300次循环加载后微动疲劳等效塑性应变极值点位置Fig.6 The location of fretting fatigue PEEQ extreme point after 300 cyclic loadings

三种晶向试样内部均存在多个等效塑性应变(equivalent plastic strain, PEEQ)极值所对应的初始起裂热点。PEEQ越大的单元,其所在的部位越可能开裂。因此,可以通过不同晶向试样的等效塑性应变分布云图预判微动疲劳初始起裂点位置[30]。不同晶向试样内部微动疲劳初始起裂点位置识别如下:试验晶向试样的初始起裂点位于微动垫-试样接触中线附近,接触面以下1.28 mm的E2608单元处;随机晶向试样的初始起裂点位于接触面以下1.41 mm的E1495单元处;(111) 晶向试样的初始起裂点位于接触面以下1.56 mm的E1837单元处。结合不同晶向试样的宏观屈服强度可知,试样弹模及屈服强度越大,其初始起裂部位越深。就PEEQ分布形态而言,三种晶向试样均表现出不同程度的各向异性。试验晶向试样PEEQ云图呈倾斜条带状分布,(111) 晶向试样PEEQ云图呈絮状分布,而随机晶向试样PEEQ云图介于条带状与絮状分布之间。可以预见,PEEQ分布形态将影响到裂缝发轫及扩展特征。由此可见,不同晶向试样,其亚表面裂缝发轫位置及扩展特征均有所不同。

3.3 晶向对塑性应变能密度及微动疲劳裂纹萌生寿命的影响

晶系滑移过程中所耗散的能量可用塑性应变能密度(PSED)ϖ衡量[31],其单次循环作用过程中所有滑移系所耗散的能量可由式(10)所计算

ϖ=∮cσ:dεp=∑∮ταdγα

(10)

式中:σ为有效应力,εp为晶系滑移产生的累积塑性应变;τα和γα分别为作用于滑移系α的剪应力和塑性剪应变。根据以往计算发现,单个循环的塑性应变能密度与其累积值成正比关系;且在循环加载过程中,单次循环作用所耗散的塑性应变能会随着加载的进展发生细微的变化,因此,本文仅考虑累积塑性应变能密度时程曲线的分布特征。

此外,初始起裂寿命N与ϖ之间的关系可由式(11)表示

(12)

式中,f(x)为能量吸收系数、最小工程可检裂纹长度、平均位错间距,以及弹性模量、泊松比等参数的函数。由此可见,塑性应变能密度与初始起裂寿命成反比关系。

不同晶向试样接触面和试样内部初始起裂点的总累积塑性应变能密度时程曲线(plastic strain energy density, PSED),如图7所示。对于接触面初始起裂部位和模型内部初始起裂部位而言,试验晶向试样的PSED曲线斜率及所对应的单次循环塑性应变能密度最大,(111)晶向试样的PSED曲线斜率及所对应的单次循环塑性应变能密度最小,因此,试验晶向试样的初始起裂寿命最小,(111) 晶向试样的初始起裂寿命最大,随机晶向试样的初始起裂寿命居中。此外,对于同一种晶向试样,其接触面初始起裂部位PSED均小于模型内部初始起裂点PSED,表明在当前采用的荷载和接触条件下,模型内部更易于起裂。综上所述,晶体试样弹模和屈服强度越大,其微动疲劳初始起裂寿命越大,即试样初始起裂寿命与材料强度成正相关关系。

(a) 接触面初始起裂点

(b) 内部初期起裂点图7 晶向试样累积塑性应变能密度时程曲线Fig.7 Time history of accumulated plastic strain energy density

4 结 论

(1) 晶向分布将对晶体材料的弹性模量和屈服及硬化强度造成明显影响。(111)晶向试样的弹性模量最大和强度最大,随机晶向试样的弹性模量和强度居中,而一般试验晶向试样的弹性模量和强度均最小。

(2) 晶向分布对微动疲劳接触面起裂部位的影响不明显,而对亚表面起裂部位影响比较显著。(111)晶向试样的亚表面初始起裂点最深,随机晶向试样次之,而一般试验晶向试样的亚表面初始起裂部位最浅。

(3) 晶向分布对微动疲劳初始起裂寿命的影响比较明显。(111) 晶向试样的初始起裂寿命最大,随机晶向试样的初始起裂寿命居中,而试验晶向试样的初始起裂寿命最小。

上述结论表明,晶向分布这一细观织构影响了晶体材料宏观尺度的弹性强度特征和塑性变形特征,并对金属构件微动疲劳初始起裂寿命、亚表面裂缝发轫位置及扩展特征均造成了影响。利用材料加工方法改变金属晶体试样浅表层的晶向分布状态,将可以抑制微动疲劳裂缝的发轫及扩展。