基于因子分析和貝葉斯估計(jì)的機(jī)床熱誤差建模*

代貴松，楊建國，朱小龍

(上海交通大學(xué) 機(jī)械與動(dòng)力工程學(xué)院，上海 200240)

0 引言

隨著科學(xué)技術(shù)的日新月異，制造業(yè)對于機(jī)床精度的要求不斷提高［1］。大量研究表明，在高精、高速加工過程中，機(jī)床熱變形引起的熱誤差占總誤差的70%［2］。因此，機(jī)床熱誤差補(bǔ)償能夠很有效地提高機(jī)床的加工精度。在熱誤差補(bǔ)償過程中，溫度測點(diǎn)的布置、選擇和熱誤差模型的建立是難點(diǎn)［3］。從經(jīng)濟(jì)上考慮，需要盡可能地減少傳感器的數(shù)目，但太少的傳感器勢必影響熱誤差模型的精度［4］。因此，對溫度測點(diǎn)進(jìn)行優(yōu)化布置，確定傳感器的最優(yōu)數(shù)目和最佳位置具有重要的實(shí)用價(jià)值。

本文采用因子分析對機(jī)床的溫度變量進(jìn)行分組［5］，再從各分組中選擇出與熱變形相關(guān)系數(shù)最大的溫度變量作為建模的典型變量，以減小溫度變量之間的耦合性對于建模精度的影響，從而提高熱模型的魯棒性［6］。最后，根據(jù)優(yōu)選出的典型溫度變量，采用因子分析和貝葉斯估計(jì)建立機(jī)床熱誤差模型。結(jié)果表明，該方法大大提高了熱誤差模型的精度和魯棒性，實(shí)現(xiàn)了熱誤差建模領(lǐng)域一次新的嘗試，具有重大的理論意義和工程實(shí)踐意義。

1 溫度測點(diǎn)在機(jī)床上的優(yōu)化布置

1.1 因子分析基本原理

因子分析是一種多變量統(tǒng)計(jì)分析方法［7］，目的在于研究原始變量的內(nèi)部關(guān)系，通過尋找變量的共同因子來簡化和分析變量中存在的復(fù)雜關(guān)系。它把每個(gè)變量分解為兩部分因子:一部分是由這些變量內(nèi)含的公共因子構(gòu)成，即所謂公共因子部分;另一部分是每個(gè)變量各自獨(dú)有的因子，即所謂獨(dú)特因子部分。從本質(zhì)上來說，因子分析根據(jù)相關(guān)性的大小將變量分組，使得同組內(nèi)的變量之間相關(guān)性較高，但不同組的變量相關(guān)性較低。由此可以認(rèn)為，每組變量代表一個(gè)因子，它們可以反映問題的一個(gè)方面。

設(shè)F = (fij)= (f(1)，…，f(k))為一個(gè)n × k 階矩陣，U=(u(1)，…，u(p))為一個(gè)n×p 階矩陣，Λ=(λij)為k ×p 階矩陣，且滿足X =FΛ +U，寫成向量表達(dá)式，有:

上式即為因子分析的數(shù)學(xué)模型，式中:k ＜p，f(1)，f(2)，…，f(k)為公共因子，u(1)，u(2)，…，u(p)為獨(dú)特因子，λij表示第i 個(gè)變量x(i)在第j 個(gè)因子f(j)上的載荷。

1.2 利用因子分析選擇關(guān)鍵溫度變量

為了檢測機(jī)床的溫度場，按圖1 所示，在一臺(tái)車削中心上安裝了16 個(gè)溫度傳感器，可分為五組:①兩個(gè)傳感器(0 和1)用于測量主軸箱的溫度;②四個(gè)傳感器(2 ～5)用于測量絲桿螺母的溫度;③兩個(gè)傳感器(6 和7)用于測量冷卻液溫度;④一個(gè)傳感器(8)用于測量室溫;⑤七個(gè)傳感器(9 ～15)用于測量床身溫度。

圖1 溫度傳感器的布置

考慮到本文的篇幅，以主軸徑向?yàn)槔?，采用?shù)字式位移傳感器測量主軸X 方向熱漂移誤差，位移傳感器固定在刀架上，其布置如圖2 所示。

圖2 位移傳感器在主軸徑向的布置

定義編號(hào)為i 的溫度傳感器所測得的溫度變量為Ti，位移傳感器測得的位移誤差為δ。在檢測過程中，每隔450s 采集一次實(shí)驗(yàn)數(shù)據(jù)，共測得71 組數(shù)據(jù)，根據(jù)實(shí)驗(yàn)結(jié)果得到的機(jī)床主軸徑向熱誤差和溫度變化曲線如圖3、圖4 所示。

圖3 機(jī)床主軸徑向熱誤差曲線

圖4 機(jī)床溫度曲線

對16 個(gè)溫度變量進(jìn)行因子分析:第一步計(jì)算出16 個(gè)溫度變量的相關(guān)系數(shù)矩陣R;第二步利用雅克比方法計(jì)算出矩陣R 的特征向量和特征值;第三步計(jì)算各因子對總變差的貢獻(xiàn)率和累計(jì)貢獻(xiàn)率;第四步計(jì)算溫度變量在各主因子上的載荷。根據(jù)上述步驟可獲得前四個(gè)因子的特征值、貢獻(xiàn)率、累積貢獻(xiàn)率，如表1 所示。前四個(gè)因子的累積貢獻(xiàn)率為96.701%，故可采用4 因子模型。最后，計(jì)算各溫度變量在這四個(gè)因子上的載荷，如表2 所示。

表1 總變差分解表

根據(jù)表2，將溫度變量分為四組:選取在因子1上載荷較大的溫度變量T3、T4、T6、T7、T10、T13、T15為第一組;選取在因子2 上載荷較大的溫度變量T4、T6、T7、T9為第二組;選取在因子3 上載荷較大的溫度變量T10、T15為第三組;選取在因子4 上載荷較大的溫度變量T0、T1、T8為第四組。再根據(jù)各組溫度變量與機(jī)床主軸熱誤差的相關(guān)系數(shù)大小，選擇出絲桿螺母溫度T4、冷卻液溫度T6、床身溫度T15、主軸箱溫度T1四個(gè)關(guān)鍵溫度變量，該結(jié)果與文獻(xiàn)［8］中通過灰色理論得到的溫度關(guān)鍵點(diǎn)一致。

表2 載荷分布表

2 熱誤差建模

2.1 熱誤差模型推導(dǎo)

將機(jī)床徑向熱誤差δ 與優(yōu)選的四個(gè)關(guān)鍵溫度變量T4、T6、T15、T1組成模型輸入數(shù)據(jù)集X =［x0，x1，x2，x3，x4］=［δ，T4，T6，T15，T1］，X 是n ×p 維矩陣。其中n=71，p=5，δ 為71 組熱誤差組成的列向量，Ti(i=4，6，15，1)為71 組溫度值組成的列向量，由此可推導(dǎo)出機(jī)床的熱誤差模型。

首先，建立X 的因子分析模型，如下式:

式中，F(xiàn) 是X 的n×m 維主因子陣，Λ 是X 的m×p 維因子載荷陣，U 是X 的n×p 維特殊因子陣，m ＜5。

其次，對式(2)進(jìn)行代數(shù)變換，兩邊同乘以C =Λ'(ΛΛ')－1可得:

令X0=x0，X1=［x1，x2，x3，x4］，則X =［X0，X1］;令C0為C 的第一行組成的行矩陣，C1為C 其余行構(gòu)成的矩陣，則C=［C0，C1］T。

由此，式(3)可變?yōu)?

再次，對式(4)進(jìn)行代數(shù)變換，兩邊同乘以C0'

(C0C0')－1可得:即得到最終的機(jī)床熱誤差模型如下式:

在上式中，F(xiàn)、U、C 均是未知參數(shù)，本文采用貝葉斯估計(jì)的方法予以求出。

2.2 模型參數(shù)的貝葉斯估計(jì)

在因子分析模型(1)中，u(i)(i =1，2，…p)相互獨(dú)立，且服從多元正態(tài)分布N(0，Ψ)，Ψ 為p×p 維協(xié)方差陣，并由此可推知:

上式是因子分析模型的極大似然函數(shù)，下面討論參數(shù)的貝葉斯估計(jì)問題［9］。

模型參數(shù)先驗(yàn)分布的設(shè)置是貝葉斯估計(jì)的基礎(chǔ)，對于因子分析模型參數(shù)，Λ 的先驗(yàn)分布為矩陣正態(tài)分布，Ψ 的先驗(yàn)分布為Wishart 分布，F(xiàn) 的先驗(yàn)分布為矩陣正態(tài)分布，假設(shè)參數(shù)(Λ，Ψ)與F 相互獨(dú)立，即:

其中，

式中:B 是對稱正定陣，H ＞0 為對稱正定陣，(v，Q，B，Λ0)為超參數(shù)，具體設(shè)置方法見文獻(xiàn)［10］。

模型參數(shù)后驗(yàn)估計(jì)采用最大后驗(yàn)估計(jì)(MAP)，是通過模型參數(shù)的后驗(yàn)分布進(jìn)行迭代得到的，因此要先計(jì)算模型參數(shù)的最大后驗(yàn)分布。

根據(jù)貝葉斯定理可推導(dǎo)出因子載荷陣Λ 的最大后驗(yàn)分布為:

在式(12)中，～Λ = (X'F + Λ0D－1)(D－1+ F'F)－1，顯然當(dāng)Λ =～Λ 時(shí)，π(Λ|F，Ψ，X)取最大，故Λ的最大后驗(yàn)估計(jì)是:

同理可得，Ψ 的最大后驗(yàn)估計(jì)是:

F 的最大后驗(yàn)估計(jì)是:

根據(jù)式(13)、(14)、(15)，利用迭代法可以得到用于熱誤差模型的的參數(shù)估計(jì)，即

在式(16)、(17)、(18)中，k 的取值范圍是0，1，…n，整個(gè)迭代過程持續(xù)到(Λk+1，Ψk+1，F(xiàn)k+1)收斂為止。

綜合式(2)、(16)可得:

2.3 模型驗(yàn)證

機(jī)床主軸徑向熱誤差δ 與四個(gè)關(guān)鍵溫度變量T4、T6、T15、T1組成的數(shù)據(jù)集為模型輸入，即2.1 所推導(dǎo)出的熱誤差預(yù)測模型，為模型1;機(jī)床徑向熱誤差δ 與全溫度變量T0、T1、…、T15組成的數(shù)據(jù)集為模型輸入，按照2.1 中同樣的方法，推導(dǎo)出類似的熱誤差預(yù)測模型，為模型2;以T4、T6、T15、T1為自變量，δ為因變量，采用最小二乘法得到機(jī)床的熱誤差模型，為模型3。模型3 的表達(dá)式為:

圖5 為三個(gè)熱誤差預(yù)測模型對原始測量熱誤差的擬合情況。將模型預(yù)測熱誤差與原始測量誤差做差值運(yùn)算，可得到三個(gè)熱誤差預(yù)測模型的殘差，如圖6 所示。通過圖5、圖6 不難看出:對于原始測量熱誤差的擬合情況:模型1 比模型3 好，模型3 比模型2好;殘差變化范圍:模型1 比模型3 小，模型3 比模型2 小。綜合比較，模型1 對熱誤差的預(yù)測結(jié)果最好。

圖5 預(yù)測模型的擬合情況

圖6 預(yù)測模型的殘差

通過因子分析和貝葉斯估計(jì)的建模補(bǔ)償，使機(jī)床的徑向熱誤差從22μm 減小到了4μm 以內(nèi)，補(bǔ)償效果達(dá)到了82%，大大提高了機(jī)床的加工精度。

3 結(jié)論

因子分析是較好的溫度測點(diǎn)優(yōu)化策略，極大地減少了傳感器的數(shù)量，有效地避免了熱誤差模型中的耦合現(xiàn)象，提高了熱誤差的辨識(shí)精度。

通過因子分析模型推導(dǎo)出熱誤差預(yù)測模型，然后對模型中的參數(shù)進(jìn)行貝葉斯估計(jì)，實(shí)驗(yàn)驗(yàn)證了這一方法可用于熱誤差建模，模型精度比最小二乘法高，而且魯棒性也比最小二乘法好。

［1］W.S. Yun，S.K. Kim，D.W. Cho，Thermal error analys is for a CNC lathe feed drive system，International Journal of Machine Tools and Manufacture，1999(39):1087－1101.

［2］Zhao Haitao，Yang Jianguo，Shen Jinhua，Thermal error optimization modeling and real-time compensation on a CNC turning center. International Journal of Machine Tools and Manufacture，2007，207(2008):172－179.

［3］蘇鐵明，葉三排，孫偉. 基于FCM 聚類和RBF 神經(jīng)網(wǎng)絡(luò)的機(jī)床熱誤差補(bǔ)償建模［J］. 組合機(jī)床與自動(dòng)化加工技術(shù)，2011，2011(10):1－4.

［4］林偉青，傅建中. 數(shù)控機(jī)床熱誤差建模中的溫度傳感器優(yōu)化研究［J］. 成組技術(shù)與生產(chǎn)現(xiàn)代化，2007，24(3):5－8.

［5］嚴(yán)樂，劉紅忠，盧秉恒，等. 高精度壓印機(jī)熱誤差補(bǔ)償中溫度變量的辨識(shí)［J］. 西安交通大學(xué)學(xué)報(bào)，2006，40(7):827－830.

［6］楊建國，薛秉源. CNC 車削中心熱誤差模態(tài)分析及魯棒建模［J］. 中國機(jī)械工程，1998，9(5):31－35.

［7］張潤楚. 多元統(tǒng)計(jì)分析［M］. 北京:科學(xué)出版社，2003.

［8］李永祥，童恒超，楊建國. 灰色系統(tǒng)理論在機(jī)床熱誤差測點(diǎn)優(yōu)化中的應(yīng)用［J］. 機(jī)械設(shè)計(jì)與研究，2006，22(3):78－80.

［9］吳喜之，科茨(Kotz，Samuel). 現(xiàn)代貝葉斯統(tǒng)計(jì)學(xué)［M］.北京:中國統(tǒng)計(jì)出版社，2000.

［10］Liu Y-K and Liu B.Fuzzy Random Variables:A Scalar Expected Value Operator［J］. Fuzzy Optimization And Decision Making，2003，11(2):143－160.