吳 磊，于慎波，于言明，竇汝桐，翟鳳晨

(1.沈陽(yáng)工業(yè)大學(xué)機(jī)械工程學(xué)院, 沈陽(yáng) 110870；2.中國(guó)郵政集團(tuán)公司遼寧省機(jī)要通信局, 沈陽(yáng) 110000)

0 引言

永磁同步電主軸作為高速機(jī)床的核心部件，將電機(jī)和主軸組成一個(gè)整體，實(shí)現(xiàn)了機(jī)床主軸系統(tǒng)的高效傳動(dòng)[1]。隨著永磁材料的不斷創(chuàng)新與發(fā)展，使得永磁同步電主軸在航空航天、機(jī)械加工和精密制造等領(lǐng)域得到了越來(lái)普遍的應(yīng)用[2]。永磁同步電主軸轉(zhuǎn)子細(xì)長(zhǎng)，且轉(zhuǎn)速較高，和其它轉(zhuǎn)子相比，其轉(zhuǎn)子系統(tǒng)接近臨界轉(zhuǎn)速的機(jī)率更大。另外，由于外部環(huán)境、負(fù)載的變化、加工誤差等原因會(huì)導(dǎo)致轉(zhuǎn)子發(fā)生動(dòng)態(tài)偏心[3]，進(jìn)而產(chǎn)生不平衡磁拉力。在正常運(yùn)行過(guò)程中，電主軸轉(zhuǎn)子固有頻率和轉(zhuǎn)子不平衡力導(dǎo)致的強(qiáng)迫振動(dòng)頻率相耦合，當(dāng)轉(zhuǎn)速接近該階模態(tài)下的臨界轉(zhuǎn)速時(shí)，主軸將會(huì)發(fā)生多重共振，導(dǎo)致噪聲過(guò)大，嚴(yán)重時(shí)，甚至導(dǎo)致轉(zhuǎn)子不穩(wěn)定或損壞[4]。為了避免轉(zhuǎn)子在運(yùn)轉(zhuǎn)中產(chǎn)生共振，影響電主軸正常運(yùn)行，需要對(duì)轉(zhuǎn)子偏心引起的不平衡磁拉力以及轉(zhuǎn)子臨界轉(zhuǎn)速和穩(wěn)定性進(jìn)行相關(guān)的動(dòng)力學(xué)分析[5-6]。

文獻(xiàn)[7-9]在轉(zhuǎn)子系統(tǒng)的不平衡響應(yīng)分析、臨界轉(zhuǎn)速和穩(wěn)定性等領(lǐng)域做了許多的研究，為旋轉(zhuǎn)機(jī)械的發(fā)展奠定了理論基礎(chǔ)，但主要集中在水、汽輪機(jī)，對(duì)永磁同步電主軸的研究涉及很少。目前求解臨界轉(zhuǎn)速常采用傳統(tǒng)Prohl傳遞矩陣法、Riccati傳遞矩陣法和有限元法等[10]，其中Prohl傳遞矩陣法應(yīng)用在各項(xiàng)同性轉(zhuǎn)子系統(tǒng)，Riccati應(yīng)用在各項(xiàng)異性轉(zhuǎn)子系統(tǒng)，但由于復(fù)雜的外部環(huán)境以及算法本身的缺陷，計(jì)算時(shí)會(huì)出現(xiàn)根的溢出和丟根的現(xiàn)象，導(dǎo)致臨界轉(zhuǎn)速計(jì)算不精確。有限元法計(jì)算轉(zhuǎn)子系統(tǒng)臨界轉(zhuǎn)速時(shí)，對(duì)于復(fù)雜轉(zhuǎn)子由于建模復(fù)雜，占用計(jì)算機(jī)內(nèi)存較大，計(jì)算速度較慢，在一定程度上限制了其應(yīng)用范圍。因此，如何提高解析法計(jì)算轉(zhuǎn)子臨界轉(zhuǎn)速的精確度以及進(jìn)行機(jī)電耦合的分析，成為了轉(zhuǎn)子動(dòng)力學(xué)研究熱點(diǎn)。

為了解決上述問(wèn)題，本文中提出了一種用加時(shí)間因t的方法構(gòu)建轉(zhuǎn)子傳遞矩陣，運(yùn)用了雙重QR法[11]求解傳遞矩陣特征值，進(jìn)而可求得轉(zhuǎn)子系統(tǒng)的臨界轉(zhuǎn)速，既提高了計(jì)算速度，又解決了根的上下溢出問(wèn)題。同時(shí)運(yùn)用編程計(jì)算了轉(zhuǎn)子運(yùn)行中由于動(dòng)態(tài)偏心而導(dǎo)致的不平衡磁拉力，使其與機(jī)械系統(tǒng)進(jìn)行耦合，分析其對(duì)動(dòng)力學(xué)特性的影響，對(duì)研究轉(zhuǎn)子系統(tǒng)動(dòng)力學(xué)特性有著重要的參考價(jià)值。

1 永磁同步電主軸臨界轉(zhuǎn)速計(jì)算

1.1 轉(zhuǎn)子離散建模與傳遞矩陣的建立

根據(jù)文獻(xiàn)[12]中節(jié)點(diǎn)劃分原則，應(yīng)用集總質(zhì)量法對(duì)所求轉(zhuǎn)子進(jìn)行離散處理[8]。本文中，永磁同步電主軸轉(zhuǎn)子如圖1所示，將其簡(jiǎn)化為25個(gè)軸段，26個(gè)節(jié)點(diǎn)(節(jié)點(diǎn)7和19為軸承支撐點(diǎn))并在對(duì)應(yīng)節(jié)點(diǎn)附加上相應(yīng)的質(zhì)量和轉(zhuǎn)動(dòng)慣量。

圖1 電主軸轉(zhuǎn)子總成模型

對(duì)于離散后的第i個(gè)單元，將i和i+1分別定義為其左右兩端的截面編號(hào)，根據(jù)轉(zhuǎn)子受力特點(diǎn)，考慮軸承支承為各項(xiàng)異性，將截面的x方向和y方向彎矩、剪力、位移和轉(zhuǎn)角組合成該截面的狀態(tài)向量，如式(1)所示：

(1)

式中，Mx、My、θx、θy、x、y、Qx、Qy分別為水平和垂直方向的彎矩、轉(zhuǎn)角、撓度和剪力。

引入以復(fù)變量表示的振動(dòng)量w=τ+σi，其中，w為復(fù)渦動(dòng)頻率；實(shí)部τ為衰減指數(shù)；虛部σ為阻尼圓周率。根據(jù)圓盤的受力分析，可得出加時(shí)間因子t的剛性薄圓盤的傳遞矩陣，如式(2)所示：

(2)

根據(jù)軸段左右截面狀態(tài)向量關(guān)系，可以得出等截面彈性軸的傳遞矩陣，如式(3)所示：

(3)

將狀態(tài)向量{Z}i分別寫成力向量{Z}g和位移向量{Z}h，根據(jù)邊界條件以及轉(zhuǎn)子系統(tǒng)各截面狀態(tài)向量的關(guān)系，可得出方程(4)：

(4)

其中，[Q] =[B]K[D]K[B]K-1[D]K-1…[B]1[D]1。

由方程(4)可得：

[Q12]{Z}h,1=0

(5)

由于{Z}h,1≠0時(shí)，上式有非零解的條件為：

|Q12|=0

(6)

由于式(6)中的元素為多項(xiàng)式，運(yùn)用代數(shù)余子式法將其展開，頻率方程實(shí)際上是關(guān)于w的8K階一元高次方程，即：

e0+e1w+e2w2+···+e8K-1w8K-1+w8K=0

(7)

求方程(7)的根可以轉(zhuǎn)化為求其對(duì)應(yīng)友矩陣全部特征值的問(wèn)題。

電主軸轉(zhuǎn)子為細(xì)長(zhǎng)階梯軸，為了準(zhǔn)確模擬實(shí)際轉(zhuǎn)子工況，需要離散足夠的節(jié)點(diǎn)數(shù)，本文中離散節(jié)點(diǎn)數(shù)K=26，方程最高階會(huì)達(dá)到208，求解過(guò)程中會(huì)出現(xiàn)根的溢出，求解難度特別大。文中，為了縮小求解過(guò)程中大數(shù)與小數(shù)之間的差異，防止根的溢出，提出加時(shí)間因子t的方法建立轉(zhuǎn)子系統(tǒng)傳遞矩陣。計(jì)算過(guò)程中其他參數(shù)與時(shí)間因子相乘，轉(zhuǎn)速與時(shí)間因子相除，運(yùn)用雙重步QR法編程求解出復(fù)頻率的實(shí)部與虛部，計(jì)算完成后將虛部乘以時(shí)間因子得到臨界轉(zhuǎn)速。

1.2 解析法計(jì)算永磁同步電主軸臨界轉(zhuǎn)速

本文中以額定功率9.5 kW永磁同步電主軸轉(zhuǎn)子為例，取時(shí)間因子t等于15 000構(gòu)建傳遞矩陣，取3000 r/min為轉(zhuǎn)速步長(zhǎng)，運(yùn)用編程的方法求解出復(fù)頻率并畫出坎貝爾圖。同頻線與各階渦動(dòng)速度曲線的交點(diǎn)即是轉(zhuǎn)子系統(tǒng)臨界轉(zhuǎn)速，如圖2所示。從圖中可以看出，隨著轉(zhuǎn)子自轉(zhuǎn)速度逐步升高，一階正向和反向渦動(dòng)速度相近且變化很小；二階反向渦動(dòng)速度逐步減小，正向渦動(dòng)速度逐步增大。

圖2 解析法計(jì)算轉(zhuǎn)子坎貝爾圖

1.3 有限元法計(jì)算永磁同步電主軸臨界轉(zhuǎn)速

利用ANSYS計(jì)算臨界轉(zhuǎn)速時(shí)，運(yùn)用實(shí)體模型建模復(fù)雜，且計(jì)算結(jié)果精度達(dá)不到要求。因此本文運(yùn)用有限元法時(shí)，以MASS21、BEAM188和COMBI214分別模擬質(zhì)量點(diǎn)、彈性軸和軸承單元。同時(shí)質(zhì)量點(diǎn)單元、軸承單元分別施加旋轉(zhuǎn)約束和全約束，并充分考慮了陀螺效應(yīng)的影響，通過(guò)編寫ANSYS APDL命令流的方式建立二維有限元模型如圖3所示。

圖3 ANSYS二維質(zhì)量點(diǎn)單元模型

本文用QR阻尼法進(jìn)行模態(tài)提取，利用PLCAMP命令畫出軸承-轉(zhuǎn)子系統(tǒng)的CAMP-BELL圖，如圖4所示。

圖4 ANSYS有限元計(jì)算轉(zhuǎn)子坎貝爾圖

1.4 解析法和有限元法計(jì)算臨界轉(zhuǎn)速結(jié)果對(duì)照

各階臨界轉(zhuǎn)速為同頻線與正向渦動(dòng)速度的交點(diǎn)。轉(zhuǎn)子臨界轉(zhuǎn)速計(jì)算結(jié)果對(duì)比如表1所示。由表1可知，解析法與有限元計(jì)算結(jié)果對(duì)照，誤差均在3.2%以內(nèi)，驗(yàn)證了本文中解析法計(jì)算臨界轉(zhuǎn)速的有效性。

表1 臨界轉(zhuǎn)速結(jié)果對(duì)比 (r/min)

2 穩(wěn)定性分析

穩(wěn)定性裕度體現(xiàn)在系統(tǒng)自身的修復(fù)能力，也就是能夠回到平穩(wěn)狀態(tài)的能力，度量方法主要有對(duì)數(shù)衰減率、系統(tǒng)阻尼和系統(tǒng)抗定常干擾界限值法，其中對(duì)數(shù)衰減率是以振幅衰減的快慢程度為指標(biāo)，對(duì)系統(tǒng)受干擾后恢復(fù)到平衡狀態(tài)過(guò)程的度量[8]。

本文以對(duì)數(shù)衰減率為依據(jù)判定軸系的穩(wěn)定性裕度。由本文第一部分解析法得到復(fù)頻率w=τ+σi，則第i階特征值對(duì)應(yīng)的對(duì)數(shù)衰減率：

(8)

γi越大表明系統(tǒng)穩(wěn)定程度越高，當(dāng)γi<0時(shí)，表明系統(tǒng)處于失穩(wěn)狀態(tài)。

計(jì)算得到的軸系對(duì)數(shù)衰減率曲線如圖5所示。圖中表明，一階反向、二階正向和三階反向?qū)?shù)衰減率呈逐漸下降趨勢(shì)，其中三階反向?qū)?shù)衰減率最小且有加速下降的趨勢(shì)，最易引起系統(tǒng)失穩(wěn)，設(shè)計(jì)時(shí)應(yīng)盡量避開三階反向渦動(dòng)臨界轉(zhuǎn)速。

圖5 永磁同步電主軸轉(zhuǎn)子對(duì)數(shù)衰減率

3 不平衡磁拉力

3.1 不平衡磁拉力的計(jì)算

通常氣隙磁場(chǎng)分布均勻，轉(zhuǎn)子不承受不平衡磁拉力，但由于轉(zhuǎn)子偏心導(dǎo)致徑向氣隙磁場(chǎng)畸變，會(huì)產(chǎn)生不平衡磁拉力，使得轉(zhuǎn)子受力不均勻。傳統(tǒng)公式計(jì)算不平衡磁拉力時(shí)，將模型簡(jiǎn)化為理想狀態(tài)，計(jì)算結(jié)果不夠精確。因此本文基于Ansoft有限元軟件對(duì)某型號(hào)8極36槽永磁同步電主軸轉(zhuǎn)子動(dòng)態(tài)偏心引起的不平衡磁拉力進(jìn)行分析。電主軸的基本參數(shù)如表2所示。

表2 永磁同步電主軸部分參數(shù)

電主軸有限元模型建立完成并求解后可提取氣隙磁密，然后由式(9)求出不平衡磁拉力[13-15]。

(9)

其中，F(xiàn)x、Fy分別是不平衡磁拉力在X和Y方向上的分量；l是轉(zhuǎn)子硅鋼片的軸向長(zhǎng)度；rp是不平衡磁拉力積分半徑。

轉(zhuǎn)子偏心后，徑向畸變氣隙磁場(chǎng)可看作由多個(gè)同心機(jī)磁場(chǎng)疊加而成，且切向氣隙磁密對(duì)徑向不平衡磁拉力影響極小，可忽略不計(jì)。因此可以將氣隙積分圓周離散成N個(gè)相等的子區(qū)域，求得各區(qū)域徑向畸變氣隙磁場(chǎng)產(chǎn)生的不平衡磁拉力然后進(jìn)行疊加，所以式(9)可以表示成式(10)的形式。

(10)

式中，δiω、δiρ分別為第i個(gè)區(qū)域的起始和結(jié)束角度。

則不平衡磁拉力的合力為：

(11)

為了防止轉(zhuǎn)子和定子出現(xiàn)掃膛現(xiàn)象，轉(zhuǎn)子偏心距不應(yīng)大于平均氣隙的10%[16]。本文中電主軸氣隙總長(zhǎng)為1.6 mm，因此偏心距不應(yīng)超過(guò)0.16 mm。根據(jù)此原則，以無(wú)偏心和動(dòng)態(tài)偏心0.16 mm為例，將模型積分圓周劃分為1000個(gè)區(qū)域，分別提取其徑向氣隙磁密，如圖6、圖7所示。

圖6 無(wú)偏心徑向氣隙磁密

圖7 偏心距0.16 mm徑向氣隙磁密

應(yīng)用式(10)，利用編程求出不平衡磁拉力在X和Y軸的分力，如圖8所示。圖中可看出，當(dāng)無(wú)偏心時(shí)，X和Y方向不平衡電磁力均在0 N附近上下波動(dòng)，且浮動(dòng)很??；當(dāng)偏心距為0.16 mm時(shí)，X和Y方向上均產(chǎn)生了隨時(shí)間周期性變化的幅值為230 N的不平衡磁拉力。

圖8 不平衡磁拉力在X和Y軸的分力

應(yīng)用上述方法，可求出偏心距分別為0.032 mm、0.064 mm、0.096 mm和0.128 mm時(shí)在X和Y方向的不平衡磁拉力，由式(10)、式(11)可求得不同偏心距下不平衡磁拉力在Y方向分量Fy和在X方向分量Fx隨圓周的變化以及不平衡磁拉力的合力，如圖9、圖10所示。

圖9 Fy隨Fx在圓周內(nèi)的變化

圖10 不平衡磁拉力合力隨時(shí)間變化

由圖9可看出，相同偏心距下，不同時(shí)刻的X和Y方向的不平衡磁拉力在圓周上組成了一個(gè)同心圓，且圓的半徑即為該偏心狀態(tài)下的不平衡磁拉力，這是由于動(dòng)態(tài)偏心時(shí)最小氣隙隨著時(shí)間在變化。由圖9和圖10可得到不同偏心距下不平衡磁拉力的數(shù)值如表3所示，且不平衡磁拉力隨著偏心距的增大呈線性增加。

表3 不同偏心距下的不平衡磁拉力

3.2 機(jī)電耦合

當(dāng)電主軸轉(zhuǎn)子在高速旋轉(zhuǎn)時(shí)，不平衡磁拉力被引入機(jī)械系統(tǒng)引發(fā)共振[17]，嚴(yán)重情況下會(huì)造成定轉(zhuǎn)子接觸，刮蹭永磁體。因此在動(dòng)力學(xué)分析時(shí)，考慮不平衡磁拉力的影響十分必要。將不同動(dòng)態(tài)偏心下產(chǎn)生的不平衡磁拉力與轉(zhuǎn)子機(jī)械系統(tǒng)進(jìn)行耦合后，臨界轉(zhuǎn)速的變化如圖11所示。

圖11 不平衡磁拉力對(duì)臨界轉(zhuǎn)速的影響

通過(guò)分析圖11可以看出，考慮不平衡磁拉力時(shí)的臨界轉(zhuǎn)速低于無(wú)不平衡磁拉力時(shí)的臨界轉(zhuǎn)速。隨著不平衡磁拉力的增加，一階臨界轉(zhuǎn)速由23 657 r/min降低到23 339 r/min，降幅為1.4%，二階臨界轉(zhuǎn)速幾乎無(wú)明顯變化。

4 結(jié)論

(1)運(yùn)用加時(shí)間因子t的方法構(gòu)建轉(zhuǎn)子傳遞矩陣，同時(shí)運(yùn)用雙重步QR法計(jì)算轉(zhuǎn)子系統(tǒng)臨界轉(zhuǎn)速，既提高了計(jì)算精度又解決了根的上下溢出問(wèn)題。同時(shí)將解析法與有限元法計(jì)算的各階臨界轉(zhuǎn)速進(jìn)行對(duì)比，誤差均小于3.2%，驗(yàn)證了解析法的有效性。

(2)通過(guò)分析對(duì)數(shù)衰減率曲線，發(fā)現(xiàn)三階反向渦動(dòng)衰減率曲線最低，且下降速度有逐漸加快的趨勢(shì)。說(shuō)明三階反向渦動(dòng)穩(wěn)定性裕度最小，最易引起轉(zhuǎn)子系統(tǒng)失穩(wěn)。

(3)分析了動(dòng)態(tài)偏心工況下的不平衡磁拉力，結(jié)果表明，不平衡磁拉力隨著偏心距的增大呈線性增加。隨著不平衡磁拉力增加，一階臨界轉(zhuǎn)速將下降，二階臨界轉(zhuǎn)速幾乎不受影響。