基于代理模型的加工中心滑鞍有限元模型修正

陸原超，周俊榮，王瑞超，李會(huì)軍

（五邑大學(xué) 智能制造學(xué)部，廣東 江門 529020）

0 引言

伴隨著機(jī)械制造加工業(yè)的飛速發(fā)展，對(duì)加工中心的加工精度和加工效率等要求不斷上升，對(duì)加工中心的設(shè)計(jì)需求提出了新的挑戰(zhàn)，傳統(tǒng)的經(jīng)驗(yàn)設(shè)計(jì)方法已經(jīng)不能滿足新的設(shè)計(jì)需求，有限元分析越來(lái)越多地應(yīng)用到加工中心的設(shè)計(jì)過(guò)程之中[1]。有限元分析的第一步，就是建立一個(gè)準(zhǔn)確的有限元模型，但由于加工中心結(jié)構(gòu)的復(fù)雜性[2]，模型修正技術(shù)應(yīng)運(yùn)而生[3]。有限元計(jì)算結(jié)構(gòu)輸出與輸入為隱函數(shù)關(guān)系[4]，修正過(guò)程中頻繁調(diào)用有限元模型，需要大量的計(jì)算量，通過(guò)代理模型將有限元模型用顯函數(shù)近似可以大幅減少計(jì)算量，提高模型修正效率[5]。常用的代理模型有多項(xiàng)式響應(yīng)面[6]、支持向量回歸[7]和Kriging模型[8]等。

1 理論基礎(chǔ)

1.1 模型修正理論

有限元模型修正起源于20世紀(jì)60年代，可以分為矩陣法和參數(shù)法。矩陣法通過(guò)修改結(jié)構(gòu)剛度矩陣與質(zhì)量矩陣等減小與真實(shí)結(jié)構(gòu)間的誤差，但由于矩陣法缺少物理意義且不具備工程意義，已經(jīng)逐漸被淘汰。參數(shù)法通過(guò)修改密度、彈性模量、連接剛度等物理參數(shù)，保證了修正過(guò)程的物理意義，逐漸成為研究重點(diǎn)[2]。

通過(guò)修正有限元模型的設(shè)計(jì)參數(shù)，參數(shù)型模型修正方法逐步減小響應(yīng)計(jì)算值與實(shí)測(cè)值之間的誤差，將模型修正這一參數(shù)辨識(shí)問(wèn)題轉(zhuǎn)變?yōu)閮?yōu)化問(wèn)題，可將優(yōu)化目標(biāo)函數(shù)F(x)寫為

其中：x是待修正參數(shù)的向量；ri是第i個(gè)響應(yīng)的殘差函數(shù)；

ωi是第i個(gè)殘差函數(shù)的權(quán)重系數(shù)。

1.2 Kriging模型

傳統(tǒng)參數(shù)型模型修正方法需要不斷地調(diào)用有限元模型參與計(jì)算，由于有限元模型計(jì)算的黑箱特性，通過(guò)建立一個(gè)可以模擬輸入與輸出之間函數(shù)關(guān)系的顯式代理模型可以大幅減少計(jì)算量，Kriging模型就是一種常用的代理模型[4，6]。Kriging模型可以分為回歸部分與隨機(jī)響應(yīng)部分，寫為：

其中：y(x)是Kriging模型的響應(yīng)預(yù)估值；f(x)是Kriging模型的回歸部分，是關(guān)于x的多項(xiàng)式函數(shù)，用來(lái)擬合整體趨勢(shì)；Z(x)是Kriging模型的隨機(jī)響應(yīng)部分，用來(lái)擬合局部變化。Kriging模型不僅能得到預(yù)估值，還能得到方差，并利用方差在設(shè)計(jì)空間中重新取樣，對(duì)自身進(jìn)行更新從而獲得較高的模型預(yù)測(cè)精度。重新取樣依據(jù)的規(guī)則被稱為加點(diǎn)準(zhǔn)則，最大期望改善（Expected Improvement，EI）準(zhǔn)則是使用較為廣泛的加點(diǎn)準(zhǔn)則，其將EI函數(shù)定為評(píng)價(jià)指標(biāo)來(lái)表征改善程度I(x)，并將EI函數(shù)最大值的設(shè)計(jì)點(diǎn)作為新設(shè)計(jì)點(diǎn)添加到設(shè)計(jì)點(diǎn)集中更新Kriging模型，如此循環(huán)直至滿足收斂準(zhǔn)則。改善程度I(x)可以寫為

其中：Φ是標(biāo)準(zhǔn)正態(tài)累計(jì)分布函數(shù)；ψ是標(biāo)準(zhǔn)正態(tài)分布概率密度函數(shù)。通過(guò)不斷迭代添加新設(shè)計(jì)點(diǎn)更新Kriging模型，可以使得模型自適應(yīng)進(jìn)化，從而滿足工程需要。

2 滑鞍模態(tài)試驗(yàn)

滑鞍是加工中心主要移動(dòng)部件，直接與電主軸相連，提供加工中心Z軸移動(dòng)。為了獲得滑鞍的模態(tài)參數(shù)，對(duì)其進(jìn)行了力錘敲擊模態(tài)試驗(yàn)。實(shí)驗(yàn)過(guò)程中選用了東方所INV9832-50三向加速度傳感器與奇石樂(lè)9722A2000力錘，在滑鞍表面均布測(cè)點(diǎn)，共布置了20個(gè)測(cè)點(diǎn)。測(cè)試過(guò)程采用移動(dòng)傳感器法，保持力錘在同一敲擊點(diǎn)不變，移動(dòng)傳感器遍歷所有測(cè)點(diǎn)，數(shù)據(jù)采集依靠基鈦克Impaq頻譜分析儀，采樣頻率為5000 Hz，最終結(jié)果導(dǎo)入到計(jì)算機(jī)中，運(yùn)用專業(yè)軟件ME’SCOPE進(jìn)行模態(tài)分析，得到滑板模態(tài)參數(shù)，如表1所示。

表1 滑鞍模態(tài)試驗(yàn)實(shí)測(cè)值

3 構(gòu)建Kriging代理模型

3.1 有限元模型建立

滑鞍整體尺寸約為320 mm×220 mm×110 mm，根據(jù)滑鞍幾何尺寸在SoildWorks中建立起三維模型，將三維模型導(dǎo)入到Hypermesh 中劃分結(jié)構(gòu)化網(wǎng)格，網(wǎng)格模型如圖1所示。

圖1 滑鞍有限元網(wǎng)格模型

滑鞍依靠導(dǎo)軌滑塊和螺母絲桿與滑板相連，用彈簧阻尼單元對(duì)滑鞍邊界條件進(jìn)行模擬，用彈簧的剛度值代替接合面剛度值，其中導(dǎo)軌滑塊可以簡(jiǎn)化為法向彈簧與軸向彈簧，絲杠螺母可以簡(jiǎn)化為徑向彈簧和軸向彈簧，其接觸剛度值可以通過(guò)廠家手冊(cè)查得，滑塊的法向剛度k1為8.62×108N/m，切向剛度k2為5.64×108N/m；絲桿的徑向剛度k3為6.32×108N/m，軸向剛度k4為7.32×108N/m?；安牧蠟榛诣T鐵，材料參數(shù)取自ANSYS材料庫(kù)，彈性模量E為110 GPa，密度ρ為7200 kg/m3，泊松比μ為0.28。邊界條件設(shè)定模擬模態(tài)試驗(yàn)時(shí)的外部條件，將有限元模型導(dǎo)入到ANSYS Workbench中進(jìn)行有限元模態(tài)分析，與模態(tài)試驗(yàn)的結(jié)果進(jìn)行對(duì)比，如表2所示。

表2 滑鞍固有頻率實(shí)測(cè)值與計(jì)算值的對(duì)比

從對(duì)比結(jié)果可以看出，計(jì)算固有頻率和實(shí)測(cè)固有頻率最大相對(duì)誤差為19.37%，發(fā)生在第六階，前六階的平均誤差為11%左右，表明有限元模型雖然能一定程度上模擬真實(shí)模型，但仍需要進(jìn)一步修正以滿足工程需求。

3.2 參數(shù)選擇與靈敏度分析

由于滑鞍材料參數(shù)選用的是ANSYS材料庫(kù)中灰鑄鐵的材料參數(shù)，與滑鞍實(shí)際材料可能有所偏差，故將滑鞍材料參數(shù)納入修正參數(shù)的選擇當(dāng)中。導(dǎo)軌滑塊和絲桿螺母的剛度參數(shù)是從廠家手冊(cè)上查得，但在有限元模型中將其簡(jiǎn)化為彈簧阻尼單元，因此將彈簧剛度也納入修正參數(shù)當(dāng)中。綜上所述，初步選擇了滑鞍密度參數(shù)ρ、滑鞍彈性模量E、滑鞍泊松比μ、法向彈簧剛度k1、切向彈簧剛度k2、徑向彈簧剛度k3和軸向彈簧剛度k4等7個(gè)結(jié)構(gòu)參數(shù)作為設(shè)計(jì)參數(shù)。

對(duì)上述7個(gè)待修正參數(shù)進(jìn)行模態(tài)頻率靈敏度計(jì)算，計(jì)算結(jié)果如圖2所示?？梢钥闯龌懊芏圈?、彈性模量E、法向彈簧剛度k1、徑向彈簧剛度k3和軸向彈簧剛度k4對(duì)滑鞍模態(tài)頻率的靈敏度較高，因此，有限元模型修正的修正參數(shù)選 擇 ρ、E、k1、k3和k4等5 個(gè)設(shè)計(jì)參數(shù)。結(jié)合現(xiàn)實(shí)考慮，修正參數(shù)的取值范圍在其初值的0.8～1.2倍之間。

圖2 設(shè)計(jì)參數(shù)靈敏度分析

3.3 自適應(yīng)Kriging模型構(gòu)建

實(shí)驗(yàn)設(shè)計(jì)是一種基于數(shù)理統(tǒng)計(jì)和概率論的科學(xué)方法，可以通過(guò)較少的試驗(yàn)次數(shù)得到輸入與輸出間的關(guān)系。拉丁超立方實(shí)驗(yàn)設(shè)計(jì)是蒙特卡羅方法的一種，具有較高的采樣效率，在實(shí)際應(yīng)用中使用較廣。五設(shè)計(jì)參數(shù)的拉丁超立方實(shí)驗(yàn)設(shè)計(jì)，取樣方式選為中心復(fù)合設(shè)計(jì)，初始設(shè)計(jì)點(diǎn)共有27個(gè)。將選定的參數(shù)提交到有限元模型中參與計(jì)算，得到對(duì)應(yīng)的前六階固有頻率數(shù)值，根據(jù)此輸入與輸出關(guān)系可以構(gòu)建出初始的Kriging模型，然后可以根據(jù)EI加點(diǎn)準(zhǔn)則增加細(xì)化點(diǎn)，以0.01為收斂值，最終新增了4個(gè)細(xì)化點(diǎn)，構(gòu)建成立自適應(yīng)Kriging模型。以法向彈簧剛度k1和徑向彈簧剛度k3與第一階固有頻率f1的關(guān)系為例，標(biāo)準(zhǔn)和自適應(yīng)Kriging響應(yīng)面模型如圖3所示。

從圖3可以看出，自適應(yīng)Kriging模型在函數(shù)的極大值點(diǎn)附近插入了新設(shè)計(jì)點(diǎn)，從單點(diǎn)極大值變成了兩點(diǎn)極大值。從擬合優(yōu)度來(lái)看，標(biāo)準(zhǔn)Kriging模型的擬合優(yōu)度為96.52%，自適應(yīng)Kriging 模型的擬合優(yōu)度為99.24%，相較于標(biāo)準(zhǔn)Kriging 模型擬合精度得到了提高，是對(duì)有限元模型更高精度的近似。

圖3 標(biāo)準(zhǔn)Kriging和自適應(yīng)Kriging一階固有頻率響應(yīng)面

4 模型修正

基于自適應(yīng)Kriging模型，采用遺傳算法搜尋設(shè)計(jì)空間內(nèi)目標(biāo)函數(shù)的最小值，遺傳算法初代種群大小設(shè)定為1000，子代種群大小設(shè)定為500，最大迭代代數(shù)為20，允許的最大帕累托百分比設(shè)定為70%，穩(wěn)定收斂閾設(shè)定為2%，修正前后各參數(shù)與結(jié)果的對(duì)比列在表3。從中可以看出，變動(dòng)量較大的為滑鞍的彈性模量、滑塊法向接觸剛度與絲桿軸向接觸剛度，但均未接近設(shè)計(jì)空間邊緣，有較高的置信度。

表3 修正參數(shù)在模型修正前后的對(duì)比

表4給出了修正之后的有限元模型前六階模態(tài)計(jì)算頻率與實(shí)測(cè)頻率的對(duì)比結(jié)果。從中可以看出，與未修正的有限元模型相比，修正后的模型與實(shí)際模型的相似度更高，最大相對(duì)誤差從原先的19.37%降低為4.52%，最大相對(duì)誤差在5%以內(nèi)，可以滿足工程計(jì)算需求。

表4 模型修正后計(jì)算頻率與實(shí)測(cè)頻率的對(duì)比

5 結(jié)語(yǔ)

基于自適應(yīng)Kriging模型，本文對(duì)有限元模型進(jìn)行了參數(shù)修正，較好地平衡了計(jì)算效率與計(jì)算精度，最終修正后的有限元模型與實(shí)測(cè)模型相對(duì)誤差降至5%以下，可以滿足工程計(jì)算需求。