彈丸空氣阻力定律全定義域解析函數(shù)經(jīng)驗(yàn)公式*

倪慶樂(lè),王雨時(shí),聞 泉,張志彪

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

彈丸空氣阻力定律全定義域解析函數(shù)經(jīng)驗(yàn)公式*

倪慶樂(lè),王雨時(shí),聞 泉,張志彪

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

為找出1943年彈丸空氣阻力定律、航彈阻力定律在全定義域內(nèi)的有理解析函數(shù)經(jīng)驗(yàn)公式,從而便于通過(guò)計(jì)算機(jī)編程進(jìn)行外彈道數(shù)值解算和求解析解,運(yùn)用1stOpt數(shù)學(xué)軟件對(duì)1943年彈丸空氣阻力定律和航彈阻力定律在其全定義域內(nèi)進(jìn)行了曲線擬合,對(duì)應(yīng)各阻力定律得到了通用于全定義域的有理解析函數(shù)經(jīng)驗(yàn)公式。這些經(jīng)驗(yàn)公式均不再分段。其中擬合1943年彈丸空氣阻力定律的經(jīng)驗(yàn)公式最大相對(duì)誤差為3.271%,擬合高阻和低阻航彈阻力定律的分別為1.914%和3.328%。

外彈道學(xué);經(jīng)驗(yàn)公式;空氣阻力;阻力定律;阻力系數(shù)

0 引言

在外彈道計(jì)算中,空氣阻力系數(shù)是重要的參數(shù)。目前工程設(shè)計(jì)中的外彈道解算多數(shù)仍采用彈丸彈形系數(shù)概念和彈丸空氣阻力定律,并且其精度一般均能滿足工程設(shè)計(jì)要求[1]。目前身管武器彈丸常用的阻力定律主要是1943年阻力定律,但1943年彈丸空氣阻力定律是以表值形式給出的。為了便于空氣外彈道解算計(jì)算機(jī)程序編制,并作進(jìn)一步的解析分析,需要對(duì)該表值定律進(jìn)行曲線擬合,轉(zhuǎn)換成解析形式。文獻(xiàn)[2]給出的1943年彈丸空氣阻力定律的阻力函數(shù)經(jīng)驗(yàn)公式,以250 m/s、400 m/s、1 400 m/s為速度節(jié)點(diǎn)將阻力函數(shù)曲線分為4段,同樣應(yīng)用不方便。文獻(xiàn)[3]曾用三次拋物線分段擬合方法給出了1943年阻力定律的解析函數(shù)表達(dá)式,其誤差雖然不大于4%,但分段過(guò)多,在外彈道編程應(yīng)用中不方便。文獻(xiàn)[4]對(duì)文獻(xiàn)[3]的曲線分段擬合方法進(jìn)行了改進(jìn),以跨音速段的最大值為界分亞音速段和超音速段兩段,利用Logistic曲線進(jìn)行擬合處理,給出了1943年阻力定律的分段解析函數(shù)表達(dá)式,其最大誤差為5.05%。文獻(xiàn)[5]采用三次B樣條曲線對(duì)1943年阻力定律進(jìn)行擬合,最大誤差為2.2%,但由于是根據(jù)遞推定義所做的曲線,對(duì)曲線中間點(diǎn)的值只能通過(guò)前面點(diǎn)推導(dǎo)得出,不能由一個(gè)解析式表達(dá)出整條曲線,過(guò)程繁瑣,不便于推廣使用。航彈阻力定律也有與此類似的問(wèn)題。文中為找出1943年阻力定律、航彈阻力定律在全定義域內(nèi)的直觀的有理解析函數(shù)經(jīng)驗(yàn)公式,擬運(yùn)用1stOpt軟件對(duì)其進(jìn)行擬合。

1 曲線回歸理論基礎(chǔ)與1stOpt軟件簡(jiǎn)介

1.1 曲線回歸理論基礎(chǔ)

曲線回歸是以最小二乘法分析變量間在數(shù)量變化上的特征和規(guī)律的方法,主要用來(lái)確定兩個(gè)變量間數(shù)量變化的某種特定的規(guī)律,估計(jì)表示該種曲線關(guān)系特點(diǎn)的一些重要參數(shù)(極大值、極小值、漸近值等),進(jìn)行數(shù)值內(nèi)插或理論上的外推。

最小二乘法通過(guò)最小化誤差的平方和來(lái)尋找數(shù)據(jù)的最佳函數(shù)匹配,將實(shí)測(cè)值yi與擬合函數(shù)對(duì)應(yīng)值yj離差的平方和∑(yi-yj)2最小作為優(yōu)化判據(jù)。

相關(guān)系數(shù)是用以反映變量之間關(guān)系密切程度的統(tǒng)計(jì)指標(biāo),是按積差方法計(jì)算,以兩變量與各自平均值的離差為基礎(chǔ),通過(guò)兩個(gè)離差相乘來(lái)反映兩變量之間的相關(guān)程度。相關(guān)系數(shù)恒小于1,其值越接近于1表明兩者相關(guān)程度越高。

方差用來(lái)度量隨機(jī)變量與其數(shù)學(xué)期望之間的偏離程度,也可以反映出隨機(jī)變量或一組數(shù)據(jù)的離散程度,方差越大,離散程度越大。設(shè)x是一個(gè)隨機(jī)變量,若E{[x-E(x)]2}存在,則稱E{[x-E(x)]2}為x的方差,其中E(x)表示x的數(shù)學(xué)期望。

1.2 1stOpt軟件簡(jiǎn)介

1stOpt是我國(guó)七維高科有限公司獨(dú)立開(kāi)發(fā)的一套數(shù)學(xué)優(yōu)化分析綜合工具軟件包,主要用于曲線擬合、非線性復(fù)雜模型參數(shù)估算求解、線性/非線性規(guī)劃等領(lǐng)域。1stOpt軟件不僅操作界面簡(jiǎn)單,而且其核心的通用全局優(yōu)化算法克服了當(dāng)今世界上優(yōu)化計(jì)算領(lǐng)域使用迭代法必須給出合適初始值的難題,即不需要用戶給出參數(shù)初始值,在絕大多數(shù)情況下,僅依靠自身的全局搜索能力,從任意隨機(jī)值出發(fā),即可求得最優(yōu)解。

1stOpt擁有強(qiáng)大的非線性回歸功能,可自動(dòng)搜索最匹配的模型公式。對(duì)常用的兩變量及三變量問(wèn)題,1stOpt模型公式分別預(yù)置有超過(guò)3 700和近1 000個(gè)函數(shù)式,且可根據(jù)需要自行隨意添加修改。

2 1943年彈丸空氣阻力定律的擬合

2.1 1stOpt軟件自動(dòng)擬合函數(shù)

文獻(xiàn)[6]給出了1943年彈丸空氣阻力定律表值函數(shù),包含了Ma從0～4.9之間的共169組數(shù)據(jù)(Ma在0.7～2.1之間間隔0.01,Ma在2.1～4.9之間間隔0.1)。在1stOpt軟件代碼框內(nèi)輸入相應(yīng)的命令代碼和1943年阻力定律的表值數(shù)據(jù),設(shè)置好算法后進(jìn)行擬合,擬合過(guò)程采用的算法為麥夸特法+通用全局優(yōu)化法。1stOpt軟件會(huì)自動(dòng)將數(shù)據(jù)和預(yù)置函數(shù)進(jìn)行匹配擬合,并按相關(guān)系數(shù)由大到小依次列出。

擬合所得的相關(guān)系數(shù)最高的曲線方程為:

(1)

式中:y為空氣阻力系數(shù);x為馬赫數(shù);p1、p2、p3、p4、p5、p6、p7、p8、p9為待定系數(shù)。

擬合所得的待定系數(shù)共有15位有效數(shù)字,分別是:

p1=0.154 944 442 142 445

p2=-2.645 771 025 713 59

p3=-0.389 894 040 604 623

p4=2.742 885 185 889 92

p5=0.381 002 068 711 428

p6=-1.499 032 921 978 55

p7=-0.237 775 531 351 884

p8=0.471 410 603 010 32

p9=0.114 120 378 908 769

其最大相對(duì)誤差為3.355%,方差為2.047 8×10-3,最大相對(duì)誤差點(diǎn)在Ma=0.91處。擬合曲線如圖1所示。

圖1 1943年阻力定律擬合曲線

2.2 數(shù)據(jù)均勻化處理

在2.1節(jié)擬合所依據(jù)的數(shù)據(jù)中,Ma在0.7～2.1之間時(shí)間隔0.01,在2.1～4.9之間時(shí)間隔0.1,這就造成了2.1～4.9之間的數(shù)據(jù)在擬合中的權(quán)重偏小。為了能使各數(shù)值點(diǎn)在擬合中的權(quán)重相等,對(duì)2.1～4.9之間的數(shù)據(jù)進(jìn)行線性插值,使其間隔為0.01,得到共490組數(shù)據(jù)。利用1stOpt軟件對(duì)490組數(shù)據(jù)進(jìn)行擬合,得到相關(guān)系數(shù)最大的方程形式仍同式(1),但擬合所得的待定系數(shù)值與式(1)的略有不同:

p1=0.157 182 242 292 442

p2=-2.697 961 914 717 5

p3=-0.422 569 444 115 854

p4=2.834 248 695 049 56

p5=0.448 824 524 761 816

p6=-1.495 221 097 605 56

p7=-0.266 930 620 397 897

p8=0.412 922 615 928 505

p9=0.100 836 885 531 747

其最大誤差為3.906%,方差為2.242×10-3,最大相對(duì)誤差點(diǎn)在Ma=0.91處。與2.1節(jié)擬合所得結(jié)果比較可知:對(duì)數(shù)據(jù)進(jìn)行均勻化處理后擬合得到的最大誤差和方差與未進(jìn)行處理的結(jié)果基本相同,因此推薦使用2.1節(jié)所得結(jié)果,在后續(xù)擬合中不再進(jìn)行數(shù)據(jù)均勻化處理。

2.3 待定系數(shù)有效數(shù)字位數(shù)的選取

由上面所得結(jié)果可以看出,擬合所得待定系數(shù)的位數(shù)較多,不便于計(jì)算和應(yīng)用,因此需要討論待定系數(shù)的位數(shù)對(duì)擬合曲線精度的影響。

不同有效數(shù)字位數(shù)的待定系數(shù)擬合效果的比較列于表1中。從表1可以看出,當(dāng)待定系數(shù)有效數(shù)字位數(shù)為3～5位時(shí),其擬合效果與最優(yōu)解擬合效果基本上接近,所以在應(yīng)用中可以將待定系數(shù)取3位有效數(shù)字。

表1 不同有效數(shù)字位數(shù)的待定系數(shù)擬合效果的比較

2.4 待定系數(shù)的優(yōu)化

將所有待定系數(shù)按四舍五入法則處理后取5位有效數(shù)字所得的擬合曲線與最優(yōu)解擬合曲線雖然很接近,但是仍存在一定的偏差,為此,可依次將其中一個(gè)待定系數(shù)作為待優(yōu)化項(xiàng),將其余待定系數(shù)賦予5位有效數(shù)字的近似值,確定新的擬合方程進(jìn)行二次擬合,依此來(lái)對(duì)待定系數(shù)進(jìn)行優(yōu)化。

從表2可以看出,對(duì)待定系數(shù)進(jìn)行二次擬合與將

表2 二次擬合與最優(yōu)解取5位有效數(shù)字的擬合效果比較

最優(yōu)解取3位有效數(shù)字的擬合效果基本相同,因此,1943年彈丸空氣阻力定律的有理式擬合結(jié)果為:

(2)

或

(3)

3 航彈阻力定律的曲線擬合

航彈阻力定律分為高阻航彈阻力定律和低阻航彈阻力定律。文獻(xiàn)[7]給出了高阻航彈阻力定律、低阻航彈阻力定律的表值函數(shù),包含了Ma從0～1.5之間的各20組數(shù)據(jù)。利用相同方法對(duì)高阻彈、低阻彈阻力定律進(jìn)行曲線擬合。擬合所得的相關(guān)系數(shù)最高的高阻彈阻力定律擬合方程為:y=p1+p2x2+p3x4+p5x8+p6x10+p7x12+

p8x14+p9x16+p10x18+p11x20

(4)

(5)

高阻彈、低阻彈阻力定律擬合方程待定系數(shù)取不同位數(shù)時(shí)擬合效果的比較如表3所列。

由表3可看出,利用此方法可以很好的對(duì)高阻航彈阻力定律、低阻航彈阻力定律進(jìn)行曲線擬合,擬合結(jié)果為擬合方程中的待定系數(shù)均取5位有效數(shù)字,即對(duì)于高阻航彈:

y=0.581 98+0.397 84x2-5.088 6x4- 53.419x8+49.707x10-0.020 752x12- 37.437x14+30.065x16-9.889 8x18+ 1.216 4x20

(6)

而對(duì)于低阻航彈:

(7)

高阻航彈阻力定律、低阻航彈阻力定律的擬合曲線分別如圖2和圖3所示。其最大相對(duì)誤差分別在Ma=0.9和Ma=0.8處。

圖2 高阻航彈阻力定律擬合曲線

圖3 低阻航彈阻力定律擬合曲線

4 結(jié)束語(yǔ)

文中利用1stOpt軟件對(duì)1943年阻力定律、高阻航彈阻力定律、低阻航彈阻力定律進(jìn)行了曲線擬合,結(jié)果表明:該方法可以很好的擬合各阻力定律,可表

示為一個(gè)全定義域內(nèi)的通用解析函數(shù)經(jīng)驗(yàn)公式,不再是分段函數(shù),便于計(jì)算程序編制和做進(jìn)一步解析求解;該方法擬合精度也更高,對(duì)1943年阻力定律擬合的最大誤差為3.271%,低于兩段Logistic曲線擬合時(shí)的最大相對(duì)誤差5.05%,對(duì)高阻航彈阻力定律擬合的最大誤差為1.914%,對(duì)低阻航彈阻力定律擬合的最大誤差為3.328%,均可以滿足彈道解算精度要求。而對(duì)于通過(guò)風(fēng)洞試驗(yàn)和計(jì)算流體力學(xué)仿真所得的彈丸阻力系數(shù)與馬赫數(shù)之間的對(duì)應(yīng)數(shù)值關(guān)系,也可按此方法給出類似的全定義域內(nèi)的解析函數(shù)經(jīng)驗(yàn)公式。

[1] 楊翔, 王雨時(shí), 聞泉. 應(yīng)用阻力系數(shù)擬合曲線解析式數(shù)值解算外彈道諸元 [J]. 彈箭與制導(dǎo)學(xué)報(bào), 2014, 34(5): 151-155.

[2] 浦發(fā), 芮筱亭. 外彈道學(xué) [M]. 修訂本. 北京: 國(guó)防工業(yè)出版社, 1989: 39-40.

[3] 王雨時(shí). 外彈道學(xué)阻力定律的三次拋物線分段擬合 [C]∥遼寧省兵工學(xué)會(huì)第5屆學(xué)術(shù)年會(huì)論文, 1994.

[4] 王雨時(shí). 彈丸戰(zhàn)斗部及其破片空氣阻力系數(shù)的Logistic曲線分段擬合 [J]. 彈箭與制導(dǎo)學(xué)報(bào), 2006, 26(1): 242-244.

[5] 趙靜, 杜忠華, 趙永平, 等. 基于三次B樣條曲線的一維彈道修正彈空氣阻力系數(shù)擬合 [J]. 火力與指揮控制, 2015, 40(4): 123-126.

[6] 徐明友. 火箭外彈道學(xué) [M]. 北京: 兵器工業(yè)出版社, 1989: 352.

[7] 浦發(fā), 薛曉中, 程豫生. 航空炸彈彈道學(xué) [M]. 北京: 中國(guó)人民解放軍空軍后勤部軍械部, 1986: 41.

Empirical Formulas of Projectile Air Resistance Law in Whole Definition Domain in Analytic Function

NI Qingle,WANG Yushi,WEN Quan,ZHANG Zhibiao

(School of Mechanical Engineering, Nanjing University of Science and Technology, Nanjing 210094, China)

In order to find out rational analytic function empirical formulas of projectile air resistance law and aerial bombs resistance law in 1943 in the whole definition domain to proceed numerical calculation of exterio trajectory and analytical solution by computer programming, using 1stOp mathematical software was used to fit curves of projectile air resistance law in 1943 and aerial bombs resistance law in the whole definition domain. Rational analytic function empirical formulas common to the whole definition domain was got corresponding to the resistance laws. These empirical formulas are no longer segmented. The maximum relative error of empirical formula fitting projectile air resistance law in 1943 was 3.271%, and the maximum relative error of empirical formula fitting high resistance aerial bomb and low resistance aerial bomb air resistance law were 1.914% and 3.328% respectively.

exterior ballistics; empirical formula; air resistance; resistance law; resistance coefficient

2016-03-08

江蘇省2015年度普通高校研究生科研創(chuàng)新(實(shí)踐)計(jì)劃項(xiàng)目(KYLX15_0334)資助

倪慶樂(lè)(1991-),男,河北衡水人,碩士研究生,研究方向:引信設(shè)計(jì)及其動(dòng)態(tài)特性的研究。

TJ012.3

A