嚴(yán)國軍, 楊 彥, 顧建華, 姚月琴,祁 淼,4
(1. 鹽城工業(yè)職業(yè)技術(shù)學(xué)院 機(jī)電工程學(xué)院,江蘇 鹽城 224000; 2. 東南大學(xué) 機(jī)械工程學(xué)院,江蘇 南京 210000; 3. 南京理工大學(xué) 控制工程學(xué)院,江蘇 南京 210000;4. 江蘇大學(xué) 汽車工程學(xué)院,江蘇 鎮(zhèn)江 212000)
在我國,電動(dòng)農(nóng)業(yè)機(jī)械還處于研究階段。無論從市場(chǎng)需求,還是能源發(fā)展方向等方面考慮,電動(dòng)農(nóng)業(yè)機(jī)械都是我國新農(nóng)村建設(shè)和現(xiàn)代農(nóng)業(yè)發(fā)展的需要。我國“十二五”規(guī)劃曾明確提出:“提高設(shè)施裝備水平,提升產(chǎn)業(yè)發(fā)展質(zhì)量,促進(jìn)設(shè)施農(nóng)業(yè)可持續(xù)發(fā)展”。隨著新能源車輛的發(fā)展,作為農(nóng)業(yè)大國,純電動(dòng)拖拉機(jī)正在慢慢替代傳統(tǒng)燃油拖拉機(jī),并且應(yīng)用在農(nóng)業(yè)的各個(gè)方面。
純電動(dòng)車輛動(dòng)力傳動(dòng)系統(tǒng)中普遍存在機(jī)電耦合現(xiàn)象,研究機(jī)電系統(tǒng)耦合振動(dòng),可減小系統(tǒng)非線性特性的影響。純電動(dòng)車輛動(dòng)力傳動(dòng)系統(tǒng)的扭轉(zhuǎn)振動(dòng)問題在電磁激勵(lì)和外部激勵(lì)的作用下變得復(fù)雜與突出。在機(jī)電耦合方面,林利紅等人對(duì)純電動(dòng)車輛動(dòng)力傳動(dòng)系統(tǒng)進(jìn)行機(jī)電耦合動(dòng)力學(xué)分析和實(shí)驗(yàn)研究,建立了系統(tǒng)的機(jī)電耦合振動(dòng)仿真模型,仿真分析了由于電流調(diào)節(jié)器參數(shù)、阻尼、諧波擾動(dòng)、間隙以及負(fù)載擾動(dòng)等因素引起的機(jī)電耦合振動(dòng)動(dòng)態(tài)過程,并將物理實(shí)驗(yàn)研究與理論分析結(jié)果進(jìn)行對(duì)比和驗(yàn)證[1]。在車用電機(jī)非線性扭轉(zhuǎn)振動(dòng)方面,陳星等人采用多尺度法研究了純電動(dòng)車輛動(dòng)力傳動(dòng)系統(tǒng)的扭轉(zhuǎn)振動(dòng)系統(tǒng)特性,仿真結(jié)果表明合理的設(shè)計(jì)和控制電機(jī)的電磁參數(shù)和機(jī)械參數(shù),可以解決電機(jī)扭振系統(tǒng)出現(xiàn)的跳躍分岔等非線性特性現(xiàn)象[2-3]。
針對(duì)低速重載起步工況下純電動(dòng)拖拉機(jī)機(jī)電耦合傳動(dòng)系統(tǒng)受到非周期沖擊干擾力作用,分析扭轉(zhuǎn)角與電磁轉(zhuǎn)矩的耦合影響。為此筆者通過建立純電動(dòng)拖拉機(jī)機(jī)電耦合傳動(dòng)系統(tǒng)的非線性扭轉(zhuǎn)振動(dòng)模型,求解無擾動(dòng)Hamilton系統(tǒng)的平衡點(diǎn);利用Melnikov法求解分岔閾值曲線,分析電機(jī)控制參數(shù)對(duì)純電動(dòng)拖拉機(jī)機(jī)電耦合傳動(dòng)系統(tǒng)非線性動(dòng)力學(xué)行為的影響;基于Melnikov理論分析研究傳動(dòng)系統(tǒng)的分岔和混沌特性,從而求解系統(tǒng)穩(wěn)定性運(yùn)行邊界條件。
不同情況下,永磁同步電機(jī)的穩(wěn)定運(yùn)行時(shí)的相量如圖1。
圖1 永磁同步電機(jī)相量Fig. 1 Phasor diagram of permanent magnet synchronous motor
由圖1得到如下關(guān)系式:
(1)
通過求解式(1)可得到電機(jī)定子電流的直、交軸分量:
(2)
定子相電流:
(3)
而電機(jī)的輸入功率:
P1=mUI1cosφ
(4)
當(dāng)忽略定子電阻,由式(4)可得到電機(jī)的電磁功率:
(5)
則電機(jī)的電磁轉(zhuǎn)矩:
(6)
式中:m為電機(jī)相數(shù);Ω為機(jī)械角速度,且ω=pΩ,其中,ω為電機(jī)的電角速度,p為電機(jī)的極對(duì)數(shù);E0為空載反電動(dòng)勢(shì),可由式(7)計(jì)算得到。
(7)
式中:f為電源頻率;Kdp為基波繞組系數(shù);N為每相串聯(lián)匝數(shù);Φ10為永磁體基波磁通;Am為永磁體面積;Br為永磁體剩磁;σ0為空載漏磁系數(shù);bm0為永磁體空載工作點(diǎn)假設(shè)值;KΦ為氣隙磁通波形系數(shù);αi為計(jì)極弧系數(shù);αp為極弧系數(shù),δ為氣隙長度,D為定子內(nèi)徑。
Te=b1sinθ+b2sin 2θ
(8)
由圖1可知,θ=φ+ψ,因此可得到sinθ和sin 2θ三角函數(shù)展開式。
利用sinx≈x-1/6x3,cosx≈1-1/2x2,泰勒展開的近似式將含φ的三角函數(shù)展開并化簡,即可得到作用于轉(zhuǎn)子的電磁轉(zhuǎn)矩:
Te=k0+k1φ-k2φ2-k3φ3
(9)
式中:k0=b1sinψ+b2cos 2ψ;k1=b1cosψ+2b2sinψ;k2=1/2b1sinψ+2b2sin 2ψ;k3=1/6b1cosψ+4/3b2cos 2ψ。
通過分析車載永磁同步電機(jī)機(jī)電耦合扭轉(zhuǎn)振動(dòng)特性,將轉(zhuǎn)子系統(tǒng)和傳動(dòng)系統(tǒng)進(jìn)行簡化,可以等效為二質(zhì)體機(jī)電耦合轉(zhuǎn)子模型[4-5],如圖2。圖中永磁同步電機(jī)轉(zhuǎn)子與機(jī)械轉(zhuǎn)子之間彈性聯(lián)軸器聯(lián)接,電機(jī)機(jī)座可視為剛性。圖中Te是電磁轉(zhuǎn)矩,TL是負(fù)載轉(zhuǎn)矩。
圖2 二質(zhì)體機(jī)電耦合轉(zhuǎn)子系統(tǒng)Fig. 2 Two plastid electromechanical coupled rotor system
根據(jù)牛頓定律,二質(zhì)體轉(zhuǎn)子系統(tǒng)扭轉(zhuǎn)振動(dòng)動(dòng)力學(xué)方程為:
(10)
(11)
當(dāng)轉(zhuǎn)子系統(tǒng)處于穩(wěn)態(tài)運(yùn)行時(shí),式(10)通過方程變換并消除得到公式中的靜態(tài)項(xiàng):
(12)
則機(jī)電耦合扭轉(zhuǎn)振動(dòng)動(dòng)力學(xué)方程可轉(zhuǎn)化為無量綱形式:
(13)
式中:v=J1/J2;r=1/(1+v);μ=(1/J1+1/J2)C;β=k2r2/J1;γ=k3r3/J1;η2=ω0-k1r/J1,ω0為不考慮機(jī)電耦合時(shí)轉(zhuǎn)子系統(tǒng)的機(jī)械固有頻率,可由式(14)進(jìn)行計(jì)算:
(14)
從式(14)可看出,動(dòng)力學(xué)方程具有平方項(xiàng)和立方項(xiàng),方程為典型的非線性動(dòng)力學(xué)方程。因此,考慮到機(jī)電耦合效應(yīng),傳動(dòng)系統(tǒng)轉(zhuǎn)子扭轉(zhuǎn)振動(dòng)方程將有明顯的非線性特征。
選取的拖拉機(jī)電機(jī)額定功率為7.5 kW,額定電壓為96 V,其他相關(guān)電機(jī)參數(shù)如表1。
表1 電機(jī)參數(shù)表Table 1 Motor parameter list
式中:ε為小參數(shù),且0<ε?1。
無擾動(dòng)系統(tǒng)的勢(shì)函數(shù)在研究非線性動(dòng)力學(xué)行為中具有十分重要的作用,因此需要討論無擾動(dòng)系統(tǒng)的函數(shù)。當(dāng)ε=0時(shí),排除阻尼項(xiàng)和擾動(dòng)項(xiàng),函數(shù)方程變?yōu)闊o擾動(dòng)系統(tǒng):
(16)
在排除阻尼項(xiàng)和擾動(dòng)項(xiàng)后,轉(zhuǎn)子動(dòng)力學(xué)系統(tǒng)為Hamilton系統(tǒng),其Hamilton方程為:
(17)
式(17)對(duì)應(yīng)的勢(shì)函數(shù)為:
(18)
Hamilton系統(tǒng)平衡點(diǎn)穩(wěn)定性分為如下幾種情況:
1)情況A
當(dāng)β2≤4γ(1-κ),Hamilton系統(tǒng)只有一個(gè)平衡點(diǎn)P0(0,0):
①情況A1:當(dāng)κ<1時(shí),P0是中心點(diǎn)。
②情況A2:當(dāng)κ>1時(shí),P0是鞍點(diǎn)。
2)情況B
當(dāng)β2>4γ(1-κ)和4γ(1-κ)>0時(shí),有3個(gè)平衡點(diǎn),分別為P0(0,0),P1(x1,0),P2(x2,0):
①情況B1:當(dāng)γ>0時(shí),P1(x1,0)是鞍點(diǎn),P0(0,0)和P2(x2,0)是中心點(diǎn)。當(dāng)β=0.8,κ=0.35,γ=0.2時(shí),可以得到勢(shì)函數(shù)和相軌跡如圖3(a)和圖3(b),從圖3(a)和圖3(b)中可看出存在橫截同宿軌道。
②情況B2:當(dāng)γ<0時(shí),P1(x1,0)是中心點(diǎn),P0(0,0)和P2(x2,0)是鞍點(diǎn)。當(dāng)β=0.8,κ=1.65,γ=-0.2時(shí),可得到勢(shì)函數(shù)和相軌跡如圖3(c)和3圖(d),從圖3(c)和3圖(d)中可以看出存在橫截異宿軌道。
圖3 B1和B2情況下Hamilton系統(tǒng)勢(shì)函數(shù)和相軌跡Fig. 3 Potential functions and phase trajectories of Hamilton systems in the case of B1 and B2
3)情況C
當(dāng)β2>4γ(1-κ)和4γ(1-κ)<0時(shí),有3個(gè)平衡0點(diǎn),分別為P0(0,0),P1(x1,0),P2(x2,0):
①情況C1:當(dāng)γ>0時(shí),P0(0,0)是鞍點(diǎn),P1(x1,0)和P2(x2,0)是中心點(diǎn)。當(dāng)β=0,κ=1.5,γ=0.2時(shí),可得到勢(shì)函數(shù)和相軌跡如圖4(a)和圖4(b),從圖4(a)和圖4(b)中可以看出存在橫截同宿軌道。
②情況C2:當(dāng)γ<0時(shí),P0(0,0)是中心點(diǎn),P1(x1,0)和P2(x2,0)是鞍點(diǎn)。當(dāng)β=0.8,κ=1.65,γ=-0.2時(shí),可得到勢(shì)函數(shù)和相軌跡如圖4(c)和圖4(d),從圖中可以看出存在橫截異宿軌道。
圖4 C1和C2情況下Hamilton系統(tǒng)勢(shì)函數(shù)和相軌跡Fig. 4 Potential functions and phase trajectories of Hamilton systems in the case of C1 and C2
綜合以上分析可知:情況A中只存在一個(gè)平衡點(diǎn),將不做討論;情況B和C均具有3個(gè)平衡點(diǎn),而在情況B1和C1中,對(duì)較小的擾動(dòng)幅值,Hamilton系統(tǒng)具有封閉軌道,此時(shí)將產(chǎn)生橫截同宿軌道,可能導(dǎo)致同宿混沌;在情況B2和C2中,增加激勵(lì)幅值,將會(huì)產(chǎn)生橫截異宿軌道,可能導(dǎo)致異宿混沌。
從圖3和圖4中可以看出,考慮機(jī)電耦合作用的Hamilton系統(tǒng)是雙勢(shì)對(duì)稱或不對(duì)稱波動(dòng)系統(tǒng)原型,同時(shí)由于機(jī)電耦合動(dòng)力學(xué)方程系數(shù)的復(fù)雜多變,所以機(jī)電系統(tǒng)扭轉(zhuǎn)振動(dòng)將會(huì)表現(xiàn)出復(fù)雜的非線性動(dòng)力學(xué)行為。車用永磁同步電機(jī)傳動(dòng)系統(tǒng)為獲得更高的傳動(dòng)范圍和高功率密度,通常設(shè)計(jì)和控制于弱磁效應(yīng),因此,筆者只討論ψ>0的情況。當(dāng)內(nèi)功率因數(shù)角ψ在[0,π]內(nèi)變化時(shí),Hamilton系統(tǒng)的情況如圖5。分析圖可知,系統(tǒng)只存在B1和C1兩種情況,筆者將以C1為例進(jìn)行進(jìn)一步討論(B1可采用相似的過程進(jìn)行計(jì)算討論)。
圖5 機(jī)電傳動(dòng)系統(tǒng)的平衡點(diǎn)的穩(wěn)定性分類Fig. 5 Stability classification of the equilibrium points ofelectromechanical transmission system
已知例C1中,P0(0,0)是鞍點(diǎn),P1(x1,0)和P2(x2,0)是中心點(diǎn)。系統(tǒng)Hamilton量為:
(19)
當(dāng)H(0,0)=0時(shí),存在兩條連接鞍點(diǎn)(0,0)的同宿軌道{q0(t)}。同宿軌道求解過程如式(20):
(20)
由式(20)可得:
(21)
積分得:
(22)
當(dāng)t=0,x0=±1時(shí),c=0。整理后得兩條同宿軌的參數(shù)方程為:
(23)
當(dāng)擾動(dòng)項(xiàng)增加到Hamilton系統(tǒng)時(shí),封閉的同宿軌將被破壞,此時(shí)可能會(huì)導(dǎo)致分岔,全局同宿分岔預(yù)示著混沌運(yùn)動(dòng)。Melnikov方法可以估算存在同宿分岔的參數(shù)空間,并可以轉(zhuǎn)化為混沌,從而確定機(jī)電耦合扭轉(zhuǎn)振動(dòng)混沌運(yùn)動(dòng)的閾值[6-7]。因此,在分析同宿分岔時(shí),將系統(tǒng)(15)寫成狀態(tài)方程的形式為:
(24)
Melnikov方法實(shí)際上是通過計(jì)算雙曲不動(dòng)點(diǎn)的穩(wěn)定和不穩(wěn)定流形之間的距離來確定系統(tǒng)是否存在橫截同宿點(diǎn)。根據(jù)Smale-birkhoff同宿理論,當(dāng)μ≠0,f≠0,且足夠小時(shí),系統(tǒng)(15)可能存在Smale馬蹄存在意義下的混沌。根據(jù)Melinkov理論,Melnikov方程為:
(25)
式中:q0(τ)=(x±,y±),代表Hamilton系統(tǒng)的同宿軌道;F∧G=f1g2-f2g1。
將式(24)代入式(25)中,初始條件為t=0,x0=±1,y0=0時(shí),y±(τ)是時(shí)間τ的奇函數(shù)。則可得到:
(26)
式中,I1和I2分別為部分積分:
所以,對(duì)于一定頻率,當(dāng)滿足式(27)時(shí),系統(tǒng)可能產(chǎn)生Smale馬蹄變換意義下的混沌。
(27)
式(26)第一項(xiàng)積分I1可通過數(shù)學(xué)軟件計(jì)算直接求得I1=0.372 7;第二項(xiàng)積分I2在求解過程相當(dāng)復(fù)雜,所以采用數(shù)值積分法,通過采用數(shù)值積分法畫出混沌閾值曲線,然后根據(jù)閾值曲線同樣能夠?qū)ο到y(tǒng)做定性和定量分析。Melnikov函數(shù)數(shù)值積分法求解過程如下:
(28)
對(duì)式(28)分離變量,并對(duì)兩邊同時(shí)積分,得到:
(29)
將式(28)和式(29)代入I2中得到:
(30)
式(30)可以計(jì)算得到對(duì)應(yīng)ω值的I2值,代入式(31)即可得出隨ω變化的Melnikov閾值曲線。以f作為控制參數(shù)時(shí),Melnikov閾值曲線如圖6。在一定條件下,在f-ω平面上存在混沌區(qū)域。當(dāng)f在該區(qū)域取值時(shí),系統(tǒng)存在Smale馬蹄變換意義下的混沌。
圖6 同宿分岔的閾值曲線(μ=2.5)Fig. 6 Threshold curve of homoclinic bifurcation (μ=2.5)
利用Runge-Kutta法對(duì)擾動(dòng)系統(tǒng)(15)進(jìn)行數(shù)值仿真,來進(jìn)行驗(yàn)證上述用Melnikov方法做出的判斷是否正確,從而對(duì)得到對(duì)應(yīng)的時(shí)域波形、相軌跡、龐加萊截面圖和幅值譜進(jìn)行分析。圖7為μ=2.5,ω=1.0,ε=0.1,時(shí),F(xiàn)=εf,以f為控制參數(shù)時(shí)的系統(tǒng)分岔圖,初始條件為x=0,y=0。
圖7 f為控制參數(shù)時(shí)的系統(tǒng)分叉Fig. 7 Bifurcation diagram of the system when f is thecontrol parameter
從圖7中可以看出,當(dāng)f為控制參數(shù)時(shí),機(jī)電耦合非線性扭振系統(tǒng)表現(xiàn)出復(fù)雜的非線性振動(dòng)現(xiàn)象[8]:當(dāng)激勵(lì)值f較小時(shí)(f<3.6),系統(tǒng)為1周期振蕩;當(dāng)激勵(lì)值f增大到3.6~3.8區(qū)間,2周期振蕩取代1周期振蕩;當(dāng)f>3.8時(shí),控制參數(shù)f進(jìn)一步增加,系統(tǒng)在倍周期分岔后進(jìn)入混沌運(yùn)動(dòng);當(dāng)f>3.8時(shí),機(jī)電耦合系統(tǒng)表現(xiàn)為周期運(yùn)動(dòng)與混沌運(yùn)動(dòng)交替存在。
便于更好的理解,選取3組數(shù)值對(duì)周期運(yùn)動(dòng)變?yōu)榛煦邕\(yùn)動(dòng)的過程進(jìn)行具體分析系統(tǒng),分別選取控制參數(shù)f=2,f=3.2,f=4,如圖6(b)。當(dāng)ω=1.0時(shí),控制參數(shù)f由0逐漸增加,橫截同宿點(diǎn)橫截交叉發(fā)生于f≥3.6處。圖8~圖10表示了3種情況下時(shí)域波形、相軌跡、龐加萊截面和幅值譜。f=2對(duì)應(yīng)的點(diǎn)處于閾值曲線之下,沒有達(dá)到混沌運(yùn)動(dòng)的前提條件,如圖6。
圖8 f=2時(shí)1周期運(yùn)動(dòng)Fig. 8 1 period motion graph when f=2
圖9 f=3.2時(shí)2周期運(yùn)動(dòng)Fig. 9 2 period motion graph when f=3.2
圖10 f=4時(shí)混沌態(tài)Fig. 10 Chaotic state graph when f=4
由圖8可以看出,f=2時(shí),機(jī)電耦合扭轉(zhuǎn)系統(tǒng)為1周期運(yùn)動(dòng);當(dāng)f開始增大并超過臨界閾值,系統(tǒng)開始失去穩(wěn)定性,并存在發(fā)生混沌運(yùn)動(dòng)可能性。從圖9可以看出,f=3.2時(shí),此時(shí)運(yùn)行工況點(diǎn)處于閾值曲線之上,系統(tǒng)運(yùn)動(dòng)為2周期運(yùn)動(dòng),系統(tǒng)發(fā)生倍周期分岔,但并未發(fā)生混沌運(yùn)動(dòng)。從圖10可以看出:f=4時(shí),運(yùn)行工況點(diǎn)處于閾值曲線之上,進(jìn)入如圖7中的混沌區(qū)域,其計(jì)算結(jié)果如圖10;時(shí)域波形并無周期性,相軌跡為多個(gè)交叉的環(huán)形圈,龐加萊截面圖具有雜亂無章的多個(gè)點(diǎn),頻譜分析結(jié)果具有能量集中區(qū)域,這種現(xiàn)象表明該情況下系統(tǒng)發(fā)生混沌運(yùn)動(dòng)[9-10],這與我們之前Melnikov方法解析得到的結(jié)果是一致的。
以永磁同步電機(jī)的機(jī)電復(fù)合傳動(dòng)系統(tǒng)為研究對(duì)象,考慮機(jī)電耦合關(guān)系,建立永磁同步電機(jī)的電磁轉(zhuǎn)矩模型,并通過非線性動(dòng)力學(xué)理論得到電磁參數(shù)作用下的機(jī)電耦合參數(shù)振動(dòng)界限以及規(guī)律。分析得到如下結(jié)論:
1)產(chǎn)生非線性振動(dòng)的主要原因是由于車用永磁同步電機(jī)扭轉(zhuǎn)角引起的電磁激勵(lì)。
2)考慮機(jī)電耦合作用,扭轉(zhuǎn)振動(dòng)系統(tǒng)的固有特性與電機(jī)的運(yùn)行狀態(tài)以及結(jié)構(gòu)參數(shù)相關(guān),且在某些情況下,固有頻率有較大的降幅。
3)考慮到車用電機(jī)采用弱磁控制,并且固有頻率受到內(nèi)功率因數(shù)角的影響,因此車用永磁同步電機(jī)需要合理控制內(nèi)功率因數(shù)角的范圍。
4)純電動(dòng)拖拉機(jī)機(jī)電耦合傳動(dòng)系統(tǒng)的控制參數(shù)激勵(lì)值f逐漸增大時(shí),系統(tǒng)具有周期、擬周期和混沌等復(fù)雜的動(dòng)力學(xué)行為。當(dāng)出現(xiàn)非線性振動(dòng)時(shí),機(jī)電耦合傳動(dòng)系統(tǒng)受到非周期沖擊干擾力作用,產(chǎn)生較大強(qiáng)度動(dòng)載荷并引起非線性扭振,嚴(yán)重影響純電動(dòng)拖拉機(jī)機(jī)電耦合傳動(dòng)系統(tǒng)可靠性和耐久性。因此,應(yīng)避免純電動(dòng)拖拉機(jī)在低速重載工況下長時(shí)間運(yùn)行。