小行星(99942)Apophis的軌道改進與Yarkovsky效應測定

田偉

(上海歷元信息科技有限公司,上海 200241)

1 引言

近地小天體 (99942)Apophis于2004年6月19日首先被R.A.Tucker、D.J.Tholen和F.Bernardi發現,因最早期的軌道分析顯示該顆小行星在2029年4月13日撞擊地球的概率超過了2%,得到了廣泛關注[1]。隨著對Apophis觀測資料的不斷積累,尤其是高精度的雷達數據,研究人員發現該顆小行星撞擊地球的概率遠小于早期預報[2-4]。在2021年3月Apophis再次飛掠地球時,研究人員利用這次飛掠期間的觀測資料以及歷史資料已基本排除2068年之前Apophis撞擊地球的可能性[5]。

類似于Apophis,很多近地小天體會多次飛掠地球。由于受地月系統的影響,天體軌道和自轉狀態都會受到不同程度的改變。Brozovic等人[4]的研究結果顯示2029年Apophis飛掠地球期間,其日心半長徑將從0.92AU變化為1.10AU(相應地,其軌道類型將從Aten型變化為Apollo型)。Souchay等人[6]的研究結果表明Apophis在2029年飛掠地球期間,在地月系統的引力作用下,其自轉狀態也將發生顯著改變。因此,在利用實測數據確定Apophis軌道之前,本文將定量分析2004-2022年期間幾次飛掠地球過程中太陽系歷表的誤差和相關動力學參數的不確定度對Apophis軌道的影響。圖1列出了觀測弧段期間Apophis飛掠地球、月球和金星的概況。

圖1 1994-2022期間,小行星 (99942)Apophis飛掠地球、月球和金星的概況Fig.1 Overview of encounters of asteroid(99942)Apophis with the Earth,Moon and Venus from 1994 to 2022

2 觀測資料

國際小行星中心 (MPC)收集了關于Apophis的地面光學觀測和無線電雷達測量數據。文中采用了截至2021年5月20日該網站公布的所有光學數據,以及2005-2021年期間由美國Arecibo和Goldstone雷達站獲取的高精度測距和多普勒頻移觀測數據[4]。表1給出了部分觀測數據的統計信息。

表1 小行星 (99942)Apophis(部分)地面光學 (赤經/赤緯)和雷達 (測距/多普勒頻移)觀測數據的統計Table 1 Statistics of partial ground-based optical observations(right ascension and declination)and radar observations(range and doppler frequency shift)of asteroid(99942)Apophis

3 動力學模型

文中考慮的天體包括太陽、八大行星、月球和冥王星,以及48顆大質量的主帶小行星。動力學模型考慮的因素包括:(1)點質量相互作用:EIH方程; (2)非球形引力場:太陽、地球和月球的低階引力場;(3)主帶小行星的攝動效應; (4)Yarkovsky效應對Apophis的攝動影響等。

3.1 點質量相互作用

因中心天體太陽與Apophis的相互作用中相對論效應十分顯著,需要考慮點質量相互作用中的一階后牛頓近似下的相對論效應。太陽、月球、八大行星、冥王星對小行星 (99942)Apophis的引力作用由Einstein-Infeld-Hoffmann(EIH)方程給出[7]:

式中,下標A表示Apophis,GMB和GMC分別是天體B和C的質量與引力常數G的乘積,rAB=rB-rA和rAB=|rAB|分別是天體A和B之間的相對位置矢量和相對距離,rA(或rB)和vA(或vB)分別為對應天體A(或B)的位置和速度矢量。相對位置矢量rAC和rBC的定義類似于rAB,而相對距離rAC和rBC的定義類似于rAB,以上各量均在太陽系質心廣域參考系中給出。c=299792458m·s-1為光子在真空中的傳播速度。一階后牛頓參數β=γ=1。

3.2 非球形引力場

除考慮主要天體 (作為點質量)的引力作用外,文中也考慮了地球、月球和太陽的非球形引力場對Apophis軌道的影響。非球形引力場的引力作用取決于Apophis與非球形天體之間距離。通常近地小行星會 (近似)周期性地與地球或金星交會 (見圖1)。以地球為例,其引力場高階項 (J2,Cnm,Snm)對某外部天體的影響可由式 (2)表示[8]

式中,在ξ-η-ζ坐標系中,ξ坐標軸代表從地球到Apophis的指向,ξ-η平面包含了地球自轉軸,ζ軸與ξ和η軸形成右手坐標系。GME和R分別代表了地球的質量與引力常數乘積和平均半徑。表示n階m次的締合Legendre函數。n1和n2分別是Jn和 (Cnm,Snm)的最高階數。(r,φ,λ)是從地球到Apophis構成的位置向量在地固坐標系中的球面坐標。需要說明的是,在式 (2)的基礎上計算太陽系質心坐標系下的相互作用力時,需要實施一系列坐標轉換。太陽、地球和月球的非球形引力場分別取J2項、3階以下項和2階以下項。限于篇幅,文中未考慮金星非球形引力場對Apophis的影響。

3.3 主帶小行星攝動模型

某一特定小行星對Apophis的攝動效應取決于攝動小行星的質量和相對距離。目前被精確測量質量參數的小行星在已發現的小行星總數中占比較小,因而小行星攝動效應是動力學模型中一個難以確定的因素。Giorgini等人[2]統計了約37.3萬顆已知小行星對Apophis軌道的攝動效應,并給出了4顆攝動影響最大的小行星:Ceres、Vesta、Pallas和Hygiea,其中Ceres、Pallas和Vesta是主帶小行星中質量最大的3顆。文中未對單顆小行星的攝動效應做逐個篩選,而是選取了48顆大質量小行星作為攝動小天體。

3.4 Yarkovsky效應

由于目前對Apophis的觀測主要依賴于地面設備,其形狀、自轉軸指向、平均密度、反照率、熱慣量等基本幾何和物理參數確定精度偏低,在計算Apophis受到的Yarkovsky效應時,文中采用了簡化的非保守力線性模型[3,9-10],

式中,標準參數A1、A2和A3分別代表與徑向、橫切線方向和軌道面法線方向加速度分量相關的參數,r和r分別為小行星相對太陽的位置向量和距離,h=r×r˙。文中僅考慮了橫切線方向的作用力,即A1=A3=0,A2≠0。線性模型中的非零參數A2將與軌道參數一起解算。

4 部分太陽系歷表參數對Apophis軌道的影響

作為太陽系內主要天體,太陽、地球和月球的質量不確定度對Apophis軌道有著不同程度的影響。目前太陽、地球和月球的質量參數可以通過不同技術手段精確確定[11-14],其數值和不確定度見表2。圖2給出了相應天體質量不確定度對Apophis軌道的攝動效應。

表2 太陽、地球和月球的質量參數 (單位:km3/s2)Table 2 The mass parameters of the sun,Earth and Moon(Unit:km3/s2)

圖2 由太陽、地球和月球的質量參數不確定度引起的小行星 (99942)Apophis軌道變化 (單位:km)Fig.2 Orbit changes of asteroid(99942)Apophis induced by uncertainties of mass parameters for the sun,Earth and Moon(Unit:km)

在飛掠大行星過程中,被飛掠天體的歷表位置誤差是Apophis軌道確定的一個誤差項。其影響程度取決于最近飛掠距離,圖3給出了2004-2022年期間,地球的歷表位置誤差 (假設位置誤差為1km)對Apophis軌道的攝動效應。相比于與地球歷表誤差的攝動效應,月球的影響較小。

圖3 地球歷表誤差對Apophis軌道的影響 (假設位置誤差為1km)(單位:km)Fig.3 Influences on orbit of asteroid(99942)Apophis by ephemeris error of the Earth(the supposed position error is 1 km)(Unit:km)

文中考慮了主要天體非球形引力場的影響,其中太陽、地球和月球的J2項是主要影響項。與被飛掠天體的歷表位置誤差類似,其非球形引力場的引力攝動效應與最近飛掠距離直接相關。在2004-2022年期間飛掠地球的最近距離約為0.1AU,地球和月球的J2項對Apophis軌道的影響可完全忽略 (其量級均小于10m)。太陽J2項的影響與飛掠事件無關,在百年尺度上該效應可以對Apophis軌道產生公里量級的攝動。

值得注意的是,當Apophis的軌道穿越地月系統時,比如2029年Apophis飛掠地球的最近距離約為5個地球半徑,地球和月球的歷表誤差和非球形引力場的影響將有量級上的變化,需要進一步研究太陽系歷表誤差對Apophis飛掠前后運動狀態的影響。

5 數據分析與定軌結果

在分析光學和雷達觀測資料時,我們采用了經典最小二乘原理下的批處理方法[15,16]。觀測方程和法方程的建立步驟參閱文獻[17]。考慮到部分數據存在由不同原因引起的系統性偏差,我們采用了如下的定權方法:首先將全部觀測數據根據臺站進行分組,計算每組擬合后的均值a0和均方根σRMS;然后將作為該組每個觀測量的測量中誤差;考慮到觀測量N＜10的分組,我們將作為相應觀測量的測量中誤差。雷達觀測資料的權重由觀測文件中的測量誤差信息給出。考慮到觀測數據中由不同原因引起的粗差,本文采用了3σ0的原則對粗差進行剔除,具體剔除數量見表1。

初始軌道在歷元T0=JD2454466.5處的狀態向量由JPL#199提供的軌道根數給出。由于初始軌道較為準確,通過3次迭代便得到收斂數值解(驗前和驗后中誤差的比值趨于1)。解算參數包括Apophis的軌道參數和Yarkovsky效應參數A2。我們采用2004-2021年期間的所有光學和雷達數據得到解 (Sol-A)。圖4和圖5分別給出了部分光學觀測 (在赤經和赤緯方向上)和雷達測距數據的擬合后觀測殘差分布。部分擬合后觀測殘差的加權均方根誤差 (WRMS)在表1中給出。圖6給出了2004-2022年期間Apophis軌道誤差的傳播情況 (9倍形式誤差或9σ)。

圖4 部分光學觀測 (赤經和赤緯方向)的擬合后殘差分布Fig.4 Distribution of the post-fit residuals for partial optical observations(RA and Dec.)

圖5 無線電雷達測距的擬合后殘差分布Fig.5 Distributionof the post-fit residuals for the radio radar ranging measurements

圖6 2004-2022年期間小行星 (99942)Apophis的軌道誤差 (9σ)Fig.6 Orbit errors(9σ)during the period 2004-2022 for asteroid(99942)Apophis

表3給出了不同參考文獻的Yarkovsky效應參數A2估計值。與Pérez-Hernández等人[19]和JPL官方網站[20]給出的結果相比,Sol-A給出的Yarkovsky效應參數值A2更接近Yarkvosky效應的線性模型理論計算值。值得注意的是,解Sol-A 的形式誤差 (1σ) 與 Pérez-Hernández 等人[19]和JPL官方網站[20]給出的結果略有不同,這可能是因采用了不同于后兩者的定權策略而造成的。

表3 不同文獻對小行星 (99942)Apophis的Yarkovsky參數A2的估計 (單位:10-14AU/d2)Table 3 Estimations of Yarkovsky parameter A2from different references for asteroid(99942)Apophis(Unit:10-14AU/d2)

為了考察2021年期間觀測數據對Apophis軌道確定和Yarkovsky效應測定的影響,我們利用2004-2020年期間的所有數據得到另一組解 (Sol-B)。比較表2中解Sol-A和Sol-B的Yarkovsky效應參數值不難發現2021年的觀測數據對Yarkovsky效應約束較強,被估計參數A2的信噪比有明顯提高。

為了檢驗Yarkovsky效應參數A2對數據的依賴性。我們將所有光學和雷達數據分割為四組,每組觀測數據的數量和類型近似等同。利用每組數據進行單獨解算Apophis軌道參數和Yarkovsky效應參數值A2,其中Yarkovsky效應參數值A2的結果見表4。表4中四組解的加權平均值與解Sol-A保持一致。與解 Sol-A中給出的形式誤差 (1σ)相比,表4中加權平均值的WRMS(約9σ)更能反映Yarkovsky效應參數A2的真實不確定度。相應地,我們在圖6中取9σ來評估2004-2022年期間的Apophis軌道誤差。

表4 數據分組后,四組不同解中的Yarkovsky參數A2估計值 (單位:10-14AU/d2)Table 4 Results of Yarkovsky parameter A2from four solutions with sub-grouped observations(Unit:10-14AU/d2)

6 總結與討論

文中初步分析了大行星歷表誤差,以及部分歷表動力學參數誤差對小行星Apophis軌道的影響,發現觀測期間2004-2021年的幾次飛掠,大行星歷表誤差對Apophis軌道的影響較小,不會對軌道參數和Yarkovsky效應參數A2的確定產生顯著影響。大行星歷表誤差的影響程度取決于小行星與大行星的最近飛掠距離。2029年Apophis將以5倍地球半徑的距離飛掠地球,期間大行星歷表誤差的影響將會顯著增強,需進一步研究。

Yarkovsky效應是Apophis軌道確定和預報的主要誤差來源。本文利用2004-2021年期間的數據同時對Apophis的軌道和Yarkovsky效應參數A2進行確定。通過比較Sol-A和Sol-B兩組解,不難發現2021年的觀測數據顯著地提高了Apophis軌道的確定精度和Yarkovsky效應的測量信噪比。考慮到解算結果對觀測數據的依賴性,我們給出了在2004-2022年期間Apophis軌道誤差小于10km(9σ)和Yarkovsky效應參數值A2=-2.73±0.12(對應信噪比約為22.8)的結果。