基于最優運輸無網格法的Whipple 屏超高速撞擊數值模擬*

樊 江，袁 圓，廖祜明，袁慶浩，陳高翔，黎 波

（1. 北京航空航天大學能源與動力工程學院，北京 100191；2. 凱斯西儲大學機械與航空航天工程學院，美國 克利夫蘭，OH44106）

太空中微流星體和空間碎片雖然體形很小，但由于通常都具有超高的相對速度（2～15 km/s），會對航天器的安全造成巨大的威脅。例如，長時間暴露的航天器外部超過30 000 個直徑大于0.3 mm 的隕石坑都是由微流星體或軌道碎片撞擊形成[1]。航天防護中常采用Whipple 屏對特征尺寸1 cm 以下的碎片進行破碎[2-3]。Whipple 屏在超高速撞擊后的損傷狀態對于Whipple 屏的優化設計至關重要。

超高速撞擊的物理過程非常復雜。撞擊過程中，材料的慣性、可壓縮性效應或相變效應比結構效應更顯著，伴隨著大變形、熱流固耦合、相變（液化、氣化、凝固等）、碎裂等現象，采用數值方法進行準確模擬是一項巨大的挑戰[4]。基于網格的超高速撞擊模擬存在著網格畸變或者需要不斷重新劃分網格[5]的缺點；而無網格法由于不需要進行網格離散及采用高階插值函數，有利于解決大變形和流固耦合問題，在超高速撞擊問題上應用更為廣泛。目前常用的無網格方法有光滑粒子流體動力學法（smoothed particle hydrodynamics, SPH）和再生核粒子法（reproducing kernel particle method, RKPM）以及質量點法（material point method, MPM）等。有大量的研究應用以上的方法[6-7]及其改進方法[8]對Whipple 屏的超高速撞擊進行模擬，雖然取得了豐富研究成果，但由于算法本身的固有缺陷，仍然無法解決準確設置位移邊界條件[9]、計算效率不高[10]、拉應力不穩定[11]及有效計算帶摩擦的動態接觸等問題，同時也缺乏嚴密的收斂性與誤差理論分析。

最優運輸無網格方法（optimal transportation meshfree, OTM）是由Li 等[12]提出的一種基于最優運輸理論和局部最大熵插值函數的無網格方法，該方法結合基于能量釋放率的物質點失效方式，能夠很好地解決現有Whipple 屏超高速撞擊模擬方法存在的問題。

本文首先采用OTM 法模擬鋁球超高速撞擊單層鋁板，通過與實驗結果以及其他數值方法的計算結果的比較，驗證OTM 法在超高速撞擊問題上的適用性；然后采用OTM 法對Whipple 屏超高速撞擊問題進行模擬，研究不同速度不同撞擊角度下的碎片云形狀、緩沖墻彈孔尺寸以及后墻剝落穿透的損傷情況，并與文獻[13]的實驗結果進行對比。

1 OTM 理論

OTM 法的主要特點是采用局部最大熵無網格插值函數，克服了一般無網格法中插值函數不滿足Kronecker delta 屬性的本質缺陷，解決了傳統無網格法難以準確施加位移邊界條件的問題；采用物質點充當積分點，有效避免了計算結果在拉伸載荷下的不穩定性；采用最優運輸理論對時間離散，保證了哈密頓作用量的時間離散形式滿足動量守恒條件，且收斂性能得到嚴格的數學證明[14]；采用基于能量釋放率的物質點失效方式，能夠很好地模擬材料的損傷情況，并且該方式已被證明可收斂到 Griffith 斷裂準則[15]。

1.1 時間離散

1.2 空間離散

OTM 法將計算域分為兩類點的集合。一類為節點，用下標a 表示，包含了位置、速度、加速度的信息；一類為物質點，用下標p 表示，包含了材料參數、質量、應力應變等信息，相當于一般有限元法中的積分點。定義在時刻節點的坐標為，物質點的坐標為如圖1 所示。

圖1 空間離散示意圖[17]Fig. 1 Spatial discrete diagram[17]

1.3 失效準則

圖2 評估能量釋放率的局部鄰域[15]Fig. 2 The local neighborhood used to estimate the energy-release rate[15]

2 鋁球超高速撞擊鋁板算例分析

文獻[8]采用基于擬流體模型的SPH 新方法（擬流體SPH 法）對鋁球超高速撞擊鋁板進行了數值模擬，并與Hiermaier 等的實驗結果、傳統SPH 法和自適應光滑粒子流體動力學法（adaptive smoothed particle hydrodynamics, ASPH）的計算結果進行了對比。本節將采用OTM 法對相同的算例進行模擬，并與文獻[8]中的結果進行比較。

圖3 鋁球撞擊鋁板離散模型Fig. 3 Discrete model of aluminum ball impacting single aluminum plate

2.1 數值模擬模型

2.2 材料模型

超高速撞擊需要考慮彈丸和靶材的應變率硬化和熱軟化等問題，因此選擇能較好地模擬鋁合金材料塑性響應的J2 粘塑性模型（J2-viscoplasticity）。

J2 粘塑性模型有效屈服應力為

式中：J=V/V0，V 為當前體積，V0為初始體積；為參考熔化溫度；為參考體積下的Grüneisen 參數，為常數。鋁合金型號為LY12，材料參數如表1 所示。本算例中依據文獻[8]提供的材料信息進行了部分參數修改。式(15)和式(16)中提到的各參考系數的取值如表2 所示。高溫高壓下材料的變形與溫度、壓力的關系采用SESAME 狀態數據庫描述[18]。

表1 LY12 材料參數Table 1 Material parameters of LY12

表2 J2 黏塑性模型參數Table 2 Parameters of J2 viscoplasticity model

2.3 結果分析

OTM 法能夠很好地模擬內核碎片云的位置、外泡碎片云的形態以及反濺碎片云的形態等特征信息。特別是反濺碎片云的膨脹距離和寬度，與實驗很好吻合。一般的SPH 法由于采用Johnson-Cook 損傷模型而導致薄板屈服應力小于真實的屈服應力，外濺碎片云的反濺程度過大，與實驗偏差較大（如圖4 所示）。

圖4 OTM 法與各類SPH 方法計算結果對比Fig. 4 Comparison of OTM and various SPH methods’ simulation results

表3 鋁球超高速撞擊鋁板結果對比Table 3 Comparison of high-velocity impact results between aluminum projectile and plate

3 Whipple 屏超高速撞擊模擬

3.1 對比實驗簡介

文獻[13]進行了一系列Whipple 屏超高速撞擊實驗（如圖5 所示），球形彈丸直徑0.4～0.5 cm，撞擊速度4.47～6.15 km/s，撞擊角度分為0°和45°兩種。靶材間距為10 cm，厚度為0.192 cm。實驗得到不同撞擊速度和撞擊角度下的彈孔尺寸、后墻損傷情況和碎片云激光陰影照片等結果。

3.2 數值模擬模型

本文按照對比實驗建立模型，模擬了撞擊角度為0°和45°兩種情況。在劃分網格時細化了緩沖墻和后墻的中心區域，最小網格尺寸為0.2 mm。在正撞模擬中，物質點共有340 382 個，節點共有67 389 個（如圖6 所示）；在斜撞模擬中，物質點共有340 473 個，節點共有67 401 個。

彈丸和靶材的材料型號均為LY12，材料模型仍然使用J2 粘塑性模型（材料參數見表1、表2）和SESAME 狀態方程。

3.3 結果分析

文獻[13]進行了一系列不同撞擊速度的Whipple 屏正撞與斜撞實驗。參數設置如表4 所示。

OTM 法采用相同的參數進行對應的數值模擬，實驗與仿真中的緩沖墻彈孔尺寸對比結果如表5 所示。其中后四組是斜撞實驗，得到的彈孔呈現橢圓形，彈孔尺寸指的是橢圓的長軸和短軸。

圖5 實驗模型示意圖Fig. 5 Schematic diagram of experimental model

圖6 Whipple 屏超高速撞擊數值模擬模型Fig. 6 The numerical simulation model of Whipple shield hypervelocity impact

表4 實驗參數設置Table 4 Parameters in experiments

表5 緩沖墻彈孔尺寸對比Table 5 Bullethole size comparison of outer bumper

緩沖墻彈孔尺寸的模擬結果與實驗吻合得較好，如在撞擊速度5.29 km/s 的實驗中，彈孔直徑為1.15 cm，而仿真的結果為1.05 cm，相對誤差為8.69%（如圖7 所示）。

圖7 緩沖墻損傷對比圖(撞擊速度5.29 km/s)Fig. 7 Damage characteristics comparison chart of outer bumper (Impact velocity 5.29 km/s)

后墻的損傷形式一般有成坑、產生鼓包、層裂、剝落和穿孔[19]。文獻[13]只關注剝落和穿透，如圖8(a)所示。圖8(b)為OTM 仿真中的剝落和穿透。

對比仿真與實驗中后墻損傷情況（如表6 所示），可見正撞仿真中，撞擊速度越大，后墻的損傷越小。這是由于速度大的彈丸被緩沖墻破碎得更充分，形成了更小的碎片，減輕了對后墻的破壞作用。這與文獻[13]中的結論一致。

圖8 實驗與仿真中的剝落和穿透Fig. 8 Definitions of spalling and penetration in experiments and simulations

表6 后墻損傷情況對比Table 6 Damagecomparison of spacecraft wall

仿真顯示，彈丸碎片撞擊到后墻上，殘余應力約為150 MPa，而撞擊較嚴重的區域殘余應力達到300 MPa 以上，如圖9 所示。

圖9 后墻損傷圖(撞擊速度5.29 km/s)Fig. 9 Damage characteristics of spacecraft wall (impact velocity is 5.29 km/s)

目前超高速彈丸與Whipple 防護屏撞擊的數值模擬，很難精確模擬出后墻的損傷情況。大部分后墻的損傷數據都來自于實驗。文獻[8]用擬流體SPH 法對Whipple 屏撞擊的研究中，也只給出了后墻的中心損傷區域，缺乏對于后墻剝落與穿透等損傷狀態的探究。由此可看出，OTM 法基于能量釋放率的物質點失效的方式，在模擬材料的斷裂損傷方面有著明顯的優勢。

OTM 法對Whipple 屏的超高速撞擊模擬不僅能模擬出彈丸穿透緩沖墻形成碎片這一過程，還能很好地模擬出碎片云的形態，包括碎片云呈現出頭部橢圓形、尾跡扇形并有擴散趨勢、碎片大多集中在頭部和緩沖墻穿孔處等細節特征，如圖10～11 所示。由于實驗中的碎片云照片是采用激光陰影技術得到的，細微的碎片無法在照片中顯示，因而實驗中的碎片云輪廓較為清晰；而OTM 計算得到的碎片云對比圖中，所有碎片均有顯示，輪廓有一定發散性。

圖10 正撞碎片云對比圖Fig. 10 Fragment cloud comparison chart of vertical impact

圖11 斜撞碎片云對比圖Fig. 11 Fragment cloud comparison chart of oblique impact

4 結 論

本文采用OTM 法對鋁球超高速撞擊鋁板和Whipple 防護屏超高速撞擊進行了數值模擬。通過鋁球超高速撞擊鋁板這一驗證算例，可看出OTM 法能夠為超高速撞擊問題提供有效的數值模擬手段。在Whipple 屏的超高速撞擊模擬中，OTM 法能夠很好地預測Whipple 防護屏與彈丸撞擊時緩沖墻和后墻的損傷情況，尤其在模擬后墻的剝落、穿透等損傷狀態和破碎過程中碎片云的形態變化方面有著顯著的優勢，驗證了其在超高速撞擊數值模擬方面的適用性，為Whipple 防護屏在航天防護方面的相關探究提供了有效的數值模擬手段。