侯東曉,徐 良,時培明
(1.東北大學秦皇島分校 控制工程學院,河北 秦皇島 066004; 2.燕山大學 電氣工程學院,河北 秦皇島 066004)
軋機振動問題一直是困擾鋼鐵行業的一個基礎性難題,振動的頻繁發生嚴重影響著產品質量和生產效率[1-3]。因此如何減少或抑制振動的發生,對穩定軋制生產有著重要的作用[4-5]。
針對長久存在的軋機振動問題,國內外專家學者對此進行了廣泛的研究。劉浩然等[6]建立了軋件變形滯后非線性作用下輥系垂直振動動力學模型。Shi等[7]考慮軋機主傳動系統接軸傾角和輥隙摩擦力的影響,建立了含參數激勵剛度和非線性摩擦阻尼的主傳動系統非線性振動模型。Amer等[8]采用多尺度攝動法求解了軋機的扭轉振動系統,運用數值方法研究了系統的穩定性和穩態響應。范小彬等[9]考慮軋制界面非線性黏滑摩擦特性, 建立了軋輥水平方向摩擦顫振模型和工作輥水平方向“跳振”模型。孫建亮等[10]基于軋制理論和動力學理論建立了軋機水平振動與軋制過程參數的關系模型。Zeng等[11]建立了非線性摩擦下軋機垂直-水平-扭轉耦合動力學模型,利用Hopf分岔代數準則對系統的穩定域進行了分析。Liu等[12]研究了軋機輥系水平-扭轉耦合振動的靜態和動態分岔特性,得到了系統的分岔和穩定性條件。
在軋制潤滑方面,付括等[13]建立了穩態和非穩態軋制變形區油膜厚度分布模型,并應用卡爾曼方程分析其應力分布;Wang等[14]建立了不穩定潤滑條件下的界面膜約束多重耦合模型,定量分析了一些參數對振動速度和振動位移的影響;Azushima等[15]提出了一種入口油膜厚度的估算方法;Noriki等[16]對軋件表面油膜和入口處乳化液的成膜機理的相互作用進行了研究。
本文考慮混合潤滑摩擦的影響,建立了與入口油膜厚度有關的動態軋制力模型。考慮軋機輥系結構振動的影響,建立了軋機垂直非線性振動模型,并得到相應的動力學方程。通過數值仿真,分析了以軋制力參數和外擾力為變量的主共振幅頻特性和分岔特性,為軋機振動問題的研究提供參考。
考慮軋輥在垂直方向上的振動,建立了軋機輥系幾何形狀圖,如圖1所示。

圖1 軋機輥系幾何形狀
圖1中:y0為帶鋼軋制前的厚度;y1為帶鋼軋制后的厚度;R為軋輥半徑;l為變形區長度;s為后張力;ve為軋件的入口速度;vd為帶鋼的出口速度;vr為軋輥的軋制速度;xe為帶鋼的入口位置;xd為帶鋼的出口位置。
假設輥縫形狀呈二次曲線分布[17],則變形區中任意位置的軋件厚度可以表示為
(1)
根據式(1)可以得到入口位置xe(yc)。
(2)
根據金屬秒流量相等的原則,在考慮軋輥垂直振動時,金屬秒流量公式應修正為
(3)
在變形區的中性點處,軋件的速度與軋輥的水平分速度相等。
(4)
當軋輥發生振動時,根據幾何關系可得
(5)
聯立式(1)和式(5),并代入式(3),得到方程
(6)
(7)
在軋件的變形區任意位置選取一微元體,如圖2所示。

圖2 變形區任意位置的應力關系示意圖
根據力的平衡關系建立平衡方程。
(σx+dσx)(y(x)+dy(x))-σx·dy(x)+
2pdxtanφ?2τsdx=0
(8)
式中:負號為從入口處到中性點處;正號為從中性點到出口處。根據幾何關系可知tanφ=dy(x)/2dx,代入式(8),化簡方程為
(9)
式中的屈服準則采用均勻變形下的von Mises屈服準則。
(10)
式中:σ為軋件材料的屈服強度;σx為正應力。
代入可得正應力的微分方程為
(11)
混合潤滑摩擦狀態下,軋件的摩擦切應力τs由干摩擦產生的摩擦切應力τa和全膜潤滑產生的摩擦切應力τb按比例組成。
τs=λτa+(1-λ)τb
(12)
式中:τa為干摩擦產生的摩擦切應力,由黏著摩擦理論式(13)求得;τb為全膜潤滑產生的摩擦切應力(由式(14)求得);λ為界面間的接觸比。
(13)
(14)
式中:η0為乳化液的黏度;z變形區內潤滑油膜的厚度,假設乳化液是不可壓縮的,變形區內的潤滑油膜厚度可由式(15)求得。
(15)
式中,z0為入口油膜厚度。在忽略表面形貌因素影響時,z0可由Wilson入口膜厚計算公式[18]求得。
(16)
式中:γ為乳化液的黏壓系數;η0為標準大氣壓下乳化液的黏度。
將式(12)、式(13)、式(14)代入正應力微分方程式(11),并對式(11)進行積分可分別得到入口處、中性點處和出口處的應力。
(17)
(18)
(19)
將入口位置、中性點和出口位置的應力表達式分別代入式(10),可以分別得到入口位置、中性點和出口位置的壓力pe,pn,pd。
假設軋制過程中軋件變形區的軋制壓力近似線性分布[19],采用如圖3所示的軋制力模型。根據代數關系式可以得到軋制力的表達式為

圖3 變形區軋制壓力分布示意圖
(20)
將壓力表達式和位置表達式代入其中并將軋制力進行Taylor展開,得到如下形式的動態軋制力。
F=F0+ΔF
(21)
式中:F0為穩態軋制力;ΔF為軋制力動態變化部分。
進一步考慮軋機機械結構影響,可建立了軋機輥系垂直振動模型,如圖4所示。

圖4 軋機垂直振動模型
考慮到支撐輥質量遠大于工作輥質量[20],可將其考慮為一個質量體,m1和m2分別為上下部輥系的等效質量;k1和c1為上輥系的等效剛度和等效阻尼;k2和c2為下輥系的等效剛度和阻尼;P1和P2分別為上下輥系受到的外擾力。
(22)
考慮到軋機機構具有對稱性,有m1=m2,k1=k2,c1=c2,F1=F2,P1=P2,式(22)可以寫成
(23)

(24)
(25)
當軋機入口板帶厚度或張力出現波動時,可導致周期性的外部擾動P=P0cos(Ωt)[21]。設式(25)是弱非線性系統,在式(25)右側冠以小參數ε。

(26)
采用多尺度法對非線性系統求解,令Tn=εnt。
設式(26)的解為
yc(t)=y0(T0,T1)+εy1(T0,T1)
(27)
將式(27)代入式(26),令ε冪指數相同項相等,得
(28)


(29)
式中,Di=?/?Ti。
設式(28)的解為
(30)
主共振時外擾力頻率Ω=ω0+εσ,將Ω表達式和式(30)代入式(29),消除方程的久期項,得
(31)
式中,σ為協調頻率。
幅值A引入極坐標形式
(32)
將幅值的極坐標形式式(32)代入式(31),分離實部、虛部,并設ψ=σT1-φ1,得
(33)
(34)
采用某廠1780軋機參數進行仿真。軋機和軋件的參數:m=1.44×105kg,k=1.08×1011N/m,c=2.2×107N·s/m,y0=0.014 1 m,y1=0.008 2 m,R=0.42 m,s=3.6 MPa,σ=150.6 MPa;乳化液參數:γ=2×10-8Pa-1,η0=0.42 Pa·s,λ=0.8,z0=1.8×10-6m。動態軋制力的系數:a1=1.39×1010,a2=1.26×108,a3=4.81×1012,a4=8.77×1010,a5=4.58×108,a6=1.22×1016,a7=3.39×1014,a8=3.18×1012,a9=9.97×109。
動態軋制力參數a6變化時系統的幅頻曲線,如圖5所示。從圖5可知,較大的三次剛度系數會使系統產生跳躍現象,但不會影響主共振的幅值。

圖5 參數a6變化時的主共振幅頻曲線
外部擾動力幅值P0變化時系統的幅頻響應曲線,如圖6所示。從圖6可知,外擾力幅值不僅會使系統產生跳躍現象,還會影響主共振幅值的大小。

圖6 外擾力幅值P0變化時的主共振幅頻曲線
外部擾動力幅值P0變化時系統在垂直方向的全局分岔曲線,如圖7所示。從圖7可知,外擾力幅值的變化會使系統在倍周期運動和混沌運動之間相互轉換。

圖7 外擾力幅值P0變化時的分岔圖
入口油膜厚度z0變化時系統在垂直方向的全局分岔曲線,如圖8所示。從圖8可知:當k=1.3×1011N/m時油膜厚度的變化會使系統的運動狀態更加復雜、不穩定;當k=1.08×1011N/m時,系統為2倍周期運動,油膜厚度的變化不會影響系統的運動狀態。比較圖8(a)和圖8(b)可知,合適的剛度會使系統對變化的油膜厚度不再敏感。

圖8 入口油膜厚度z0變化時的分岔圖
系統的相圖和Poincare截面圖,如圖9~圖11所示。

圖9 P0=1.6×106 N時垂直方向的周期運動

圖10 P0=4.1×106 N時垂直方向的周期運動

圖11 P0=6.96×106 N時垂直方向的混沌運動
在圖9中,相圖為2條封閉的曲線,Poincare截面為2個點,說明系統是二倍周期運動;在圖10中,相圖分別為4條封閉的曲線,Poincare截面分別為4個點,說明系統是4倍周期運動;在圖11的Poincare截面圖可看到許多雜亂無章的點,說明系統處在混沌運動狀態。
(1)通過將Karman微分方程中的干摩擦模型替換為混合潤滑摩擦模型,建立了混合潤滑摩擦狀態下軋機輥系垂直振動模型,并得到了與入口油膜厚度有關的動態軋制力表達式。
(2)動態軋制力參數a6增大,系統會出現跳躍現象;外擾幅值P0的變化同時影響著系統的主共振幅值大小和系統的跳躍性。因此,可以通過降低外部擾動力來減小軋機振動幅值,減小振動帶來的危害。
(3)改變外部擾動力幅值,系統會在周期運動和混沌運動之間相互轉換;當k=1.3×1011N/m時,入口油膜厚度的變化會影響系統的運動狀態;當k=1.08×1011N/m時,系統的運動狀態不受入口油膜厚度變化的影響,系統更加穩定。