基于ANSYS 隨機粗糙表面接觸分析

周朝朋，楊照

（200093 上海市 上海理工大學 機械工程學院）

0 引言

隨著國民經濟不斷發展，高新技術快速更新迭代，涌現了以光刻機、數控設備、機器人、航空航天設備等為代表的現代裝備，為了滿足人民的需求，設備也越加復雜化，所有的裝備都離不開接觸面，而接觸面的特性對于機械設備的正常運行有著極大的影響。1966 年，Williamson 和Greenwood[1]提出了著名的GW 接觸模型，模型假設粗糙表面的微凸體的分布符合高斯分布，并將結合部等效為光滑的剛性面和粗糙的表面；在國內，張學良[2]采用螺栓固定的懸臂梁模型，通過敲擊法獲得了模型的固有頻率，利用MATLAB 軟件識別隨著面壓載荷的變化法向接觸剛度的變化規律，經過計算獲得不同面壓載荷下的單位結合面上的法向接觸剛度的數值。但由于加工工藝的多樣性和復雜性，到目前為止還沒有對接觸面有明確的理論基礎，本文以MATLAB 生成模擬隨機粗糙表面關鍵點數據并導入仿真軟件中進行有限元分析，就接觸力、接觸面積隨著位移的關系進行研究。研究表明，接觸力、接觸面積隨位移呈三次函數關系，并對之后研究接觸面特性有積極的參考價值。

1 隨機粗糙表面生成

1.1 生成粗糙表面數據函數

為了模擬復雜的接觸表面，經研究得知，加工表面凸起峰值符合高斯分布，現階段生成表面數據方式有2 種：第1 種利用Rand（）函數，y=sqrt（0.09）×rand（1 000，1）即可生成1 000個點方差為0.09 的隨機數據；第2 種利用wgn（）函數，y=wgn（1000，1，10×log10（0.09））即生成1 000 個點方差為0.09 的符合高斯分布的隨機數據。本文采用第2 種方式。

1.2 生成隨機表面數據

使用MATLAB 軟件，寫M 文件并運行生成隨機表面數據[3]。圖1、圖2 所示分別為粗糙度為Ra=3.2、Ra=6.3 的輪廓數據圖，表1 為局部數據表。從圖中可以看出隨機數據的復雜性。為了研究生成數據的分布規律，繪制數據分布規律直方圖，如圖3、圖4 所示。由統計學來看，數據分布具有規律性，完全可以模擬符合高斯分布的復雜曲面，分析結果也具有代表性。

圖1 粗糙度為3.2 輪廓數據圖Fig.1 Contour data diagram with roughness of 3.2

圖2 粗糙度為6.3 輪廓數據圖Fig.2 Contour data with roughness of 6.3

表1 MATLAB 生成的粗糙表面數據點Tab.1 Rough surface data points generated by MATLAB

圖3 第1 組白噪聲分布直方圖Fig.3 Histogram of white noise distribution

圖4 第1 組白噪聲分布箱線圖Fig.4 Boxed diagram of white noise distribution

2 Ansys 接觸特性分析

2.1 模型的建立

將MATLAB 生成的為粗糙表面的關鍵數據點直接導入ANSYS 中，在軟件中利用Spline 函數生成光滑的曲線。圖5 所示為粗糙表面局部圖。光滑的曲線可以避免在之后的網格劃分中畸形網格對分析結果產生影響。在建好光滑曲線之后，通過K，L 命令建立模型其他的輪廓線，由于直線的建立較為簡單，此處不再贅述。如圖6 為建立好的關鍵線輪廓。在粗糙接觸面曲線和關鍵輪廓線建好之后，利用AL 命令，建立分析所需要的面，如圖7 所示。

圖5 ANSYS 生成的光滑曲線Fig.5 Smooth curve generated by ANSYS

圖6 關鍵線輪廓線圖Fig.6 Contour diagram of key lines

圖7 模型面圖Fig.7 Model plane diagram

2.2 定義模型屬性及網格劃分

本次分析采用45#，其彈性模量：E=207 GPa，屈服強度：σs=0.355 GPa，泊松比：υ=0.3，密度：ρ=7 850 kg/mm3。分析類型為平面應變，使用平面182 單元。由于此次分析為接觸分析，分析的關鍵在接觸面[4]，為了節約計算空間且對計算結果不產生影響，所以在Y 方向采用分層并網格逐步加密的方法。圖8 所示為生成的網格圖，網格數：313 446；網格節點：306 011。

圖8 網格圖Fig.8 Grid diagram

2.3 接觸與約束

分析的塑性選項采用雙線性等向強化（BISO），此屈服模式下的 Von Mises 屈服準則也是適用的。隨機粗糙表面與剛性平面的接觸采用拉格朗日法則建立接觸關系。在粗糙模型Y 方向定義位移約束，為了清晰地得到接觸力、接觸面積隨位移的變化趨勢，采用300 個等步長加載，并關閉自動判斷選項，增加步長的同時也使求解結果更易收斂，且結果采用Von Mises 判斷依據。

2.4 結果分析

以粗糙度為Ra=3.2 為例，圖9、圖10 為接觸分析結果位移云圖和局部接觸應力圖。

圖9 位移云圖Fig.9 Displacement cloud map

圖10 局部接觸應力圖Fig.10 Local contact stress diagram

圖11、圖12 為應力應變曲線。從曲線可以看出，當應力達到0.3 GPa 時，零件由彈性階段轉變為彈塑性階段，即隨機粗糙表面接觸前期處于純彈性階段，達到屈服應力時，此時零件既有彈性變形，也有塑性變形[5]。

圖11 應力應變曲線Fig.11 Stress-strain curves

圖12 應力應變曲線局部放大圖Fig.12 Local enlargement of stress-strain curves

以壓入深度為橫坐標，分別以接觸力、接觸面積為縱坐標，做壓入深度-接觸力關系圖和壓入深度-接觸面積關系圖，如圖13、圖14 所示。

圖13 壓入深度和接觸力曲線Fig.13 Pressing depth and contact force curve

圖14 壓入深度和接觸面積曲線Fig.14 Pressing depth and contact area curve

以式（1）為擬合函數擬合[6]曲線，得到壓入深度接觸力擬合參數，見表2。

表2 壓入深度接觸力擬合參數Tab.2 Fitting parameters of pressing depth and contact force

由表2 可知，擬合精度完全符合要求。

表2 中數據分別代入式（1），可得：

3 結論

通過ANSYS 仿真分析，得到粗糙度分別為Ra=3.2，Ra=6.3，Ra=12.5 的隨機粗糙表面下的接觸應力、接觸面積隨著位移的變化關系，并根據變化曲線以三次函數擬合得到擬合曲線、擬合方程。由擬合誤差來看，擬合符合要求，可為今后研究接觸特性提供數據支持和參考。