Monte-Carlo 方法在尺寸鏈方程組計(jì)算中的應(yīng)用

蔡偉妹，仲其伐

(南京理工大學(xué) 機(jī)械工程學(xué)院，江蘇 南京 210094)

0 引言

引信作為一種微型化和精密化的軍工產(chǎn)品，技術(shù)和尺寸精度要求都比較高，因此尺寸鏈計(jì)算是引信產(chǎn)品設(shè)計(jì)中的一個(gè)重要環(huán)節(jié)，是引信設(shè)計(jì)的特點(diǎn)之一。引信產(chǎn)品的尺寸鏈計(jì)算主要為零件設(shè)計(jì)、產(chǎn)品裝配、制造工藝、發(fā)火性能以及自毀性能等的計(jì)算和校核提供數(shù)據(jù)，通過(guò)理論計(jì)算及時(shí)發(fā)現(xiàn)產(chǎn)品設(shè)計(jì)的缺陷，改善產(chǎn)品的可生產(chǎn)性、可裝配性以及產(chǎn)品使用的安全性與可靠性。

引信產(chǎn)品設(shè)計(jì)一般采用極大值極小值法計(jì)算尺寸鏈。文獻(xiàn)［1］從時(shí)間發(fā)生的概率大小和生產(chǎn)經(jīng)濟(jì)性以及可生產(chǎn)性方面出發(fā)，提出用Monte-Carlo 方法解尺寸鏈的問(wèn)題，將Monte-Carlo 方法應(yīng)用到一般尺寸鏈的計(jì)算中;文獻(xiàn)［2］系統(tǒng)的介紹了求解由方程組組成的尺寸鏈問(wèn)題的方法;文獻(xiàn)［3］通過(guò)Monte-Carlo 法求解多環(huán)尺寸鏈估算裝配的成功率問(wèn)題;文獻(xiàn)［4］和文獻(xiàn)［5］則應(yīng)用了Monte-Carlo 方法研究球轉(zhuǎn)子旋轉(zhuǎn)偏心情況，針對(duì)相關(guān)原則的同軸度誤差的處理問(wèn)題，通過(guò)對(duì)計(jì)算條件的假設(shè)與抽樣檢驗(yàn)，得出了精度比較高的計(jì)算仿真結(jié)果。文獻(xiàn)［7］解決了用常規(guī)方法對(duì)尺寸鏈組成環(huán)的尺寸分布服從非正態(tài)分布而無(wú)法求解的問(wèn)題。

本文所要做的工作，是在文獻(xiàn)［1］～［6］所做工作的基礎(chǔ)上，將Monte-Carlo 方法具體應(yīng)用到由方程組構(gòu)成的尺寸鏈計(jì)算的實(shí)際問(wèn)題中，進(jìn)行分析求解，得出所求尺寸鏈封閉環(huán)的尺寸分布情況，進(jìn)而求出封閉環(huán)的極值。

1 尺寸鏈計(jì)算方法

一般來(lái)說(shuō)，進(jìn)行尺寸鏈計(jì)算的主要目的是求解極值。所使用的計(jì)算方法一般有三類:極大值極小值法、概率統(tǒng)計(jì)法和Monte-Carlo 法，另外還有高等數(shù)學(xué)中的函數(shù)微分求極值法。

1.1 極大值極小值法

極大值極小值法是尺寸鏈計(jì)算中應(yīng)用得較為廣泛的一種計(jì)算方法。應(yīng)用該方法進(jìn)行計(jì)算的基本原則［7］是:

1)計(jì)算封閉環(huán)最大尺寸時(shí)，所有增環(huán)都取最大極限尺寸，所有減環(huán)都取最小極限尺寸。

2)計(jì)算封閉環(huán)最小尺寸時(shí)，所有增環(huán)都取最小極限尺寸，所有減環(huán)都取最大極限尺寸。

3)計(jì)算封閉環(huán)基本尺寸時(shí)，所有組成環(huán)都取基本尺寸。

4)封閉環(huán)的最大最小尺寸的差值即為封閉環(huán)的公差。

極大值極小值法考慮了極端不利的情況，一般情況下極端條件出現(xiàn)的可能性很小，把這種極端情況考慮進(jìn)去所計(jì)算出來(lái)的結(jié)果使得尺寸公差范圍小，生產(chǎn)的難度變大。一般來(lái)說(shuō)，極大值極小值法適用于精度要求比較高，使用環(huán)境比較極端多變并要求產(chǎn)品裝配的100%互換的產(chǎn)品設(shè)計(jì)上。從引信設(shè)計(jì)角度和使用環(huán)境來(lái)說(shuō)，引信的尺寸鏈計(jì)算一般采用極大值極小值法。

1.2 概率統(tǒng)計(jì)法

概率統(tǒng)計(jì)法從概率原理出發(fā)，考慮尺寸鏈中各組成環(huán)的實(shí)際分布情況，能客觀的反映加工精度和尺寸分布的本質(zhì)。概率統(tǒng)計(jì)法的計(jì)算方法是，將各組成環(huán)的尺寸分布看作某一特定形態(tài)的分布，然后根據(jù)這一形態(tài)分布進(jìn)行計(jì)算［7-8］。

一般情況下，當(dāng)各組成環(huán)的分布規(guī)律按正態(tài)分布時(shí)，其封閉環(huán)尺寸也必然符合正態(tài)分布規(guī)律;當(dāng)各組成環(huán)分布不按正態(tài)分布，在組成環(huán)數(shù)不太少，各組成環(huán)變化范圍大小相差不大的情況下，封閉環(huán)的尺寸仍趨于正態(tài)分布;當(dāng)尺寸鏈組成環(huán)數(shù)較少而各組成環(huán)分布又偏離正態(tài)分布較大時(shí)，封閉環(huán)也將偏離正態(tài)分布。

一般而言，應(yīng)用概率統(tǒng)計(jì)法計(jì)算引信零部件的尺寸鏈時(shí)，必須明確知道各組成環(huán)的尺寸分布形態(tài)。但是要明確知道各組成環(huán)的分布律是很困難的，需要長(zhǎng)時(shí)間的統(tǒng)計(jì)資料積累，因此，概率統(tǒng)計(jì)法只適用于尺寸分布律已知的情況下的尺寸鏈計(jì)算。

1.3 微分求導(dǎo)法

在尺寸鏈計(jì)算中，由函數(shù)式表示的尺寸鏈可利用微分法來(lái)求解，該方法的本質(zhì)是通過(guò)微分來(lái)求解函數(shù)的極值。在使用微分法進(jìn)行尺寸鏈計(jì)算之前，必須建立封閉環(huán)與各組成環(huán)之間的函數(shù)關(guān)系。

設(shè)封閉環(huán)與各組成環(huán)的函數(shù)關(guān)系為:

式中:f—封閉環(huán);

n—組成環(huán)個(gè)數(shù)。

當(dāng)各個(gè)組成環(huán)xi之間相互獨(dú)立時(shí)，對(duì)函數(shù)全微分，得:

假設(shè)x1，…，xi－1，xi+1，…xn為已知常數(shù)，xi為未知變量，求出xi在其取值范圍內(nèi)的增減性，從而判斷其取何值時(shí)使得f 取極值。

當(dāng)各組成環(huán)相關(guān)時(shí)，對(duì)式(3)進(jìn)行求偏微分，即:

假設(shè)x1，…，xi－1，xi+1，…，xn為已知常數(shù)，xi為未知變量，求出xi在其取值范圍內(nèi)的增減性，從而判斷其取何值時(shí)使得f 取極值。

由以上可知，在封閉環(huán)與各組成環(huán)之間的函數(shù)關(guān)系比較明確的情況下，可通過(guò)微分求導(dǎo)法求出封閉環(huán)的極值，從原理上看，該方法是通用的。但是，對(duì)于由方程組組成的尺寸鏈，對(duì)方程組進(jìn)行求導(dǎo)不僅繁瑣，并且各組成環(huán)的增減性不容易判斷，因此該方法就存在著一定的局限性。

1.4 Monte-Carlo 方法

a)Monte-Carlo 方法的特點(diǎn)

Monte-Carlo 方法被稱為隨機(jī)模擬方法(random simulation)，也被稱為隨機(jī)抽樣(random sampling)方法。該方法主要用于解決確定性的數(shù)學(xué)問(wèn)題和隨機(jī)性問(wèn)題［9-10］，是一種獨(dú)具風(fēng)格的數(shù)值計(jì)算方法。它的理論基礎(chǔ)來(lái)源于概率的大數(shù)定理和伯努利定理，其優(yōu)點(diǎn)以及與其他方法的不同點(diǎn)可歸納如下:

1)Monte-Carlo 方法及其程序結(jié)構(gòu)簡(jiǎn)單，只需要產(chǎn)生符合要求的隨機(jī)數(shù)，通過(guò)重復(fù)抽樣，求得平均值即可。

2)收斂的概率性和收斂速度與問(wèn)題的維數(shù)無(wú)關(guān)，Monte-Carlo 方法可適用于多維問(wèn)題的求解。Monte-Carlo方法的收斂是概率意義下的收斂，其收斂速度比一般數(shù)值方法的收斂速度要慢得多。

3)Monte-Carlo 方法通用性強(qiáng)，適用范圍廣，在求解問(wèn)題時(shí)受條件限制的影響小。

b)計(jì)算機(jī)偽隨機(jī)數(shù)的產(chǎn)生和抽樣方法

在應(yīng)用Monte-Carlo 方法模擬某問(wèn)題的求解過(guò)程時(shí)，需要產(chǎn)生各種概率分布的隨機(jī)變量。服從［0，1］分布的隨機(jī)變量是最簡(jiǎn)單、最基本并且是最重要的隨機(jī)變量也成為隨機(jī)數(shù)。其他分布的隨機(jī)變量的抽樣是通過(guò)隨機(jī)數(shù)來(lái)實(shí)現(xiàn)的。

引信零部件的尺寸分布可看作服從以下三類分布［5］:

1)均勻分布。服從［a，b］分布的偽隨機(jī)數(shù)可由式(4)生成:

2)正態(tài)分布。當(dāng)影響因素具有確定的公稱值和公差范圍時(shí)，便可認(rèn)為該影響因素服從正態(tài)分布。生成該類分布的表達(dá)式為:

其中:tN01，j=)，i=1，2，…，n。

3)瑞利分布。服從瑞利分布的隨機(jī)數(shù)可由式(6)生成:

隨機(jī)抽樣的方法一般有直接抽樣法、舍選抽樣法、復(fù)合抽樣法、復(fù)合舍選抽樣法、近似抽樣法和變換抽樣法六種。針對(duì)不同的影響因素分別進(jìn)行隨機(jī)抽樣，通過(guò)計(jì)算便可獲得相應(yīng)的尺寸分布結(jié)果。

需要補(bǔ)充說(shuō)明的是，在生成具有正態(tài)分布特性的隨機(jī)數(shù)后，需要對(duì)不在公差范圍內(nèi)的數(shù)值進(jìn)行剔除，以保證計(jì)算結(jié)果的正確性。

2 應(yīng)用Monte-Carlo 方法計(jì)算尺寸鏈的實(shí)例

現(xiàn)在以兩個(gè)例子具體說(shuō)明Monte-Carlo 方法在尺寸鏈方程組中的實(shí)際應(yīng)用。

2.1 求某引信體上端退刀槽斜面的寬度

a)問(wèn)題的引出

某引信體及其相關(guān)尺寸如圖1 所示。在該引信體零件中，需要計(jì)算引信體上端退刀槽斜面的寬度，為檢查風(fēng)帽與引信體安裝正確性提供數(shù)據(jù)。該斜面的寬度是通過(guò)一組方程組來(lái)表示的，無(wú)法通過(guò)簡(jiǎn)單的極大值極小值法來(lái)進(jìn)行計(jì)算，也無(wú)法通過(guò)微分求導(dǎo)法和概率統(tǒng)計(jì)法得出結(jié)果。

圖1 引信體及相關(guān)尺寸

為此，首先通過(guò)數(shù)學(xué)上的幾何關(guān)系，以引信體的中心軸線為y 軸來(lái)建立一個(gè)直角坐標(biāo)系，并分別將以上各個(gè)相關(guān)尺寸用a，b，c，d，e，f，g，h，α 來(lái)表示，引信體坐標(biāo)系建立如圖2 所示。

圖2 引信體零件坐標(biāo)系

在該坐標(biāo)系中，圖中標(biāo)注的A，B，C，D 四點(diǎn)的坐標(biāo)分別表示為(－a/2，e－d)，(－c/2，h)，(－b/2，e－df+g)，(x，y)。根據(jù)A，B 兩點(diǎn)坐標(biāo)與半徑r1，r2之間的數(shù)學(xué)關(guān)系以及C，D 兩點(diǎn)的關(guān)系，通過(guò)建立關(guān)系式，得到如下方程組:

方程組中的l 即為所要求解的退刀槽斜面的寬度。

通過(guò)上述所列方程組可以發(fā)現(xiàn)，一個(gè)看似很簡(jiǎn)單的問(wèn)題，已經(jīng)轉(zhuǎn)化為一系列的解釋式，并且該方程組無(wú)解析解。在此，將Monte-Carlo 方法應(yīng)用到該問(wèn)題的求解過(guò)程中。

b)計(jì)算的基本假設(shè)

在應(yīng)用Monte-Carlo 方法進(jìn)行尺寸鏈計(jì)算前，需要對(duì)引信體零件相關(guān)尺寸作以下假設(shè):

1)不考慮有關(guān)零件表面形狀誤差影響;

2)不考慮有關(guān)零件同軸度誤差的影響;

3)假設(shè)零件的相關(guān)尺寸均服從正態(tài)分布，取其公差帶中心為散布中心，公差帶寬度取6σ，這相當(dāng)于取工藝能力系數(shù)為1(引信生產(chǎn)工藝能力系數(shù)一般為1.2～1.5);

4)所有尺寸均以抽樣105次作為整體，以保證所得結(jié)果與實(shí)際相符合。

c)Monte-Carlo 法求解

通過(guò)上述假設(shè)，將Monte-Carlo 方法與MATLAB 編程相結(jié)合對(duì)該問(wèn)題進(jìn)行求解，程序流程圖如圖3 所示。

圖3 程序流程圖

在模擬的過(guò)程中，先隨機(jī)生成1.2×105個(gè)隨機(jī)數(shù)，將不在尺寸公差范圍內(nèi)的隨機(jī)數(shù)剔除，然后再?gòu)姆弦蟮碾S機(jī)數(shù)中隨機(jī)抽取105個(gè)隨機(jī)數(shù)代表整體，得出如圖4 所示頻數(shù)直方圖。

圖4 退刀槽斜面寬度尺寸頻數(shù)分布直方圖

注:正態(tài)分布抽樣的置信度為99.73%。

d)計(jì)算結(jié)果分析

通過(guò)上面的抽樣分析可以知道:

1)在置信度為99.73%的情況下，引信體上端退刀槽斜面寬度的最大值為1.238 mm，最小值為1.01 mm，平均值為1.12 mm，其分布曲線基本上呈正態(tài)分布。

2)在方程組沒(méi)有解析解的情況下，Monte-Carlo 方法給出了較為精確的計(jì)算結(jié)果，并且可以看出，在隨機(jī)抽樣得出的結(jié)果中，該尺寸值的取值是隨機(jī)的，并且大量分布在1.05 mm～1.2 mm 之間。

2.2 風(fēng)帽收口前的軸向間隙

a)問(wèn)題分析

某引信安裝風(fēng)帽后的裝配簡(jiǎn)圖如圖5 所示，在本文中只給出引信體與風(fēng)帽兩零件裝配的示意圖。

圖5 某引信安裝風(fēng)帽后的裝配簡(jiǎn)圖

在求解該問(wèn)題時(shí)，可應(yīng)用極大值極小值法;但通過(guò)計(jì)算發(fā)現(xiàn)，極大值極小值法計(jì)算結(jié)果的最小值為負(fù)值，與實(shí)際不符合。究其原因，是尺寸鏈計(jì)算中連續(xù)用了兩個(gè)通過(guò)其他尺寸鏈計(jì)算出來(lái)的引用尺寸。根據(jù)誤差傳遞可知，計(jì)算的結(jié)果是不準(zhǔn)確的。下面根據(jù)2.1 條中所用方法來(lái)建立幾何關(guān)系方程，坐標(biāo)系的建立與2.1 條相同，如圖6所示。

圖6 裝配圖坐標(biāo)系

方程的建立過(guò)程在此不再贅述。在該問(wèn)題中，要解決的是求出圖中B，C 兩點(diǎn)的坐標(biāo)，其中已知B，C 兩點(diǎn)的坐標(biāo)分別表示為(－i/2，y2)，(x1，e－d－f)。

將兩點(diǎn)的坐標(biāo)分別代入方程 (x－a0)2+(y－b0)2=，即可求得兩點(diǎn)的坐標(biāo)值，問(wèn)題的解即為:

b)基本假設(shè)

根據(jù)引信體與風(fēng)帽的裝配情況，現(xiàn)作以下假設(shè):

1)不考慮有關(guān)零件表面形狀誤差影響;

2)不考慮有關(guān)零件同軸度誤差的影響;

3)假設(shè)風(fēng)帽軸線與引信體軸線重合;

4)假設(shè)零件的相關(guān)尺寸均服從正態(tài)分布，取其公差帶中心為散布中心，公差帶寬度取6σ;

5)所有尺寸均以抽樣105次作為整體，以保證所得結(jié)果與實(shí)際相符合。

c)Monte-Carlo 法求解

應(yīng)用Monte-Carlo 方法求解得到的結(jié)果如圖7 所示的頻數(shù)直方圖。

圖7 風(fēng)帽收口前軸向間隙尺寸頻數(shù)分布直方圖

d)計(jì)算結(jié)果分析

通過(guò)抽樣分析可知:

1)風(fēng)帽收口前軸向間隙的最大值為0.989 6 mm，最小值為0.435 5 mm，平均值為0.684 2 mm，其分布曲線基本上呈正態(tài)分布。

2)在極大值極小值法計(jì)算中，最小值為負(fù)值，并且最大值與最小值的差值比較大，而該結(jié)果中的最小值都大于0.4 mm，因此是與實(shí)際結(jié)果是相符合的。

3 結(jié)論

通過(guò)上述計(jì)算可知:

1)Monte-Carlo 方法通過(guò)隨機(jī)抽樣模擬得出的值，可準(zhǔn)確的反映尺寸鏈封閉環(huán)的尺寸分布情況，并得出在尺寸鏈封閉環(huán)的極值。

2)Monte-Carlo 方法可解決尺寸鏈計(jì)算中遇到的復(fù)雜的尺寸鏈計(jì)算問(wèn)題，通過(guò)大量的模擬運(yùn)算，得出問(wèn)題的統(tǒng)計(jì)特征值。

3)Monte-Carlo 方法結(jié)合MATLAB 程序設(shè)計(jì)可計(jì)算組成環(huán)除正態(tài)分布外的其他分布規(guī)律的尺寸鏈，程序簡(jiǎn)單，節(jié)省時(shí)間。

［1］王克沖.用蒙特卡羅法解尺寸鏈問(wèn)題［J］.華東工學(xué)院學(xué)報(bào)，1985，(01):53-59.

［2］魏傳禮，殷仁龍，劉之生.尺寸鏈的計(jì)算及數(shù)值解法［J］.華東工學(xué)院學(xué)報(bào)，1989，(01):52-59.

［3］王磊，丁紅宇，楊東拜，等.基于Monte Carlo 方法求解多環(huán)尺寸鏈的分析研究［J］.制造業(yè)信息化，2002，(01):28-30.

［4］聞泉，王雨時(shí)，陳會(huì)光.引信機(jī)構(gòu)旋轉(zhuǎn)偏心的蒙特卡羅模擬［J］.探測(cè)與控制學(xué)報(bào)，2007，29(3):72-75，79.WEN Quan，WANG Yu-shi，CHEN Hui-guang.Simulation for Spin Eccentricity of Devicein Fuze with Monte Carlo Method［J］.Journal of Detection ＆Control，2007，29(3):72-75，79.

［5］陳會(huì)光.引信球轉(zhuǎn)子機(jī)構(gòu)研究［D］.南京:南京理工大學(xué)，2007.

［6］陳松濤.蒙特卡洛法在尺寸鏈求解中的應(yīng)用研究［C］.江蘇省計(jì)量測(cè)試學(xué)術(shù)論文集，2005:374-379.

［7］引信設(shè)計(jì)手冊(cè)編寫(xiě)組.引信設(shè)計(jì)手冊(cè)［M］.北京:國(guó)防工業(yè)出版社，1987.

［8］李柱主編.互換性與測(cè)量技術(shù)基礎(chǔ)［M］.北京:計(jì)量出版社，1985.

［9］徐鐘濟(jì).蒙特卡羅方法［M］.上海:上海科學(xué)技術(shù)出版社，1985.

［10］朱本仁編.蒙特卡羅方法引論［M］.濟(jì)南:山東大學(xué)出版社，1986.