強洪夫，范樹佳，2，陳福振，劉 虎

基于擬流體模型的SPH新方法及其在彈丸超高速碰撞薄板中的應用＊

強洪夫1，范樹佳1，2，陳福振1，劉 虎1

（1.火箭軍工程大學動力工程系，陜西 西安710025；2.中國人民解放軍96656部隊，北京100085）

引入顆粒動力學理論（擬流體模型）建立了適用于超高速碰撞的SPH新方法。將超高速碰撞中處于損傷狀態(tài)的碎片等效為擬流體，在描述其運動過程中引入了碎片間相互作用和氣體相對碎片的作用。采用該方法對球形彈丸超高速碰撞薄板形成碎片云的過程進行了數(shù)值模擬，得到了彈坑直徑、外泡碎片云和內(nèi)核碎片云的形狀、分布，并與使用傳統(tǒng)SPH方法、自適應光滑粒子流體動力學（ASPH）方法的模擬結果進行對比，結果顯示：新方法在內(nèi)核碎片云形狀和分布上計算結果更加準確。同時對Whipple屏超高速碰撞問題進行了研究，分析了不同撞擊速度下防護屏彈坑尺寸及艙壁損傷特性等特性，計算結果與實驗吻合較好且符合Whipple防護結構的典型撞擊極限曲線。

SPH；碎片云；超高速碰撞；Whipple屏

隨著航天技術的高速發(fā)展，大量空間碎片遺留在地球軌道并處于高速運動狀態(tài)，時刻威脅著在軌航天器的安全［1］。在航天器艙壁外加裝防護屏的被動防護技術作為一種有效的防護方式越來越受到重視，這種結構最先由F.L.Whipple［2］提出，防護機理為空間碎片超高速碰撞防護屏后，生成大量更為細小的碎片和顆粒，少部分碎片向相反方向運動，大部分碎片向運動前方拋出。目前，研究碎片運動的方法主要包括實驗，理論計算以及數(shù)值模擬方法。其中，A.J.Piekutowski［3］和管公順等［4］進行了大量的實驗，得出了彈丸直徑、薄板厚度以及碎片云各特征點速度的變化規(guī)律，但研究存在著成本高、周期長、得到的實驗數(shù)據(jù)有限等問題；H.F.Swift等［5］在實驗的基礎上，作出若干適當假設，結合經(jīng)驗公式對碎片云特性進行研究，提出了球形彈丸高速碰撞薄板的碎片云模型，思想具有開創(chuàng)性意義，模型雖然不斷進行細節(jié)改進，但仍存在著對碎片云描述過于簡化，與實際情況差距較大的局限性；數(shù)值模擬具有快速、直觀的特點，不僅能彌補實驗研究的不足并為實驗提供指導，還可以充分了解超高速碰撞全過程和內(nèi)在機理。目前，基于網(wǎng)格的計算方法對彈塑性材料中低速變形比較可靠，但超高速碰撞問題涉及了移動的邊界和運動的物質(zhì)交界面，當采用Lagrange有限元方法計算超高速碰撞問題時，需要引入網(wǎng)格重分和界面重構等算法進行處理，將引入難以忽略的數(shù)值誤差，且算法復雜，計算量大。雖然有限差分（FD）、有限體積（FVM）等Euler方法不存在計算大畸變問題的困難，但是難以準確描述各類界面，并且無法追蹤單碎片的運動軌跡。

光滑粒子流體動力學方法（smoothed particle hydrodynamics，SPH）［6－7］作為一種完全Lagrange無網(wǎng)格粒子法，相對于傳統(tǒng)網(wǎng)格方法具有精確追蹤自由表面、變形邊界和材料交界面等優(yōu)勢，非常適合模擬穿透、層裂、剝落及碎片等現(xiàn)象，在處理多介質(zhì)、大變形破壞問題上具有網(wǎng)格方法無法比擬的優(yōu)勢，被廣泛應用于高速沖擊，爆炸和斷裂等問題。其中，L.D.Libersky等［8］率先將材料強度效應引入SPH方法，成功地開展了高速碰撞數(shù)值模擬計算。G.R.Johnson等［9］將SPH方法和有限元方法結合開展侵徹方面的數(shù)值模擬，取得一些進步。S.Hiermaier等［10］在對鋁球超高速碰撞不同材料的薄靶實驗的基礎上，使用傳統(tǒng)SPH方法對其作了數(shù)值模擬，得到的彈坑直徑、碎片云形狀等結果與實驗結果相似。徐金宏等［11］將改進的SPH方法應用于超高速碰撞問題，研究了靶板厚度、彈丸速度和形狀對碎片云的影響。C.E.Zhou等［12］使用傳統(tǒng)SPH方法對鋁球超高速撞擊鋁薄板的實驗進行了三維數(shù)值模擬。G.R.Liu等［13］使用ASPH方法求解彈丸超高速碰撞薄板問題，相比傳統(tǒng)SPH方法增加了空間求解能力。綜上所述，SPH算法經(jīng)過不斷完善與發(fā)展，成為近年來超高速碰撞數(shù)值計算中一種非常有吸引力和有前途的新的數(shù)值方法。但目前國內(nèi)外學者在求解超高速碰撞問題時，對碰撞后的物質(zhì)大多采用Grüneisen和Tillotson狀態(tài)方程進行描述，通過將“碎片云”中顆粒的運動類比于“流體”的行為進行相態(tài)劃分，仿真結果僅能達到形態(tài)上相似、定量相近求解精度較低，無法得到碰撞后碎片運動理想效果。同時，由于碎片云在膨脹過程中碎片與碎片之間以及碎片與空氣間的相互作用起著至關重要的作用，碎片表現(xiàn)出典型的離散顆粒的性質(zhì)。而傳統(tǒng)SPH方法大多僅用于離散求解連續(xù)介質(zhì)，對于與連續(xù)介質(zhì)完全不同的離散顆粒來說傳統(tǒng)求解連續(xù)性物質(zhì)的數(shù)值方法不再有效。因此對于彈丸超高速碰撞形成碎片云的數(shù)值模擬需要采用適用于離散顆粒相求解的理論進行補充研究。

目前基于擬流體模型的用于描述氣體－顆粒兩相流動問題的雙流體模型理論逐漸成熟，陳福振等［14］已成功將顆粒動理學理論引入到傳統(tǒng)SPH方法之中，建立了適用于氣粒兩相流問題求解的光滑離散顆粒流體動力學方法（smoothed discrete particle hydrodynamics，簡稱SDPH）。該方法不僅突破了傳統(tǒng)SPH方法材料屬性固定、粒子僅作為幾何質(zhì)點的局限，將SPH粒子拓展至具有物理質(zhì)點特性的范疇，還拓展了SPH方法的應用范圍，現(xiàn)已成功應用于風沙運動過程［14］、氣－粒兩相間傳熱蒸發(fā)［15］和燃料拋撒成霧及其燃燒爆炸［16］等問題。本文中在此基礎上，將新方法引入到超高速碰撞領域，建立基于擬流體模型的SPH新方法，將超高速碰撞后產(chǎn)生的物質(zhì)離散為SPH粒子，將處于損傷狀態(tài)的顆粒相等效于擬流體，同時在描述其運動過程中引入氣體相的作用，通過施加空氣外力和計算碎片阻力的方式用來表征空氣對碎片云的影響；對球形彈丸超高速碰撞薄靶板和Whipple防護結構進行三維數(shù)值模擬，分析靶板開坑直徑、碎片云膨脹過程、形態(tài)及不同撞擊速度下Whipple結構的損傷特性。

1 數(shù)學模型

在球形彈丸超高速碰撞鋁薄板數(shù)值模擬過程中，未處于失效狀態(tài)的鋁球和鋁薄板屬于SPH方法中連續(xù)相的離散求解問題。超高速碰撞后處于失效狀態(tài)的碎片屬于離散顆粒相求解問題，因此本文在Johnson－Cook損傷模型的基礎上引入SDPH模型，采用適用于離散顆粒相求解的方法求解碎片的運動過程，將碎片等效為擬流體，引入碎片間相互作用和氣體相對碎片的作用。

1.1 SPH基本方程

在SPH方法中，采用粒子近似方法對核函數(shù)進行插值離散并對粒子i處場函數(shù)進行粒子估計有：

式中：i、j為里編號，N 為緊支域內(nèi)的粒子總數(shù)，f為任意的連續(xù)光滑函數(shù)，〈〉為近似符號，mi、ρi和ri分別表示粒子i的質(zhì)量、密度和位置矢量，Wij＝W（ri－rj，h）為核函數(shù)，h為定義核函數(shù)影響區(qū)域的光滑長度。對式（1）求導可得：

未處于失效狀態(tài)的鋁球和鋁薄板屬于SPH方法中連續(xù)相的離散求解問題，不考慮熱傳導過程，Lagrange框架下的控制方程組表述為：

式中：ρ、v、P、x分別為流體的密度、速度、壓力、位置，τ為偏應力張量，g為質(zhì)量體積力。

最后，考慮人工粘性的影響，對式（3）～（5）進行SPH離散，得到適用于超高速碰撞問題的SPH方程組：

式中：vij＝vi－vj，σ＝－!P＋!·τ，Π 為人工黏性，α和β表示坐標方向。

1.2 彈丸和靶板的本構模型及狀態(tài)方程

為描述靶板材料的屈服應力及損傷演化，這里引入Johnson－Cook損傷模型［17］，該模型中將材料的屈服強度表示為損傷變量、等效應變、等效應變率和溫度的函數(shù)：

式中：A、B、C、n、m是材料常數(shù)，D為損傷變量，D＝0表示材料沒有損傷，D＝1表示材料完全失效。r是累積損傷塑性應變，r＝（1－D）ε，rr＊＝（1－D）rε＊，rε＊＝rε／rε0，ε 是累積塑性應變，rε0是自定義參考應變率。溫度T＊＝（T－T0）／（Tm－T0），T0是室溫，Tm是材料熔點。

損傷變D是累積塑性應變ε的函數(shù)，當D＝1時發(fā)生損傷破壞：

式中：εf是斷裂塑性應變，與材料的應力三軸度、應變率和溫度相關。本構模型中的剪切損傷演化模型將εf描述如下：

式中：D1～D5為材料常數(shù)，σ＊＝σm／σeq為應力三軸度，σm＝（σx＋σy＋σz）／3為平均正應力。

超高速碰撞問題中最常用的兩個狀態(tài)方程為Tillotson狀態(tài)方程和Grüneisen狀態(tài)方程，Tillotson狀態(tài)方程形式如下：

式中：P 為壓力；a、b、AT為擬合常數(shù)；BT為調(diào)節(jié)參數(shù)；αT、βT為材料常數(shù)；E 為比內(nèi)能；E0為初始比內(nèi)能；V 為比容；μ＝ρ／ρ0－1，η＝ρ／ρ0，ρ為材料的實時密度，ρ0為材料的初始密度。

Grüneisen狀態(tài)方程形式如下：

式中：系數(shù)Γ＝1.99，CS＝3 940，SS＝1.489；ρ為材料的實時密度；ρ0為材料的初始密度；e為材料的內(nèi)能。

1.3 SDPH模型方程

超高速碰撞后處于失效狀態(tài)的碎片屬于離散顆粒相求解問題。本節(jié)中列出了SDPH模型中顆粒相質(zhì)量、動量、擬溫度守恒方程以及顆粒相的SPH離散方程等，具體方法及各參數(shù)表達式參見文獻［14］。

同時，碎片區(qū)別于具有不同粒徑、任意分布的礦石、風沙等顆粒，相對于風沙運動過程等二維模型中單位面積的顆粒數(shù)，球形彈丸超高速碰撞鋁薄板產(chǎn)生的碎片數(shù)目相對較少、計算量較小，本文SDPH單粒子不必為減少計算量而表征具有一系列粒徑分布的顆粒群。由于SDPH單粒子與單個碎片屬于幾何近似關系，因此SDPH單粒子代表一個真實碎片，攜帶質(zhì)量、密度、速度、擬溫度以及壓力等參量。但此時的SDPH粒子不等同于傳統(tǒng)SPH方法中用于求解連續(xù)相流體問題的連續(xù)性介質(zhì)的一個質(zhì)點粒子，而是作為離散顆粒的形式用來求解擬流體模型和描述顆粒性質(zhì)的變化歷程。

固體顆粒相的控制方程中，連續(xù)性方程、動量方程、擬溫度守恒方程分別為：

式中：ρs和vs分別為顆粒相密度和速度，P 為連續(xù)相壓力，Ps為離散相壓力，τs為顆粒相粘性應力張量，ρsg為外部體積力，Rsg為相間相互作用力；Ts為顆粒相擬溫度，I為單位矩陣，（－psI＋τs）∶!vs為由顆粒相應力產(chǎn)生的能量，kTs!Ts為能量耗散項，kTs為能量耗散系數(shù)，γTs為顆粒間碰撞產(chǎn)生的能量耗散項，φgs為連續(xù)相與顆粒相間的能量交換，具體參數(shù)表達式詳見文獻［14］。

固體顆粒相SPH離散方程組中，連續(xù)性方程、動量方程、擬溫度守恒方程分別為：

式中：應力σ＝－psI＋τs，ρi為SDPH粒子i的密度（即顆粒相有效密度）；ρs為顆粒的實際密度；速度矢量vij＝vi－vj，R′sphgs為作用于SPH粒子上的單位質(zhì)量曳力，fbpi為作用于粒子i上的壁面力；擬溫度梯度!Ts的SPH離散公式為：

式中：ρij＝ （ρi＋ ρj）／2。

2 算例驗證分析

本章算例中室溫及周圍環(huán)境的初始溫度均設為273K，時間積分采用蛙跳格式，時間步長為0.1μs。假設碰撞過程中，靶板不受側面邊界約束，光滑長度采用1.5倍粒子間距。

2.1 鋁球超高速碰撞鋁薄板

因為航天器的防護結構多為鋁合金，所以本文參照Hiermaier等［10］的實驗中鋁球超高速碰撞鋁薄板的過程進行數(shù)值計算，并結合Zhou等［12］和Liu等［13］采用傳統(tǒng)SPH方法及ASPH方法的計算結果進行對比驗證。

如圖1所示，鋁球直徑為10mm，鋁靶板尺寸為40mm×40mm×4mm。鋁球的撞擊速度為6.18km／s。SPH 粒子直徑為0.67mm，整個模型共有23 391個粒子，其中鋁球離散為1 791個粒子，鋁薄板離散為21 600個粒子。人工粘性參數(shù)同Hiermaier［10］和Zhou保持一致。由于Hiermaier并沒有列出所用鋁合金型號，本文根據(jù)密度、屈服應力等強度模型參數(shù)確定其類別，其他參數(shù)參考同類材料［18－19］，具體材料參數(shù)參見表1，其中：CV為定容比熱。狀態(tài)方程選用同 Hiermaier等［10］和Zhou等［12］相同的Tillotson狀態(tài)方程，具體參數(shù)見表2。

圖1 算例模型結構Fig.1Experimental model

表1 本構模型參數(shù)Table 1Parameters of constitutive model

表2 狀態(tài)方程參數(shù)Table 2Parameters of state equation

圖2為采用SPH新方法計算得到的碎片分布與Hiermaier實驗［10］、Hiermaier二維數(shù)值模擬結果［10］、Zhou等采用傳統(tǒng)SPH方法的三維數(shù)值模擬結果［12］和Liu等采用ASPH方法的二維數(shù)值模擬結果［13］對比圖。從圖2中可以看出，采用SPH新方法和傳統(tǒng)SPH方法中的內(nèi)核碎片云均位于外泡碎片云的前方且占了整個碎片云質(zhì)量的大部分。但通過對比，可以發(fā)現(xiàn)存在兩處不同：

第一處不同為采用傳統(tǒng)SPH方法的內(nèi)核碎片云成圓錐型，隨著外泡碎片云的不斷膨脹錐面不斷彎曲，分布比較集中。而采用SPH新方法的將處于失效狀態(tài)的碎片等效于擬流體，充分考慮碎片間相互作用及數(shù)值模擬中經(jīng)常被忽略的氣體相對碎片的影響，故內(nèi)核碎片云呈“圓球”型，分布比較分散，與實驗吻合較好；

第二處是采用SPH新方法的計算結果中反濺碎片云比采用傳統(tǒng)SPH方法的反濺碎片云膨脹的距離和寬度都小，這是因為SPH方法引入Johnson－Cook損傷模型后，薄板屈服應力要小于未引入損傷模型的屈服應力，如式（9）所示；根據(jù)反濺碎片云形成過程研究［20］，當鋁球與薄板接觸點后移速度高于薄板材料反濺速度向后的分量時，鋁球阻擋薄板該部分材料繼續(xù)反濺，反之無法阻止薄板材料反濺；而薄板的屈服應力直接決定了其反濺速度，所以采用SPH新方法的計算結果中反濺碎片云比采用傳統(tǒng)SPH方法中的反濺碎片云膨脹的距離和寬度都??；但從反濺碎片云形態(tài)來看，新方法計算得到的反濺碎片云運動形態(tài)與實驗更為接近，而傳統(tǒng)SPH方法中反濺碎片云的計算結果明顯過大。

表3給出了采用SPH新方法得到的彈坑直徑、外泡碎片云和內(nèi)核碎片云的形狀、分布及其與使用傳統(tǒng)SPH方法［10］、ASPH方法的模擬結果［13］的對比，其中：d1和d2分別為不包括和包括彈坑邊緣的彈坑直徑，l為碎片云的膨脹距離，w 為碎片云寬度，Δ為l／w 的相對誤差。從表3可以看出，采用SPH新方法得到的結果與實驗吻合較好。與使用傳統(tǒng)SPH方法、ASPH方法的模擬結果進行對比，新方法在內(nèi)核碎片云形狀和分布上計算結果更加準確，與實驗更加符合。表明SPH新方法適合模擬鋁球沖擊鋁薄板等超高速碰撞問題。

圖2 不同SPH方法數(shù)值模擬結果與實驗結果對比Fig.2Comparison of different SPH algorithms‘simulation results and the experimental result

表3 超高速碰撞數(shù)值模擬結果對比Table 3Comparison of high－velocity impact simulation results

為形象展示了碎片云形成、膨脹的過程，圖3～4分別以俯視、左視兩種視角展示了不同時刻鋁球超高速撞擊鋁薄板三維數(shù)值模擬過程，其中圖4以左視角展示了“碎片云”形成、發(fā)展及最終形態(tài)，并與Zhou采用傳統(tǒng)SPH方法的三維數(shù)值模擬結果［12］進行對比。

圖3 鋁彈超高速碰撞鋁薄板俯視圖Fig.3Top view of hypervelocity impact of projectile on thin plates

圖4 SPH新方法與傳統(tǒng)SPH方法計算結果對比Fig.4Comparison of the new SPH algorithm and traditional SPH algorithm‘s simulation results

2.2 Whipple防護屏超高速碰撞數(shù)值模擬研究

柳森等［21］進行了一系列有關 Whipple防護屏的超高速碰撞實驗，本文中選取正碰撞04－0080，04－0090兩組典型過程［21］對其進行數(shù)值模擬研究。圖5為Whipple屏結構示意圖。為了減小計算量，在不影響防護屏彈孔尺寸及艙壁彈坑分布的前提下，分別將防護屏和艙壁長寬減少至50和80mm。SPH粒子直徑為0.3mm，彈丸直徑為5mm，共離散為2 440個粒子，撞擊速度為5.29和6.15km／s。防護屏和艙壁厚度均為2mm，分別離散為169 344和427 734個粒子。彈丸、防護屏及艙壁材料為LY12鋁合金。Johnson－Cook損傷模型的具體參數(shù)［22］見表4。由于本算例的撞擊速度處于超高速撞擊范疇中相對較低的撞擊速度范圍，碰撞過程中所有材料處于固體狀態(tài)，且Grüneisen狀態(tài)方程膨脹態(tài)實際上是固相材料的線性膨脹，且沖擊絕熱線關系與高壓固體狀態(tài)方程之間聯(lián)系密切，沖擊絕熱關系的實驗數(shù)據(jù)較多［23］。因此本算例中狀態(tài)方程采用 Grüneisen 狀 態(tài) 方 程，LY12 鋁 合 金 式 中 各 系 數(shù) 為ρ＝2 790kg／m3，Γ ＝2.0，CS＝5 328，SS＝1.338。

圖5 Whipple防護屏結構Fig.5Construction of Whipple shield

表4 LY12鋁合金本構模型參數(shù)Table 4Parameters of LY12aluminium alloy of constitutive model

圖6給出了撞擊速度為5.29km／s時Whipple屏損傷的實驗和數(shù)值模擬結果，圖7給出了撞擊速度為6.15km／s時Whipple屏損傷的數(shù)值模擬結果?？梢钥闯?，艙壁的損傷分布集中在撞擊中心，損傷有彈坑、斑點和侵蝕三種類型。撞擊中心處的中心損傷區(qū)域彈坑密集重疊，損傷程度較重，主要由內(nèi)核碎片云撞擊形成，中心損傷區(qū)外圍出現(xiàn)了由外泡碎片云撞擊形成若干個環(huán)形分布的圓形撞擊坑群，損傷程度較輕。將撞擊速度5.29和6.15km／s的艙壁損傷情況進行對比，撞擊速度值越大，彈丸破碎程度更大，彈坑及斑點數(shù)量越多，撞擊中心所在的中心損傷區(qū)的損傷值越小。這表明在該速度范圍內(nèi)，撞擊速度值越大，對艙壁的損傷破壞越小，Whipple防護屏防護效果越好，這符合Whipple防護結構的典型撞擊極限曲線中粉碎段（中速段）的變化規(guī)律。由計算結果與實驗結果對比可見，防護屏彈孔尺寸、艙壁中心撞擊坑損傷區(qū)的尺寸和艙壁彈坑形狀、分布與實驗吻合較好，驗證了新方法對超高速碰撞問題數(shù)值模擬的有效性，具體數(shù)值參見表5。

表5 采用SPH新方法的數(shù)值模擬結果與實驗結果比較Table 5Comparison of the new SPH algorithm and experimental results

圖6 撞擊速度為5.29km／s時 Whipple屏損傷情況Fig.6Damage characteristics of the Whipple shield at impact velocity of 5.29km／s

圖7 撞擊速度為6.15km／s時 Whipple屏損傷情況Fig.7Damage characteristics of the Whipple shield at impact velocity of 6.15km／s

由于實驗中存在靶板安裝和撞擊點散布誤差引起的碰撞角度誤差 ，防護屏彈孔直徑的計算結果與實驗結果的誤差處在合理區(qū)間；此外，實驗中并沒有直接給出中心損傷區(qū)域?qū)挾?，本文中對實驗結果中進行短尺測量估計，計算與測量結果在容許誤差內(nèi)；因為數(shù)值模擬層裂破壞存在一定困難，在計算過程中只考慮了不穿孔和穿孔兩種損傷狀態(tài)，實驗中的剝落等層裂損傷狀態(tài)的數(shù)值模擬有待進一步研究。

3 結 論

引入顆粒動力學理論（擬流體模型）建立適用于超高速碰撞的SPH新方法，將SDPH方法的應用領域拓展至超高速碰撞領域，對球形彈丸超高速碰撞薄板形成碎片云的過程和Whipple防護屏超高速碰撞問題進行了數(shù)值模擬，得到以下結論：

（1）對球形彈丸超高速碰撞薄板形成碎片云的基礎算例進行了數(shù)值模擬，得到的彈坑直徑、碎片云形態(tài)及內(nèi)核碎片云的形狀和分布均與實驗吻合較好。與使用傳統(tǒng)SPH方法、自適應光滑粒子流體動力學（ASPH）方法的仿真結果進行對比，新方法的計算結果在內(nèi)核碎片云形狀和分布上更加準確；驗證了新方法在求解超高速碰撞問題的可行性和準確性；

（2）通過Whipple屏超高速碰撞的算例，模擬了不同撞擊速度下防護屏彈坑尺寸及艙壁損傷情況等內(nèi)容；計算結果與實驗結果對比吻合較好，符合Whipple防護結構的典型撞擊極限曲線，驗證了新方法對球形彈丸超高速碰撞Whipple防護結構問題的有效性；為下一步研究彈丸斜撞擊Whipple屏全過程以及撞擊損傷的內(nèi)在機理等打下了工作基礎。

［1］ Klinkrad H.Collision risk analysis for LEO［J］.Advances in Space Research，1993，13（2）：105－112.

［2］ Whipple F L.Meteorites and space travel［J］.Astronomical Journal，1947，52（13）：132－137.

［3］ Piekutowski A J.Debris clouds produced by the hypervelocity impact of nonspherical projectiles［J］.International Journal of Impact Engineering，2001，26（1－10）：613－624.

［4］ 管公順，龐寶君，哈躍，等.鋁合金 Whipple防護結構高速撞擊實驗研究［J］.爆炸與沖擊，2005，25（5）：461－466.Guan Gongshun，Pang Baojun，Ha Yue，et al.Experimantal investigation of high－velocity impact on aluminum alloy Whipple shield［J］.Explosion and Shock Waves，2005，25（5）：461－466.

［5］ Swift H F，Bamford R，Chen R.Designing space vehicle shields for meteoroid protection a new analysis［J］.Advances in Space Research，1983，2（12）：219－234.

［6］ Lucy L B.A numerical approach to the testing of the fission hypothesis［J］.Astronomical Journal，1977，82：1013－1024.

［7］ Gingold R A，Monaghan J J.Smoothed particle hydrodynamics：theory and application to non－spherical stars［J］.Monthly Notices of the Royal Astronomical Society，1977，181（3）：375－389.

［8］ Libersky L D，Petscheck A G.Smoothed particle hydrodynamics with strength of materials［C］∥ Trease H，F(xiàn)ritts J，Crowley W.Proceedings of The Next Free Lagrange Conference.Springer Verlag，New York，1991，395：248－257.

［9］ Johnson G R，Petersen E H，Stryk R A.Incorporation of an SPH option into the EPIC code for a wide range of high velocity impact computations［J］.International Journal of Impact Engineering，1993，14（1－4）：385－394.

［10］ Hiermaier S，K?nke D，Stilp A J，et al.Computational simulation of the hypervelocity impact of Al－spheres on thin plates of different materials［J］.International Journal of Impact Engineering，1997，20（1－5）：363－374.

［11］ 徐金宏，湯文輝，羅永.光滑粒子模擬法在超高速碰撞現(xiàn)象中的應用［J］.爆炸與沖擊，2006，26（1）：53－58.Xu Jinhong，Tang Wenhui，Luo Yong.Applications of the smoothed particle hydrodynamics method to hypervelocity impact simulations［J］.Explosion and Shock Waves，2006，26（1）：53－58.

［12］ Zhou C E，Liu G R，Lou K Y.Three－dimensional penetration simulation using smoothed particle hydrodynamics［J］.International Journal of Computational Methods，2007，4（4）：671－691.

［13］ Liu G R，Liu M B.Smoothed particle hydrodynamics：a meshfree particle method［M］.Singapore：Word Scientific Publishing Co.Pte.Ltd.，2003.

［14］ 陳福振，強洪夫，高巍然.風沙運動問題的SPH－FVM 耦合方法數(shù)值模擬［J］.物理學報，2014，63（13）：130202.Chen Fuzhen，Qiang Hongfu，Gao Weiran.Simulation of aeolian sand transport with SPH－FVM coupled method［J］.Acta Physica Sinica，2014，63（13）：130202.

［15］ 陳福振，強洪夫，高巍然.氣－粒兩相流傳熱問題的光滑離散顆粒流體動力學方法數(shù)值模擬［J］.物理學報，2014，63（23）：230206.Chen Fuzhen，Qiang Hongfu，Gao Weiran.Numerical simulation of heat transfer in gas－particle two－phase flow with smoothed discrete particle hydrodynamics［J］.Acta Physica Sinica，2014，63（23）：230206.

［16］ 陳福振，強洪夫，高巍然.燃料爆炸拋撒成霧及其爆轟過程的SDPH方法數(shù)值模擬研究［J］.物理學報，2015，64（11）：110202.Chen Fuzhen，Qiang Hongfu，Gao Weiran.Numerical simulation of fuel dispersal into cloud and its combustion and explosion with smmoothed discrete particle hydrodymamics［J］.Acta Physica Sinica，2015，64（11）：110202.

［17］ Johnson G R，Cook W H.A constitutive model and data for metals subjected to large strains，high strain rates and high temperatures［C］∥Proceedings of the Seventh International Symposium on Ballistics.Hague，Netherlands，1983：541－547

［18］ 李春雷.2A12鋁合金本構關系實驗研究［M］.哈爾濱：哈爾濱工業(yè)大學，2006.

［19］ 費思聰，孫秦.鋁合金的J－C失效參數(shù)標定與仿真分析［J］.計算機仿真，2013，30（9）：46－50.Fei Sicong，Sun Qin.Investigation on parameters calibration for the J－C failure model of aluminum alloy［J］.Computer Simulation，2013，30（9）：46－50.

［20］ 遲潤強.彈丸超高速撞擊薄板碎片云建模研究［D］.哈爾濱：哈爾濱工業(yè)大學，2010.

［21］ 柳森，李毅，黃潔，等.用于驗證數(shù)值模擬仿真的 Whipple屏超高速撞擊試驗結果［J］.宇航學報，2005，26（4）：505－508.Liu Sen，Li Yi，Huang Jie，et al.Hypervelocity impact test results of Whipple shield for the validation of numerical simulations［J］Journal of Astronautics，2005，26（4）：505－508.

［22］ 侯日立，周平，彭建祥 沖擊波作用下LY12鋁合金結構損傷的數(shù)值模擬［J］.爆炸與沖擊，2012，32（5）：470－474.Hou Rili，Zhou Ping，Peng Jianxiang，Numerical simulation of shock damage of LY12aluminium alloy structure［J］.Explosion and Shock Waves，2012，32（5）：470－474.

［23］ 徐金中.基于SPH方法的空間碎片超高速碰撞特性及其防護結構設計研究［D］.長沙：國防科學技術大學，2008.

A new smoothed particle hydrodynamics method based on the pseudo－fluid model and its application in hypervelocity impact of
aprojectile on a thin plate

Qiang Hongfu1，F(xiàn)an Shujia1，2，Chen Fuzhen1，Liu Hu1
（1.Power Engineering Department，Rocket Force University of Engineering，Xi’an 710025，Shaanxi，China；2.Unit 96656，PLA，Beijing100085，China）

Based on the kinetic theory of granular flow （pseudo－fluid model），a new Smoothed Particle Hydrodynamics（SPH）algorithm suited for hyper velocity collision was presented in this paper.The damaged debris of the hyper velocity impact was equated with the pseudo－fluid and the effects of the debris’interaction and the effects of the gas on the formation process of the debris cloud were investigated.The new SPH algorithm was employed to simulate the 3Dhyper velocity impact of an alloy projectile on thin target plates，and the numerical results of crater diameters，the structure and morphology characteristics of the debris cloud and the core debris cloud’s shape and distribution were in good agreement with the experimental results.Compared with the simulations of the standard SPH and ASPH，the simulation of the new algorithm is more accurate in the core debris cloud’s characteristics.Meanwhile，the hyper velocity impact of the Whipple shield problem was also simulated at different impact velocities.It was found that the crater diameters and the damage characteristics of the rear walls agree well with the experimental results，and that the simulation results are consistent with the typical ballistic limit curve of a Whipple shielding structure.

smoothed particle hydrodynamics；debris cloud；hypervelocity impact；Whipple shield

O389 國標學科代碼：1303530

A

10.11883／1001－1455（2017）06－0990－11

2016－07－09；

2016－12－12

國家自然科學基金項目（51276192）

強洪夫（1963－ ），男，博士，教授，博士生導師；通信作者：范樹佳，fan＿shu＿jia＠163.com。

（責任編輯 王小飛）