一種基于平穩(wěn)滑翔彈道解析解的快速彈道規(guī)劃方法
【專利摘要】本發(fā)明公開了一種基于平穩(wěn)滑行彈道解析解的快速彈道規(guī)劃方法，包括以下幾個(gè)步驟：步驟1：滑翔段彈道規(guī)劃問題建模；步驟2：滑翔段彈道規(guī)劃變量設(shè)計(jì)；步驟3：滑翔段彈道解析解求解；步驟4：滑翔段彈道終端速度控制方案；步驟5：滑翔段彈道再入走廊調(diào)整方案；步驟6：滑翔段彈道規(guī)劃初值生成；步驟7：滑翔段彈道規(guī)劃流程設(shè)計(jì)。本發(fā)明提出以縱向機(jī)動(dòng)加速度比例系數(shù)和橫向機(jī)動(dòng)加速度比例系數(shù)為滑翔段彈道規(guī)劃變量，使得運(yùn)動(dòng)方程中彈道傾角、彈道偏角、高度、經(jīng)度與緯度的微分方程中不含速度項(xiàng)；本發(fā)明獲得了固定縱向機(jī)動(dòng)加速度比例系數(shù)和橫向機(jī)動(dòng)加速度比例系數(shù)對(duì)應(yīng)的滑翔段彈道解析解。
【專利說明】一種基于平穩(wěn)滑翔彈道解析解的快速彈道規(guī)劃方法

【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及一種基于平穩(wěn)滑行彈道解析解的快速彈道規(guī)劃方法，屬于航天技術(shù)、 武器【技術(shù)領(lǐng)域】。

【背景技術(shù)】
[0002] 隨著高超聲速技術(shù)的快速發(fā)展，再入彈道規(guī)劃成為了研究熱點(diǎn)?；瓒问窃偃霃?道的主要構(gòu)成部分，決定了再入飛行的航程和機(jī)動(dòng)范圍。滑翔段彈道快速規(guī)劃不僅可以分 析高超聲速飛行器的性能，還可用于在線彈道規(guī)劃和預(yù)測(cè)制導(dǎo)，具有較高的研究價(jià)值。
[0003] 由于滑翔段彈道具有飛行時(shí)間長，對(duì)規(guī)劃變量高度敏感，且再入過程非線性約束 強(qiáng)，使得滑翔彈道的可行域較窄，采用傳統(tǒng)優(yōu)化算法求解時(shí)往往需要很長的計(jì)算時(shí)間。盡管 偽譜法、SQP等優(yōu)化算法在計(jì)算效率上有較大的提高，其計(jì)算耗時(shí)依然在秒級(jí)以上，并且當(dāng) 過程約束碰邊界時(shí)，其計(jì)算效率會(huì)大大降低，甚至影響到結(jié)果的收斂性。為了提高滑翔段彈 道規(guī)劃的速度，一些特殊的約束被引入滑翔段彈道規(guī)劃，如阻力曲線、等熱流曲線、等動(dòng)壓 曲線、平衡滑翔條件等。這些約束的共同特點(diǎn)是將滑翔段彈道限定在某種特殊的平穩(wěn)滑翔 狀態(tài)，從而降低了彈道規(guī)劃問題的敏度。但由于在彈道規(guī)劃時(shí)還需要進(jìn)行數(shù)值積分，限制了 規(guī)劃速度的提高。為了擺脫對(duì)彈道的積分的依賴，動(dòng)態(tài)逆被引入滑翔段彈道規(guī)劃。該方法 直接規(guī)劃幾何彈道形狀，并通過調(diào)整彈道形狀來滿足約束要求。但由于滑翔段前后的環(huán)境 變化劇烈，使得滑翔段全程彈道難以用低階曲線擬合。


【發(fā)明內(nèi)容】

[0004] 本發(fā)明的目的是通過求解平穩(wěn)滑翔彈道高精度解析解，獲得不依賴彈道積分的滑 翔段彈道快速規(guī)劃方法，為分析高超聲速飛行器性能和預(yù)測(cè)制導(dǎo)提供技術(shù)支持。
[0005] -種基于平穩(wěn)滑行彈道解析解的快速彈道規(guī)劃方法，包括以下幾個(gè)步驟：
[0006] 步驟1 :滑翔段彈道規(guī)劃問題建模；
[0007] 步驟2 :滑翔段彈道規(guī)劃變量設(shè)計(jì)；
[0008] 步驟3 :滑翔段彈道解析解求解；
[0009] 步驟4 :滑翔段彈道終端速度控制方案；
[0010] 步驟5 :滑翔段彈道再入走廊調(diào)整方案；
[0011] 步驟6 :滑翔段彈道規(guī)劃初值生成；
[0012] 步驟7 :滑翔段彈道規(guī)劃流程設(shè)計(jì)。
[0013] 本發(fā)明的優(yōu)點(diǎn)在于：
[0014] (1)提出以縱向機(jī)動(dòng)加速度比例系數(shù)和橫向機(jī)動(dòng)加速度比例系數(shù)為滑翔段彈道 規(guī)劃變量，使得運(yùn)動(dòng)方程中彈道傾角、彈道偏角、高度、經(jīng)度與緯度的微分方程中不含速度 項(xiàng)；
[0015] (2)獲得了固定縱向機(jī)動(dòng)加速度比例系數(shù)和橫向機(jī)動(dòng)加速度比例系數(shù)對(duì)應(yīng)的滑翔 段彈道解析解，包括高度解析解、飛行距離解析解、彈道偏角解析解、經(jīng)度解析解和緯度解 析解；
[0016] (3)給出了縱向機(jī)動(dòng)加速度比例系數(shù)和橫向機(jī)動(dòng)加速度比例系數(shù)的初值生成方 法；
[0017] (4)提出以縱向機(jī)動(dòng)加速度比例系數(shù)、橫向機(jī)動(dòng)加速度比例系數(shù)、橫向反轉(zhuǎn)彈道 傾角和初始彈道傾角增量為滑翔段彈道規(guī)劃變量，其中：縱向機(jī)動(dòng)加速度比例系數(shù)對(duì)應(yīng)滑 翔段飛行距離；橫向機(jī)動(dòng)加速度比例系數(shù)對(duì)應(yīng)終端速度；橫向反轉(zhuǎn)彈道傾角對(duì)應(yīng)終端經(jīng)緯 度；初始彈道傾角增量對(duì)應(yīng)過程約束，上述四者之間相互獨(dú)立。
[0018] (5)本發(fā)明采用滑翔彈道解析解規(guī)確定彈道，而積分僅僅用于終端速度校正，使得 彈道規(guī)劃具有非常快的速度。

【專利附圖】

【附圖說明】
[0019] 圖1是平穩(wěn)滑翔彈道規(guī)劃建模流程圖；
[0020] 圖2是大圓坐標(biāo)系；
[0021] 圖3是4與乂{的關(guān)系；
[0022] 圖4是A Y對(duì)過程約束的影響；
[0023] 圖5是球面飛行距離估計(jì)；
[0024] 圖6是1與8{的關(guān)系；
[0025] 圖7是滑翔段彈道快速規(guī)劃方法流程；
[0026] 圖8是縱向彈道解析解精度驗(yàn)證；
[0027] 圖9是橫向彈道解析解精度驗(yàn)證；
[0028] 圖10是Bell解析解；
[0029] 圖11是可達(dá)區(qū)域邊界彈道攻角曲線族；
[0030] 圖12是可達(dá)區(qū)域邊界彈道傾側(cè)角曲線族；
[0031] 圖13是可達(dá)區(qū)域邊界橫向彈道曲線族；
[0032] 圖14是可達(dá)區(qū)域邊界縱向彈道曲線族；
[0033] 圖15是可達(dá)區(qū)域邊界彈道速度曲線族；
[0034] 圖16是可達(dá)區(qū)域邊界彈道熱流密度曲線族；
[0035] 圖17是可達(dá)區(qū)域邊界彈道動(dòng)壓曲線族；
[0036] 圖18是可達(dá)區(qū)域邊界彈道過載曲線族；
[0037] 圖19是滑翔段彈道規(guī)劃可達(dá)區(qū)域；
[0038] 圖20是平面規(guī)劃彈道；
[0039] 圖21是規(guī)劃彈道攻角曲線；
[0040] 圖22是規(guī)劃彈道傾側(cè)角曲線。

【具體實(shí)施方式】
[0041] 下面將結(jié)合附圖和實(shí)施例對(duì)本發(fā)明作進(jìn)一步的詳細(xì)說明。
[0042] 本發(fā)明在分析運(yùn)動(dòng)方程的耦合條件的基礎(chǔ)上，提出了一種以固定縱向機(jī)動(dòng)加速度 比例系數(shù)和橫向機(jī)動(dòng)加速度比例系為控制變量的彈道模式。該彈道模式下，速度與其他狀 態(tài)變量之間解耦，從而獲得了解耦的彈道形狀微分方程，并利用該微分方程求解得到了高 精度的三維彈道形狀解析解。此基礎(chǔ)上，構(gòu)造了滑翔段彈道的快速規(guī)劃方法：首先利用彈 道形狀解析解規(guī)劃出滿足位置約束的彈道；然后通過調(diào)整橫向機(jī)動(dòng)大小來滿足終端速度約 束；最后通過調(diào)整初始彈道傾角來滿足過程約束要求。
[0043] 本發(fā)明是一種基于平穩(wěn)滑行彈道解析解的快速彈道規(guī)劃方法，規(guī)劃方法建模流程 如圖1所示，規(guī)劃流程如圖7所示，包括以下幾個(gè)步驟：
[0044] 步驟1 :滑翔段彈道規(guī)劃問題建模
[0045] 為了便于理論分析，假設(shè)地球?yàn)閳A球體，不考慮地球自轉(zhuǎn)的影響，則半速度坐標(biāo)系 下的三自由度質(zhì)點(diǎn)運(yùn)動(dòng)方程如下所示，

【權(quán)利要求】
1. 一種基于平穩(wěn)滑行彈道解析解的快速彈道規(guī)劃方法，包括以下幾個(gè)步驟： 步驟1 :滑翔段彈道規(guī)劃問題建模 設(shè)地球?yàn)閳A球體，不考慮地球自轉(zhuǎn)，則半速度坐標(biāo)系下的三自由度質(zhì)點(diǎn)運(yùn)動(dòng)方程如下 所示，
式中，h為飛行器的高度，r為從地心到飛行器質(zhì)心的距離，h = I-Rtl，其中Rtl為地球半 徑；Θ、Φ分別為經(jīng)度和緯度；s為飛行距離；V為相對(duì)地球的速度；Y為當(dāng)?shù)貜椀纼A角；V 為彈道偏角；A、i、6、一、戶、夕和P分別為高度、飛行距離、經(jīng)度、緯度、彈道傾角、彈道 偏角和速度對(duì)時(shí)間的導(dǎo)數(shù)；σ為傾側(cè)角；g = μ /r2為當(dāng)?shù)刂亓铀俣?，μ為地球重力加?度相關(guān)的常數(shù)；L和D分別為升力加速度和阻力加速度，其表達(dá)式為， L = pV2SrefCL/(2m) (2) D = P V2SrefCn/(2m) (3) 式中，m為飛行器質(zhì)量；Sref為氣動(dòng)參考面積；Q和Cd分別為升力系數(shù)和阻力系數(shù)，與馬 赫數(shù)Ma和攻角α相關(guān)；P為大氣密度，如下所示， f = f seae ⑷ 式中，Psm為海平面大氣密度；β為指數(shù)大氣模型常數(shù)； 過程約束如下所示
式中，為熱流密度，為最大熱流密度，k為常系數(shù)；q為動(dòng)壓，qmax為最大動(dòng)壓；η 為總過載，nmax為最大總過載，gQ = 9. 81m/s2 ; 步驟2 :滑翔段彈道規(guī)劃變量設(shè)計(jì) 引入縱向機(jī)動(dòng)加速度aE和橫向機(jī)動(dòng)加速度ae，表達(dá)式如下， aE = Lcos 〇 + (V2/r-g) cos y (6) ae = Lsin 〇/cos γ + (V2/r) cos Y sin Ψ tan Φ (7) 式中，aE和&0分別為縱向機(jī)動(dòng)加速度和橫向機(jī)動(dòng)加速度；將式（6)和式（7)分別帶入 式（1)，則彈道傾角微分方程和彈道偏角微分方程為，
假設(shè)aE <0,則滑翔段的Y單調(diào)遞減，將運(yùn)動(dòng)方程中自變量替換為Y ;將式（1)除以 式（8)中的彈道傾角微分方程，得
式中
分別為高度、飛行距離、彈道偏角、經(jīng)度、緯度 和速度的微分方程；其中，第一行為滑翔段縱向運(yùn)動(dòng)微分方程，第二行為滑翔段橫向運(yùn)動(dòng)微 分方程，第三行為速度微分方程； 設(shè)滑翔段的aE與&0形式如下，
式中，kE和1^分別為縱向機(jī)動(dòng)加速度比例系數(shù)和橫向機(jī)動(dòng)加速度比例系數(shù)；將式 (10)帶入式（9)，并取 sinY = y、cosy = 1，則有
以kE和ke作為規(guī)劃變量，將滑翔段彈道規(guī)劃問題轉(zhuǎn)化為求解非線性方程的問題； 步驟3 :滑翔段彈道解析解求解； 設(shè)滑翔段的初始彈道傾角為γd，初始高度為Iv初始彈道偏角為Vci、初始經(jīng)度和緯度 為分別為Qtl和（K;滑翔段的終點(diǎn)高度為hf，對(duì)于給定的kE和1^，對(duì)式（11)至（15)式進(jìn) 行解析求解； 第一、高度解析解 將式⑷帶入式（11)可得，
對(duì)式（16)積分可得，
式中，為初始彈道傾角；P ^為滑翔段初始大氣密度；對(duì)式（17)進(jìn)一步整理可得，
式中，C1為一常數(shù)，取值為 第二、飛行距離解析解
將式（18)帶入式（12)可得，
設(shè)初始飛行距離Stl = 0,則式（19)積分可得， s = fs(y)-fs(y〇) (20) 式中，fs( Y)為由式（19)不定積分獲得的函數(shù)，其表達(dá)式與C1的正負(fù)相關(guān)，如下所示，
式中，X為自變量； 第三、道偏角解析解 進(jìn)一步，由式（13)積分可得， Ψ = ￥〇+(k0/kE) (y-y〇) (21) 式中，Ψο為初始彈道偏角； 第四、緯度解析解 將cos Ψ在π /2附近5階泰勒展開，可得，
取y= π/2-Ψ，其中y為積分中間變量，并將式（18)、式（21)和式（22)帶入式（15) 可得，
式中，C2 為一常量,滿足 C2 = Y0+kE 〇 /2-￥。)/1^ ;取， Ca〇 = k0/(60r β kE)Cal = 2C2k0/kE (24)
式中，Ca(l、Cal和Ca2均為常量；將式（24)帶入式（23)可得，
式中，Cbl、Cb2、Cb3和Cb4均為常量；利用式（26)可將式（25)化為，
對(duì)式（27)積分可得， Φ = Φ ο+f φ (y) -f φ (y〇) (28) 式中，Φci為初始緯度；y〇與初始彈道偏角相關(guān)，取值為y〇= π/2-(y)為與y相 關(guān)的函數(shù)，如下所示
式中，gjy)為與y相關(guān)的函數(shù)，如下所示
第五、經(jīng)度解析解 式（18)、式（21)和式（29)分別給出了滑翔段彈道的高度、彈道偏角和緯度隨彈道傾角 變化規(guī)律的解析解，將它們分別代入式（14)可得，
利用Guass-Legendre求積公式直接獲得結(jié)果，具體表達(dá)式如下，
式中：
；Y i為[Y，Y〇]中的第i個(gè)高斯節(jié)點(diǎn)；Ai為求積系數(shù)；η為高斯節(jié)點(diǎn) 總數(shù)； 步驟4 :滑翔段彈道終端速度控制方案； 設(shè)負(fù)比能量e = μ/Γ-ν2/2,則由式（1)可得， e = VD (32) 式中，?為負(fù)比能量的導(dǎo)數(shù)；由式（8)、式（10)和式（32)可得，
由式（33)看出，Cd是速度求解的關(guān)鍵；通過情況下，Cd可寫成Q和e的函數(shù)，如下所 示， Cd = fCD (Cl, e) (34) 式中，fCD為升力系數(shù)與阻力系數(shù)的關(guān)系函數(shù)；求解式（33)積分的關(guān)鍵為求解Q，由式 (6)、（7)式和式（10)可得，
另外，由式（4)和式（18)可得，
由（35)、式（36)和式（37)可求得Q，將Q帶入式（34)可的Cd的大?。? 最終，將r和Cd的表達(dá)式帶入式（33)，并采用龍格-庫塔法進(jìn)行積分，獲得終端負(fù)比能 量，從而解得終端的速度大??； 當(dāng)規(guī)劃的彈道包含傾側(cè)反轉(zhuǎn)點(diǎn)時(shí)，則反轉(zhuǎn)點(diǎn)前后的ke需要滿足如下連接條件，
式中,1^_和ke+分別為反轉(zhuǎn)點(diǎn)前后的橫向機(jī)動(dòng)加速度比例系數(shù)；P。、Ψ。和Φ。分別為 Y。處的大氣密度、彈道偏角和緯度； 步驟5 :滑翔段彈道再入走廊調(diào)整方案； 初始彈道傾角取值如下所示， Y 0 = Y *+ Δ γ (39) 式中，Λ Y為調(diào)整彈道再入走廊位置的初始彈道傾角分量；滿足攻角平滑性要求 的初始彈道傾角分量，滿足，
式中Kn = Lcos 〇 /D,為縱向升阻比； 步驟6 :滑翔段彈道規(guī)劃初值生成； 在進(jìn)行滑翔段彈道規(guī)劃前，需要對(duì)sf、Y'kE和ke的初值進(jìn)行估算； 第一、Sf估算 球面坐標(biāo)系下，（Θ。，Φ。）和（Θ f，Cj5f)的距離為，
式中，Ss為（θ〇, φ〇)至（0f，Cj5f)的球面距離；F0和$分別滿足，
則終端飛行距離為，
式中，Ko =<-K，其中以為滑翔段起點(diǎn)子午線與射面所在大圓的夾角； 第二、^估算 在滑翔段，縱向升阻比與終端飛行距離存在如下關(guān)系，
將式（42)帶入式（40)，求得^的估值為，
第二、kE估算 設(shè)Λ Y = 0,則Ytl = ^ ;由式（17)和式（20)可知，Yf、kE和Sf滿足如下關(guān)系，
Sf = fs ( Y f) -fs ( Y 〇) 式中，。和Pf分別為終端彈道傾角和終端大氣密度；由式（44)可知，給定sf，則存在 唯一的kE與之對(duì)應(yīng)；將式（41)得到的Sf帶入式（44)求解即可獲得kE的估值； 第四、ke估算 假設(shè)Vf-Vtl僅用于修正V1O的影響，則根據(jù)平面圓弧彈道假設(shè)有， ￥f-￥〇 = 2 Ψ10 ke的初始估值為， k0 = 2￥10CN1/(y-y0) (45) 步驟7 :滑翔段彈道規(guī)劃流程設(shè)計(jì)； 具體為： 第一、根據(jù)步驟6,獲取sf、和1^的規(guī)劃初值，并假設(shè)Λ Y = 〇和Y。= Yq; 第二、根據(jù)給定的Λ Y、Yc;、kE和1^，利用步驟3計(jì)算Θ Uf)和φ(Υ?);其中， Y。為反轉(zhuǎn)點(diǎn)彈道傾角；Θ ( Y f)和Φ ( Y f)為規(guī)劃彈道終端經(jīng)度和緯度； 第三、修正Y。，使得終端緯度滿足要求，修正的方法為牛頓迭代法，迭代的終止條件為 Φ ( I < Φ1?ηι--;式中Φ1?πιη為諱度規(guī)劃精度；若I Φ ( I > Φ1?πιη則轉(zhuǎn)向步 驟7中的第二步； 第四、修正Sf或者kE使得終端經(jīng)度滿足要求，修正方法為牛頓迭代法，迭代終止條件 為I θ (Yf)-0f| < 0limit;式中0limit為經(jīng)度規(guī)劃精度；若I θ (Yf)-0f|彡0limit則轉(zhuǎn)向 步驟7中的第二步； 第五、根據(jù)Λ Y、Yc;、kE和1^，利用步驟4積分活動(dòng)終端速度V(Yf)，并 采用牛頓迭代法修正ke，使得|VUf)-Vf| < Vlimit ;式中Vlimit為速度規(guī)劃精度；若 V(Yf)-VfI彡Vlimit則轉(zhuǎn)向步驟7中的第二步； 第六，根據(jù)上一步的積分彈道判斷是否滿足過程約束，并對(duì)Λ Y進(jìn)行如下修正，
式中，h為實(shí)際飛行高度；hlimit為由最大熱流密度約束、最大動(dòng)壓約束和最大過載約束 確定的再入走廊高度下邊界；min(h-hlimit) < O表示滑翔段彈道超出再入走廊，需要增大 Λ Y，并轉(zhuǎn)向步驟7中的第二步；min (h-hlimit)彡O表示滑翔段彈道滿足再入走廊要求；Ylim 為min(h_hlimit)取極值時(shí)的彈道傾角；Λ Y(°ld)和Λ Y(nM)分別為修正前和修正后的初始彈 道傾角修正量；h(Ylini)、S(Y lini)和hliDlit( Y1J分別為彈道傾角取Ylini時(shí)的高度、飛行距 尚和聞度下邊界。
【文檔編號(hào)】G06F17/50GK104392047SQ201410691412
【公開日】2015年3月4日 申請(qǐng)日期:2014年11月25日 優(yōu)先權(quán)日:2014年11月25日
【發(fā)明者】陳萬春, 胡錦川 申請(qǐng)人:北京航空航天大學(xué)