劉 濤，劉 娟

（1.呼和浩特鐵路局 包頭車(chē)輛段，內(nèi)蒙古 包頭 150200；2.中國(guó)石油大學(xué)（北京）機(jī)械與儲(chǔ)運(yùn)工程學(xué)院，北京 102249）

旋轉(zhuǎn)機(jī)械在我們的生活中起著非常重要的作用，伴隨著現(xiàn)代工業(yè)的發(fā)展，旋轉(zhuǎn)設(shè)備向著高精度，高速化和智能化的方向發(fā)展，從而導(dǎo)致轉(zhuǎn)子和定子間的間隙越來(lái)越小，結(jié)構(gòu)越來(lái)越緊密，轉(zhuǎn)子碰摩故障發(fā)生得越來(lái)越頻繁[1]。

轉(zhuǎn)子與定子碰摩是一類(lèi)非常復(fù)雜的動(dòng)力學(xué)問(wèn)題。采用剛?cè)狁詈系霓D(zhuǎn)子有限元模型，研究一個(gè)或少數(shù)幾個(gè)參數(shù)對(duì)碰摩響應(yīng)行為的影響，為常用的研究路徑[2]。如何建立一個(gè)精度較高的含有碰摩故障的轉(zhuǎn)子模型，對(duì)轉(zhuǎn)子碰摩故障的研究，是很重要的。以模態(tài)試驗(yàn)為依據(jù)，進(jìn)行模型修正，可以有效地解決模型精度較低的問(wèn)題。從1968年Fox和Kapoor首次提出特征向量靈敏度的計(jì)算方法開(kāi)始，基于靈敏度的模型修正方法得到廣泛應(yīng)用[3]。徐張明，高天明和Amar Bouazzouni，Mohammad等都對(duì)靈敏度模型修正方法進(jìn)行了研究，且取得了一定的成果[4–6]。現(xiàn)在應(yīng)用較為廣泛的方法為基于有限元模型的設(shè)計(jì)參數(shù)型修正方法，其優(yōu)點(diǎn)是物理意義清晰，工程實(shí)用性更強(qiáng)。其中，基于靈敏度分析的設(shè)計(jì)參數(shù)型模型修正方法應(yīng)用最廣泛[7]。隨著工程需求的日趨提高，人們對(duì)模型修正技術(shù)的計(jì)算精度和計(jì)算效率也有了更高的要求。

為了提高剛?cè)狁詈限D(zhuǎn)子模型柔性轉(zhuǎn)軸的精度，以及模型修正技術(shù)的計(jì)算精度與計(jì)算效率。本文基于模態(tài)試驗(yàn)數(shù)據(jù)，將二分法–PSO作為模型修正的優(yōu)化算法，實(shí)現(xiàn)有限元模型修正，從而提高了有限元模型修正的速度與精度，且避免了在SPO迭代尋優(yōu)過(guò)程中陷入局部最優(yōu)解的問(wèn)題?；谀Ｐ托拚碚?，利用轉(zhuǎn)子系統(tǒng)模態(tài)試驗(yàn)獲取的轉(zhuǎn)軸固有頻率，以及WorkBench建立的有限元模型，進(jìn)行模型修正。研究結(jié)果表明，修正后的有限元模型精度有了大幅度的提高，可以較為準(zhǔn)確地反應(yīng)試驗(yàn)?zāi)Ｐ偷慕Y(jié)構(gòu)特征。最后利用修正后的模型進(jìn)行轉(zhuǎn)子碰摩故障動(dòng)力學(xué)特性的研究，考慮轉(zhuǎn)速變化對(duì)碰摩故障的影響，通過(guò)分析結(jié)果可以知道修正后的轉(zhuǎn)子模型可以較為準(zhǔn)確地反映試驗(yàn)?zāi)Ｐ偷倪\(yùn)行結(jié)果。

1 模型修正方法理論

1.1 靈敏度模型修正

以靈敏度為基礎(chǔ)的迭代修正法與有限元模型的離散本質(zhì)一致，只需要修正系統(tǒng)中的某幾個(gè)單元或節(jié)點(diǎn)的參數(shù)，就可以使有限元模型迅速收斂于與試驗(yàn)數(shù)據(jù)良好相關(guān)的修正模型，因此備受青睞[8–10]。目標(biāo)函數(shù)是最優(yōu)的設(shè)計(jì)性能，它是有關(guān)設(shè)計(jì)變量的函數(shù)

其中p代表設(shè)計(jì)參數(shù)向量，{fe}，{fp(p)}分別代表結(jié)構(gòu)動(dòng)態(tài)特性的實(shí)驗(yàn)值和分析值；r(p)為殘差項(xiàng)；VLB、VUB分別代表結(jié)構(gòu)設(shè)計(jì)參數(shù)變化的上下限；Wf代表結(jié)構(gòu)各個(gè)特征量之間的加權(quán)矩陣。

靈敏度的計(jì)算公式為

則修正公式可以寫(xiě)成

其中S?為靈敏度矩陣S的廣義逆，因此只要得到系統(tǒng)的靈敏度矩陣，確定了特征量的變化量就可以通過(guò)優(yōu)化迭代法確定參數(shù)的變化量。

1.2 粒子群優(yōu)化算法

PSO是一種隨機(jī)搜索全局優(yōu)化新方法，相較于其他算法，具有實(shí)現(xiàn)容易、精度高、收斂快等特點(diǎn)[11]。

對(duì)于有N個(gè)粒子群Xi（xi1，xi2，…，xid）的速度，位置更新公式為

其中Vi(vi1，vi2，…，vid)為始速度，Pbesti為粒子X(jué)i找到的最優(yōu)點(diǎn)，Gbest為所有粒子找到的最優(yōu)點(diǎn)，C1，C2為加速度因子，r1，r2是均值，用以增加粒子的多樣性，ω為慣性權(quán)。

1.3 二分法–PSO

PSO算法在迭代過(guò)程中隨機(jī)性較強(qiáng)，要達(dá)到高精度需要較長(zhǎng)時(shí)間。為此將二分法引進(jìn)迭代中，由盲目搜索變成計(jì)劃搜索。

因此，按照?qǐng)D1中的流程，完成對(duì)轉(zhuǎn)子系統(tǒng)的模型修正。

圖1 基于二分法–SPO的模型修正總體思路

具體算法：

(1)初始化X，V，Gbest，Pbest；

(2)更新X，V；

(3)給定的搜索區(qū)間：

計(jì)算函數(shù)Fj=f(X1min，X2min，…，Xjmax，…Xnmax)，

其中j=1，2…n；選取最小值Fi=min(F1，F(xiàn)2，…Fn)，即Fi=f(X1min，X2min，…，Ximax)；

(4)按式(5)或者式(6)更新參數(shù)X；

(5)判斷：設(shè)定誤差e，若Fi＜e，進(jìn)行步驟6，若若Fi＞e，返回步驟3；

(6)更新Pbest，Gbest；

(7) 判 斷 如 果 iter＞iterMax，或 者 Gbest＞GbestMax，則退出計(jì)算，否則返回步驟2。

2 分析模態(tài)

2.1 試驗(yàn)?zāi)B(tài)分析

為了獲得試驗(yàn)數(shù)據(jù)，通過(guò)對(duì)自由-自由狀態(tài)下的轉(zhuǎn)子結(jié)構(gòu)進(jìn)行激振，測(cè)量激振力和響應(yīng)數(shù)據(jù)，進(jìn)行數(shù)據(jù)采集和信號(hào)分析，從而獲得模態(tài)數(shù)據(jù)[12]。試驗(yàn)時(shí)用充氣的輪胎墊住轉(zhuǎn)子結(jié)構(gòu)以模擬自由-自由狀態(tài)。如圖3所示。

圖2 模態(tài)試驗(yàn)設(shè)計(jì)圖

圖3 轉(zhuǎn)子有限元模型

選擇轉(zhuǎn)軸的8個(gè)節(jié)點(diǎn)作為測(cè)點(diǎn)，采用單點(diǎn)激勵(lì)輸入、多點(diǎn)輸出方式，對(duì)轉(zhuǎn)子系統(tǒng)的每一個(gè)測(cè)點(diǎn)施加一定的激勵(lì)，用加速度傳感器采集激勵(lì)力信號(hào)，采集節(jié)點(diǎn)上的振動(dòng)響應(yīng)信號(hào)。其中設(shè)備包括：轉(zhuǎn)子振動(dòng)實(shí)驗(yàn)臺(tái)，傳感器，力錘，等。根據(jù)采集到的激勵(lì)信號(hào)和響應(yīng)信號(hào)，識(shí)別出轉(zhuǎn)軸前4階固有頻率，如表1所示。

表1 轉(zhuǎn)軸前4階固有頻率

2.2 計(jì)算模態(tài)分析

在SolidWorks軟件中建立轉(zhuǎn)子系統(tǒng)的三維實(shí)體模型，模型構(gòu)建過(guò)程中忽略該結(jié)構(gòu)件中某些局部結(jié)構(gòu)，如圓角、倒角等細(xì)節(jié)以減少有限元計(jì)算時(shí)間，提高有限單元的協(xié)調(diào)性。將適當(dāng)簡(jiǎn)化后的三維模型導(dǎo)入ANSYS Workbench進(jìn)行有限元分析。由于模型中的軸承為滑動(dòng)軸承。因此，通過(guò)彈性阻尼連接方式的Bushing連接進(jìn)行模擬，由于阻尼較小忽略不計(jì)，而軸承剛度通過(guò)K值來(lái)計(jì)算。

邊界條件為：在轉(zhuǎn)軸的兩端以及中間部分添加圓柱約束，同時(shí)，在轉(zhuǎn)軸中安裝轉(zhuǎn)子的部分添加應(yīng)力。網(wǎng)格劃分：采用Workbench自動(dòng)網(wǎng)格劃分技術(shù)，智能地控制網(wǎng)格大小，使得計(jì)算精度不會(huì)太低，又可以保證計(jì)算速度不會(huì)太低。最后進(jìn)行受力分析，建立有限元模型，同時(shí)進(jìn)行模態(tài)分析計(jì)算，得到相應(yīng)的模態(tài)參數(shù)。轉(zhuǎn)子系統(tǒng)有限元模型如圖3所示，圖4為轉(zhuǎn)軸前4階模態(tài)的振型圖。

3 模型修正結(jié)果

參數(shù)的選擇與試驗(yàn)中靈敏度有很大的關(guān)系。首先考慮可能影響有限元模型的幾個(gè)特征參量：彈性模量E，密度M，質(zhì)量系數(shù)M、C，轉(zhuǎn)軸的截面半徑R，轉(zhuǎn)軸的長(zhǎng)度L，轉(zhuǎn)軸的支承剛度K。而各個(gè)參數(shù)對(duì)模型特性的影響程度不同，鑒于參數(shù)靈敏度是選擇優(yōu)化參數(shù)的一個(gè)重要指標(biāo)，表示各參數(shù)對(duì)模型特性的影響程度。因此，可以通過(guò)靈敏度的分析確定需要優(yōu)化的參量。

圖4 轉(zhuǎn)軸前4階模態(tài)特性

如圖5所示，各參數(shù)在各階次模態(tài)下的靈敏度是不同的。為了提高迭代優(yōu)化的效率，防止在迭代過(guò)程中發(fā)生病理現(xiàn)象，在參數(shù)選擇中，需選擇靈敏度較高的參量。質(zhì)量系數(shù)MC與長(zhǎng)度L的靈敏度相對(duì)于其他參量較小可以忽略不計(jì)，轉(zhuǎn)軸的截面半徑R，彈性模量E，軸承剛度K靈敏度相對(duì)較大，可以選擇為待修正的設(shè)計(jì)參量。同時(shí)，密度M的靈敏度雖不是很大，但對(duì)計(jì)算結(jié)果影響較大，因此也作為優(yōu)化參量。

圖5 優(yōu)化參數(shù)及其靈敏度分析

在模型修正的過(guò)程中，各階次的固有頻率與試驗(yàn)結(jié)果間的誤差是不斷變化的。如圖6所示。

這是分別采用PSO算法與二分法–PSO算法對(duì)目標(biāo)函數(shù)進(jìn)行優(yōu)化迭代的過(guò)程。PSO算法與二分法–PSO算法都是伴隨著迭代次數(shù)的增加，目標(biāo)函數(shù)值不斷變小，直到最后趨于穩(wěn)定。然而PSO算法對(duì)目標(biāo)函數(shù)的優(yōu)化迭代在20次以上誤差逐漸穩(wěn)定且趨于0，而二分法–PSO方法對(duì)目標(biāo)函數(shù)的優(yōu)化迭代在10次左右，誤差就逐漸穩(wěn)定且趨于0。同時(shí)，二分法–PSO方法每次迭代計(jì)算所用時(shí)間比PSO算法每次迭代計(jì)算所用時(shí)間少。經(jīng)過(guò)簡(jiǎn)單計(jì)算可以看出，改進(jìn)后的方法迭代計(jì)算所用時(shí)間少于傳統(tǒng)方法所用時(shí)間。

由此可以證明采用二分法–PSO算法對(duì)目標(biāo)函數(shù)的優(yōu)化迭代的速度更快。同時(shí)得到表2所示的修正前后的設(shè)計(jì)參數(shù)值。

圖6 目標(biāo)函數(shù)收斂過(guò)程

從表3的結(jié)果可以看出，修正前固有頻率與試驗(yàn)分析獲取的固有頻率間的差異相對(duì)較大，說(shuō)明仿真前的模型在反映試驗(yàn)?zāi)Ｐ偷慕Y(jié)構(gòu)特性時(shí)，誤差相對(duì)較大。而修正后的有限元模型的計(jì)算結(jié)果與試驗(yàn)測(cè)試結(jié)果更相近，即采用PSO方法與基于二分法–PSO方法進(jìn)行修正后的結(jié)果都比初始模型分析結(jié)果的準(zhǔn)確度要高，而且基于二分法–PSO方法修正結(jié)果的誤差比PSO方法修正的誤差更小。說(shuō)明了修正方法及設(shè)計(jì)參數(shù)優(yōu)化結(jié)果的正確性，同時(shí)可以表明基于二分法–PSO的修正方法可以得到更高的計(jì)算精度。

4 模型應(yīng)用

圖7 轉(zhuǎn)子碰摩試驗(yàn)

以圖7所示的試驗(yàn)臺(tái)為仿真研究的原型，進(jìn)行轉(zhuǎn)子系統(tǒng)碰摩故障的研究，其中碰撞位置采用塑料塊與轉(zhuǎn)子圓盤(pán)摩擦模擬碰摩故障。

將修正后的有限元模型導(dǎo)入ADAMS中，構(gòu)建剛?cè)崤己系霓D(zhuǎn)子模型，然后進(jìn)行轉(zhuǎn)子系統(tǒng)碰摩故障的動(dòng)力學(xué)仿真研究。由于試驗(yàn)轉(zhuǎn)子的1階臨界轉(zhuǎn)速為5 160 r/min，偏心量為m?r=2.5×10-4kg?m，碰摩間隙設(shè)置為0.15 mm。在ADAMS中，碰摩故障的研究主要利用IMPACT函數(shù)來(lái)模擬兩個(gè)構(gòu)件之間的碰撞力，IMPACT函數(shù)構(gòu)造的碰撞力由彈性力與阻尼力組成的，由于在試驗(yàn)中，用塑料塊與轉(zhuǎn)子圓盤(pán)摩擦模擬碰摩故障，因此根據(jù)ADAMS推薦，碰摩剛度為1 150 N/m，阻比系數(shù)為0.68，滲透深度為0.1，碰撞指數(shù)為2。然后添加約束條件，即轉(zhuǎn)子與轉(zhuǎn)軸間的固定副，軸與底座間的旋轉(zhuǎn)副。添加一個(gè)恒定轉(zhuǎn)速，大小為9 000(度每秒)，即1 500 r/min來(lái)模擬電機(jī)的驅(qū)動(dòng)，在轉(zhuǎn)子系統(tǒng)的另一邊添加一個(gè)負(fù)載轉(zhuǎn)矩0，為了降低突變帶來(lái)的誤差，因此使用STEP函數(shù)控制，使得速度穩(wěn)定增加到9 000(度每秒)。

圖8為修正前，修正后，實(shí)驗(yàn)測(cè)量獲取的振動(dòng)頻率圖。

從圖中可以看出，修正后模型的振動(dòng)頻率與試驗(yàn)測(cè)量的振動(dòng)頻率都在26 Hz左右，而修正前的振動(dòng)頻率在35 Hz附近，對(duì)比試驗(yàn)結(jié)果，誤差較大，因此可以看出修正后的結(jié)果更為接近試驗(yàn)結(jié)果。

圖10是轉(zhuǎn)速為2 760 r/min時(shí)的頻譜圖與軸心軌跡圖，通過(guò)圖9與圖10可以看出，在轉(zhuǎn)速低于臨界轉(zhuǎn)速時(shí)，存在大于基頻的高倍頻成分，軸心軌跡有內(nèi)凹現(xiàn)象，伴隨著轉(zhuǎn)速的增長(zhǎng)，二倍頻增加的速度遠(yuǎn)大于一倍頻增長(zhǎng)的速度。

表2 修正前后的模型參數(shù)

表3 模型修正前后固有頻率與誤差的變化情況

圖8 碰摩轉(zhuǎn)子修正前后及實(shí)驗(yàn)測(cè)量頻率圖

圖9 臨界轉(zhuǎn)速以?xún)?nèi)的仿真結(jié)果

圖10 臨界轉(zhuǎn)速以?xún)?nèi)的仿真結(jié)果

圖11為臨界轉(zhuǎn)速情況下的振動(dòng)情況，此時(shí)轉(zhuǎn)速在5 160 r/min附近，出現(xiàn)小于一倍頻的成分，基頻幅值遠(yuǎn)大于之前，同時(shí)，軸心軌跡為多環(huán)嵌套；圖12中，轉(zhuǎn)速在6 190 r/min附近，是臨界轉(zhuǎn)速的1.2倍，因此出現(xiàn)了1X/1.2的低頻成分，同時(shí)軸心軌跡為多環(huán)嵌套。

圖11 臨界轉(zhuǎn)速時(shí)仿真結(jié)果

圖12 臨界轉(zhuǎn)速以外仿真結(jié)果

5 結(jié)語(yǔ)

為了提高剛?cè)狁詈限D(zhuǎn)子模型柔性轉(zhuǎn)軸的精度，以及模型修正技術(shù)的計(jì)算精度與計(jì)算效率，本文利用模態(tài)試驗(yàn)數(shù)據(jù)，應(yīng)用二分法–PSO算法實(shí)現(xiàn)有限元模型修正的優(yōu)化計(jì)算，從而提高了有限元模型修正的速度與精度，且避免了在SPO迭代尋優(yōu)過(guò)程中陷入局部最優(yōu)解的問(wèn)題。對(duì)比PSO模型修正的結(jié)果與二分法–PSO模型修正的結(jié)果，發(fā)現(xiàn)基于二分法–PSO的模型修正方法在計(jì)算速度，計(jì)算精度方面相都有所提高。

將修正后的模型應(yīng)用在轉(zhuǎn)子碰摩故障的動(dòng)力學(xué)特性研究中，考慮轉(zhuǎn)速變化對(duì)轉(zhuǎn)子碰摩動(dòng)力學(xué)特性的影響。研究結(jié)果表明，當(dāng)轉(zhuǎn)速小于臨界轉(zhuǎn)速時(shí)，伴隨著轉(zhuǎn)速的增加二倍頻幅值的增長(zhǎng)迅速；在臨界轉(zhuǎn)速時(shí)，振動(dòng)幅值遠(yuǎn)大于先前；當(dāng)轉(zhuǎn)速大于臨界轉(zhuǎn)速時(shí)，出現(xiàn)低于一倍頻的成分。仿真結(jié)果與試驗(yàn)結(jié)果一致性較高，為實(shí)際工程中轉(zhuǎn)子系統(tǒng)的安全穩(wěn)定運(yùn)行提供一定的實(shí)際參考。

[1]MING ZENG,YU YANG,JINDE ZHENG.Normalized complex Teager energy operator demodulation method and its application to fault diagnosis in a rubbing rotor system[J].Mechanical Systems and Signal Processing,2015(50-51):380-399.

[2]FERNANDO L PONTA,ALEJANDRO D OTERO,LUCAS I LAGO,Effects of rotor deformation in windturbine performance:The dynamic rotor deformation bladeelementmomentum model(DRD-BEM)[J].Renewable Energy,2016(92):157-170.

[3]FOX R L,KAPOOR M P.Rate of change of eigenvectors and eigenvalues[J].AIAA Journal,1968,12(6):2426-2429.

[4]徐張明，高天明，沈榮瀛，等.一種改進(jìn)的利用頻響函數(shù)進(jìn)行有限元模型修正的方法[J].振動(dòng)與沖擊，2002，21(3)：43-45.

[5]ASMA F,BOUAZZOUNI A.Finite element model updating using variable separation[J].European Journal of Mechanics-A/Solids,2007,26(4):728-735.

[6]MOHAMMAD R D,JAVAD V A,SIROUS G,et al.Determination of nonlinear behavior of a ball joint system by model updating[J].Journal of Constructional Steel Research,2012,71:52-62.

[7]GUILLAUME ROSENZVEIG,FRAN?OIS LOUF,A FE model updating method for the simulation of the assembly process of large and lightweight aeronautical structures[J].Finite Elements in Analysis and Design,2016(111):56-63.

[8]BAIYONG GOU,WEIJIE ZHANG,QIUHAI LU.A successive selection method for finite element model updating[J].Mechanical Systems and Signal Processing,2016(70-71):320-333.

[9]NEHETE D V,MODAK S V,GUPTA K.Coupled vibroacoustic modelupdating using frequency response functions[J].Mechanical Systems and Signal Processing,2016(70-71):308-319.

[10]NYSSEN F,GOLINVAL J C.Identification of mistuning and model updating of an academic blisk based on geometry and vibration measurements[J].Mechanical Systems and Signal Processing,2016(68-69):252-264.

[11]ZHANG C,SUN J.An alternate two phases particle swarm optimization algorithm for flow shop scheduling problem[J].Expert Systems with Applications,2009,36(3):5162-5167.

[12]慶光蔚.基于復(fù)模態(tài)和頻響修正方法的磁軸承-轉(zhuǎn)子系統(tǒng)參數(shù)識(shí)別[D].南京：南京航空航天大學(xué)，2010.