劉 軍，孫致遠(yuǎn)，張鳳國(guó)，王 裴,2

（1. 北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所，北京 100094；2. 北京大學(xué)應(yīng)用物理與技術(shù)研究中心，北京 100871）

激波在金屬自由面卸載后，金屬材料內(nèi)部可能會(huì)產(chǎn)生層裂（或微層裂）現(xiàn)象。根據(jù)激波加載強(qiáng)度、波形、演化時(shí)間以及金屬材料物性等的不同，金屬層裂區(qū)內(nèi)部會(huì)形成層裂片、孔洞等不同的拉伸稀疏狀態(tài)。當(dāng)金屬層裂區(qū)后方再次受到加載波或慣性聚心等作用時(shí)，金屬層裂區(qū)會(huì)被密實(shí)基體持續(xù)追趕并再次壓實(shí)為密實(shí)狀態(tài)，直至處于拉伸狀態(tài)的物質(zhì)全部被壓實(shí)后金屬層裂區(qū)消失。該過(guò)程即為本文中所述的“金屬層裂區(qū)的再壓實(shí)過(guò)程”。

在實(shí)驗(yàn)觀測(cè)上，金屬層裂區(qū)的再壓實(shí)過(guò)程體現(xiàn)為自由面測(cè)速曲線在首次沖擊起跳后經(jīng)過(guò)一定時(shí)間間隔的二次起跳行為。而且由于不同位置上二次加載波到達(dá)時(shí)刻與多種因素相關(guān)，實(shí)驗(yàn)上難以確定層裂區(qū)進(jìn)入再壓實(shí)的時(shí)刻，僅能確定層裂區(qū)再壓實(shí)的完成時(shí)刻（即速度二次起跳時(shí)刻）。目前實(shí)驗(yàn)在層裂區(qū)拉伸狀態(tài)及再壓實(shí)后物質(zhì)狀態(tài)的精確測(cè)量上仍存在較大難度。在理論研究中，金屬層裂過(guò)程一般認(rèn)為是材料內(nèi)部孔洞產(chǎn)生、演化發(fā)展及匯合的過(guò)程。當(dāng)孔洞演化到一定時(shí)間后，相互貫通形成層裂片。相關(guān)模擬研究表明，大的層裂片內(nèi)部也可能存在未貫通的小孔洞。由此來(lái)看，金屬層裂區(qū)內(nèi)部的拉伸狀態(tài)可能是孔洞、層裂片或兩者共存的狀態(tài)，其是隨時(shí)間變化的。在層裂區(qū)進(jìn)入再壓實(shí)過(guò)程的初始時(shí)刻難以確定情況下，很難對(duì)金屬層裂區(qū)的初始拉伸狀態(tài)進(jìn)行準(zhǔn)確預(yù)測(cè)。正是由于進(jìn)入再壓實(shí)階段的初始拉伸狀態(tài)的復(fù)雜性及再壓實(shí)后狀態(tài)的不確定性，宏觀模擬結(jié)果與少數(shù)實(shí)驗(yàn)結(jié)果能夠較好比對(duì)，尤其難以驗(yàn)證很多復(fù)雜加載下宏觀模擬該過(guò)程的可靠性，制約了金屬層裂區(qū)再壓實(shí)過(guò)程的深入研究。

在無(wú)法確定層裂區(qū)初始拉伸狀態(tài)及再壓實(shí)狀態(tài)的情況下，本文中首先將層裂區(qū)初始含有不同拉伸狀態(tài)的因素引入該問題，并通過(guò)直接模擬得到不同拉伸狀態(tài)下的再壓實(shí)狀態(tài)。然后，開展與直接模擬中工況完全對(duì)應(yīng)的宏觀模擬，從而保證直接模擬與宏觀模擬具有可比性。最后，通過(guò)與對(duì)應(yīng)工況下直接模擬結(jié)果對(duì)比，對(duì)宏觀精細(xì)模擬此類問題的正確性和適應(yīng)性進(jìn)行討論。

1 金屬斷裂再壓實(shí)的動(dòng)力學(xué)過(guò)程

基于宏觀實(shí)驗(yàn)診斷結(jié)果及微介觀分子動(dòng)力學(xué)模擬總結(jié)，本文中將強(qiáng)沖擊條件下金屬層裂再壓實(shí)的力學(xué)行為分為如下6 個(gè)典型過(guò)程（見圖1）：①材料在稀疏波作用下的拉伸稀疏過(guò)程；②材料內(nèi)部微孔洞增長(zhǎng)過(guò)程；③大孔洞貫通后，含真空間隙層裂片形成的過(guò)程；④二次加載作用下大層裂片的閉合過(guò)程；⑤層裂片閉合后形成的含孔洞物質(zhì)再壓實(shí)過(guò)程；⑥壓縮成密實(shí)材料后的力學(xué)行為演化過(guò)程。需要說(shuō)明的是，雖然本文中將材料層裂再壓實(shí)過(guò)程的時(shí)間序列簡(jiǎn)化為6 個(gè)過(guò)程，但真實(shí)情況下上述幾個(gè)過(guò)程間的轉(zhuǎn)換并不是序列化完成的，過(guò)程之間存在弛豫、交疊，例如：若無(wú)二次加載，則上述過(guò)程簡(jiǎn)化為①→② 或①→②→③；若二次加載時(shí)間間隔較短，則材料內(nèi)部孔洞匯合形成層裂片過(guò)程可能無(wú)法完成，斷裂再壓實(shí)過(guò)程變?yōu)棰佟凇荨蓿壳澳承┛锥丛鲩L(zhǎng)模型中考慮了孔洞縮小因素，則該過(guò)程亦可認(rèn)為是①→②）；若二次加載時(shí)間間隔較長(zhǎng)，層裂區(qū)發(fā)展較充分，在受到一個(gè)較強(qiáng)的二次加載后上述6 個(gè)過(guò)程全部發(fā)生；若金屬在熔化狀態(tài)下形成微層裂現(xiàn)象，拉伸區(qū)充分發(fā)展后主要以金屬小液滴形式存在，上述描述過(guò)程不適用。

圖1 金屬層斷裂再壓實(shí)力學(xué)過(guò)程示意圖Fig. 1 Schematic diagram of the mechanical process of metal spallation and recompression

典型的①→②過(guò)程在實(shí)驗(yàn)測(cè)試診斷中表現(xiàn)為自由面速度由最高點(diǎn)回落及回落速率變化等。宏觀模擬①→②過(guò)程的計(jì)算，可使用損傷模型計(jì)算孔隙度 α ?？紤]其對(duì)材料狀態(tài)方程、強(qiáng)度的影響，重新計(jì)算得到該過(guò)程材料壓力、溫度等信息。具體形式如下：

典型的①→②過(guò)程在實(shí)驗(yàn)測(cè)試診斷中表現(xiàn)為自由面速度由最高點(diǎn)回落到一定高度后再次起跳（層裂片信號(hào)）。宏觀模擬②→③過(guò)程的計(jì)算，則涉及損傷模型參數(shù)之外的另一個(gè)重要參數(shù)：截?cái)鄵p傷度。隨著材料損傷度的增長(zhǎng)，當(dāng)＞情況下，認(rèn)為網(wǎng)格內(nèi)部大孔洞貫通導(dǎo)致含真空間隙層裂片形成。該過(guò)程計(jì)算方式可采用如下形式：

上述過(guò)程在單次沖擊金屬層裂行為模擬研究中已取得較好效果，能夠?qū)崿F(xiàn)對(duì)單次沖擊層裂現(xiàn)象的模擬預(yù)測(cè)（如自由面速度、層裂片厚度等）。

典型的③→⑥過(guò)程在實(shí)驗(yàn)中表現(xiàn)為界面速度平臺(tái)及二次起跳過(guò)程。由于相關(guān)物理量的實(shí)驗(yàn)測(cè)試診斷存在困難，目前很少有針對(duì)③→⑥全過(guò)程的深入研究。這就使金屬層裂再壓實(shí)過(guò)程的宏觀模擬存在一些難以驗(yàn)證的問題。首先，目前針對(duì)該過(guò)程的壓力計(jì)算處理方法種類較多：直接使用金屬負(fù)壓段EOS 處理該過(guò)程；有些則在壓實(shí)前全程設(shè)置零壓；還有的借鑒使用多孔介質(zhì)的相應(yīng)物理建模方法。但對(duì)于在特定工況下哪種處理更合理并沒有較明確的認(rèn)識(shí)。而且，針對(duì)該過(guò)程網(wǎng)格溫度的計(jì)算處理方法也是沒有相關(guān)文獻(xiàn)進(jìn)行闡述。其次，在缺乏有效驗(yàn)證手段情況下，難以對(duì)該過(guò)程的密度變化、再壓實(shí)狀態(tài)等是否正確做出明確的判斷。下面通過(guò)對(duì)幾種典型拉伸狀態(tài)下層裂區(qū)再壓實(shí)的直接模擬，對(duì)上述兩點(diǎn)在宏觀模擬此類問題中的正確性及適應(yīng)性進(jìn)行討論。

2 金屬層裂再壓實(shí)過(guò)程的直接數(shù)值模擬

考慮到金屬層裂區(qū)初始進(jìn)入再壓實(shí)過(guò)程（③→④過(guò)程）真實(shí)拉伸狀態(tài)的復(fù)雜性，在假設(shè)層裂片/孔洞的大小、間隙相等情況下，本節(jié)將該初始拉伸狀態(tài)分為僅含層裂片、僅含孔洞、同時(shí)存在孔洞與層裂片3 類典型情況（見圖2），使用歐拉方法進(jìn)行建模和直接模擬。圖2 中 α為層裂區(qū)初始孔隙度，為層裂區(qū)單位長(zhǎng)度上的貫通層裂片數(shù)量，為層裂區(qū)單位長(zhǎng)度上的孔洞數(shù)量。

圖2 相同初始孔隙度 α0=1.3 下不同層裂區(qū)建模Fig. 2 Different models of spallation zones with the same initial porosity of α0=1.3

直接數(shù)值模擬使用自主開發(fā)的三維多介質(zhì)彈塑性流體力學(xué)歐拉方法計(jì)算軟件MEPH進(jìn)行計(jì)算，其中金屬使用Mie-Grüneisen 狀態(tài)方程及Steinberg-Guinan 本構(gòu)模型進(jìn)行計(jì)算。

對(duì)于金屬銅的計(jì)算，首先，需要驗(yàn)證的是此種細(xì)觀建模情況下金屬再壓實(shí)問題直接模擬的正確性。在圖2(b)孔洞數(shù)量=12情況下 ，圖3 給出了三維多孔銅1/4 建模（該問題為純?nèi)S問題，二維軸對(duì)稱模型無(wú)法反映真實(shí)的邊界效應(yīng)），計(jì)算范圍[0,1] mm×[0,0.05] mm×[0,0.05] mm，計(jì)算網(wǎng)格長(zhǎng)度為 0.00125 mm，軸正方向設(shè)置為連續(xù)邊界條件、其他均設(shè)置對(duì)稱邊界條件。改變孔隙度 α、再壓實(shí)初始速率情況下，圖4 給出了8 組由MEPH 直接模擬得到的壓實(shí)密度與對(duì)應(yīng)壓力（實(shí)心點(diǎn)）。通過(guò)與實(shí)驗(yàn)結(jié)果（空心點(diǎn)）的對(duì)比來(lái)看，本文的直接建模模擬能夠適應(yīng)考察范圍附近的含孔隙物質(zhì)再壓實(shí)問題中壓實(shí)密度與壓力的模擬估計(jì)。

圖3 孔洞數(shù)量 nb=12 下的三維多孔銅1/4 模型Fig. 3 Three-dimensional 1/4 model of porous copper with nb=12

圖4 不同孔隙度下多孔銅三維直接模擬與實(shí)驗(yàn)結(jié)果對(duì)比Fig. 4 Comparison of 3D direct simulation and experimental results of porous copper at different porosity

表1 金屬銅的狀態(tài)方程參數(shù)[26]Table 1 The EOS parameters of Cu samples[26]

表2 金屬銅的SG 本構(gòu)模型參數(shù)[26]Table 2 The SG parameters of Cu samples[26]

回到金屬層裂區(qū)再壓實(shí)問題當(dāng)中，來(lái)考察圖2 中的幾種不同層裂區(qū)初始拉伸狀態(tài)對(duì)再壓實(shí)后物質(zhì)狀態(tài)帶來(lái)的影響。

在僅含層裂片（圖2(a)情況下，設(shè)置初始孔隙度 α=2 ，初始銅層裂片以- 1.0 km/s 速度向=0 固壁面運(yùn)動(dòng)（即再壓實(shí)速度=1 km/s ），層裂片數(shù)取不同值的模擬結(jié)果參見圖5～6。其中，圖5 給出了=8,20下層裂區(qū)在再壓實(shí)過(guò)程中典型時(shí)刻的材料內(nèi)部密度狀態(tài)；圖6(a)給出了不同層裂片數(shù)量下層裂區(qū)平均密度 ρˉ 隨再壓實(shí)時(shí)間變化曲線，曲線峰值即為層裂區(qū)再壓實(shí)密度；圖6(b)給出了層裂區(qū)再壓實(shí)密度隨建模中層裂片數(shù)的變化情況?？梢钥吹?，僅含層裂片情況下，隨著層裂片數(shù)的增加，層裂區(qū)再壓實(shí)狀態(tài)逐漸收斂，收斂情況下的再壓實(shí)狀態(tài)可作為宏觀模擬的比對(duì)參考值。

圖5 層裂片數(shù) nspall=8,20 下典型時(shí)刻的層裂區(qū)壓縮狀態(tài)Fig. 5 The compression state of the spallation zone at different times with nspall = 8, 20

圖6 相同初始孔隙度情況下平均密度的變化Fig. 6 Variation of the average density under the same initial porosity

在僅含孔洞（圖2(b)的情況下，設(shè)置初始孔隙度 α=1.3 ，同時(shí)在孔洞大小、間距相等假設(shè)下設(shè)置孔洞數(shù)=12。圖7～9 給出了再壓實(shí)速度=1,0.5,0.3 km/s 下的三維直接模擬結(jié)果，及對(duì)應(yīng)速度下層裂區(qū)平均密度隨時(shí)間變化情況。受平均密度統(tǒng)計(jì)方法等因素影響，孔洞塌縮射流擊穿層裂區(qū)外界面后的層裂區(qū)平均密度、再壓實(shí)密度已沒有實(shí)質(zhì)意義。由于層裂區(qū)再壓實(shí)后形成了界面向外的高速物質(zhì)噴射，本文中認(rèn)為此時(shí)的金屬層裂再壓實(shí)后已難以作為單一問題研究，應(yīng)分為孔洞閉合后的基體區(qū)和界面噴射區(qū)兩個(gè)子區(qū)域分別進(jìn)行考慮。

圖7 初始孔隙度1.3、再壓實(shí)速率1.0 km/s， nb=12 情況下的密度圖和層裂區(qū)平均密度隨時(shí)間的變化Fig. 7 Results of the density distribution and average density over time by direct simulations with α0=1.3, v0=1 km/s ,nb=12

圖8 初始孔隙度1.3、再壓實(shí)速率0.5 km/s， nb=12 情況下的密度圖和層裂區(qū)平均密度隨時(shí)間的變化Fig. 8 Results of the density distribution and average density over time by direct simulations with α0=1.3, v0=0.5 km/s ,nb=12

圖9 初始孔隙度1.3、再壓實(shí)速率0.3 km/s， nb=12 情況下的密度圖和層裂區(qū)平均密度隨時(shí)間的變化Fig. 9 Results of the density distribution and average density over time by direct simulations with α0=1.3, v0=0.3 km/s ,nb=12

對(duì)于同時(shí)含有孔洞、層裂片的層裂區(qū)（圖2(c)），保持初始孔隙度 α=1.3不變 ，設(shè)置孔洞數(shù)=10,層裂片數(shù)=10。需要指出的是=10,=10 情況下的建模并不唯一。圖10～12 給出了再壓實(shí)速率分別為=1.0,0.5,0.3 km/s 情況下的直接數(shù)值模擬典型時(shí)刻密度圖，及對(duì)應(yīng)情況下層裂區(qū)平均密度隨時(shí)間變化情況。可以看到，同時(shí)含有孔洞、層裂片的層裂區(qū)再壓實(shí)后的表面噴射較少、界面速度差異不大，并沒有出現(xiàn)純孔洞情況下那樣嚴(yán)重的表面噴射現(xiàn)象。此時(shí)統(tǒng)計(jì)得到的再壓實(shí)密度具有一定的參考比對(duì)價(jià)值。

圖10 初始孔隙度1.3、再壓實(shí)速率1.0 km/s， nb=10 ， nspall=10 情況下的密度圖和層裂區(qū)平均密度隨時(shí)間的變化Fig. 10 Results of the density distribution and average density over time by direct simulations with α0=1.3, v0=1 km/s ,nb=10, nspall=10

圖11 初始孔隙度1.3、再壓實(shí)速率0.5 km/s， nb=10 ， nspall=10 情況下的密度圖和層裂區(qū)平均密度隨時(shí)間的變化Fig. 11 Results of the density distribution and average density over time by direct simulations with α0=1.3, v0=0.5 km/s ,nb=10, nspall=10

圖12 初始孔隙度1.3、再壓實(shí)速率0.3 km/s， nb=10 ， nspall=10 情況下的密度圖和層裂區(qū)平均密度隨時(shí)間的變化Fig. 12 Results of the density distribution and average density over time by direct simulations with α0=1.3, v0=0.3 km/s ,nb=10, nspall=10

3 金屬層裂區(qū)再壓實(shí)過(guò)程的宏觀有限元模擬

在層裂區(qū)內(nèi)部層裂片/孔洞的大小、間隙相等假設(shè)下，初始時(shí)刻所有宏觀網(wǎng)格密度相等，即有初始時(shí)刻 ρ=ρ(=1,···,) 。圖13 給出了將層裂區(qū)初始均勻剖分為=20 個(gè)網(wǎng)格的宏觀計(jì)算模型。本節(jié)宏觀模擬使用基于四邊形網(wǎng)格的二維有限元方法程序，該程序在宏觀彈塑性流體力學(xué)模擬中具有較好的模擬效果。其中，宏觀模擬中使用的邊界條件、材料本構(gòu)及EOS 均與直接數(shù)值模擬相同。

圖13 宏觀上將層裂區(qū)初始均勻剖分為nelem=20Fig. 13 Finite element method is used to simulate the spallation zone described by nelem=20

網(wǎng)格斷裂后處理算法是宏觀層裂再壓實(shí)問題模擬準(zhǔn)確性的基礎(chǔ)，本文處理如下：對(duì)于網(wǎng)格斷裂后的應(yīng)力計(jì)算，結(jié)合式(2)、(4)可知，在網(wǎng)格斷裂后（＞）應(yīng)設(shè)置網(wǎng)格應(yīng)力為零；對(duì)于網(wǎng)格斷裂后的溫度計(jì)算，以最簡(jiǎn)單的保持?jǐn)嗔褧r(shí)刻溫度不變進(jìn)行處理，將網(wǎng)格斷裂狀態(tài)下塑性功、沙漏力等造成的溫升全部置零。

圖14 中給出了宏觀模擬平均密度 ρ（黑色線）與不同初始拉伸狀態(tài)下直接模擬得到的平均密度ρˉ（虛線）隨時(shí)間變化的對(duì)比，其中平均密度曲線峰值及對(duì)應(yīng)時(shí)刻即為再壓實(shí)密度及再壓實(shí)所用的間隔時(shí)間。由圖14 的對(duì)比看到，宏觀模擬所設(shè)置的斷裂后處理算法能夠較好比對(duì)層裂區(qū)僅含層裂片及同時(shí)含有層裂片及孔洞的情況，再壓實(shí)密度、再壓實(shí)所用的時(shí)間間隔及卸載歷程均基本相符。

圖14 不同拉伸狀態(tài)下的層裂區(qū)再壓實(shí)平均密度演化的對(duì)比Fig. 14 Comparison of the average recompression density under different tensile conditions

為了進(jìn)一步驗(yàn)證上述認(rèn)識(shí)，圖15 額外給出了3 組不同孔隙度、不同再壓實(shí)速度下的宏觀模擬平均密度與僅含層裂片的直接模擬結(jié)果（虛線）對(duì)比看到，在較大的孔隙度、較高的再壓實(shí)速度下=50 的宏觀模擬結(jié)果仍然能較好比對(duì)直接模擬得到的再壓實(shí)密度及時(shí)間間隔。=12的宏觀模擬結(jié)果在再壓實(shí)狀態(tài)上存在一定差異，說(shuō)明宏觀模擬中層裂區(qū)劃分網(wǎng)格數(shù)也會(huì)對(duì)該現(xiàn)象的模擬結(jié)果造成一定影響。即宏觀上較好模擬層裂再壓實(shí)問題的前提是層裂區(qū)的網(wǎng)格數(shù)需要超過(guò)一定數(shù)量，該網(wǎng)格數(shù)與宏觀模擬程序相關(guān)（經(jīng)驗(yàn)證對(duì)本文程序而言＞20 即可）。

圖15 不同孔隙度、不同再壓實(shí)速度下的宏觀模擬與直接模擬的對(duì)比Fig. 15 Comparison of the average density obtained from the macro- and direct simulations at different α0 and v0

與圖14 的對(duì)比同時(shí)看到，宏觀模擬不能較好應(yīng)對(duì)的是拉伸狀態(tài)為純孔洞的層裂區(qū)再壓實(shí)過(guò)程。這已在第3 節(jié)進(jìn)行了初步說(shuō)明：純孔洞再壓實(shí)情況下產(chǎn)生的界面噴射行為是導(dǎo)致宏觀模擬無(wú)法較好比對(duì)的原因。其直接原因是：在孔隙度不大情況下僅含孔洞層裂區(qū)可能是連通的（如圖3 所示），層裂區(qū)后界面追趕形成的彈塑性波可以通過(guò)連續(xù)體傳導(dǎo)至層裂區(qū)前界面從而造成層裂區(qū)前界面提前減速，而宏觀模擬中在斷裂網(wǎng)格再壓實(shí)前彈塑性波不會(huì)跨過(guò)斷裂網(wǎng)格向前傳導(dǎo)。

4 總 結(jié)

本文中考慮金屬層裂區(qū)初始進(jìn)入再壓實(shí)過(guò)程時(shí)真實(shí)拉伸狀態(tài)的復(fù)雜性，將金屬層裂區(qū)初始拉伸狀態(tài)設(shè)置為僅含層裂片、僅含孔洞、同時(shí)含有孔洞與層裂片3 類情況，對(duì)不同孔隙度、再壓實(shí)速率、層裂片數(shù)及孔洞數(shù)下的層裂再壓實(shí)過(guò)程開展了直接數(shù)值模擬研究。在驗(yàn)證含孔直接建模模擬正確性的基礎(chǔ)上，統(tǒng)計(jì)得到了多種工況下的層裂區(qū)再壓實(shí)密度、再壓實(shí)時(shí)間間隔等物理量，為對(duì)應(yīng)工況下的宏觀模擬提供比對(duì)依據(jù)。在直接模擬中還觀察到：拉伸狀態(tài)為純孔洞且再壓實(shí)速率較高情況下，即使層裂區(qū)外界面平整光潔，在再壓實(shí)后仍然會(huì)產(chǎn)生一定程度的外界面物質(zhì)噴射現(xiàn)象。

由于宏觀模擬不能區(qū)分網(wǎng)格內(nèi)部拉伸狀態(tài)，同時(shí)探討了宏觀模擬金屬層裂再壓實(shí)過(guò)程的可靠性問題。得到了如下初步認(rèn)識(shí)：斷裂后處理算法使用了應(yīng)力置零和溫度不變情況下，宏觀模擬可以較好模擬含層裂片的層裂區(qū)再壓實(shí)過(guò)程，再壓實(shí)密度、時(shí)間間隔、卸載路徑均能夠與直接模擬結(jié)果較好符合，其較好比對(duì)的前提是層裂區(qū)所剖分的宏觀網(wǎng)格數(shù)超過(guò)一定數(shù)量；若層裂區(qū)內(nèi)部?jī)H含孔洞，則無(wú)論孔洞塌縮是否形成界面噴射宏觀模擬均無(wú)法較好反映該層裂再壓實(shí)過(guò)程。上述認(rèn)識(shí)僅為金屬層裂再壓實(shí)問題的宏觀精細(xì)化模擬提供支撐，在模擬精度要求不高或再壓實(shí)時(shí)間歷程較短的情況下對(duì)斷裂后處理算法和層裂區(qū)網(wǎng)格剖分可不做要求。