李鑫冉 王 歆

(1中國科學(xué)院紫金山天文臺(tái)南京210008)(2中國科學(xué)院空間目標(biāo)與碎片觀測(cè)重點(diǎn)實(shí)驗(yàn)室南京210008)(3中國科學(xué)院大學(xué)北京100049)

基于遺傳算法的極短弧定軌?

李鑫冉1,2,3?王 歆1,2?

(1中國科學(xué)院紫金山天文臺(tái)南京210008)
(2中國科學(xué)院空間目標(biāo)與碎片觀測(cè)重點(diǎn)實(shí)驗(yàn)室南京210008)
(3中國科學(xué)院大學(xué)北京100049)

空間目標(biāo)的巡天觀測(cè)獲取了海量的極短弧觀測(cè)數(shù)據(jù),而經(jīng)典初軌計(jì)算方法對(duì)于極短弧幾乎不能獲得合理的結(jié)果.將初軌計(jì)算問題轉(zhuǎn)換為兩個(gè)三變量的分層優(yōu)化問題,采用遺傳算法,針對(duì)具體問題選擇了優(yōu)化變量以及相應(yīng)的遺傳操作,建立了一種極短弧初軌計(jì)算方法.基于實(shí)測(cè)資料的數(shù)值實(shí)驗(yàn)表明,方法可為后續(xù)工作提供有效的初值.

航天器,天體力學(xué),方法:數(shù)值

1 引言

短弧定軌,也稱初軌計(jì)算或軌道計(jì)算,是指沒有先驗(yàn)信息的情況下,根據(jù)單站單圈的少量觀測(cè)數(shù)據(jù)計(jì)算二體意義下的軌道.測(cè)角資料的軌道計(jì)算是天體力學(xué)中的經(jīng)典問題,兩百多年來實(shí)測(cè)數(shù)據(jù)的軌道計(jì)算方法研究從未停止過,并且隨著新的天體以及新的觀測(cè)手段的出現(xiàn),不斷涌現(xiàn)出新的需求和困難[1].

經(jīng)典的軌道計(jì)算方法是由Gauss和Laplace針對(duì)小行星的軌道計(jì)算問題建立的,人造衛(wèi)星上天以后有了較大發(fā)展,Escobal對(duì)其進(jìn)行了詳細(xì)的介紹[2].Ta ff對(duì)各種方法進(jìn)行了分類比較[3].我國學(xué)者也做了較多的研究,主要集中在改進(jìn)Laplace方法方面,并結(jié)合我國自主的實(shí)測(cè)資料情況做了一些改進(jìn)[4-10].

初軌計(jì)算問題始終沒有得到滿意的解決,根本上是由其求解方程的本征病態(tài)引起的,短弧定軌問題是不適定的.面對(duì)這種困難,為了能夠獲得較好的結(jié)果,只能依賴采集盡量長的弧段,對(duì)于光學(xué)觀測(cè)5～6 min的弧段才足以獲得較好的結(jié)果[11],在實(shí)踐中一般要求最少連續(xù)采集3 min以上的弧段.

隨著空間目標(biāo)數(shù)量的不斷增加,近年來空間目標(biāo)觀測(cè)中新出現(xiàn)了針對(duì)空域的觀測(cè)模式,通過對(duì)空域的持續(xù)觀測(cè),記錄下經(jīng)過該空域的所有目標(biāo)的方位信息.這種方式改變了原來一臺(tái)設(shè)備只能同時(shí)對(duì)單個(gè)目標(biāo)進(jìn)行跟蹤的情況,對(duì)同時(shí)經(jīng)過空域的目標(biāo)可以全部進(jìn)行監(jiān)視.這種方式提高了觀測(cè)效率,大幅度提高了觀測(cè)數(shù)量,但與此同時(shí)也帶來了新的問題,即使采用大視場(chǎng)望遠(yuǎn)鏡,對(duì)于低軌道空間目標(biāo)能夠采集到的弧長也是十分有限的,例如中國科學(xué)院空間目標(biāo)與碎片光學(xué)觀測(cè)網(wǎng)采用的超大視場(chǎng)望遠(yuǎn)鏡采集到的弧段大多不超過30 s,相當(dāng)多的僅有10余秒.近年來國際上許多大型巡天項(xiàng)目也都將空間目標(biāo)作為其觀測(cè)對(duì)象,這些項(xiàng)目采集了海量數(shù)據(jù),極短弧定軌問題就成了有效利用數(shù)據(jù)的關(guān)鍵[12].

經(jīng)典方法對(duì)于如此短弧段的初軌計(jì)算幾乎無法適用,主要體現(xiàn)在計(jì)算不收斂或者得到完全不合理(unrealistic)的解(軌道半長徑小于地球赤道半徑).這樣短的弧段被稱為“極短弧”(TSA,Too Short Arc),以區(qū)別傳統(tǒng)意義上的短弧[13?14].極短弧的具體弧長尚無明確定義,文獻(xiàn)[11]認(rèn)為2 min以下弧段為極短弧,而文獻(xiàn)[14]則采用了1 min的弧段,一般認(rèn)為極短弧指采用經(jīng)典方法很難獲得合理解的弧段,而具體弧長可以根據(jù)具體對(duì)象確定.

除了經(jīng)典計(jì)算方法外,軌道計(jì)算中還有一類方法為優(yōu)選法,采用優(yōu)選法可克服經(jīng)典方法中迭代不收斂的現(xiàn)象.但常應(yīng)用在一些特殊情況,特別是可轉(zhuǎn)換為一維優(yōu)選問題的情況,如圓軌道的初軌計(jì)算等.對(duì)于多維情況,計(jì)算方法復(fù)雜不實(shí)用,應(yīng)用極少.另一方面這些方法通常需要一個(gè)較為準(zhǔn)確的初值,而初值選取本身就是一個(gè)初軌計(jì)算問題,對(duì)于極短弧軌道計(jì)算問題也不適用[1].

Ansalone等[14]提出將遺傳算法用于解決極短弧定軌問題,由此可以避免迭代不收斂以及傳統(tǒng)優(yōu)選法需要較為準(zhǔn)確的初值和多維問題求解復(fù)雜的問題.他們針對(duì)天基模擬資料采用雙r優(yōu)選法,構(gòu)造了遺傳算法的一個(gè)具體實(shí)現(xiàn),但其采用的遺傳運(yùn)算與一般規(guī)律很不一致.Hinagawa等[15]將該方法用于同步衛(wèi)星的一般初軌計(jì)算,并對(duì)部分遺傳運(yùn)算環(huán)節(jié)做了調(diào)整.

本文也采用遺傳算法用于極短弧定軌問題,采用了完全不同于現(xiàn)有方法的變量和運(yùn)算構(gòu)造方法,實(shí)現(xiàn)了一種具體方法,并結(jié)合中國科學(xué)院空間目標(biāo)與碎片光學(xué)觀測(cè)網(wǎng)的實(shí)際數(shù)據(jù)情況做了參數(shù)選擇和計(jì)算驗(yàn)證.

2 遺傳算法

遺傳算法(Genetic Algorithm,GA)是模擬達(dá)爾文生物進(jìn)化論的自然選擇的計(jì)算模型,通過模擬自然進(jìn)化優(yōu)勝劣汰的過程搜索最優(yōu)解,是一種啟發(fā)式方法.遺傳算法自上世紀(jì)70年代由Holland和De Jong提出以來,得到了蓬勃發(fā)展,應(yīng)用廣泛,特別是算法的數(shù)學(xué)原理已得到充分研究,其收斂性和全局搜索能力已得到證明[16].

遺傳算法將問題的解向量編碼為染色體,解向量中的每個(gè)元素稱為染色體上的基因.算法模擬生物進(jìn)化的過程,在解的取值范圍內(nèi)隨機(jī)產(chǎn)生G個(gè)染色體作為初始種群(population),對(duì)每個(gè)染色體進(jìn)行適應(yīng)度( fitness)評(píng)估,以衡量解的優(yōu)劣性.通過選擇(selection),交叉(crossover),變異(mutation)遺傳運(yùn)算,使染色體得到進(jìn)化,獲得適應(yīng)度更好的新種群,通過不斷進(jìn)化從而求解出問題的最佳近似解.其計(jì)算流程見圖1.

遺傳算法中種群中個(gè)體的數(shù)目稱為種群數(shù),種群數(shù)越大其搜索能力也越強(qiáng),但計(jì)算效率也越低.用來進(jìn)行適應(yīng)度評(píng)估的函數(shù)稱為適值函數(shù),函數(shù)值稱為適值,即適應(yīng)度越好,越接近最佳近似解.選擇操作是為了從當(dāng)代種群中選出優(yōu)秀的個(gè)體,使其可以作為父代繁衍子代,將優(yōu)秀的基因傳遞下去.選擇依據(jù)適應(yīng)度的大小,適應(yīng)度越高,被選中的概率越大.交叉和變異是遺傳算法中最重要的兩個(gè)運(yùn)算步驟,交叉操作通過對(duì)兩個(gè)父代染色體進(jìn)行組合從而產(chǎn)生新的子代個(gè)體,交叉操作使得優(yōu)秀的基因有機(jī)會(huì)組合在一起,從而獲得更優(yōu)的個(gè)體.對(duì)于通過選擇操作得到的子代,以一定的交叉概率Pc再進(jìn)行交叉操作,交叉操作決定了遺傳算法的全局搜索能力,一般取值大多在0.6～0.9之間.單純通過交叉操作無法產(chǎn)生新的基因,遺傳算法通過變異操作來產(chǎn)生新的基因,以實(shí)現(xiàn)較強(qiáng)的局部搜索能力.對(duì)于經(jīng)過選擇和交叉后的子代各染色體的基因以變異概率Pm進(jìn)行變異操作,與生物進(jìn)化類似,變異發(fā)生概率較低,一般取值在0.01～0.1之間.經(jīng)過選擇、交叉和變異后完成一代的進(jìn)化.如此迭代進(jìn)化,遺傳算法將收斂到全局最優(yōu)解.對(duì)于實(shí)際計(jì)算總需要一個(gè)終止條件,一般當(dāng)最優(yōu)解的適應(yīng)度達(dá)到給定閾值或者進(jìn)化代數(shù)(迭代次數(shù))達(dá)到指定的最大代數(shù)作為終止條件.

對(duì)于不同問題,遺傳算法流程是一致的,重點(diǎn)在于結(jié)合具體問題特點(diǎn)設(shè)計(jì)各種運(yùn)算,關(guān)鍵在于避免種群出現(xiàn)早衰,即使得解陷入局部極值,而不是全局最優(yōu)解.

圖1 遺傳算法的流程圖Fig.1 The flowchart of genetic algorithm

3 初軌計(jì)算的遺傳算法

3.1 變量和編碼的選擇

對(duì)于優(yōu)選問題,首先需要明確優(yōu)化變量.對(duì)于軌道計(jì)算問題最容易想到的是以軌道根數(shù)作為優(yōu)選變量,同時(shí)優(yōu)選6個(gè)變量使得參數(shù)空間維數(shù)太高,不利于求解.基于這個(gè)原因,文獻(xiàn)[14]將問題轉(zhuǎn)換為兩個(gè)變量的優(yōu)化問題,優(yōu)選變量選擇為觀測(cè)首末時(shí)刻的斜距.這種選擇必須將優(yōu)選變量結(jié)合首末觀測(cè)資料共同得到軌道,計(jì)算受到首末觀測(cè)資料精度的制約,特別當(dāng)首末資料出現(xiàn)較大誤差時(shí),無法得到合理的結(jié)果.因此我們選擇和觀測(cè)量不相關(guān)的歷元時(shí)刻t0的3個(gè)Kepler根數(shù)(a,e,M0)作為優(yōu)選變量,這樣參數(shù)空間維數(shù)不太高,同時(shí)解不和具體觀測(cè)資料相關(guān),便于后續(xù)資料處理.

遺傳算法常采用二進(jìn)制編碼,對(duì)于二進(jìn)制編碼,交叉、變異算子比較容易設(shè)計(jì),但對(duì)于軌道計(jì)算這種對(duì)于精度要求較高的連續(xù)問題,采用二進(jìn)制編碼將使編碼長度過長,因此本文采用了實(shí)數(shù)編碼,這樣無需再做轉(zhuǎn)換[17].

3.2 初始種群的生成和終止條件

為了生成初始種群,首先根據(jù)一些已知信息定義3個(gè)優(yōu)選變量的值域,有a∈[al,au], e∈[el,eu],M0∈[Ml,Mu].不同于傳統(tǒng)優(yōu)選方法需要較為準(zhǔn)確的初值,取值范圍并無特殊要求,可根據(jù)已知信息確定,沒有太多信息時(shí),可將范圍取得大一些,這是容易做到的.初始種群采用均勻隨機(jī)采樣方式建立,對(duì)于每個(gè)優(yōu)化變量,由

可得到一個(gè)取值,其中X為任一優(yōu)化變量,λ為[0,1)上均勻分布的隨機(jī)數(shù).對(duì)于3個(gè)優(yōu)化變量分別重復(fù)G次上述過程,即可得到種群數(shù)為G的初始種群{xi,i=1,2,···,G}.

為了簡化計(jì)算流程,選擇迭代次數(shù)達(dá)到最大進(jìn)化代數(shù)K為終止條件.

3.3 適值函數(shù)

令已知一組觀測(cè)量{ti,αi,δi,i=1,2,···,N},其中(αi,δi)為赤經(jīng)和赤緯.根據(jù)歷元時(shí)刻t0的(a,e,M0),對(duì)于任意時(shí)刻ti,可得到Mi:

其中Zi為ti時(shí)刻目標(biāo)的天頂距.對(duì)任意一對(duì)觀測(cè)時(shí)間(tk,tj),tk＞tj,即可得到(rk,rj)和相應(yīng)的(fk,fj).定義適值函數(shù):

對(duì)于每個(gè)個(gè)體xi,其適值為:

顯然,對(duì)于軌道計(jì)算問題,適值越小越優(yōu).

適值除了表明解的優(yōu)劣,在遺傳算法中最主要的作用是決定每個(gè)個(gè)體被選擇進(jìn)入子代的概率,因此需要對(duì)適值進(jìn)行一定的變換.對(duì)于討論的問題,使得適值越小的個(gè)體有越大的選擇概率,這個(gè)過程稱為適值的標(biāo)定(Scaling)[18]7-21.

遺傳算法中全局搜索能力主要取決于種群的多樣性,基于適值F的選擇方式,容易在進(jìn)化前期使得一些相對(duì)較好的解成為超強(qiáng)個(gè)體,被大量選擇進(jìn)入子代,降低種群的多樣性;而在進(jìn)化后期由于種群相對(duì)集中,個(gè)體間適值差異較小,從而變?yōu)殡S機(jī)搜索.本文選擇基于適值的秩作為選擇操作的依據(jù),即標(biāo)定后的適值和原適值取值無關(guān),只與其排序有關(guān).對(duì)于每代種群,對(duì)每個(gè)個(gè)體計(jì)算其適值Fi后,從小到大排序,令個(gè)體xi的排名為ri.顯然最優(yōu)解的ri=1,最劣解的ri=G.標(biāo)定后的適值取

其中c為常數(shù),用于維持種群的多樣性,可根據(jù)具體問題調(diào)整.c＜1時(shí)可減小各適值間的差異,增加樣本多樣性,使種群中較差個(gè)體的選擇概率不會(huì)太低;而c＞1時(shí)則增大了各適值間的差異,使種群中較差個(gè)體的選擇概率更低.由于標(biāo)定后的適值和具體取值無關(guān),顯然在進(jìn)化的各階段都有助于保持種群的多樣性.

3.4 選擇操作

選擇操作的原則是適值越大的被選擇的概率越大,對(duì)于個(gè)體xi,選擇概率Pi取為:

可見經(jīng)過適值標(biāo)定后,最優(yōu)個(gè)體的選擇概率是確定的,從而可使不太優(yōu)的解也有一定的被選擇的概率.選擇操作根據(jù)選擇概率可采用輪盤賭(Roullete wheel)算法實(shí)現(xiàn),考慮采用一般輪盤賭算法對(duì)于一次具體選擇操作,隨機(jī)因素會(huì)使得較優(yōu)個(gè)體沒有被選擇,為了避免這種情況,采用隨機(jī)均勻采樣方法(Stochastic uniform sampling)實(shí)現(xiàn)選擇操作.隨機(jī)均勻采樣方法等價(jià)于只隨機(jī)轉(zhuǎn)動(dòng)一次輪盤,之后的轉(zhuǎn)動(dòng)按照預(yù)設(shè)的規(guī)則進(jìn)行.

初軌計(jì)算問題是要獲得最優(yōu)解,并不需要種群中所有個(gè)體都達(dá)到最優(yōu),因此在選擇操作后引入精英保留策略,即對(duì)于進(jìn)化中最優(yōu)的Ne個(gè)個(gè)體總是選擇進(jìn)入子代,采用精英保留策略后,每代的最優(yōu)個(gè)體不會(huì)發(fā)生衰退.

3.5 交叉操作

由于采用了實(shí)數(shù)編碼,無法直接使用類似二進(jìn)制編碼的交叉操作,文獻(xiàn)[2]仿照二進(jìn)制單點(diǎn)交叉方法,通過直接交換兩個(gè)父代相關(guān)解向量元素的方法實(shí)現(xiàn)交叉.本文采用了算術(shù)交叉的方法[18]42-44,令和分別是父代(k?1代)選為交叉的兩個(gè)個(gè)體,則子代(k代)的兩個(gè)個(gè)體為:

其中,α為[0,1)上均勻分布的隨機(jī)數(shù).對(duì)于M0∈[0,2π],為角度量,因此可有兩個(gè)方向進(jìn)行交叉,約定總是在劣弧(＜π)上進(jìn)行交叉.從上式可看出,采用的交叉運(yùn)算為凸變換,即

交叉新生成的子代一定滿足約束條件,避免了子代個(gè)體有效性的判斷.

3.6 變異操作

遺傳算法中變異算子決定了解的局部搜索能力,文獻(xiàn)[14]采用了高斯變異.高斯變異的變異范圍較難確定,文獻(xiàn)[14]是直接給定了10 km作為變異范圍.本文為了減少參數(shù),采用了動(dòng)態(tài)變異方法[18]42-44,令xi,j表示第i個(gè)個(gè)體的第j個(gè)基因,變異操作為:

其中β為[0,1)上均勻分布的隨機(jī)數(shù),決定變異發(fā)生的方向.和分別表示第j個(gè)基因的上下界.Δ函數(shù)具體形式為:

γ為[0,1)上均勻分布的隨機(jī)數(shù),k和K分別是當(dāng)前進(jìn)化代數(shù)和最大進(jìn)化代數(shù),b為常數(shù),主要調(diào)整全局搜索能力.從式中可看出Δ(k,x)隨k增大逐漸趨于0,這樣在進(jìn)化早期有較大的搜索范圍,而到后期搜索范圍較小,增強(qiáng)細(xì)調(diào)能力.顯然當(dāng)b＜1時(shí)傾向于全局開拓能力,而b＞1時(shí)則更傾向于盡快進(jìn)入局部搜索.同樣變異算子也可以確保變異后代符合約束條件.

3.7 (i,?,ω)的求解

上述各節(jié)實(shí)現(xiàn)了一個(gè)完整的遺傳算法,經(jīng)過求解后可得到t0時(shí)刻的(a,e,M0),由(5)式可得到ti時(shí)刻的地心位置矢量ri,在此基礎(chǔ)上求解完整的6個(gè)軌道根數(shù)已無任何困難．

仍可采用和上述遺傳算法相同過程求解(i,?,ω).求解中(a,e,M0)作為已知量, (i,?,ω)作為優(yōu)化變量,仍是三變量優(yōu)化問題.由(a,e,M0,i,?,ω)按照二體問題模型可計(jì)算ti時(shí)刻的,只需將適值函數(shù)變?yōu)?

由于已有一定信息,初始種群也無需在全區(qū)域搜索,根據(jù)每對(duì)(ri,rj)可得到一組(i,?),取其平均可得到參考的作為初值范圍選取依據(jù),值域可取

4 數(shù)值驗(yàn)證

根據(jù)上述過程,采用MATLAB1http://www.mathworks.com編寫了求解程序,并選取中國科學(xué)院空間目標(biāo)與碎片光學(xué)觀測(cè)網(wǎng)中的一圈實(shí)測(cè)資料進(jìn)行了數(shù)值驗(yàn)證,資料精度優(yōu)于5′′,采樣頻率1 Hz.為了便于比對(duì),選取弧段開始的不同弧長進(jìn)行了計(jì)算.歷元選擇為第1個(gè)觀測(cè)時(shí)刻,考慮到空間目標(biāo)探測(cè)中絕大多數(shù)目標(biāo)都是近圓軌道(e＜0.003)[1],值域選擇為a∈[1.03,4.15], e∈[0,0.003]以及M0∈[0,2π],對(duì)于a并無利用任何先驗(yàn)信息,取值范圍已覆蓋了整個(gè)中低軌道區(qū)域;初軌計(jì)算中一般軌道面相對(duì)較準(zhǔn),取δi=0.1?,δ?=1.5?.

算法中的參數(shù)采用了常見取值,G=50,K=50,Ne=2,Pc=0.8,Pm=0.05,對(duì)于初軌計(jì)算更側(cè)重全局搜索能力,因此取(9)式中c=0.5,(14)式中b=0.5.分別計(jì)算了10～60 s弧長的結(jié)果,圖2和圖3分別給出了一組計(jì)算過程中從初始種群開始到進(jìn)化結(jié)束,每代最優(yōu)個(gè)體的適值F和軌道半長徑a的結(jié)果.從圖中可看出收斂速度非?？?從第10代起,適值F已經(jīng)很接近最小值;a的變化也類似,第10代以后也基本穩(wěn)定在準(zhǔn)確解附近的較小區(qū)域內(nèi),算法搜索效率較高.

圖2 適值F的收斂過程Fig.2 The convergence process of the fitness value

表1給出了不同弧長下的計(jì)算結(jié)果,其中POD表示由多天多站精密定軌得到的結(jié)果,作為參考標(biāo)準(zhǔn),Laplace為改進(jìn)Laplace方法得到的結(jié)果,GA為本文方法得到的結(jié)果.從表中可見,對(duì)于Laplace方法,40 s弧長時(shí)雖能獲得結(jié)果,但半長徑偏差已近100 km,當(dāng)弧長縮短到20 s時(shí),已無法得到合理的結(jié)果.采用本文方法對(duì)于各種弧長均得到了合理的結(jié)果,在10 s弧長下結(jié)果依然較好.

圖3 軌道半長徑a的收斂過程Fig.3 The convergence process of the semi-major axis a

從遺傳算法的計(jì)算過程可見,計(jì)算結(jié)果依賴隨機(jī)數(shù)的生成.為了驗(yàn)證隨機(jī)數(shù)對(duì)計(jì)算的影響,避免隨機(jī)數(shù)對(duì)極短弧定軌結(jié)果產(chǎn)生普適的影響,采用不同的隨機(jī)數(shù)種子進(jìn)行了多組計(jì)算.圖4為10 s弧長采用不同隨機(jī)數(shù)種子計(jì)算100次得到的軌道半長徑直方圖.從結(jié)果看,采用不同隨機(jī)數(shù)得到的結(jié)果有一定差異,但這并不是遺傳算法不收斂的體現(xiàn),而是初軌計(jì)算問題本身病態(tài)的體現(xiàn).圖4中a仍在合理區(qū)間內(nèi),遠(yuǎn)優(yōu)于極短弧定軌本身能夠達(dá)到的精度.文獻(xiàn)[19-20]給出的35 s弧長的定軌精度只有79 km,而70 s弧段的定軌精度也僅為16 km.表2和表3分別給出了20 s和30 s弧長計(jì)算中的10組結(jié)果,也說明了這點(diǎn). 10組計(jì)算結(jié)果有一定差異,但最終殘差都非常小,相互差異小于0.1′′,在精度范圍內(nèi)這些解應(yīng)當(dāng)是等價(jià)的.對(duì)實(shí)測(cè)資料的短弧定軌問題而言解并不唯一,極短弧情況下更為突出,遺傳算法的全局搜索能力使得不同組只是得到了不同的最優(yōu)解.在實(shí)踐中可進(jìn)行多組計(jì)算,選擇適值相對(duì)較小或者多組平均結(jié)果作為定軌結(jié)果.

為了評(píng)估初軌計(jì)算精度,文獻(xiàn)[19-20]引入Bootstrap方法實(shí)現(xiàn)了不依賴其他信息的精度估計(jì)方法,方法結(jié)合經(jīng)典迭代過程的特點(diǎn),通過對(duì)殘差的重采樣構(gòu)建偽觀測(cè)量進(jìn)行精度評(píng)估.而遺傳算法本身的特點(diǎn)使得殘差重采樣方法不再有效,因此采用參數(shù)化的Bootstrap方法對(duì)計(jì)算精度進(jìn)行估計(jì)[21],在原始觀測(cè)資料上根據(jù)測(cè)量誤差分布來構(gòu)建偽觀測(cè)量用于精度評(píng)估.對(duì)觀測(cè)資料增加5′′的隨機(jī)誤差后采用本文方法定軌,采用不同隨機(jī)數(shù)種子對(duì)上述過程重復(fù)100次.由于原始觀測(cè)資料中已有誤差,增加5′′隨機(jī)誤差后得到的偽觀測(cè)量誤差應(yīng)大于5′′,因此結(jié)果可作為5′′誤差資料定軌精度的上界.由于不同隨機(jī)數(shù)種子的使用,結(jié)果中包含了資料誤差和隨機(jī)數(shù)對(duì)求解過程的兩重影響因素,可全面反映求解的實(shí)際精度.由于定軌精度與弧長密切相關(guān),選擇了最短的10 s弧長進(jìn)行數(shù)值計(jì)算,得到的均值和標(biāo)準(zhǔn)差列于表4.結(jié)果表明采用本文方法定軌結(jié)果的精度特征與經(jīng)典方法相仿,軌道面精度略優(yōu),沿跡量精度略差.對(duì)于10 s弧長的定軌結(jié)果,半長徑a精度優(yōu)于10 km,該精度已可為后續(xù)軌道改進(jìn)以及軌道識(shí)別提供一個(gè)有效的初值.

表1 各弧長下的定軌結(jié)果Table 1 The results of orbit determination for various arc lengths

5 結(jié)論與討論

綜上所述,本文方法克服了經(jīng)典方法對(duì)于極短弧定軌的一些不足,能夠應(yīng)用于極短弧定軌,計(jì)算無需先驗(yàn)初值,具有普適性.本文采用的遺傳運(yùn)算和參數(shù)選擇與遺傳算法的常規(guī)選取保持一致,保證了其求解的穩(wěn)定性.

遺傳算法是基于種群的算法,其計(jì)算效率較經(jīng)典方法大大降低,當(dāng)測(cè)量數(shù)據(jù)弧段較長時(shí),經(jīng)典方法已可獲得很好的結(jié)果,遺傳算法毫無優(yōu)勢(shì).但對(duì)于極短弧而言,計(jì)算效率上的差異就當(dāng)前計(jì)算能力而言已可忽略,在經(jīng)典方法失效時(shí),仍可獲得有效的定軌結(jié)果,為之后的工作提供初值.

本文將定軌問題轉(zhuǎn)換為分層優(yōu)化問題,分解為兩個(gè)三參數(shù)優(yōu)化問題,在優(yōu)化變量僅增加一個(gè)的基礎(chǔ)上,解不再依賴具體的觀測(cè)量,提高了方法對(duì)于資料誤差和野值的容錯(cuò)能力,也便于各種資料處理方法的應(yīng)用.

此外,在確定(a,e,M0)時(shí)放棄了對(duì)軌道面的約束,可使得對(duì)殘差較為敏感的沿跡方向相關(guān)量的計(jì)算精度有所提高,由于初軌計(jì)算中得到較為準(zhǔn)確的a是關(guān)鍵,這種處理對(duì)于定軌結(jié)果的后續(xù)應(yīng)用顯然是有幫助的.

圖4 不同隨機(jī)數(shù)種子下的半長徑的分布Fig.4 The histogram of the semi-major axis with di ff erent random seeds

表2 采用不同隨機(jī)數(shù)種子下20 s弧長的軌道根數(shù)Table 2 The results of 20-second arc length with di ff erent random seeds

在初軌實(shí)踐中固定根數(shù)等處理都是十分有效的手段[22],在原有方法中引入這些約束需要在前幾次迭代中得到合理的結(jié)果或者較好的初值,這些對(duì)于極短弧定軌是很難滿足的,而采用遺傳算法后這些處理都可以很自然地融入其中,沒有任何困難,大大方便了實(shí)際使用.

表3 采用不同隨機(jī)種子下30 s弧長的軌道根數(shù)Table 3 The results of 30-second arc length with di ff erent random seeds

表4 軌道根數(shù)的均值和標(biāo)準(zhǔn)差Table 4 The mean and standard errors of elements

[1]吳連大.人造衛(wèi)星和空間碎片的軌道和探測(cè).北京:中國科學(xué)技術(shù)出版社,2011:158-191

[2]Escobal P R.Methods of Orbit Determination.New York:John Wiley&Sons,1965:239-289

[3]Ta ffL G.AJ,1984,89:1426

[4]賈沛璋,吳連大.天文學(xué)報(bào),1997,38:353

[5]陸本魁,戎鵬志,吳建民,等.宇航學(xué)報(bào),1997,18:1

[6]賈沛璋,吳連大.天文學(xué)報(bào),1998,39:337

[7]Jia P Z,Wu L D.ChA&A,1999,23:384

[8]陸本魁,李劍峰,馬靜遠(yuǎn),等.宇航學(xué)報(bào),1999,20:14

[9]劉林,王歆.天文學(xué)報(bào),2003,47:175

[10]Liu L,Wang X.ChA&A,2003,27:335

[11]Vallado D A,Carter S S.JAnSc,1998,46:195

[12]Milani A,Knezevic Z.CeMDA,2005,92:1

[13]Tommei G,Milani A,Rossi A.CeMDA,2007,97:289

[14]Ansalone L,Curti F.AdSpR,2013,52:477

[15]Hinagawa H,Yamaoka H,Hanada T.AdSpR,2014,53:532

[16]De Jong K A.Evolutionary Computation—A Uni fied Approach.Cambridge:MIT Press,2006:26-27

[17]Michalewicz Z.Genetic Algorithm+Data Structures=Evolution Programs.Berlin:Springer,1992: 97-105

[18]玄光男,程潤偉.遺傳算法與工程優(yōu)化.北京:科學(xué)出版社,2000

[19]王歆.天文學(xué)報(bào),2013,54:73

[20]Wang X.ChA&A,2013,37:338

[21]Efron B,Tibshirani R.An Introduction to the Bootstrap.New York:Chapman&Hall,1993:53-56

[22]賈沛璋,吳連大.紫金山天文臺(tái)臺(tái)刊,1997,16:175

Genetic Algorithm for Initial Orbit Determination with Too Short Arc

LI Xin-ran1,2,3WANG Xin1,2
(1 Purple Mountain Observatory,Chinese Academy of Sciences,Nanjing 210008)
(2 Key Laboratory for Space Object and Debris Observation,Purple Mountain Observatory,Chinese Academy of Sciences,Nanjing 210008)
(3 University of Chinese Academy of Sciences,Beijing 100049)

The sky surveys of space objects have obtained a huge quantity of tooshort-arc(TSA)observation data.However,the classical method of initial orbit determination(IOD)can hardly get reasonable results for the TSAs.The IOD is reduced to a two-stage hierarchical optimization problem containing three variables for each stage. Using the genetic algorithm,a new method of the IOD for TSAs is established,through the selection of optimizing variables as well as the corresponding genetic operator for speci fic problems.Numerical experiments based on the real measurements show that the method can provide valid initial values for the follow-up work.

space vehicles,celestial mechanics,methods:numerical

P135;

:A

10.15940/j.cnki.0001-5245.2016.01.007

2015-07-16收到原稿,2015-08-21收到修改稿?國家自然科學(xué)基金項(xiàng)目(11373072)資助

?lixr@pmo.ac.cn

?wangxin@pmo.ac.cn