改進(jìn)粒子群算法與fmincon函數(shù)混合尋優(yōu)的平面度、垂直度誤差評定

許家赫 陳岳坪

摘? ? 要：平面度、垂直度誤差是檢驗(yàn)幾何產(chǎn)品互換性的重要參考指標(biāo)，對使用性能檢測及產(chǎn)品品質(zhì)的評價(jià)有著重大影響.最小二乘法雖然使用廣泛，卻只能適用于對測量精度需求度不高的應(yīng)用場合.為使平面度、垂直度誤差評定值精準(zhǔn)度更高，依據(jù)新一代產(chǎn)品幾何技術(shù)規(guī)范與驗(yàn)證為框架，采用最小區(qū)域法，建立平面度、垂直度誤差評定數(shù)學(xué)模型.運(yùn)用一種根據(jù)粒子早熟收斂程度進(jìn)行自適應(yīng)調(diào)整權(quán)重的PSO算法與fmincon函數(shù)混合尋優(yōu)的方法實(shí)現(xiàn)對平面度、垂直度誤差的優(yōu)化與評定.通過測量箱體尺寸數(shù)據(jù)分析計(jì)算，對比文獻(xiàn)數(shù)據(jù)結(jié)果進(jìn)行驗(yàn)證，證明了該評定方法的有效性及評定結(jié)果精確性.

關(guān)鍵詞：平面度誤差;垂直度誤差;fmincon函數(shù);自適應(yīng)調(diào)整權(quán)重;混合尋優(yōu)

中圖分類號：TG83? ? ? ? ? DOI：10.16375/j.cnki.cn45-1395/t.2019.04.016

0? ? 引言

隨著社會生產(chǎn)對產(chǎn)品使用要求日趨增高，平面度和垂直度作為評價(jià)幾何產(chǎn)品位置精度的重要影響因素，是產(chǎn)品應(yīng)用性能的重要參考指標(biāo).國外研究起步早，提出了引入垂直度偏轉(zhuǎn)參數(shù)，進(jìn)而校正補(bǔ)償[1]等方法.國內(nèi)則多用最小二乘法（Least Square Method，LSM）構(gòu)建基準(zhǔn)直線[2]的方法.最小二乘法由于應(yīng)用便捷，目前廣泛應(yīng)用，缺點(diǎn)是誤差精度不高，結(jié)果不穩(wěn)定.考慮到產(chǎn)品生產(chǎn)和社會應(yīng)用的實(shí)際情況，對平面度、垂直度誤差精準(zhǔn)而有效計(jì)算的需求更為迫切.

利用最小區(qū)域法[3]，根據(jù)新一代產(chǎn)品幾何技術(shù)規(guī)范（Geometrical Product Specification and verification，GPS）體系標(biāo)準(zhǔn)，計(jì)算平面度誤差和線對面的垂直度誤差，把一種以粒子的早熟收斂程度為調(diào)整權(quán)重依據(jù)[4]的自適應(yīng)性PSO算法與fmincon函數(shù)進(jìn)行混合尋優(yōu)，并且運(yùn)用計(jì)算基準(zhǔn)要素與被測要素之間的最大變化范圍為最小的原則，用以對測量數(shù)據(jù)進(jìn)行有效評定.

1? ? 構(gòu)建目標(biāo)函數(shù)模型

1.1? ?建立基準(zhǔn)平面

因加工誤差無法消除，導(dǎo)致實(shí)際工作表面非理想狀態(tài)，實(shí)際表面不能有效且完整表達(dá)，以新一代GPS誤差體系標(biāo)準(zhǔn)為依據(jù)，應(yīng)用前需要對實(shí)際表面進(jìn)行檢測、擬合、操作、認(rèn)定[5]，并用認(rèn)證表面代替實(shí)際工作面.平面評定過程中運(yùn)用此方法時(shí)，因基準(zhǔn)無法直接獲取理想要素，所以誤差基準(zhǔn)要素評定的精準(zhǔn)度最為關(guān)鍵，故可使用分離、提取等操作手法[6]，用以高效分析和研究誤差基準(zhǔn)選取的正確與否，將基準(zhǔn)理想化，從而達(dá)到使用要求.

按照平面度的定義，用兩平行平面將被測區(qū)域中所有點(diǎn)完全包圍，兩平面的最小間距即稱作平面度誤差.若被測平面內(nèi)所有距離的最大值與最小值的絕對值之差為最小，則根據(jù)定義該平面就是要被創(chuàng)建的基準(zhǔn)平面，上述計(jì)算所得差值的大小便是所求平面度的誤差值.

設(shè)空間坐標(biāo)oxyz的基準(zhǔn)擬合平面方程為：

[z=ax+by+c]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?   ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（1）

式中，a、b、c為待定系數(shù).

被測輪廓上的點(diǎn)則設(shè)為[Pi（xi，? yi，? zi）]，[（i=1，? 2，? 3，? …，? n）]，接著需要計(jì)算該點(diǎn)到基準(zhǔn)平面距離d：

[di=z-axi-byi-ci1+a2+b2]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （2）

最后得到滿足最小區(qū)域法的平面的目標(biāo)函數(shù)：

[f（a，? b，? c）=min（dmax-dmin）]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （3）

使目標(biāo)函數(shù)[fa，? b，? c]的函數(shù)值尋優(yōu)為最小，即將計(jì)算擬合平面問題轉(zhuǎn)化為優(yōu)化[a、b、c]值的問題.

1.2? ?基準(zhǔn)平面坐標(biāo)變換

任取基準(zhǔn)平面上某一點(diǎn)[P0（x0，? y0，? z0）]作原點(diǎn)，創(chuàng)建新坐標(biāo)系[OXYZ]，且令[ox//OX]，設(shè)[XOY]新坐標(biāo)平面同[xoy]坐標(biāo)平面的夾角α，坐標(biāo)系轉(zhuǎn)換如圖1所示.

將原坐標(biāo)系[o，? x，? y，? z]與新坐標(biāo)系（[O]，X，Y，Z）按照轉(zhuǎn)換矩陣作如下變換[7]：

[XYZ=-ac1-c2 ? ? -bc1-c2 ? 1-c2 b1-c2 -a1-c2 0 ? ? ? ? ? a ? ? ? ? b cx-x0y-y0z-z0]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（4）

其中，[x0、? y0、? z0]為原坐標(biāo)系原點(diǎn)坐標(biāo).實(shí)測點(diǎn)[Pj（xj， yj， zj）]，[（j=1， 2， 3，…， n）]則通過矩陣轉(zhuǎn)換得到坐標(biāo)系理論新測點(diǎn)[Pj（Xj，? Yj，? Zj）]，[（j=1， 2， 3，…， n）].

1.3? ?關(guān)于軸線相對平面的任意方向上的垂直度誤差評定

設(shè)向量[S={p，? q，? r}]作擬合直線方向向量.參數(shù)[p、q、r]是其分量表現(xiàn)形式，利用最小區(qū)域法作擬合依據(jù)，計(jì)算被測直線的測點(diǎn)集.可知其理想直線就是以被測直線的測點(diǎn)集算術(shù)平均中心點(diǎn)為基點(diǎn)的最小區(qū)域擬合直線.根據(jù)上述已知，運(yùn)用被測直線算術(shù)平均中心點(diǎn)P0（X0，Y0，Z0）為原點(diǎn)，計(jì)算出最小區(qū)域擬合直線方程：

[X-X0p=Y-Y0q=Z-Z0r]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （5）

任意測點(diǎn)[PjXj，? Yj，? Zj]與擬合直線間距離D：

[Dj=qXj-X0-pYj-Y02+rXj-X0-pZj-Z02+r（Yj-Y0）-q（Zj-Z0）2p2+q2+r2]? ? ? （6）

最后計(jì)算直線度誤差：

[f（p，? q，? r）=Dmax]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（7）

要求得滿足最小區(qū)域擬合直線的p、q、r值，需要將式（6）優(yōu)化計(jì)算取得最大值后再代入式（7）尋求最小值.

所求任意方向的傾斜度誤差值可以根據(jù)定向誤差評定理論知，是以擬合直線為軸線的，滿足最小條件的圓柱面直徑值.通過對理想正確角度α進(jìn)行設(shè)置，例如α=0°或90°時(shí)，線對面傾斜度誤差則會轉(zhuǎn)化成為平行度誤差及垂直度誤差，因此垂直度誤差為：

[F（p，? q，? r）=max（2Dj）]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （8）

2? ?自適應(yīng)粒子群優(yōu)化與fmincon函數(shù)混合尋優(yōu)操作步驟

2.1? ?粒子群算法和fmincon優(yōu)化函數(shù)

粒子群優(yōu)化算法（Particle Swarm Optimization，PSO），其特點(diǎn)是將解的值初始隨機(jī)化，利用適應(yīng)度評價(jià)解的品質(zhì)，迭代尋求最優(yōu)解.每個(gè)粒子均“配置”有適應(yīng)值和速度矢量，通過跟隨當(dāng)前空間里搜索到的最優(yōu)值來搜尋全局最優(yōu)[8].

將隨機(jī)產(chǎn)生N個(gè)粒子當(dāng)做一個(gè)部落，每次迭代粒子會通過跟蹤兩個(gè)極值自我更新，分別是個(gè)體極值pbest，和全局極值gbest.搜索到最優(yōu)值后，粒子依照如下公式更新速度和位置[9]：

[vβd=ωvβd+c1r1（pβd-xβd）+c2r2（pgd-xβd）]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（9）

[xβd=xβd+vβd]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（10）

其中，[vβd]表示粒子[β]在d維時(shí)的速度，[pβd]表示粒子[β]在d維時(shí)的位置，[ω]為慣性權(quán)重. fmincon函數(shù)通過MATLAB優(yōu)化系統(tǒng)工具箱提供，運(yùn)行時(shí)首先檢查有無梯度表達(dá)提供，從而對目標(biāo)函數(shù)進(jìn)行不同的算法迭代運(yùn)算.因此對非線性有約束優(yōu)化問題可以快速且高效的解決.

2.2? ?以早熟收斂程度作權(quán)重調(diào)節(jié)依據(jù)的自適應(yīng)粒子群算法[10]

該算法根據(jù)種群里早熟收斂情況和個(gè)體適應(yīng)值為參考依據(jù)，動(dòng)態(tài)調(diào)節(jié)權(quán)重變化.設(shè)定粒子[pi]的適應(yīng)值為[fi]，最優(yōu)粒子適應(yīng)度[fm]，則平均適應(yīng)值[favg=1ni=1nfi]，計(jì)算中將優(yōu)于平均適應(yīng)值粒子的適應(yīng)值求平均記為（[favg]），定義[Δ=fm-favg].根據(jù)適應(yīng)值不同可將群體分成3部分，可將慣性權(quán)重[ω]作如下調(diào)節(jié)：

1）若[fi]優(yōu)于[favg]，則

[ω=ω-（ω-ωmin）·fi-favgfm-favg]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （11）

2）若[fi]優(yōu)于[favg]小于[fm]，慣性權(quán)重不變.

3）若[fi]小于[favg]，則

[ω=1.5-11+k1·exp（-k2·Δ）]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （12）

其中，k為控制參數(shù)，k1控制[ω]的上限，k2控制權(quán)重的調(diào)節(jié)能力.若粒子處于疏散狀態(tài)，[Δ]則較大，[ω]變小，此時(shí)算法局部搜索能力加強(qiáng)，使群體趨于收斂;若粒子表現(xiàn)為高度聚集狀態(tài)，[Δ]則較小，[ω]變大，從而快速有效跳出局部最優(yōu)解的情況.

該方法能使粒子具有更強(qiáng)適應(yīng)性，提高了局部搜索能力，通過與fmincon優(yōu)化函數(shù)進(jìn)行混合尋優(yōu)，更顯著提高了收斂速度與計(jì)算精度，增強(qiáng)了粒子多樣性，提供了自我調(diào)節(jié)能力，并且優(yōu)化了粒子全局搜索能力.

2.3? ?基本實(shí)現(xiàn)步驟[11]

Step1? 隨機(jī)化種群中粒子位置和速度.

Step2? 計(jì)算各粒子適應(yīng)值的大小，并對每個(gè)粒子進(jìn)行評價(jià).同時(shí)將相關(guān)粒子的位置和適應(yīng)值的信息存儲于個(gè)體極值pbest中，將最優(yōu)的個(gè)體位置和適應(yīng)值保存在全局極值gbest中.

Step3? 根據(jù)式（9）、式（10）更新位置和速度.

Step4  根據(jù)式（11）、式（12）更新權(quán)重.

Step5? 比較當(dāng)前所有的pbest和gbest，并更新gbest.

Step6? ?若達(dá)到停止條件，停止搜索并輸出結(jié)果，否則返回Step3繼續(xù)搜索，直到滿足停止準(zhǔn)則為止.

Step7? ?為防止搜索的結(jié)果有局部最優(yōu)解的情況出現(xiàn)，將粒子群優(yōu)化結(jié)果值輸入fmincon函數(shù)進(jìn)行處理，進(jìn)行最小值尋優(yōu)運(yùn)算.

Step8? ?通過output指令輸出最終優(yōu)化值.

3? ? 實(shí)例驗(yàn)證

實(shí)驗(yàn)采用三坐標(biāo)測量機(jī)作為檢測工具，將箱體任意放置于工作臺上，進(jìn)行測量前，先創(chuàng)立基于工件位置的測量坐標(biāo)系，錄入坐標(biāo)值時(shí)是用手動(dòng)采點(diǎn)，最后通過數(shù)據(jù)處理系統(tǒng)提取測點(diǎn)坐標(biāo).箱體尺寸數(shù)據(jù)來源文獻(xiàn)[11].

平面的優(yōu)化目標(biāo)函數(shù)為：

[f（a，b，c）=max（z-axi-byi-c）-min（z-axi-byi-c）1+a2+b2]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（13）

設(shè)平面[z=a0x+b0y+c0]為初始平面，[Δa、Δb、Δc]為對應(yīng)參數(shù)增量，優(yōu)化目標(biāo)函數(shù)[f]則又可表示為[f（a0+Δa， b0+Δb， c0+Δc）]，令粒子編碼表示為[T（Δa， Δb， Δc）]，問題搜索維數(shù)[D=3]，可將該優(yōu)化目標(biāo)函數(shù)轉(zhuǎn)換成[dΔa， Δb， Δc]計(jì)算，目標(biāo)優(yōu)化函數(shù)采用適應(yīng)度函數(shù)：

[f（T）=d（Δa， Δb， Δc）=min（dmax-dmin）]? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（14）

初始化種群粒子個(gè)數(shù)[N=50]，其中粒子維數(shù)D=3，令學(xué)習(xí)因子c1=c2=2，隨機(jī)產(chǎn)生初始加速權(quán)重，最大迭代次數(shù)M不超過算法要求的最大值.輸入初始基準(zhǔn)平面數(shù)值，a=﹣0.000 09，b=﹣0.000 49，由于平面度誤差跟基準(zhǔn)位置無關(guān)，c值只是起偏置作用，c值可以任意選取，這里設(shè)c=﹣1.經(jīng)過混合尋優(yōu)計(jì)算，優(yōu)化后基準(zhǔn)平面參數(shù)為：a=﹣0.000 091 335 529，b=﹣0.000 381 215 095，c=﹣1.將上述結(jié)果代入可得

f =0.058 950 905 044 mm.

根據(jù)文獻(xiàn)[5]和文獻(xiàn)[11]中幾種不同的算法評定所得結(jié)果值進(jìn)行分析比較，結(jié)果如表1所示.

根據(jù)數(shù)據(jù)對比得出：對于擬合基準(zhǔn)平面，運(yùn)用混合尋優(yōu)算法，使所得基準(zhǔn)平面計(jì)算精度更高更準(zhǔn)確，誤差值更小，創(chuàng)建的基準(zhǔn)平面擬合效果更好，尋優(yōu)值更加穩(wěn)定.

圖2為運(yùn)行算法時(shí)，平面度誤差收斂曲線.觀察收斂數(shù)值發(fā)現(xiàn)，從第3代到第10代誤差值均為23.844 406 099 134 mm，短暫陷入了局部極值的情況.當(dāng)誤差計(jì)算到22代0.092 960 135 913 mm時(shí)，達(dá)到穩(wěn)定狀態(tài)，直到運(yùn)行到72代，得到最小誤差? ? ? ? ?0.058 950 846 983 mm，算法結(jié)束.圖中也明顯可以看出，這種混合尋優(yōu)快速穩(wěn)定的特點(diǎn)，能夠有效擺脫局部最小值困境且沒有大范圍震蕩現(xiàn)象發(fā)生，不僅適應(yīng)性強(qiáng)而且相比文獻(xiàn)精度高.

接著對垂直度誤差進(jìn)行計(jì)算，垂直度誤差值尋優(yōu)迭代過程的情況如圖3所示.通過觀察迭代計(jì)算值，可以發(fā)現(xiàn)當(dāng)運(yùn)行到第2代到第9代時(shí)誤差值均為0.804 020 151 392 mm，陷入了局部最小值情況，當(dāng)計(jì)算到第14代，誤差為0.506 726 534 095 mm時(shí)基本穩(wěn)定，直到55代，達(dá)到最優(yōu)值0.499 839 967 830 mm.從圖中也能明顯看出當(dāng)尋優(yōu)迭代到15代時(shí)，垂直度誤差值基本穩(wěn)定在0.5? mm左右，等到了20代以后迭代后續(xù)過程平穩(wěn)，沒有突變值的發(fā)生，直到算法結(jié)束.

4? ? 總結(jié)

本文建立了滿足最小區(qū)域法的平面度和垂直度誤差評定模型，通過使用粒子群算法與fmincon函數(shù)混合尋優(yōu)的方法，有效地解決了非線性約束優(yōu)化問題.該方法能精準(zhǔn)地計(jì)算出基準(zhǔn)平面擬合所需各項(xiàng)參數(shù)值，同時(shí)提高了計(jì)算結(jié)果的精度，并擬合出了精度更高的基準(zhǔn)平面，最后實(shí)現(xiàn)了平面度和線對面垂直度誤差的評定.

參考文獻(xiàn)

[1]? ? ?ZHAO G Y，WANG S J，SUI W T.Compensation approach to error of perpendicularity and mismatch error of position loop gains in CNC machine[J].Applied Mechanics and Materials， 2012（121﹣126）：3190﹣3194.

[2]? ?  郭誼.直線垂直度誤差的一種評定方法[J].計(jì)量與測試技術(shù)，2008（5）：21﹣22.

[3]? ? ?黃富貴.平面度誤差各種評定方法的比較[J].工具技術(shù)，2007（8）：107﹣109.

[4]? ? ?溫正，孫克華.MATLAB智能算法[M].北京：清華大學(xué)出版社，2017.

[5]? ? ?李向前. 基于新一代GPS的形狀公差規(guī)范設(shè)計(jì)方法研究[D].桂林：桂林電子科技大學(xué)，2010.

[6]? ? ?馬利民. 新一代產(chǎn)品幾何量技術(shù)規(guī)范（GPS）理論框架體系及關(guān)鍵技術(shù)研究[D].武漢：華中科技大學(xué)，2006.

[7]? ? ?周東棟，樊軍.基于改進(jìn)粒子群算法的垂直度誤差評定[J].現(xiàn)代制造工程，2017（2）：25﹣29.

[8]? ? ?趙乃剛，鄧景順.粒子群優(yōu)化算法綜述[J].科技創(chuàng)新導(dǎo)報(bào)，2015 （26）：216﹣217.

[9]? ? ?張頂學(xué)，關(guān)治洪，劉新芝.一種動(dòng)態(tài)改變慣性權(quán)重的自適應(yīng)粒子群算法[J].控制與決策，2008，23（11）：1253﹣1257.

[10]? ?杜振鑫，王兆青.一種改進(jìn)的動(dòng)態(tài)改變慣性權(quán)重的粒子群算法[J].微電子學(xué)與計(jì)算機(jī)，2011，28（3）：85﹣88.

[11]? ?高明，鮑家定，周晴倫.基于粒子群算法的線對面的垂直度誤差評定[J].中小企業(yè)管理與科技，2012（6）：203﹣205.

Evaluation of planarity and verticality errors for mixed optimization of improved particle swarm optimization and fmincon function

XU Jiahe， CHEN Yueping*

（School of Mechanical and Traffic Engineering， Guangxi University of Science and Technology，

Liuzhou 545006，China）

Abstract： Flatness and verticality errors are important reference indicators for testing the interchangeability of geometric products， which have a significant impact on performance testing and product? ?quality evaluation. Although the least squares method is widely used， it can only be applied to? ? ? ? ? ? ?applications where the measurement accuracy is not required. In order to improve the accuracy of? ? ?flatness and verticality error evaluation， a mathematical model for flatness and verticality error? ? ? ? evaluation was established based on the framework of new generation product geometry specification and verification by using the minimum area method. In this paper， a mixed optimization method of PSO algorithm and fmincon function is used to optimize and evaluate the flatness and verticality errors? ? ? ?according to the degree of premature convergence of particles. The validity and accuracy of the? ? ? ? evaluation method are proved by analyzing and calculating the measurement box size data and? ? ? ? comparing the results of literature data.

Key words： flatness error; verticality error; fmincon function; adaptive weight adjustment; hybrid? ? ? optimization

（責(zé)任編輯：黎? ?婭）