高階線性比例制導系統脫靶量冪級數解

赫泰龍，陳萬春，周浩

北京航空航天大學 宇航學院，北京 100083

比例導引是最經典的制導律，由于其簡潔、有效和易于物理實現，目前世界上戰術導彈幾乎都采用比例導引制導[1-3]。脫靶量是設計分析制導系統的關鍵性能指標，通過研究線性化比例制導系統，Zarchan[1]提出了零控脫靶量的概念，用于解釋比例導引，還用于推導和理解其他先進的制導律[3-6]。直接對制導系統微分方程進行數值仿真是求解脫靶量的通用方法，但是伴隨法才是戰術導彈制導系統最主要的分析設計工具[1,7-9]。對于一階線性比例制導系統，當有效導引比為正整數時，利用伴隨方程可直接得到脫靶量的解析解[1,10]；但是對于一般的高階制導系統，并不存在脫靶量的解析解。另外，比例導引伴隨系統的初始時刻是該微分方程的正則奇點[11-12]，在伴隨仿真時需要在伴隨時間上增加一個小量來避免奇異，不過大量的應用算例表明伴隨仿真結果是可靠的[1]。

冪級數法在求解特殊線性微分方程、非線性微分方程和實際工程問題中都有重要應用[13-17]。但是利用冪級數法來研究比例導引制導系統的工作還很少，文獻[11]給出了一階制導系統冪級數解的簡單示例；Holt[18]研究了一個特殊的三階比例制導系統，該系統由3個不同帶寬的單延遲環節表示，并且得到了脫靶量的冪級數解。文獻[11,18]的工作都是針對特殊的制導系統，缺乏通用性,都沒有給出冪級數解的收斂性證明；也沒有研究冪級數部分和逼近精確解的速度[15,19]、精度和適用區間。

本文針對一般的高階線性比例導引制導系統，推導得到了冪級數和指數函數e-kt乘積的形式的脫靶量的解，并證明了該冪級數解的收斂性。然后給出了一階系統和高階二項式系統冪級數解簡化遞推關系。考慮實際應用，本文還重點分析了參數k對冪級數解收斂速度的影響，同時對比伴隨仿真結果，給出了參數k的選擇方法。冪級數解為比例制導系統的設計評估提供了一種新的手段。

1 線性比例導引制導系統及脫靶量

考慮一般的線化比例導引制導系統，研究目標階躍機動和導彈初始瞄準誤差角引起的脫靶量。圖1給出了原線化比例導引制導回路以及伴隨系統框圖，脫靶量則可以由伴隨系統一次仿真得到[1]。N為有效導引比，δnT為目標階躍機動轉化而來得脈沖輸入，δHE為導彈瞄準誤差角初始狀態轉化而來的脈沖輸入，VC為導彈和目標接近速度，nL為導彈活度的加速度，伴隨系統中t為剩余飛行時間或總的飛行時間，nT為目標機動水平大小，VM為導彈速度，θHE為導彈初始瞄準誤差角，MHE為導彈初始瞄準誤差角引起的脫靶量，MnT為目標階躍機動引起的脫靶量；記伴隨系統狀態分別為z1、z2、zu和ζ，則可得到伴隨系統微分方程為

(1)

圖1 線性比例導引制導系統及其伴隨系統Fig.1 Linear proportional navigation guidance system and adjoint system

式中:z1的初值為1,是由伴隨系統的脈沖輸入轉化而來。該伴隨系統輸出得到脫靶量為

(2)

G(s)表示一般的穩定傳遞函數，可能包含導引頭動力學、噪聲濾波、飛控系統等環節；通常G(s)可以表示為

(3)

即G(s)是由Q1個一階環節和Q2個二階環節組成，αi、ξj和βj為各環節特征參數系數，無量綱的正數；不失一般性，令

(4)

T為參考時間常數或總制導系統時間常數，具有時間的量綱。為了推導一般形式的脫靶量的解，將G(s)關于s分母多項式展開可得到

(5)

其中:Q=Q1+2Q2、λ0，λ1，…，λQ為多項式系數，由式(3)中的αi、ξj和βj唯一確定，結合式(4)可以得到

λ0=1,λ1=1

(6)

進一步，便于應用，將伴隨式(1)、式(2)進行無量綱化：

(7)

(8)

經過簡單求導運算可以得出，無量綱的伴隨狀態關于無量綱的時間的微分方程與式(1)相同，只需將傳遞函數替換為

(9)

因此，為了簡化表達式，如無特別說明，后文中無量綱或歸一化的變量仍使用原變量符號，方程的求解和結果討論都是針對無量綱化的變量。

2 伴隨方程的冪級數解

2.1 脫靶量冪級數的系數的遞推關系

一般地，伴隨系統微分式(1)不存在由有限項初等函數構成的解析解，但是當G(s)為一階環節且有效導引比N為正整數時，利用Laplace變換可以得到式(1)解析解[1]，N等于3、4和5的結果(未進行歸一化)如表1所示。

表1 N為正整數時一階制導系統伴隨方程的解析解

注意表1中解析解都是e指數與多項式函數乘積的形式，受此啟發，當N不是正整數或者制導系統的階數是高階時，可以求解和探討無量綱式(1)如下形式的冪級數解：

(10)

(11)

(12)

(13)

式中:an、bn、cn、dn分別為各級數的待定系數，參數k為指數項衰減常數，可用來調節冪級數解整體收斂速度，這些冪級數都假設為收斂的。注意ζ的級數解中含有tQ項，這是因為關于ζ動態是由傳遞函數式(9)來描述的，等價于如下微分方程：

(14)

式中:ζ(q)為變量ζ的q階導數。對z1、ζ、z2和zu級數求導可得

(15)

式中:A和C及其上下標分別表示排列數和組合數。

將式(10)～式(13)及其相應的導數式(15)代入微分式(1)和式(14)，等號兩側的指數部分消掉，利用關于時間多項式各次冪的系數相等，并結合伴隨方程的狀態初值，可以得到如下遞推關系：

(16)

當n≥1時

(17)

式中:

Bn和Pn都是中間變量，用于簡化表達式書寫；A和C及其上下標仍然表示排列數和組合數;cn和dn分別為歸一化脫靶量MHE和MnT冪級數的系數。

2.2 脫靶量冪級數的收斂半徑為無窮大

定理對任意參數k，脫靶量冪級數的收斂半徑為無窮大。

證明由于e-kt在t=0處泰勒展開冪級數的收斂半徑是無窮大，容易證明如果參數k=0時，由遞推關系式(16)和式(17)所確定的冪級數式(10)～式(13)是收斂的，則k為任意數時，由同樣遞推關系所確定的冪級數仍然收斂，且收斂半徑相同。

記狀態向量

(18)

則微分方程式(1)和式(14)可以寫成如下狀態空間描述：

(19)

式中:

(20)

A(t)=

(21)

以及狀態初值

(22)

顯然，矩陣R的特征值只有0，任何正整數都不是R的特征值；函數矩陣A(t)是常值矩陣，在t=0處是解析的，而且其冪級數展開收斂半徑為無窮大，則由文獻[14]中定理6.1可得，微分方程式(19)的解可以表示為收斂半徑為無窮大的冪級數：

(23)

式中:Xn為向量值級數系數，維數與X相同。將式(23)代入式(19)，利用冪級數恒等條件，可以得到遞推關系為

(24)

事實上，t=0是微分式(1)的正則奇點，在t=0處式(1)的解是解析的，也就是在t=0的鄰域內，式(1)的解可以表示為收斂的冪級數展開式，而且由于R和A(t)的特殊性，或者從遞推式(24)也可以直接得到，冪級數解的收斂半徑為無窮大。

為了便于研究脫靶量冪級數解的性質，本文還會用到級數的前n+1項部分和Sn與級數的n+1項以后的余式Rn。例如，當研究目標階躍機動引起的脫靶量冪級數解式(13)時,有前面已經證明冪級數解收斂半徑為無窮大，這意味著在我們關心的歸一化的時間區間[0,tf]里，冪級數是一致收斂的，可以用部分和Sn來一致地逼近脫靶量真實的解，而且理論上可以以任意精度逼近，只要n足夠大。

(25)

3 脫靶量冪級數解討論分析

本節將針對一階環節、高階二項式環節以及一般的高階環節等制導系統，研究其在特定參數k時的脫靶量冪級數解，并對比伴隨仿真結果進行討論分析。

3.1 一階環節制導系統

考慮如下一階滯后環節的制導系統：

(26)

這里G(s)是歸一化后傳遞函數，可以用來代表簡化的飛控系統或噪聲濾波等環節；相應地，Q=1,λ0=1,λ1=1，代入一般的高階系統冪級數解待定系數的遞推關系式(16)和式(17)，可得如下遞推關系:

(27)

當n≥1時

(28)

注意到bn遞推中含有(1-k)，如果取參數k=1，則遞推關系可以得到很大程度簡化，而且經過簡單的數列技巧變換，可以得出級數z2和zu的待定系數cn和dn的遞推關系:

(29)

(30)

以及cn和dn的通項表達式分別為

(31)

(32)

顯然，當N為大于2的正整數時成立

(33)

所以當N為正整數時，歸一化脫靶量z2和zu的冪級數部分就轉化為有限項的多項式，從而得到脫靶量的解析解；特別地，當取N=3,4,5時，就可以得到表1中所列出的結果，也就是用冪級數法同樣可以得到解析解。

當N為正整數時，脫靶量冪級數解就是有限項的解析解，與伴隨仿真結果重合，當N為非正整數時，得不到脫靶量的解析解。所以作為算例，這里給出N=3.5和N=4.5時脫靶量冪級數解部分和Sn與伴隨仿真結果的對比，如圖2所示。

圖中橫軸和縱軸分別表示歸一化的時間和脫靶量。由于脫靶量解中指數部分隨著時間快速衰減，而且主要研究導彈末段自尋的制導回路，所以重點關心歸一化的時間在10以內的脫靶量情況。從圖中可以看出，一階制導系統脫靶量冪級數解收斂速度很快，取級數的前7項部分和就可以很好的逼近精確脫靶量，特別是在總飛行時間較小的時候。

參數k表示級數中指數部分衰減常數，對于一階制導系統傳遞函數式(26)，其齊次解或脈沖響應輸出的衰減常數正是1，仿真結果顯示也符合脫靶量隨飛行時間衰減的速率，因此選擇參數k=1是合理的，而且還直接得到了待定系數cn和dn的通項，以及N為正整數時脫靶量的有限項解析表達式。這里只給出目標階躍機動引起的脫靶量結果，導彈初始瞄準角誤差引起脫靶量冪級數解性質類似，后文算例也只討論目標階躍機動脫靶量。

圖2 一階制導系統目標階躍機動脫靶量的伴隨仿真結果與冪級數解部分和Fig.2 Adjoint simulation results and partial sums of power series solutions for single-lag guidance system in step target maneuver

3.2 高階二項式環節制導系統

考慮如下歸一化的Q階二項式制導系統傳遞函數：

(34)

是由Q個相同的一階環節乘積得到，形式簡單，常用于高階制導系統性能的初步分析[1]。對比式(9)可以得到

參考3.1節一階制導系統參數k的選取，考慮到Q階二項式傳遞函數式(34)齊次響應指數衰減常數為Q，選取k=Q符合二項式制導系統脫靶量隨飛行時間衰減的速率，而且經過排列組合等運算，可以將通用的待定系數的遞推關系式(17)簡化，即當n≥1時

(35)

(36)

(37)

進一步利用一些技巧變換，可以得到n≥1時，待定系數cn和dn的等效遞推關系為

(38)

(39)

當然，在求解cn的前提下，利用式(35)中的遞推關系求解dn計算量更少。

3.3 一般的高階制導系統

對于一般的高階制導系統，可以利用通用的待定系數的遞推關系式(16)和式(17)，得到脫靶量冪級數；以五階制導系統為例，傳遞函數G(s)分別選取二項式形式：

(40)

各時間常數不等的一階環節乘積形式：

(41)

以及三個一階環節和一個二階環節乘積形式：

(42)

式中:G2(s)和G3(s)中各參數取值如表2所示。

參數k=9時，3種時間常數分布不同的五階制導系統，目標階躍機動引起的脫靶量冪級數解部分和Sn(t)與伴隨仿真結果的對比如圖3所示，算例中N=4。從圖3中可以看出，隨著n增加，部分和Sn都逐漸逼近脫靶量伴隨仿真結果(或者精確解)，S70已經與伴隨結果相差無幾，而且直觀上3個五階制導系統逼近速度(級數收斂速度)差不多。另外，不同的時間常數分布的五階制導系統脫靶量結果很相近[1]。

表2 G2(s)和G3(s)中各參數取值Table 2 Values of parameters in G2(s) and G3(s)

圖3 不同時間常數分布的五階制導系統目標階躍機動脫靶量的伴隨仿真結果與冪級數解部分和逼近Fig.3 Adjoint simulation results and partial sums approximation of normalized miss distance in target maneuver for fifth-order guidance systems with different configurations

4 參數k對冪級數解收斂速度的影響

本節將進一步分析參數k對冪級數解收斂速度的影響，給出選取參數k的方案。

首先來看一下不同參數k對同一制導系統冪級數解的影響，這里選取前面含有一個二階環節的五階制導系統G(s)=G3(s)。圖4給出了k分別取3、7和10時脫靶量冪級數解(部分和Sn)與伴隨仿真結果，算例中N=4。

圖4 含有一個二階環節的五階制導系統目標階躍機動脫靶量的伴隨仿真結果與冪級數解部分和逼近Fig.4 Adjoint simulation results and partial sums approximation of normalized miss distance in target maneuver for fifth-order guidance system with a quadratic distribution

從圖4中可以看出，k=7時大約需要60項部分和(即S60)就可以很好地逼近精確解，而k=10時大約80項，k=3時則要超過150項，k=7時冪級數解收斂速度要快于k=3和k=10時。這可以從冪級數解zu的假設形式(13)來給出解釋，具體地，zu是由常規冪級數和指數函數e-kt乘積組成。當k較大時，指數函數隨著t的衰減速度很快，部分和Sn(t)中指數部分占主導(在n較小時)，比如k=10時S60(t)的曲線在t=8處就衰減到0附近，導致偏離脫靶量精確解。當k較小時，指數函數隨著t的衰減速度要小于脫靶量精確解的實際衰減速度，而多項式在t較大時隨著t的增大是發散的，此時部分和Sn(t)中多項式部分在t較大處可能占主導，比如k=3時S80(t)、S100(t)、S120(t)的曲線隨著t的增加嚴重偏離脫靶量精確解。

至此，討論冪級數部分和Sn來逼近脫靶量精確解的誤差都是指截斷誤差，也就是余式Rn，理論上，Rn隨著n的增加趨于0。如圖4(a)所示，隨著n的增加，Sn逼近誤差較小的區間逐漸擴大，這與普通函數的冪級數展開性質類似；但是當n增大到約150時，就無法通過增加部分和項數n來進一步改善逼近效果，甚至在t=8附近S150(t)的曲線呈現出近似隨機振蕩的特征，這是由于數值計算過程中舍入誤差[20]的累積導致的，本文所有數值算例默認使用雙精度浮點數，對于S150(t)中高階多項式部分在t較大時求值運算(Horner方法[21])雙精度是不夠的。

為了進一步準確客觀分析，引入用來衡量冪級數收斂速度指標變量ncr，其意義是使得部分和逼近誤差余式Rn小于指定精度ε的索引變量n的最小值，換句話說，至少需要ncr+1項部分和，才會使得逼近精確解的誤差小于ε；ncr與k有關，將ncr表示成關于k的函數形式[19]

ncr(k)?minnRn(t)≤ε;k

(43)

從定義可以看出ncr越小越好，意味著收斂速度越快。使得ncr取值最小的k，是最優的，記為

(44)

注意ncr和kopt與整個制導系統的參數和分析條件有關，包括制導系統階數、時間常數分布、有效導引比N和分析時間t、誤差精度ε等。圖5給出了不同參數和條件下冪級數解收斂速度指標ncr(k)的曲線，為了保證數值精度，采用四精度甚至更高精度浮點數運算，并且用級數前1 001項(足夠精確了)部分和S1000(t)作為脫靶量精確解來計算余式，即Rn(t)=S1000(t)-Sn(t)。

圖5(a)考慮高階二項式制導系統，N=4，t=10，ε=10-4，針對不同的系統階數Q得到ncr(k)的曲線，可以看出最佳參數kopt與階數Q有關，但都在Q附近，選取參數k=Q是合理的。當k較小時，二項式系統階數越高，冪級數解收斂速度越慢；當k較大時，不同階數制導系統冪級數解收斂速度幾乎相同。

圖5(b)考慮一般的五階制導系統，N=4，t=10，ε=10-4，針對時間常數分布不同的傳遞函數G(s)得到ncr(k)的曲線，可以看出最佳參數kopt與系統時間常數分布有關，分布范圍越大(如G2(s))kopt越大，二項式分布kopt最小。當k較小時，時間常數散布越大(如G2(s))，ncr越大，冪級數解收斂速度越慢；當k較大時，不同時間常數分布的制導系統冪級數解收斂速度幾乎相同。

圖5(c)針對一般的含有一個二階環節的五階制導系統G(s)=G3(s)，N=4，t=10，得到不同誤差精度ε時ncr(k)曲線，顯然精度越高(即ε越小)，ncr越大，逼近所需要的部分和項數越多。不同誤差精度ε，ncr(k)曲線形式相同，kopt也幾乎不變。

圖5 不同參數和條件下冪級數解收斂速度指標ncr(k)的曲線Fig.5 Plots of convergence rate index ncr(k) as a function of k with different system parameters and conditions

圖5(d)針對一般的含有一個二階環節的五階制導系統G(s)=G3(s),N=4,ε=10-4，在不同時間t處得到ncr(k)的曲線，可以看出t越大，逼近到同樣精度ε所需要的部分和項數越多；ncr(k)曲線形式相同，kopt也幾乎不變。

綜上分析，同時考慮冪級數收斂速度和遞推關系的簡化，可以按照如下方案選取參數k：一階制導系統選取k=1，Q階二項式系統選取k=Q，對于一般的高階系統式(3)或式(9)可選取

(45)

而相應的部分和項數可由ncr(k)確定。通常k取值大一些更好，這時不同制導系統收斂速度幾乎相同，而且數值運算穩定性更好，避免雙精度浮點數運算舍入誤差累積發散；考慮到收斂速度，k取值一般不超過10。

5 結 論

本文研究了線性比例導引制導系統脫靶量的冪級數解，通過理論論證和大量數值算例表明冪級數解是分析比例導引制導系統新的有效手段。

1) 推導了一般高階比例制導系統脫靶量的冪級數解系數的遞推關系，通用性強。

2) 嚴格證明了脫靶量冪級數解的收斂半徑為無窮大。

3) 針對一階環節和高階二項式環節等制導系統，通過選取適當的參數k，得到了冪級數解系數簡化的遞推關系，并研究了不同參數和條件下冪級數解部分和的逼近性質。

4) 通過數值方法討論分析了參數k對脫靶量冪級數解收斂速度的影響，給出了選取參數k以及部分和項數的方法，為冪級數解的實際應用奠定了基礎。