巖石中爆炸成坑的DDA模擬

孫淑美 陳 潔 范一鍇

(1.國家知識產權局專利局專利審查協作北京中心材料工程發明審查部，北京 100190;2.中國林業科學研究院，北京 100091;3.空軍第一建筑安裝工程總隊，北京 100076)

物理試驗是研究爆炸沖擊問題的重要手段，但由于爆炸荷載的特殊性，物理試驗有著固有的局限。隨著計算機技術和數學科學的發展，數值模擬能力越來越強大，逐漸成為一個與物理試驗互補的最佳研究手段。在物理試驗的邊界條件及提供的關鍵參數基礎上，基于數學的推導和計算機大規模運算，可以模擬整個物理試驗的過程，獲取研究者所關注的技術指標［1］。

非連續變形分析(DDA，Discontinuous Deformation Analysis)是近年來在非連續大變形分析領域的后起之秀，對于爆破這種伴隨大變形、開裂、材料塊體的分離運動效應的物理過程的模擬，具有獨特的優勢［2］。趙根、劉紅巖、朱傳云等［3-5］對 DDA 應用于爆破領域進行探討，證實其可行性，但在模擬過程中，爆炸荷載的取值和添加都采用經驗估計方法，本研究擬對非連續變形分析在爆破成坑過程的模擬上做出一些改進，即在建模和參數設置時考慮了巖體中爆炸破壞特性以及炸藥的狀態方程，規范其爆炸荷載的取值和添加方法，使之能更好地應用于爆破分析。

1 總體流程

DDA方法已具備塊體系統動力學分析的能力，但是用于爆炸荷載效應的計算，還需要考慮爆炸荷載的相關細節，簡言之就是爆炸荷載的計算和加載問題。炸藥引爆后除了剛開始的電磁、熱輻射、沖擊波，作用于介質時間最長的還是爆炸后生成的氣體。本模擬中所需要的爆炸荷載即爆轟產物(氣體)的壓力p，這一參數可以通過炸藥的狀態方程求出。狀態方程是用來描述爆轟產物的p—V—T關系，可以精確地反映爆轟產物的膨脹驅動做功過程。在圍巖中爆破時，爆炸近區發生急劇的變化，包括壓縮波的傳播、塑性變形的發生、介質運動的破壞以及爆炸空腔的形成。從炮孔向外依次可分為汽化區、液化區、破碎區、裂紋區和彈性區。炮孔和裂紋區之間的區域變化過程十分復雜，要對這部分進行模擬耗費巨大，目前DDA方法也不能勝任這一工作，但本研究重點是彈坑形成過程分析，關注這一區域也沒有必要，因此，可進行合理的簡化。在模擬中，不考慮裂紋區以內的區域，參考相關文獻，在炮孔周圍預設一定大小的“空腔區”，使荷載直接作用于“空腔區”外圍的塊體上。解決爆炸荷載計算和加載的問題后，結合DDA程序，可進行爆炸成坑的模擬，模擬的總體流程如圖1所示。本模擬用到的DDA程序由三個子程序組成:線條生成程序DL.exe、塊體切割程序DC.exe和分析計算主程序 DF.exe。其中DL是DC的前處理程序，它可以通過輸入少量的數據即可自動生成大量的線條，包括邊界線、節理線、洞室邊界線等在塊體生成中需要用到的各種屬性的線，生成的線的信息供DC使用;同樣的，DC是DF的前處理程序，通過DL提供的各種屬性的線的信息，生成帶外部邊界、內部洞室和材料屬性的模型，并對模型進行切割剖分，形成塊體系統，塊體系統信息供DF使用;塊體系統信息和外部荷載信息以及時步控制的數據文件讀入分析計算主程序DF進行計算，計算過程以圖形實時顯示。本文所用DDA程序為石根華博士提供。

圖1 總體流程圖

2 爆轟產物狀態方程

爆轟產物的狀態方程分為兩類:一類與化學反應有關，稱為顯含化學反應狀態方程，廣泛用于爆轟的熱力學計算，又叫熱力學狀態方程;另一類不處理化學反應，稱為不顯含化學反應狀態方程，廣泛應用于爆炸力學數值計算，在爆炸效應數值模擬方面取得比較滿意的效果，又叫動力學狀態方程［6］。典型的動力學狀態方程是JWL(Jones-Wilkins-Lee)狀態方程，它能描述爆轟產物的壓力—比容—內能(p—V—E)的關系，廣泛應用于數值模擬中，其表達式為:

其等熵方程為:

其中，p為壓力;E為爆轟產物單位體積的內能;V為爆轟產物的相對體積(即爆轟產物體積與裝藥體積之比);A，B，C，R1，R2，ω均為材料常數，需通過圓筒試驗標定得到。

JWL狀態方程這6個系數隨著炸藥種類和裝藥密度的變化而變化，要標定出各種炸藥在不同裝藥密度下爆轟產物的狀態方程，必須做大量的圓筒試驗，從工作量而言非常困難。針對這一問題，薛再清［7］提出了通用JWL(JWLG)狀態方程，用單位質量裝藥爆轟產物體積Vm代替爆轟產物相對體積V，適用于不同裝藥密度，其表達式為:

其等熵方程為:

其中，p為爆轟產物壓力;V為單位質量裝藥爆轟產物體積;EG為單位質量炸藥的爆熱;ρ0為裝藥密度;A，B，CG，Ra，Rb，ω 均為待定系數。

3 預設“空腔區”的計算

文獻［8］中認為空腔半徑R1與爆炸能量近似成正比，依賴于爆炸能量及介質的性質，包括其可壓縮性和強度，其計算可采用如下公式:

其中，ρ為介質的密度;c0為介質中聲波波速;fc為介質的單軸抗壓強度;Q為炸藥的能量;β=0.61。破碎區半徑Rd與空腔半徑成正比:

破碎區是指巖體在爆炸荷載作用下發生粉碎性破壞的區域，本文中爆炸成坑模擬時，直接忽略這一區域，預設“空腔區”半徑根據破碎區半徑Rd取值。

4 爆坑模擬案例

取靶體材料為以下參數的凝灰巖:密度ρ=2.83×103kg/m3，巖體中波速 c0=5 870 m/s，單軸抗壓強度 fc=114.5 MPa［9］。設藥包半徑a=0.1 m，為裝藥密度 ρTNT=1.65×103kg/m3的 TNT藥包，能量密度 Q0=6.803 MJ/kg［10］，則炸藥總能量:

將以上參數代入式(5)，式(6)，可得 R1=0.214 m，Rd=1.28 m，可見壓碎區半徑已接近藥包半徑的13倍，這一結果偏大，考慮到壓碎區也不是完全變成粉末，碎屑也可以承擔一定的壓力，模擬時所需的空腔半徑可在R1和Rd之間取值，取預設“空腔區”半徑為0.60 m。

根據文獻［7］，該TNT藥包的JWLG狀態方程的具體參數取值見表1。

表1 藥包所用TNT的JWLG狀態方程參數

將表1中的參數代入式(3)，式(4)可得爆轟產物氣體在0.6 m空腔中產生的壓力為2.464 MPa，設該荷載在20 μs衰減至0。

圖2 模型示意圖（單位：m）

建立的模型尺寸如圖2所示，巖體尺寸為20 m×10 m，藥包埋深為1 m，外部固定邊界為內部尺寸20 m×15 m，外部尺寸24 m×17 m，寬2 m的“凹”形。

爆破漏斗的形成過程如圖3所示，分別截取1 000步、2 000步和5 000步時的爆破漏斗瞬時狀態。

圖3 爆破漏斗形成過程

從該結果可以看出，用DDA進行爆炸成坑的模擬，由于預設了“空腔區”，從粉碎區外緣開始計算，此時爆炸作用相對于爆炸中心已經衰減不少，爆轟產物的壓力不足以使圍巖體破碎，而是沿著劃分單元時劃分出的節理面開裂。在接近自由面(原地面)的藥包上方，塊體被沖擊力拋起，呈輻射狀飛散。被拋出的塊體一部分塊體落回坑內，一部分落到原自由面(原地面)上，共同形成可見彈坑，可見彈坑的輪廓如圖3c)所示。將圖3c)與實際工程中觀測到的彈坑形狀對比可以發現，DDA模擬結果同樣存在明顯的回填物組成的區域以及開裂區，用DDA模擬的爆炸成坑效果與試驗結果符合較好。

5 結語

本文基于DDA程序，采用JWLG狀態方程計算爆炸荷載和引入預設“空腔區”解決爆炸荷載的計算和加載問題，模擬巖石中爆炸成坑的過程及最終可見彈坑的形態。采用實際巖石和炸藥參數進行了案例演示，模擬結果與工程實際相吻合，研究思路可行。

在用DDA模擬爆炸成坑過程的大量嘗試中發現，單元形狀、尺寸和模擬參數的選擇，對成坑過程以及最終的彈坑形狀和尺寸有較大的影響。要想將模擬結果用于指導工程建設，在材料、單元參數設計和模擬參數的選擇上還需要進一步的研究。

［1］范一鍇，黃 新.平頭彈低速侵徹黏土的數值模擬研究［J］.科學技術與工程，2012，12(22):5537-5540.

［2］范一鍇.撞擊及爆炸沖擊荷載效應的動力離心模型試驗和數值研究［D］.北京:北京航空航天大學，2012.

［3］趙 根.深水條件下圍堰拆除爆破技術研究［D］.合肥:中國科學技術大學，2008.

［4］劉紅巖，楊 軍，陳鵬萬.爆破漏斗形成過程的DDA模擬分析［J］.工程爆破，2004，10(2):17-20.

［5］朱傳云，戴 晨，姜清輝.DDA方法在臺階爆破仿真模擬中的應用［J］.巖石力學與工程學報，2002，21(S2):2461-2464.

［6］章冠人，陳大年.凝聚炸藥起爆動力學［M］.北京:國防工業出版社，1991.

［7］薛再清.爆轟產物狀態方程及含鋁炸藥的爆炸過程［D］.北京:北京理工大學，1998.

［8］錢七虎，王明洋.巖土中的沖擊爆炸效應［M］.北京:國防工業出版社，2010.

［9］水利水電科學研究院.巖石力學參數手冊［M］.北京:水利電力出版社，1991.

［10］周聽清.爆炸動力學及其應用［M］.合肥:中國科技大學出版社，2001.