超高速撞擊中靶后碎片云的內(nèi)外邊界方程研究

賈光輝，劉巧霞，張 榜

（北京航空航天大學(xué) 宇航學(xué)院，北京 100191）

0 引言

空間碎片超高速撞擊嚴(yán)重威脅著低地球軌道航天器的安全。為保障在軌航天器的安全運(yùn)行，Whipple在1947年提出在航天器艙壁外面一定距離處放置一塊防護(hù)屏，即Whipple防護(hù)，后來又發(fā)展了填充式防護(hù)或多層防護(hù)[1-2]。當(dāng)空間碎片擊穿防護(hù)屏，會(huì)在防護(hù)屏后方形成碎片云，使得空間碎片對(duì)航天器艙壁的威脅變成碎片云對(duì)艙壁的威脅[3]。因此，研究超高速撞擊防護(hù)靶形成的碎片云特性對(duì)航天器防護(hù)結(jié)構(gòu)設(shè)計(jì)具有重要意義。不同撞擊條件下形成的碎片云具有不同的特性[4]，根據(jù)撞擊角度可將碎片云分為斜撞擊碎片云和正撞擊碎片云，本文針對(duì)超高速正撞擊形成的碎片云特性進(jìn)行研究。

在已有的碎片云形狀模型中，Swift等[5]首先提出靶后碎片云均勻分布在球殼上的模型，該模型雖結(jié)果不甚精確，但建模思想具有開創(chuàng)性；Corvonato等[6]采用伯努利雙紐線建立了靶后碎片云的邊界方程，使得碎片云形狀更貼近實(shí)際；遲潤(rùn)強(qiáng)[7]在Corvonato等研究的基礎(chǔ)上，保留靶后碎片云的雙紐線邊界方程，并且利用雙曲線方程建立了反濺碎片云的邊界方程。然而，以上邊界方程模型只建立了碎片云的外邊界方程，沒有區(qū)分靶板和彈丸的內(nèi)外邊界。

本文將碎片云根據(jù)材料來源分為靶板碎片云和彈丸碎片云，基于數(shù)值仿真結(jié)果分別建立這兩種碎片云的內(nèi)外邊界方程，進(jìn)而得到整體碎片云的形狀模型，旨在為航天器防護(hù)結(jié)構(gòu)設(shè)計(jì)提供參考。

1 數(shù)值仿真的有效性驗(yàn)證

文獻(xiàn)[8]給出一組碎片超高速撞擊實(shí)驗(yàn)案例：直徑 0.953cm 鋁球（Al 2017.T4）以 6.64 km/s速度垂直撞擊厚度 0.22cm 鋁板（Al 6061.T6）。案例的實(shí)驗(yàn)照片如圖1所示，實(shí)驗(yàn)所用球/板的材料參數(shù)見表1。

圖1 撞擊發(fā)生后 6 μs 及 19.8 μs 的實(shí)驗(yàn)碎片云Fig.1 Debris clouds obtained in experiment at 6 μs and 19.8 μs

表1 超高速撞擊實(shí)驗(yàn)所用材料參數(shù)Table 1 Material parameters for the hypervelocity impact experiment

采用SPH算法對(duì)上述算例的工況建模并進(jìn)行仿真計(jì)算，建模時(shí)SPH立方體粒子的邊長(zhǎng)取0.004 4cm，粒子數(shù)為 6 773 984，材料模型的主要參數(shù)見表1。仿真結(jié)果與實(shí)驗(yàn)結(jié)果的疊加對(duì)比如圖2所示，可見數(shù)值仿真結(jié)果與對(duì)應(yīng)實(shí)驗(yàn)結(jié)果的一致性很好，說明采用SPH算法是可行的。

圖2 撞擊發(fā)生后 6 μs及 19.8 μs的 SPH 仿真碎片云與實(shí)驗(yàn)碎片云對(duì)比Fig.2 Comparison of debris clouds obtained in simulation and experiment at 6 μs and 19.8 μs

2 碎片云形狀的雙紐線邊界方程

2.1 雙紐線邊界方程

2001年，Corvonato等基于碎片云滿足自相似性的假設(shè)，采用伯努利雙紐線作為碎片云形狀的邊界函數(shù)[6]，建立了碎片云形狀模型。所謂碎片云的自相似性即在演化過程中碎片云的形狀始終保持不變。其理論基礎(chǔ)為：彈丸與靶板發(fā)生撞擊后，撞擊體的材料迅速破碎，碎片之間沒有相互作用力可以用來獲取相對(duì)加速度，故碎片做勻速直線運(yùn)動(dòng)。基于該假設(shè)，可將各時(shí)刻碎片云粒子的物理坐標(biāo)(x,y)分別沿著橫坐標(biāo)/縱坐標(biāo)進(jìn)行碎片云粒子物理坐標(biāo)的歸一化，得到自然坐標(biāo)(ξ,η)下的碎片云分布，其變換式為

式中：xmax和ymax分別表示對(duì)應(yīng)時(shí)刻碎片云粒子在正x向和正y向的最大值。

在自然坐標(biāo)下，為使η=φ(ξ)滿足伯努利雙紐線方程，需要對(duì)自然坐標(biāo)(ξ,η)進(jìn)一步變換得到新自然坐標(biāo) (ξ′,η′)，變換式為

Corvonato提出碎片云形狀的邊界方程滿足雙紐線函數(shù)[6]，并將碎片云根據(jù)材料來源分為靶板碎片云和彈丸碎片云，其中靶板碎片云的邊界為外邊界，其邊界方程為

彈丸碎片云的邊界為內(nèi)邊界，其邊界方程為

式中K為內(nèi)外邊界方程的比例系數(shù)，可由靶板碎片云最前端速度與彈丸碎片云最前端速度的比值得到，且K關(guān)于厚徑比和撞擊速度的函數(shù)關(guān)系[6]為

式中：t為鋁板厚度，cm；dp為鋁球直徑，cm；v為撞擊速度，km/s；c為材料中的聲速，km/s。

若取K=0.8，則碎片云邊界如圖3所示，圖中外邊界為靶板碎片云的外邊界，內(nèi)邊界既是靶板碎片云的內(nèi)邊界又是彈丸碎片云的外邊界。

圖3 由文獻(xiàn)[6]碎片云邊界方程得到的碎片云邊界Fig.3 Boundary curves of debris cloud for Ref.[6]

2.2 雙紐線邊界方程的適用性分析

本文針對(duì)直徑 0.6cm的鋁球以5.5～8 km/s的速度垂直撞擊厚度0.16cm鋁板的工況研究靶后碎片云的內(nèi)外邊界方程，其中算例的工況設(shè)置與文獻(xiàn)[6]相同，以便對(duì)比。

撞擊速度為6 km/s的工況下，根據(jù)式(5)確定內(nèi)外邊界方程的比例系數(shù)K為0.93，根據(jù)式(3)、式(4)得到算例靶后碎片云的內(nèi)外邊界方程，整體碎片云、靶板碎片云和彈丸碎片云在新自然坐標(biāo)系的分布分別如圖4、圖5和圖6所示。

圖4 整體碎片云分布Fig.4 Distribution of the entire debris cloud

圖5 靶板碎片云分布Fig.5 Distribution of target debris cloud

圖6 彈丸碎片云分布Fig.6 Distribution of projectile debris cloud

為了衡量邊界方程的適用性，提出邊界方程適用率概念。在新自然坐標(biāo)系下，若碎片云粒子p分布在所給的邊界內(nèi)，說明粒子p滿足該邊界方程，反之則不滿足該邊界方程。將滿足邊界方程的粒子數(shù)與碎片云總粒子數(shù)的比值稱為邊界方程適用率，簡(jiǎn)稱適用率。統(tǒng)計(jì)圖5中分布在內(nèi)外邊界間的靶板粒子數(shù)n1，n1與靶板總粒子數(shù)n2的比值即為靶板適用率；統(tǒng)計(jì)圖6中分布在內(nèi)邊界里的彈丸粒子數(shù)n3，n3與彈丸總粒子數(shù)n4的比值即為彈丸適用率；(n1+n3)與(n2+n4)的比值即為碎片云整體適用率。根據(jù)文獻(xiàn)[6]所給雙紐線邊界方程，得到算例對(duì)應(yīng)的雙紐線邊界方程的靶板、彈丸和整體適用率分別為0.615 2、0.146 4 和0.353，可見其整體適用率較低，說明雙紐線邊界方程的描述不完善，應(yīng)根據(jù)碎片云的分布特點(diǎn)修正碎片云的內(nèi)外邊界方程。

3 碎片云的內(nèi)外邊界方程修正

3.1 內(nèi)外邊界方程

由圖5和圖6可看出：對(duì)于靶板碎片云部分，頭部出現(xiàn)縮頸；對(duì)于彈丸碎片云部分，頭部明顯后凸。因此本文認(rèn)為碎片云邊界方程應(yīng)用分段函數(shù)表示，設(shè)分界線為ξ′=m。根據(jù)圖1，假設(shè)碎片云邊界的分段函數(shù)由雙紐線和圓弧曲線組成。在新自然坐標(biāo)系下，靶板碎片云粒子基本滿足等厚度（設(shè)厚度符號(hào)為tr，量綱為1）的殼體分布，根據(jù)分界點(diǎn)將靶板碎片云粒子分為靶板前段和靶板后段2部分，如圖7所示。同樣，根據(jù)分界點(diǎn)將彈丸碎片云粒子分為彈丸前段（包括彈丸碎片云頭部后凸部分）和彈丸后段2部分，如圖8所示。在新自然坐標(biāo)系下討論：1）靶板后段內(nèi)外邊界方程；2）靶板前段內(nèi)外邊界方程；3）彈丸后段內(nèi)外邊界方程；4）彈丸前段內(nèi)外邊界方程。

圖7 靶板碎片云分布特征Fig.7 Distribution characteristics of target debris cloud

圖8 彈丸碎片云分布特征Fig.8 Distribution characteristics of projectile debris cloud

1）靶板后段內(nèi)外邊界方程

根據(jù)圖7可知靶板后段外邊界基本滿足雙紐線方程

根據(jù)靶板碎片云等厚度tr殼體分布，可確定靶板后段內(nèi)邊界雙紐線方程。具體方法如下：

首先，靶板后段內(nèi)邊界曲線滿足雙紐線

式中：(ξ0′,η0′)為雙紐線的頂點(diǎn)；a為雙紐線系數(shù)，a＞0。

由靶后碎片云關(guān)于ξ′軸對(duì)稱，得η0′=0；靶板碎片云等厚度tr殼體分布，因此由式(6)可知靶板內(nèi)邊界函數(shù)在頂點(diǎn)處的法線方向與ξ′軸的夾角為45°，故得靶板碎片云內(nèi)邊界函數(shù)過頂點(diǎn)和端點(diǎn)(1?tr,0)，代入式(7)可得到靶板后段內(nèi)邊界的雙紐線方程，其中。

2）靶板前段內(nèi)外邊界方程

利用三點(diǎn)定圓原理確定靶板前段內(nèi)外邊界方程。如圖7 所示：由A(1,0)、B(m,φ(m)?tr)、B′(m,-[φ(m)?tr])三點(diǎn)可確定靶板前段外邊界圓弧方程（點(diǎn)B′是點(diǎn)B的ξ′軸對(duì)稱點(diǎn)）；由于靶板前段碎片云等厚度tr殼體分布，其內(nèi)外邊界是同心圓，所以可由靶板前段外邊界圓弧方程和殼體厚度tr直接確定靶板前段內(nèi)邊界圓弧方程。

3）彈丸后段內(nèi)外邊界方程

根據(jù)圖7和圖8的對(duì)比可知，彈丸后段和靶板后段分布區(qū)域基本重合，故彈丸后段與靶板后段有相同的內(nèi)外邊界方程。

4）彈丸前段內(nèi)外邊界方程

為了得到彈丸前段內(nèi)外邊界條件，取彈丸前段在ξ′方向的特征點(diǎn)D(n,0)和C(k,0)，利用三點(diǎn)定圓原理確定彈丸前段內(nèi)外邊界方程。如圖8所示，由C(k,0)、B(m,φ(m)?tr)、B′(m,-[φ(m)?tr])三點(diǎn)可確定彈丸前段外邊界圓弧方程；由D(n,0)、E(m?tr,φ(m)?tr)、E′(m?tr,-[φ(m)?tr])三點(diǎn)可確定彈丸前段內(nèi)邊界圓弧方程（點(diǎn)E′是點(diǎn)E的ξ′軸對(duì)稱點(diǎn)）。

雙紐線與圓弧在分界ξ′=m處不是連續(xù)函數(shù)，因此在分界處將靶板后段外邊界雙紐線對(duì)應(yīng)邊界點(diǎn)與靶板前段外邊界圓弧對(duì)應(yīng)邊界點(diǎn)用直線連接，并用同樣的方法處理靶板前后內(nèi)邊界、彈丸前后外邊界、彈丸前后內(nèi)邊界。圖5、圖6對(duì)應(yīng)算例的特征參數(shù)為：tr=0.04，m=0.886，n=0.76，k=0.98。

為研究撞擊速度對(duì)邊界方程的影響，本文開展了直徑為0.6cm的鋁球以5.5～8.0 km/s的速度垂直撞擊厚度0.16cm鋁板的數(shù)值仿真，工況見表2。根據(jù)仿真結(jié)果，用同樣的處理方法，統(tǒng)計(jì)各算例的4 個(gè)特征參數(shù)（tr、m、n、k），見表2。

表2 各算例的 4 個(gè)特征參數(shù)Table 2 Determination of parameters

進(jìn)一步分析撞擊速度對(duì)4個(gè)特征參數(shù)的影響：參數(shù)tr為常數(shù)，tr=0.04；參數(shù)k隨速度的變化較小，且接近于1；參數(shù)m隨撞擊速度增加而減小，即分界點(diǎn)后移，

參數(shù)n相對(duì)于其他3個(gè)參數(shù)變化幅度最大，其值隨撞擊速度增加而減小，即彈丸碎片云后端點(diǎn)后移，

根據(jù)m、n值的變化可知，在撞擊速度v=5.5～8.0 km/s的范圍內(nèi)，撞擊速度越大，撞擊材料的破碎程度越高，使得m、n值越小。

3.2 修正后內(nèi)外邊界方程的具體形式

以表2中算例2為例， 4個(gè)特征參數(shù)為：tr=0.04、m=0.886、n=0.76、k=0.98，其碎片云的邊界方程見式(6)和式(10)～(13)。

1）靶板后段內(nèi)外邊界方程

靶板后段外邊界方程為式(6)。

靶板后段內(nèi)邊界方程為

2）靶板前段內(nèi)外邊界方程

靶板前段外邊界方程為

靶板前段內(nèi)邊界方程為

3）彈丸后段內(nèi)外邊界方程

彈丸后段外邊界方程與靶板后段外邊界方程一致，為式(6)。彈丸后段內(nèi)邊界方程與靶板后段內(nèi)邊界方程一致，為式(10)。

4）彈丸前段內(nèi)外邊界方程

彈丸前段外邊界方程與靶板前段外邊界方程一致，為式(11)。

彈丸前段內(nèi)邊界方程為

根據(jù)式(6)和式(10)～(13)，靶板碎片云邊界和彈丸碎片云邊界分別如圖9和圖10所示。

在新自然坐標(biāo)下，根據(jù)雙紐線方程和圓弧方程組成的邊界方程，可得算例2的適用率，與雙紐線邊界方程的適用率進(jìn)行對(duì)比（見表3）可知：靶板適用率由61.52%提高到95.07%；彈丸適用率由14.64%提高到97.22%，碎片云整體適用率由35.3%提高到96.27%。表4是各算例的兩種邊界方程碎片云整體適用率的對(duì)比結(jié)果。

圖9 靶板碎片云內(nèi)外邊界Fig.9 The inner and outer boundaries of target debris cloud

圖10 彈丸碎片云內(nèi)外邊界Fig.10 The inner and outer boundaries of projectile debris cloud

表3 算例 2 兩種邊界方程適用率對(duì)比Table 3 Comparison of applicability rate for Example II

表4 各算例的兩種邊界方程碎片云整體適用率對(duì)比Table 4 Comparison of applicability rate for various Examples

最小撞擊速度（算例 1：v=5.5 km/s）和最大撞擊速度（算例 6：v=8.0 km/s）工況下碎片云的內(nèi)外邊界分布分別如圖11和圖12所示，可以看出由雙紐線方程和圓弧方程組成的邊界方程比單一的雙紐線邊界方程與碎片云實(shí)際分布情況更貼合。

圖11 算例 1（v=5.5 km/s、t/dp=0.267）碎片云的內(nèi)外邊界分布Fig.11 The inner and outer boundaries of debris cloud in Example I (v=5.5 km/s,t/dp =0.267)

圖12 算例 6（v=8.0 km/s、t/dp=0.267）碎片云的內(nèi)外邊界分布Fig.12 The inner and outer boundaries of debris cloud in Example VI (v=8 km/s,t/dp=0.267)

4 結(jié)束語

本文采用SPH算法，選取文獻(xiàn)[6]中的算例，對(duì)直徑0.6cm的鋁球以5.5～8.0 km/s的速度正撞擊厚度0.16cm的鋁板進(jìn)行了仿真計(jì)算，研究了碎片云的內(nèi)外邊界方程。

通過分析靶板碎片云、彈丸碎片云以及整體碎片云的適用率，發(fā)現(xiàn)單一的雙紐線邊界方程與碎片云實(shí)際分布情況的貼合度不夠；基于靶后碎片云的幾何分布特點(diǎn)，分別建立靶板碎片云和彈丸碎片云內(nèi)外邊界方程的分段函數(shù)（后段仍然采用雙紐線方程，前段采用圓弧方程），完整表述了靶后碎片云彈丸內(nèi)凸、靶板內(nèi)凹的形狀特點(diǎn)；將修正后的由雙紐線方程和圓弧方程組成的邊界方程與原來單一的雙紐線邊界方程的整體碎片云適用率進(jìn)行對(duì)比，發(fā)現(xiàn)前者與碎片云實(shí)際分布情況的貼合度更好。

本文僅探討了直徑為0.6cm的鋁球正撞擊厚度0.16cm鋁板的工況，修正后邊界方程在彈靶尺寸方面的普適性有待進(jìn)一步研究。

[1]賈光輝,姚光樂,張帥.填充式防護(hù)結(jié)構(gòu)彈道極限方程的差異演化優(yōu)化[J/OL].北京航空航天大學(xué)學(xué)報(bào)(2017-11-28)[2018-01-26].http://onlinelibrary.wiley.com/doi/10.13700/j.bh.1001-5965.2017.0515 JIA G H,YAO G L,ZHANG S.Differential evolution optimization for Stuffed Whipple shield ballistic limit equations[J/OL].Journal of Beijing University of Aeronautics and Astronautics (2017-11-28)[2018-01-26].http://onlinelibrary.wiley.com/doi/10.13700/j.bh.1001-5965.2017.0515

[2]CHRISTIANSEN E L.Design and performance equations for advanced meteoroid and debris shields[J].International Journal of Impact Engineering,1993,14(1):145-156

[3]賈光輝,黃海,胡震東,等.空間碎片超高速撞擊極限穿透比動(dòng)能研究[J].宇航學(xué)報(bào),2006,27(5):1000-1003 JIA G H,HUANG H,HU Z D,et al.Research on debris kinetic energy per unit area in hypervelocity impacting[J].Journal of Astronautics,2006,27(5):1000-1003

[4]唐蜜,柏勁松,李平,等.超高速撞擊中影響碎片云形狀因素分析[J].高壓物理學(xué)報(bào),2007,21(4):425-432 TANG M,BAI J S,LI P,et al.Factors analysis of debris cloud ’s shape of hypervelocity impact[J].Chinese Journal of High Pressure Physics,2007,21(4):425-432

[5]SWIFT H F,BAMFORD R,CHEN R.Designing space vehicle shields for meteoroid protection:a new analysis[J].Advances in Space Research,1982,2(12):219-234

[6]CORVONATO E,DESTEFANIS R,FARAUD M.Integral model for the description of the debris cloud structure and impact[J].International Journal of Impact Engineering,2001,26(1):115-128

[7]遲潤(rùn)強(qiáng).彈丸超高速撞擊薄板碎片云建模研究[D].哈爾濱:哈爾濱工業(yè)大學(xué),2010:104-106

[8]LACERDA D,LACOME J-L.Simulations of hypervelocity impacts with smoothed particle hydrodynamics:DYNALIS Report[R].Paris,France,2003