陶瓷巴西圓盤動態劈裂的離散元模擬

趙翰卿，任會蘭

(北京理工大學 爆炸科學與技術國家重點實驗室， 北京 100081)

陶瓷材料具有低密度、高硬度、高強度等優秀性能，因此被廣泛應用于抗中高速沖擊防護裝甲設計[1-2]，軍事應用前景十分廣闊。

侯海量等[3]研究了高速破片沖擊下輕型陶瓷/金屬復合裝甲的沖擊響應，發現陶瓷面板在高速破片的高速沖擊下，首先受沖擊壓縮誘發微裂紋，持續壓縮致使材料進一步損傷，材料強度降低，隨后微裂紋在反射拉伸波作用下發生擴展，形成宏觀裂紋而發生失效。因此，沖擊載荷作用下，陶瓷材料的裂紋擴展機理及動態拉伸強度對陶瓷裝甲性能的研究十分重要。

目前關于陶瓷材料動力學性能的理論、實驗研究有很多，但主要集中于沖擊壓縮載荷條件。寧建國等[4]基于翼型裂紋擴展過程中的能量守恒原理，建立了陶瓷類脆性材料的動態損傷本構模型。張曉晴等[5]利用改進的SHPB實驗方法對氧化鋁陶瓷的動態力學性能進行了研究。

陶瓷等脆性材料的脆性使得直接測量動態拉伸強度十分困難，因此相關的實驗研究較少。基于分離式霍普金森桿(SHPB)的巴西圓盤試驗憑借其巧妙的設計，成為一種有效的間接測量手段。然而陶瓷等脆性材料的動態裂紋擴展過程十分迅速，憑借現有的技術手段很難完整觀測，而數值模擬手段可以作為實驗研究的有效補充。

美國學者Cundall[6]提出的離散元方法，對比傳統的有限元方法，在求解大位移和非線性問題方面優勢明顯[7]。本文基于顆粒離散元法(軟件PFC2D)，建立分離式霍普金森壓桿數值實驗平臺，對陶瓷巴西圓盤的動態劈裂進行數值模擬，為開展實驗研究提供支持。

1 建立陶瓷離散元模型

顆粒離散元法中的本構模型是針對接觸的細觀本構模型，模型參數一般包含接觸的剛度和強度參數。PFC的本構模型中的參數無法直接引入宏觀參數，需要以目標材料宏觀力學性能參數為標準，通過一系列數值試驗校準模型細觀參數。當一套細觀參數能同時滿足對應材料的各種宏觀力學性能時，便可認為這套細觀參數能夠模擬目標材料[6]。

1.1 選擇黏結模型

離散元中的細觀本構模型又稱黏結模型，目前使用最廣的模型為平行黏結模型(parallel-bonded model)，該模型被廣泛應用于模擬巖石等脆性材料[8-9]。但是在材料參數標定過程中，使用該模型往往會得到偏高的單軸拉壓強度比[6]，主要原因是顆粒間的平行黏結鍵在達到破壞條件后會退化為線性接觸，無法對顆粒的旋轉施加約束，導致模型的壓縮強度相對偏低。

平節理黏結模型(Flat-joint Model)是D.O.Potyondy[10]提出的一種適用于硬質巖石的改進模型，見圖1。該模型通過將虛擬交界面分段的方式，將球形顆粒構造成多邊形顆粒，當平節理黏結被破壞后，平節理黏結鍵不再提供黏結作用，但是不會被刪除，繼續起到抑制顆粒旋轉的作用，提高裂紋合并的難度，以此相對提高模型的壓縮強度，從而使模型的拉壓強度比更接近實驗室實驗結果。

圖1 平節理黏結模型示意圖

本文將選取平節理黏結模型(flat-joint model)作為模擬陶瓷材料的細觀本構模型，陶瓷材料的宏觀力學參數取自我們對中國建材研究院提供的15%氧化鋯增韌氧化鋁陶瓷進行實驗室試驗所得的實測數據。

1.2 選擇模型初始參數

根據Potyondy和Cundall的研究[6]，離散元顆粒模型的斷裂韌性與模型顆粒半徑及顆粒間黏結鍵的拉伸強度有關，因此顆粒半徑可以暫定一個大小，在后續標定斷裂韌性時重新調節。最小顆粒半徑暫定為20 μm，最大粒徑與最小粒徑之比為1.5，模型的孔隙率為0.1。離散元中的孔隙率只是用于估算填充模型區域所需顆粒數的數學參數，并非表征真實材料的氣孔率。

顆粒間接觸模量、黏結模量、黏結強度初始值可設置為材料真實彈性模量，接觸、黏結剛度比可暫定初始值為2，在后續標定過程中進行調節。顆粒間接觸摩擦系數對宏觀強度標定影響不大[11]，可取經驗值0.577。

1.3 通過數值試驗標定參數

1.3.1單軸壓縮試驗

通過單軸壓縮試驗來校準模型的彈性模量E、泊松比ν和抗壓強度σc。壓縮試驗的幾何模型：高3.4 mm，寬1.7 mm 的矩形，包含6 184個顆粒，上下兩面剛性墻體分別以17 mm/s的速度向中間壓縮。模型破壞過程形貌及應力應變曲線如圖2所示。校準后得到的模型彈性模量為 350 GPa，泊松比0.236，抗壓強度2 843 MPa。

圖2 單軸壓縮試驗破壞形貌及應力-應變曲線

1.3.2 三點彎曲試驗

通過三點彎曲試驗來校準模型的彎曲強度σb。三點彎曲試驗的幾何模型為：長6 mm，寬1.2 mm的矩形梁，包含7 706個顆粒，梁下方為兩個相距4.8 mm的圓形剛體，梁上方為12 mm/s速度向下施壓的圓形剛體。模型破壞過程形貌及載荷-位移曲線如圖3。將模型斷裂時的最大載荷F=118.607 kN代入三點彎曲強度計算公式，得到校準后的彎曲強度σb=593.0 MPa。

圖3 三點彎曲試驗破壞形貌及載荷-位移曲線

1.3.3三點彎曲單邊切口梁試驗(斷裂韌性測試試驗)

通過帶單邊切口的三點彎曲試驗來校準模型的斷裂韌性KIC。幾何模型為：長6 mm，寬1.2 mm的矩形梁，梁的正下方有一個寬40 μm，高0.6 mm的切口，模型包含7 677個顆粒，梁下方為兩個相距4.8 mm的圓形剛體，梁上方為12 mm/s速度向下施壓的圓形剛體。模型破壞過程形貌以及載荷-位移曲線如圖4。將模型斷裂時的最大載荷F=21.116 kN代入三點彎曲斷裂韌性計算公式，得到校準后的斷裂韌性為KIC=6.498 MPa·m1/2。

圖4 斷裂韌性試驗破壞形貌及載荷-位移曲線

1.3.4巴西劈裂試驗

通過準靜態巴西劈裂試驗來校準模型的拉伸強度。模擬所采用的巴西圓盤直徑為5 mm，包含21 083個顆粒。左右兩端剛性墻體分別以5 mm/s的速度向中間壓縮。模型破壞過程形貌以及載荷-位移曲線如圖5。將圓盤劈裂時的最大載荷F=2 346.243 kN代入巴西圓盤試件抗拉強度計算公式，得到校準后的拉伸強度為299 MPa。

圖5 巴西劈裂試驗破壞形貌及載荷-位移曲線

在上述數值實驗中，涉及載荷的讀數都很大，這是因為模擬采用2D模型，模型的厚度為單位1(1 m)。表1為離散元模型的校準結果與實驗測量結果，由表1可以看出各項宏觀力學性質參數均符合得很好，可以認為此離散元材料模型能夠很好地模擬15%ZrO2-Al2O3增韌陶瓷材料。表2為校準后離散元模型中主要細觀參數的取值。

表1 15%ZrO2-Al2O3增韌陶瓷的實驗數據與數值模擬結果

表2 15%ZrO2-Al2O3增韌陶瓷材料離散元模型的主要細觀參數

實驗測量結果缺少拉伸強度，一般脆性材料的拉壓強度比在0.05～0.15，數值試驗所得拉壓強度比為0.105，與實際情況相符。

2 陶瓷巴西圓盤動態劈裂的模擬

2.1 建立SHPB數值實驗平臺

參考文獻[12-13]，基于離散元構建二維SHPB結構模型，圖6為結構示意圖。入射桿與透射桿長度均為300 mm，撞擊桿長30 mm，圓桿直徑10 mm，為了改善桿端面與陶瓷圓盤試件的接觸情況，將桿端面邊界顆粒細化為直徑25 μm的細小顆粒，撞擊桿與入射桿接觸端面顆粒也進行細化對齊處理。桿件材料采用與陶瓷圓盤試件同樣的平節理黏結模型，桿件之間以及試件與桿件之間的接觸則采用線性接觸模型，模擬接觸但不黏結的接觸狀態。

圖6 離散元SHPB模型結構示意圖

桿件材料的部分宏觀參數取自文獻[14]，見表3。由于離散元中彈塑性材料標定較為復雜，且本次模擬過程中不涉及桿件發生損壞的情況，因此將桿件材料簡化為理想彈性體，顆粒間的接觸拉伸強度與接觸剪切強度均設置為極大值(100 GPa)。綜合考慮數值模擬精確性與計算資源，校準后桿件材料的細觀參數如表4所示。

表3 SHPB桿件材料參數

表4 SHPB桿件材料離散元模型的主要細觀參數

利用SHPB對脆性材料進行高速沖擊加載時，為保證SHPB實驗的均勻性要求，需要設法延長加載波的上升前沿時間，讓應力波有足夠長的時間在試件內部反射至均勻。常用的波形整形技術有添加波形整形器、使用異形炮彈[15]等，可以將傳統SHPB的矩形波整形為三角波以及半正弦波。

本文通過控制撞擊桿的運動速度以半正弦曲線模式變化，對入射桿端面施加動應力，實現近似半正弦波的加載波。因此在本文模擬中，入射應力脈沖的寬度與撞擊桿長度無關，由設置的加載作用時間控制。運動速度曲線表達式為v=vm*sin(πt/T)，其中vm為加載速度峰值，t為持續時間，T為半正弦波的作用時間。

2.2 驗證SHPB系統模型的可靠性

在入射桿與透射桿上設置一系列半徑為3.2 mm的測量圓(見圖6)，來監測桿件中應力波的傳播情況。測量圓1、圓2分別位于入射桿與透射桿的正中間，測量圓1、3、4、5的間距為30 mm。不插入試件，在空桿狀態下，設置撞擊桿的加載速度峰值為13 m/s，作用時間32 μs，在位置1～5監測到的應力波形如圖7所示。

圖7 撞擊桿加載速度峰值13 m/s、持續時間32 μs時入射和透射波形

由圖7可見，撞擊桿的動應力加載產生了近似半正弦波的波形，上升沿和下降沿均約為16 μs。自左向右，測量圓5、4、1、3、2測得的各個位置的應力時程曲線均無明顯彌散現象，入射波(measure 1)與透射波(measure 2)形狀、幅值基本一致，由此可以認為建立的SHPB系統模型是可靠的。

2.3 陶瓷巴西圓盤動態劈裂試驗的模擬

陶瓷巴西圓盤試件模型的直徑為8 mm，包含53 935個顆粒。設置撞擊桿的加載峰值速度為13 m/s，作用時間32 μs，產生如圖7所示的入射波形對陶瓷圓盤試件進行沖擊加載，模擬基于SHPB實驗裝置的陶瓷巴西圓盤動態劈裂試驗。模擬所得的入射波、反射波及透射波的應變時域曲線如圖8所示。

圖8 SHPB數值實驗平臺記錄的入射、反射與透射波

對于傳統巴西圓盤試驗，圓盤中心位置的拉伸應力σt可計算為：

(1)

式(1)中：P(t)為徑向壓縮載荷；D為圓盤的直徑；t為圓盤的厚度(2D模型厚度為單位1，即1 m)。

基于一維應力波理論，徑向壓縮載荷P(t)可按照式(2)計算。

P(t)=EA0εT(t)

(2)

式(2)中：E為桿的彈性模量；A0為桿的橫截面積(2D模型中為桿的直徑)；εT(t)為透射波時域曲線。

根據模擬所得應變時域曲線按照式(1)計算的圓心拉伸應力，得到圓盤中央拉伸應力峰值(即動態拉伸強度)為231 MPa。

3 模擬結果與實驗結果對比分析

圖9為基于14.5 mm直徑分離式霍普金森桿實驗裝置的陶瓷巴西圓盤動態劈裂過程，圓盤試件直徑14 mm，材料為15%氧化鋯增韌氧化鋁陶瓷，撞擊桿的速度為13.082 m/s。

從圖9可以看出，陶瓷圓盤試件在17.37 μs時，在與入射桿、透射桿接觸位置出現細小裂紋；在34.74 μs時，圓盤已沿加載直徑方向被主裂紋貫穿，裂紋擴展速度超過了高速攝相的曝光速度，同時在圓盤受壓區域附近生成了次生裂紋；在52.12 μs時，更多的次生裂紋生成和匯聚，圓盤試件完全劈裂破壞。高速攝相未能完整記錄裂紋擴展的全程。

圖9 陶瓷巴西圓盤劈裂過程的高速攝相

圖10為實驗所得入射波、反射波及透射波波形數據。由圖10可以看出，反射波存在一個明顯的轉折點并緊隨一個急速上升的信號，這說明波在試件內傳播過程中，波阻抗出現了明顯的變化。這種波阻抗迅速的減小對應了試件的瞬間破壞，并將加載信號完全反射回入射桿內。試件由初始加載到完全破壞的響應時間基本和透射波的信號相等，證實了試件破壞的時間節點為圖10中的轉折點。

圖10 實驗記錄的波形

圖11與圖12分別為撞擊桿速度10.427 m/s與11.918 m/s時，高速攝相記錄的破壞過程與實驗記錄的波形，與撞擊桿速度13.082 m/s時的實驗結果具有相同的特征。

圖11 撞擊桿速度10.427 m/s時的高速攝相記錄

圖12 撞擊桿速度11.918 m/s時的高速攝相記錄

模擬結果所得應變時域曲線也具有與實驗相似的轉折點。從模擬結果中截取從反射波信號初始加載時刻到透射波完全卸載時刻的數據點，計算反射波由初始加載到轉折點的響應時間，與透射波信號的響應時間作比較，結果如圖13所示。可以發現二者的響應時間基本相等，這說明模擬結果中的轉折點也是試件破壞的節點，模擬結果與實驗結果的數據曲線具有相同的形態特征。

圖13 由初始加載到轉折點與透射波信號響應時間的反射波信號

根據圖10、圖11、圖12中實驗所得應變時域曲線，按照式(1)計算得到動態拉伸強度分別為205 MPa、208 MPa和187 MPa。同二維離散元模型相比，陶瓷巴西圓盤試件是三維的；此外，陶瓷材料是由粉體材料燒結而成，材料內部不可避免存在有微孔洞、微裂紋等初始微缺陷，實驗結果比模擬結果的動態拉伸強度稍低。

以圓盤試件端面出現明顯微裂紋的時刻為時間起點，數值模擬的陶瓷巴西圓盤動態劈裂裂紋擴展過程如圖14所示。由圖14可以看出，增韌陶瓷巴西圓盤動態劈裂的過程中，試件從端面開始產生明顯微裂紋，到裂紋從端面擴展、匯集至沿徑向貫穿圓盤。陶瓷材料內有大量的微缺陷，如晶界裂紋、微孔洞、微裂紋等，在外載荷下這些微損傷將進一步成核、擴展直至匯合而引起陶瓷材料的宏觀破壞。

圖14 數值模擬的裂紋擴展過程及破壞形貌示意圖

4 結論

1) 基于顆粒離散元法，獲得了能夠模擬15%氧化鋯增韌氧化鋁陶瓷的離散元模型細觀參數，并成功建立了可靠的SHPB數值實驗平臺；

2) 實驗與數值模擬所得數據曲線具有相同的形態特征：兩者的反射波信號均存在一個明顯的轉折點并緊隨一個急速上升的信號，該轉折點即試件破壞的時間節點；

3) 模擬結果顯示：陶瓷圓盤動態破壞是試件與彈性桿接觸處最先產生裂紋，裂紋聚集形成局部損傷區域，并沿軸線向中心擴展而導致試件劈裂破壞，是一個動態過程，圓盤試件的破壞形貌與實驗結果基本一致。