遠地點發動機推力矢量不確定度評估方法

鄭 科 張 偉

(北京航天試驗技術研究所，北京 100074)

1 引 言

遠地點發動機常被用于星、箭、彈的精確軌道控制和姿態調整，其推力矢量直接關系到衛星能否入軌和發射任務的成敗，準確測出推力矢量參數，能夠為發動機整星組裝時提供基本依據。目前，推力矢量測量有多種方案，在我國的遠地點發動機高空模擬熱標定試驗中，技術和工藝都比較成熟且經多次飛行驗證效果顯著的是北京航天試驗技術研究所1992年研制的動態軸矢量推力架和該所2010年研制的基于壓電晶體測力板技術的多分力動態矢量推力架。

遠地點發動機一般用于衛星等航天器的軌道控制和姿態調整。推力矢量測量的準確性對衛星能否入軌以及入軌成本有很大的影響，甚至影響衛星的工作壽命，所以推力矢量參數是發動機整星組裝的重要參考依據。對于推力矢量的測量，在我國遠地點發動機高空模擬熱標定試驗中，北京航天試驗技術研究所1992年研制的動態軸矢量推力架和該所2010年研制的基于壓電晶體測力板技術的多分力動態矢量推力架經過多次實際飛行驗證，技術相對成熟。

2 動態軸矢量推力架測量模型

測試系統的推力矢量模型如圖1所示，

O

-

XYZ

為轉臺坐標系，

O

’-

X

’

Y

’

Z

’為發動機坐標系，發動機熱試車時產生的發動機推力矢量

F

作用于對接面

Y

’

O

’

Z

’平面的A點上。

圖1 推力矢量模型

推力矢量

F

在推力軸線方向(

X

軸方向)的分量為

F

,其垂直于推力軸線方向的分量為

F

(

Z

’

O

’

Y

’平面內)，在轉臺上產生一個繞

Z

軸的旋轉力矩

M

,經過簡單力學推導可得-

M

=

F

(

δ

cos

β

+

R

tan

α

cos

γ

)

(1)

S

=

δ

cos

β

+

R

tan

α

cos

γ

(2)

當推力矢量裝置使

Y

’

O

’

Z

’平面(即發動機對接面)繞

OX

軸旋轉(角，公式(2)變為

S

=[

δ

cos(

β

+

θ

)+

R

tan

α

cos(

γ

+

θ

)]

(3)

其中

M

，

F

，

R

，

θ

為可測參數，

α

，

β

，

γ

，

δ

為待求參數。分別改變

R

和

θ

的值，令

i

=1,2,3,...,

n

,發動機點火

n

(

n

=6～8)次，求解方程組即可得到矢量參數

S

=[

δ

cos(

β

+

θ

)+

R

tan

α

cos(

γ

+

θ

)]=

δ

cos

β

cos

θ

-

δ

sin

β

sin

θ

+

R

tan

α

cos

γ

cos

θ

-

R

tan

α

sin

γ

sin

θ

(4)

令

(5)

公式(4)可寫為矩陣形式

在這里我們主要對帽殼進行研究分析后，決定采用如圖2亞光表面處理，作為共享產品，使用時間長，接觸人多，頭盔表面難免留有很多劃痕，而亞光表面的日常的劃痕磨損會相對沒有拋光表面易發現，所以，從造型美觀來看，我們更傾向于使用亞光處理表面。同時亞光表面的處理，將會使頭盔的造型感更加沉穩有質感。

AX

=

S

(6)

X

=[

x

x

x

x

]，

若

A

A

的逆矩陣(

A

A

)存在，則得到

X

的最小二乘解為

X

=(

A

A

)

A

S

(7)

推力矢量參數為

(8)

3 GUM不確定度評估

3.1 GUM法簡介

GUM法通過測量模型直接傳遞不確定度，基于方和根法計算合成不確定度。方法源自《測量不確定度表示指南》，假設輸出量近似為正態分布或者縮放位移t分布，通過設定包含概率、計算自由度以及查表求出對稱的包含區間，也叫做不確定度傳播率，尤其適合于測量模型為線性模型，是目前常用的測量不確定度評估方法。

3.2 GUM法評估過程

計算參數合成不確定度一般表達式為

(9)

式中：

u

(

y

)——合成不確定度；

f

——被測量

y

與直接測得量

x

的函數關系；

u

(

x

)——分量

x

的不確定度。

根據公式(8)可以得到

(10)

根據公式(9)計算得到各參數合成標準不確定度表達式為

(11)

B

=(

A

A

)

A

(12)

公式(7)可以寫為

X

=

BS

(13)

可得

(14)

(15)

從表1中計算結果可以得出，

α

合成標準不確定度在0.005°以內；

β

合成標準不確定度基本在10°以內，

γ

合成標準不確定度在10°以內，

δ

合成標準不確定度在0.04mm以內。

3.3 GUM法在該模型評估中的局限性

計算合成標準不確定度的公式(9)忽略了泰勒級數的高價項以及相關性，當模型為復雜非線性模型時，近似為線性模型存在一定的誤差。輸出量的概率密度函數(Probability Density Function, PDF)偏離正態分布或縮放位移

t

分布時，概率分布偏離對稱分布，采用GUM法得出的包含區間存在誤差。由于

R

和

θ

的測量存在誤差，導致

A

自身存在誤差，由最小二乘法解出的矩陣

B

存在誤差，最后通過GUM法得出的不確定度存在一定誤差。

4 MCM不確定度評估方法

4.1 MCM原理

蒙特卡洛法為測量不確定度評估提供了一個通用的計算方式，適合多個輸入量以及單個輸出量的測量模型，輸入量可以為任意的概率分布，相比GUM法，其對測量模型沒有限制，省去了復雜靈敏度的計算。蒙特卡洛法的核心是對輸入量的PDF重復采樣，即利用對輸入量概率分布的隨機抽樣代入測量模型進行分布傳遞，區別于GUM法的傳遞不確定度，其傳遞的是數值特性。計算出輸出量Y的PDF離散采樣值，由于離散采樣值包含輸出量Y的數值特性，因此可以從離散采樣值中得到輸出量的Y的最佳估計值和Y的標準不確定度，隨著對概率密度函數的隨機采樣數M的增加，可提高不確定度、包含區間等數值特性的可靠性。蒙特卡洛法的概率分布傳遞原理圖如圖2所示。

表1 8臺次熱標矢量參數不確定度Tab.1 Uncertaintyofthrustvectorparametersineighttimes試車代號uc(α)/°uc(β)/°uc(γ)/°uc(δ)/mm10.00202.01.70.01820.00201.75.20.01930.00302.13.40.02540.00303.71.20.01850.00200.851.30.01660.00507.51.60.03370.00303.73.90.02980.00100.681.10.0060

圖2 由多輸入量PDF得到單一輸出量PDF的示意圖

4.2 MCM評估過程

(1)定義輸出量

Y

，即推力矢量參數(

δ

,

α

,

β

,

γ

)。(2)確定與

Y

相關的輸入量

X

(

F

)、

X

(

M

)、

X

(

R

)、

X

(

θ

)。(3)建立

Y

和

X

之間的模型

Y

=f(

X

,

X

,

X

,

X

)。(4)利用可獲信息，為

X

設定PDF，如正態分布、均勻分布等。對于

X

和

X

可以采用反復動作調整推力架50次，再直接測量長度(

X

)和角度(

X

)各50次的方法，分別確定

X

和

X

的PDF。對于

X

(

F

)和

X

(

M

)，在確定

X

和

X

的某一特定狀態下，采用現場反復施加靜態力的方法，獲取

X

和

X

的分布樣本各50次，再根據樣本的概率統計評估其PDF。(5)根據蒙特卡洛法的相關國家標準，試驗樣本量的大小

M

，本次初步選取經驗值10個。(6)從輸入量

X

的PDF中隨機抽取

M

個樣本(

X

,

X

,

X

,

X

)，

r

=1,2,...,

M

。(7)對每個樣本向量(

X

,

X

,

X

,

X

)，計算相應

Y

的模型值y=f[(

x

,

x

,

x

,

x

)]。

(8)將計算得到的M個模型值按照遞增順序排序，得出輸出量Y的離散表示G。

(9)由輸出量的離散表示G計算Y的期望值以及標準差，求出估計值y和不確定度μ(y)。

(10)通過輸出量的離散表示G求出設定包含概率下的包含區間。

4.3 MCM的MATLAB實現及結果比較

極限情況下，

α

=0時，

γ

計算結果失真；

δ

=0時，

β

計算結果失真。按照機械加工的計量精測和往次試驗結果統計，實際

α

值大多分布在(0.02～0.2)°以內，合格驗收標準為

α

≤0.2°；

δ

值大多分布在(0.06～0.95)mm以內，合格驗收標準為

δ

≤0.9mm；

β

和

γ

作為

δ

和

α

的定位角，存在批次性在類似象限的散布問題。根據推力矢量計算的數學模型，給定一組真值：

F

=489；(

δ

,

α

,

β

,

γ

)=(0.08,0.05,16,200)；對應

R

=[418

.

9 418

.

9 418

.

9 418

.

9 496

.

72 496

.

72 496

.

72 496

.

72]，

θ

=[0 30 60 90 90 60 30 0]。計算得到

M

=[130

.

37 87

.

73 21

.

58-50

.

36-61

.

71 27

.

34 109

.

07 161

.

58]。通過MATLAB編寫軟件，固化GUM法流程，導入數據，不確定計算結果為：

u

(

δ

)=0.0115mm，

u

(

α

)=0.0014°，

u

(

β

)=4.9699°，

u

(

γ

)=1.0514°，如圖3所示。依據概率統計相關理論，擴展因子

k

=2時，

p

=0.95；

k

=3時，

p

=0.99。故可將當前的不確定評估結果根據需要，加權擴展因子

k

。

圖3 GUM法計算程序

依據

F

、

M

、

R

、

θ

各自的PDF(

F

±0.5%隨機發散、

M

±1%隨機發散、

R

±0.05%隨機發散、

θ

±0.05%隨機發散)，進行MCM不確定度評估。通過MATLAB實現MCM的軟件實現，采用4.2節的評估過程，如圖4所示，針對該真值的不確定度評估結果為：

u

(

δ

)=0.0103mm，

u

(

α

)=0.0013°，

u

(

β

)=2.7089°，

u

(

γ

)=0.5135°，10個樣本散布如圖5所示(橫坐標為數值，縱坐標為出現次數)。相比GUM法不確定度評估，MCM可以計算任何概率下的包含區間。MCM與GUM法結果對比如表2。

圖4 MCM計算程序

圖5 結果散布圖

表2 MCM與GUM不確定度評估對比Tab.2 ComparisonofMCMandGUMuncertaintyevaluation(δ,α,β,γ)評估方法u(δ)/mmu(α)/°u(β)/°u(γ)/°(0.08,0.05,16,200)GUM0.0120.00145.01.1MCM0.0100.00132.70.51

4.3.1 其他數據的對比驗證

不改變

F

，

R

，

θ

，改變(

δ

,

α

,

β

,

γ

)進行驗證，初始值

F

=489，

R

=[418

.

9 418

.

9 418

.

9 418

.

9 496

.

72 496

.

72 496

.

72 496

.

72]，

θ

=[0 30 60 90 90 60 30 0]，計算

M

，結果如表3，不確定度的評估對比如表4。

通過數據對比得出，MCM評估的不確定度數值基本小于GUM法計算的不確定度數值，GUM法相對保守。

5 GUM法與MCM的區別和聯系

從該遠地點發動機推力矢量不確定度評估計算過程來看，GUM法是對多個輸入變量的某一種特定組合進行計算，得到一個輸出結果；而MCM在各個輸入變量PDF的約束下，隨機進行多個輸入變量的

M

次采樣組合，計算得出

M

個輸出結果，再將

M

個輸出結果作為統計樣本，計算樣本的期望和方差，

表3 Mz計算結果Tab.3 CalculationresultofMZ(δ,α,β,γ)Mz/(N·m)(0.06,0.02,16,200)[-38.99-25.58-5.3216.3720.91-7.63-34.12-51.47](0.9,0.02,16,200)[355.86259.7694.05-96.85-92.3191.75251.22343.38](0.06,0.2,16,200)[-643.71-439.24-117.07236.47281.90-140.13-524.62-768.53](0.9,0.2,16,200)[-248.86-153.90-17.70123.25168.68-40.76-239.28-373.69]

表4 MCM與GUM不確定度評估對比Tab.4 ComparisonofMCMandGUMuncertaintyevaluation(δ,α,β,γ)評估方法u(δ)/mmu(α)/°u(β)/°u(γ)/°(0.06,0.02,16,200)GUM0.00360.000402.10.83MCM0.00320.000401.10.40(0.9,0.02,16,200)GUM0.0280.00341.05.9MCM0.0250.00320.522.9(0.06,0.2,16,200)GUM0.0560.0069321.3MCM0.0410.0064530.61(0.9,0.2,16,200)GUM0.0250.00311.00.60MCM0.0220.00290.550.30