晶粒梯度對AA04-T35合金微動疲勞裂紋萌生的影響研究

中圖分類號：0346.0；TH117 文獻標志碼：A

Abstract：Accurateidentification of influencing factors of freting fatigue is crucial for itsassessment and control.To addressthe fracture problem of AA2O24-T351 alloy under freting fatigue，this study utilizes a crystal plasticity finite element modelthatincorporates microstructuralsensitivityto investigate the impactof graingradient representation on freting fatigue in the AA2024-T351 alloy.Additionally，a submodel approach is employed to calculatethecontact regionofthefreting fatigue specimens，andthecumulativeplastic slip is analyzed todetermine thecrack initiation location and predict crack initiation life.The study reveals that the grain gradient results in an uneven distribution of stress fields within metallic materials due to variations in grainsizesand grain boundary conditions.Compared tothefretingfatiguecontactsurface，thegrain gradient exhibitsasignificant change in crack initiationlocations in thesubsurface.Moreover，thegrain gradient evidentlyaffects thecrack initiation lifeoffreting fatigue.Theresearch findings indicate that thedistribution of graingradient alters theplasticdeformation characteristics between grains.As thegrain gradient increases，the subsurfacecrack initiation life gradually surpasses that of the contact surface.The variation in grain gradient significantly influences the subsurface crack nucleation locations，thereby determining the tendencies of crack propagation.Building upon theseresearch conclusions，modifying the grain gradient distribution through material processng methods effectively mitigates the crack initiation and propagation of fretting fatigue.

Key Words：fretting fatigue;crack initiation;grain gradient;crystal plasticity;finite element method

微動疲勞指在循環載荷作用下，構件中的兩個不同部件在接觸表面發生微小的相對位移，從而引發材料損傷的破壞過程1.在航空航天、能源、交通等工程領域，微動疲勞是造成構件失效的重要原因之一[2-5].在工程結構中，例如輪轂連接、大跨度橋梁的鉸鏈和螺栓連接件等，都會發生微動疲勞現象.已有研究表明，在微動疲勞中，構件失效前的相對位移幅度非常小，裂紋的形成往往發生在材料未達到疲勞極限的情況下.此外，裂紋的萌生位置往往難以察覺，一旦裂紋形成，幾乎沒有有效的手段來阻正其進一步擴展，因此，微動疲勞會降低構件的疲勞強度，加速構件的疲勞失效，從而引起災難性后果.

AA2024-T351合金是一種具有良好的強度、韌性和抗疲勞性能的鋁合金.這些優異的性能使得它在航空航天、汽車、軌道交通以及結構工程等領域得到廣泛應用.它常被用于制造飛機結構件、飛機翼梁、車身構件、車輪、橋梁支撐等.在這些高強度和輕量化的應用中，AA2024-T351合金構件面臨著微動疲勞帶來的結構性能下降和斷裂失效的風險[7-9].微動疲勞全過程壽命共經歷裂紋萌生、裂紋擴展、構件斷裂三個階段.其中，疲勞主要取決于裂紋萌生，因此，在AA2024-T351合金微動疲勞研究中，準確識別裂紋萌生形核點位置并預測裂紋萌生壽命對于提高結構可靠性和延長使用壽命至關重要.

從損傷機理的角度來看，微動疲勞的影響因素可以分為材料因素、力學因素和磨損狀況因素[10].在保持力學因素和磨損狀況因素條件不變的情況下，材料因素成為微動疲勞的唯一可控因素，引起了國內外研究者的廣泛關注.通過相關研究發現，采用合適的材料加工方法和表面處理技術能夠對接觸構件表層材料的性質產生積極影響，進而提高其微動疲勞抗力[1].因此，我們需要建立微動疲勞與材料細觀結構（晶粒尺寸、晶粒形態和晶粒取向等）之間的關系，以便更好地進行抗微動疲勞設計.

材料的微觀結構演變是一個涉及多種機制的復雜過程，利用晶體塑性有限元法來預測材料的結構演變，可以在工業生產之前預知材料在變形過程中微觀組織的變化，從而實現優化和控制工藝參數的目的，并最終獲得所需的理想內部結構[12].目前，已對材料微觀結構中晶粒尺寸和取向對構件微動疲勞性能的影響展開了研究.McCarthy等[13]采用晶體塑性有限元進行了關于晶體織構對微動疲勞裂紋萌生壽命影響的研究.Zhang等[14]利用晶體塑性有限元研究了Ti-6A1-4V合金中晶粒形貌、晶粒尺寸和晶粒分布等微觀結構特性對微動疲勞行為的影響.國內學者林恩強等[15]、丁俊等[6]和趙杰江等[7]分別利用數值模擬、室內試驗對構件微動疲勞裂紋萌生、微動機制及磨損速率進行了相關探索.晶粒梯度作為晶粒尺寸和取向的空間變化，在材料微觀結構上呈現復雜性，直接影響材料的整體力學性能和疲勞行為.晶粒梯度的不均勻分布導致局部應力集中和裂紋萌生位置的變化，顯著影響微動疲勞裂紋的萌生和擴展過程[18]，因此，深入研究晶粒梯度對微動疲勞的影響，有助于全面理解微動疲旁機理，提高構件的抗疲勞性能，并為工程實踐提供材料加工和設計指導.然而，在現有的研究中，對晶粒梯度在材料微動疲勞的整體性能影響方面尚未有系統的研究.

本文采用晶體塑性有限元模擬方法，對AA2024-T351合金在循環載荷下的微動疲勞過程進行了分析，研究探討了晶粒梯度對AA2024-T351合金試樣的應力場分布、裂紋萌生位置和裂紋萌生壽命的影響，通過累積塑性滑移確定微動疲勞裂紋的形核位置，并進行裂紋萌生壽命的預測.

1微動疲勞試驗模型

本研究基于Hojjati-Talemi等[19]的微動疲勞實驗數據，采用有限元軟件ABAQUS進行分析.圖1是微動疲勞試驗裝置示意圖.由于試驗裝置和加載方式具有軸對稱性，可將裝置簡化一半進行計算.微動疲勞有限元模型如圖2所示.微動墊和試樣均由AA2024-T351合金構成，彈性模量 E 為 74.1GPa ，泊松比 v 為0.3.墊片兩側水平位移 U_x=0 ，試樣底部豎向位移 U_y=0. 墊片上表面施加法向載荷 F ，采用多點約束（MPC）確保各點之間位移響應均勻.在試樣兩側施加周期性變化的遠程振動載荷、軸向作用力 σ_axial 和軸向反作用力 σ_reaction 根據法向載荷 F 、切向作用力 Q 和接觸表面之間的摩擦系數 μ ，疲勞磨損的滑動狀態可分為三種類型：整體滑移、部分滑移和混合滑移[20].由文獻[19]可知，本文涉及的微動疲勞屬于部分滑移，摩擦系數 μ 為0.65.采用壓力過閉合方法捕捉接觸表面法向行為，Lagrange乘子定義切向行為.將微動墊底表面設置為接觸對主表面，而試樣上表面定義為從表面.微動墊和試樣均使用CPE4R單元進行網格劃分.

圖1微動疲勞試驗裝置示意圖

為提高計算效率降低成本，本文對試樣接觸區域采用子模型方法進行材料多晶建模.文獻[19]中接觸半寬 α_a 為 462μm ，子模型半寬應超過 2a ，設為1mm. 子模型深度遵循圣維南原理2]，距接觸表面盡可能遠，設為 0.5mm. 研究表明，當單元面積與晶粒面積比值為0.01時，可以保證計算精度并實現收斂[22].據此，定義單元尺寸為 6μm ，子模型共包含27639個單元.使用NEPER軟件隨機生成多晶體[23-24]，如圖3所示，子模型中的不同顏色對應不同的晶粒.

圖2微動疲勞有限元模型的邊界和單元網格尺寸

圖3試樣接觸區域子模型示意圖

Fig.3Illustrationof the submodel of the specimen at the contact region

2晶體塑性有限元方法

2.1晶體塑性模型

用于制造微動疲勞試樣的AA2024-T351合金具有面心立方（FCC）結構，共包含12組滑移系[25].本文采用基于Huang等[26開發的晶體彈塑性本構，并對

與速率相關的黏塑性硬化規律進行了調整，以滿足微動疲勞循環硬化效應，

考慮到金屬材料的大變形，用變形梯度來描述晶體變形方程.將變形梯度張量 F 分解為彈性分量F^e 和塑性分量 F^p ，可由式（1）表示：

F=F^eF^p

速度梯度 L^p 是由晶體位錯滑移引起的，它是所有滑移系剪切速率之和，由式（2）表示：

式中： n 為滑移系的總數； 是第 α 滑移系的剪切應變率； s 和 m 分別為滑移平面內的滑移方向和法向方向.

剪切應變率表示如式（3）：

式中： 是初始剪切應變； τ^α 是滑移系臨界剪切應力； n 是應變速率敏感系數[27-28]；為了考慮疲勞循環加載下的運動硬化效應，引入了背應力 χ^α 采用Armstrong-Frederick運動硬化定律來描述背應力的演化方程[29]，其定義如式（4）：

式中 b 和 r 分別是硬化常數和動態恢復系數.

g^α 是滑移系當前強度，可表示為：

式中： h 為硬化模量，反映了變形過程中由位錯增殖相互作用引起的滑移阻力的增加情況.

潛硬化 h_αβ 表示不同滑移系之間交互硬化行為，而自硬化 h_αα 表示單個滑移系硬化情況.分別見式（6和式（7）：

式中： τ_o 和 τ_s 分別表示初始和最大的滑移系強度； h₀

是初始硬化強度； q 是潛硬化系數，取值為1.

是所有滑移系累積滑移量，即：

2.2晶體材料參數校準

為獲得準確的微動疲勞宏一微觀耦合多尺度計算結果，需要建立與宏觀應力-應變關系一致的晶體材料參數.本文采用“試錯法”，對單位代表體積元（RVE）模型進行參數擬合和校正.為了解決循環加載下多晶應力場分布不均勻的問題，對RVE模型進行均質應力 σ_ij 和應變 ε_ij 處理[30]，見式（9）和式（10）：

式中： V 為RVE的體積.

本文使用的AA2024-T351合金晶體材料參數源于先前的研究工作[31]，主要參數見表1.

表1AA2024-T351合金晶體材料參數Tab.1CrystallinematerialparametersforAA2024-T351 alloy

2.3晶粒梯度子模型

為了驗證子模型方法的合理性，本文對微動疲勞全局模型接觸區域的VonMises分布與子模型進行了比較.如圖4所示，全局模型接觸區域與子模型應力分布一致，表明子模型方法的準確性.

圖4全局模型和子模型的VonMises分布Fig.4Distribution of Von Mises stress for global modeland submodel

根據文獻[32]中AA2024-T351合金晶粒尺寸分布映射，確定平均晶粒尺寸 d 為 40μm 為了研究晶粒梯度對AA2024-T351合金微動疲勞的影響，本文建立了三組晶粒梯度子模型，具體參數如表2所示.

表2晶粒梯度子模型尺寸

Tab.2Thesizeforgraingradientof submodels

圖5表示三組晶粒梯度子模型，晶粒尺寸沿子模型上表面逐漸增大到下表面，在中間存在一個粒徑分布的過渡區，該區域同時存在小粒徑和大粒徑.

圖5晶粒梯度子模型

Fig.5Realisations for grain distribution of the submodel

3模擬結果及分析

3.1晶粒梯度對應力場的影響

三組晶粒梯度子模型的應力場分布如圖6所示，圖4（a是與子模型對應的微動疲勞全局模型應力分布.三組子模型的應力分布趨勢與全局模型接觸區域一致，沿中軸線對稱，最大應力位于子模型中心靠近下表面，并向外逐漸擴展.通過對比三組子模型的應力集中區域，觀察到應力分布不均勻.不同晶粒尺寸對子模型應力場產生不同的影響，隨著晶粒梯度增大，高應力點的分布變得更加離散，同時子模型應力場梯度也隨之增大.這個現象歸因于晶體各向異性，且塑性變形主要發生在大尺寸晶粒的晶界處.在子模型的高應力梯度工況下，疲勞過程區域內的應力分量大小和相對貢獻將迅速變化，增加了裂紋萌生和擴展風險.

圖6接觸區域VonMises分布

Fig.6Distribution of Von Mises stress for the contact region

如圖7（a所示，三組子模型的上表面應力分布與全局模型基本一致.然而，在接觸后緣附近，因晶粒織構而存在輕微差異，其中 Φ_a 為微動疲勞接觸半徑，取值為 0.477mm 總體來說，晶粒梯度對微動疲勞接觸表面應力分布影響較小.圖7（b）是子模型在不同深度處的應力分布，從子模型上表面中點 （x=0 y=0 ）開始，沿中軸線一直延伸至下表面.結果顯示，在較淺的深度范圍（ ylt;0.15mm 內，各組晶粒梯度子模型應力分布基本一致；在深度為 0.15～0.45mm 之間，各子模型應力分布出現明顯差異，這表明晶粒梯度對子模型近下表面深度范圍內的應力分布有顯著影響.微觀結構中晶粒的各向異性導致應力的波動，且通常在晶界處應力值最大.在子模型No.1中，晶粒數量較其他子模型更多，晶界分布密集，因此在晶界附近的應力波動最為活躍[33].相反，子模型No.3中的晶粒數量最少，晶界分布離散，因此晶界附近的應力波動不明顯.

圖7接觸區域應力分布

Fig.7Stress distribution for the contact region

3.2晶粒梯度對微動疲勞裂紋萌生位置的影響

微動疲勞裂紋萌生位置可分為兩種類型：接觸表面和亞表面34.然而，由于在試驗過程中觀察亞表面裂紋萌生形核點存在困難，往往會忽視這種類型的裂紋.文獻[19]研究表明，微動疲勞裂紋萌生位置在接觸界面后緣，對應等效多軸損傷應力 σ^* 最大值，因此，可以通過晶粒梯度子模型的 σ^* 預測微動疲勞裂紋萌生位置.等效多軸損傷應力 σ^* 可由三軸因子R_v 和VonMises等效應力 σ_eq 表示為：

其中， σ_H 是體積應力.

圖8是三組晶粒梯度子模型等效多軸損傷應力分布，其中， σ_y 為材料的屈服強度，取值為 383MPa^[19] 如圖8（a）所示，子模型的微動疲勞裂紋萌生位置幾乎相同，均位于接觸表面后緣，與文獻[19]結論一致；圖8（c）～（e）分別是不同晶粒梯度子模型中等效多軸損傷應力云圖.結果表明，各模型 σ^* 最大值分布在接觸表面后緣以及亞表面區域.然而，亞表面區域的裂紋萌生形核點并不明顯，這表明在微動疲勞晶體塑性有限元模型中，僅僅使用等效多軸損傷應力來識別裂紋萌生位置是不夠全面的.

為了準確識別微動疲勞接觸表面和亞表面的裂紋萌生位置，本文采用累積塑性滑移 p 作為裂紋萌生位置的判定標準.累積塑性滑移 p 基于微觀結構物理框架，能夠捕捉晶體位錯滑移，并適用于評估疲勞損傷和壽命預測[35]，表示如式（13）：

圖9是三組子模型累積塑性滑移 p 分布圖.對于接觸表面的裂紋萌生位置而言，各組子模型的 p 最大值均出現在接觸后緣，與前述結論一致，說明晶粒梯度對微動疲勞接觸表面裂紋萌生位置影響并不顯著.與圖8亞表面區域裂紋萌生形核點不明顯相比，圖9累積塑性滑移清晰地顯示了裂紋萌生位置.在不同晶粒梯度子模型中，對于微動疲勞亞表面裂紋萌生位置的識別如下：隨著晶粒梯度的增大，亞表面裂紋萌生形核點的水平位置從子模型中心軸左側向右側移動，豎向位置都在接觸面以下約 0.3mm 處，由此可見，晶粒梯度越大，亞表面裂紋萌生形核點的水平位置距離中心軸線右側最遠.就累積塑性滑移云圖形態而言，三組晶粒梯度均表現出不同程度的各向異性，這將影響后期裂紋擴展走向.子模型No.1的累積塑性滑移云圖呈團狀分布，子模型No.3呈傾斜條帶狀分布，而子模型No.2介于兩者之間.在不同晶粒梯度的子模型中，裂紋主要在材料的晶界處萌生，僅少量裂紋出現在晶粒內部，這表明晶界易受損傷，促使裂紋的萌生和擴展，這一結果與文獻[36]中的結論一致.因此，在金屬的設計和成形過程中，應重視提升晶界的強度和抗斷裂性能.

圖9接觸區域累積塑性滑移p分布 Fig.9Distribution of accumulated plastic slip p （2 for the contact region

3.3晶粒梯度對微動疲勞裂紋萌生壽命的影響

在本研究中，使用累積塑性滑移 p 預測微動疲勞裂紋萌生壽命.單次循環累積塑性滑移 p_cyc 在疲勞循環加載穩定階段應保持恒定，累積塑性滑移 p 與裂紋萌生壽命 N_i 成反比，裂紋萌生壽命由式（14）表示：

式中： p_crit 是疲勞過程中裂紋萌生累積塑性滑移臨 界值.

圖10是三組子模型裂紋萌生形核點的累積塑性滑移隨微動疲勞循環周期的演化過程.結果表明，在循環初期，子模型No.1亞表面裂紋萌生的累積塑性滑移小于接觸面，但隨著循環次數增加，亞表面累積塑性滑移逐漸超過接觸面，最終導致裂紋在亞表面處萌生，屬于文獻[34]中描述的微動疲勞裂紋萌生位置之一.相比之下，子模型No.2接觸面和亞表面裂紋萌生的累積塑性滑移趨勢與子模型No.1相反，即由亞表面累積塑性滑移大于接觸面逐漸過渡到小于接觸面.而子模型No.3接觸面裂紋萌生的累積塑性滑移始終大于亞表面.隨著晶粒梯度的增大，這種趨勢變化更為顯著.這是因為在同一子模型中，大尺寸晶粒內部取向更加一致，使得位錯難以穿越晶界，從而減少了塑性滑移的發生，導致材料亞表面的累積塑性滑移減少，進而影響了裂紋萌生壽命，使亞表面的裂紋萌生壽命明顯大于接觸表面.由此可見，晶粒梯度對微動疲勞兩種類型的裂紋萌生累積塑性滑移具有顯著影響.通過對比所有裂紋萌生點，發現子模型No.1累積塑性滑移最小，裂紋萌生壽命最大.這是因為當晶粒梯度較小時，晶粒間應力分布更均勻，能夠減少塑性滑移的集中程度，從而降低裂紋萌生傾向，延長裂紋萌生壽命，這符合子模型低應力梯度在疲勞過程中變化不明顯、不易發生裂紋萌生和擴展的特征.相反，當晶粒梯度較大時，子模型內的應力梯度也會增加，這種應力梯度的增加會促進位錯滑移，導致累積塑性滑移增大，從而增加裂紋的萌生和擴展風險.

圖10子模型裂紋萌生點累積塑性滑移演化過程 Fig.10Evolutionsofaccumulated plasticslip p forthecrack initiationofthesubmodel

4結論

1）晶粒梯度對微動疲勞接觸區域應力場分布產生顯著影響.在微動疲勞中，變形主要發生在大尺寸晶粒中，晶粒梯度越大，意味著晶粒間尺寸變化越大，從而導致應力空間梯度增大.

2）晶粒梯度對微動疲勞接觸面上裂紋萌生位置影響較小，但對亞表面裂紋萌生影響很顯著.當晶粒梯度d/d值為1時，亞表面裂紋萌生位置位于子模型中心軸左側較遠處；隨著晶粒梯度增大，即dd值為1.75時，亞表面裂紋萌生位置逐漸向子模型中心軸靠近，但仍然在左側；當晶粒梯度d/d值進一步增加至2.5時，亞表面裂紋萌生位置則移至子模型中心軸右側.

3）晶粒梯度對微動疲勞裂紋萌生壽命影響顯著.晶粒梯度 d_c/d 值為1時，微動疲勞裂紋萌生壽命最長.隨著晶粒梯度增大，微動疲勞亞表面的累積塑性滑移小于接觸表面，即亞表面的裂紋萌生壽命大于接觸表面.

以上結論表明，晶粒梯度影響晶體材料應力場分布，并對金屬構件的微動疲勞裂紋萌生位置和壽命產生影響.通過金屬材料加工方法改變晶粒梯度，將有效改善微動疲勞裂紋萌生位置和裂紋擴展方向.

參考文獻

[1] SZOLWINSKI M P，FARRIS T N.Observation，analysis and

prediction of frettingfatigue in2024-T351 aluminumalloy[J]. Wear，1998，221（1）：24-36.

[2]SU Y，HAN Q N，ZHANG C C，et al．Effects of secondary orientation and temperature on the fretting fatigue behaviors of Nibased single crystal superallys ［J].Tribology International， 2019，130：9-18.

[3]GURER G，GUR C H.Failure analysis of fretting fatigue initiation and growth on railway axle press-fits[J]．Engineering Failure Analysis，2018，84：151-166.

[4]POURHEIDAR A， REGAZZI D，CERVELLO S，et al. Fretting fatigue analysis offull-scale railway axles in presence of artificial micro-notches[J].Tribology International，2020，150：106383.

[5]HELMI A M.Freting fatigue and wear damage of structural components in nuclear power stations：Fitness for service and life management perspective［J]. Tribology International，2006， 39（10）：1294-1304.

[6］趙華，周仲榮．數值方法在微動疲勞研究中的應用進展[J]. 摩擦學學報，2000，20（4）：317-320. ZHAO H，ZHOU Z R.Application and progress of numerical methodsinanalysisof fretting fatigue[J].Tribology，2000， 20（4）：317-320.（in Chinese）

[7]LI J，ZHOU JZ，FENG A X，et al. Analysis of microstructure and tensile properties produced by cryogenic laser peening on 2024- T351 aluminum alloy[J].Vacuum，2018，158：141-145.

[8]ZHOU J Z，LI J， XUSQ，et al. Improvement in fatigue properties of 2024-T351 aluminum ally subjected to cryogenic treatment andlaser peening[J]. Surface and Coatings Technology，2018， 345：31-39.

[9]WANG JT，ZHANG Y K，CHEN JF，et al.Effect of laser shock peening on the high-temperature fatigue performance of 7075 aluminum ally[J]. Materials Science and Engineering：A，2017， 704：459-468.

[10］劉道新，何家文．微動疲勞影響因素及鈦合金微動疲勞行為 [J].航空學報，2001，22（5）：454-457. LIU DX，HEJW.Review of factors that influence fretting fatigue（ff）and investigation on ffbehavior of Ti-alloy[J]．Acta Aeronautica et Astronautica Sinica，2001，22（5）： 454-457.（in Chinese）

[11］張天楊，陳思源，TOURSANGSARAKIM，等．2024-T351鋁合 金激光沖擊微觀組織強化機制分析[J]．航空制造技術，2023， 66（16） ：111-116. ZHANG T Y，CHEN S Y，TOURSANGSARAKI M，et al. Analysis of microstructure strengthening mechanisms of laser peening on2024-T351 aluminumalloy[J］．Aeronautical Manufacturing Technology，2023，66（16）：111-116.（inChinese）

[12］秦芳誠，李永堂，巨麗，等．環件輾擴成形過程微觀組織及性 能控制研究進展[J]．機械工程學報，2016，52（16）：42-56. QIN FC，LIY T，JUL，et al. Advances in control of microstructure and properties during rolling forming of ring parts[J]. Journal of MechanicalEngineering，2016，52（16）：42-56.（inChinese）

[13]MCCARTHY OJ，MCGARRY JP，LEEN SB. The effect of grain orientation on fretting fatigue plasticity and lifeprediction[J]. Tribology International，2014，76：100-115.

[14]ZHANG M，MCDOWELL D L，NEU R W. Microstructure sensitivity of freting fatigue based on computational crystal plasticity[J].TribologyInternational，2009，42（9）：1286-1296.

[15］林恩強，郭然，羅吉祥，等.螺栓緊固鋁板疲勞裂紋萌生的有 限元參數分析[J].工程力學，2010，27（6）：245-251. LIN EQ，GUO R，LUO JX，et al. Analysis of the parameter effectonfatiguecrack inriveted componentswith finite element method[J]. Engineering Mechanics，2010，27（6）：245-251.（in Chinese）

[16]丁俊，康國政，彭金方，等．調質42CrMo鋼的彎曲微動疲勞實 驗研究[J]．工程力學，2014，31（8）：230-234. DING J，KANG G Z，PENG JF，et al．Experimental study on bending fretting fatigue of tempered42CrMo steel [J]. Engineering Mechanics，2014，31（8）：230-234.（in Chinese）

[17］趙杰江，尹東，唐力晨，等.微動接觸狀態對鋯合金磨損速率 影響的數值計算與分析[J].工程力學，2016，33（增刊1）： 275-282. ZHAOJJ，YIN D，TANG L C，et al. Numerical calculation and analysis of the effect offretting contact conditions on the wear rate of zircaloy[J]. Engineering Mechanics，2016，33（Sup.1） ：275- 282.（in Chinese）

[18]LIULY，YANGQS，LIU X，et al.Crystal cracking of graingradient aluminum by a combined CPFEM-CZM method [J]. EngineeringFracture Mechanics，2021，242：107507.

[19]HOJJATI-TALEMIR，ABDEL WAHABM，DE PAUWJ，et al. Prediction of fretting fatigue crack initiation and propagation lifetimefor cylindricalcontact configuration[J].Tribology International，2014，76：73-91.

[20］劉道新，劉軍，劉元鏞．微動疲勞裂紋萌生位置及形成方式研 究[J]．工程力學，2007，24（3）：42-47. LIUDX，LIUJ，LIUYY.Studyon nucleatinglocationand formationof frettingfatiguecracks[J].EngineeringMechanics， 2007，24（3）：42-47.（in Chinese）

[21]LAI Y S，WANG TH.Verification of submodeling technique in thermomechanical reliability assessment of flip-chip package assembly[J]．MicroelectronicsReliability，2005，45（3/4）： 575-582.

[22]HAREWOODFJ，MCHUGH PE.Investigationof fnite element meshindependencein ratedependentmaterials[J] Computational Materials Science，2006，37（4）：442-453.

[23]WANG J，ZHOU C Z.A meso-scale fretting fatigue simulation method based on submodelling technique[C]//Proceedingsof the ASME 2021 International Mechanical Engineering Congress and Exposition. Virtual Online：American Society of Mechanical Engineers，2021.

[24]CORMIER N G，SMALLWOOD B S，SINCLAIRG B，et al. Aggressivesubmodellingofstressconcentrations[J] International Journal forNumericalMethodsinEngineering，1999， 46（6）：889-909.

[25]ZHANGKS，SHIYK，JUJW.Grain-level statistical plasticity analysison straincycle fatigue ofaFCC metal[J].Mechanics of Materials，2013，64：76-90.

[26]HUANG Y G.A user-material subroutine incroporating single crystal plasticity inthe ABAQUS finite element program[M]. Cambridge，UK：HarvardUniversity，1991.

[27]PEIRCE D，ASARO R J，NEEDLEMAN A.An analysis of nonuniformand localized deformation in ductile single crystals[J]. ActaMetallurgica，1982，，30（6）：1087-1119.

[28]HORSTEMEYER MF，MCDOWELL DL. Modeling effectsof dislocation substructure in polycrystal elastoviscoplasticity [J]. MechanicsofMaterials，1998，27（3）：145-163.

[29］BENNETT VP，MCDOWELLD L.Polycrystal orientation distribution effects on microslip inhigh cycle fatigue[J]. InternationalJournal ofFatigue，2003，25（1）：27-39.

[30]DIARDO，LECLERCQ S，ROUSSELIERG，et al. Evaluation of finiteelement based analysis of 3D multicrystalline aggregates plasticityapplicationtocrystalplasticitymodel identificationand thestudy of stress and strain fields near grain boundaries[J]. InternationalJournalofPlasticity，2005，21（4）：691-722.

[31]CHENQQ，XUYZ，LINXH，etal．Acrystal plasticitybased methodology forcrack initiation lifetimeprediction of fretting fatiguein AA2024-T351 alloy[J].Construction and Building Materials，2023，408：133578.

[32]LI J，ZHOU J Z，SUN Y H，et al． Study on mechanical properties and microstructure of 2024-T351 aluminum alloy treated by cryogenic laserpeening[J].Opticsamp;Laser Technology，2019， 120：105670.

[33]ZHOU|H，WEI PT，LIU HJ，et al.Crystal elasticity analysis of contact fatigue behavior ofa wind turbine gear[J].Journal of Mechanical Science and Technology，2019，33（10）：4791- 4802.

[34]WANGJ，ZHOU CZ. Analysis of crack initiation location and its influencing factors offreting fatigue inaluminum alloy components[J].Chinese Journal of Aeronautics，2022，35（6）： 420-436.

[35]MANONUKUL A，DUNNEFPE.High-and low-cycle fatigue crackinitiation usingpolycrystal plasticity[J]. Proceedings of the RoyalSocietyofLondon SeriesA：Mathematical，Physical and EngineeringSciences，2004，460（2047）：1881-1903.

[36]YALCINKAYA TiO，FIRATA O. Inter-granular cracking throughstrain gradient crystal plasticity and cohesive zone modeling approaches[J].Theoretical and Applied Fracture Mechanics，2019，103：102306.