基于射擊誤差狀態方程的高炮毀傷概率建模

王向民，王 軍，謝杰濤，郭 治

(1.南京理工大學 自動化學院，南京 210094;2.中國人民解放軍 32200 部隊，遼寧 錦州 121000)

為應對來自空中的威脅，防空反導系統作為最直接有效的反制手段，受到了各國的高度重視.以高炮為主的近程末端防空武器系統，其具有較高的費效比，常作為防御體系的最后屏障，因而它對來襲目標的毀傷關系到整個防空作戰的成敗.為獲取高炮武器系統的毀傷概率：一種方法是通過大量的實彈射擊，統計出射彈數目和實際毀傷目標的數量，二者比值即為武器的毀傷概率,該方法可信度高，但檢驗費用昂貴，且不能適用于武器系統的設計與論證階段；另一種方法是建立基于高炮射擊誤差的毀傷概率數學模型，利用數值計算的方法獲得高炮武器系統的毀傷概率，該方法通過對高炮的射擊誤差進行適當分解，利用誤差模型轉換法構建毀傷概率的計算模型.在我國現行的國家軍用標準(簡稱國軍標)中，把誤差轉換為重復誤差和不可重復誤差，雖然降低了計算復雜度，但存在較大的模型近似[1-4].射擊誤差中的弱相關誤差源對于毀傷概率計算具有顯著影響，陶德敬等[5]提出將射擊誤差分解為共有分量和非共有分量，運用誤差序列預測理論給出了確定強相關誤差和弱相關誤差下射擊過程的平均毀傷概率表達式，但該模型只能滿足誤差序列為2階狀態方程以下的情況；姚志軍等[6]利用弱相關射擊誤差狀態方程給出了遞推形式的毀傷概率計算方法，但同樣對于高階自相關誤差序列的相關性計算沒有明確給出解析式.隨著現代高炮武器系統射頻的提高，射擊過程中前一發彈藥的毀傷概率對其后各發都有影響，因而每次發射的彈藥對毀傷概率的貢獻都是有差別的，而將誤差序列簡單地近似為1階或2階系統的模型簡化方法[7]，有悖于高炮射擊的真實模型.高炮的誤差序列從本質上來說，可以視為武器身管的隨機振動、跟蹤誤差以及火控解算誤差等綜合而成的平穩時間序列，有關時間序列相關性的分析方法，包括線性和非線性回歸方法進行預測[8]、基于數據驅動的時間序列預測等[9-10]，需要有足夠的實驗數據才能獲得較高的預測精度[11-12]，并且隨著時間序列相關性階數的增加，建模的復雜性和計算量也會大幅增加[13].

為了建立更精確的毀傷概率計算模型，本文以射擊誤差序列所構建的狀態方程為出發點，運用卡爾曼濾波的思想，以遞推方式給出連續脫靶條件下射擊諸元誤差的密度函數，將包括由射擊沖擊載荷影響在內的射擊誤差的弱相關分量分解為完全可預測與完全不可預測的分量，再分別與完全可預測的強相關分量和完全不可預測的不相關分量合并，構成2個隨機分量，根據高炮進行點射時的時空特性，構建射擊誤差的隨機狀態模型，進而推導出相應的毀傷概率遞推計算模型.

1 射擊誤差序列狀態方程構建

x(k)=xb(k)+xr(k)+xq(k)

(1)

式中：xr(k)與xq(k)之和稱為射擊諸元誤差xs(k).射擊誤差序列的狀態方程可進一步表示為

(2)

作為隨機狀態方程，還要已知它的初始狀態方差才能是完備的.設它的初始狀態方差為

(3)

考慮到射擊誤差為平穩高斯過程，因而可得[5]

(4)

式中：Rn為n階預測系數；w(k)為白噪聲.此時隨機誤差序列可分解成可預測分量與不可預測分量.

2 射擊誤差序列狀態方程構建

射擊誤差序列的一步預測特性是影響毀傷概率的一個決定性因素.對一個正常的射擊過程而言，其射擊誤差是由兩部分構成：射擊諸元誤差，它是經歷了足夠長的時間后，暫態分量可以忽略了的相關平穩序列；射彈散布誤差，它雖然也是平穩序列，但是它是在擊發后才加入到射擊誤差序列之中，并且在有限次數的發射后，即行終止的平穩序列.

2.1 射擊誤差序列的統計特性

對射擊誤差而言，它的不相關、弱相關和強相關3個平穩隨機序列中的每一個都是由更多的獨立的子序列之和構成.對誤差序列的一步預測方程,可表示為

(5)

(6)

并且有

(7)

(8)

兩個獨立序列和的均值函數是其k步預測值，而協方差函數是其k步預測方差.兩個隨機序列之和存在一個過渡過程.當0<|r|<1且k→∞時，此過渡過程才能結束，兩個獨立序列的和的均值函數與方差函數才能如同隨機變量一樣求和與方差.顯然，只有兩種情況是例外，即兩個獨立序列均為強相關序列或均為不相關序列，這是因為在這兩種情況下，其任何一步的預測均值與方差都是不變的，且與其初值相同.

如果射擊誤差序列中的弱相關分量的過渡過程可以在射擊準備時間內衰減到可以忽略的條件，射擊過程中的所有弱相關分量與強相關分量才可以依隨機變量求和規則予以處理.例如：射擊諸元誤差xs(k)有3個獨立序列構成，即

xs(k)=xg(k)+xw(k)+xq(k)

(9)

(10)

且有xs(k)的一步預測系數表達式為

(11)

2.2 連續脫靶條件下射擊諸元誤差的密度函數

f(xs(M),xs(M-1),…,xs(n),…,xs(1))=

(12)

Xs(k-1)=

[xs(k-1)xs(k-2) …xs(k-n)]T

(13)

若射擊諸元誤差的最高階次為n，根據式(12)給出的處理方法，xs(k)在當k>n時，則預測樣本空間為

{xs(k)|X(k-1)}∈

(14)

此時xs(k)為一種雙重正態分布，它作為射擊誤差的動態均值，與前面的各個動態均值都有關.

由于k∈[1,n]期間，高炮武器控制系統完成了振動的過渡過程，當k>n，xs(k)進入平穩序列，表明高炮隨時可以實施射擊.

對首次擊發的第一枚彈藥，由于射彈散布誤差xb(1)使x(1)=xb(1)+xs(1)的方差產生突變，此時，x(1)的條件密度函數可表述為

f(x(1))=f(x(1)|xs(1)f(xs(1)))

(15)

因而在首發彈藥脫靶條件下，射擊諸元誤差的密度函數為

g1(xs(1))=

(16)

式中:s1為當k=1時刻靶標的迎彈面，其補集為R1.因此首發彈藥脫靶的概率為

(17)

根據條件概率密度的定義，連續兩次發射均脫靶條件下，射擊諸元誤差的密度函數為

g1(xs(1))dxs(1)

(18)

而連續兩發均脫靶的概率為

P(x(1)?s1,x(2)?s2)=

(19)

依據上述推導方式，可以得到前N發彈藥均脫靶條件下，射擊諸元誤差的密度函數為

gN-1(xs(N-1))dxs(N-1)=

gN-1[xs(N-1)]dxs(N-1)

(20)

上述密度函數的積分即為連續發射N發彈藥均脫靶的概率.對于滿足0-1毀傷律的高炮武器系統，在連續發射N發彈藥的條件下，至少命中一發的概率可表示為

(21)

如果將毀傷定義為命中，則式(21)為射擊諸元誤差存在相關性時，高炮武器進行N發連續射擊的毀傷概率.

3 單管高炮武器系統的毀傷概率

單管高炮武器系統是指在同一瞬時僅能發射一發彈藥.對于轉管火炮，雖然在外形上具有多個發射管，但每個火炮身管只能在轉到特定位置后，依單管規則順序發射，因此在計算毀傷概率時，仍然可按照單管方式加以計算.

(22)

再考慮到射彈散布誤差序列是在開始時與射擊準備誤差序列求和的，故必須以已知射擊諸元Zs(k)條件下的分布特性來表述.基于上述分析可知第一發彈藥射擊諸元誤差的密度函數為

(23)

從式(23)出發，可遞推出前k發均不毀傷條件下，射擊諸元誤差的密度函數為

gk(Zs(k))=

gk-1(Zs(k-1))dZs(k-1)

(24)

則單管高炮武器一次N發的點射的毀傷概率為

(25)

式中：Σb和Σs為協方差矩陣；det為行列式；I為單位矩陣.

4 數值仿真與分析

當射擊誤差中的弱相關分量為1階平穩序列模型時，分別采用國軍標、Monte Carlo(MC)模擬、文獻[6]遞推法以及本文方法進行計算，仿真計算結果如圖1所示.當射擊誤差中的弱相關分量采用3階平穩序列模型時，誤差序列的各階自相關系數為：ri,i=0.81,ri,i+1=0.72,ri,i+2=0.66 (i=1,2,…,21)，仿真計算結果如表1所示.當射擊誤差中的弱相關分量采用4階平穩序列模型時，誤差序列的各階自相關系數為：ri,i=0.81,ri,i+1=0.72,ri,i+2=0.66，ri,i+3=0.60 (i=1,2,…,20)，仿真計算結果如表2所示.

圖1 相關系數與點射毀傷概率關系對比圖Fig.1 Comparison of bust firing damage probability with correlation coefficient

表1 具有3階自相關性的射擊誤差序列毀傷概率Tab.1 Damage probabilities of firing error sequence with third-order autocorrelation

表2 具有4階自相關性的射擊誤差序列毀傷概率Tab.2 Damage probabilities of firing error sequence with fourth-order autocorrelation

上述仿真結果表明，當射擊誤差序列考慮為1階平穩序列模型時，4種計算方法的結果毀傷概率結果差異性很小.但隨著誤差相關性階數的增加，國軍標的計算結果與其他3種存在較為明顯的差距，這表明將高階自相關性誤差序列模型近似為1階模型，會產生一定的模型計算誤差.而運用本文給出的預測系數計算方法可以較好地適應各類高階模型，計算結果更加符合射擊誤差的實際情況.由于仿真中的各弱相關誤差序列的隨機特性已經給定，采用MC方法模擬射擊過程中的各弱相關誤差源，可以忽略模型誤差，通過大量的模擬抽樣所得到的數值計算結果是可信的.因此可以將MC方法的計算結果作為比較基準進行準確性驗證，通過表1和表2計算結果對比可看出，國軍標方法有20%以上的偏差，文獻[6]方法的偏差為10%以上，而本文計算結果的偏差小于3%.

對于計算過于復雜而難以求得解析解的隨機過程問題，雖然MC模擬可以通過構造符合一定規則的隨機數來進行數值求解，但正如前文所述，考慮到高炮武器系統驗證的費效比問題，高炮武器系統從設計到部隊列裝，不可能通過大量的試驗進行毀傷概率的驗證，在試驗數據較少的情況下，使用MC模擬射擊打靶過程，存在模型誤差的問題，有可能產生較大的計算誤差.為進一步說明此問題，利用GPML Toolbox version 4.1工具箱在MATLAB R2013a下，生成一個均值為0，方差為1的高斯隨機過程，并產生 10 000 組樣本，每個樣本均勻采樣20個點，統計采樣點落在某一區間的占比，用于模擬高炮打靶發射彈藥數與命中目標數的比例，即為命中概率.根據MC模擬的原則，分別抽取不同數量的樣本m1，每組總采樣點個數為20m1，記錄采樣點落在[-0.5，0.5]之間的個數為n1，記p1=n1/(20m1)，稱為抽樣落點占比.而總的采樣點個數為 200 000 個，統計得到全體采樣點落在[-0.5，0.5]之間的個數為 93 662 個，記p2=93 662/200 000=0.468，稱為全體落點占比.不同采樣數量下的落點占比與全體落點占比的統計結果如圖2所示.

圖2 MC算法抽樣數量與落點占比的關系Fig.2 Number of samples versus proportion of dropping points in MC algorithm

從圖2中可知，當抽樣數量較少時，MC模擬的結果具有較大的不確定性.模擬過程中還發現，即使抽樣數量相同的不同樣本，MC模擬的數值統計結果也存在明顯差異，這對于試驗數據較少，且要求毀傷概率計算具有嚴格數學推理的武器研制與驗證單位是不能接受的.而本文的遞推估計模型具有嚴格的數學推導，只需確定弱相關誤差的均值、方差和相關系數，就可以準確計算出毀傷概率.

5 結語

提出了一種利用高炮武器射擊誤差的狀態方程估計毀傷概率的方法，在將射擊誤差序列中的弱相關部分，分解為可預測誤差和不可預測誤差的，給出了可預測系數的計算方法，并利用遞推方法建立了高炮武器連續點射的毀傷概率計算模型，可以較好地適用于射擊誤差序列具有多階自相關性時高炮武器系統的毀傷概率計算.本方法也可以作為各類速射火炮進行連續射擊時毀傷概率計算的一種參考.