劉成穎，張 潔，劉 巍，張 智,3

(1.清華大學 機械工程系，北京 100084；2.精密超精密制造裝備及控制北京市重點實驗室，北京 100084；3.海軍航空工程學院 飛行器工程系，山東 煙臺 264001)

基于步長控制的銑削穩(wěn)定性預測全離散法*

劉成穎1,2，張 潔1，2，劉 巍1,2，張 智1,2,3

(1.清華大學 機械工程系，北京 100084；2.精密超精密制造裝備及控制北京市重點實驗室，北京 100084；3.海軍航空工程學院 飛行器工程系，山東 煙臺 264001)

基于控制離散步長來保證預測精度的思想，提出了一種改進的用于預測銑削穩(wěn)定性葉瓣圖的時域全離散法。一直以來，銑削工藝系統(tǒng)的動力學過程被看作帶有單個延時反饋量的周期性線性系統(tǒng)，其數(shù)學模型是無限維時滯微分方程組。全離散法在一個周期內(nèi)將時滯微分方程離散為有限個常微分方程，求解前后兩步的狀態(tài)轉移矩陣，繼而得到一個完整周期內(nèi)的狀態(tài)轉移矩陣。在選定的轉速范圍內(nèi)循環(huán)求解過程中，通過控制離散步長小于設定值來達到提高低轉速條件下全離散法的收斂率的目的。通過仿真和切削實驗，證明方法具有較高的預測精度和計算效率。

銑削穩(wěn)定性；Floquet理論；步長控制；收斂率

0 引言

銑削加工過程中的顫振是一種嚴重的系統(tǒng)失穩(wěn)現(xiàn)象。顫振迫使工藝人員選取保守的切削參數(shù)，降低材料去除率，而且顫振會大大降低加工質量甚至使工件報廢。再生理論較好地解釋了顫振的機理并獲得了大量的實驗驗證[1]，該理論認為上一個刀齒切削時留下的表面波紋會影響當前刀齒的瞬時切削厚度，即切削時存在延時反饋，該反饋可能引起系統(tǒng)發(fā)散。工藝系統(tǒng)的動力學模型被歸結為包含周期系數(shù)矩陣的時滯微分方程組，由于方程中存在時滯項，其單值算子和相空間均為無限維，因此不存在閉形解[2]。

對于該方程組的近似求解方法，可歸納為頻域法和時域法兩種。Altintas等將方程中的周期系數(shù)矩陣展開成傅立葉級數(shù)，只保留其平均量，得到了可在頻域中解析計算穩(wěn)定性邊界的單頻率法[3]。單頻率法具有極高的計算效率，但在小徑向切深時不能預測出附加不穩(wěn)定區(qū)域。之后Altintas等考慮周期系數(shù)矩陣的前r階級數(shù)，推廣得到了多頻率法解決了這一問題[4]，但該方法計算效率較低，應用不多。Insperger等提出了半離散法[5]，該方法在時域上將時滯項離散，利用插值法將時滯微分方程化簡為一系列的常微分方程，首先獲得相鄰兩步之間的狀態(tài)轉移矩陣，繼而獲得一個完整時滯周期上的狀態(tài)轉移矩陣，并應用Floquet定理判定銑削過程的穩(wěn)定性。Ding等將狀態(tài)項全部離散，基于直接積分格式推導得到了計算效率更高的全離散法[6]。時域離散方法兼顧了計算精度和計算效率。使用穩(wěn)定性判斷方法，掃描選定的工藝參數(shù)范圍(通常是轉速和切深)，便可得到穩(wěn)定性葉瓣圖。穩(wěn)定性葉瓣圖直觀地給出了工藝系統(tǒng)穩(wěn)定區(qū)域和不穩(wěn)定區(qū)域的邊界，對于通過合理選取工藝參數(shù)來避免顫振并獲得較大的材料去除率具有極大的指導意義。

在掃描轉速區(qū)間時，模型的時滯項發(fā)生變化，在低轉速條件下時滯項要大于高轉速條件下的時滯項。雖然通過不同的插值方法和插值階次，發(fā)展出了多種具有更低階局部離散誤差的離散方法[7-9]，但通常的通過固定離散步數(shù)來進行穩(wěn)定性判斷的做法，可能存在低轉速下收斂率過低的問題。通過設定初始時間步長來確定每個轉速較為合理的離散步數(shù)，不僅可以提高低轉速區(qū)間的收斂率，而且可以提高整個計算流程的計算效率。

1 動力學模型及全離散法判穩(wěn)

不失一般性，考慮再生效應并略去切削力中的靜態(tài)成分，剛性工件柔性刀具雙自由度銑削工藝系統(tǒng)的動力學模型為:

(1)

(2)

其中，

設定初始時間步長τ0，則一個刀齒切削周期等距離散為m步，且：

(3)

則實際時間步長τ=T/m。

其中，ceil(·)表示向上取整，以控制實際時間步長不大于設定值τ0，mlim是人工設定的最小離散步數(shù)，以避免時滯量較小時因離散步數(shù)過小而導致預測精度不足。離散時間點為：

tk=t1+(k-1)τ,k=1,2,…,m+1

在第k步和第k+1步之間，即kτ≤t≤(k+1)τ時，狀態(tài)空間方程的解為：

(4)

將周期項A(kt+t-x)和B(kt+t-x)，狀態(tài)項x(kτ+τ-ξ)，時滯項x(kτ+τ-ξ-T)進行插值處理并積分可得xk+1的顯式表達式：

xk+1=[I-Fk+1]-1(F0+F0,k)xk+[I-Fk+1]-1Fm-1xk+1-m+
[I-Fk+1]-1Fmxk-m

(5)

其中，

構造離散映射yk+1=Dkyk

(6)

其中，

yk=col(xkxk-1…xk+1-mxk-m)

(7)

(8)

可構造一個完整周期上的狀態(tài)轉移矩陣

φ=Dm-1Dm-2…D1D0

(9)

最后可應用Floquet定理判斷系統(tǒng)的穩(wěn)定性，即若轉移矩陣的所有特征值的模均小于等于1，則系統(tǒng)穩(wěn)定，否則系統(tǒng)發(fā)散。

2 葉瓣圖仿真分析

為檢驗方法的精度和計算效率，選取文獻中的雙自由度算例進行葉瓣圖仿真分析：

(10)

其中，φj(t)是第j個刀齒的角位置；g是用來判斷刀齒切削狀態(tài)的窗函數(shù)，當?shù)洱X切削時g=1，刀齒退出切削時g=0；N是刀齒數(shù)目。刀具為兩齒刀，刀具x和y方向的模態(tài)對稱，模態(tài)質量mt=0.03993kg，阻尼比ζ=0.011，固有頻率fn=922Hz；切削力系數(shù)分別為Kt=6×108N/m2和Kr=4×108N/m2。刀具浸入比為1，順銑。

Insperger等[10]證明當離散步數(shù)m足夠大時，時域離散法獲取的穩(wěn)定性界限將收斂到真實值。因此，應用1階全離散法，取m=200繪制的葉瓣圖作為真實參考值，用于對比仿真曲線的精度。所有仿真曲線的葉瓣圖均將轉速-切深二維平面劃分為200×100的網(wǎng)格。仿真結果如圖1所示，點劃線表示準確的穩(wěn)定性邊界。結果表明，基于步長控制的方法在計算效率占優(yōu)的情況下依然具有更快的收斂速率，尤其是在低轉速即大時滯情況下對比更為明顯。所有運算均基于MATLAB?2016a平臺并在同一臺電腦上(intel core(TM) i7-6700HQ, 2.6 GHz, 8 GB)完成。

圖1 1階全離散法和控制步長全離散法仿真對比

3 實驗驗證

驗證試驗在VMP-23加工中心上展開，機床最高轉速10000rpm。實驗所使用的刀具為普通硬質合金兩齒立銑刀，直徑為12mm，安裝長度60mm，工件材料是7075-T6航空鋁合金。通過模態(tài)試驗和切削力系數(shù)辨識實驗獲取的模態(tài)參數(shù)及切削力系數(shù)如表1所示。刀具浸入比為1，逆銑。

表1 實驗刀具部分參數(shù)

使用步長控制方法得到的葉瓣圖如圖2所示，初始步長為1×10-6s。為驗證方法的有效性和葉瓣圖的準確性，在穩(wěn)定區(qū)域和顫振區(qū)域選取不同的轉速切深組合進行切削實驗。實驗過程中，使用Kistler 9257B測力儀和Kistler 5070A功率放大器來采集切削力，實驗過程中均沒有使用切削液。圖中，圓圈代表穩(wěn)定點，叉代表顫振點，圓圈和叉重疊代表不確定點。從實驗結果來看，葉瓣圖很好地預測出了工藝系統(tǒng)的穩(wěn)定區(qū)域。

圖2 實驗刀具穩(wěn)定性葉瓣圖

選取葉瓣圖中典型的穩(wěn)定點和顫振點進行對比。圖2中參數(shù)點A(7500rpm，0.5mm)和參數(shù)點B(8000rpm,0.5mm)的Y向切削力頻譜分別如圖3a和圖3b所示。工藝系統(tǒng)在A點發(fā)生顫振，顫振頻率為2172Hz，接近于系統(tǒng)固有頻率，而在B點穩(wěn)定，其頻譜的峰值均為刀齒通過頻率的分頻和倍頻。

4 結論

(1)通過步長控制，可避免在使用時域離散方法判定大時滯量時滯微分方程斂散性時因離散步數(shù)不足而導致的收斂率過低的問題，即可在保證計算效率的前提下提高葉瓣圖計算過程中低轉速區(qū)間的收斂率。

(2)穩(wěn)定性葉瓣圖可以準確預測工藝系統(tǒng)的穩(wěn)定區(qū)域。在兩個葉瓣之間存在穩(wěn)定切深較大的穩(wěn)定區(qū)域，通過合理選取穩(wěn)定區(qū)間內(nèi)的工藝參數(shù)，可以在穩(wěn)定工況下獲取較大的加工效率。

(3)在切削非薄壁件工況下，刀具往往是工藝系統(tǒng)最薄弱的環(huán)節(jié)，而不同的刀具以及相同刀具的不同安裝方式均具有不同的穩(wěn)定性葉瓣圖。怎樣快速獲取穩(wěn)定性葉瓣圖以及如何通過合理選取刀具和安裝方式來獲取更優(yōu)的穩(wěn)定區(qū)域需要進一步的研究。

圖3 穩(wěn)定點與顫振點切削力頻譜圖

[1] Quintana G, Ciurana J. Chatter in machining processes: a review[J]. International Journal of Machine Tools and Manufacture, 2011, 51(5): 363-376.

[2] 丁漢, 丁燁, 朱利民. 銑削過程穩(wěn)定性分析的時域法研究進展[J]. 科學通報, 2012, 57(31): 2922-2932.

[4] Merdol S D, Altintas Y. Multi Frequency Solution of Chatter Stability for Low Immersion Milling[J]. Journal of Manufacturing Science & Engineering, 2004, 126(3):459-466.

[5] Insperger T, Stépán G. Semi-discretization method for delayed systems[J]. International Journal for numerical methods in engineering, 2002, 55(5): 503-518.

[6] Ding Y, Zhu L M, Zhang X J, et al. A full-discretization method for prediction of milling stability[J]. International Journal of Machine Tools and Manufacture, 2010, 50(5): 502-509.

[7] Insperger T, Stépán G. Updated semi-discretization method for periodic delay-differential equations with discrete delay[J]. International Journal for Numerical Methods in Engineering, 2004, 61(1): 117-141.

[8] Ding Y, Zhu L M, Zhang X J, et al. Second-order full-discretization method for milling stability prediction[J]. International Journal of Machine Tools and Manufacture, 2010, 50(10): 926-932.

[9] Qiang Q, Sun Y, Jiang Y. On the accurate calculation of milling stability limits using third-order full-discretization method[J]. International Journal of Machine Tools & Manufacture, 2012, 62(1):61-66.

[10] Insperger T. Full-discretization and semi-discretization for milling stability prediction: Some comments[J]. International Journal of Machine Tools & Manufacture, 2010, 50(7):658-662.

AFull-discretizationMethodforMillingStabilityPredicationBasedonTime-stepControl

LIU Cheng-ying1,2, ZHANG Jie1,2, LIU Wei1,2, ZHANG Zhi1,2,3

(1. Department of Mechanical Engineering, Tsinghua University, Beijing 100084, China;2. Beijing Key Lab of precision/Ultra-precision Manufacturing Equipment and Control, Beijing 100084, China)

A modified full-discretization method is proposed for milling stability prediction based on time-step control. The dynamic process of milling system is regarded as a linear periodic system with a single time delay and its mathematical model is formulated as delay-differential equations(DDEs). By dividing the delay period equally into a set of time steps, the state transition matrixes between the adjacent two steps could be deduced by the full-discretization method. Then the state transition matrix over one delay period is structured. By controlling the time steps not exceeding the set value when scanning the objective range of rotating speed, the full-discretization method converges faster than setting a constant time interval value. At last the modified method is verified by simulation and experiments.

milling stability; Floquet theory; time-step control; rate of convergence

TH161；TG506

A

1001-2265(2017)12-0059-03

10.13462/j.cnki.mmtamt.2017.12.014

20147-03-15

國家04科技重大專項課題(2013ZX04001-021)

劉成穎(1960—)，女，遼寧大連人，清華大學副教授，研究方向為超精密加工工藝與裝備、精密驅動與控制、CAD/CAM/CAPP，(E-mail)liucy@tsinghua.edu.cn。

(編輯李秀敏)