程胡華，王益柏，蔡其發(fā)，趙 亮

(1. 63729部隊，太原030027；2. 61741部隊，北京 100094；3. 中國科學院大氣物理研究所大氣科學和地球流體力學數(shù)值模擬國家重點實驗室，北京 100029)

0 引 言

近年來，中國在航天領(lǐng)域取得一系列豐碩成果，運載火箭作為航天技術(shù)領(lǐng)域里的核心之一，其內(nèi)部結(jié)構(gòu)非常復雜，并匯聚了無數(shù)尖端技術(shù)，為持續(xù)提高其完成任務的成功率，大量科學家對運載火箭結(jié)構(gòu)特征、飛行控制等方面開展了一系列研究[1-4]，并取得了豐富的成果；在發(fā)射之前，除確保運載火箭本身滿足發(fā)射、飛行條件外，還需考慮大氣環(huán)境的影響，目前主要關(guān)注高空風對運載火箭飛行的影響[5-13]，例如，程鎮(zhèn)煌[5]通過全面敘述大型火箭風載試驗的各種問題，指出火箭設計必須經(jīng)過靶場風載試驗和高空風載試驗，測出火箭各艙段載荷和儀器艙處風激振幅等，準確確定各種設備的適應能力；余夢倫[7]指出發(fā)現(xiàn)高空風引起的氣流攻角對火箭飛行中的氣動載荷和飛行環(huán)境有較大的影響，這種影響有時可造成火箭飛行失敗。

目前，在定量分析大氣環(huán)境對運載火箭飛行影響(即qαmax)的計算程序中，除了需要大氣高空風參數(shù)外，還需要大氣密度參數(shù)，但程序中的大氣密度參數(shù)通常采用參考大氣密度，而非發(fā)射場零日的真實大氣密度，針對大氣密度對運載火箭或其它飛行器飛行的影響，廣大科學工作者已開展了一系列工作，并取得大量有意義成果[14-23]，黃世勇和聞悅[14]指出大氣密度對飛行器的安全及返回有非常重要影響，李健等[20]研究發(fā)現(xiàn)大氣密度是建立飛行動力學模型的重要參數(shù)，不同的大氣密度模型對彈道終點高度、最大過載、最大熱流等參數(shù)影響明顯；程旋等[21]指出飛行器的氣動力和氣動力矩都是大氣密度的函數(shù)，大氣密度的變化直接影響到飛行狀態(tài)下飛行器的受力，對飛行器的姿態(tài)影響顯著，引起攻角、側(cè)滑角、偏航角等角度的變化，導致再入飛行器偏離預定軌跡；因此，大量科技工作者針對大氣密度開展一系列研究工作[24-27]。

獲取參考大氣密度的參考大氣模型有CIRA系列、MSIS系列等，由于大氣密度受來自太陽輻射及大氣波動的共同影響，變化十分復雜，人們至今尚未完全掌握其變化的物理機制，導致已有的參考大氣模型精度不高[19,24,28]，因此，在采用參考大氣密度定量計算大氣環(huán)境對運載火箭飛行的qαmax時，由于參考大氣密度與真實大氣密度之間存在偏差，但針對其偏差對運載火箭飛行的qαmax精度特征研究尚未見相關(guān)文獻，考慮到目前計算程序仍采用參考大氣密度，本文利用2017年8月1日-2019年7月31日晉西北地區(qū)探空資料，分別利用MSIS00模式參考大氣密度、真實大氣密度計算某型運載火箭飛行時的qαmax，一方面對MSIS00模式參考大氣密度與真實大氣密度之間的差異特征進行研究，另一方面對它們之間的qαmax差異進行統(tǒng)計，最后，利用多元線性回歸、BP神經(jīng)網(wǎng)絡方法對基于MSIS00模式參考大氣密度得到的qαmax建立訂正模型，并對訂正效果進行檢驗。

1 資料及方法

1.1 探空資料

探空資料為晉西北地區(qū)每日08時(北京時，下同)，時間長度為2017年8月1日-2019年7月31日，所需氣象要素為高度、氣壓、風速、風向、溫度和相對濕度，單位分別為m，hPa，m·s-1，(°)，℃和%，考慮到qαmax值常出現(xiàn)在16 km高度以下，故選取探空資料的幾何高度范圍為0～20000 m，其中，0～10000 m高度范圍內(nèi)間隔為200 m，10000～20000 m高度范圍內(nèi)間隔為500 m。

1.2 計算真實大氣密度

利用高垂直分辨率探空資料計算真實大氣密度，對于濕空氣，利用溫度、相對濕度、氣壓可計算大氣密度，計算步驟如下[29]：

1)計算飽和水汽壓

(1)

式中：T為溫度。

2)計算水汽壓

e=10-2RH·es

(2)

式中：RH為相對濕度。

3)計算比濕

(3)

式中：p為氣壓。

4)計算濕空氣的比氣體常數(shù)

Rv=Rd(1+0.608q)

(4)

式中：Rd為干空氣的比氣體常數(shù)，數(shù)值為287.05，單位為J/(kg·K)。

5)利用濕空氣的狀態(tài)方程，可計算濕空氣大氣密度

(5)

1.3 MSIS參考大氣模式

MSIS(Mass Spectrometer and Incoherent Scatter empirical model)大氣模式最初版本是由美國國家航空航天局空間飛行中心編制的經(jīng)驗大氣預報模式，通過利用多顆衛(wèi)星質(zhì)譜計所測得的成分資料和非相干散射雷達探測資料而建立，目前MSIS00是該系列模式的較新版本[30]，其利用了更新、更多的雷達、衛(wèi)星和火箭資料，是當今國際上廣泛使用的參考大氣模式[31-35]，MSIS00模式的適用高度范圍為地面至高層大氣，可給出不同高度層的大氣溫度、大氣密度及成分結(jié)構(gòu)。

1.4 統(tǒng)計方法

以真實大氣密度為基準，通過計算參考大氣密度與真實大氣密度之間的偏差、相對偏差和偏差區(qū)間占有率，反映參考大氣密度誤差特征，其相應的計算公式如下：

(6)

(7)

式中：Bias為偏差，Rerr為相對偏差，n為樣本數(shù)，x和y分別代表參考大氣密度與真實大氣密度。

文中的偏差區(qū)間占有率計算式定義如下：

(8)

式中：若xi-yi<-N結(jié)果成立，則定義xi-yi<-N的值為1，否則值為0;該定義同樣適用于表達式-N≤xi-yi≤N和xi-yi>N。

2 密度差異特征

2.1 平均差異特征

在2017年8月1日-2019年7月31日，參考大氣密度與真實大氣密度之間的偏差、相對偏差隨高度變化特征如圖1所示，偏差隨高度增加呈現(xiàn)正值、負值的交替變化特征(見圖1(a))，其正偏差值主要位于高度2.00～7.60 km，負偏差值主要位于高度11.40～19.90 km，整個高度范圍內(nèi)的值范圍為-0.0083 kg·m-3(14.90 km)～0.0065 kg·m-3(4.20 km)，平均值為5.2843×10-4kg·m-3；類似偏差(見圖1(a))，相對偏差隨高度增加同樣呈現(xiàn)正值、負值的交替變化特征(見圖1(b))，但高層的相對偏差值明顯大于低層，其正相對偏差值主要位于高度2.00～7.60 km，負相對偏差主要位于高度11.40～19.90 km，整個高度范圍內(nèi)的值范圍為-4.2142%(15.90 km)～5.2210%(21.40 km)，平均值為-0.3063%；對照圖1(a)與圖1(b)可知，雖然偏差、相對偏差隨高度變化趨勢基本一致，但大值區(qū)所處的高度層存在一定差異，例如，偏差極大值(0.0065 kg·m-3)、極小值(-0.0083 kg·m-3)分別位于4.20 km，14.90 km，而相對偏差極大值(5.2210%)、極小值(-4.2142%)分別位于21.40 km和15.90 km。

圖1 在2017年8月1日-2019年7月31日，參考大氣密度與真實大氣密度之間的密度偏差、密度相對偏差隨高度變化特征Fig.1 The altitude variations of deviations and relative deviations between the reference atmospheric density and the true atmospheric density from August 1, 2017 to July 31, 2019

2.2 年內(nèi)差異特征

在2017年8月1日-2019年7月31日，參考大氣密度與真實大氣密度之間偏差、相對偏差的年內(nèi)—高度變化特征如圖2所示；對于年內(nèi)變化特征，在低層，偏差值呈現(xiàn)負值、正值、負值變化特征，其正偏差極大值出現(xiàn)在7月，而在高層則表現(xiàn)為正值、負值、正值變化特征，其負偏差極大值出現(xiàn)在8月(見圖2(a))；對于高度變化特征，1-3月、11-12月偏差隨高度增加呈現(xiàn)出由負偏差向正偏差轉(zhuǎn)變過程，而在4-10月則表現(xiàn)為正偏差向負偏差轉(zhuǎn)變過程；在圖2(a)中，低層(2.0～4.0 km)的6-8月存在一個正偏差大值區(qū)，高層(13～17 km)的7-9月存在一個負偏差大值區(qū)，最小值為-0.0538 kg·m-3(12月的2.00 km)，最大值為0.0558 kg·m-3(7月的2.0 km)，平均值為5.3219×10-4kg·m-3；類似偏差(見圖2(a))，相對偏差同樣存在明顯的年內(nèi)變化及隨高度變化特征(見圖2(b))，對于年內(nèi)變化特征，低層(10 km以下)呈現(xiàn)出負值、正值、負值的變化特征，而高層(10 km以上)則表現(xiàn)為正值、負值、正值的變化特征；對于高度變化特征，1-3月、11-12月隨高度增加呈現(xiàn)出負相對偏差向正相對偏差轉(zhuǎn)變特征，而4-10月隨高度增加呈現(xiàn)出正相對偏差向負相對偏差轉(zhuǎn)變，圖2(b)中的負相對偏差最大值為-13.2699%(8月的15.90 km)，正相對偏差最大值為10.1708%(12月的21.40 km)，平均值為-0.3046%。

圖2 2017年8月1日-2019月7月31日，參考大氣密度與真實大氣密度之間的密度偏差、密度相對偏差隨年內(nèi)—高度變化特征Fig.2 The time-altitude cross section of deviations and relative deviations between the reference atmospheric density and the true atmospheric density from August 1, 2017 to July 31, 2019

2.3 季節(jié)差異特征

在2017年8月1日-2019年7月31日，不同季節(jié)參考大氣密度與真實大氣密度之間的偏差、相對偏差隨高度變化特征如圖3所示，偏差、相對偏差隨高度變化特征與季節(jié)有關(guān)，春、夏、秋季的大氣密度偏差在低層以正值為主，偏差值由大到小分別為夏季、春季、秋季，在高層以負值為主，偏差值由大到小分別為夏季、秋季、春季；而冬季的大氣密度偏差在低層以負值為主，在高層以正值為主(見圖3(a))，在整個高度范圍內(nèi)，春、夏、秋、冬季的大氣密度偏差范圍分別為-0.0055 kg·m-3(11.90 km)～0.0120 kg·m-3(2.00 km)，-0.0234 kg·m-3(14.90 km)～0.0515 kg·m-3(2.00 km)，-0.0097 kg·m-3(14.90 km)～0.0047 kg·m-3(4.20 km)，-0.0505 kg·m-3(2.00 km)～0.0134 kg·m-3(11.40 km)，平均值分別為4.8702×10-4kg·m-3，0.0083 kg·m-3，-0.0013 kg·m-3，-0.0054 kg·m-3；類似圖3(a)，相對偏差隨高度變化特征同樣與季節(jié)密切相關(guān)(見圖3(b))，其中，春、夏、秋季的大氣密度相對偏差在低層以正值為主，相對偏差值由大到小分別為夏季、春季、秋季，在高層以負值為主，相對偏差值由大到小分別為夏季、秋季、春季；而冬季的大氣密度相對偏差在低層以負值為主、高層以正值為主，在整個高度范圍內(nèi)，春、夏、秋、冬季的大氣密度相對偏差范圍分別為-3.0923%(17.90 km)～5.9970%(21.40 km)，-11.6516%(15.90 km)～5.4388%(2.00 km)，-4.7972%(15.90 km)～5.2188%(21.40 km)，-4.7440%(2.00 km)～10.0290%(21.40 km),平均值分別為-0.2805%，-0.7126%，-0.7011%，0.4637%，對照圖3(a)和圖3(b)及上述分析可知，雖然大氣密度偏差與相對偏差隨高度變化趨勢基本一致，但出現(xiàn)極大值、極小值的高度層存在差異，且在整個高度層，春季和夏季的平均偏差為正值、秋季和冬季的平均偏差為負值；而春、夏、秋季的相對偏差均為負值，冬季為正值。

圖3 2017年8月1日-2019年7月31日，不同季節(jié)參考大氣密度與真實大氣密度之間的密度偏差、密度相對偏差隨高度變化特征Fig.3 The altitude variations of deviations and relative deviations between the reference atmospheric density and the true atmospheric density from August 1, 2017 to July 31, 2019 in different seasons

2.4 平均年內(nèi)變化特征

在2017年8月1日-2019年7月31日的1-12月，整個高度層內(nèi)的平均大氣密度偏差呈現(xiàn)“∧”字符變化特征(見圖4(a))，其中，1-3月、10-12月為負偏差，4-9月為正偏差，圖4(a)中的值范圍為-0.0060 kg·m-3(11月)～0.0099 kg·m-3(7月)，平均值為5.3219×10-4kg·m-3；整個高度層內(nèi)的平均大氣密度相對偏差年內(nèi)變化特征如圖4(b)所示，與圖4(a)基本相反(從圖2(a)與圖2(b)之間的高層、低層差異可以解釋)，圖4(b)中的1-3月和12月為正相對偏差、4-11月為負相對偏差，其值范圍為-0.9569%(8月)～0.7225%(2月)，平均值為-0.3046%。

圖4 2017年8月1日-2019年7月31日，整個高度層內(nèi)參考大氣密度與真實大氣密度之間的密度偏差、密度相對偏差年內(nèi)變化特征Fig.4 The annual variation of deviations and relative deviations between the reference atmospheric density and the true atmospheric density within the entire height from August 1, 2017 to July 31, 2019

3 qαmax差異特征

3.1 qαmax差異統(tǒng)計

通過上述研究發(fā)現(xiàn)，參考大氣密度與真實大氣密度之間的差異特征與高度、季節(jié)均有密切關(guān)系，在2017年8月1日-2019年7月31日，參考大氣密度與真實大氣密度之間的qαmax偏差特征如圖5所示，qαmax偏差呈現(xiàn)“M”字型變化特征(見圖5(a)，圖中粗虛線分別對應+100 Pa·rad，-100 Pa·rad，圖5(a)中的值范圍為-217.1721～219.0648 Pa·rad，平均值為-9.4684 Pa·rad，其中，正偏差超過100 Pa·rad的樣本數(shù)有24個，負偏差超過-100 Pa·rad的樣本數(shù)有49個；qαmax偏差在不同區(qū)間范圍內(nèi)的占有率特征如圖5(b)所示，沿橫坐標軸方向呈現(xiàn)增大、減小的變化趨勢，其中，qαmax偏差在(-20, 0]范圍內(nèi)占有率最大(20.4110%)，在[0，20)范圍內(nèi)占有率為次大值(16.0274%)，在(-∞， -100]、(100， +∞)范圍內(nèi)的占有率分別為6.7123%，3.2877%；類似平均大氣密度偏差存在明顯的年內(nèi)變化特征(見圖4(a))，qαmax偏差同樣具有明顯的年內(nèi)變化特征(見圖5(c))，但變化特征與圖4(a)相反，圖5(c)中的qαmax偏差在1-3月、11-12月為正偏差，4-10月為負偏差，在1-12月的值范圍為-51.3659 Pa·rad(6月)～46.3820 Pa·rad(1月)，平均值為-9.5690 Pa·rad，qαmax偏差同樣與季節(jié)有關(guān)(見圖5(d))，其中，春、夏、秋季均為負偏差，冬季為正偏差，負偏差最大值出現(xiàn)在夏季(-42.5615 Pa·rad)，正偏差最大值出現(xiàn)在冬季(37.8531 Pa·rad)。

3.2 qαmax偏差超過100 Pa·rad

參考大氣密度與真實大氣密度之間的qαmax偏差超過100 Pa·rad時，其分布特征如圖6所示，其中，有24個樣本的qαmax偏差超過100 Pa·rad(見圖6(a))，有49個樣本的qαmax偏差超過-100 Pa·rad(見圖6(b))，在圖6(a)中，qαmax偏差值范圍為101.1910～219.0648 Pa·rad，其中，基于參考大氣密度得到的qαmax值范圍為1.8522×103～4.7914×103Pa·rad，而基于真實大氣密度得到的qαmax值范圍為1.7337×103～4.6090×103Pa·rad；在圖6(b)中，qαmax偏差值范圍為-101.4770～-217.1721Pa·rad，其中，基于參考大氣密度得到的qαmax值范圍為873.2646～2.8181×103Pa·rad，而基于真實大氣密度得到的qαmax值范圍為982.6712～3.0027×103Pa·rad。

在2017年8月1日-2019年7月31日，qαmax偏差值超過100 Pa·rad時，對應的風場、大氣密度偏差隨高度變化特征如圖7所示。從圖7(a)可以看出，雖然qαmax偏差≤-100、qαmax偏差≥100對應的風場隨高度變化趨勢基本一致，且極大值出現(xiàn)在基本相同的高度層(12 km左右)、極小值均出現(xiàn)在2 km，但qαmax偏差≥100對應的風速在各高度層均大于相應的qαmax偏差≤-100，在整個高度范圍內(nèi)，qαmax偏差≤-100，qαmax偏差≥100對應的風速范圍分別為5.9592 m·s-1(2.8 km)～40.0294 m·s-1(12.90 km)，7.2083 m·s-1(2.0 km)～58.8748 m·s-1(11.90 km)，平均值分別為18.0332 m·s-1，33.8605 m·s-1，qαmax偏差≤-100、qαmax偏差≥100對應的大氣密度偏差隨高度變化特征存在明顯差異(見圖7(b))，其中，qαmax偏差≤-100對應大氣密度偏差隨高度增加呈現(xiàn)正偏差向負偏差轉(zhuǎn)變，而相應的qαmax偏差≥100為負偏差向正偏差轉(zhuǎn)變，即兩者表現(xiàn)為相反的變化特征，在整個高度范圍內(nèi)，qαmax偏差≤-100，qαmax偏差≥100對應的大氣密度偏差范圍分別為-0.0217 kg·m-3(13.90 km)～0.0386 kg·m-3(2.0 km)，-0.0539 kg·m-3(2.0 km)～0.0178 kg·m-3(11.0 km)，平均值分別為0.0058 kg·m-3，-0.0052 kg·m-3。

圖5 2017年8月1日-2019年7月31日，qαmax偏差變化及其在不同區(qū)間范圍內(nèi)的占有率、年內(nèi)變化和季節(jié)變化特征Fig.5 The change characteristics, percentages in different ranges, annual variation and seasonal variation of qαmax deviation from August 1, 2017 to July 31, 2019

圖6 在qαmax偏差≥100 Pa·rad、qαmax偏差≤-100 Pa·rad時，基于真實大氣密度與參考大氣密度計算得到的qαmax值變化特征Fig.6 Under the qαmaxdeviation≥100 Pa·rad, qαmax deviation≤-100 Pa·rad conditions, the variation characteristics of qαmax value calculated based on true atmospheric density and reference atmospheric density

圖7 qαmax偏差值超過100 Pa·rad時，風速、密度偏差隨高度變化特征Fig.7 The altitude variations of upper wind speed, atmospheric density deviation when the value of qαmaxdeviation exceeds 100 Pa·rad

4 建 模

基于目前計算大氣環(huán)境對運載火箭飛行的qαmax程序包中，大氣密度為參考大氣密度、而非真實大氣密度，通過上述研究可以看出，由于參考大氣密度與真實大氣密度之間存在差異，導致該程序包計算得到的qαmax與真實值之間存在差異，偏差最大值超過200 Pa·rad，為減小運載火箭飛行的風險，一方面，強烈建議采用真實大氣密度代替程度包中的參考大氣密度，另一方面，可考慮采用多元線性回歸、BP神經(jīng)網(wǎng)絡等建模方法建立訂正模型，本文采用常用的多元線性回歸方法、BP神經(jīng)網(wǎng)絡方法建立訂正模型，并對訂正模型的效果進行檢驗。

4.1 建模思想

基于2017年8月1日-2019年7月31日，一共730個樣本，其中:1)前600個樣本作為訓練樣本，后130個樣本作為檢驗樣本；2)基于參考大氣密度計算得到qαmax對應的運載火箭飛行高度h，獲取0～h高度范圍內(nèi)的緯向風、經(jīng)向風、參考大氣密度平均值，分別記為umean，vmean，ρmean；3)日期以天數(shù)表示，其中，1月1日為1，1月31日為31，…，記為d；4)計算參考大氣密度與真實大氣密度之間的qαmax偏差；5)對2017年8月1日-2019年7月31日(共730個樣本)的h，umean，vmean，ρmean，d歸一化處理，即將數(shù)值轉(zhuǎn)換為[-1,1]范圍內(nèi)；6)針對前600個樣本，利用多元線性回歸、BP神經(jīng)網(wǎng)絡方法對h，umean，vmean，ρmean，d與qαmax偏差之間建立模型；7)基于已建立好的訂正模型，利用后130個樣本的h，umean，vmean，ρmean，d參數(shù)，得到訂正模型的qαmax偏差值(記為Dqαmax)，Dqαmax與參考大氣密度計算得到的qαmax(記為Rqαmax)之和，即為所需要的結(jié)果。

基于參考大氣密度、真實大氣密度得到的qαmax分別記為Rqαmax，Tqαmax，BP神經(jīng)網(wǎng)絡、多元線性回歸方法建立訂正模型得到的qαmax，分別記為Bqαmax，Mqαmax；本文基于多元線性回歸方法建立的模型如下：

49.89umean-40.42vmean-83.18ρmean-9.94d-

103.56h-23.08=Dqαmax

(9)

Mqαmax=Rqαmax-Dqαmax

(10)

在式(9)～(10)中，umean，vmean，ρmean，d和h的變量意義與前面描述一致，為后130個樣本的已知量，Rqαmax為基于參考大氣密度得到的后130個樣本qαmax；利用式(9)可計算得到Dqαmax，Dqαmax反映了Rqαmax與基于多元線性回歸方法得到Mqαmax之間的偏差；利用式(10)得到的Mqαmax，即為基于多元線性回歸方法建立模型得到的qαmax訂正值。

本文使用的BP神經(jīng)網(wǎng)絡方法，是采用MATLAB軟件自帶的神經(jīng)網(wǎng)絡工具箱進行建模，同樣先得到Dqαmax，其中，Dqαmax反映了Rqαmax與基于BP神經(jīng)網(wǎng)絡方法得到Bqαmax之間的偏差；利用式(11)即可得到基于BP神經(jīng)網(wǎng)絡方法的Bqαmax，即

Bqαmax=Rqαmax-Dqαmax

(11)

式中：Rqαmax為基于參考大氣密度得到的后130個樣本qαmax。

記Rqαmax，Bqαmax，Mqαmax與Tqαmax之間的偏差分別為RqαmaxDev，BqαmaxDev，MqαmaxDev，下面利用數(shù)理統(tǒng)計方法對BP神經(jīng)網(wǎng)絡、多元線性回歸方法建立的訂正模型效果進行檢驗。

4.2 模型檢驗

在后130個樣本中，RqαmaxDev，BqαmaxDev，MqαmaxDev的變化特征及在不同區(qū)間范圍占有率分布特征如圖8所示，從圖8(a)可以看出，通過對參考大氣密度的結(jié)果進行BP神經(jīng)網(wǎng)絡、多元線性回歸建模訂正處理，所得到的qαmax偏差明顯減小，在圖8(a)中，RqαmaxDev，BqαmaxDev，MqαmaxDev的值范圍分別為-201.6939～122.1358 Pa·rad，-75.9311～39.4346 Pa·rad，-52.6029～60.1846 Pa·rad；RqαmaxDev，BqαmaxDev，MqαmaxDev在不同區(qū)間范圍內(nèi)的占有率存在明顯差異(見圖8(b))，其中，RqαmaxDev在區(qū)間(-20，0]占有率最高(19.2308%)，在區(qū)間(0，20]占有率次高(15.3846%)，BqαmaxDev在區(qū)間(-20，0]占有率最高(31.5385%)，在區(qū)間(-40，-20]占有率次高(27.6923%)，MqαmaxDev在區(qū)間(-40，-20]占有率最高(27.6923%)，在區(qū)間(0，20]占有率次高(26.9231%)，從圖8(b)還可以看出，RqαmaxDev，BqαmaxDev，MqαmaxDev占有率為0.00%的區(qū)間范圍也存在明顯差異，RqαmaxDev在區(qū)間(60,80]占有率為0.00%，BqαmaxDev在區(qū)間(-∞，-100]，(-100，-80]，(40，60]，(60，80]，(80，100]，(100，∞]內(nèi)的占有率均為0.00%，MqαmaxDev在區(qū)間(-∞，-100]，(-100，-80]，(-80，-60]，(80，100]，(100，∞]內(nèi)的占有率均為0.00%。Rqαmax，Bqαmax，Mqαmax與Tqαmax之間的平均偏差、平均絕對差、平均相對偏差、平均標準差及相關(guān)系數(shù)統(tǒng)計結(jié)果見表1。

圖8 RqαmaxDev，BqαmaxDev，MqαmaxDev的變化及其在不同區(qū)間范圍占有率分布特征Fig.8 The variation characteristics, and percentages in different ranges of RqαmaxDev，BqαmaxDev，MqαmaxDev

表1 Rqαmax，Bqαmax，Mqαmax與Tqαmax之間的數(shù)理統(tǒng)計特征Table 1 The mathematical statistics results between Rqαmax，Bqαmax，Mqαmax and Tqαmax

5 結(jié) 論

運載火箭飛行前，需確定發(fā)射場大氣環(huán)境對運載火箭飛行的影響，若計算得到的qαmax超過閾值，則發(fā)射推遲；在計算qαmax的程序中，所用的大氣密度通常是參考大氣密度，而非發(fā)射場零日的真實大氣密度，目前針對參考大氣密度對qαmax精度的影響研究較少，本文以2017年8月1日-2019年7月31日高垂直分辨率探空資料為例，分析參考大氣密度與真實大氣密度以及它們計算得到某型運載火箭飛行qαmax之間的差異特征，最后建立訂正模型，得到主要結(jié)論如下：

1)密度偏差、相對偏差隨高度增加而呈現(xiàn)正值、負值的交替變化特征，整個高度層內(nèi)的平均大氣密度偏差、相對偏差具有明顯年內(nèi)變化特征，其中，偏差值在1-3月、10-12月為負值，4-9月為正值，值范圍為-0.0060 kg·m-3(11月)～0.0099 kg·m-3(7月)，相對偏差在4-11月為負值，1-3月、12月為正值，值范圍為-0.9569%(8月)～0.7225%(2月)。

2)參考大氣密度與真實大氣密度之間的qαmax偏差隨時間變化呈現(xiàn)“M”字型特征，其值范圍為-217.1721～219.0648 Pa·rad，平均值為-9.4684 Pa·rad；qαmax偏差值在不同區(qū)間范圍內(nèi)占有率存在較明顯差異，其最大占有率(20.4110%)位于(-20,0]區(qū)間范圍。

3)參考大氣密度與真實大氣密度之間的qαmax偏差，具有明顯的年內(nèi)、季節(jié)變化特征；其中，1-3月、11-12月為正偏差，4-10月為負偏差，值范圍為-51.3659 Pa·rad(6月)～46.3820 Pa·rad(1月)，春、夏、秋季為負偏差，冬季為正偏差，值范圍為-42.5615 Pa·rad(夏季)～37.8531 Pa·rad(冬季)。

4)通過對qαmax偏差≥100 Pa·rad、qαmax偏差≤-100 Pa·rad對應的風場、大氣密度偏差進行對比發(fā)現(xiàn)，qαmax偏差≥100 Pa·rad對應的風速在各高度層均明顯偏大、大氣密度偏差隨高度增加呈現(xiàn)負值向正值轉(zhuǎn)變，而qαmax偏差≤-100 Pa·rad對應的風速在各高度層均偏小、且大氣密度偏差為正值向負值轉(zhuǎn)變。

5)基于參考大氣密度得到的qαmax偏差，通過利用多元線性回歸、BP神經(jīng)網(wǎng)絡方法建立模型進行訂正，結(jié)果表明均能夠較好地提高qαmax精度，且多元線性回歸建立的模型改進效果更好；與參考大氣密度得到的qαmax相比，多元線性回歸、BP神經(jīng)網(wǎng)絡方法建立模型得到的qαmax偏差由-36.1569 Pa·rad減小到-7.7012 Pa·rad，-17.2267 Pa·rad，絕對差由49.3254 Pa·rad減小到20.3244 Pa·rad，23.9357 Pa·rad，相對偏差由-2.0201%減小到-0.7741%，-1.2917%，標準差由55.5491 Pa·rad減小到23.2982 Pa·rad，23.8563 Pa·rad，相關(guān)系數(shù)由0.9952提高到0.9960，0.9960。

通過上述研究結(jié)果可以看出，由于參考大氣密度與真實大氣密度之間的偏差存在年內(nèi)、季節(jié)變化特征，導致所計算得到的qαmax偏差同樣存在較明顯年內(nèi)、季節(jié)變化，且正、負偏差最大值均超過200 Pa·rad；考慮到目前計算程序中所用的大氣密度通常是參考大氣密度，而非發(fā)射場零日真實大氣密度，若程序計算得到的qαmax接近閾值，則此時真實的qαmax可能超過閾值，導致運載火箭飛行存在重大風險。雖然采用多元線性回歸或BP神經(jīng)網(wǎng)絡建立訂正模型可以提高qαmax精度，但同樣存在一定風險，因此，在分析大氣環(huán)境對運載火箭飛行的影響時，最好采用發(fā)射場零日真實大氣密度代替參考大氣密度，從而得到運載火箭飛行時真實的qαmax，進一步提高運載火箭飛行的保障能力，避免出現(xiàn)重大事故；考慮到火箭發(fā)射前，需提前幾天確定發(fā)射場區(qū)發(fā)射零日是否滿足發(fā)射條件，然而，此時沒有未來幾天的大氣密度實況資料、且考慮大氣密度與多種因素有關(guān)(溫度、相對濕度等)，因此，通過建立大氣密度預報模型，獲取更接近真實大氣密度的大氣密度預報值，代替計算程序中的參考大氣密度，可能是后期進一步提高精度的研究方向。