爆轟加載下金屬錫層裂破碎數(shù)值模擬*

賀年豐，任國(guó)武，陳永濤，郭昭亮

（中國(guó)工程物理研究院流體物理研究所，四川 綿陽(yáng) 621999）

在爆轟等強(qiáng)加載中，金屬樣品在強(qiáng)沖擊波壓縮后又經(jīng)歷持續(xù)稀疏波拉伸作用，可能發(fā)生熔化狀態(tài)下的層裂破碎，形成不斷膨脹的準(zhǔn)連續(xù)空間分布碎片云團(tuán)[1-3]。目前對(duì)材料在固體狀態(tài)下傳統(tǒng)層裂問(wèn)題的研究已經(jīng)比較成熟[4]。但若金屬樣品在加-卸載過(guò)程中進(jìn)入熔化狀態(tài)，由于材料強(qiáng)度急劇下降，在稀疏波拉伸作用下的破碎霧化現(xiàn)象與傳統(tǒng)層裂現(xiàn)象明顯不同[3]，且目前該方面的研究和認(rèn)識(shí)相對(duì)不足。從基礎(chǔ)科學(xué)研究角度，層裂破碎現(xiàn)象是沖擊加-卸載波與材料性能耦合的復(fù)雜動(dòng)力學(xué)問(wèn)題；在武器研究領(lǐng)域，層裂破碎的形成將改變金屬飛層的動(dòng)力學(xué)響應(yīng)及行為，進(jìn)而影響武器內(nèi)爆性能。由于強(qiáng)沖擊加-卸載下的層裂破碎問(wèn)題在沖擊動(dòng)力學(xué)和武器物理研究中的重要理論研究?jī)r(jià)值和強(qiáng)烈需求背景，該問(wèn)題引起了研究人員的密切關(guān)注。

早在1914年，Hopkinson[5]就研究了強(qiáng)動(dòng)載下鋼板的層裂破壞現(xiàn)象。對(duì)于錫、鉛等一類低強(qiáng)度、低熔點(diǎn)金屬，強(qiáng)沖擊加載下的層裂破碎和霧化噴射現(xiàn)象更明顯。1983年Andriot等[1]發(fā)現(xiàn)了爆轟加載錫熔化后不同于固體層裂的破碎現(xiàn)象。Holtkamp等[3]利用多幅質(zhì)子照相技術(shù)給出了爆轟加載下不同厚度金屬錫對(duì)應(yīng)的從熔化破碎到固體層裂破碎清晰對(duì)比圖像。Signor等[6]和Resseguier等[7]利用激光加載產(chǎn)生的三角波對(duì)錫的熔化破碎開(kāi)展了研究，并通過(guò)回收技術(shù)、X光斷層攝影技術(shù)等手段給出了破碎顆粒的尺寸分布規(guī)律，此外也驗(yàn)證了Grady層裂模型的適用性。陳永濤等[8-9]、Chen等[10]也利用中能X光照相和DPS（Doppler pin system）測(cè)速給出了錫樣品的質(zhì)量、密度空間分布和自由面速度剖面。張林等[11]則通過(guò)改進(jìn)Asay窗技術(shù)來(lái)診斷主體層裂破碎中不同特征區(qū)信息。雖然目前實(shí)驗(yàn)診斷技術(shù)的發(fā)展大大促進(jìn)了對(duì)低熔點(diǎn)金屬層裂破碎的認(rèn)識(shí)，但在診斷早中期樣品內(nèi)部應(yīng)力波系相互作用以及與材料失效耦合過(guò)程信息上面臨較大的難度。微觀機(jī)制研究方面，Soulard等[12]利用分子動(dòng)力學(xué)模擬研究了加-卸載中熔化狀態(tài)下金屬的層裂破碎，發(fā)現(xiàn)熔化狀態(tài)下材料層裂強(qiáng)度將急劇下降，并指出了表面層裂片的存在。Xiang等[13-14]模擬了鉛的層裂破碎，指出早期損傷失效是由空穴化主導(dǎo)的。但由于時(shí)-空尺度的巨大差異，以及分子動(dòng)力學(xué)模型與真實(shí)材料在微介觀結(jié)構(gòu)上的差異，微觀模擬分析的置信度以及與宏觀物理現(xiàn)象間的關(guān)聯(lián)性仍有待進(jìn)一步研究。由于爆轟加載下低熔點(diǎn)金屬的層裂破碎涉及相變、層裂碎裂等應(yīng)力波與材料物性耦合復(fù)雜問(wèn)題，宏觀數(shù)值模擬研究難點(diǎn)較多，目前公開(kāi)文獻(xiàn)資料較少。比較詳細(xì)的有Signor等[6]采用多相物態(tài)方程和拉伸層裂失效模型對(duì)錫的沖擊熔化問(wèn)題開(kāi)展了一維模擬以預(yù)測(cè)碎片尺寸分布；曹結(jié)東等[15]基于提出的熔化度損傷模型開(kāi)展了二維軸對(duì)稱爆轟驅(qū)動(dòng)錫熔化破碎模擬，得到了和實(shí)驗(yàn)結(jié)果可對(duì)照的面密度分布圖像?？梢?jiàn)，目前的宏觀模擬集中于一、二維模擬，且偏重于碎片尺寸及密度空間分布等特征分析；而三維狀態(tài)層裂破碎過(guò)程中，復(fù)雜波系與材料耦合作用過(guò)程及機(jī)制仍缺乏相關(guān)研究分析。

本文中，采用光滑粒子流體動(dòng)力學(xué)(smoothed particle hydrodynamics, SPH)方法[16-17]，并結(jié)合Steinberg-Guinan高壓本構(gòu)[18]以及Grady層裂模型[19]針對(duì)爆轟加載低熔點(diǎn)錫樣品實(shí)驗(yàn)開(kāi)展相應(yīng)的三維宏觀數(shù)值模擬，計(jì)算結(jié)果中自由面速度和中后期的碎片質(zhì)量分布演化均與實(shí)驗(yàn)結(jié)果基本吻合。進(jìn)一步，通過(guò)對(duì)不同厚度樣品中內(nèi)部壓力、速度的時(shí)空演化以及破碎物質(zhì)各特征區(qū)對(duì)比分析，剖析了爆轟加載層裂破碎問(wèn)題中加載三角波及自由面反射稀疏波、側(cè)面稀疏波相互作用下層裂破碎演化機(jī)制。

1 計(jì)算方法和模型

1.1 SPH方法

SPH方法是目前應(yīng)用較廣泛的一種無(wú)網(wǎng)格粒子類方法，其采用一組拉格朗日粒子來(lái)表示連續(xù)體，粒子自身攜帶所有物理量。由于采用核函數(shù)插值來(lái)構(gòu)建空間任一點(diǎn)上物理量的近似值，SPH方法不需要有限元的網(wǎng)格信息，避免了拉格朗日方法中網(wǎng)格畸變導(dǎo)致的數(shù)值求解困難，因此SPH方法能方便地處理斷裂、層裂、破碎等材料破壞問(wèn)題[17]。

1.2 材料本構(gòu)和狀態(tài)方程

針對(duì)所研究問(wèn)題的強(qiáng)動(dòng)載屬性，采用適用于高應(yīng)變率和高壓情況下的Steinberg-Guinan（SG）本構(gòu)模型以及Mie-Grüneisen狀態(tài)方程。其中SG本構(gòu)模型中剪切模量G和屈服強(qiáng)度Y如下：

式中：η為壓縮比，εp為等效塑性應(yīng)變，T為材料溫度(K)，p為壓強(qiáng)，G0和Y0分別為初始狀態(tài)下的剪切模量和屈服強(qiáng)度；屈服強(qiáng)度還滿足為最大屈服強(qiáng)度；上述式子中其他變量皆為材料參數(shù)。另外，假設(shè)熔化狀態(tài)下的材料強(qiáng)度可忽略。

Mie-Grüneisen狀態(tài)方程：

能描述絕大部分金屬固體在沖擊載荷下的熱力學(xué)行為，其中μ=ρ/ρ0-1，ρ和ρ0分別為材料密度和初始密度，c0為初始波速，S1為材料參數(shù)，γ為Mie-Grüneisen系數(shù)，e為質(zhì)量?jī)?nèi)能。

采用JWL狀態(tài)方程來(lái)描述炸藥爆轟產(chǎn)物壓力pe、內(nèi)能ee、 相對(duì)體積ηe之間的關(guān)系：

式中：A、B、ω、R1、R2為炸藥參數(shù)。

1.3 Grady層裂模型

為了避免對(duì)損傷演化過(guò)程的復(fù)雜描述，Grady[19]從能量平衡的角度對(duì)材料損傷破碎問(wèn)題進(jìn)行理論分析，并推導(dǎo)了層裂強(qiáng)度、層裂時(shí)間、碎片平均尺寸的解析表達(dá)式。假設(shè)材料動(dòng)態(tài)破碎過(guò)程中消耗的能量W由彈性能U和局部動(dòng)能E提供，且有約束條件s≤2c0t，其中s為碎片尺寸，能量平衡可表示為：

對(duì)于延性材料，塑性斷裂耗能可近似為W=Yεc， εc為臨界孔洞分?jǐn)?shù)，則有：

式中為應(yīng)變率。將式(6)代入約束條件解得層裂強(qiáng)度：

對(duì)于液體材料，則耗能項(xiàng)由新增表面能和黏性耗能構(gòu)成，能量平衡可表示為：

式中：τ為表面能，η為黏性系數(shù)。

由于局部動(dòng)能項(xiàng)和黏性耗散項(xiàng)相對(duì)較小，通?？珊雎?，推導(dǎo)得到液體材料層裂強(qiáng)度：

需要指出的是，Grady模型從能量平衡出發(fā)，由能量準(zhǔn)則判斷該狀態(tài)是否發(fā)生層裂失效，而未顯式描述早期微損傷萌生及生長(zhǎng)過(guò)程。由于早期微損傷演化時(shí)間相對(duì)本文中所關(guān)注的層裂破碎過(guò)程是較小量，且微損傷影響范圍為其周邊小區(qū)域，其對(duì)整體波系演化的影響較小。因此，針對(duì)本文中波系與材料宏觀相互作用過(guò)程的模擬研究，Grady模型是可適用的。

2 爆轟加載平面錫樣品數(shù)值模擬

2.1 模型和參數(shù)

本研究對(duì)應(yīng)的實(shí)驗(yàn)裝置和簡(jiǎn)化后的計(jì)算模型如圖1所示，根據(jù)對(duì)稱性，實(shí)際計(jì)算采用四分之一模型，炸藥底面中心一點(diǎn)起爆。高能炸藥(high explosive,HE)半徑為12.5 mm，厚度為30 mm。錫（Sn）樣品半徑為12.5 mm，厚度和實(shí)驗(yàn)一致，分別為4、6、8 mm。鋁（Al）支座內(nèi)徑為12.5 mm，厚度為3 mm。詳細(xì)參數(shù)材料見(jiàn)表1～3。

根據(jù)陳永濤等[8-10]和Holtkamp等[3]的實(shí)驗(yàn)結(jié)果，爆轟加載下錫樣品已發(fā)生了卸載熔化。實(shí)際上，熔化相變和層裂失效都發(fā)生在反射稀疏波卸載過(guò)程中很短時(shí)間內(nèi)，且兩者相互耦合。為簡(jiǎn)化模擬計(jì)算復(fù)雜度，可假設(shè)卸載過(guò)程中樣品立即進(jìn)入卸載熔化，再發(fā)生熔化狀態(tài)下的層裂失效。本爆轟實(shí)驗(yàn)中應(yīng)變率約為105s-1量級(jí)，而金屬錫 τ=0.54 J/m2，可由式（9）得到熔化狀態(tài)錫的層裂強(qiáng)度為65.2 MPa。鋁支座在爆轟加載中將保持固相，可依照一般延性材料取εc=0.15[19]，利用式（7）計(jì)算得到其層裂強(qiáng)度為2.48 GPa。以上層裂強(qiáng)度和文獻(xiàn)[19]中實(shí)驗(yàn)數(shù)據(jù)基本一致。

圖1 實(shí)驗(yàn)裝置示意圖以及計(jì)算模型圖Fig.1 Schematics of experimental device and simulation

表1 Sn和Al的SG本構(gòu)參數(shù)Table 1 The material parameters in SG constitutive relation for Sn and Al

表2 Sn和Al的Mie-Grüneisen狀態(tài)方程參數(shù)Table 2 The material parameters in Mie-Grüneisen equation of state for Sn and Al

表3 高能炸藥JWL狀態(tài)方程參數(shù)Table 3 The parameters in JWL equation of state for high explosive

2.2 數(shù)值模擬結(jié)果及分析

對(duì)不同厚度的錫樣品開(kāi)展數(shù)值模擬。在柱形高能炸藥底面中心一點(diǎn)起爆爆轟加載中，由于爆轟波在傳播過(guò)程中受到邊側(cè)稀疏的影響，加載沖擊波壓力呈類倒錐形分布。以6 mm厚（h=6 mm）的錫樣品為例，如圖2所示，波陣面下方和兩側(cè)分別為入射稀疏和邊側(cè)稀疏區(qū)。沖擊波到達(dá)自由面之前，沿厚度方向的壓力剖面演化如圖3所示?？梢?jiàn)加載沖擊波為典型三角波形，且在向樣品內(nèi)部推進(jìn)中峰值逐漸減小。

圖2 沖擊波到達(dá)自由表面時(shí)刻壓力分布Fig.2 Pressure distribution in the Sn specimen when the shock wave arriving at the free surface

圖3 樣品中心軸線位置沿厚度方向的壓力剖面演化Fig.3 Evolution of pressure distribution along the central axis of the Sn specimen

由于樣品徑向尺寸大于厚度，在加載的早期階段，可認(rèn)為樣品中心軸附近區(qū)域未受到邊側(cè)稀疏波的影響。以6 mm厚的錫樣品為例，中間點(diǎn)的壓力和速度的演化如圖4所示，其中包含入射沖擊加載段（a）、入射稀疏段（b）、反射稀疏段（c）。圖4中pmax、vmax為該質(zhì)點(diǎn)受到的沖擊加載壓力和速度峰值，其大小表征了入射沖擊波的強(qiáng)度；Δp、Δv為自由面反射稀疏作用下壓力和速度下降值，其表征了所受反射稀疏波的強(qiáng)度。樣品中心軸上不同（初始）厚度位置的壓力和速度演化見(jiàn)圖5，可見(jiàn)入射沖擊波在向內(nèi)部推進(jìn)中逐漸減弱，同時(shí)反射稀疏波從自由界面向內(nèi)部傳播中由于多層裂發(fā)生迅速減弱。因此，隨著質(zhì)點(diǎn)靠近自由表面，其入射稀疏段縮短，而所受反射稀疏波的卸載強(qiáng)度增大；反之則相反。

圖4 錫樣品中間點(diǎn)的壓力和速度的演化Fig.4 Evolutions of pressure and velocity at the centre of the Sn specimen

圖5 錫樣品中心軸上不同厚度處的壓力和速度的演化Fig.5 Evolution of pressure and velocity at different sites along the central axis of the Sn specimen

不同厚度錫樣品的模擬結(jié)果中碎裂物質(zhì)分布以及與實(shí)驗(yàn)結(jié)果真實(shí)比例對(duì)比見(jiàn)圖6，樣品自由表面速度比較見(jiàn)表4。總體上，模擬結(jié)果和實(shí)驗(yàn)結(jié)果基本吻合，這也驗(yàn)證了本文所采用的層裂模型和熔化狀態(tài)下相關(guān)參數(shù)的合理性。在爆轟加載下，4 mm厚的錫樣品殘留主體部分出現(xiàn)明顯凸起；6 mm厚的樣品中間輕微凸起；8 mm厚的樣品則為下凹。限于X光照相的分辨能力，自由表面附近的較低密度碎片云無(wú)法被捕捉。此外，實(shí)驗(yàn)圖像中碎片呈現(xiàn)隨機(jī)分布特征，而模擬結(jié)果為分層破碎特征，主要原因在于真實(shí)材料和加載具有隨機(jī)性和不均勻性。受限于求解規(guī)模，模擬中錫樣品的空間離散尺度為0.1 mm，而實(shí)際破碎顆粒能達(dá)到微米量級(jí)，因此模擬結(jié)果中的每個(gè)分散層裂片代表的是周邊一個(gè)破碎物質(zhì)帶。進(jìn)一步觀察碎片云分布可見(jiàn)其中包含3個(gè)不同區(qū)域。以6 mm厚的錫樣品為例，見(jiàn)圖6(f)，其中A為反射稀疏波（主導(dǎo)）區(qū)；B為入射稀疏波（主導(dǎo)）區(qū)；C為邊側(cè)稀疏波（主導(dǎo)）區(qū)。對(duì)于4 mm厚的錫樣品，A區(qū)直接連通至底部而無(wú)B區(qū)。隨著厚度的增大，A區(qū)逐漸上移，并與B區(qū)交匯于樣品內(nèi)部，同時(shí)C區(qū)增大。

圖6 爆轟加載下不同厚度錫樣品碎片云分布的數(shù)值模擬結(jié)果與實(shí)驗(yàn)結(jié)果[9]的比較Fig.6 Comparisons of fragment distributions between simulation results and experimental results[9] for the tin specimens with different thicknesses under explosive loading

表4 不同方法得到的不同厚度樣品自由面速度Table 4 Free-surface velocities by different methods for the tin specimens with different thicknesses

不同厚度的錫樣品中心軸上離底面0.5 mm處的速度時(shí)間歷程見(jiàn)圖7。由圖7可見(jiàn)，當(dāng)樣品厚度較薄時(shí)，雖然在沖擊波之后短時(shí)間內(nèi)入射稀疏波將使得質(zhì)點(diǎn)速度下降，但很快自由表面反射稀疏波到達(dá)樣品底面并導(dǎo)致其速度上升。同時(shí)由于邊側(cè)稀疏波的影響，兩側(cè)速度上升較弱。因此，4 mm厚的錫樣品主體底部呈現(xiàn)上凸形狀。隨著厚度增加，反射稀疏波作用迅速減弱。同時(shí)反射稀疏波到達(dá)樣品底部時(shí)間延后，導(dǎo)致底部經(jīng)歷的入射稀疏作用時(shí)間更長(zhǎng)，速度下降幅度更大。因此，隨著厚度增加，樣品底部變?yōu)橄掳?。反射卸載波到達(dá)樣品底面后同一時(shí)刻不同厚度樣品的內(nèi)部速度分布見(jiàn)圖8，在入射、反射和邊側(cè)稀疏波共同作用下，4 mm厚的樣品在中心線兩側(cè)形成低速區(qū)；6 mm厚的樣品在底面中間形成一片較平整的低速區(qū)；而8 mm厚的樣品的低速區(qū)更大，且靠近中心軸速度更低。可見(jiàn)隨著樣品厚度增加，邊側(cè)和入射稀疏波的影響增強(qiáng)，同時(shí)反射稀疏波伴隨層裂推進(jìn)而迅速變?nèi)?，底部的低速區(qū)向內(nèi)部擴(kuò)大。這也可解釋圖6(f)中A、B、C區(qū)域在不同厚度錫樣品中物理量以及樣品底部形貌的演化。在自由面附件區(qū)域，沖擊波在自由面反射稀疏波并與入射稀疏波相互作用而引發(fā)層裂，并向內(nèi)部推進(jìn)。由于類倒錐形的入射波壓力分布（圖2），碎片云中反射稀疏波主導(dǎo)區(qū)也呈類倒錐形。

圖7 不同厚度錫樣品中心軸上離底面0.5 mm處速度演化Fig.7 The velocity evolution at x0=0.5 mm in the Sn specimens with different thicknesses

圖8 起爆后6.4 μs時(shí)刻不同厚度Sn樣品內(nèi)部速度分布Fig.8 Velocity distributions in the Sn specimens with different thicknesses at 6.4 μs after detonation

3 結(jié) 論

對(duì)爆轟加載下低熔點(diǎn)金屬錫的層裂破碎問(wèn)題開(kāi)展了數(shù)值模擬。計(jì)算結(jié)果中自由面速度和中后期的碎片質(zhì)量分布演化均與實(shí)驗(yàn)結(jié)果基本吻合。模擬結(jié)果也反映了自由面附近低密度物質(zhì)的存在。錫樣品的層裂破碎是爆炸沖擊波后反射稀疏波、邊側(cè)稀疏波、入射稀疏波與材料性能共同作用的結(jié)果。碎片云分布中也明顯可見(jiàn)各波影響主導(dǎo)區(qū)域。對(duì)于薄樣品，反射稀疏波短時(shí)間內(nèi)到達(dá)樣品底部，而兩側(cè)由于邊側(cè)稀疏波的干擾，速度上升較弱，因此主體底部呈現(xiàn)上凸形貌。隨著樣品厚度的增大，在向內(nèi)傳播過(guò)程中，伴隨層裂發(fā)生，反射稀疏波強(qiáng)度迅速減弱，且樣品底部經(jīng)歷的入射稀疏作用過(guò)程更長(zhǎng)，因此樣品底部將逐漸變?yōu)橄掳?。另外，隨樣品厚度的增大，反射稀疏波的影響主導(dǎo)區(qū)將上移，而邊側(cè)和入射稀疏波影響區(qū)擴(kuò)大。樣品的層裂破碎是自由面反射稀疏波與入射稀疏波共同作用的結(jié)果，且由自由面附近向內(nèi)推進(jìn)。由于入射波的類倒錐形壓力分布，破碎區(qū)也將出現(xiàn)了類倒錐形碎片分布。需要指出的是，本文中認(rèn)為樣品在稀疏波到達(dá)時(shí)即刻進(jìn)入了卸載熔化狀態(tài)，而實(shí)際上沖擊加-卸載中金屬的相變?yōu)閯?dòng)態(tài)演化過(guò)程，該問(wèn)題也將是本文后續(xù)研究的關(guān)注點(diǎn)之一。