王 濤，汪 兵，林健宇，柏勁松，李 平，鐘 敏，陶 鋼

（1.中國(guó)工程物理研究院流體物理研究所，四川 綿陽(yáng) 621999；2.南京理工大學(xué)能源與動(dòng)力工程學(xué)院，江蘇 南京 210094）

當(dāng)?shù)兔芏攘黧w在慣性力作用下加速高密度流體時(shí)，界面上的預(yù)制初始擾動(dòng)會(huì)增長(zhǎng)，這種現(xiàn)象稱(chēng)為Rayleigh-Taylor（RT）不穩(wěn)定性[1-2]。同樣，在金屬材料界面也會(huì)發(fā)生RT不穩(wěn)定性，產(chǎn)生與流體界面不穩(wěn)定性類(lèi)似的尖釘和氣泡結(jié)構(gòu)，但是與流體界面不穩(wěn)定性的最大不同之處在于金屬材料具有強(qiáng)度效應(yīng)，而研究表明金屬材料強(qiáng)度能夠抑制擾動(dòng)增長(zhǎng)[3-6]。并且，金屬的強(qiáng)度效應(yīng)也會(huì)使金屬界面不穩(wěn)定性問(wèn)題比流體界面不穩(wěn)定性更復(fù)雜，需要考慮的因素更多，比如金屬材料的本構(gòu)方程、損傷斷裂、相變、熔化、微結(jié)構(gòu)特征等。金屬界面不穩(wěn)定性的發(fā)生常常預(yù)示著材料經(jīng)歷高溫、高壓、高應(yīng)變率等極端環(huán)境，因而常存在于慣性約束聚變[7-9]、超新星爆炸[10]、小行星撞擊[11-12]、地球內(nèi)核運(yùn)動(dòng)和板塊構(gòu)造[13-14]以及爆炸焊接等領(lǐng)域，因此具有重要的研究意義。

目前，金屬材料界面不穩(wěn)定性的理論研究主要利用能量平衡[15-17]或力平衡[18-20]原理，推導(dǎo)擾動(dòng)增長(zhǎng)所滿足的色散關(guān)系，即線性分析；由于這些線性分析方法均采用理想塑性本構(gòu)關(guān)系，因此預(yù)測(cè)擾動(dòng)增長(zhǎng)受到較大的限制。Bai等[21]基于能量平衡原理，采用Steinberg-Guinan（SG）和Johnson-Cook（JC）本構(gòu)模型及變壓力的準(zhǔn)等熵加載歷程，針對(duì)有限厚度金屬平板問(wèn)題，推導(dǎo)了相應(yīng)的擾動(dòng)增長(zhǎng)方程，對(duì)線性分析方法進(jìn)行了拓展，使其能夠用于平面幾何中金屬材料不穩(wěn)定性線性至非線性段增長(zhǎng)的預(yù)測(cè)，并對(duì)爆轟驅(qū)動(dòng)和等離子體驅(qū)動(dòng)的金屬界面不穩(wěn)定性實(shí)驗(yàn)[6,20,22-23]進(jìn)行了分析，結(jié)果和實(shí)驗(yàn)吻合較好，但是這種線性分析方法仍然存在較大的局限性，如目前還無(wú)法考慮加載過(guò)程等。

金屬材料RT不穩(wěn)定性的實(shí)驗(yàn)研究始于20世紀(jì)70年代，首次由Barnes等[22,24]通過(guò)炸藥爆轟準(zhǔn)等熵驅(qū)動(dòng)鋁平板實(shí)現(xiàn)，并利用高能X射線照相觀察界面擾動(dòng)增長(zhǎng)。之后的平面爆轟驅(qū)動(dòng)金屬材料RT不穩(wěn)定性實(shí)驗(yàn)[25-27]均采用相似的裝置。美國(guó)和俄羅斯開(kāi)展了大量有關(guān)金屬材料界面不穩(wěn)定性的研究工作，包括界面擾動(dòng)增長(zhǎng)規(guī)律及影響因素等。Olson等[28]通過(guò)爆轟驅(qū)動(dòng)實(shí)驗(yàn)，研究了銅材料的微結(jié)構(gòu)和機(jī)加工對(duì)RT不穩(wěn)定性增長(zhǎng)的影響，發(fā)現(xiàn)單晶晶向和材料加工導(dǎo)致的應(yīng)變硬化會(huì)影響擾動(dòng)增長(zhǎng)，而多晶材料的晶粒尺寸和晶界強(qiáng)化在所研究的加載條件下對(duì)擾動(dòng)增長(zhǎng)無(wú)影響。Wang等[27]通過(guò)爆轟驅(qū)動(dòng)實(shí)驗(yàn)和數(shù)值模擬，研究了純鋁的RT不穩(wěn)定性增長(zhǎng)規(guī)律和強(qiáng)度因素（基于SG本構(gòu)模型的初始剪切模量和初始屈服強(qiáng)度）對(duì)擾動(dòng)增長(zhǎng)的影響。何長(zhǎng)江等[29]通過(guò)數(shù)值模擬研究了初始擾動(dòng)波長(zhǎng)對(duì)平面爆轟驅(qū)動(dòng)鋁RT不穩(wěn)定性增長(zhǎng)的影響。郝鵬程等[30]和劉軍等[31]通過(guò)對(duì)內(nèi)爆加載金屬鋼殼的數(shù)值模擬，研究了材料強(qiáng)度（剪切模量和基于理想彈塑性本構(gòu)的屈服強(qiáng)度以及是否考慮彈塑性）和初始波長(zhǎng)對(duì)界面擾動(dòng)增長(zhǎng)的影響。平面爆轟準(zhǔn)等熵加載的壓力一般在30 GPa左右，為了提高加載壓力，F(xiàn)rahan等[26]采用兩級(jí)炸藥爆轟驅(qū)動(dòng)技術(shù)，即第一級(jí)炸藥驅(qū)動(dòng)鐵飛片撞擊起爆第二級(jí)炸藥，再加載樣品，使加載壓力提高到50 GPa，由此研究了該加載壓力下鈹?shù)腞T不穩(wěn)定性問(wèn)題，并結(jié)合數(shù)值模擬對(duì)材料不同本構(gòu)模型的適用性進(jìn)行了檢驗(yàn)。為了研究極端加載壓力和應(yīng)變率條件下的金屬界面不穩(wěn)定性及材料混合行為，人們又發(fā)展了電磁加載[32-33]和激光加載技術(shù)[5-6,23]，并開(kāi)展了大量研究工作，如Park等[5-6,34]利用Omega激光器，結(jié)合數(shù)值模擬，研究了釩和鉭分別在100 GPa和350 GPa的準(zhǔn)等熵加載壓力下的擾動(dòng)增長(zhǎng)行為。目前，美國(guó)國(guó)家點(diǎn)火裝置NIF已經(jīng)能夠?qū)崿F(xiàn)超過(guò)1 TPa的準(zhǔn)等熵加載[35]。

由于材料強(qiáng)度對(duì)金屬界面不穩(wěn)定性發(fā)展有抑制作用，所以近年來(lái)人們又將其作為一種手段來(lái)研究金屬材料的強(qiáng)度—擾動(dòng)增長(zhǎng)法[4,25,36]。Lorenz等[23]利用基于激光準(zhǔn)等熵加載技術(shù)的RT擾動(dòng)增長(zhǎng)法研究了峰值壓力為20 GPa時(shí)6061-T6鋁的屈服強(qiáng)度，結(jié)果顯示該條件下6061-T6鋁的屈服強(qiáng)度比正常條件下高出3.6倍。Park等[5-6]利用激光準(zhǔn)等熵加載實(shí)驗(yàn)和數(shù)值模擬研究發(fā)現(xiàn)了釩在100 GPa加載條件下具有一種新的RT不穩(wěn)定性致穩(wěn)機(jī)制—聲子拖曳導(dǎo)致有效晶格黏性增大，進(jìn)而導(dǎo)致釩的強(qiáng)度增大，對(duì)RT不穩(wěn)定性的致穩(wěn)作用增強(qiáng)；而且研究顯示釩的強(qiáng)度增大了3.5倍，而350 GPa加載壓力下鉭的強(qiáng)度增大了13倍[34]。Belof等[37]通過(guò)爆轟驅(qū)動(dòng)鐵的RT不穩(wěn)定性實(shí)驗(yàn)和數(shù)值模擬，確認(rèn)了鐵在經(jīng)歷α→ε→α′相變過(guò)程中強(qiáng)度大大提高。

本研究利用自研的爆轟與沖擊動(dòng)力學(xué)歐拉計(jì)算程序和SG本構(gòu)模型，針對(duì)爆轟驅(qū)動(dòng)金屬錫RT不穩(wěn)定性增長(zhǎng)對(duì)兩類(lèi)模型初始參數(shù)的敏感性進(jìn)行數(shù)值模擬分析，一類(lèi)是樣品的初始參數(shù)，包括初始擾動(dòng)振幅、初始擾動(dòng)波長(zhǎng)和樣品初始厚度，另一類(lèi)是SG本構(gòu)模型的初始參數(shù)。

1 計(jì)算方法

對(duì)于考慮爆轟和材料彈塑性行為的多物質(zhì)、大變形及強(qiáng)沖擊波物理問(wèn)題，發(fā)展了高保真度的爆轟與沖擊動(dòng)力學(xué)歐拉有限體積計(jì)算程序，其守恒型控制方程組為

式中：下標(biāo)i、j分別代表x、y、z3個(gè)方向，遵循張量運(yùn)算法則，V為控制體體積，S為控制體表面積，ni為外法向單位矢量，ρ、uk(k=i,j)、p、E分別為密度、速度、壓強(qiáng)和單位質(zhì)量總能量，sij為偏應(yīng)力張量分量。

采用維數(shù)分裂技術(shù)，將式(1)描述的物理問(wèn)題分解為多個(gè)一維問(wèn)題。對(duì)于每個(gè)一維問(wèn)題，采用PPM（Piecewise parabolic method）方法對(duì)單元內(nèi)物理量的分布進(jìn)行插值和重構(gòu)，然后通過(guò)歐拉型兩步算法進(jìn)行計(jì)算，即物理量的Lagrange推進(jìn)求解，再將拉氏網(wǎng)格上的物理量映射到靜止的歐拉網(wǎng)格上。材料強(qiáng)度效應(yīng)、炸藥爆轟過(guò)程和人工黏性在Lagrange步中實(shí)現(xiàn)。多物質(zhì)界面采用體積分?jǐn)?shù)方法進(jìn)行捕捉。

數(shù)值模擬中，炸藥和金屬材料分別采用JWL和Mie-Grüneisen狀態(tài)方程（EOS）描述

式中：v=ρ0/ρ為相對(duì)比容，T為溫度，α、σ、R1、R2、ω為常數(shù)，為單位體積內(nèi)能。

金屬材料采用SG本構(gòu)模型描述。SG本構(gòu)模型是在彈塑性本構(gòu)方程中引入壓力、溫度和應(yīng)變率項(xiàng)，并且壓力與應(yīng)變率對(duì)流動(dòng)應(yīng)力的耦合效應(yīng)具有可分離變量特性。因?yàn)镾G本構(gòu)模型中流動(dòng)應(yīng)力依賴于壓力，所以材料的本構(gòu)方程與狀態(tài)方程之間存在某種耦合關(guān)系。這一耦合關(guān)系反映了高壓下金屬材料的壓力硬化特性。SG本構(gòu)模型給出的動(dòng)態(tài)屈服強(qiáng)度YSG和剪切模量G分別為

式中：Y0和G0分別為初始屈服強(qiáng)度和剪切模量，β和n分別為材料應(yīng)變硬化系數(shù)和硬化指數(shù)，A為壓力硬化系數(shù)，η=ρ/ρ0為材料壓縮比，B為熱軟化系數(shù)。

2 計(jì)算模型和條件

簡(jiǎn)化的二維計(jì)算模型見(jiàn)圖1，炸藥為JO-9159，樣品為金屬錫。炸藥高度為30 mm；樣品直徑為32 mm，厚度3 mm，樣品面向炸藥一側(cè)的界面上預(yù)置了初始振幅a0= 0.1 mm、初始波長(zhǎng)λ0= 3 mm的正弦形式初始擾動(dòng)。炸藥和樣品之間的真空間隙尺寸為10 mm。該真空間隙保證了爆轟產(chǎn)物經(jīng)過(guò)真空間隙加載樣品的過(guò)程中壓力連續(xù)增大，形成輕介質(zhì)對(duì)重介質(zhì)的準(zhǔn)等熵加載，進(jìn)而誘發(fā)RT不穩(wěn)定性；它的另一個(gè)作用是確保在高壓下樣品溫度依然低于熔化溫度[23]。

圖1 二維計(jì)算模型Fig.1 Two dimensional computational model

JO-9159炸藥的JWL狀態(tài)方程參數(shù)見(jiàn)表1，錫的Mie-Grüneisen狀態(tài)方程參數(shù)和SG本構(gòu)模型標(biāo)準(zhǔn)參數(shù)見(jiàn)表2和表3。其中：pCJ和DCJ分別為炸藥的爆壓和爆速，Ymax為材料最大屈服強(qiáng)度。

表1 JO-9159炸藥JWL狀態(tài)方程參數(shù)Table 1 EOS parameters of JO-9159 explosive

表2 錫的Mie-Grüneisen狀態(tài)方程參數(shù)Table 2 Mie-Grüneisen EOS parameters of Sn

表3 錫的SG本構(gòu)模型參數(shù)Table 3 SG constitutive model parameters of Sn

3 結(jié)果分析

3.1 程序確認(rèn)

在開(kāi)展金屬錫RT不穩(wěn)定性數(shù)值模擬研究之前，先通過(guò)對(duì)Lindquist等[38]的爆轟驅(qū)動(dòng)T-6061鋁的RT不穩(wěn)定性實(shí)驗(yàn)進(jìn)行模擬，確認(rèn)本計(jì)算程序的可靠性。該實(shí)驗(yàn)采用HMX炸藥爆轟加載1.5 mm厚的鋁樣品，加載壓力峰值為30 GPa。初始擾動(dòng)波長(zhǎng)λ0= 2 mm，初始擾動(dòng)振幅a0= 0.15 mm, 0.11 mm。圖2給出了本研究通過(guò)自研歐拉程序計(jì)算得到的擾動(dòng)振幅a和文獻(xiàn)中采用ARES程序計(jì)算的結(jié)果，及其與實(shí)驗(yàn)結(jié)果的比較（橫軸y為自由面運(yùn)動(dòng)位移）。可以看出，吻合度較好，表明自研的歐拉程序模擬這類(lèi)爆轟驅(qū)動(dòng)金屬材料界面不穩(wěn)定性是準(zhǔn)確的。

圖2 Lindquist等爆轟驅(qū)動(dòng)鋁實(shí)驗(yàn)的擾動(dòng)振幅比較Fig.2 Comparison of perturbation amplitudes of Lindquist et al.'s experiments driven by explosion

3.2 樣品初始參數(shù)敏感性分析

首先，考察模擬結(jié)果對(duì)計(jì)算網(wǎng)格的收斂性問(wèn)題。圖3和圖4分別給出了網(wǎng)格尺寸為10.0、5.0和2.5μm時(shí)加載壓力剖面和擾動(dòng)振幅增長(zhǎng)曲線。可以看出，2.5 μ m時(shí)模擬結(jié)果已經(jīng)達(dá)到收斂，所以后續(xù)的所有計(jì)算均采用2.5 μ m的計(jì)算網(wǎng)格。從圖3還可以看出：加載壓力p是連續(xù)光滑增大的，近似形成對(duì)界面的準(zhǔn)等熵加載；之后由于爆轟產(chǎn)物膨脹壓力持續(xù)減小及由自由面多次反射的稀疏波抵達(dá)擾動(dòng)界面形成持續(xù)的卸載作用，加載壓力逐漸減小。

圖3 不同網(wǎng)格尺寸時(shí)的加載壓力剖面Fig.3 Loading pressure profiles for different grid size

圖4 不同網(wǎng)格尺寸時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.4 Perturbation amplitude growth for different grid size

接下來(lái)，通過(guò)數(shù)值模擬研究金屬錫樣品初始參數(shù)，包括初始擾動(dòng)振幅、初始擾動(dòng)波長(zhǎng)、樣品初始厚度對(duì)爆轟驅(qū)動(dòng)金屬錫RT不穩(wěn)定性增長(zhǎng)的影響。在研究初始振幅的影響時(shí)，保持初始波長(zhǎng)λ0為3 mm，樣品初始厚度h為3 mm，初始振幅a0分別取0.3、0.1、0.05、0.025和0.012 5 mm。圖5顯示了不同初始振幅時(shí)的擾動(dòng)振幅增長(zhǎng)曲線?？梢钥闯觯弘S著初始振幅的減小，擾動(dòng)振幅增長(zhǎng)也減弱；當(dāng)初始振幅減小到一定值后，擾動(dòng)振幅不再增長(zhǎng)，即爆轟驅(qū)動(dòng)金屬錫RT不穩(wěn)定性增長(zhǎng)存在一個(gè)截止初始振幅。

圖5 不同初始振幅時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.5 Perturbation amplitude growth for different initial amplitude

在研究初始波長(zhǎng)的影響時(shí)，保持初始振幅為0.1 mm，樣品初始厚度為3 mm，而初始波長(zhǎng)λ0分別取9、6、3、1和0.5 mm。圖6給出了不同初始波長(zhǎng)時(shí)的擾動(dòng)振幅增長(zhǎng)曲線。從初始波長(zhǎng)為9 mm開(kāi)始，隨著初始波長(zhǎng)的減小，擾動(dòng)振幅逐漸增大；但是，當(dāng)初始波長(zhǎng)減小到一定程度后，擾動(dòng)振幅增長(zhǎng)又隨著初始波長(zhǎng)的減小而減小。所以爆轟驅(qū)動(dòng)金屬錫RT不穩(wěn)定性增長(zhǎng)存在一個(gè)最不穩(wěn)定的波長(zhǎng)，或者增長(zhǎng)速度最快的模態(tài)。

在研究樣品初始厚度的影響時(shí)，保持初始振幅為0.1 mm，初始波長(zhǎng)為3 mm，而樣品初始厚度h分別取3、5、7和9 mm。圖7顯示了樣品初始厚度不同時(shí)的擾動(dòng)振幅增長(zhǎng)曲線。隨著樣品初始厚度的增大，擾動(dòng)增長(zhǎng)減弱，表明樣品厚度增大可以抑制擾動(dòng)增長(zhǎng)；而且當(dāng)樣品厚度增大到一定程度時(shí)，擾動(dòng)增長(zhǎng)停止，即爆轟驅(qū)動(dòng)金屬錫RT不穩(wěn)定性增長(zhǎng)存在一個(gè)樣品截止厚度。

圖6 不同初始波長(zhǎng)時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.6 Perturbation amplitude growth for different initial wavelength

圖7 不同樣品初始厚度時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.7 Perturbation amplitude growth for different initial thickness of sample

3.3 SG本構(gòu)模型參數(shù)敏感性分析

研究金屬錫RT不穩(wěn)定性增長(zhǎng)對(duì)SG本構(gòu)模型參數(shù)的敏感性時(shí)，計(jì)算模型的初始厚度取3 mm，初始振幅取0.1 mm，初始波長(zhǎng)取3 mm。SG本構(gòu)模型的可調(diào)參數(shù)包括應(yīng)變硬化系數(shù)β、應(yīng)變硬化指數(shù)n、壓力硬化系數(shù)A和熱軟化系數(shù)B，對(duì)于錫材料，它們的標(biāo)準(zhǔn)值分別為2 000、0.06、0.086 6 GPa-1、2 120 K-1（見(jiàn)表3）。在研究其中一個(gè)參數(shù)的影響時(shí)，其余參數(shù)均取標(biāo)準(zhǔn)值。材料應(yīng)變硬化系數(shù)β依次取1 000、1 500、2 000、2 500、3 000；應(yīng)變硬化指數(shù)n依次取 0.01、0.03、0.06、0.12、0.24；壓力硬化系數(shù)A依次取 0.021 7、0.043 3、0.086 6、0.173 2、0.346 4 GPa-1；熱軟化系數(shù)B依次取 0.53×10?3、1.06×10?3、2.12×10?3、4.24×10?3、8.48×10?3K-1。

圖8～圖11分別給出了改變錫SG本構(gòu)模型的應(yīng)變硬化系數(shù)β、應(yīng)變硬化指數(shù)n、壓力硬化系數(shù)A和熱軟化系數(shù)B時(shí)的擾動(dòng)振幅增長(zhǎng)曲線?？梢钥闯觯涸诒Z準(zhǔn)等熵驅(qū)動(dòng)下金屬材料發(fā)生復(fù)雜流動(dòng)，此時(shí)材料應(yīng)變硬化系數(shù)和應(yīng)變硬化指數(shù)的變化對(duì)擾動(dòng)振幅增長(zhǎng)的影響很小，即這兩個(gè)SG本構(gòu)模型參數(shù)的改變對(duì)材料屈服強(qiáng)度的影響較小；而壓力硬化系數(shù)和熱軟化系數(shù)的變化對(duì)擾動(dòng)增長(zhǎng)的影響顯著，隨著壓力硬化系數(shù)的增大，擾動(dòng)增長(zhǎng)受到的抑制作用越顯著，而隨著熱軟化系數(shù)的減小，擾動(dòng)增長(zhǎng)受到抑制效果的變化卻逐漸減小，即熱軟化系數(shù)小到一定程度時(shí)其變化將很難改變材料屈服強(qiáng)度。根據(jù)文獻(xiàn)分析可知，高壓高應(yīng)變率加載條件下，金屬材料的強(qiáng)度一般都會(huì)增大，而標(biāo)準(zhǔn)參數(shù)的本構(gòu)模型一般都低估了材料強(qiáng)度。因此基于SG本構(gòu)模型，采用擾動(dòng)增長(zhǎng)法預(yù)估高壓高應(yīng)變率下錫的材料強(qiáng)度，修正壓力硬化系數(shù)以獲得合理的強(qiáng)度是一條合理的途徑。

圖8 應(yīng)變硬化系數(shù)不同時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.8 Perturbation amplitude growth for different strain hardening coefficient

圖9 應(yīng)變硬化指數(shù)不同時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.9 Perturbation amplitude growth for different strain hardening exponent

圖10 壓力硬化系數(shù)不同時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.10 Perturbation amplitude growth for different pressure hardening coefficient

圖11 熱軟化系數(shù)不同時(shí)的擾動(dòng)振幅增長(zhǎng)曲線Fig.11 Perturbation amplitude growth for different thermal softening coefficient

4 結(jié) 論

開(kāi)展了爆轟驅(qū)動(dòng)金屬錫RT不穩(wěn)定性增長(zhǎng)對(duì)樣品初始參數(shù)和SG本構(gòu)模型初始參數(shù)敏感性的數(shù)值模擬分析，得出以下結(jié)論。

（1）樣品初始參數(shù)（初始振幅、初始波長(zhǎng)、樣品初始厚度）對(duì)爆轟驅(qū)動(dòng)錫RT不穩(wěn)定性增長(zhǎng)有重要影響。RT不穩(wěn)定性增長(zhǎng)隨著初始振幅的減小而減小，且存在一個(gè)截止初始振幅；存在一個(gè)最不穩(wěn)定的模態(tài)（波長(zhǎng)），當(dāng)初始波長(zhǎng)大于該波長(zhǎng)時(shí)，RT不穩(wěn)定性增長(zhǎng)隨著初始波長(zhǎng)的減小而增大，反之，RT不穩(wěn)定性增長(zhǎng)隨著初始波長(zhǎng)的減小而減??；樣品厚度的增大可以抑制RT不穩(wěn)定性增長(zhǎng)，而且存在一個(gè)樣品截止厚度。

（2）金屬錫的RT不穩(wěn)定性增長(zhǎng)對(duì)其SG本構(gòu)模型應(yīng)變硬化系數(shù)和應(yīng)變硬化指數(shù)的變化不敏感，而對(duì)壓力硬化系數(shù)和熱軟化系數(shù)比較敏感，但是從采用擾動(dòng)增長(zhǎng)法預(yù)估材料強(qiáng)度的角度來(lái)說(shuō)，修正壓力硬化系數(shù)以獲得合理的錫材料強(qiáng)度是合理的途徑。