DOI: 10.11817/j.issn.1672-7207.2018.05.010
接触载荷作用下非均质材料表层应力集中分析
杨万友，周青华，王家序，杨勇，苗强
(四川大学 空天科学与工程学院，四川 成都，610065)
摘要：利用数值化等效夹杂方法分析接触载荷作用下规则分布杂质单元影响材料表层最大von Mises 应力大小和位置，并与相同工况下均质材料分析结果进行对比。在此基础上，研究不同分布参数下正态分布杂质单元对材料应力集中行为的影响。研究结果表明：在特定接触载荷作用下，均质材料表层最大von Mises 应力随摩擦因数增大而增大，但杂质单元的存在将改变最大von Mises 应力位置变化规律。非均质材料基体表层应力集中随分布杂质单元x方向平均坐标值增大而逐渐增大；随y方向平均坐标值增大而逐渐减小；随杂质单元半径增大而先增后减；分布杂质单元体积分数增大将导致基体最大von Mises应力上升。
关键词：接触载荷；非均质材料；应力集中；失效
中图分类号：TB330 文献标志码：A 文章编号：1672-7207(2018)05-1095-08
Subsurface stress concentration of heterogeneous material under contact loading
YANG Wanyou, ZHOU Qinghua, WANG Jiaxu, YANG Yong, MIAO Qiang
(School of Aeronautics and Astronautics, Sichuan University, Chengdu 610065, China)
Abstract: A numerical equivalent inclusion method was used to analyze the maximum subsurface von Mises stress and its location for heterogeneous materials with regularly distributed heterogeneity elements under contact loading. The obtained results were compared with that of a homogeneous material under the same load. Further, parametric analyses were conducted to investigate the stress concentration caused by normally distributed heterogeneities in heterogeneous materials. The results show that the maximum subsurface von Mises stress of homogeneous material under a given load increases when friction coefficient becomes larger. Location of the maximum von Mises stress changes due to the existence of distributed heterogeneity elements. Subsurface stress concentration of the matrix becomes severer if the average x coordinate value of the normally distributed heterogeneities increases, or the average y coordinate value decreases. The maximum von Mises stress of materials rises firstly and then declines as the radius of inhomogeneity increases. The increase of volume fraction in heterogeneity results in a severer stress concentration in the matrix.
Key words: contact loading; heterogeneous material; stress concentration; failure
传递运动和动力的界面往往直接决定航空、航天及机械等领域重要装备关键零部件传动性能，其中材料表层质量是影响零部件界面性能的关键因素。然而，工程材料中不可避免地存在如杂质、裂纹等缺陷[1-2]。在接触载荷作用下，这些缺陷的存在将引起材料表层甚至表面应力集中，诱发零部件破坏失效[3-4]。前期针对材料接触应力行为的研究建立在所分析材料为均质材料的理想条件之上[5]，然而材料中分布的不同程度缺陷容易导致局部应力集中[3]。针对该类问题，国内外研究者利用有限元等数值模拟方法开展了大量研究。CHOI等[6]利用考虑晶体滑移和形变孪晶的晶体塑性有限元法来模拟面内压缩过程中AZ31镁合金的空间应力集中表现，其结果表明形变孪晶在非孪晶晶粒与孪晶晶粒的晶界之间会加剧局部应力集中。冯磊等[3]利用有限元法分析非均质材料受单调拉伸和循环拉压载荷条件下夹杂物周围的局部力学行为，研究夹杂性能、位置、多夹杂、循环塑性对局部应力集中的影响。KUBAIR等[7]用数值化方法研究了材料的非均匀特性对由功能梯度材料板上圆孔产生的应力集中因子的影响，其结果表明当材料杨氏模量从圆孔处向外逐渐增大时，应力集中因子减小。CERIT等[8-10]也开展了类似研究。有限元等方法在研究非均质材料应力集中问题时，由于材料内部微结构的复杂性，往往建模困难；且随着分析网格数目的增加，有限元等方法的时间和内存消耗急剧增大，导致无法进行精确求解[11]。部分研究者尝试用解析方法来替代数值模拟方法进行该类问题研究。YU等[12]通过改进等效夹杂方法，对在外载作用下2个理想固结在一起的半无限体中嵌入椭球形杂质的弹性场进行了求解。除等效夹杂方法外，MUSKHELISHVILI[13]提出复变方法来处理 二维杂质问题的求解。该方法已被很多研究者采用[14-15]。然而，上述针对非均质材料应力场求解的解析方法只能处理简单形状微观单元问题。本文作者基于前期研究提出用于求解非均质材料弹性场的半解析方法[16]，对在线接触载荷作用下非均质材料表层应力集中问题进行研究。
1 非均质材料应力集中分析
1.1 均质材料应力集中分析
图1所示为在线接触载荷作用下的均质材料接触模型及其特性(其中σvm,max为最大von Mises应力)。首先利用接触模型(见图1(a))对均质材料应力集中行为进行研究。半径为r的载荷作用于弹性模量Em为210 GPa，泊松比νm为0.3的半无限均质材料基体上，基体表面与x轴平行。半圆形载荷在直角坐标系下沿y轴对称分布，最大正向载荷为P，最大切向载荷S可通过摩擦因数f求得，即S=f P。利用SMITH等[17]提出的算法，将摩擦因数f由0变化至1，计算可得下表面最大应力和位置变化，如图1(b)所示。
图1 线接触载荷作用下均质材料接触模型及其特性
Fig. 1 Contact model and performance of homogeneous material under line contact loading
从图1(b)可以看出：当摩擦因数为0~0.273时，均质材料下表面最大von Mises应力增长较为缓慢，但当摩擦因数大于0.273之后却呈线性增长。当f小于0.273时，最大von Mises应力位置保持在下表面，但逐步向材料表面靠近；然而当f=0.273时，最大应力位置发生突变，从下表面直接上移到表面。最大von Mises应力位置往往也称为材料的第一屈服点[18]，对于材料疲劳失效行为起着决定性作用。故当f大于0.273时，线载荷作用下均质材料的最大von Mises应力将在表面产生，极易导致表面磨损与剥落，造成接触零部件的破坏失效。
1.2 非均质材料应力集中分析
利用含多个杂质单元的非均质材料模型来研究材料非均质特性导致的材料应力集中行为及其与均质材料之间的区别。分别采用规则分布的100个圆形、正方形杂质单元模拟材料内部非均质特性，如图2(a)~(b)所示。圆形杂质半径R=0.05r，正方形杂质边长
(以保持与圆形单元面积一致)，呈10×10阵列分布，表层杂质离材料表面之间的距离h=0.1r，相邻单元之间的中心距离d=0.2r。Ei和νi分别为杂质单元弹性模量和泊松比，其余参数设置与图1(a)模型相同。利用前期研究中所提出的算法[16]，分别设定杂质单元材料参数如下：Ei=2Em，νm=νi(硬杂质)；Ei=0.5Em，νm=νi(软杂质)。分别计算2种不同形状分布单元在不同摩擦因数(0~1)下最大von Mises应力及其位置的变化，结果如图2(c)~(d)所示。
从图2(c)~(d)可以看出：与均质基体情形类似，当含规律分布杂质单元时，基体最大von Mises 应力同样随摩擦因数f增大而先缓慢增大，之后逐渐呈近似线性增长。但由于杂质单元的加入，在相同摩擦因数作用下，基体下表面最大von Mises应力都比均质材料的要大，且软杂质单元较硬杂质单元增大更显著。由于规律分布杂质单元的存在，最大von Mises 应力位置随着摩擦因数f增大，从下表面逐渐分段阶跃上移至接近表面位置(该位置由杂质阵列中表层单元位置决定)。2种圆形杂质引起的最大应力位置阶跃上移最明显的阶段分别发生在f=0.31和f=0.24处；2种正方形杂质引起的最大应力位置阶跃上移最明显的阶段分别发生在f=0.27和f=0.22处。由图2(c)与图1(b)对比可知：圆形硬质单元可相对延缓材料第一屈服点跳跃至近表面位置，而软杂质的存在则使得当摩擦因数较小时材料第一屈服点就已上移到近表面，增大材料发生破坏的风险。由图2(d)与图1(b)对比可知：基体中不管是存在正方形硬质单元还是软质单元，都使得摩擦因数较小时材料第一屈服点就已上移到近表面，加剧材料破坏的发生。
当摩擦因数f=0.4时，均质材料和含规则分布杂质单元的非均质基体应力场如图3所示(其中σvm为von Mises应力)。从图3(a)可以看出：均质材料在摩擦力作用下，主要应力分布集中于表面。而在含规则分布异质单元的非均质材料内部应力则被极大地扰动，并且在杂质单元附近易产生应力集中。从图3(b)~(e)可看出：硬、软2种杂质都会导致基体产生应力集中，且软杂质更为显著。圆形硬质单元周围的应力集中程度要明显比正方形硬质单元的低，而圆形与正方形软质单元的应力集中程度则相当。图3所示结果有效地解释了图2中的模型计算结果。
图2 规则分布杂质接触模型及其最大von Mises应力和位置
Fig.2 The maximum von Mises stresses and its locations for contact model involving regularly distributed heterogeneities
进一步利用规则分布杂质接触模型(见图2(a)~(b))研究单元材料属性变化对基体应力集中行为的影响。分别设定杂质单元与基体弹性模量比λ=Ei/Em(或λ=Em/Ei)为1，2，4和8，这2种材料泊松比相同，即νm=νi=0.3。经计算分析得到在不同摩擦因数作用下，因分布杂质弹性模量对基体最大von Mises应力影响如图4所示。由图4可知：在相同特定载荷作用下含有规则分布杂质单元的非均质材料中产生的最大von Mises 应力比均质材料的高。异质单元与基体材料属性差异越大，导致下表面应力集中更为严重。杂质形状的不同，导致下表面应力场扰动相异，使得最大von Mises应力随摩擦因数的变化幅度有所差别。
图3 不同材料在摩擦因数f =0.4时的von Mises 应力场
Fig. 3 Von Mises stress field of different materials with friction coefficient f=0.4
图4 不同摩擦因数下规则分布圆形、正方形杂质对基体最大von应力影响
Fig. 4 Influence of regularly distributed circle and square heterogeneities with different friction coefficients on the maximum von Mises stresses of matrix
2 正态分布杂质应力集中参数化研究
材料中杂质、裂纹等缺陷往往呈随机分布，因此，有必要对随机分布杂质单元影响材料应力集中行为进行研究，以期分析结果更贴近工程实际。基于上述分析结果，设定摩擦因数f=0.1，以保证材料内部最大应力位于表面之下。在此基础上，对接触载荷作用下含正态分布杂质单元的非均质材料应力响应进行参数化分析。
含正态分布杂质单元的非均质材料线接触模型如图5(a)所示。令杂质单元中心位置(x, y)分别服从
和
(
和
分别为分布杂质在x和y方向平均坐标，dx和dy分别为对应的标准差)，单元总数同样设定为100，其弹性模量与基体弹性模量比λ=Ei/Em分别为1，2，4，6和8，其余参数与图2中模型参数一致。当dx=dy=r，f=0.1时，模型对应的应力场如图5(b)~(c)所示。显然，正态分布杂质单元对基体应力场的影响与杂质单元的分布及弹性模量密切相关。因此，下面将对正态分布杂质位置、尺寸和体积分数等进行参数化分析，以掌握杂质分布参数影响材料应力集中行为的规律。
图5 正态分布杂质接触模型及其计算结果
Fig. 5 Contact model with normally distributed heterogeneities and its computation results
2.1 正态分布杂质x方向平均坐标
设定杂质单元深度位置=r、半径R=0.05r，改变，使其在-0.9r~0.9r之间变化，计算得到的基体最大von Mises应力结果如图6所示。图6(a)和6(b)分别对应Ei/Em≥1及Ei/Em≤1这2种情形，且图中应力为5组计算算例的平均值。从图6(a)可以看出：随着硬质单元x方向横坐标的增大，对基体最大von Mises 应力也随之增长，但变化幅度较小。硬质杂质弹性模量与基体弹性模量差异越大，其对基体最大von Mises 应力影响越大。由图6(b)可得到与图6(a)类似的结论，但软质单元计算结果显然比硬质单元的结果更为显著。由于本文设定摩擦力沿x轴为正方向，因此，当增大时，位于接触载荷影响区域内的异质单元数目将增加，加剧基体应力集中。
图6 正态分布圆形杂质x方向位置对基体最大von Mises应力影响
Fig. 6 Influence of x locations of normally distributed circle heterogeneities on the maximum von Mises stresses of matrix
2.2 正态分布杂质y方向平均坐标
设定圆形杂质单元横坐标位置=0、半径不变，研究在0.1r~1.9r之间变化时基体应力集中变化规律，结果如图7所示。由图7(a)可知：随着正态分布杂质沿y方向平均坐标的增大，杂质单元引起的基体最大von Mises 应力逐渐减小。增大，意味着大部分单元将远离接触表面，将削弱材料内部应力集中现象。
图7 正态分布圆形杂质y方向位置对基体最大von应力影响
Fig. 7 Influence of y locations of normally distributed circle heterogeneities on the maximum von Mises stresses of matrix
2.3 正态分布杂质单元半径
设定正态分布单元=0，=r，并保持其体积分数不变，改变其半径R(0~0.10r之间)，计算得到的基体最大von Mises应力结果如图8所示。从图8(a)可以看出：随着分布单元半径的增大，导致的基体最大von Mises 应力呈先增后减的趋势。当硬质单元半径约为0.05r时，基体产生的最大von Mises 应力达到最大，最易发生应力集中；而软杂质单元则在R=0.06r左右时应力集中最为明显。
2.4 正态分布杂质单元体积分数
设定正态分布单元=0，=r，R=0.05r，Ei/Em分别为1，2，4，6和8，其体积分数由0逐渐增大至50%。计算得到相应的基体材料内部最大应力变化如图9所示。由图9可知：在计算范围内，随着正态分布杂质体积分数增大，基体内部最大von Mises 应力逐渐增大。体积分数增大说明材料内部杂质单元数目增多，这将导致两单元之间距离减小，相互作用增强，引起更大程度的应力集中。但需注意的是，在体积分数较小时，基体内部最大应力增长较快。随着体积分数继续增加，最大应力增长幅度有所减缓。
图8 相同体积分数圆形正态分布圆形杂质半径对基体最大von Mises应力影响
Fig. 8 Influence of radius of normally distributed circle heterogeneities with the same volume fraction on the maximum von Mises stresses of matrix
3 结论
1) 在接触载荷作用下，当均质材料的摩擦因数超过一定值时，表面下最大von Mises 应力位置将跳跃到基体表面，从而使得表面磨损加剧。而下表面存在分布杂质单元则会显著影响非均质材料中最大应力位置的分布规律，即最大应力位置随着摩擦因数增大，从下表面逐渐分段阶跃至近表面杂质阵列最表层单元位置。
图9 相同半径正态分布圆形杂质体积分数对基体最大von Mises应力影响
Fig. 9 Influence of volume fraction of normally distributed circle heterogeneities with the same radius on the maximum von Mises stresses of matrix
2) 不同形状均布杂质单元对基体下表面应力扰动程度由杂质单元形状及材料的硬质或软质决定。当杂质单元为软质材料时，材料形状对基体下表面应力扰动影响区别并不明显。
3) 正态分布杂质单元与基体材料属性差异越大，更易引起材料内部应力集中，且相应的软质单元引起的基体应力增大较硬质单元的更为显著。
4) 基体最大von Mises 应力将随圆形分布颗粒x方向平均坐标的增大而增大，随y坐标增大而减小；杂质单元半径增大将引起基体最大应力先增后减；而体积分数增加将加剧材料下表面应力集中。
参考文献：
[1] 王雪明, 谢富原, 李敏, 等. 热压罐成型加筋板L形筋条纤维密实影响因素研究[J]. 航空材料学报, 2011, 31(4): 64-68.
WANG Xueming, XIE Fuyuan, LI Min, et al. Experimental research on fiber compaction of L-shaped stiffeners in stiffened skins by autoclave process[J]. Journal of Aeronautica Materials, 2011, 31(4): 64-68.
[2] 顾轶卓, 张佐光, 李敏, 等. 复合材料变厚层板热压成型缺陷类型与成因实验研究[J]. 复合材料学报, 2008, 25(2): 41-46.
GU Yizhuo, ZHANG Zuoguang, LI Min, et al. Experimental study on type and cause of defects in variable thickness composite laminates during hot pressing process[J]. Acta Materiae Compositae Sinica, 2008, 25(2): 41-46.
[3] 冯磊, 轩福贞. 非金属夹杂物对材料内局部应力集中的影响[J]. 机械工程学报, 2013, 49(8): 41-48.
FENG Lei, XUAN Fuzhen. Effect of non-metallic inclusions on the local stress concentration within materials[J]. Journal of Mechanical Engineering, 2013, 49(8):41-48.
[4] 陈康华, 方玲, 李侠, 等. 颗粒失效对SiCp/Al复合材料强度的影响[J]. 中南大学学报(自然科学版), 2008, 39(3): 493-499.
CHEN Kanghua, FANG Ling, LI Xia, et al. Influence of particle failure on strength of SiCp/Al composites[J]. Journal of Central South University(Science and Technology), 2008, 39(3): 493-499.
[5] SADEGHI F, JALALAHMADI B, SLACK T S, et al. A review of rolling contact fatigue[J]. Journal of Tribology, 2009, 131(4): 220.
[6] CHOI S H, KIM D H, PARK S S, et al. Simulation of stress concentration in Mg alloys using the crystal plasticity finite element method[J]. Acta Materialia, 2010, 58(1): 320-329.
[7] KUBAIR D V, BHANU-CHANDAR B. Stress concentration factor due to a circular hole in functionally graded panels under uniaxial tension[J]. International Journal of Mechanical Sciences, 2008, 50(4): 732-742.
[8] CERIT M, GENEL K, EKSI S. Numerical investigation on stress concentration of corrosion pit[J]. Engineering Failure Analysis, 2009, 16(7): 2467-2472.
[9] ACHOUR T, BOUIADJRA B B, SERIER B. Numerical analysis of the performances of the bonded composite patch for reducing stress concentration and repairing cracks at notch[J]. Computational Materials Science, 2003, 28(1): 41-48.
[10] 徐慧, 伍晓赞, 程仕平, 等. 复合裂纹的应力强度因子有限元分析[J]. 中南大学学报(自然科学版), 2007, 38(1): 79-83.
XU Hui, WU Xiaozan, CHENG Shiping, et al. Finite element analysis of stress intensity factor in composite mode crack[J]. Journal of Central South University(Science and Technology), 2007, 38(1): 79-83.
[11] ZHOU Qinghua, JIN Xiaoqing, WANG Zhanjiang, et al. Numerical modeling of distributed inhomogeneities and their effect on rolling contact fatigue life[J]. Journal of Tribology, 2012, 137(1): 313-315.
[12] YU Jiahuan, KUANG Zhenbang. The stress analysis of an ellipsoidal inhomogeneity in dissimilar media[J]. Composites Science and Technology, 2003, 63(7): 955-966.
[13] MUSKHELISHVILI N I. Some basic problems of the mathematical theory of elasticity[J]. The Mathematical Gazette, 1953, 48(365): 351.
[14] KUSHCH V I, SHMEGERA S V, BURYACHENKO V A. Interacting elliptic inclusions by the method of complex potentials[J]. International Journal of Solids and Structures, 2005, 42(20): 5491-5512.
[15] BURYACHENKO V, KUSHCH V. Effective transverse elastic moduli of composites at non-dilute concentration of a random field of aligned fibers[J]. Journal of Applied Mathematics and Physics, 2006, 57(3): 491-505.
[16] 杨万友, 周青华, 王家序, 等. 考虑异质材料的线接触性能建模与分析[J]. 复合材料学报, 2016, 33(8): 1848-1858.
YANG Wanyou, ZHOU Qinghua, WANG Jiaxu, et al. Modeling and analysis on line contact performance considering heterogeneous properties of material[J]. Acta Materiae Compositae Sinica, 2016, 33(8): 1848-1858.
[17] SMITH J O, LIU C K. Stresses due to tangential and normal loads on an elastic solid with application to some contact stress problems[J]. Journal of Applied Mechanics-Transactions of the ASME, 1953, 20(2): 157-166.
[18] LI Xiaoyong, WANG Shilong, WANG Zhanjiang, et al. Location of the first yield point and wear mechanism in torsional fretting[J]. Tribology International, 2013, 66(7): 265-273.
(编辑 伍锦花)
收稿日期：2017-05-20；修回日期：2017-06-30
基金项目(Foundation item)：国家自然科学基金资助项目(51435001，51405316)；中央高校基本科研业务费资助项目(2017SCU12021)；航空科学基金资助项目(20150219001) (Projects(51435001, 51405316) supported by the National Natural Science Foundation of China; Project (2017SCU12021) supported by the Fundamental Research Funds for the Central Universities of China; Project(20150219001) supported by Aeronautical Science Foundation of China)
通信作者：周青华，博士，副教授，从事空间摩擦学与可靠性工程、多尺度优化设计理论研究；E-mail: qh.zhou@foxmail.com
