艾延廷，翟 學(xué)，王 志，喬永利

(沈陽航空航天大學(xué) 遼寧省數(shù)字化工藝仿真與試驗(yàn)技術(shù)重點(diǎn)實(shí)驗(yàn)室，沈陽 110136)

在工程領(lǐng)域，結(jié)構(gòu)分析的研究對(duì)象往往是由眾多零件、組件和部件經(jīng)過裝配而形成的機(jī)械系統(tǒng)。工程結(jié)構(gòu)裝配形式繁多，包括配合連接、螺栓連接、銷連接、焊接、膠結(jié)、嚙合連接等。在對(duì)裝配體進(jìn)行模態(tài)仿真和動(dòng)力學(xué)分析時(shí)，由于結(jié)合面處理的技術(shù)問題，不得不對(duì)分析模型進(jìn)行簡(jiǎn)化。如不考慮零件之間的連接配合關(guān)系，而是將裝配體中的各個(gè)零件直接合并為一個(gè)整體，認(rèn)為裝配件之間是剛性連接。由于有預(yù)載荷的機(jī)械裝配體(如螺栓法蘭連接)的接觸表面存在接觸應(yīng)力，不同的接觸應(yīng)力和接觸面粗糙度，使接觸面間形成分布不均勻的法向接觸剛度，會(huì)對(duì)整個(gè)裝配體的振動(dòng)模態(tài)產(chǎn)生影響。剛性連結(jié)假設(shè)不能正確反映裝配結(jié)構(gòu)間的連接剛度和阻尼，其動(dòng)力學(xué)特性計(jì)算結(jié)果往往與實(shí)際情況相差很大，無法實(shí)際應(yīng)用。本文的研究起源于對(duì)航空發(fā)動(dòng)機(jī)盤式拉桿轉(zhuǎn)子各輪盤之間以及盤軸間的連接剛度對(duì)發(fā)動(dòng)機(jī)轉(zhuǎn)子的軸向和橫向振動(dòng)模態(tài)影響的分析。該方面的研究主要方法有兩種，一是模態(tài)實(shí)驗(yàn)［1］，二是對(duì)部件做集中參數(shù)處理后通過鍵合圖法［2-3］或有限元法［4-5］進(jìn)行數(shù)值計(jì)算。隨著發(fā)動(dòng)機(jī)轉(zhuǎn)子以及整機(jī)結(jié)構(gòu)設(shè)計(jì)技術(shù)的發(fā)展，要求在設(shè)計(jì)初期即進(jìn)行轉(zhuǎn)子和整機(jī)的有限元建模，進(jìn)而通過仿真實(shí)現(xiàn)其動(dòng)力學(xué)優(yōu)化設(shè)計(jì)。文獻(xiàn)［1］提出的方法是對(duì)設(shè)計(jì)并加工出的實(shí)物模型進(jìn)行實(shí)驗(yàn)?zāi)B(tài)分析，因而無法在設(shè)計(jì)之初進(jìn)行模態(tài)計(jì)算與分析;文獻(xiàn)［2-3］使用的鍵合圖法與輪盤等結(jié)構(gòu)的有限元模型無法聯(lián)合使用。文獻(xiàn)［4］采用有限元法中的接觸單元進(jìn)行接觸剛度分析與模態(tài)計(jì)算。在一個(gè)裝配體中，常常有多處采用螺栓聯(lián)接結(jié)構(gòu)，如在某型航空發(fā)動(dòng)機(jī)轉(zhuǎn)子兩級(jí)盤-盤聯(lián)接中有50余螺栓，如果對(duì)每個(gè)螺栓都進(jìn)行有限元建模，則需要大量的計(jì)算機(jī)資源，造成大量的時(shí)間消耗，甚至計(jì)算機(jī)無法運(yùn)行，因此，該方法的實(shí)際應(yīng)用受到限制。文獻(xiàn)［5］通過對(duì)傳統(tǒng)的有限元方法進(jìn)行改進(jìn)，并編程計(jì)算，實(shí)現(xiàn)接觸面接觸剛度對(duì)轉(zhuǎn)子動(dòng)力特性影響的研究與分析。本文試圖對(duì)帶有預(yù)緊力的螺栓聯(lián)接結(jié)構(gòu)的接觸剛度進(jìn)行有限元等效計(jì)算，以便為工程中的螺栓聯(lián)接有限元建模與仿真問題提供一個(gè)有效方法。

本文在靜態(tài)接觸分析的基礎(chǔ)上，依據(jù)接觸力學(xué)相關(guān)理論，根據(jù)接觸面的粗糙度參數(shù)確定出接觸剛度［6-7］;然后在有限元模型中加入實(shí)體層單元，通過改變實(shí)體層單元的彈性模量來等效接觸剛度;最后再通過有限元方法對(duì)裝配體進(jìn)行模態(tài)分析。

1 接觸剛度理論

取表面的平均水平線為基準(zhǔn)線，稱剛性表面到基準(zhǔn)距離為間隔d，稱頂點(diǎn)到基準(zhǔn)的距離為頂點(diǎn)高度zs。zs的平均值為，其概率密度函數(shù)為φ(zs)。如果在名義表面積A0中有N個(gè)頂點(diǎn)，則在間隔d處接觸的頂點(diǎn)數(shù)為:

如果頂點(diǎn)高度超出了間隔，那么它將被壓縮δ=zs-d，并且將在半徑為a的小圓形區(qū)域中與平面構(gòu)成接觸。因此，第i個(gè)頂點(diǎn)具有的接觸區(qū)域?yàn)?

壓縮該頂點(diǎn)所需要的力可寫成:

式中f(δ)和g(δ)與接觸表面的材料性質(zhì)有關(guān)。如果變形完全是在彈性極限之內(nèi)，由Hertz方程有:

為了求出真實(shí)的總接觸面積A和總的名義應(yīng)力=P/A0，必須對(duì)高度超出間隔的所有凹凸起伏求和。于是有:

將式(4)代入式(5)，并進(jìn)行歸一化處理可得:

由式(6)可得接觸剛度kc:

2 模型建立

2.1 接觸剛度計(jì)算

許多實(shí)際的表面，特別是未處理過的基表面，其高度分布服從Gauss分布。對(duì)于標(biāo)準(zhǔn)偏差為σ的Gauss分布，頂點(diǎn)高度的分布非常接近Gauss分布，其標(biāo)準(zhǔn)偏差σs≈σ，頂點(diǎn)的高度位于表面的平均水平之上0.5σ～1.5σ之間;頂點(diǎn)曲率與該表面的均方根曲率量級(jí)相同，即=σk;如果波狀表面是隨機(jī)的各向同性表面，采用趨近于零的抽樣間隔，則單位面積的頂點(diǎn)數(shù)為ηs=1.209，其中ηp為跡線單位長(zhǎng)度上的峰點(diǎn)數(shù)，通過分析輪廓跡線中的峰點(diǎn)數(shù)確定。均值為zs，方差為σs的Gauss分布的概率密度函數(shù)為:

由式(6)和式(7)可得:

基于式(9)得到真實(shí)接觸面積和歸一化間隔隨無量綱載荷的變化關(guān)系如圖1所示。

圖1 真實(shí)接觸面積和歸一化間隔隨無量綱載荷的變化Fig.1 The variation of real contact area and normal step as a function of dimensionless load

根據(jù)式(7)和圖1，可以確定接觸剛度。

上面接觸剛度理論是基于假設(shè)每一個(gè)凹凸的變形與其相鄰的凹凸無關(guān)。當(dāng)實(shí)際接觸面積與名義接觸面積相比不再很小時(shí)，例如在歸一化間隔小于0.5的地方，誤差將加大。在另一端，如果歸一化間隔大約超過3.0時(shí)，則使用統(tǒng)計(jì)方法來研究接觸問題變得不精確。

2.2 有限元等效剛度建立

通用有限元計(jì)算軟件能夠?qū)佑|問題進(jìn)行靜力學(xué)分析，但是這些分析假定了接觸面連續(xù)，對(duì)于包含粗糙度的接觸面不能夠進(jìn)行分析。本文利用通用有限元軟件計(jì)算出接觸表面的接觸壓力;然后使用這些接觸壓力數(shù)據(jù)和接觸表面的粗糙度以及接觸部件的材料屬性，根據(jù)接觸力學(xué)的知識(shí)計(jì)算出接觸表面每個(gè)單元的接觸剛度;然后使用層單元，將接觸剛度折算為層單元的彈性模量來修改模型。

折算的具體過程如圖2所示。當(dāng)具有名義面積A0的接觸單元和下面的目標(biāo)單元接觸時(shí)。該法向接觸剛度可等效為層單元的剛度，如圖2所示。建模時(shí)比較方便的做法是將原來靠在一起的接觸單元和目標(biāo)單元移開一定距離h，可以取接觸時(shí)的最大間隔值。等效公式為:

通過調(diào)節(jié)單元的彈性模量E'就可以等效接觸剛度kc了。

圖2 使用層單元來模擬接觸剛度Fig.2 Contact stiffness simulated with solid elements

由于不考慮切向接觸剛度，材料的密度和泊松比都可以看做為0。接觸面用層單元模擬，接觸面的剛度用層單元彈性模量等效，層單元的厚度為h，層單元兩端與原來實(shí)體的連結(jié)采用Nastran軟件中的多點(diǎn)約束技術(shù)完成。

3 算例

3.1 接觸剛度計(jì)算

本文研究的物理模型長(zhǎng)l=10 mm，寬b=10 mm，截面面積A=100 mm2，高h(yuǎn)=50 mm。材料參數(shù)為:彈性模量E=2.1e11 Pa，泊松比 μ =0.3，密度 ρ=7 850 kg/m3。接觸面取在立方體的中間位置。立方體的下端面支撐在剛性水平面上，立方體的上端面施加壓強(qiáng)108 Pa。

中間接觸面的間隔即層單元的厚度h取1 mm。文獻(xiàn)［4，7］使用靜態(tài)實(shí)驗(yàn)和模態(tài)實(shí)驗(yàn)，并通過參數(shù)識(shí)別對(duì)類似結(jié)構(gòu)的法向和切向接觸剛度進(jìn)行了測(cè)量。

3.2 基于解析解的模型固有頻率計(jì)算

接觸面間的接觸剛度可以通過接觸面間的接觸參數(shù)的實(shí)驗(yàn)測(cè)量及其分布特征分析按上述公式計(jì)算。而接觸剛度一旦確定，即可按層單元法和多點(diǎn)約束技術(shù)進(jìn)行有限元分析與計(jì)算。計(jì)算的準(zhǔn)確性，一般采用實(shí)驗(yàn)的方法對(duì)最終的模態(tài)分析結(jié)果進(jìn)行驗(yàn)證，實(shí)際包含對(duì)接觸剛度計(jì)算的準(zhǔn)確性和有限元模型的準(zhǔn)確性兩方面的驗(yàn)證。本文假設(shè)接觸剛度計(jì)算是準(zhǔn)確的，只進(jìn)行有限元模型的正確性驗(yàn)證。研究?jī)蓚€(gè)長(zhǎng)立方體，其中間接觸，接觸剛度用彈性系數(shù)為kn的彈簧模擬，研究其縱向振動(dòng)問題，建立圖3所示力學(xué)模型。

圖3 計(jì)算的力學(xué)模型Fig.3 Mechanical model for calculation

根據(jù)有關(guān)知識(shí)，可以推導(dǎo)出:

或:

由式(13)，可以求出各階固有頻率的解析解。

3.3 基于有限元等效的振動(dòng)模態(tài)計(jì)算

本文只分析結(jié)構(gòu)在豎直方向的振動(dòng)，所有節(jié)點(diǎn)水平方向約束。圖4中的h取為1 mm。由式(10)計(jì)算出等效實(shí)體單元的彈性模量為E'=1.04×1010N/m2，這比上下立方體的材料的彈性模量要小，所以其固有頻率比由上下立方體所構(gòu)成的整體部件要低。表1列出了考慮法向接觸剛度的等效模型解、整體有限元模型解與固有頻率的解析解對(duì)比，可以看出，有限元等效模型解與解析解非常接近，而不考慮接觸問題的整體模型解較有限元模型解誤差大。

圖5給出了二者前3階振型的對(duì)比，不同的顏色代表不同的位移量，如在第一階振型圖中，顏色從藍(lán)到紅，對(duì)應(yīng)振動(dòng)位移從0到最大值1。由圖5可以看出，等效模型與整體模型二者在振型上也存在差異。

圖4 物理模型(左)、有限元接觸模型(右)Fig.4 Physical model(left)，element contact model(right)

表1 考慮法向接觸剛度的等效模型、整體有限元模型與固有頻率解析解對(duì)比(Hz)Tab.1 Contrast of the frequencies for the equivalent model considering contact stiffness，the solid model and theoretical solution(Hz)

圖5 等效模型(右)同整體模型(左)的前3階模態(tài)振型對(duì)比Fig.5 Comparisons of the first three vibration models for the equivalent model(right)and solid model(left)

4 結(jié)論

本文探討了在接觸面間添加一層單元，通過改變實(shí)體單元彈性模量來建立考慮法向接觸剛度對(duì)裝配體振動(dòng)模態(tài)影響的有限元模型。研究表明，不考慮法向接觸剛度的有限元整體模型解較解析解誤差很大，一般超過10%以上，且頻率越低差異越大;而考慮接觸剛度的有限元等效模型解與固有頻率解析解差異很小，即只要接觸剛度計(jì)算準(zhǔn)確，本文提出的方法就能得到理想的計(jì)算結(jié)果，且計(jì)算簡(jiǎn)單、實(shí)用。另外，考慮與不考慮接觸剛度，其各階模態(tài)振型也存在一定差異，這主要是接觸剛度單元的存在引起的。本文分析中假設(shè)了通過預(yù)載所形成的法向接觸剛度在振動(dòng)中一直保持不變，事實(shí)上，由于在振動(dòng)中接觸面間的接觸間隔一直在改變，所以接觸剛度也會(huì)跟著改變，這屬于非線性振動(dòng)模態(tài)范疇，下一步將要進(jìn)行深入研究。本文的研究，有助于分析像發(fā)動(dòng)機(jī)盤式拉桿轉(zhuǎn)子這種連接剛度對(duì)系統(tǒng)振動(dòng)模態(tài)影響較大的情況。

［1］施麗銘，張艷春.拉桿式模型轉(zhuǎn)子固有頻率的實(shí)驗(yàn)與計(jì)算研究［J］.振動(dòng)與沖擊，2008，27(S):47-49.

［2］王艾倫，駱 舟.拉桿轉(zhuǎn)子軸向振動(dòng)的動(dòng)力學(xué)模型［J］.中國(guó)機(jī)械工程，2009，20(13):1524-1527.

［3］章圣聰，王艾倫.盤式拉桿轉(zhuǎn)子的振動(dòng)特性研究［J］.振動(dòng)與沖擊，2009，28(4):117-120.

［4］ Zheng Y，Rong Y，Hou Z.A finite element analysis for stiffness of fixture units［J］.Journal of Manufacturing Science and Engineering，2005，127:429-432.

［5］高 銳，袁 奇，高 進(jìn).燃?xì)廨啓C(jī)拉桿轉(zhuǎn)子有限元模型研究及臨界轉(zhuǎn)速計(jì)算［J］.熱能動(dòng)力工程，2009，24(3):305-308.

［6］饒柱石，夏松波.粗糙平面接觸剛度的研究［J］.機(jī)械強(qiáng)度，1994，16(2):71-75.

［7］約翰遜K J.接觸力學(xué)［M］.北京:高等教育出版社，1992.

［8］饒柱石.欄桿組合式特種轉(zhuǎn)子動(dòng)力學(xué)特性及其接觸剛度的研究［D］.哈爾濱:哈爾濱工業(yè)大學(xué)，1992.

［9］Greenwood J A，Tripp J H.The contact of two nominally flat rough surface［J］.Pro.IMechE，1970，185:102-108.