呂財(cái),王廣軍,2,陳紅,2,章廣祥,陳澤弘

(1.重慶大學(xué)動(dòng)力工程學(xué)院,400044,重慶;2.重慶大學(xué)低品位能源利用技術(shù)及系統(tǒng)教育部重點(diǎn)實(shí)驗(yàn)室,400044,重慶)

孔鉆削是機(jī)械加工過程常見的金屬去除成型工藝[1-2],與其他切削工藝相比,鉆削是半封閉狀態(tài)加工,摩擦嚴(yán)重、散熱困難。在鉆孔過程中,工件的溫度分布是一個(gè)重要的參數(shù),它不僅會(huì)影響鉆具的壽命,還直接影響到工件的性能指標(biāo),如殘余應(yīng)力、尺寸誤差和加工表面的硬度等。因此,全面了解孔鉆削過程中工件的溫度分布,有助于提高加工效率和產(chǎn)品工件的加工質(zhì)量[3-5]。

通常,實(shí)驗(yàn)測(cè)量只能獲得部分的溫度信息,為了獲得鉆孔工件完整的瞬態(tài)溫度場(chǎng),必須明確進(jìn)入工件的瞬時(shí)熱負(fù)荷。由于進(jìn)入工件的瞬時(shí)熱負(fù)荷難以直接準(zhǔn)確測(cè)量,借助導(dǎo)熱反問題方法,利用工件內(nèi)部或表面的部分溫度測(cè)量信息反求進(jìn)入工件的瞬時(shí)熱負(fù)荷,進(jìn)而重構(gòu)工件的瞬態(tài)溫度場(chǎng),是一種方便有效的解決方案。

目前,求解導(dǎo)熱反問題的方法主要有梯度類優(yōu)化方法[6-7]、正則化技術(shù)[8]、分散式模糊推理[9-10]、順序函數(shù)法(SFSM)[11-12]、隨機(jī)優(yōu)化方法[13-14]等。其中,SFSM利用未來一段時(shí)間內(nèi)的測(cè)量溫度信息反演當(dāng)前時(shí)刻的未知參數(shù),在非穩(wěn)態(tài)導(dǎo)熱反問題中得到了廣泛的應(yīng)用。然而,它需要假定未來一段時(shí)間內(nèi)未知變量的變化形式,這會(huì)影響反演結(jié)果的準(zhǔn)確性[15]。針對(duì)這一問題,Wang等結(jié)合系統(tǒng)預(yù)測(cè)控制思想,建立了一種基于動(dòng)態(tài)矩陣控制(DMC)的反演方法,實(shí)現(xiàn)了非穩(wěn)態(tài)導(dǎo)熱系統(tǒng)多個(gè)熱邊界條件的同時(shí)反演[16-17]。與SFSM相比,DMC反演方案無需假定待反演參數(shù)在未來時(shí)間段的變化形式,降低了反演結(jié)果對(duì)于未來觀測(cè)信息的依賴性,提高了反演結(jié)果的可靠性。作為一種順序反演方案,DMC為鉆孔工件的溫度場(chǎng)在線重構(gòu)提供了可能。

關(guān)于導(dǎo)熱反問題方法在鉆孔傳熱過程中的應(yīng)用研究已經(jīng)積累了一些成果。Huang等采用最速下降法反求了進(jìn)入鉆頭的熱流量[18]。Tai等基于序列二次優(yōu)化方法研究了深孔鉆削時(shí)鉆孔壁面和底面產(chǎn)生的熱量,發(fā)現(xiàn)深孔鉆削時(shí)預(yù)測(cè)工件溫度需要考慮孔壁面熱量的影響[19]。De Sousa等建立了一種基于格林函數(shù)和動(dòng)態(tài)觀測(cè)器的反演方法,并將其用于估算鉆孔過程中進(jìn)入工件的瞬時(shí)熱流量[20]。需要說明的是,在工件的鉆孔加工過程中,隨著鉆具的持續(xù)運(yùn)動(dòng),工件的受熱邊界不斷發(fā)生移動(dòng),系統(tǒng)輸入(瞬態(tài)熱負(fù)荷)與系統(tǒng)輸出(測(cè)點(diǎn)處溫度)之間的動(dòng)力學(xué)映射關(guān)系可能發(fā)生明顯變化。本質(zhì)上,鉆孔傳熱系統(tǒng)是一類非線性傳熱系統(tǒng),相對(duì)于線性系統(tǒng),非線性傳熱系統(tǒng)反問題的不適定性更加突出[7]。

本文基于模型預(yù)測(cè)控制思想研究了鉆孔工件移動(dòng)邊界熱流反演及溫度場(chǎng)重構(gòu)問題,以鉆孔工件傳熱模型為基礎(chǔ),建立了工件瞬態(tài)溫度場(chǎng)階躍響應(yīng)函數(shù)模型以及工件溫度預(yù)測(cè)模型;利用工件的部分溫度測(cè)量信息,通過滾動(dòng)優(yōu)化實(shí)現(xiàn)工件移動(dòng)傳熱邊界條件的反演及工件瞬態(tài)溫度場(chǎng)的在線重構(gòu)。文中利用數(shù)值仿真試驗(yàn)驗(yàn)證了上述方法的有效性,并詳細(xì)討論了鉆具進(jìn)給速度、溫度測(cè)量誤差、測(cè)點(diǎn)數(shù)等對(duì)反演結(jié)果的影響。

1 鉆孔工件溫度分布數(shù)學(xué)模型

鉆孔工件的簡(jiǎn)化物理模型如圖1所示。工件在加工過程中,一方面受到來自鉆削產(chǎn)生的部分摩擦熱q(t),另一方面與環(huán)境進(jìn)行對(duì)流換熱。對(duì)工件的傳熱過程引入如下假設(shè)[5,18,20]:

①由于工件的對(duì)稱性,只考慮工件沿徑向和軸向的傳熱,將其簡(jiǎn)化為二維傳熱問題;

②工件的對(duì)稱面視為絕熱,其他邊界則與環(huán)境發(fā)生對(duì)流換熱;

③進(jìn)入工件的熱負(fù)荷沿孔徑方向均勻分布,只隨時(shí)間發(fā)生變化;

④工件材料的熱物性各向同性且恒定。

根據(jù)上述假設(shè),建立圓柱坐標(biāo)下鉆孔工件溫度分布數(shù)學(xué)模型

圖1 鉆孔工件的簡(jiǎn)化物理模型

(1)

T(r,z,0)=T0(r,z)

(2)

(3)

(4)

(5)

式中:λ、ρ和c分別為工件的導(dǎo)熱系數(shù)、密度和比熱容;T0(r,z)為工件的初始溫度;q(t)為進(jìn)入到工件的瞬時(shí)熱流密度;T∞(t)為環(huán)境溫度;h為對(duì)流換熱系數(shù);n為界面的法線方向。

2 模型預(yù)測(cè)控制方法求解反問題

2.1 反問題目標(biāo)函數(shù)

如果方程組(1)～(5)中的其他條件均已知,移動(dòng)邊界熱流q(t)未知時(shí),根據(jù)溫度測(cè)點(diǎn)m(m=1,2,…,M)在時(shí)刻k,k+1,…,k+R-1的測(cè)量值來估計(jì)qk,構(gòu)成了對(duì)應(yīng)的導(dǎo)熱反問題。這里的R為預(yù)定的未來時(shí)間步數(shù),qk為q(t)在k時(shí)刻的估計(jì)值。

上述的導(dǎo)熱反問題可以通過使以下目標(biāo)函數(shù)最小化求解

minJ(Q)=(Y-T)T(Y-T)+αQTQ

(6)

2.2 預(yù)測(cè)模型

(7)

(8)

將上述方程轉(zhuǎn)化為矩陣方程

(9)

定義工件傳熱模型的階躍響應(yīng)函數(shù)

X(r,z,t)=?T(r,z,t)/?qk

(10)

X(r,z,t)反映了qk對(duì)溫度場(chǎng)T(r,z,t)的影響程度。將方程組(1)～(5)對(duì)qk求導(dǎo),則在[tk,tk+R-1]時(shí)間段內(nèi)X(r,z,t)應(yīng)滿足下列的階躍響應(yīng)方程組

(11)

X(r,z,tk-1)=0

(12)

(13)

(14)

(15)

由方程組(11)～(15)可以看出,X(r,z,t)隨著邊界AB的移動(dòng)而發(fā)生變化。因此,可以根據(jù)有限容積法和高斯賽德爾迭代法在每一時(shí)刻求解上述的階躍響應(yīng)方程組,得到階躍響應(yīng)函數(shù)的離散值,進(jìn)而按下式來確定測(cè)點(diǎn)m處的階躍響應(yīng)系數(shù)

(16)

2.3 滾動(dòng)優(yōu)化

模型預(yù)測(cè)反演方法通過采用滾動(dòng)優(yōu)化獲得邊界熱流向量Q=[qk,qk+1,…,qk+R-1]T。

將式(6)對(duì)Q求導(dǎo),并令dJ(Q)/dQ=0,得到移動(dòng)邊界熱流向量的最優(yōu)估計(jì)值

(17)

式中E為R階的單位矩陣。

在當(dāng)前時(shí)刻k,只需根據(jù)式(17)確定即時(shí)熱流量qk,下一時(shí)刻按同樣的優(yōu)化方法求得qk+1,此即滾動(dòng)優(yōu)化策略。

由式(17)得即時(shí)熱流量的最優(yōu)估計(jì)值為

(18)

d=[1,0,…,0]1×R

3 仿真試驗(yàn)結(jié)果及討論

3.1 仿真試驗(yàn)條件

對(duì)于圖1所示的鉆孔工件傳熱系統(tǒng),其結(jié)構(gòu)和物理參數(shù)為:外半徑r2=2 cm,內(nèi)半徑r1=0.5 cm,高度H=5 cm;工件的材料為AISI 1018碳鋼,導(dǎo)熱系數(shù)λ=51.9 W·m-1·K-1,密度ρ=7 870 kg·m-3,比熱容c=486 J·kg-1·K-1。

取工件初始溫度T0(r,z)=20 ℃,對(duì)流換熱系數(shù)h=30 W·m-2·K-1,環(huán)境溫度T∞(t)=20 ℃。

按式(19)和式(20)給定兩組不同變化形式的熱流密度

q1(t)=6+20sin(πt/24)e-0.2t(MW·m-2)

(19)

(20)

根據(jù)網(wǎng)格無關(guān)性驗(yàn)證結(jié)果,在區(qū)域離散時(shí),取空間步長(zhǎng)Δr=Δz=0.5 mm,時(shí)間步長(zhǎng)Δt=0.5 s。

仿真試驗(yàn)過程中,假定溫度測(cè)點(diǎn)均勻布置于工件的外表面,具體位置見表1。

表1 溫度測(cè)點(diǎn)位置

考慮到實(shí)際溫度測(cè)量結(jié)果不可避免地存在測(cè)量誤差,在數(shù)值仿真試驗(yàn)中,按下式模擬測(cè)點(diǎn)處的實(shí)際測(cè)量溫度

(21)

為了評(píng)估反演方案的準(zhǔn)確性,引入均方根誤差Sq及平均相對(duì)誤差η

(22)

(23)

式中:K為待反演熱流總的時(shí)間離散點(diǎn)數(shù);(qk)exa為熱流密度q(t)在k時(shí)刻的真實(shí)值。

3.2 移動(dòng)熱邊界條件反演

取未來時(shí)間步數(shù)R=3,鉆具進(jìn)給速度v=0.002 m/s,測(cè)量誤差的標(biāo)準(zhǔn)差σ=0,利用測(cè)點(diǎn)T2和T4提供的溫度信息反演未知熱邊界條件q(t)。圖2和圖3分別給出了邊界熱流q1(t)和q2(t)的反演結(jié)果。

圖2 熱流q1(t)的反演結(jié)果

圖3 熱流q2(t)的反演結(jié)果

由圖2和圖3可知,對(duì)于不同變化形式的移動(dòng)邊界熱流q(t),本文方法均可以獲得良好的反演結(jié)果,證明了本文反演方案的有效性。

3.3 進(jìn)給速度的影響

采用測(cè)點(diǎn)T2和T4的溫度信息,在R=3、σ=0的條件下,分別取v=0.001,0.002,0.003 m/s來考察鉆具進(jìn)給速度對(duì)反演結(jié)果的影響,如圖4和圖5所示,對(duì)應(yīng)的反演誤差見表2。

圖4 不同進(jìn)給速度下q1(t)的反演結(jié)果

圖5 不同進(jìn)給速度下q2(t)的反演結(jié)果

v/m·s-1η/%q1q2Sq/W·m-2q1q20.0010.010.021.04×1031.77×1030.0020.380.503.14×1044.08×1040.0033.236.232.92×1055.35×105

由圖4、圖5及表2可知,隨著鉆具進(jìn)給速度的增大,導(dǎo)致傳熱過程非線性特征進(jìn)一步加強(qiáng),當(dāng)R一定時(shí),對(duì)應(yīng)的熱邊界條件反演精度有所降低,但在仿真試驗(yàn)所討論的進(jìn)給速度范圍內(nèi),仍然可以獲得較好的反演結(jié)果。

3.4 測(cè)量誤差的影響

采用測(cè)點(diǎn)T2和T4的溫度信息,在R=3、v=0.002 m/s的條件下,分別取σ=0.01,0.03,0.05來考察溫度測(cè)量誤差對(duì)反演結(jié)果的影響,如圖6和圖7所示,對(duì)應(yīng)的反演誤差見表3。

圖6 不同測(cè)量誤差下q1(t)的反演結(jié)果

圖7 不同測(cè)量誤差下q2(t)的反演結(jié)果

ση/%q1q2Sq/W·m-2q1q20.012.275.591.86×1054.58×1050.036.008.034.92×1056.59×1050.059.9712.948.18×1051.06×106

根據(jù)σ=0.05條件下的反演結(jié)果q1(t)和工件的傳熱模型,重構(gòu)工件的瞬態(tài)溫度場(chǎng)如圖8所示。作為對(duì)比,圖9給出了空間點(diǎn)1(0.5 cm,4.5 cm)、2(0.5 cm,3.5 cm)、3(0.5 cm,2.5 cm)、4(0.5 cm,1.5 cm)和5(0.5 cm,0.5 cm)處工件溫度的重構(gòu)值和根據(jù)實(shí)際熱邊界條件得到的模擬結(jié)果。

(a)t=2 s(b)t=7 s(c)t=12 s(d)t=17 s(e)t=22 s 圖8 工件瞬態(tài)溫度場(chǎng)的重構(gòu)結(jié)果

圖9 σ=0.05時(shí)工件內(nèi)部溫度重構(gòu)結(jié)果

由以上結(jié)果可知,雖然隨著溫度測(cè)量誤差的增大,移動(dòng)邊界熱流的反演精度有所降低,但仍然可以獲得較好的反演結(jié)果,表明本文的模型預(yù)測(cè)反演方案具有一定的抗干擾能力。

3.5 測(cè)點(diǎn)數(shù)的影響

在R=3、σ=0.01、v=0.002 m/s的條件下,分別取測(cè)點(diǎn)數(shù)M=1,2,5進(jìn)行數(shù)值試驗(yàn),考察測(cè)點(diǎn)數(shù)對(duì)熱邊界條件反演結(jié)果的影響。其中,當(dāng)M=1時(shí),溫度測(cè)點(diǎn)取T3;當(dāng)M=2時(shí),溫度測(cè)點(diǎn)取T2和T4;當(dāng)M=5時(shí),溫度測(cè)點(diǎn)取T1～T5。對(duì)應(yīng)的結(jié)果見圖10、圖11和表4。

圖10 不同測(cè)點(diǎn)數(shù)下q1(t)的反演結(jié)果

圖11 不同測(cè)點(diǎn)數(shù)下q2(t)的反演結(jié)果

當(dāng)M=1時(shí),測(cè)點(diǎn)位置離初始和結(jié)束時(shí)的熱邊界較遠(yuǎn),對(duì)應(yīng)的階躍響應(yīng)系數(shù)較小,使得在鉆削起始段和結(jié)束段的熱流反演結(jié)果較差。在鉆削方向上增加溫度測(cè)點(diǎn)數(shù),可以明顯改善熱流的反演結(jié)果。對(duì)于本文工件的鉆削深度而言,M=2,5時(shí)的反演結(jié)果相差不大。

表4 不同測(cè)點(diǎn)數(shù)下的反演誤差

3.6 未來時(shí)間步數(shù)的影響

采用測(cè)點(diǎn)T2和T4的溫度信息,在σ=0.01、v=0.002 m/s的條件下,分別取R=3,6,9來討論未來時(shí)間步數(shù)的選取對(duì)反演結(jié)果的影響,對(duì)應(yīng)的結(jié)果見圖12、圖13及表5。

圖12 不同未來時(shí)間步數(shù)下q1(t)的反演結(jié)果

圖13 不同未來時(shí)間步數(shù)下q2(t)的反演結(jié)果

Rη/%q1q2Sq/W·m-2q1q232.275.591.86×1054.58×10563.916.073.26×1055.03×10592.796.032.37×1055.05×105

由圖12、圖13及表5可知,采用本文的模型預(yù)測(cè)反演方法,在不同的未來時(shí)間步數(shù)下均能夠獲得較為滿意的結(jié)果,降低了移動(dòng)邊界熱流反演結(jié)果對(duì)于未來時(shí)間步數(shù)的依賴性。

4 結(jié) 論

本文采用模型預(yù)測(cè)控制方法研究了鉆孔工件移動(dòng)邊界熱流的反演以及溫度場(chǎng)重構(gòu)的問題,并討論了鉆具進(jìn)給速度、溫度測(cè)量誤差、測(cè)點(diǎn)數(shù)等對(duì)邊界熱流反演結(jié)果的影響,得到的主要結(jié)論如下。

(1)本文建立的反演方案能夠有效反演鉆孔工件的移動(dòng)熱邊界條件,為鉆孔工件瞬態(tài)溫度場(chǎng)的在線重構(gòu)提供了一種解決方案。

(2)當(dāng)測(cè)量溫度數(shù)據(jù)存在一定的誤差時(shí),本文的模型預(yù)測(cè)控制方案仍然可以獲得較為滿意的反演結(jié)果,表明該方案具有一定的抗干擾能力。

(3)由于不需要對(duì)未來一段時(shí)間內(nèi)移動(dòng)熱邊界條件的變化進(jìn)行假設(shè),使得本文的模型預(yù)測(cè)反演方案降低了反演結(jié)果對(duì)于未來時(shí)間步數(shù)的依賴性。