衛(wèi)星解體碎片生成數(shù)值模擬

張曉天, 賈光輝

(北京航空航天大學 宇航學院, 北京 100191)

0 引 言

航天器解體碎片生成研究是空間碎片環(huán)境建模的關鍵。目前的主要研究手段是對解體碎片數(shù)據進行統(tǒng)計分析，建立經驗公式。解體碎片數(shù)據來源于兩個方面：軌道碎片雷達觀測和地面解體實驗。

美國太空觀測網(SSN，Space Surveillance Network)對地球軌道的物體進行跟蹤和編目。以目前的觀測能力，能夠對低地球軌道特征長度10cm以上的物體和高地球軌道特征長度1m以上的物體進行可靠跟蹤，難以獲得小尺寸碎片的數(shù)據。

在空間碎片群體中，小碎片的數(shù)目要遠大于大碎片。為了研究小碎片的特性，NASA 開展了SOCIT(Satellite Orbital Debris Characterization Impact Test)地面解體實驗[1]。解體實驗中使用了OSCAR 小衛(wèi)星，通過加速鋁彈丸對其碰撞造成解體。日本T. Hanada 研究小組與NASA 約翰遜研究中心合作，開展了復合材料小衛(wèi)星地面解體實驗研究[2]。柳森、黃潔、蘭勝威和李毅等開展了小型衛(wèi)星模型的解體實驗研究[3-5]，進行了三種衛(wèi)星模型的撞擊解體實驗，并提出了新的解體閾值。王海福、余慶波和金學科等開展了圓柱殼結構爆炸解體實驗和數(shù)值仿真研究[6-7]。地面解體實驗受到加速設備能力限制，難以構造常規(guī)尺寸衛(wèi)星解體的碰撞條件，地面解體實驗難度很大。

本文通過超高速碰撞數(shù)值模擬技術對航天器解體過程進行研究，從而獲得解體碎片的特性。數(shù)值模擬技術是一種相對成本低的技術，一旦建立了有效的模擬算法，就可以大量開展。由于目前廣泛使用的SPH方法在單個碎片的識別和置信度方面存在不足，本文通過一種SPH和有限元結合的方法來實現(xiàn)對碎片的識別和統(tǒng)計，并與衛(wèi)星模型地面解體實驗結果進行對比。

1 衛(wèi)星解體仿真方法

目前的超高速撞擊數(shù)值模擬研究主要針對航天器防護結構分析，廣泛使用的數(shù)值模擬方法是SPH(Smoothed Particle Hydrodynamics)方法[8]。但在航天器防護結構仿真中，由于撞擊體較小，航天器一般不會發(fā)生解體；另外這種分析的關注點為防護結構是否被碎片穿透，并不關注撞擊產生的二次碎片云中單個碎片的信息(如碎片數(shù)目、大小和速度矢量等)。而解體模型就是要根據解體條件給出解體碎片特征信息，解體碎片特性是解體動力學數(shù)值模擬最關心的結果。因此，將現(xiàn)有數(shù)值模擬算法成果應用于航天器解體動力學數(shù)值模擬，首先需要對單個碎片的識別性和置信度等方面進行分析。

對于超高速碰撞問題，SPH方法的仿真結果常表現(xiàn)為離散的粒子云，難以從中識別出哪些粒子構成一個碎片，也難以分辨破碎材料和未破碎材料的界面。雖然可以通過搜索算法直接將結果中鄰近的粒子劃分為一個碎片，但是這種方法存在幾個方面的不足。第一，置信度問題。在初始模型中，粒子被以一定間距δ填充入物體所占據的空間。此時粒子作為1個群體代表了1個物體，其可信性較高。但是單個粒子代表了1塊δ立方的物質微元，其作為孤立個體可信性較低。因為該物質微元有可能進一步分裂破碎，而一個粒子不可能再被分裂。所以仿真結果中的孤立粒子，甚至是少量幾個粒子都不能作為置信的單個碎片存在。第二，流動性問題。超高速碰撞過程中，物質流動范圍很大，初始不相鄰的物質微元完全有可能在撞擊過程中被“擠”到一起。因此僅僅根據某個時間點上粒子的空間坐標是否相鄰來判斷這兩個物質微元是否屬于同一個碎片是不夠可靠的。第三，鄰近閾值問題。判斷粒子是否鄰近，一般是判斷兩者間距是否小于給定閾值。這個閾值的給定缺乏力學依據，難以制定統(tǒng)一的規(guī)則，這會影響不同問題之間的一致性。

SPH方法在航天器解體碎片仿真問題中存在上述不足；另一方面有限元方法材料界面清晰，碎片效果直觀，單個碎片置信度高，但是在大變形問題中會發(fā)生網格畸變，使得計算無法進行。將SPH與有限元方法相結合的方法，稱為有限元重構方法，提高了對單個碎片的識別性和置信度。在有限元重構方法中首先通過立方體填充，將任意形狀的幾何體離散為粒子-形心重合模型；進而通過顯式積分迭代步中的SPH積分、有限元重構與失效分析三個步驟，獲得包含F(xiàn)E單元與SPH粒子的混合仿真結果。

1.1粒子-形心重合模型

圖1是一個離散化的實例。A代表一個任意形狀的物體，B表示使用立方體FE單元對物體區(qū)域進行填充，C表示在每個單元的形心處創(chuàng)建一個SPH粒子，相應的粒子間距等于FE的單元邊長。此時的模型稱為“粒子-形心重合模型”，包括FE網格與SPH粒子兩部分。

A B C D E

1.2顯式積分迭代步

對于粒子-形心重合模型中的粒子部分，按照SPH控制方程進行顯示時間積分迭代，在每個時間步后能夠獲得所有粒子的位置、速度信息。

AEi={EH|Ni∈EH}

(1)

在立方體填充模型中，根據節(jié)點被多少個單元所共有——也就是AEi集合中單元的個數(shù)——共分為8種情況。如果1個節(jié)點被8個單元共有，則稱其為內點，否則稱其為邊界點。有限元重構分為3個步驟依次執(zhí)行：內點重構、第一類邊界點重構，第二類邊界點重構。

(1)內點重構

用IN表示立方體填充網格中的內點集合，它包含了所有被8個單元共享的節(jié)點。內點的重構坐標定義為這8個單元對應的SPH粒子當前坐標的算術平均：

(2)

其中，pH表示EH形心對應的SPH粒子的坐標向量，ri表示Ni的坐標向量。

(2)第一類邊界點重構

定義一個與Ni關聯(lián)的單元標識集合：

(3)

第一類邊界點集合定義為：

(4)

第一類邊界點的重構坐標按下式進行計算：

(5)

(3)第二類邊界點重構

定義一個與Ni關聯(lián)的單元標識集合：

(6)

第二類邊界點集合定義為：

(7)

第二類邊界點的坐標按下式進行計算：

(8)

如果一個邊界點所有的背點都是第二類邊界點，將陷入一個死循環(huán)，這些點無法被更新。為避免這種極端情況，在幾何體離散的過程中需要對最薄的部位采取兩層以上的單元來填充。

在大變形問題中，SPH粒子流動范圍很大，重構后的單元不一定能維持常規(guī)的形狀，很多單元會發(fā)生畸變，此時失效單元將被刪除，單元形心對應的SPH粒子被保留；未失效單元被保留，對應的SPH粒子被刪除。本文采用Johnson-Cook失效準則。FE重構分析和失效分析按照SPH計算的時間步逐步進行，第n+1個時間步的分析繼承第n個時間步的結果。在結果中FE單元代表未失效材料，SPH粒子代表失效材料。

2 兩體碰撞解體算例

使用上述方法對實心球體撞擊實心立方體解體問題進行仿真，兩物體相對位置關系如圖2所示，從右上角局部放大圖中可以看出粒子和單元形心的重合關系。球體的直徑為1m，立方體邊長1m。兩物體的材料密度均為0.2g/cm3，球體質量為104.72kg，立方體質量為200kg。圖中的撞擊速度v=7.6km/s，相對撞擊速度為15.2km/s。建模中有限元單元邊長定為2cm，因此在仿真分析中可以得到邊長大于2cm的置信碎片結果。

圖2 球體與立方體碰撞解體仿真模型

圖3給出了350μs時刻碎片云的對比。左圖給出了SPH仿真結果，從圖中可以看出，SPH的結果表現(xiàn)為粒子云。由局部放大圖中難以分辨出哪些粒子構成一個碎片，破碎材料和未破碎材料的界面難以清晰定義。中圖為有限元重構的結果，結果由SPH粒子和有限元單元混合表達。在局部放大圖中，由單元表達的單個碎片清晰可見。右圖為有限元重構結果中的單元部分，也就是置信碎片部分。SPH粒子代表了充分破碎的材料，作為整體表達了碎片云的形狀、演化與材料分布，但是作為單個碎片是不夠可信的。

圖3 球體與立方體碰撞解體數(shù)值模擬結果

固定單元是否失效的狀態(tài)不變，將節(jié)點和粒子坐標反推回0時刻，就得到了反推圖。圖4所示為350μs的反推圖。其中由單元表示的部分在撞擊過程結束后會轉變?yōu)橹眯潘槠蛘哒f是較大的碎片；由粒子表示的部分在撞擊結束后會轉變?yōu)榉侵眯潘槠?，這些部分會充分破碎為小碎片云。由圖中可見，撞擊區(qū)域內的材料大部分都充分破碎了，大碎片主要由遠離撞擊區(qū)的材料形成。

圖4 球體與立方體碰撞結果反推圖

有限元重構方法識別了SPH仿真結果中的碎片，確定了清晰的破碎與未破碎材料的界面，對SPH仿真結果中單個碎片的置信部分和非置信部分進行了清晰的劃分。首先，在解體過程中如果單元失效，它將被轉換為粒子；而粒子繼續(xù)受變形破裂，就無法描述，因為其不可再分。因此在解體數(shù)值模擬結果中，以單元形式存在的物質就構成了置信碎片；而以粒子形式存在的失效單元則屬于非置信碎片。對于單個碎片的分析應該基于置信碎片或者說基于單元部分進行。第二，由于引入了單元，刻畫了單元之間的物質連通性，初始相鄰的單元如果解體結束后仍然相鄰則屬于同一個碎片；如果初始不相鄰，而在解體結束后相鄰，則說明它們的鄰近是由于材料的運動引起的，因此不構成同一個碎片。這樣，哪些單元構成一個碎片將由單元部分的拓撲結構決定。第三，單元轉換粒子判據是材料失效或相變，而不是人為指定幾何閾值，因此置信與否的分界線從根本上是由力學機理決定的。

3 碎片統(tǒng)計方法

在有限元重構的仿真結果中，置信碎片由單元部分表達，碎片的總數(shù)以及碎片的單元構成是由單元部分的拓撲結構決定的。碎片統(tǒng)計就是要設計一種自動的算法，能夠統(tǒng)計出碎片的數(shù)目、碎片的單元構成以及每個碎片的尺寸、質量和速度矢量等信息。

把置信部分的單元看作頂點，如果兩個單元之間存在一個公共節(jié)點，則在兩單元之間定義一條邊，這樣頂點和邊的集合可以表示為：

(9)

〈Vtx,Edg〉構成一個無向圖。如果該圖中存在一個首尾相接的邊序列，能夠把圖中兩個頂點連接起來，則這兩個頂點是連通的。如果〈Vtx,Edg〉的一個子集中任意兩個頂點是連通的，則稱這個子集為該無向圖的一個連通分支。每一個連通分支就對應了碎片云中的一個碎片。因此，碎片的統(tǒng)計與特征量計算就是對連通分支的分析。

本文使用了基于廣度優(yōu)先的連通分支遍歷算法，時間復雜度為o(n2)。對圖4中的碎片進行分析，圖5給出了連通分支搜尋的結果，其中每一種顏色代表了一個連通分支，也就是一個碎片。該算例中共產生了大于2cm的碎片562個，其中178來自球體，384個來自立方體。

在獲得了碎片的單元構成之后，可以對每個碎片的特征長度進行計算：

(10)

圖5 碎片統(tǒng)計結果

4 衛(wèi)星模型解體碎片數(shù)值模擬

利用上述方法對文獻[3]開展的地面衛(wèi)星模型撞擊解體實驗進行數(shù)值模擬。

圖6給出了衛(wèi)星解體仿真模型，圖7為Case B的仿真結果。從圖7中可以看出，雖然彈丸尺寸較小，但是造成了很大的破壞。由于衛(wèi)星模型在撞擊法平面內的兩個軸向并不對稱，所以結構中的開孔也不對稱，而是呈橢圓形。仿真結果中，衛(wèi)星結構并沒有像實驗中一樣板結構之間互相脫離。這是因為在實驗結構中，相鄰板結構采用了L型鈑金配合螺釘連接，連接強度較弱；在仿真中，連接形式進行了近似處理，直接采用了相鄰邊焊接的方式。

圖6 衛(wèi)星模型撞擊解體仿真建模

圖8給出了Case A和Case B的累計碎片數(shù)目對碎片尺寸的分布曲線對比?？傮w來說仿真結果與實驗數(shù)據比較一致。仿真曲線橫坐標值在0.001以下的缺失是由模型分辨率引起的，為了獲得更小碎片的特性，需要進一步提高模型分辨率。另外，在小尺寸碎片區(qū)域，NASA模型的計算值明顯高于實驗值，這可能是由于NASA模型的建模數(shù)據源與文獻[3]實驗的差異造成的。NASA模型小碎片特性部分的建模數(shù)據主要來源于SOCITE4解體實驗，該實驗中用到的OSCAR衛(wèi)星的內部結構與材料和文獻[3]實驗中用的衛(wèi)星模型相差較大。這也體現(xiàn)了衛(wèi)星的結構和材料對解體碎片特性的影響，表明了基于解體動力學過程分析、計算碎片特性的必要性。

圖7 Case B數(shù)值模擬結果

(a) Case A

(b) Case B

5 結 論

本文基于超高速碰撞數(shù)值模擬，建立了一種用于模擬衛(wèi)星解體碎片生成的有限元重構方法，對衛(wèi)星模型解體問題進行了仿真研究，得到如下結論：

(1) 有限元重構方法能夠對衛(wèi)星解體生成的碎片進行模擬，結合無向圖連通域分析方法，能夠獲得碎片分布以及單個碎片的數(shù)據；

(2) 對地面衛(wèi)星模型解體實驗進行了數(shù)值模擬分析。對于Case A和Case B兩種工況，模擬獲得的累計碎片數(shù)對碎片尺寸的分布曲線與實驗曲線基本一致。表明了本文方法的有效性。另外，通過和NASA模型計算結果對比，表明解體衛(wèi)星的結構形式和材料對解體碎片特性有較大影響，通過對解體動力學過程的研究，可以更加有效的計算解體碎片特性。

參考文獻:

[1]Krisko P H, Horstman M, Fudge M L. SOCIT4 collisional-breakup test data analysis: With shape and materials characterization[J]. Advances in Space Research, 2008, 41(7): 1138-1146.

[2]Hanada T, Liou J C, Nakajima T, et al. Outcome of recent satellite impact experiments[J]. Advances in Space Research, 2009, 44(5): 558-567.

[3]Lan S, Li Y, Huang J, et al. Fragment characteristic of simulated spacecraft under hypervelociy impact[C]//Proceeding of 62nd International Astronautical Congress. Cape Town, SA, 2011.

[4]蘭勝威, 李毅, 黃潔, 等. 模擬航天器超高速撞擊解體碎片特性研究[C]//第六屆全國空間碎片學術交流會論文集, 四川成都, 2011.

Lan Shengwei, Li Yi, Huang Jie, et al. Investigation on the characteristics of fragments generated by the simulated spacecraft hypervelocity impact test[C]//Proceeding of 6thNational conference of space debris. Chengdu, Sichuan, 2011.

[5]李毅, 黃潔, 馬兆俠, 等. 一種新的衛(wèi)星超高速撞擊解體閾值模型研究[C]//第六屆全國空間碎片學術交流會論文集. 四川成都, 2011.

Li Yi, Huang Jie, Ma Zhaoxia, et al. Investigation on new

breakup threshold model of satellite breakup under hypervelocity impact[C]//Proceeding of 6thNational conference of space debris. Chengdu, Sichuan, 2011.

[6]余慶波, 金學科, 王勤智, 等. 結構強度對衛(wèi)星爆炸解體影響研究[C]//第六屆全國空間碎片學術交流會論文集. 四川成都, 2011.

Yu Qingbo, Jin Xueke, Wang Qinzhi, et al. The influence of structure strength on satellite breakup caused by explosion[C]//Proceeding of 6thNational conference of space debris. Chengdu, Sichuan, 2011.

[7]金學科, 余慶波, 王勤智, 等. 燃料量對衛(wèi)星爆炸解體影響研究[C]//第六屆全國空間碎片學術交流會論文集. 四川成都, 2011.

Jin Xueke, Yu Qingbo, Wang Qinzhi, et al. The influence of the amount of fuel on satellite breakup caused by explosion[C]//Proceeding of 6thNational conference of space debris. Chengdu, Sichuan, 2011.

[8]Johnson G R, Stryk R A, Beissel S R. SPH for high velocity impact computations[J]. Computer Methods in Applied Mechanics and Engineering, 1996, 57(4): 573-580.

作者簡介:

張曉天(1984-)，男，山東濟南人，博士，講師。研究方向：高速碰撞動力學，航天器防護結構設計，空間碎片與空間環(huán)境。通信地址：北京航空航天大學宇航學院航天飛行器技術系(100191)。E-mail: zhangxiaotian@buaa.edu.cn