拉桿結構中彈簧剛度和有限元模型剛度融合修正方法研究

張子陽，謝壽生，錢征文，龍 門，彭靖波

（1.空軍工程大學 工程學院，西安 710038；2.空軍工程大學 電訊工程學院，西安 710077）

建立航空發動機高壓轉子的有限元模型對于分析結構動態特性、提高裝配水平具有十分重要的意義［1］。某型航空發動機高壓轉子采用的拉桿結構增強了結構的剛度，但拉桿轉子的接觸面會使轉子局部剛度降低，而高壓轉子盤與盤之間的聯接采用的是過盈聯接的方式，處于彈性或彈塑性接觸狀態。因此，復雜的接觸情況使得拉桿轉子的結構變化復雜，如何在有限元建模過程中考慮這些因素的影響是進一步提高模型精度需要解決的問題［2，3］。

基于梁模型的高壓轉子整體有限元建模方法，不僅沒有考慮接觸部位縱向纖維的不連續性，而且忽視了橫向擠壓與松弛的可能性，因而夸大了接觸處的剛度。文獻［4］提出了刻畫熱套連接處性狀的力學模型，開始對整體轉子進行修正。文獻［5］分析了拉桿轉子的受力情況，然后考慮接觸面接觸剛度對轉子動力特性的影響，將拉桿和接觸面等效為一個鉸鏈和一個抗彎彈簧，對傳統的有限元方法進行了改進。文獻［6］利用ANSYS分析非線性接觸的模塊，對拉桿與圓盤以及圓盤之間的接觸行為進行研究，在一定程度上反映了拉桿結構的接觸狀態。

本文對盤與盤、拉桿與盤之間的接觸力學模型進行研究，對復雜接觸部位的有限元模型進行重構，提出將彈簧剛度矩陣的修正和有限元模型剛度矩陣的修正融合的方法，構建了考慮復雜接觸的有限元剛度整體優化模型，并且運用拉直算子將優化模型中需要修正的參數分離出來，結合實驗模態分析結果進行了修正。最后對某型航空發動機高壓轉子拉桿結構運用本文提出的方法進行了算法驗證，并且基于MAC對算法的有效性進行了檢驗。

1 拉桿結構接觸狀態的力學模型

根據對拉桿結構的受力分析，該結構主要包含兩個影響有限元分析結果的接觸方式：盤與盤之間的接觸以及盤與拉桿之間的接觸。

由于采用熱套工藝，并且盤與盤之間是過盈聯接的方式，裝配后的盤面摩擦力能夠使相鄰盤面的橫向位移保持一致，但不能使兩個接觸面的轉角相等。這將造成相鄰兩個擠壓面上一部分松弛，而另一部分受到再擠壓作用。因此，盤與盤之間的接觸可以用一組沿擠壓面均布的抗壓彈簧和鉸接頭來表示，如圖1所示。

拉桿通過螺母和圓盤發生的擠壓作用與圓盤間的擠壓作用相類似，可用鉸接頭和抗彎彈簧來表征，因此該接觸單元沒有質量，只有等效抗彎剛度。

高壓轉子的軸和拉桿均采用有限元的梁模型表示。

圖1 盤與盤之間接觸的力學模型Fig.1 Mechanical modal of contact between wheel disks

圖2 拉桿與盤之間接觸的力學模型Fig.2 Mechanical modal of contact between rod and wheel disks

2 彈簧剛度和有限元剛度的融合修正

2.1 修正參數的選擇

根據受力分析可以確定，拉桿結構中的復雜接觸是對結構動力特性影響顯著而且建模誤差集中的部位。為提高有限元建模精度，對這些局部的建模進行重構。但是在ANSYS模態計算中僅使用面接觸單元的原始參數是不能對這種復雜接觸進行完全描述的，因此需要根據結構的實驗模態參數，結合靈敏度分析，對參數進行修正。

抗壓彈簧和抗彎彈簧的剛度值反映了接觸的動態特性，經過模型結構調整后，初始估計的彈簧剛度值需要根據試驗數據進一步修正，因此選為需要修正的參數。

模型剛度陣的修正屬于矩陣型修正方法，而彈簧剛度的修正屬于元素型的修正方法［7］。由于拉桿轉子的剛度變化復雜，這兩個剛度的修正結果必然相互影響，且修正后的剛度陣改變了原矩陣的帶狀和稀疏性，因此本文將兩種修正方法進行融合，建立統一的優化模型，并運用拉直算子將需要修正的參數分離出來單獨進行修正。

2.2 融合修正模型的建立

通過對拉桿結構的理論分析給出初始的有限元模型，其質量矩陣為Ma，剛度矩陣為Ka，并且Ma∈Rn×n與Ka∈Rn×n均為實對稱矩陣。

假設通過實驗模態分析得到前s階的精確模態參數準確地反映了結構的動態特性，記為：

一方面，運用 Smith［8，9］提出的元素型修正方法，建立特征頻率誤差和彈簧剛度的關系式，通過求解使特征頻率誤差最小的剛度修正因子建立彈簧剛度的迭代計算式。第i階特征頻率ωi的誤差Δωi對抗壓彈簧和抗彎彈簧剛度的修正量Δg的1階Taylor展式為：

寫出矩陣形式為：

定義δ=Aα－b為方程的誤差。

另一方面，運用實驗模態參數對模型剛度陣進行修正時，需要滿足的特征方程及正交性條件為：

當方程（5）有解時，滿足方程（5）的矩陣K有無窮多個，而修正的目標就是這無窮多個K中找出最符合實際結構的剛度Ka的矩陣，即矩陣應滿足：

因此，通過引入系數λ建立統一的優化模型為：

2.3 基于拉直算法的模型優化求解

運用拉直算法［10］對優化模型進行求解。該方法能夠將需要修正參數分離出來直接對其修正。在保證計算模型的帶狀稀疏特性基礎上，將實驗模態數據融于修正模型。

令剛度修正矩陣為ΔK，首先對特征向量φ進行歸一化，令ΔKφi=1。滿足式（5）的條件下，可得：

進行拉直運算可得：

對式（9）左端變量進行調整：

（1）把對角線上的變量按對角線方向依次放在列向量vec（ΔK）的最前端；

（2）把下三角陣中非對角線上的變量按對角線方向依次放在列向量vec（ΔK）的中間部分；

（3）把上三角陣上非對角線上的變量按對角線方向依次放在列向量vec（ΔK）的最后部分。

記調整后ΦT?I為E，而調整后的vec（ΔK）仍用原符號表示，式（9）表示為：

對vec（ΔK）進行半拉直變換：

令F=E×V則：

在svec（ΔK）中僅有前qK個變量為所求變量。因此，刪除多余變量，取svec（ΔK）的前qk列組成新矩陣Y，得到：

對彈簧剛度修正的拉直算法中，取α中的前qK個變量，滿足δ=0的最優解可以通過使下述函數最小化得到：

其中μ為比例因子。

拉直后的優化模型為：

建立彈簧剛度修正的迭代方程為：

運用梯度法對優化模型進行求解即可得到優化后的剛度矩陣和適宜的彈簧剛度。

2.4 基于MAC的模型修正檢驗

為檢驗修正后的有限元模型的準確性，采用模態確認標準 MAC（modal assurance criterion）［11］對修正后的結果進行檢驗。

因此，彈簧剛度和有限元剛度的融合修正過程如圖3所示。

圖3 剛度融合修正過程Fig.3 Fusion stiffness modification of spring and FE modal

3 算例

為驗證本文提出的有限元修正算法的有效性，在拉桿轉子力學模型的基礎上，對某型航空發動機高壓轉子的拉桿結構進行建模分析。

對盤與盤之間以及盤與拉桿之間的接觸進行了重構，有限元模型的接觸面采用面接觸單元（CONTA174），重構后的矩陣形式與原來的矩陣形式相比，在ANSYS建模過程中由于包含了接觸單元，矩陣的規模擴大了，而且和接觸單元相關的單元的參數也發生了變化。按裝配工藝進行正常預緊，接觸部位的有限元模型如圖3所示。

以本實驗室研制的高壓轉子裝配檢測儀為實驗平臺，獲取高壓轉子拉桿結構的實驗模態參數。其中，實驗模態分析中的高壓轉子拉桿結構參數與仿真模型的參數是一致的。

通過求解初始有限元模型，可以得到修正前的模態參數。與實驗模態分析的數據進行比較，前4階特征頻率如表1所示。

圖4 某型航空發動機拉桿結構中接觸部位的有限元模型Fig.4 Contact surface in FE modal of high pressure spool

表1 實驗模態分析結果與修正前頻率對比HzTab.1 Frequency comparison between EMA and modal before updated Hz

運用本文提出的優化剛度修正方法，結合實驗模態分析得到的模態參數，對初始有限元模型進行修正。圖5以1號螺栓位置的4個接觸單元的剛度值為代表，給出了修正前后的單元剛度的取值和修正過程中參數的變化過程。修正后的頻率與實驗模態分析結果進行對比，如表2所示。

圖5 接觸單元剛度值的變化過程Fig.5 Evolutions of stiffness in contact element

表2 實驗模態分析結果與修正后頻率對比HzTab.2 Frequency comparison between EMA and modal after updated Hz

表3 模型修正前后MAC的對比Tab.3 MAC comparisons before and after modal updating

計算修正前后的MAC，對修正后的有限元模型的準確性進行檢驗，模型修正前后MAC的對比如表3所示，圖6、7表示了修正前后的模型與實驗模態分析結果的模態相關性。

圖6 修正前模型和實驗模型的模態相關性Fig.6 MAC between modal before updated and EMA

圖7 修正后模型和實驗模型的模態相關性Fig.7 MAC between modal after updated and EMA

4 結論

通過與航空發動機高壓轉子拉桿結構的實驗結果進行對比，說明修正后的基于復雜接觸的有限元模型真實反映了拉桿聯接的圓盤之間的接觸狀態。

本文提出的方法從拉桿結構的力學模型出發，在明確的物理模型基礎上，結合實驗模態分析結果，運用拉直算子對挑選過的參數進行修正，保證了修正后的質量陣與剛度陣之間的關聯性，避免了虛元與負剛度值的出現。

［1］章圣聰，王艾倫.盤式拉桿轉子的振動特性研究［J］.振動與沖擊，2009，28（4）：117 －120.

［2］ Bin X，Nan C.Frictional contact analysis and experiment study of complex mechanism with pinned connections［J］.Finite Elem.Anal.Des，2010，47（2）：142 －150.

［3］Alkatan F，Stephan P，Daidie A.Equivalent axial stiffness of various components in bolted joints subjected to axial loading［J］.Finite Elements in Analysis and Design，2007，43（8）：589－598.

［4］夏松波，武新華，汪光明.熱套轉子力學模型研究［J］.航空學報，1987，8（9）：449 －455.

［5］汪光明，饒柱石，夏松波.拉桿轉子力學模型的研究［J］.航空學報，1993，14（8）：419 －423.

［6］高 銳，袁 奇，高 進.燃氣輪機拉桿轉子有限元模型研究及臨界轉速計算［J］.熱能動力工程，2009，，24（3）：305－308.

［7］袁永新，蔣家尚.一種利用模態測量數據修正剛度矩陣的新方法［J］.振動與沖擊，2009，28（11）：34－37.

［8］ Smith M J，Hutton S G.Frequency modification using Newton’s method and inverse iteration eigenvector updating［J］.AIAA Journal，1992，30（7）：1886 －1891.

［9］李 輝，丁 樺.結構動力模型修正方法研究進展［J］.力學進展，2005，35（2）：170 －180.

［10］周星德，林小國，徐艷紅.有限元模型修正的虛擬拉直算法［J］.振動工程學報，2008，21（2）：152－156.

［11］Ewins D J.Modal testing：theory and practice［M］.Research Studies Press，2000.