航空榫連結構接觸有限元分析及影響因素探究

任浩杰，何法江，陳志雄

（上海工程技術大學 航空運輸學院，上海 201620）

0 引 言

接觸問題是指兩物體之間由較小的接觸面傳遞作用力，普遍存在于工程應用中，如：過盈裝配、傳動齒輪和連桿齒形配合等［1］。航空發動機榫連接結構就是典型的接觸工作條件下的結構，被廣泛應用在發動機風扇、壓氣機等部件的葉盤之間。研究表明，榫頭和榫槽在接觸面上的相互作用復雜，局部存在較大的接觸應力，容易導致損傷發生破壞，需要計算接觸應力，確定其危險區域［2－3］。但是接觸問題是一種復雜的邊界條件非線性問題，存在各種因素影響，工程應用中很難通過解析獲得接觸應力。目前主要采用有限元法求解復雜接觸問題，有限元法將一個復雜的系統分解為大量的小單元，通過對這些相互作用的單元進行求解，去逼近無限未知量的真實系統。當前國內對榫連接結構的數值計算集中在沿接觸面的應力平均值，對高應力梯度區域分析較少且精度偏低，不能精準確定危險失效區域［4］。然而，直接提高精度會導致有限元計算消耗大量的時間與資源，因此需要確定合適建模方式，使之能夠高效、高質量求解，并且要嚴格保證計算收斂，才能得到準確的接觸應力分布，進而確定危險失效區域。

本文首先從典型的赫茲接觸問題著手，對比有限元解和解析解，以此確定接觸問題的有限元法分析流程；其次，對壓氣機燕尾榫連接結構的接觸問題進行分析，探究網格密度對接觸應力求解精度的影響，確定危險失效區域，及進一步分析接觸應力的影響因素。

1 赫茲接觸問題

由于接觸問題的非線性，通常采用有限元法進行求解，但在處理和計算中存在計算量過大、收斂困難和精度控制難等問題［5］。赫茲接觸理論是接觸應力分析的基礎，可以通過有限元法對典型的赫茲接觸問題進行求解，對比理論解析解，判斷有限元分析結果是否合理，并以此確定采用有限元法求解接觸問題的流程［6］。以兩半徑為50 mm、寬100 mm、泊松比0.3、彈性模量200 GPa 的平行圓柱體的接觸為例，假設載荷20 KN，分別基于赫茲公式和ANSYS軟件計算接觸半寬及最大接觸應力。算例的接觸形式及應力分布，如圖1 所示。根據赫茲接觸中兩個圓柱體的接觸應力計算公式可以推導接觸半寬a公式（1）和最大接觸應力σHmax公式（2）。

圖1 接觸形式及應力分布Fig.1 Contact pattern and stress distribution

由式（1）和式（2）計算得到a＝0.240 7 mm，σHmax＝528.992 3 MPa。考慮到模型的對稱性，對模型進行簡化，取模型的對稱部分進行有限元仿真分析如圖2 所示，得到有限元解a＝0.247 3 mm，σHmax＝528.37 MPa，與根據赫茲接觸理論所得的解析解吻合。

圖2 兩圓柱體接觸的有限元模型Fig.2 Finite element model of contact between two cylinders

算例中接觸區域網格尺寸設置為0.01 mm。然而，在后續的零部件有限元分析中過小的網格尺寸會導致計算量巨大，因此需要確定合適的網格尺寸范圍。本文以接觸半寬為參照，接觸部分網格尺寸／接觸半寬（d／a）設置為0.02、0.05、0.1、0.2、0.5、1、2、3、4、5，得到10 種密度的網格模型，如圖3 所示，求解得到最大接觸應力隨網格密度的變化，如圖4所示。結果表明要得到可靠的最大接觸應力，d／a應在0.5 以下。

圖3 不同網格密度的有限元模型Fig.3 Finite element models with different mesh densities

圖4 最大接觸應力隨網格密度變化Fig.4 Maximum contact stress varies with mesh density

2 燕尾榫連接結構的有限元分析

在明確了有限元法求解接觸問題的流程和難點后，對某型號航空發動機燕尾榫連接結構進行了有限元分析。

2.1 結構模型及材料

燕尾榫連接結構模型如圖5 所示，榫槽和榫頭材料均為TC11，化學成分為Ti 基，常溫下機械性能為：σ－1＝540 MPa，σb＝1 130 MPa，σs＝866 MPa，E＝120 GPa，μ＝0.33；TC11 鈦合金材料拉伸的應力—應變曲線，如圖6 所示［7］。

圖5 燕尾榫連接結構模型Fig.5 Dovetail joint structure model

圖6 TC11 拉伸的應力—應變曲線Fig.6 Stress－strain curve of TC11 tension

2.2 有限元計算

接觸區域的應力應變非常復雜，為了更準確分析接觸應力，本文在有限元模型的接觸區域細化網格，遠離接觸區域網格設置較粗糙，以減少計算量如圖7 所示。根據燕尾榫連接結構的實際工況設置載荷和約束，輪盤上施加15 KN 均布拉力模擬周向力，在榫頭上端面施加10 KN 均布拉力模擬葉片產生的離心力，輪盤施加位移約束，摩擦系數0.5［8］。初始狀態下接觸面間距為0，接觸力為0。計算結果表明，接觸區域存在高應力梯度，該工況下接觸應力最大值位于接觸區域下邊緣點，如圖8、圖9 所示。

圖7 燕尾榫連接結構有限元模型Fig.7 Finite element model of dovetail joint structure

圖8 應力分布云圖Fig.8 Stress distribution nephogram

圖9 接觸區域應力分布Fig.9 Stress distribution in the contact area

2.3 計算結果的收斂性分析

根據計算結果已初步確定接觸區域下邊緣存在高應力梯度，且網格密度對有限元解產生影響。為保證所得到的有限元解有效，對接觸區域下邊緣依次分別細化為0.048 mm、0.024 mm、0.012 mm、0.008 mm、0.004 mm。求解得到不同網格密度下應力最大值，并計算各種應力最大值之間的誤差，并根據公式（3）判斷結果得收斂性［9］。即相鄰兩次模型的解誤差在5%內即可認為收斂，不必再繼續細化網格。各種應力最大值隨網格密度的變化趨勢如圖10 所示，可以看出最大等效Mises 應力隨網格加密而明顯增大，但到一定尺寸后變化趨勢減小；切向應力σxx、法向應力σyy以及剪切應力τxy同等效Mises應力變化趨勢相同。不同模型計算誤差分析結果見表1，表明隨著網格密度的增大，最大應力值隨之增大，計算結果逐漸收斂。網格1 同網格2 以及網格2 同網格3 的計算結果相差較大，且都超過5%，表明該網格密度下的計算結果沒有收斂。網格3 和網格4 的計算結果相差4.3%，表明網格3 的計算結果已經達到網格收斂性要求。力分即若要得到準確的接觸應布，求解模型的網格密度需大于或等于網格3。

圖10 應力最大值隨網格密度的變化趨勢Fig.10 Trend of maximum stress with grid density

表1 不同尺寸網格模型計算誤差分析Tab.1 Calculation error analysis of mesh models of different sizes

2.4 接觸區域網格的優化

通過上文的分析計算，已經確定在榫連接接觸區域下邊緣存在高應力梯度，計算結果的精度受網格密度影響大。高密度網格模型直接導致計算量大幅度增加，需要對接觸區域的網格進行優化。由于高應力梯度區域只存在于接觸區域下邊緣，因此只對接觸區域下邊緣進行網格加密。改進后結果相差不大，但是計算時間明顯減少，對于后續處理更復雜的接觸結構具有很大幫助。

3 接觸應力影響因素的進一步分析

3.1 載荷對接觸應力的影響

發動機在實際工作運轉中存在多種工況，與之對應的榫連接結構也存在著不同的工況。不同載荷下的接觸應力也會產生變化，為了分析不同工況下接觸應力的變化，取10 KN、11 KN、13 KN、15 KN 4種工況載荷，接觸面初始間距為0，摩擦系數取0.5，建立相同的模型進行分析。4 種工況下，接觸應力都存在明顯的峰值，最大接觸應力值接近且都位于接觸下邊緣區，最大值所在區域相差不大。4 種工況下接觸應力的分布和變化趨勢一致，如圖11 所示，接觸應力隨載荷增大而增大，沿接觸面長度急劇增大達到峰值，隨后急劇減小最后趨近穩定，4 種載荷下危險失效區域一致。

圖11 不同工況下接觸應力分布圖Fig.11 Contact stress distribution under different working conditions

3.2 摩擦系數對接觸應力的影響

摩擦系數也是影響接觸應力的重要因素之一。為分析摩擦系數對接觸應力的影響，在接觸面設置不同摩擦系數的模型進行數值分析。鈦合金材料在常溫時摩擦系數一般在0.5 左右，因此摩擦系數取0.2、0.4、0.5、0.6、0.8，采用相同網格密度的模型，載荷為10 KN，初始狀態接觸間距為0。5 種不同的接觸面摩擦系數下，接觸應力的分布及變化趨勢一致，如圖12 所示，隨著摩擦系數增大，接觸面上接觸應力反而減小，同時越靠近接觸邊緣，摩擦系數的影響越明顯；隨著摩擦系數的改變，最大接觸應力出現區域變化不大，危險失效區域依舊在接觸面下邊緣。

圖12 不同摩擦系數下接觸應力分布圖Fig.12 Contact stress distribution under different friction coefficients

4 結束語

本文對赫茲接觸問題進行了有限元分析，明確了有限元分析接觸問題的流程和難點，結果表明網格密度影響數值求解精度。要得到可靠結果，網格尺寸和接觸半寬之比應在0.5 以下。對航空發動機燕尾榫連接結構接觸應力進行分析，著重探究了網格密度對計算結果收斂性的影響及載荷和摩擦系數對接觸應力的影響，得出結論：在燕尾榫連接結構接觸面存在較高的應力梯度；常規的網格密度在高應力梯度求解結果的精度不夠，需要優化高應力梯度位置網格；最大接觸應力點位于接觸區域下邊緣，接觸區域下邊緣為危險失效區域；隨著載荷和接觸面摩擦系數的變化，接觸應力分布趨勢不變，數值隨載荷增大而增大，隨摩擦系數增大而減小。上述的高應力梯度區域的網格優化、網格收斂性的分析，對有限元法接觸分析建模、平衡計算精度和計算量具有參考意義；榫槽在接觸區域下邊緣位置存在危險，對榫連接結構的設計和應用具有實際意義。