運載火箭動力學天地差異性分析技術

程 興，高 晨，賈大玲，李 君，王紫揚

（北京宇航系統工程研究所，北京，100076）

0 引 言

近年來，空間發射任務逐漸多樣化、復雜化，因此對運載火箭提出了更高的要求[1～4]。彈道、姿控和制導等專業的設計和集成對于運載火箭非常重要，吸引了眾多學者的關注并且取得了許多成果[5～7]。

運載火箭的工程化設計一般采用“獨立設計+聯合分析仿真”的模式，即彈道、姿控、制導、載荷、結構、動力等專業根據標準彈道各自獨立設計，各個專業根據自身需求進行偏差條件下的包絡設計。這樣的方法可以有效保證設計效率，但是由于各個專業進行設計時使用的動力學模型不同，參數偏差的取值和考慮因素各有側重，難以準確分析多項偏差組合對運載火箭姿態、軌跡、載荷等參數的影響。因此，如何準確、高效進行總體回路的集成建模，量化分析各種因素對運載火箭質量特性以及飛行動力學的影響，是現階段運載火箭設計中急需解決的問題。

另一方面，運載火箭是復雜的系統工程，影響因素多且相互間耦合，尤其是一級飛行段（含助推飛行段）經歷惡劣氣動環境，缺乏對其飛行力學特征進行量化分析的手段，難以開展動力學關鍵參數的天地一致性確認及非故障狀態的偏差辨識。從飛行結果來看，盡管飛行過程中姿態穩定并高精度地將有效載荷送入預定軌道，但各子級的姿態及軌道參數乃至飛行時間均較標稱值存在偏差，尤其是一級飛行段的參數偏差更為突出。若能有效辨識飛行偏差背后的物理因素，據此開展設計改進與優化，提升飛行與設計的一致性，將能夠提高設計精細化水平，從而有利于提高飛行可靠性及潛力挖掘，需采用有效的辨識算法進行參數辨識和一致性分析[8～10]。

目前，常用的辨識算法包括卡爾曼濾波、極大似然估計、分割算法等。但是卡爾曼濾波算法對模型精確程度要求較高；基于極大似然估計的辨識方法一般要求解雅克比矩陣，算法計算量較大，工程實現較為困難，而且算法在系統非線性較強時可能失效，需要克服矩陣奇異的問題。運載火箭飛行環境干擾因素復雜繁多，因此所設計算法的魯棒性尤為重要。粒子群算法、遺傳算法等算法不依賴于初始值，具有較好的尋優能力，而且計算量較小，具有較好的工程應用價值[11～14]。

針對速度偏差、姿態偏差及關機時間偏差幅值較大的狀況，以實現火箭動力學天地一致性分析，本文首先提出源變量動力學建模技術，并據此實現總體回路集成仿真；在此基礎上，利用粒子群算法具有較快的收斂速度，遺傳算法具有較強全局搜索能力的特點，采用分層結構設計混合優化算法進行參數辨識，保證系統具有較快的搜索速度和較強的全局搜索能力，從而實現運載火箭動力學參數的天地一致性分析。最后通過仿真驗證了所提方法的有效性。

1 基于源變量的總體回路集成仿真技術

液體運載火箭普遍采用彈道、姿控、制導、原始數據計算、分離等多專業獨立設計的模式[15,16]。該設計模式下各專業以標準彈道為起點各自獨立設計并按專業準則進行偏差包絡設計，既保證了設計效率，也保證了系統適應有限偏差的能力。該模式下，各專業使用的動力學數學模型各有差異、參數偏差取值的考慮因素各有側重，這使得針對特定偏差的分析具有周期長、匹配性低的特點，難以量化分析單項或多項偏差對飛行姿態、飛行軌跡、飛行載荷的綜合影響，譬如推進劑流量下降后將影響推進劑質量及晃動特性，影響全箭質量分布及質心、轉動慣量，過載及速度等。另外，彈道設計采用瞬時平衡準則，忽略姿態控制的動態響應過程，姿態動力學建模以標準彈道為基礎，通常采用質量瞬時凝固法等，即專業間的模型缺乏一致性，以提高各專業間參數的匹配性、量化分析參數偏差對系統影響為目標，提出了基于源變量的動力學集成建模技術。

為了精確量化分析包括發動機推力異常、推進劑流量異常等偏差或故障對全箭質量特性及飛行動力學的影響，提出源變量建模思路：以發動機安裝角、方位角、擺角、推進劑秒耗量、加注量、全箭分布質量等參數為基本變量（即源變量），通過在線計算推進劑流量及剩余量、全箭質量特性（質量、質心、轉動慣量）、發動機推力及推進劑晃動特性、結構彈性振動及載荷分布，時時迭代計算全箭飛行軌跡和姿態動力學特性，實現總體-彈道-姿控-制導-動力-推進劑晃動-載荷計算的集成建模，并以此推動一體化仿真。

相對文獻[15]和文獻[16]中提出的運載火箭六自由度仿真方法不同，本仿真模型具有3個特點：

a）發動機推力P為安裝角、安裝象限、擺角的矢量化模型。

式中ckε為第k臺發動機的安裝角；ckμ為第k臺發動機的象限角；ckδ為第k臺發動機的擺動指令角。

該建模方式下，通過模型變換便能簡潔、高效地消除剛-晃-彈動力學模型中的代數環問題而顯著提高仿真效率及精度。這里不作動力學模型的逐項詳細分析，而僅列出耦合模型形式。動力學的完整表達形式及其物理含義詳見文獻[1]。

即有：

進一步變換有：

通過式（5）消除了動力學模型中的代數環問題。

相對于將質心動力學方程建立在發射慣性系或發射系的傳統建模方法，基于箭體坐標系的質心方程（2）能大大簡化變換矩陣T，并顯著提升計算效率及精度。

c）總體回路專業深度耦合。

相對基于事先計算推進劑晃動特性、彈性振動特性、全量質量特性的傳統運載火箭六自由度仿真，本仿真通過離線-在線組合、參數間耦合機理的再梳理等方法研究，實現基于分布質量的全箭質量特性參數（質心、轉動慣量、液位高度、增壓壓力）在線計算，推進劑晃動特性參數的在線計算（晃動質量、晃動頻率、晃動質心高度），以及飛行軌跡、飛行姿態、晃動響應、彈性響應等響應參數的在線計算。

源變量建模及仿真方法優化的詳細過程可參見文獻[17]。

源變量仿真（六自由度全量模型）與標準彈道（三自由度質心模型）的速度、位置偏差對比見圖1。

圖1 源變量仿真與彈道計算參數的對比Fig.1 Flight Parameter between Svarsim and Trajectory Design

圖1 的對比結果表明考慮控制過程的動態響應與僅考慮瞬時平衡的飛行速度略有差異，該偏差在系統設計允許的誤差范圍內，再次驗證彈道計算采用瞬時平衡處理的工程有效性，也驗證了基于源變量的多專業耦合仿真模型的正確性。

源變量集成仿真實現了多專業間的參數緊耦合，能量化分析一個或多個參數偏差對包括飛行軌跡、飛行姿態在內的飛行品質參數的影響。因此以該模型為數學基礎，開展運載火箭動力學的天地差異性辨識。

2 基于源變量仿真的動力學天地差異性辨識

下面以全程復現實際飛行中的軌跡、姿態、推進劑液位高度等關鍵飛行參數為目標，研究運載火箭參數的天地差異性。

影響飛行姿態、飛行軌跡及關機時間的因素多，其中包括各發動機流量偏差、比沖偏差，推力線橫移，各推進劑貯箱的加注量偏差、結構質量偏差，氣動法向力系數偏差、阻力系數偏差、壓心偏差，全箭質心橫移，發動機零位及伺服機構零位漂移等。除加注量偏差、結構質量偏差外，其余參數偏差均隨時間變化，總計216個待辨識量。

為適應待辨識偏差參數多的狀況，采用粒子群算法-遺傳算法組合的新型搜索算法。

目標函數為仿真與飛行之間的速度偏差、姿態角偏差、角速度偏差、控制擺角偏差及助推器質量（液位高度）偏差均方和最小，具體為

式中 ki=0,…, 4 分別為速度偏差、姿態偏差、角速度偏差、控制擺角及液位對應質量偏差的權重系數；為仿真中的發動機關機時間；為對應飛行值；，分別為各助推器氧化劑、燃燒劑質量仿真質量與飛行質量之間的偏差。

姿態角偏差、角速度及擺角偏差遠小于速度偏差及質量偏差，因此通過 ki來調節實現參數對結果的均衡化，本次分析中具體取值為 k1=10，k2=200，k3=400，

2.1 遺傳算法通過模擬生物進化論實現參數的辨識和尋優

遺傳算法包括 4個基本元素，即編碼、選擇、交叉和變異，具有較強的全局搜索能力。粒子群算法源于對鳥類捕食行為的研究，具有較快的收斂速度。本文基于遺傳算法和粒子群算法，采用分層設計的思路，設計參數辨識混合算法。算法的底層采用遺傳算法，對系統參數進行全局搜索和辨識；頂層設計采用粒子群算法，對遺傳算法得到的最優解作為初始種群，從而進行快速的局部搜索，得到參數辨識的結果。同時，為了提高粒子群算法的快速性，本文通過設計動態慣性權值對算法進行設計。定義待辨識系統參數為 xi，中間代個體定義為 yi，最大迭代次數定義為 G1，變異概率記為Pm。

參數辨識底層遺傳算法的流程為：

a）步驟1：通過確定待辨識參數的可能取值范圍，采用二進制的方法對各個待辨識參數進行編碼，即：

結合所有待辨識參數的可能取值范圍，通過二進制的方法對參數進行編碼。例如，可以將ix的范圍通過相應位數的二進制串進行描述：

b）步驟2：如果迭代至最大次數，則算法終止；反之，執行步驟3。

c）步驟3：通過式（9）、式（10）選出一定數目的中間個體，并且進行交叉運算。

d）步驟4：對待優化的所有個體以給定概率進行變異操作。

e）步驟5：計算辨識目標函數。

f）步驟 6：將當代最優個體與歷史最優個體進行比較并更新。

g）步驟7：判斷算法得到的最優解是否滿足停止條件，若滿足，則系統輸出最優解，反之則繼續進行尋優。

當采用遺傳算法迭代至代數G1后，將系統的最優解取出，作為粒子群算法的初始值。采用粒子群優化算法對遺傳算法得到的結果進行優化。采用粒子群算法迭代至代數G2之后，若滿足停止條件，則將粒子群優化后的結果作為系統辨識的最優結果輸出；若不滿足停止條件，則用粒子群種群中的 q個粒子隨機替換掉遺傳算法子群中的 q個粒子，對遺傳算法子群再一次進行優化求解，不斷循環，系統滿足停止條件，從而輸出最優解。

2.2 采用粒子群算法對系統參數進行辨識

采用粒子群算法對系統參數進行辨識時，慣性權值和參數收斂速度密切相關，較大的慣性權值有利于全局搜索；較小的慣性權值有利于算法的局部開發，加速算法的收斂。本文通過設計動態慣性權值來提高粒子收斂速度，通過式（11）、式（12）來更新粒子的位置和速度：

式中 r1，r2為服從均勻分布且相互獨立的隨機數；g1，g2為給定常數；vi(j)為待辨識參數粒子i在第j代的速度；為粒子i在第j代的位置； Pi( j)為粒子i在第j代的最優位置為第j代全局最優位置。設動態慣性權值?為

可以看出，隨著采用粒子群算法迭代次數的增多，系統慣性權值減小，而較小的慣性權值可以加快算法的收斂速度。

粒子群算法的流程為：

a）步驟1：對粒子群算法的種群進行初始化，即采用遺傳算法的最優解作為粒子群算法種群的初始位置，并隨機產生每個粒子的速度。

b）步驟2：判斷系統是否到達最大迭代數，若到達，則算法終止，反之則執行步驟3。

c）步驟3：計算粒子群算法種群每個粒子的適應度函數值。

d）步驟4：更新粒子群算法種群的個體極值和全局極值。

e）步驟5：更新粒子群算法中每個粒子的速度和位置。

f）步驟 6：判斷優化得到的最優解是否滿足停止準則，若滿足，則輸出最優解；反之，跳轉到步驟 2繼續進行參數優化。

2.3 提高辨識效率和有效性處理方式

a）包括增壓壓力在內的部分參數直接采用飛行參數而非設計參數；用于氣動插值的飛行速度、高度參數也直接用飛行參數而非在線計算參數；

b）氣動數據偏差包括法向力系數偏差、壓心系數偏差，采用乘法形式的修正系數，為時變參數，時間及修正值均通過仿真確定，為提高效率，約束時間最小間隔；

c）辨識基于Matlab/Simlink環境，充分利用高版本支持整數辨識的功能，而對部分參數進行整數規格化處理，譬如設定壓心補償系數=0.9～1.1，遺傳尋優中設定其上限為1100，下限為900，分層值為1，這樣便大幅提升尋優效率。

尋優辨識之后的飛行速度偏差變化歷程見圖2，姿態偏差變化歷程見圖3。另外，關機偏差從飛行相對于標準彈道提前4.73 s降低到當前的0.07 s，即可認為關機時間復現飛行狀態（圖中數據已歸一化處理）。

圖2 源變量仿真與飛行之間的速度偏差（辨識結果）Fig.2 Velocity Bias between Svarsim and Flight (by Identification)

圖3 源變量仿真與飛行之間的姿態偏差Fig.3 Attitude Bias between Svarsim and Flight (by Identification)

2.4 部分參數辨識結果

a）部分助推器流量偏差較大，其中表1為3#助推器的氧化劑、燃燒劑流量偏差結果（經過歸一化處理）。

表1 3#助推器氧化劑、燃燒劑流量偏差Tab.1 3# Booster Flow Deviation of Oxidizer and Fuel

該數據表明實際飛行中的流量偏大，對應的推力亦偏大，因此提前耗盡；同時關機前的飛行速度也一直大于彈道設計值，與實際飛行情況吻合。

b）針對圖3中60 s前后X向速度偏差較大的狀況，開展氣動阻力系數偏差辨識，結果表明氣動阻力修正難以消除該項偏差。

c）100 s后的X向速度偏差持續增大，則可能為推進劑溫升所致，分析中的火箭采用推進劑加溫后自身增壓方案，存在增壓氣體與推進劑換熱的狀況，缺乏推進劑溫度測量數據，難以進一步精細化分析。

d）對彈道計算中存在的“底部阻力”項進行了對比性分析，結果表明，不考慮“底部阻力”項，辨識后的彈道與飛行彈道更接近，否則100 s后的X向、Y向速度偏差將增加，其中X向速度偏差將增加到辨識結果的兩倍以上。

e）氣動壓心及法向力系數亦存在偏差，部分時刻的壓心偏差已經接近甚至超過設計要求考慮的±5%偏差范圍，即存在工程設計不包絡（已經出現的飛行偏差）風險。

3 結束語

基于源變量的動力學集成仿真技術有機耦合多專業間的參數，為運載火箭姿態動力學天地一致性分析提供了數學基礎，本次包括遺傳算法在內的新型搜索算法很好地解決了辨識參數多的障礙（本次辨識中超過210個參數），辨識出的“不考慮底部阻力更接近飛行”、“部分時刻氣動壓心偏差已經接近甚至超過±5%的設計考慮偏差邊界”等結果，將對后續優化設計提供支撐；同時為火箭飛行測量參數優化設置提供參考。

在本次任務分析的基礎上，快速開展了多次不同型號的飛行數據分析，顯示出本方法的任務適應能力。

分析中還發現，待辨識的偏差參數過多使得分析效率偏低（本次任務中的210個參數，大約需要自動執行3天時間）；也存在阻力修正系數、壓心系數等部分參數容易陷入局部最優點的狀況，后續將探索針對性更強的辨識尋優算法。