低速濺射小行星表面采樣的動力學仿真*

程 彬，寶音賀西

(清華大學 航天航空學院·北京·100084)

0 引 言

人類自古以來就對浩瀚的宇宙充滿好奇，渴望向未知領域探索求知。自20世紀50年代以來，人類探索宇宙的腳步不斷深入，深空探測項目在近幾十年來已成為各航天大國的重要發展方向之一，而小行星采樣探測更是其中不可或缺的一部分。小行星的表面保留著太陽系形成初期的原始成分，是研究太陽系起源和演化歷史的“活化石”，其表面的撞擊坑更是研究太陽系爆炸時期環境變化的最佳材料[1]；同時，其上有豐富的稀有礦物資源。例如，某顆直徑僅為500m的小行星蘊含的鉑元素就相當于迄今為止人類鉑開采量的總和，而其表面的水資源還可作為人類深空探測旅行的補給站[2]；另外，部分近地小行星有撞擊地球的可能，通過深空探測了解其結構組成，對于未來可能的小行星防御計劃的成功具有非常重要的作用[3]。

目前，各航天大國先后提出并實施了小行星探測規劃，以占據航天技術戰略制高點，從而在未來的太空資源爭奪中保持優勢。日本宇航局(JAXA)于2003年發射了隼鳥號(Hayabusa)，在近10年間針對S-類小行星Itokawa進行了近距離抵近探測與采樣返回，成為了世界上第一個對小行星成功進行采樣返回的國家[4]，而隼鳥二號(Hayabusa II)也已經于2018年6月到達了小行星Ryugu，并成功實施了表面巡游與撞擊采樣。美國宇航局(NASA)于2016年9月發射了OSIRIS-REx深空探測器，探訪了富含碳質生命元素的C-類小行星Bennu[5]。OSIRIS-REx深空探測器在2018年10月到達了小行星Bennu，將在2019年實施表面附著采樣。歐洲航天局(ESA)計劃于2022年實施Hera小行星探測任務。我國小行星環繞與采樣返回任務也正在規劃論證中，國內相關研究部門已經積極開展了針對小行星探測的相關研究工作。

考慮到小行星表面引力非常低，且表面覆蓋著顆粒狀風化層，不易于實現探測器的著陸、錨定，因此在微重力條件下，通常采用與傳統采樣方案不同的策略，即接觸-分離(Touch-And-Go)式采樣，利用探測器接觸小行星表層的短暫時間，通過收集濺射物質實施采樣。接觸時間往往僅有幾秒，之后需立刻啟動反沖發動機以離開小行星表面，從而不需要著陸、錨定機構。如探訪小行星Itokawa的隼鳥號，在采樣裝置接觸到小行星表層時發射高速射彈，撞擊其浮土層，向上飛濺的巖石碎片在擋板的作用下反彈入返回艙，整個收集過程僅持續了1s，并成功收集到了1534個風化層顆粒[6]。對于小行星采樣方案而言，最重要的是收集效率。在隼鳥號發射前，JAXA曾進行了大量實驗，以研究射彈形狀、射彈速度等參數的影響[7]。但受限于實驗條件，JAXA并沒有進行更接近小行星真實環境的超微重力實驗。得益于數值方法的特性，我們可以更為方便地模擬小行星的微重力環境，研究更為廣泛的參數空間的采樣過程，同時跟蹤風化層每個顆粒的運動，更深入地理解濺射過程，從而優化采樣收集效率，為高效采樣儀器的設計提供理論基礎。針對這一問題，Wada等[8]首先將軟球離散元法引入小行星表面風化層的撞擊仿真中，研究了撞擊導致的碎片速度、角度分布及撞擊坑尺寸與撞擊速度的關系；Nakamura等[9]在此工作的基礎上進行了微重力實驗，研究了射彈在低速穿過顆粒層時的速度衰減過程，發現離散元法具有較高的準確性；Schwartz等[10]首先使用離散元法進行了采樣過程仿真，研究了濺射質量與射彈形狀的關系，并將其與地球環境下的實驗進行了對比，發現離散元方法可以較為準確地復現實驗結果。但是，由于微重力環境下的仿真耗時較長，Schwartz等并沒有進行接近小行星環境的采樣過程仿真，也沒有考慮采樣器的復雜幾何結構。

本文采用數值方法研究了小行星表面的撞擊采樣過程，考慮了小行星的微重力環境與顆粒風化層等特殊地質條件，基于初步的濺射采樣裝置原型機，研究了射彈形狀對撞擊坑形狀的影響，以及其與采樣效率的依賴關系，為我國未來的小行星采樣裝置設計提供了理論基礎。

1 數值模型

現有天文觀測與采樣結果表明，大部分小行星表面覆蓋著顆粒狀的風化層，在微重力環境下其動力學行為復雜。目前針對顆粒物質的研究主要采用離散單元法，課題組基于此方法開發了三維大規模顆粒模擬軟件DEMBody[11-12]。該軟件以牛頓第二定律為基礎，所研究的對象是大量具有特定形狀的離散顆粒，其采用顆粒領域通用的軟球模型來描述散體的動力學行為，如圖1所示。其包括彈簧阻尼、Hertz接觸力、Cundall切向力等模型，根據粒子交疊的程度刻畫顆粒之間的正向力和摩擦力，即有

FN=-knxn+Cnun
FT=min{μ|FN|t;ktS+Ct|ut|t}

(1)

式(1)中，FN和FT分別為顆粒所受到的正向彈性力與摩擦力；kn和kt為彈性系數，Cn和Ct為阻尼系數，與顆粒碰撞恢復系數相關，在顆粒系統仿真中一般取經驗值0.5；無量綱參數μ為顆粒間的摩擦系數；x為顆粒的正向侵入量，S為接觸過程中顆粒的切向形變矢量；n和t分別為正向與切向的單位矢量；正向力與摩擦力中均包含阻尼力項，其分別與正向和切向相對速度(即un和ut)成正比，用于描述顆粒接觸過程中的黏滯耗散行為。

在之前的顆粒系統仿真過程中，為了計算簡便，往往會忽略滾動摩擦力矩的作用。已有研究[13]表明，這種簡化將會導致顆粒系統行為與真實行為出現較大差異。因此，引入了滾動摩擦系數來代表該動能損耗。當滾動摩擦系數不為0時，將會抑制顆粒對在接觸點的滾動相對速度，從而阻礙顆粒對的相對滾動，即有

MR=-μr|FN|Riωi

(2)

式(2)中，MR表示顆粒所受到的滾動力矩；Ri為顆粒半徑，ωi為顆粒角速度；無量綱參數μr為顆粒間的滾動摩擦系數。

(a)顆粒接觸正向模型

(b)顆粒接觸切向模型

(c)顆粒接觸滾動模型圖1 顆粒接觸的離散元模型Fig.1 The contact model in Discrete Element Model simulation

為了模擬數百萬尺度的顆粒動力學問題，采用并行計算的方式解決大規模計算問題。DEMBody采用OpenMP與MPI進行多核并行，可在天河二號超級計算機上進行大規模計算，在實際測試過程中可數百倍地提高效率，基本滿足本研究的需求。

2 數值程序測試

本部分通過與典型的顆粒動力學實驗進行對比，確定了程序的正確性與適用性。LI等[14]通過實驗，研究了在飛船與天體表面接觸時，產生濺射物質的質量與撞擊速度的關系。通過仿真重現了該實驗的過程，即使用直徑為19.71mm的玻璃射彈撞擊由丙烯酸顆粒填充的顆粒床，撞擊速度范圍為0.9m/s～3.6m/s，顆粒床的厚度范圍為46.1mm～138.3mm，直徑保持為75.6mm。通過記錄撞擊前后顆粒床的質量之差，可確定由撞擊引起的濺射質量。

圖2所示為仿真與實驗結果的比較，其中圓點表示仿真結果，叉號表示實驗結果，誤差線表示由多次重復實驗得到的濺射質量范圍；藍色數據表示46.1mm顆粒床，黃色數據表示92.2mm顆粒床，紅色數據表示138.3mm顆粒床。仿真結果表明，濺射質量隨著撞擊速度的增加而逐漸增加，這是由于撞擊動能的增加導致了顆粒床所獲得的動能增加；同時，濺射質量隨著顆粒床厚度的增加而逐漸減小，這是由于顆粒床厚度的增加使得顆粒系統的摩擦能量損耗增加，從而導致顆粒所獲得的動能減小。值得注意的是，在每一組仿真參數中，仿真得到的濺射質量都基本在實驗誤差范圍內，這表明程序可以很好地復現顆粒床的低速撞擊過程，具有較高的準確性。

圖2 仿真濺射質量與實驗結果的對比Fig.2 Comparison of ejected mass between simulation and experiment

3 仿真試驗及結果分析

3.1 仿真設置

日本隼鳥號對小行星Itokawa風化層的觀測表明，小行星風化層的力學性質接近碎石，顆粒半徑滿足0.2mm～0.25mm的立方冪律分布。由于隼鳥號的探測數據是目前小行星風化層最為精細的數據，因此仿真中采用的類似S-類小行星風化層的參數，如表1所示。

對于射彈，采用了與隼鳥號類似的SUS304鐵的物理參數，即密度為8.0g/cm3；同時，為了研究射彈形狀對濺射采樣過程的影響，采用了相同質量不同形狀的射彈，包括球形、60°錐形、90°錐形、120°錐形、半球形、1∶2圓柱、1∶1圓柱、2∶1圓柱，具體幾何參數如表2所示。其中，R表示射彈的半徑，H表示射彈的總厚度，M表示射彈的質量。

表1 濺射采樣仿真中所使用的接觸參數

表2 濺射采樣仿真中所使用的射彈形狀

對于采樣裝置，采用了與隼鳥號類似的結構，包括可伸縮采樣桶(Extendable Sampler Horn)、錐形角(Conical Horn)、捕集器(Sampler Catcher)與收集倉(Container)，分別使用了圓柱、圓錐、圓柱、平板邊界模擬。如圖3為采樣裝置的尺寸說明。

圖3 采樣裝置、部分射彈形狀及采樣原理示意圖Fig.3 Sampling equipment, several projectile shapes and sampling principle

在仿真中，首先將13萬顆粒隨機撒入直徑為28cm的圓柱容器中，在小行星重力作用下形成模擬風化層。隨后，將采樣裝置放置于顆粒層上方，同時將射彈以25m/s的速度垂直向下射出，設定采樣時間為1s。

3.2 仿真結果

不同形狀的射彈侵入風化層后，均在其后形成不斷擴大的撞擊空腔，同時撞擊物將攜帶的動能通過力鏈的方式傳遞給顆粒層，導致部分顆粒從風化層表面濺起，形成了倒圓錐狀的濺射云。如圖3(右)所示，這些濺射物質將撞擊采樣桶，隨后反彈并沿采樣桶不斷上升，直至到達收集倉。需要指出的是，由于采樣時間設置為1s，即若顆粒在1s內到達收集倉，才能認為該顆粒被成功采樣。

如圖4所示，對于頭部尖銳的射彈，其形成的撞擊坑近似為錐形，撞擊坑更深，但寬度小，如60°錐形、1：2圓柱形；對于頭部圓鈍的射彈，其形成的撞擊坑近似為U形，撞擊坑更淺，但寬度大，如120°錐形、2：1圓柱形；對于介于其中的射彈，其形成的撞擊坑近似為碗形，撞擊坑的深淺與寬度也介于兩者之間，如90°錐形。其原因在于，不同形狀射彈的截面積與接觸風化層面積均不同。侵入時，所受接觸力不同，對風化層的沖擊作用也不同，60°錐形射彈所受阻力較小，侵入深度大，易形成錐形撞擊坑；120°錐形射彈所受阻力較大，侵入深度小，易形成U形撞擊坑；90°錐形射彈所受的阻力介于兩者之間，因而易形成碗形撞擊坑。

圖4 不同形狀射彈撞擊得到的撞擊空腔。從上至下、從左至右分別為球形、60°錐形、90°錐形、120°錐形、半球形、1∶2圓柱、1∶1圓柱、2∶1圓柱Fig.4 The impact cavity generated by projectiles with different shapes. From left to right then top to bottom: sphere, 60° cone, 90° cone, 120° cone, hemisphere, 1∶2 cylinder, 1∶1 cylinder, 2∶1 cylinder

不同形狀射彈在撞擊后所形成的撞擊坑形狀、顆粒速度分布、濺射區域均有不同，從而影響了收集倉采樣的效率。如圖5所示，對于圓柱形射彈而言，其采樣質量與射彈的尖銳程度近似成反比，即射彈越尖銳(如1：2圓柱形射彈)，采樣質量越小；但對于錐形射彈而言，其采樣質量與射彈的尖銳程度近似成正比，即射彈越圓鈍(如120°錐形射彈)，采樣質量越小。在所有被研究的射彈形狀中，90°錐形射彈的采樣質量最高，60°錐形射彈相比略低。因此，我國未來的小行星采樣任務可以將90°錐形射彈作為首選。

圖5 不同形狀射彈撞擊收集到的風化層質量 Fig.5 The sampling mass collected by projectiles with different shapes

4 結 論

本文將散體動力學仿真方法應用到了適合于微重力小行星環境的射彈濺射采樣優化過程中，重點研究了同等質量的不同形狀射彈對彈坑形狀、采樣效率的影響，涉及到了不同角度的圓錐形、半球形、圓柱形等形狀，比較得出的結論是90°錐形射彈兼顧了提高濺射速度與濺射質量兩方面的因素，采樣效率最高。因此，我國在未來的小行星探測器中可首先考慮采用90°錐形射彈。但由于本文采用的離散元法未能考慮顆粒破碎、熱量傳遞能問題，在后續的研究中需結合光滑粒子動力學法(SPH)等沖擊動力學模型對這一問題繼續進行研究。