基于分子動(dòng)力學(xué)模擬的單晶硅沖擊壓縮相變研究*

劉夢(mèng)婷，李旺輝，奉蘭西，張曉晴，姚小虎

（華南理工大學(xué)土木與交通學(xué)院, 廣東 廣州 510641）

晶體硅具有純度高、晶體取向良好等特點(diǎn)，在電子等領(lǐng)域應(yīng)用廣泛，同時(shí)由于其復(fù)雜的結(jié)構(gòu)相變機(jī)制，在相圖研究中也被廣泛關(guān)注。

在準(zhǔn)靜態(tài)荷載下，單晶硅表現(xiàn)出顯著的多態(tài)性，目前已觀察到13 種相結(jié)構(gòu)。Mujica 等指出，單晶硅在10～11 GPa 的載荷下發(fā)生了從立方金剛石（cubic diamond, CD）結(jié)構(gòu)到β-白錫（β-tin）結(jié)構(gòu)的相變，并通過(guò)第一性原理計(jì)算出存在中間暫穩(wěn)態(tài) Imma 相。而增加壓力時(shí)，晶體硅繼續(xù)發(fā)生向簡(jiǎn)單六角形（simple hexagonal，SH）結(jié)構(gòu)的轉(zhuǎn)變；再進(jìn)一步增加壓力，變?yōu)镃mca 空間群結(jié)構(gòu)，而在41 GPa 的高壓下，晶體硅又從Cmca 結(jié)構(gòu)轉(zhuǎn)變?yōu)榱矫芏逊e結(jié)構(gòu)（hexagonal close pack，HCP）。在動(dòng)態(tài)荷載實(shí)驗(yàn)中，Gilev 等研究了在沖擊壓縮下單晶硅的金屬化現(xiàn)象，發(fā)現(xiàn)在沖擊壓力達(dá)到23 GPa 時(shí)，晶體結(jié)構(gòu)的缺陷導(dǎo)致電導(dǎo)率發(fā)生了偏離。Loveridge-Smith 等用X 射線衍射技術(shù)測(cè)量了沖擊加載下單晶硅的正交平面晶格參數(shù)，發(fā)現(xiàn)沿壓縮方向單晶硅的晶格間距減少了近11%，而與壓縮方向正交平面的晶格間距沒(méi)有顯著變化。Turneaure 等研究了單晶硅[100]、 [111]晶向在15.9、21.7 GPa 沖擊壓力下的動(dòng)力學(xué)響應(yīng)，發(fā)現(xiàn)其彈性波和非彈性波晶向依賴性很大，相變波晶向依賴性很小。Zhao 等對(duì)單晶硅[001]晶向進(jìn)行了激光沖擊壓縮和回收實(shí)驗(yàn)，觀察到兩個(gè)塑性響應(yīng)：靠近沖擊表面的大塊非晶層和沿{111}平面產(chǎn)生的滑移帶，該結(jié)果揭示了單晶硅在一維應(yīng)變條件下的塑性響應(yīng)。Smith 等對(duì)單晶硅[001]、[110]和[111]晶向進(jìn)行了沖擊實(shí)驗(yàn)，揭示了單晶硅材料沖擊響應(yīng)的各向異性，發(fā)現(xiàn)了其在沖擊壓力高于13 GPa 時(shí)轉(zhuǎn)變?yōu)楹?jiǎn)單六角形結(jié)構(gòu)。Liu 等研究了單晶硅沿[100]晶向沖擊壓縮的塑性變形，發(fā)現(xiàn)其塑性變形能力隨著溫度增加而增加。Kishimura 等對(duì)在38 GPa 沖擊壓力下的單晶硅進(jìn)行回收并進(jìn)行了X 射線分析和拉曼光譜分析。XRD 分析顯示，單晶硅在11 GPa 沖擊壓力下有少量亞穩(wěn)態(tài)相產(chǎn)生，在38 GPa 沖擊壓力下有CuSi 這一物質(zhì)生成，而在其余加載壓力下并沒(méi)有明顯的高壓相和亞穩(wěn)相存在。拉曼光譜的分析顯示，其峰值的整體偏移是由晶體尺寸減少引起的，該結(jié)果說(shuō)明了單晶硅沖擊壓縮下的晶?；^(guò)程。Turneaure 等用實(shí)時(shí)XRD 檢測(cè)了單晶硅在沖擊壓力低于54 GPa 下的晶體結(jié)構(gòu)，發(fā)現(xiàn)單晶硅在沖擊壓力為31～33 GPa 時(shí)產(chǎn)生沖擊熔化，并從熔融邊界再結(jié)晶到密堆積六角形結(jié)構(gòu)。由于硅和鍺屬于同族元素，具有相似的物理力學(xué)特性，因此Renganathan 等用原位XRD 檢測(cè)了鍺在[100]晶向的相變，發(fā)現(xiàn)鍺在沖擊壓力為15.7 GPa 時(shí)轉(zhuǎn)變?yōu)榘族a結(jié)構(gòu)，在沖擊壓力為31.5 GPa 時(shí)轉(zhuǎn)變?yōu)槿廴趹B(tài)；應(yīng)力卸載后，又轉(zhuǎn)變?yōu)榱⒎浇饎偸Y(jié)構(gòu)，該結(jié)果表明鍺從白錫結(jié)構(gòu)到立方金剛石結(jié)構(gòu)的相變是可逆的。McBride 等用實(shí)時(shí)X 射線檢測(cè)技術(shù)揭示了多晶硅的沖擊相變行為，觀察到不同壓力下的多種固-固相變以及固-液相變。Paul 等則利用密度泛函理論、進(jìn)化算法和晶格動(dòng)力學(xué)等多種方法，構(gòu)建了4 TPa 和6 000 K 下硅的相圖，為極高壓相圖的研究提供了參考。

盡管以上實(shí)驗(yàn)與計(jì)算研究結(jié)果豐富，但實(shí)驗(yàn)之間也存在較顯著的差異。一方面，實(shí)驗(yàn)直接獲取的物理量過(guò)少，難以直接表征結(jié)果；另一方面，各個(gè)實(shí)驗(yàn)設(shè)備和實(shí)驗(yàn)過(guò)程也存在差異，變量如應(yīng)變率、加載時(shí)間等不一致，而這些因素都會(huì)對(duì)單晶硅的相變有顯著影響。同時(shí)，單晶硅作為一種脆性材料，不僅加工困難，在沖擊壓縮時(shí)樣品也很難回收，而且單晶硅的相圖比較復(fù)雜，存在穩(wěn)態(tài)相和非穩(wěn)態(tài)中間相等多個(gè)相，因此它在沖擊壓縮下的物理和力學(xué)特性尚無(wú)完全統(tǒng)一的結(jié)論。此外，單晶硅沖擊實(shí)驗(yàn)中加載時(shí)間較短，很難實(shí)時(shí)測(cè)量其相變過(guò)程中的晶格結(jié)構(gòu)。所以，盡管檢測(cè)微觀結(jié)構(gòu)變化的原位實(shí)時(shí)診斷技術(shù)取得了很大的進(jìn)展，但對(duì)于可逆相變往往還是求助于把宏觀力學(xué)量的實(shí)驗(yàn)測(cè)試結(jié)果與靜高壓數(shù)據(jù)進(jìn)行比較，以確認(rèn)新相的結(jié)構(gòu)、相變的類型及性質(zhì)，這也會(huì)導(dǎo)致研究結(jié)果的差異。

由于分子動(dòng)力學(xué)具有超高的時(shí)間和空間分辨率，在微納米尺度研究材料的物理和力學(xué)特性時(shí)具有獨(dú)特的優(yōu)勢(shì)。而在當(dāng)前的超算能力下，分子動(dòng)力學(xué)模擬規(guī)?？梢赃_(dá)到千萬(wàn)甚至數(shù)十億原子，并能同時(shí)獲取各個(gè)原子運(yùn)動(dòng)的軌跡，進(jìn)而從原子尺度到微納米尺度揭示材料的變形與力學(xué)參量的關(guān)系，使得其能夠?qū)碚撚?jì)算和實(shí)驗(yàn)研究提供重要的補(bǔ)充。在分子動(dòng)力學(xué)模擬單晶硅方面也有諸多發(fā)現(xiàn)。Swift 等用第一性原理計(jì)算得到單晶硅的物態(tài)方程，這個(gè)物態(tài)方程可以反應(yīng)單晶硅從金剛石結(jié)構(gòu)到體心四方結(jié)構(gòu)的轉(zhuǎn)變。Demkowicz 等用SW 勢(shì)函數(shù)模擬非晶硅在施加恒定壓力下的變形情況，發(fā)現(xiàn)壓縮過(guò)程中存在類晶體和類液體部分，后者更容易發(fā)生塑性流動(dòng)。Kumagai 等改進(jìn)了Tersoff 勢(shì)函數(shù)，使之既能模擬彈性常數(shù)又能準(zhǔn)確估測(cè)熔點(diǎn)。Higginbotham等通過(guò)分子動(dòng)力學(xué)模擬，表明壓力誘導(dǎo)相變是彈性波異常的原因。Mogni 等則基于Tersoff-like 勢(shì)函數(shù)模擬[001]晶向單晶硅在激光沖擊下的相變，發(fā)現(xiàn)了新的Imma 相。Zhao 等研究了基于分子動(dòng)力學(xué)模擬的沖擊壓縮下單晶硅的力學(xué)響應(yīng)，發(fā)現(xiàn)硅在沖擊壓力低于10 GPa 時(shí)為彈性響應(yīng)，在稍高壓力下，由于位錯(cuò)和缺陷沿著特定的晶體學(xué)方向擴(kuò)展，會(huì)進(jìn)一步產(chǎn)生新的缺陷，此外在更高壓力下會(huì)產(chǎn)生交叉的非晶帶。

盡管基于分子動(dòng)力學(xué)的單晶硅變形與相變的模擬研究成果豐富，但仍有若干關(guān)鍵問(wèn)題尚未明確，包括單晶硅的動(dòng)態(tài)壓縮各向異性，以及單晶硅復(fù)雜相變機(jī)制的物理描述和不同因素對(duì)其相變的影響規(guī)律和機(jī)制等。深入研究這些關(guān)鍵問(wèn)題有利于更加全面地認(rèn)清單晶硅在動(dòng)態(tài)壓縮下的變形和相變的物理和力學(xué)特性。本文中將基于大規(guī)模分子動(dòng)力學(xué)方法模擬研究單晶硅的沖擊壓縮響應(yīng)，從微納米尺度認(rèn)識(shí)晶體硅的沖擊壓縮特性，揭示其在動(dòng)加載下的各向異性和結(jié)構(gòu)相變等特征。

1 研究方法

1.1 模型的構(gòu)建與沖擊壓縮模擬方法

基于廣泛使用的LAMMPS平臺(tái)，采用分子動(dòng)力學(xué)方法模擬研究單晶硅的沖擊壓縮響應(yīng)，主要對(duì)其非彈性變形和相變行為展開(kāi)分析。晶體硅通常情況下具有立方金剛石結(jié)構(gòu)，在300 K 溫度下晶格常數(shù)為0.543 nm。模型的建立是由單個(gè)晶胞沿、、晶向分別復(fù)制30、30 和390 個(gè)晶胞得到的超胞結(jié)構(gòu)。因此，模型尺寸約為16.3 nm×16.3 nm×212 nm，見(jiàn)表1。

表1 單晶硅計(jì)算模型詳細(xì)參數(shù)Table 1 Parameters of single crystal Si sample for MD simulation

晶體硅中Si-Si 原子間的相互作用采用改進(jìn)的半經(jīng)驗(yàn)Tersoff 勢(shì)函數(shù)描述，并使用Erhart 等優(yōu)化的參數(shù)。與改進(jìn)前的勢(shì)函數(shù)相比，該勢(shì)函數(shù)可以成功模擬晶體硅的多種性質(zhì)，包括彈性常數(shù)和一般熱力學(xué)性質(zhì)，以及在靜水壓縮和沖擊壓縮條件下立方金剛石相到β-白錫相的結(jié)構(gòu)特征。

在施加沖擊壓縮前，先對(duì)構(gòu)建的初始模型進(jìn)行弛豫，弛豫過(guò)程中保持三維方向均為周期性邊界條件，使整個(gè)晶體模型在NPT 系綜下弛豫了20 ps。從能量和殘余應(yīng)力的曲線上確認(rèn)了材料模型為能量趨于穩(wěn)定且殘余應(yīng)力在誤差允許范圍內(nèi)趨于零的狀態(tài)，此時(shí)獲得了合理穩(wěn)定的初始結(jié)構(gòu)。在沖擊壓縮加載階段，為了消除邊界對(duì)應(yīng)力波的影響，在沿軸施加自由邊界條件的同時(shí)，沿、方向施加周期性邊界條件。模擬平板沖擊采用的是活塞法，活塞對(duì)自由表面產(chǎn)生沖擊波。單晶硅沖擊加載的示意圖如圖1所示。沖擊壓縮模擬過(guò)程中系統(tǒng)處于NVE 系綜，并且不施加溫度控制，時(shí)間步長(zhǎng)為1 fs。這樣的加載方式和效果，與典型的平面沖擊實(shí)驗(yàn)相似，均為一維應(yīng)變問(wèn)題。為揭示單晶硅沖擊響應(yīng)的各向異性，尤其是沖擊結(jié)構(gòu)的演化，分別對(duì)[001]、[110]和[111]晶向進(jìn)行系列沖擊壓縮研究，沖擊粒子速度從0.3 km/s逐漸增加至3.2 km/s。

圖1 分子動(dòng)力學(xué)模擬中單晶硅的沖擊壓縮示意圖Fig. 1 Schematic of shock compression of single crystal Si in molecular dynamic simulations

1.2 主要物理量的計(jì)算和可視化分析方法

為了深入分析單晶硅的沖擊壓縮響應(yīng)，利用OVITO進(jìn)行可視化分析，并采用配位數(shù)、中心對(duì)稱參數(shù)、鍵角分析等多種分析技術(shù)，對(duì)晶體硅的變形和結(jié)構(gòu)變化進(jìn)行識(shí)別和分析。

2 結(jié)果與討論

2.1 沖擊響應(yīng)與各向異性

分子動(dòng)力學(xué)模擬計(jì)算結(jié)果與不同文獻(xiàn)中實(shí)驗(yàn)所得的晶體硅沖擊Hugoniot 應(yīng)力和體積變化的關(guān)系如圖2 所示。從圖2 中可以看出，在體積壓縮較小時(shí)，所有結(jié)果高度吻合；壓縮較大時(shí)，實(shí)驗(yàn)的應(yīng)力-體積關(guān)系呈現(xiàn)兩大分歧：一是以Goto 等和Gust 等的結(jié)果為代表，這兩個(gè)獨(dú)立研究團(tuán)隊(duì)的結(jié)果比較接近，二是以Turmeaure 等和McBride 等的實(shí)驗(yàn)結(jié)果為代表。因此，對(duì)晶體硅不同沖擊實(shí)驗(yàn)研究的差異值得關(guān)注，這可能與實(shí)驗(yàn)手段、技術(shù)和材料樣品差異有關(guān)。一方面，Gust 等早期采用爆炸驅(qū)動(dòng)的沖擊波壓縮實(shí)驗(yàn)，而Goto 等和Turmeaure 等采用輕氣炮驅(qū)動(dòng)的平板沖擊裝置對(duì)不同晶向的單晶硅進(jìn)行沖擊實(shí)驗(yàn)，McBride 等則采用激光驅(qū)動(dòng)的應(yīng)力脈沖加載對(duì)多晶硅進(jìn)行沖擊壓縮實(shí)驗(yàn)；另一方面，盡管所有實(shí)驗(yàn)都是運(yùn)用Rankine-Hugoniot 關(guān)系推算試樣的內(nèi)部壓力，但壓力推算中所需的粒子速度或相變波速度測(cè)定的差異顯著，例如Gust 等提及其結(jié)果與相關(guān)歷史數(shù)據(jù)的對(duì)比時(shí)，對(duì)應(yīng)于HEL 和相變點(diǎn)的臨界體積（或密度）是很吻合的，但壓力數(shù)值卻存在差異，他們將此歸因于粒子速度的測(cè)定差異。值得注意的是，兩大存在分歧的壓力-體積關(guān)系均有相互獨(dú)立的研究團(tuán)隊(duì)的結(jié)果支撐，可能預(yù)示著深層次的原因在于速度測(cè)定的技術(shù)差異與硅材料響應(yīng)有關(guān)。Goto 等的輕氣炮驅(qū)動(dòng)的平面沖擊實(shí)驗(yàn)中采用了與Gust 等相似的測(cè)量技術(shù)（斜鏡與條紋相機(jī)）對(duì)硅試件的自由面測(cè)速，其Hugoniot 壓力推算值也與Gust 等非常吻合；然而，Turneaure 等的輕氣炮驅(qū)動(dòng)平面沖擊實(shí)驗(yàn)與McBride 等的激光沖擊實(shí)驗(yàn)則采用了VISAR 測(cè)試系統(tǒng)對(duì)硅/LiF 窗口或者硅自由面速度測(cè)定，其結(jié)果相接近。此外，Smith 等指出硅自由表面端存在的應(yīng)力松弛和初始的速度濺射下的自由表面測(cè)量結(jié)果可能會(huì)影響應(yīng)力和體積的關(guān)系測(cè)定。因此，這些結(jié)果在壓力數(shù)值上的差異深刻反映了實(shí)驗(yàn)技術(shù)差異和數(shù)據(jù)處理方法對(duì)實(shí)驗(yàn)結(jié)果的重要影響。

圖2 不同實(shí)驗(yàn)[13,20,32-33]與模擬中晶體硅的沖擊壓縮Hugoniot 應(yīng)力和體積變化的關(guān)系Fig. 2 Shock Hugoniot stress as a function of volume change in various experiments[13,20,32-33] and simulations of Si crystals

與實(shí)驗(yàn)結(jié)果相比，本文計(jì)算結(jié)果與Goto 等和Gust 等的結(jié)果更接近，但也在一定程度上高于實(shí)驗(yàn)值。盡管如此，基于該原子間勢(shì)函數(shù)模擬計(jì)算的應(yīng)力-體積整體變化趨勢(shì)與實(shí)驗(yàn)結(jié)果相似。例如，隨著體積壓縮的增加，可以明顯看出在體積比低于0.8 時(shí)單晶硅已經(jīng)處于非彈性階段；而進(jìn)一步壓縮至體積比為約0.67 時(shí)，材料壓力急劇增加顯示出典型的高壓相變特征。因此，盡管該原子間勢(shì)函數(shù)對(duì)單晶硅的Hugoniot 壓力存在高估的不足，但其可以成功反映晶體硅在沖擊壓縮下的若干重要變形響應(yīng)特征，對(duì)本文研究動(dòng)高壓下硅的變形和相變?nèi)杂袃r(jià)值。

在不同沖擊壓縮下，單晶硅的剪切應(yīng)力隨粒子速度變化如圖3 所示。在[001]、[110]和[111]晶向中，其剪切應(yīng)力在一定范圍內(nèi)均隨著粒子速度的增加而逐漸增加，達(dá)到一定臨界粒子速度后開(kāi)始急劇下降。但是，該臨界粒子速度在不同晶向差異顯著，在[001]晶向最小，為1.4 km/s，在[110]晶向最大，為1.8 km/s，而在[111] 晶向時(shí)介于兩者之間，為1.5 km/s。在這3 個(gè)晶向中，最大剪切應(yīng)力分別為5.0、7.8 和6.2 GPa。剪切應(yīng)力的大幅度下降意味著單晶硅在這些沖擊加載下發(fā)生了重要的變形（如塑性或相變等），使得剪切應(yīng)力釋放。可以觀察到 [001]和[111]晶向中剪切應(yīng)力經(jīng)歷了由正值轉(zhuǎn)為負(fù)值的階段：分別在1.8、2.4 km/s 達(dá)到最大負(fù)剪切應(yīng)力，而后逐漸減小。負(fù)剪切應(yīng)力現(xiàn)象存在于單晶碳化硅和單晶金屬鐵中，學(xué)者們普遍認(rèn)為這是由于材料過(guò)度弛豫釋放的結(jié)果，并且具有晶向依賴性，多發(fā)生在[001]晶向。從剪切應(yīng)力的計(jì)算方法角度分析，負(fù)剪應(yīng)力意味著在沖擊加載中材料內(nèi)部的橫向應(yīng)力高于縱向應(yīng)力。最終，[001]、[110]和[111]晶向的剪切應(yīng)力分別在2.0、1.9、2.6 km/s 之后趨近零，因此材料失去剪切強(qiáng)度。此外，結(jié)合圖1 中應(yīng)力和體積壓縮的關(guān)系，單晶硅在[110]晶向的奇異表現(xiàn)還體現(xiàn)在彈塑性或相變轉(zhuǎn)變的臨界狀態(tài)，相比于[001]和[111]晶向的在約0.85 體積壓縮比時(shí)的“彈塑性”轉(zhuǎn)變，[110]晶向則表現(xiàn)出奇異的彈性特征。這種沖擊壓縮下奇異的彈性響應(yīng)現(xiàn)象在金剛石等材料的分子動(dòng)力學(xué)模擬研究也有過(guò)發(fā)現(xiàn)，這些材料呈現(xiàn)出非單調(diào)的剪應(yīng)力-應(yīng)變關(guān)系。由于剪切應(yīng)力是塑性等變形的驅(qū)動(dòng)力，因此，單晶硅沖擊壓縮下奇異的彈性響應(yīng)應(yīng)該與其剪切應(yīng)力的非單調(diào)性相關(guān)。

圖3 沖擊Hugoniot 狀態(tài)的剪切應(yīng)力(τ)與粒子速度(up) 曲線Fig. 3 Shock Hugoniot shear stress as a function of particle velocity

為了進(jìn)一步探究單晶硅沖擊加載下的材料特性及其各向異性，圖4 中給出了沿[001]、[110] 和[111]晶向加載時(shí)單晶硅的縱向應(yīng)力、剪切應(yīng)力和密度波剖面。根據(jù)單晶硅的沖擊Hugoniot 曲線，單晶硅3 個(gè)晶向的彈塑性或相變閾值是不同的，[110]是3 個(gè)晶向相變閾值最大的，閾值粒子速度為1.8 km/s，其次是[111]晶向，閾值粒子速度為1.5 km/s，[001]晶向最小，為1.4 km/s。因此，波剖面分析中選取了各晶向中具有代表性的3 個(gè)粒子速度狀態(tài)。如圖4(a)所示，在[001]晶向中，粒子速度分別為1.3、1.4、1.5 km/s的波剖面體現(xiàn)了單晶硅從單一彈性波到彈性-塑性/相變波的沖擊波結(jié)構(gòu)轉(zhuǎn)變過(guò)程。在縱向應(yīng)力波剖面中，當(dāng)粒子速度為1.3 km/s 時(shí)，波剖面平整，且波陣面十分狹窄陡峭。當(dāng)粒子速度增加至1.4 km/s 時(shí)，彈性前驅(qū)波幅值增加，但同時(shí)材料出現(xiàn)軟化現(xiàn)象，導(dǎo)致彈性波之后出現(xiàn)應(yīng)力幅值下降的現(xiàn)象，同時(shí)由于塑性或相變的激發(fā)，在靠近活塞端的波剖面區(qū)域呈現(xiàn)比軟化區(qū)域更高幅值的塑性或相變波。當(dāng)粒子速度為1.5 km/s 時(shí)，類似的現(xiàn)象依然存在，但由于沖擊強(qiáng)度增加，彈性前驅(qū)波幅值增加，這可能是由于單晶硅沖擊Hugoniot 彈性極限對(duì)加載速率敏感性很強(qiáng)導(dǎo)致的；同時(shí)，更高的沖擊強(qiáng)度使單晶硅的軟化越發(fā)明顯，使得彈性前驅(qū)波之后的波幅值下降更快；同樣，塑性或相變波部分幅值增加，對(duì)應(yīng)的波速度也在增加。從[001]晶向的剪切應(yīng)力波剖面可以進(jìn)一步揭示該現(xiàn)象的變化。剪切應(yīng)力對(duì)材料的響應(yīng)非常靈敏，單晶硅在沖擊壓縮下呈現(xiàn)復(fù)雜的材料響應(yīng)，包括彈性波及其軟化和塑性或相變波。剪切應(yīng)力隨著彈性波的出現(xiàn)急劇上升，然后隨著軟化作用而略有降低，并在塑性或相變波區(qū)域迅速下降，更強(qiáng)的塑性或相變波引起更大幅值的剪切應(yīng)力釋放。密度波剖面的變化對(duì)應(yīng)了包括剪切應(yīng)力波剖面在內(nèi)的各物理量波剖面情況，即密度隨著彈性波抵達(dá)而快速上升，同時(shí)由于軟化作用略有回落，并在相變波區(qū)域由于材料的進(jìn)一步壓縮而顯著上升。本文中研究獲得的單晶硅的波剖面特征在Higginbotham 等和Stubley等的研究工作也得到了體現(xiàn)。

圖4 單晶硅中分別沿[001]、[110]和[111]晶向的沖擊波剖面Fig. 4 Shock profiles in single crystal Si along [001], [110] and [111] crystal orientation

在[110]和[111]晶向中，復(fù)雜的多波結(jié)構(gòu)依舊存在，但軟化作用更加微弱。同時(shí)，與[001]晶向不同的是，在[110]和[111]晶向中波剖面彈性波幅值幾乎保持不變，顯示出與加載速率無(wú)關(guān)，塑性或相變波更加復(fù)雜。

基于波剖面的分析，采用可視化分析技術(shù)進(jìn)一步探究沖擊引起的單晶硅非彈性變形和相變等材料響應(yīng)。圖5 中分別給出了沿著[001]、[110]和[111]晶向沖擊壓縮時(shí)的可視化結(jié)果，粒子速度分別為1.4、1.8、1.5 km/s。結(jié)合上文分析可知，這些粒子速度沖擊壓縮下單晶硅3 個(gè)晶向中的剪切應(yīng)力均開(kāi)始產(chǎn)生急劇的下降，意味著單晶硅在這些加載區(qū)域發(fā)生了重要變形響應(yīng)。中心對(duì)稱參數(shù)和配位數(shù)是反映晶體結(jié)構(gòu)特征的重要指標(biāo)，可以體現(xiàn)單晶硅的發(fā)生與結(jié)構(gòu)相關(guān)的變化。圖5(a)中單晶硅根據(jù)中心對(duì)稱參數(shù)分析著色，而圖5(b)則根據(jù)原子最小配位數(shù)著色。在中心對(duì)稱參數(shù)分析中，3 個(gè)晶向的單晶硅均發(fā)生了顯著的中心對(duì)稱參數(shù)變化，但呈現(xiàn)出不同的帶狀紋理。未干擾區(qū)域具有比壓縮區(qū)域更高的中心對(duì)稱參數(shù)，而在帶狀條紋區(qū)域該參數(shù)則更低，并且在交叉區(qū)域該參數(shù)呈現(xiàn)出趨于零的變化。與此同時(shí)，在原子最小配位數(shù)分析中也可以發(fā)現(xiàn)3 個(gè)晶向的差異顯著。[001]單晶硅存在配位數(shù)為5 和6 的帶狀和交叉區(qū)域；而[110]單晶硅中在靠近加載端一段區(qū)域幾乎全部變成具有更大配位數(shù)的原子；在[111]單晶硅中則僅有極少數(shù)原子發(fā)生最小配位數(shù)的變化。這說(shuō)明在3 個(gè)晶向中單晶硅的晶體結(jié)構(gòu)均發(fā)生了結(jié)構(gòu)性的變化，但是呈現(xiàn)明顯的晶體各向異性。從不同的配位數(shù)和中心對(duì)稱參數(shù)看，單晶硅沖擊壓縮下發(fā)生了多種復(fù)雜結(jié)構(gòu)相變，這體現(xiàn)在某些結(jié)構(gòu)相變可引起中心參數(shù)的變化但不會(huì)引起配位數(shù)的變化。

圖5 單晶硅沿[001]、[110]、[111]晶向分別在1.4、1.8、1.5 km/s 沖擊壓縮下的可視化Fig. 5 Visualizations of single crystal Si in [001], [110] and [111] crystal orientations with shock particle velocities of 1.4, 1.8, 1.5 km/s

2.2 沿[001]晶向沖擊的變形與相變

為進(jìn)一步探究單晶硅的沖擊結(jié)構(gòu)相變等特征，重點(diǎn)分析[001]晶向單晶硅的情況。圖6 中給出了不同粒子速度下的單晶硅沖擊響應(yīng)可視化分析，通過(guò)晶格原子的偏轉(zhuǎn)角可以顯著地區(qū)分出沖擊壓縮下單晶硅變形和相變的局部變化，藍(lán)色區(qū)域原子的偏轉(zhuǎn)角非常小，黃色區(qū)域原子較初始結(jié)構(gòu)偏轉(zhuǎn)約54°，綠色區(qū)域原子較初始結(jié)構(gòu)偏轉(zhuǎn)45°。圖中[001]晶向的相變帶狀區(qū)域延伸方向與加載方向呈57°夾角，與偏轉(zhuǎn)角接近，這樣的少許差異與材料處于壓縮狀態(tài)相關(guān)。這表明相變波的相變帶遵循特定的晶體學(xué)取向，是沿著靠近{111}晶面向前傳播，其中{111}晶面與{001}晶面夾角為54°。隨著粒子速度的增加，相變帶狀區(qū)域?qū)挾扔新晕⒃黾拥内厔?shì)。當(dāng)粒子速度為2.0 km/s 時(shí)，靠近沖擊加載端的區(qū)域出現(xiàn)無(wú)序的狀態(tài)，這個(gè)無(wú)序區(qū)域隨著粒子速度的增加而擴(kuò)展。當(dāng)粒子速度為2.4 km/s 時(shí)，一個(gè)顯著的無(wú)序狀態(tài)與晶體狀態(tài)同時(shí)存在于相同加載方向位置?？紤]到強(qiáng)沖擊可能引起劇烈溫升和更高的壓力狀態(tài)，因此可能引起單晶硅的熔化，這意味著在此沖擊條件下單晶硅處于固-液共存的區(qū)域。值得注意的是，本文中所發(fā)現(xiàn)的沖擊波陣面后方固-固相變和固-液相變共存的現(xiàn)象也出現(xiàn)在其他兩個(gè)晶向，該發(fā)現(xiàn)與McBride 等的激光沖擊實(shí)驗(yàn)診斷結(jié)果一致，本文中沖擊壓縮下固-液共存甚至完全熔化所對(duì)應(yīng)的臨界體積壓縮比，也與McBride 等的實(shí)驗(yàn)高度吻合。如圖6 所示，當(dāng)粒子速度高達(dá)2.8 km/s 時(shí)熔化區(qū)域已經(jīng)擴(kuò)展很大，并在3.2 km/s 時(shí)波陣面后幾乎完全熔化。

圖6 不同沖擊粒子速度下[001]單晶硅的結(jié)構(gòu)演化Fig. 6 Structural evolutions of single crystal Si in [001]orientation at various shock particle velocities

基于圖6 的代表性局部區(qū)域，進(jìn)一步通過(guò)鍵角分析和晶向分布函數(shù)分析沖擊壓縮下單晶硅的結(jié)構(gòu)相變。圖7(a)為幾個(gè)代表性的局部區(qū)域的鍵角分析。區(qū)域 Ⅰ代表沖擊波陣面處于帶狀區(qū)域之前的部分，區(qū)域Ⅱ代表黃色帶狀區(qū)域，區(qū)域Ⅲ代表綠色的交叉區(qū)域，區(qū)域Ⅳ代表中間合圍區(qū)域，區(qū)域Ⅴ代表無(wú)定型區(qū)域。初始結(jié)構(gòu)中，鍵角呈現(xiàn)109°，區(qū)域Ⅰ主要鍵角峰值在105°，區(qū)域Ⅱ和Ⅲ則依次主要呈現(xiàn)96.5°和接近90°，而區(qū)域Ⅳ則與初始結(jié)構(gòu)非常相似。這說(shuō)明除了區(qū)域Ⅳ，其余區(qū)域均出現(xiàn)了與結(jié)構(gòu)變化相關(guān)的鍵角變化。圖7(b)的徑向分布函數(shù)進(jìn)一步揭示了可能的結(jié)構(gòu)變化。從圖7 可知，與初始結(jié)構(gòu)相比，相變區(qū)域Ⅱ～Ⅲ、Ⅴ變化明顯。其中區(qū)域Ⅱ～Ⅲ中的峰值位置和數(shù)量與初始結(jié)構(gòu)差異顯著，體現(xiàn)了晶體結(jié)構(gòu)的固-固轉(zhuǎn)變，而區(qū)域Ⅴ中峰值之后趨于平滑的變化則說(shuō)明硅材料向無(wú)定型狀態(tài)的變化。

圖7 對(duì)代表性局部區(qū)域的鍵角分析和徑向分布函數(shù)分析Fig. 7 Bond-angle and radial distribution functions analysis for representative local regions

2.3 沿[110]和[111]晶向沖擊的變形與相變

利用相似的后處理技術(shù)，對(duì)單晶硅沿[110]和[111]晶向沖擊壓縮的變形和相變進(jìn)行分析。圖8～9分別為單晶硅[110]晶向中不同粒子速度下的結(jié)構(gòu)演化和代表性局部的鍵角統(tǒng)計(jì)和徑向分布函數(shù)分析。當(dāng)粒子速度為1.8 km/s 時(shí)，在靠近左端一定范圍內(nèi)出現(xiàn)了兩種變化：垂直于加載方向呈帶狀條紋區(qū)域的晶格原子有很小的晶向偏轉(zhuǎn)，緊隨其后的是一定程度的原子無(wú)序狀態(tài)區(qū)域。帶狀區(qū)域中區(qū)域Ⅰ的鍵角分布顯示了兩個(gè)主峰值，分別對(duì)應(yīng)約90°和120°的鍵角，而區(qū)域Ⅱ的主峰值仍為109°。區(qū)域Ⅲ的徑向分布函數(shù)在第1 個(gè)峰值后趨于平緩且鍵角分布沒(méi)有獨(dú)立峰值，該區(qū)域?yàn)闊o(wú)序狀態(tài)，這可能是沖擊引起的熔化。對(duì)比圖4 的波剖面可知，區(qū)域Ⅰ和區(qū)域Ⅲ這兩種變化與波陣面后方的相變波對(duì)應(yīng)，而且是分離的兩列波，對(duì)應(yīng)著材料內(nèi)部不同程度的剪應(yīng)力釋放和密度上升。隨著沖擊粒子速度的增加，兩列波波速均在增加。并且當(dāng)粒子速度高于2.4 km/s 時(shí)，固-液相變區(qū)域逐漸占據(jù)主導(dǎo)。

圖8 不同沖擊粒子速度下[110]單晶硅的結(jié)構(gòu)演化Fig. 8 Structural evolutions of single crystal Si in [110]orientation at various shock particle velocities

圖9 對(duì)圖8 代表性局部區(qū)域的鍵角分析和徑向分布函數(shù)分析Fig. 9 Bond-angle and radial distribution functions analysis for representative local regions

圖10～11 為單晶硅[111]晶向中不同沖擊粒子速度下的結(jié)構(gòu)演化和代表性局部的鍵角統(tǒng)計(jì)與徑向分布函數(shù)分析。當(dāng)粒子速度為1.6 km/s 時(shí)，在波陣面后一定區(qū)域內(nèi)也出現(xiàn)了帶狀條紋以及靠近端部的無(wú)序狀態(tài)區(qū)域。與[001]和[110]晶向不同的是，[111]晶向中帶狀條紋出現(xiàn)后，隨著沖擊粒子速度的增加，條紋寬度逐漸減小，并在粒子速度為2.4 km/s 時(shí)逐漸消失。無(wú)序狀態(tài)區(qū)域隨著粒子速度的增加也不斷擴(kuò)展，并在3.2 km/s 時(shí)擴(kuò)展速度急劇增加。這一現(xiàn)象可能是沖擊加載下溫升越發(fā)顯著引起的軟化作用和壓力引起的結(jié)構(gòu)變化共同作用的結(jié)果。在相對(duì)低的沖擊粒子速度下軟化作用相對(duì)較弱，而壓力誘發(fā)的晶格原子偏轉(zhuǎn)更顯著，當(dāng)沖擊強(qiáng)度增加時(shí)，溫升會(huì)更顯著，軟化作用影響增加，晶格偏轉(zhuǎn)受到一定的削弱。類似于壓縮下的變形孿晶，高溫會(huì)抑制變形孿晶的形成。由于硅晶體材料的變形和相變對(duì)溫度、壓力、應(yīng)變率等因素比較敏感，因此其動(dòng)態(tài)物理和力學(xué)行為會(huì)受到這些競(jìng)爭(zhēng)性作用因素的影響。當(dāng)粒子速度高達(dá)3.2 km/s 時(shí)，固-液相變波波速急劇增大，這是軟化區(qū)域大面積激活轉(zhuǎn)化為熔化狀態(tài)導(dǎo)致的。

圖10 不同沖擊粒子速度下[111]單晶硅的結(jié)構(gòu)演化Fig. 10 Structural evolutions of single crystal Si in [111]orientation at various shock particle velocities

圖11 對(duì)代表性局部區(qū)域的鍵角分析和徑向分布函數(shù)分析Fig. 11 Bond-angle and radial distribution functions analysis for representative local regions

2.4 討論

從早期實(shí)驗(yàn)與模擬研究以及本文的模擬工作中可以看到晶體硅的復(fù)雜相變。在單軸加載條件下，沿著[001]、[110]、[111]晶向的加載過(guò)程都會(huì)發(fā)生相變，但是相變特征和機(jī)制具有各向異性。實(shí)驗(yàn)中已探明晶體硅在壓力逐漸增加作用下發(fā)生一些列相變，如cd 相→β-tin 相→Imma 相→sh 相等。其中cd 相是四配位結(jié)構(gòu)，sh 相為六配位結(jié)構(gòu)。β-tin 相是介于四配位和六配位之間的一種體心四角結(jié)構(gòu)，而Imma 相又是介于β-tin 相和sh 相之間的正交相。因此，本文模擬中涉及的固-固相變可能包括β-tin 相和Imma 相。Mogni 等基于相同的原子間勢(shì)函數(shù)的模擬分析指出，[001]晶向單晶硅的沖擊壓縮固-固相變結(jié)構(gòu)是Imma 相，而且他們認(rèn)為僅[001]晶向中帶狀區(qū)域?yàn)镮mma 相。通過(guò)本文中對(duì)[001]單晶硅沖擊壓縮下的交叉區(qū)域的偏轉(zhuǎn)角和配位數(shù)以及鍵角的分析，可以看出交叉區(qū)域的結(jié)構(gòu)應(yīng)該與Imma 相有所不同。該區(qū)域的配位數(shù)顯示可達(dá)5 甚至6，與sh 相具有相同配位，但本文中尚無(wú)足夠完整的結(jié)構(gòu)信息用于確定具體結(jié)構(gòu)。從圖5 中可以發(fā)現(xiàn)，單晶硅[001]、[110]和[111]晶向具有明顯各向異性，如沿[001]晶向加載時(shí)配位數(shù)發(fā)生變化，但在[110]和[111]晶向時(shí)配位數(shù)沒(méi)有變化，該現(xiàn)象是由于不同加載方式的原子堆積方式不同導(dǎo)致的。此外，不同晶向的變形與相變受到滑移系的影響，位錯(cuò)阻力小的滑移系更容易滑移。對(duì)于單晶硅這種金剛石結(jié)構(gòu)，密排面為(111)晶面，密排方向?yàn)?110>方向，所以位錯(cuò)更容易沿著滑移系{111}<110>產(chǎn)生。所以相較于[001]晶向，[110]晶向偏轉(zhuǎn)角非常小，發(fā)生結(jié)構(gòu)變化的平面幾乎平行于(110)平面。

晶體硅的沖擊壓縮研究揭示了不同的物理和力學(xué)行為，但不同實(shí)驗(yàn)之間的確存在較大差異。在已經(jīng)報(bào)道出的晶體硅變形中除了熔化還包括位錯(cuò)滑移、無(wú)定型化等。McBride 等指出非原位探測(cè)很難識(shí)別出單晶硅沖擊壓縮階段的熔化，一些中間相變也很難通過(guò)回收試樣的方法準(zhǔn)確觀察到，因?yàn)榫w硅在高壓卸載時(shí)可能產(chǎn)生其他相結(jié)構(gòu)。對(duì)分子動(dòng)力學(xué)模擬工作，盡管本文定性地描述了單晶硅的沖擊壓縮變形和相變，觀測(cè)到彈性、非彈性、固-固相變以及固-液相變和固-液共存等物理特征，對(duì)比了晶體的各向異性，但是也應(yīng)指出，本文所采用的勢(shì)函數(shù)仍然不能準(zhǔn)確包含單晶硅的全部變形和相變特征。因此，應(yīng)對(duì)不同原子間勢(shì)函數(shù)進(jìn)行系統(tǒng)評(píng)估和對(duì)相關(guān)勢(shì)函數(shù)進(jìn)行修改。

總之，本文已經(jīng)初步明確單晶硅沖擊壓縮下可依次呈現(xiàn)彈性、非彈性到固-固相變和固-液相變的變化趨勢(shì)及其各向異性特征，對(duì)進(jìn)一步理解單晶硅的沖擊壓縮變形和相變提供了一定的支持。針對(duì)動(dòng)態(tài)加載下單晶硅相變的復(fù)雜性，后續(xù)仍需對(duì)單晶硅在動(dòng)態(tài)條件下的變形和相變研究提供更深入全面的納觀尺度信息。

3 結(jié) 論

采用大規(guī)模分子動(dòng)力學(xué)方法，模擬研究了單晶硅的沖擊壓縮響應(yīng)，獲得了其沖擊壓力和剪切應(yīng)力的變化，分析了其沖擊波傳播特性各向異性，并重點(diǎn)研究了不同沖擊強(qiáng)度下單晶硅變形的非彈性轉(zhuǎn)變和相變行為，得到主要研究結(jié)論如下：

（1）單晶硅在沖擊壓縮下呈現(xiàn)彈性到相變的變化趨勢(shì)，其剪切應(yīng)力隨著沖擊粒子速度的增加先增加，達(dá)到臨界點(diǎn)后由于結(jié)構(gòu)的變化導(dǎo)致應(yīng)力急劇下降，其剪應(yīng)力和相應(yīng)的臨界點(diǎn)存在明顯的各向異性；

（2）基于偏轉(zhuǎn)角和鍵角以及徑向分布函數(shù)的分析發(fā)現(xiàn)，沿[001]晶向沖擊壓縮下可發(fā)生多種固-固相變以及固-液相變，并觀察到與McBride 等的實(shí)驗(yàn)高度一致的固-液共存的重要現(xiàn)象。

本文研究工作可為動(dòng)加載下晶體硅的相變研究提供新的納米尺度的結(jié)果支撐，為相關(guān)實(shí)驗(yàn)提供新的補(bǔ)充。