爆轟加載下彈塑性固體Richtmyer-Meshkov流動的擾動增長規(guī)律?

殷建偉 潘昊 吳子輝 郝鵬程 胡曉棉

1)(北京理工大學(xué)機電學(xué)院,北京 100081)

2)(北京應(yīng)用物理與計算數(shù)學(xué)研究所,計算物理重點實驗室,北京 100094)

3)(中國工程物理研究院研究生部,北京 100088)

1 引 言

沖擊波加載不同密度的流體界面時,界面擾動將增長演化為渦結(jié)構(gòu)并最終形成湍流混合層,這種流體動力學(xué)現(xiàn)象稱為Richtymer-Meshkov(RM)流動不穩(wěn)定性[1?4],大量研究表明,流體的表面張力、黏性等耗散性質(zhì)將削弱RM不穩(wěn)定性的增長.對于固體介質(zhì),其RM流動過程受到固體物性的顯著影響,呈現(xiàn)出不同的物理現(xiàn)象:在沖擊波作用下,材料的強度致穩(wěn)效應(yīng)將導(dǎo)致擾動振幅增長一定程度后圍繞某均值周期性振蕩[5?9],若沖擊波壓力極高,如慣性約束聚變(ICF)內(nèi)爆中加載靶丸的壓力達到100 GPa量級[10],擾動將劇烈增長,伴隨升溫、破碎、噴射、多物質(zhì)混合等復(fù)雜物理現(xiàn)象,是武器內(nèi)爆動力學(xué)、物質(zhì)界面微噴射與混合動力學(xué)過程重點關(guān)注的問題.

目前對固體中RM流動的理論研究基于彈塑性材料均勻性假設(shè),未考慮材料的相變、斷裂等復(fù)雜物理因素.Piriz[8]提出了剛塑性材料擾動面上尖釘振幅的極大值表達式,并在Buttler等[11,12]和Dimonte等[13]的實驗中得到初步驗證.Jensen等[14]在平面撞擊實驗中利用該表達式推算了鈰在極端加載條件下的強度,結(jié)果表明不同方法得到的強度結(jié)果偏差約20%—50%.因此,如何改進極值振幅表達式,更準確地描述彈塑性固體RM流動的擾動增長規(guī)律是目前亟待解決的問題,也是理解強度致穩(wěn)機制的核心問題,有助于提高對界面失穩(wěn)早期階段的物理認識.

本文采用AFE2D程序?qū)φㄋ幈Z加載下金屬材料的界面擾動增長過程進行數(shù)值模擬研究,定性分析了初始擾動的振幅與波長對擾動增長的影響以及材料強度的致穩(wěn)效應(yīng),利用相關(guān)沖擊波關(guān)系式及流體RM不穩(wěn)定性的研究成果,進一步推導(dǎo)得到了極值振幅的理論關(guān)系式,該式能夠描述擾動增長被抑制時振幅與初始擾動、材料性質(zhì)之間的關(guān)系,通過分析關(guān)系式中各項的相關(guān)性,初步得到了線性關(guān)系式的適用范圍.

2 界面擾動增長問題的數(shù)值模擬

考慮如圖1所示的界面擾動問題,高能炸藥爆轟驅(qū)動金屬樣品,該問題是Buttler等[11]開展的銅、錫等材料RM實驗的簡化模型,其中高能炸藥厚度12 mm,金屬樣品厚6 mm,樣品中的沖擊波向介質(zhì)/真空界面?zhèn)鞑r加載擾動自由面,擾動面為正弦形,其初始波長為λ0,初始振幅為ξ0.

圖1 沖擊波加載界面擾動問題示意圖Fig.1.Illustration of the surface perturbation problem loaded by the shock wave.

炸藥爆轟采用時間起爆方式,狀態(tài)方程為JWL類型狀態(tài)方程[15]

其中R1,R2和ω是JWL模型參數(shù);A和B是JWL模型系數(shù),由高能炸藥的CJ爆速DCJ、CJ爆壓PCJ及JWL參數(shù)確定;V=ρ/ρ0是相對比容,ρ0是初始密度,ρ是材料密度;E是單位初始體積的內(nèi)能.固體材料的狀態(tài)方程為Mie-Grüneisen狀態(tài)方程[16]

其中μ=ρ/ρ0?1;γ0為Grüneisen系數(shù);b為Grüneisen系數(shù)的一階體積修正量;c0,S1,S2和S3為材料沖擊波-粒子速度關(guān)系式的系數(shù).固體材料本構(gòu)采用Prandtl-Reuss關(guān)聯(lián)塑性流動律和von Mises屈服準則描述[7?9],即

下標(biāo)采用Einstein求和約定,Dij為應(yīng)變率張量;Sij為偏應(yīng)力張量;G為材料剪切模量,Y為屈服強度,采用彈性理想塑性模型,保持G和Y不變.幾種常見金屬材料的狀態(tài)方程參數(shù)與剪切模量見表1,炸藥的爆轟參數(shù)與JWL狀態(tài)方程參數(shù)見表2.

利用AFE2D程序數(shù)值模擬PBX9501炸藥爆轟驅(qū)動銅飛片的實驗(見文獻[11]Cu-pRad0426實驗),如圖1所示,樣品的自由面上引入kξ0(k=2π/λ0,初始擾動的波數(shù))分別為0.12,0.35,0.75和1.5的初始擾動,初始波長均為550μm,銅樣品中的沖擊波壓力約36 GPa.圖2為實驗的數(shù)值模擬結(jié)果,可以看到?jīng)_擊波加載擾動面時,初始波谷位置首先啟動,形成向外的突起尖釘(Spike)結(jié)構(gòu);初始波峰位置稍后啟動,與尖釘?shù)炔课婚g的速度差導(dǎo)致其向材料內(nèi)部凹陷,形成氣泡(Bubble)結(jié)構(gòu);反相后尖釘與氣泡結(jié)構(gòu)進一步演化,初始擾動kξ0=0.12時尖釘基本不增長,kξ0=0.35時尖釘增長一段時間后被抑制,kξ0=0.75和1.5時增長不受抑制,尖釘演化為射流.

表1 金屬材料的狀態(tài)方程參數(shù)和剪切模量Table 1.Parameters of equation of state and shear modulus of several selected metals.

表2 PBX9501炸藥的爆轟參數(shù)與JWL狀態(tài)方程參數(shù)Table 2.Parameters of the detonation performance and JWL equation of state of the high explosives PBX9501.

定義尖釘?shù)恼穹蜸P為尖釘與無擾動區(qū)域(Bulk)之間在沖擊波傳播方向的相對位移,尖釘?shù)脑鲩L速度vSP為尖釘與無擾動區(qū)域之間在沖擊波傳播方向的相對速度,即

氣泡的振幅ξBU和增長速度vBU的定義類似,

模擬結(jié)果顯示,kξ0=0.12時尖釘?shù)淖畲笳穹?max(ξSP)=10.55 μm,相比初始振幅=11μm幾乎沒有增長,但相位出現(xiàn)反轉(zhuǎn);kξ0=0.35時=155.3 μm. 在Buttler對CupRad0426實驗數(shù)據(jù)的分析中[11],kξ0=0.12時=(7±4)μm,kξ0=0.35時=160 μm.表3對比了不同kξ0時尖釘增長速度的極大值=max(vSP),模擬結(jié)果與實驗數(shù)據(jù)基本一致,表明本文采用的程序、模型和參數(shù)適用于界面擾動問題,可用于分析界面擾動的增長規(guī)律.

表3 PBX9501炸藥驅(qū)動銅飛片實驗的模擬結(jié)果Table 3.Simulation results of the Cu-pRad0426 experiment where the copper fl yer plate was driven by the high explosive PBX9501.

圖2 Cu-pRad0426實驗的模擬結(jié)果 (a)初始擾動界面;(b)沖擊波首次加載界面后0.7μs時刻的模擬結(jié)果;(c)3.7μs時刻的模擬結(jié)果;(d)3.7μs時刻實驗質(zhì)子照相圖像[11]Fig.2.The simulation results of the Cu-pRad0426 experiment:(a)Initial state;(b)plot of the simulation results at 0.7μs which relative to shockwave breakout at the free surface;(c)plot of the simulation results at 3.7μs;(d)proton radiographic picture of the experiment at 3.7μs[11].

3 界面擾動的增長規(guī)律

3.1 初始形態(tài)對擾動增長的影響

首先分析擾動初始形態(tài)對后期增長的影響,利用無量綱數(shù)kξ0表征初始擾動的振幅與波長之比,考慮圖1所示的界面擾動問題,選取初始波長λ0=450—600 μm,初始振幅ξ0=15—38 μm,則kξ0的變化范圍0.20—0.40,樣品材料為OFHC銅,屈服強度Y=0.47 GPa,炸藥為PBX9501.圖3比較了沖擊波首次加載界面后3.1μs時刻不同kξ0的擾動增長情況,尖釘?shù)恼穹S著kξ0增加逐步變大,尖釘形狀在kξ0=0.20時仍接近正弦擾動形態(tài),繼續(xù)提高kξ0尖釘將逐漸拉伸、變窄,最終演化為射流.上述結(jié)果表明,對于相同的加載條件和樣品材料,若初始擾動的振幅與波長比值越小,相應(yīng)的初始界面越“平坦”,擾動增長越不明顯,界面保持穩(wěn)定;振幅與波長比值越大,相應(yīng)的初始界面越“溝壑”,擾動越易于增長,甚至形成局部的微射流,界面失穩(wěn).

圖3 不同kξ0值的界面擾動增長情況,沖擊波首次加載界面后3.1μs時刻Fig.3.Growth results of the perturbated surfaces with different values of kξ0,3.1 μs relative to shockwave breakout at the free surface.

3.2 材料強度對擾動增長的致穩(wěn)效應(yīng)

相關(guān)研究表明,金屬材料的動態(tài)屈服強度存在較強的應(yīng)變硬化效應(yīng),并且在高應(yīng)變率下表現(xiàn)出顯著的應(yīng)變率硬化效應(yīng)[17?21].本文所考慮的沖擊波加載擾動面問題中,擾動面上形成的尖釘和氣泡等特征結(jié)構(gòu)均為局部高應(yīng)變區(qū),局部應(yīng)變率達到106—107s?1.然而目前高應(yīng)變率范圍內(nèi)金屬材料的強度數(shù)據(jù)極少,對應(yīng)變率的硬化規(guī)律認識亦有待進一步完善,本文所采用的彈性理想塑性模型雖不能反映材料真實的應(yīng)變與應(yīng)變率硬化效應(yīng),但在均勻性假設(shè)下通過調(diào)節(jié)屈服強度模擬不同的硬化效果也不失為一種可行的近似方法,能夠用于近似闡述材料強度的致穩(wěn)效應(yīng).

選擇擾動的kξ0=0.251,取OFHC銅的屈服強度為Y=0.17 GPa,令其逐步增加至0.57 GPa以模擬不同的硬化效應(yīng).圖4為擾動界面的增長情況,強度較低時擾動面的尖釘演化為射流,提高材料強度后尖釘增長受到抑制,逐漸從尖銳的射流形態(tài)退化為類似正弦波形的擾動形態(tài),定性地顯示出材料強度對擾動增長的致穩(wěn)效應(yīng),表明在相同的加載和初始擾動條件下,強度越高的材料擾動增長越不明顯.

3.3 擾動增長規(guī)律的線性近似

由圖3和圖4可知,提高初始擾動的kξ0,擾動面的尖釘振幅將逐漸增大,在爆轟加載下擾動增長受到材料強度的顯著抑制作用,因此當(dāng)擾動增長被抑制時,尖釘增長幅度的控制因素應(yīng)包含初始擾動kξ0與強度Y的某種組合形式.

圖4 擾動增長情況,調(diào)節(jié)屈服強度以模擬不同程度的硬化效應(yīng),沖擊波首次加載界面后3.1μs時刻Fig.4. Growth results of the surfaces with varied yield strength for the different equivalent effects of hardening,3.1μs relative to shockwave breakout at the free surface.

利用(4)式定義的尖釘振幅、增長速度和Dimonte等[7,13]的分析方法,類似于尖釘最大振幅和最大增長速度的定義,得到氣泡最大振幅=max(ξBU)和氣泡最大增長速度=max(vBU),利用初始波數(shù)無量綱化尖釘和氣泡的最大振幅

利用材料初始密度ρ0、屈服強度Y無量綱化尖釘和氣泡的最大增長速度,分別標(biāo)記為無量綱數(shù)和

圖5 (網(wǎng)刊彩色)尖釘和氣泡的最大振幅與最大增長速度的關(guān)系(銅強度分別取0.27,0.32和0.37 GPa)Fig.5. (color online)Relations between the nondimensional maximum amplitude and growth velocity of the spikes and bubbles(yield strength of copper varied as 0.27,0.32 and 0.37 GPa).

圖5和圖6分別顯示了銅、鋁樣品中尖釘和氣泡的最大振幅與最大增長速度無量綱數(shù)之間的關(guān)系,選擇初始擾動kξ0變化范圍0.15—0.3以避免尖釘演化為射流,屈服強度分別取三組值.模擬結(jié)果顯示,無量綱化的尖釘最大振幅與最大增長速度在一定范圍內(nèi)表現(xiàn)出近似線性關(guān)系,斜率C約為0.3,即

圖6 (網(wǎng)刊彩色)尖釘和氣泡的最大振幅與最大增長速度關(guān)系(鋁強度分別取0.24,0.29和0.34 GPa)Fig.6. (color online)Relations between the nondimensional maximum amplitude and growth velocity of the spikes and bubbles(yield strength of aluminum varied as 0.24,0.29 and 0.34 GPa).

Piriz理論預(yù)測的斜率值為0.29[7],Dimonte利用PAGOSA程序擬合的數(shù)據(jù)點見圖5[13],忽略截距后斜率值為0.24,本文結(jié)果中線性擬合通過原點,斜率約為0.3,更接近Piriz的理論預(yù)測值.

隨著擾動kξ0增大,尖釘?shù)脑鲩L速度逐步增加,最大振幅與增長速度之間表現(xiàn)出非線性關(guān)系,尖釘形態(tài)從擾動增長被抑制時的穩(wěn)定態(tài)向持續(xù)增長的不穩(wěn)定態(tài)過渡,直至尖釘頭部收縮變窄形成射流,最終脫離樣品主體形成噴射物.流體RM不穩(wěn)定性分析研究表明[1,2,10,22],在線性近似假設(shè)下尖釘?shù)淖畲笤鲩L速度可近似表示為初始擾動kξ0、沖擊波經(jīng)過后的界面起跳速度Δv和界面兩側(cè)Atwood數(shù)的乘積,即

其中Atwood數(shù)A=(ρ2?ρ1)/(ρ2+ρ1)表示界面兩側(cè)物質(zhì)的密度差,對于固體/真空界面,A=?1.沖擊波在彈塑性材料中傳播時,界面起跳速度Δv可由沖擊波強度和材料的Hugoniot沖擊線確定[11,17,23?25],若介質(zhì)的初始壓力和速度為零,壓力為P的沖擊波傳播時波后粒子速度與壓力的關(guān)系為

其中η表示沖擊波前后介質(zhì)的壓縮率,η=1?ρ0/ρ,與材料的聲阻抗和沖擊波速度有關(guān);vp為波后粒子速度.對于介質(zhì)/真空界面,自由面速度與波后粒子速度存在倍數(shù)關(guān)系,Δv≈2vp,將(9)和(10)式代入(8)式中得

因此尖釘振幅的增長因子為

(12)式表明,擾動增長被抑制時尖釘?shù)淖畲笳穹鲩L因子與參數(shù)組合kξ0/Y成線性關(guān)系,并且沖擊波壓力越高,尖釘增長越顯著.

3.4 線性近似的適用范圍

分析(12)式中制約尖釘增長因子的各因素:系數(shù)C是線性關(guān)系的斜率;壓縮率η與沖擊波壓力P共同表征了沖擊波的壓縮效應(yīng),為擾動增長提供啟動能量;kξ0描述了擾動界面的初始形態(tài),為初始擾動的振幅與波長之比;屈服強度Y是材料的物性參數(shù),擾動增長受其抑制.因此后三個因素之間是相互獨立的,本文分析系數(shù)C與其他因素之間的相關(guān)性,討論線性近似的適用范圍.

考慮三種材料,分別是6061鋁、OFHC銅和304不銹鋼,狀態(tài)方程為Mie-Grüneisen狀態(tài)方程,本構(gòu)模型為彈性理想塑性模型,參數(shù)見表1,炸藥為PBX9501,相關(guān)模型參數(shù)見表2.每種材料的屈服強度取三組值以模擬不同程度的硬化效應(yīng),初始擾動kξ0的變化范圍0.15—0.35.圖7為不同材料、強度和kξ0時系數(shù)C隨參數(shù)組合kξ0/Y的變化情況,圖8為kξ0=0.223和0.279時鋁材料擾動面上尖釘、氣泡和無擾動區(qū)域的速度曲線,分別對應(yīng)擾動增長被抑制和未被抑制的情況.

圖7 (網(wǎng)刊彩色)無量綱的尖釘最大振幅與增長速度之比隨參數(shù)組合kξ0/Y 的變化情況Fig.7.(color online)Plots of the ratios between the non-dimensional maximum amplitude and growth velocity of the spikes verse the combination of parameters kξ0/Y.

由圖7可知,系數(shù)分布存在三個階段,對應(yīng)于界面擾動的穩(wěn)定、穩(wěn)定與不穩(wěn)定之間的過渡和不穩(wěn)定階段.

1)參數(shù)組合kξ0/Y≤0.8 GPa?1,系數(shù)C基本為常數(shù),擾動處于被抑制的穩(wěn)定態(tài),尖釘?shù)淖畲笳穹c增長速度之間為線性關(guān)系,增長因子由(12)式描述,并且圖7中不同材料的系數(shù)一致,說明在穩(wěn)定段系數(shù)C基本保持獨立性,不受其它因素的影響.

2)參數(shù)組合0.8 GPa?1≤kξ0/Y≤1.0 GPa?1,系數(shù)C不再為常數(shù),擾動增長情況由穩(wěn)定態(tài)向不穩(wěn)定態(tài)過渡.此階段擾動可能會形成穩(wěn)定的尖釘形態(tài),也可能進一步發(fā)展成射流,取決于兩個特征時間的競爭,即沖擊波首次加載后在樣品內(nèi)部波系多次反射形成的二次沖擊波加載界面所需時間twave,與沖擊波首次加載后尖釘振幅從開始增長至被抑制所需的穩(wěn)定時間tstable.若twave＞tstable,表明在二次沖擊波到達前,強度的致穩(wěn)效應(yīng)已抑制了擾動增長,尖釘、氣泡和未擾動區(qū)域的速度趨于一致,如圖8所示,從被抑制到二次沖擊波到達,這段時間內(nèi)尖釘、氣泡和無擾動區(qū)域之間的相對速度保持為零,意味著尖釘、氣泡的振幅不再變化,界面形狀保持不變;若twave≤tstable,表明在二次沖擊波到達自由面時,強度的致穩(wěn)效應(yīng)未能將尖釘充分減速,尖釘與其他位置相比存在速度差,沖擊波二次加載將破壞此時的尖釘和氣泡形態(tài),界面形狀進一步演化.在過渡階段,(8)式的線性近似關(guān)系與(12)式的尖釘增長規(guī)律依然可能成立,取決于擾動增長在二次沖擊波達到前是否能被抑制.

圖8 (網(wǎng)刊彩色)擾動增長穩(wěn)定與不穩(wěn)定時尖釘、氣泡和未擾動區(qū)域的速度曲線Fig.8.(color online)Velocity pro files of the spikes,bubbles and the unperturbated areas in the stable and unstable growth of the perturbations.

3)參數(shù)組合kξ0/Y≥1.0 GPa?1,擾動增長顯著,尖釘頭部迅速演化為射流,界面失穩(wěn).

綜上所述,在擾動增長的穩(wěn)定態(tài),系數(shù)C保持獨立性,線性關(guān)系(8)式適用,描述尖釘增長規(guī)律的線性近似(12)式具備一定的材料普適性;在穩(wěn)定態(tài)向不穩(wěn)定態(tài)過渡階段以及不穩(wěn)定態(tài),系數(shù)C受其他因素影響失去獨立性,(8)式和(12)式不具備材料普適性;圖7中系數(shù)C保持為常數(shù)的平臺區(qū)域揭示了界面擾動從增長被抑制的線性關(guān)系區(qū)向形成射流的非線性區(qū)過渡時存在穩(wěn)定性閾值,表明強度介質(zhì)的RM擾動增長中也存在類似Rayleigh-Taylor擾動增長的穩(wěn)定性閾值[26],初步研究結(jié)果表明閾值對應(yīng)的物理狀態(tài)kξ0/Y取值約0.8 GPa?1.

4 結(jié) 論

本文研究了爆轟加載下不同材料飛片RM流動的擾動增長規(guī)律,分析結(jié)果表明擾動增長被抑制時,尖釘振幅增長因子與初始擾動的kξ0成正比,與材料的強度Y成反比,理論關(guān)系式符合相關(guān)實驗結(jié)果,與Piriz利用牛頓第二定律近似推導(dǎo)所得的結(jié)果基本一致,該理論關(guān)系式也為測量材料動態(tài)屈服強度提供了新的實驗技術(shù)途徑.對擾動增長規(guī)律影響因素的相關(guān)性分析表明,不同金屬材料的擾動增長均存在穩(wěn)定區(qū)域,穩(wěn)定區(qū)內(nèi)描述尖釘振幅增長因子與加載條件、初始擾動和材料性質(zhì)之間關(guān)系的(12)式具備一定的材料普適性,穩(wěn)定區(qū)范圍kξ0/Y＜0.8 GPa?1,其中閾值0.8 GPa?1的物理意義有待進一步分析研究.

[1]Richtmyer R D 1960Commun.Pure Appl.Math.13297

[2]Meshkov E E 1969Sovit.Fluid Dyn.4151

[3]Mikaelian K O 2013Phys.Rev.E87031003

[4]Brouillette M 2002Annu.Rev.Fluid Mech.34445

[5]Piriz A R,Lopez Cela J J,Tahir N A,Hoff mann D H H 2006Phys.Rev.E74037301

[6]Piriz A R,Lopez Cela J J,Cortazar O D,Tahir N A,Hoff mann D H H 2005Phys.Rev.E72056313

[7]Piriz A R,Lopez Cela J J,Tahir N A,Hoff mann D H H 2008Phys.Rev.E78056401

[8]Piriz A R,Lopez Cela J J,Tahir N A 2009Phys.Rev.E80046305

[9]Lopez Ortega A,Lombardini M,Pullin D I,Meiron D I 2014Phys.Rev.E89033018

[10]Remington B A,Rudd E R,Wark J S 2015Phys.Plasmas22090501

[11]Buttler W T,Oro D M,Preston D L,Mikaelian K O,Cherne F J,Hixson R S,Mariam F G,Morris C,Stone J B,Terrones G,Tupa D 2012J.Fluid Mech.70360

[12]Buttler W T,Oro D M,Olsen R T,Cheren F J,Hammerberg J E,Hixson R S,Monfared S K,Pack C L,Rigg P A,Stone J B,Terrones G 2014J.Appl.Phys.116103519

[13]Dimonte G,Terrones G,Cheren F J,Germann T C,Dunpont V,Kadau K,Buttler W T,Oro D M,Morris C,Preston D L 2011Phys.Rev.Lett.107264502

[14]Jensen B J,Cheren F J,Prime M B,Fezzaa K,Iverson A J,Carlson C A,Yeager J D,Ramos K J,Hooks D E,Cooley J C,Dimonte G 2015J.Appl.Phys.118195903

[15]Sun Z F,Xu H,Li Q Z,Zhang C Y 2010Chin.J.High Pressure Phys.2455(in Chinese)[孫占峰,徐輝,李慶忠,張崇玉2010高壓物理學(xué)報2455]

[16]Robinson A C,Swegle J W,1989J.Appl.Phys.662838

[17]Zhu J S,Hu X M,Wang P,Chen J,Xu A G 2010Adv.Mech.40400(in Chinese)[朱建士,胡曉棉,王裴,陳軍,許愛國2010力學(xué)進展40400]

[18]Vogler T J,Chhabildas L C 2006Int.J.Impact Engng.33812

[19]Barton N R,Bernier J V,Becker R,Arsenlis A,Cavallo R,Marian J,Rhee M,Park H S,Remington B A,Olson R T 2011J.Appl.Phys.109073501

[20]Smith R F,Eggert J H,Rudd R E,Swift D C,Blome C A,Collins G W 2011J.Appl.Phys.110123515

[21]Park H S,Rudd R E,Cavallo R M,Barton N R,Arsenlis A,Belof J L,Blobaum K J M,El-dasher B S,Florando J N,Huntington C M,Maddox B R,May M J,Plechaty C,Prisbrey S T,Remington B A,Wallace R J,Wehrenberg C E,Wilson M J,Comley A J,Giraldex E,Nikroo A,Farrell M,Randall G,Gray III G T 2015Phys.Rev.Lett.114065502

[22]Wouchuk J G 2001Phys.Rev.E63056303

[23]Pan H,Wu Z H,Hu X M,Yang K 2013Chin.J.High Pressure Phys.27778(in Chinese)[潘昊,吳子輝,胡曉棉,楊堃2013高壓物理學(xué)報27778]

[24]Pan H,Hu X M,Wu Z H,Dai C D,Wu Q 2012Acta Phys.Sin.61206401(in Chinese)[潘昊,胡曉棉,吳子輝,戴誠達,吳強2012物理學(xué)報61206401]

[25]Yu Y Y,Tan H,Hu J B,Dai C D,Chen D N,Wang H R 2008Acta Phys.Sin.572352(in Chinese)[俞宇穎,譚華,胡建波,戴誠達,陳大年,王煥然 2008物理學(xué)報 572352]

[26]Colvin J D,Legrand M,Remington B A,Schurtz G,Weber S V 2003J.Appl.Phys.935287