基于混凝土細觀模型的高壓狀態方程數值模擬研究

李昊陽，孔祥振,2，方 秦，彭 永，崔 健

(1.陸軍工程大學 國防工程學院，南京 210007；2.東南大學 土木工程學院，南京 211189；3.國防科技大學 文理學院，長沙 410073；4.天津大學 建筑工程學院，天津 300072)

混凝土材料廣泛應用于軍用和民用工程結構中，當受到炸藥爆炸和彈體侵徹等強動載作用時，近區混凝土處于10 GPa量級的高圍壓狀態，混凝土呈現近似流體特性，但仍具有明顯的剪切強度[1]。在流體動力學軟件中，為較好描述上述高圍壓特性，材料動態力學行為對于偏量和球量分別處理，應力偏量和應變偏量之間的關系描述了材料的剪切行為，通過本構關系(或破壞面)表示；而應力球量和應變球量的關系描述了材料的壓縮行為，通過高壓狀態方程(equation of state，EoS)表示。因此混凝土材料高壓狀態方程的準確性對于描述其高圍壓下的壓縮行為至關重要[2]。

當前獲得混凝土高壓狀態方程實驗手段主要有3類：靜水壓縮實驗、平板撞擊實驗和全尺寸爆炸實驗。靜水壓縮實驗通常利用伺服裝置施加壓力，但受設備限制，所施加的最大靜水壓力在1 GPa以下[3-4]，且實驗加載速率在準靜態水平，與強動載下的高應變率不符。平板撞擊實驗通常利用輕氣炮發射金屬高速飛片撞擊混凝土試件并產生強間斷波，通過記錄混凝土試件中沖擊波速和粒子速度，結合Rankine-Hugoniot方程，進而得到混凝土的高壓狀態方程。Hall等[5-8]通過平板撞擊實驗，得到了不同類型混凝土的高壓狀態方程，但受輕氣炮口徑限制，飛板和混凝土試件尺寸較小(上述實驗中混凝土試件直徑約為50 mm、厚度約為15 mm)，難以反映混凝土中粗骨料(通常粒徑在10～25 mm)引起的非均質性。全尺寸爆炸實驗可采用較大尺寸混凝土試件以充分考慮混凝土的非均質性，通過測量炸藥爆炸在混凝土試件中產生沖擊波的傳播得到高壓狀態方程。Rinehart等[9-10]通過全尺寸爆炸實驗得到了混凝土高壓狀態方程，但炸藥爆炸產生的振動會對測量設備造成較大影響，使得測量穩定性不易控制。此外，平板撞擊實驗和全尺寸爆炸實驗費時耗力，難以對各組分性質(如水灰比、骨料級配和體積率等)對高壓狀態方程的影響進行精細化系統研究。

利用基于混凝土細觀模型的數值模擬方法確定高壓狀態方程，可同時解決平板撞擊實驗存在的非均質性和全尺寸爆炸實驗存在的測試穩定性問題，且易于進行參數分析，明確各組分對高壓狀態方程的影響規律。該方法通常在細觀尺度上將混凝土看作粗骨料、砂漿及其膠結面ITZ(interfacial transition zone)三相介質，通過明確各組分的物理力學特性進而對宏觀尺度的混凝土力學行為進行描述。上述基于細觀模型的數值模擬方法已廣泛應用于強動載下混凝土材料的動態拉伸破壞[11]、動態尺寸效應[12-13]和沖擊波傳播衰減規律[14-15]的研究。對于本文研究的高壓狀態方程，考慮到沖擊壓縮過程中ITZ作用有限，因此可將混凝土看作由骨料和水泥石組成的兩相介質。基于兩相介質的細觀模型，Riedel等[16]通過在混凝土試件上施加不同粒子速度，得到了常規和高強兩種不同砂漿組成的混凝土高壓狀態方程，但其將骨料形狀簡化為球形，不能反映真實的骨料形狀。通過對比基于上述細觀模型得到的高壓狀態方程與已有宏觀尺度上的高壓狀態方程(AUTODYN軟件中自動生成的P-α狀態方程)，Gebbeken等[17]對上述方法進行了驗證，進一步探討了3種不同粒徑骨料對高壓狀態方程的影響。為考慮混凝土中自由水對其高壓狀態方程的影響，Huang等[18]通過有效應力和混合法則兩種方法將自由水引入水泥石高壓狀態方程中，得到了飽和混凝土的高壓狀態方程，并發現相比干混凝土，飽和混凝土的高壓狀態方程斜率(體積模量)更大。上述研究中，最基本假定為破壞面對高壓狀態方程影響較小，但均未進行系統分析驗證該觀點，此外研究均局限于單一配合比混凝土，難以推廣應用至任意配合比混凝土。

針對已有研究不足并基于Peng等[19]的前期工作，本文首先建立了兩相介質混凝土細觀模型，通過擬合實驗數據確定了骨料的狀態方程并采用兩種方法(擬合數據和改進的Hugoniot混合法則)確定了水泥石的狀態方程；然后基于細觀模型得到了混凝土的高壓狀態方程，并得到了已有實驗數據驗證，該方法可推廣至任意配合比的混凝土材料；最后在此基礎上，詳細探討了破壞面參數和水灰比對高壓狀態方程的影響。

1 高壓狀態方程

固體高壓狀態方程描述了材料靜水壓力、密度和能量三者的變化，用3個維度的空間曲面(Hugoniot曲面)來描述，曲面上的任意一點均代表材料的某一特定狀態，如圖1(a)所示。對于本文關注的混凝土材料，強沖擊(壓力10 GPa量級)通常引起的溫度和能量變化不大，因此其高壓狀態方程可簡化為二維平面上的壓力和密度(用體積應變描述)的關系曲線，如圖1(b)所示。混凝土為典型的多孔非均質材料，孔隙變化對高壓狀態方程影響顯著，已有實驗數據表明，混凝土高壓狀態方程表現為3個明顯階段：彈性段、孔隙壓實段和密實段，對應的壓力分界點分別為PHEL和Plock，其中PHEL為一維平面應變壓縮下材料發生非彈性變形對應的壓力(通常稱為Hugoniot彈性極限)。當壓力小于PHEL時材料處于線彈性階段，超過PHEL后混凝土中的孔隙逐漸坍塌壓實且加載體積模量逐漸變大，而當壓力超過壓實壓力Plock時，混凝土孔隙完全壓實進入密實狀態且體積模量陡增。

圖1 固體高壓狀態方程示意圖Fig.1 Schematic diagram of EoS at high pressures of solid materials

固體中沖擊波傳播的示意圖，如圖2所示。由于沖擊波速Us遠超材料粒子速度Up，在相同時間dt內，沖擊波陣面傳播至3-3，而粒子僅運動至2-2。

圖2 一維平面沖擊波傳播過程Fig.2 One-dimensional plane shock wave propagation process

考慮時間dt極短，波陣面前后物理量可通過質量和動量守恒方程聯系，即

ρ0(Us-Up0)=ρ(Us-Up)

(1)

P-P0=ρ0(Us-Up0)(Up-Up0)

(2)

式中：ρ為密度；P為壓力；帶有下標0的為沖擊波陣面(3-3)前；不帶下標0的為沖擊波陣面(2-2)后，因此在已知Us和Up情況下結合式(1)和式(2)即可得到混凝土的高壓狀態方程。

2 三維細觀模型的建立

考慮到沖擊壓縮過程中ITZ作用有限，ITZ對材料的影響耦合在水泥石中，因此在細觀尺度上將混凝土視為骨料和水泥石組成的兩相介質，其中，骨料為粒徑大于1 mm的礫石，水泥石包括粒徑小于1 mm的細砂、各類孔隙結構和水泥漿。基于Peng等提出的隨機骨料生成投放程序，可建立的兩相介質幾何和有限元模型，如圖3所示。該程序主要步驟為：①采用正六面體網格對混凝土試件進行網格劃分，并將所有網格單元定義為水泥石；②生成一個隨機尺寸的骨料幾何模型(采用隨機多面體描述)并隨機投放進試件中，落入骨料幾何范圍的網格單元重新定義為骨料；③重復步驟②直至骨料級配及體積率滿足要求。

圖3 混凝土細觀模型示意圖Fig.3 Schematic diagram of concrete mesoscopic model

考慮一維應變平面波的傳播特性，混凝土細觀模型尺寸可根據骨料最大粒徑Dmax確定，X和Y(Z)方向尺寸分別為20Dmax和5Dmax。在細觀有限元模型一側施加粒子速度產生沖擊波在試件中傳播，為達到一維應變狀態，施加邊界條件約束模型側表面的Y與Z方向位移，如圖4所示。為避免邊界效應，測點布置在距Y，Z方向邊界的1/5橫向尺寸外，同時為避免前后自由表面的影響，模型X方向前后1/5縱向尺寸內不設置測點。

圖4 細觀數值模擬方法示意圖Fig.4 Schematic diagram of mesoscopic numerical simulation method

混凝土中骨料和水泥石波阻抗和材料屬性不同，使得沖擊波陣面不再為平面，且波陣面后壓力和體積應變不均勻是混凝土材料非均質性的體現，因此需在測點區域內盡可能多的選擇單元用以表征該不均勻性。對于每次施加粒子速度，取300個單元作為測點，通過平均300個單元的壓力和體積應變，得到高壓狀態方程上的一個點，這樣通過施加不同粒子速度，可得到完整的混凝土高壓狀態方程。

3 材料模型及參數確定

第2章基于細觀模型確定高壓狀態方程的方法準確性取決于骨料和水泥石材料模型的準確性。基于作者前期開發的混凝土類材料的動態損傷材料模型(Kong-Fang模型)[20]，本章分別確定骨料和水泥石的破壞面和狀態方程。

Kong-Fang材料模型最初提出用以描述沖擊爆炸荷載作用下混凝土材料的動態力學性能，該模型為損傷塑性模型，并采用偏量(破壞面)和球量(狀態方程)分開處理的方式，并已推廣應用到巖石材料[21]，以下進行簡要介紹。

當前破壞面Y定義為

Y(σij,D)=r′[D(σr-σm)+σm]

(3)

式中：σij為應力張量；r′為當前子午線與壓縮子午線的比值；σm和σr分別為壓縮子午線上的最大強度面和殘余強度面，定義為

(4)

(5)

式中：參數a1和a2通過擬合三軸壓縮實驗數據確定；fc和T分別為單軸壓縮強度和拉伸強度；P為壓力；ψ為拉伸和壓縮子午線比率；D為損傷指數，D=1-(1-Dc)(1-Dt)，對于拉伸損傷Dt和壓縮損傷Dc分別定義為

(6)

(7)

(8)

為描述沖擊爆炸近區高靜水壓力特性，Kong-Fang模型采用LS-DYNA軟件中的8號狀態方程[22]來描述壓力和體積應變的關系，通過輸入10對壓力與體積應變和卸載模量進行線性差值確定。

3.1 破壞面參數

Kong-Fang材料模型中的破壞面參數主要包括：①強度面參數a1，a2，fc，T和模量E等；②損傷參數，壓縮損傷參數d1，d2，α和斷裂應變εfrac。

表1 混凝土組分材料部分參數Tab.1 Parameters of concrete components

3.2 狀態方程參數

輝綠巖(又稱粗玄巖)常用作混凝土中的骨料，基于Tsembelis等平板撞擊實驗和LASL實驗室得到的Hugoniot實驗數據，采用分段線性擬合方式，可得到其狀態方程，如圖5所示。考慮到巖石材料沖擊壓縮實驗數據較為豐富，其他巖石材料可采用相同的方法得到其狀態方程。

圖5 骨料分段線性高壓狀態方程Fig.5 Piecewise linear EoS at high pressures of aggregate

水灰比為0.35，水泥石的高壓狀態方程數據，如圖6所示。從圖6可知，實驗數據存在一定的離散性，采用分段線性擬合方式得到其狀態方程。

圖6 水泥石分段線性高壓狀態方程Fig.6 Piecewise linear EoS at high pressures of cement stone

水泥石是一種典型的多孔材料，孔隙的坍塌對其高壓狀態方程有顯著影響，而水灰比(w/c)決定了水泥石配合比和孔隙率，因此不同水灰比條件下水泥石的高壓狀態方程存在明顯差異，且已有研究針對水泥石的平板撞擊實驗數據較少，通過分段線性擬合方法得到其高壓狀態方程并不實用。

大量實驗數據表明，固體中的波速和粒子速度存在明顯的線性關系即：Us=cb+sUp，其中：cb為體聲速；s為斜率。Riedel等和Gebbeken等分別提出了基于質量和體積分數的Hugoniot混合法，用以確定多相介質的ρ，cb和s。Gebbeken等和Huang等的研究表明，基于體積分數的混合法具有更好的預測性。基于體積混合法則，水泥石的ρ，cb和s可由式(9)確定

H=Σfi·hi

(9)

式中：H和h分別為水泥石和其組分的ρ，cb，s；f為不同組分所占的體積分數；i為水泥石的不同組分。已有研究表明，水泥石中主要存在以下5種組分：未水化水泥、水泥凝膠、凝膠孔隙、毛細孔和氣孔，其中凝膠孔隙被凝膠孔水充滿，毛細孔被毛細孔水填充[28]。根據德國水泥協會的建議[29]，水灰比為0.35的水泥石以上5種組分的體積含量，如表2所示。基于式(9)和表2可得到水灰比為0.35的水泥石ρ，cb和s。

表2 水泥石組分性質及占比Tab.2 Cement stone component properties and proportion

水泥石可采用典型的多項式狀態方程，表達式為

P=k1μ+k2μ2+k3μ3

(10)

式中：μ為體積應變；k1，k2和k3可由式(11)～式(13)確定

(11)

(12)

(13)

可得到水灰比為0.35水泥石的狀態方程為P=28.02μ+29.97μ2+4.04μ3該狀態方程與實驗數據的對比，可以看出當體積應變小于0.25時，基于體積混合法則得到的高壓狀態方程與實驗數據吻合較好，而當體積應變大于0.25時，混合法則嚴重低估了實驗數據，如圖7所示。

水泥石在孔隙被完全壓實后，水泥凝膠的晶體結構會發生變化，因此表2中各組分體積分數及相應的ρ，cb和s也會隨之變化，因此Hugoniot混合法則不再適用。注意到水泥石壓實密度為2 584 kg/m3(通過混合法則得到)，相應的壓實應變為μ=ρ/ρ0-1=0.16，見圖7，密實段狀態方程可經驗性地用線性表示，其斜率(體積模量)可擬合為8K0，K0為水泥石的線彈性體積模量，可由基于體積分數的混合法則得到。可以看出經過改進得到的高壓狀態方程與實驗數據較為一致，壓力最大誤差在7%以內，已有平板撞擊實驗數據表明，水泥石類材料的密實段體積模量近于8倍的線彈性體積模量，進一步驗證了本文經驗性取值方法的正確性。

圖7 基于原始和改進的體積混合法則預測的高壓狀態方程與實驗數據的對比Fig.7 Comparison of predicted EoS at high pressures by the original and modified volume mixing rule with test data

4 細觀模型的預測結果及分析

本章首先基于已有平板撞擊實驗數據對細觀模型進行驗證，然后在細觀尺度上對沖擊波的傳播規律進行分析和討論。

Tsembelis等[30]進行了一組混凝土試件(直徑65 mm，厚度20 mm)的平板撞擊實驗，該實驗中混凝土由體積占比45.5%的輝綠巖骨料和54.5%的水泥石(水灰比為0.35)組成，具體級配和用料比等信息，如表3所示。建立與該實驗相同骨料體積含量和級配的混凝土細觀模型，基于第3章給出的輝綠巖骨料和0.35水灰比水泥石的材料模型，采用第2章給出的方法進行數值模擬，可得到該混凝土試件的高壓狀態方程。數值模擬得到的高壓狀態方程與實驗數據的對比，如圖8所示。由圖8可知：數值預測結果與實驗數據非常接近，驗證了基于細觀模型和改進的混合法則確定混凝土高壓狀態方程的正確性。

表3 混凝土組分占比Tab.3 The proportion of concrete components

由圖8可知：混凝土材料的Hugoniot彈性極限PHEL為0.4 GPa，這與Hall等、Grady和Gebbeken等利用平板撞擊實驗得到的結果相近。

圖8 基于細觀模型預測的高壓狀態方程與實驗數據的對比Fig.8 Comparison of numerically predicted EoS at high pressures with test data

圖9(a)給出了當粒子加載速度為500 m/s時，X軸方向50 mm，75 mm和125 mm 3個截面典型單元的壓力時程曲線，可以看出明顯的彈性前驅波平臺(約為0.4 GPa)、壓力升時段和塑性波平臺，壓力時程曲線的上下波動(在已有平板撞擊實驗中同樣發現)是由于骨料和水泥石波阻抗不同造成的復雜反射引起的。

圖9(b)給出了該工況下25 μs時的壓力分布云圖，可以看出明顯的鋸齒狀壓力分界線A-A，B-B和C-C(對應于圖9(a)的3個階段)，鋸齒狀是由混凝土材料的非均質性引起的。

圖9 細觀水平壓力分布Fig.9 Pressure distribution on mesoscopic level

5 參數敏感性分析

在計算流體動力學軟件中，材料行為對于偏量和球量分別處理，注意到在一維應變狀態下，應力偏量和球量同時存在，因此第3章基于細觀模型確定混凝土高壓狀態方程的方法中基本假定為偏量不會對球量產生較大影響，本章對破壞面參數進行敏感性分析驗證該假設的正確性。骨料體積含量及形狀的影響已在作者前期工作中進行了細致探討，這里不再贅述。考慮到水灰比對水泥石的性能影響較大，因此本章進一步研究水灰比對混凝土高壓狀態方程的影響并進行了單參數敏感性分析。

5.1 破壞面參數

為探討應力偏量(用破壞面描述)對混凝土高壓狀態方程的影響，本節首先對Kong-Fang模型中的強度面參數進行單參數敏感性分析(同時變化骨料和水泥石參數如a1)，即僅變化所討論參數(按±10%和±20%幅度變化)而保持其他參數不變。

Kong-Fang模型中最大強度面和殘余強度面大小由參數a1和a2控制，a1和a2對混凝土高壓狀態方程的影響，如圖10所示。由圖10可知：在討論的參數變化范圍內(±20%內)，a1和a2對高壓狀態方程幾乎沒有影響。

圖10 強度面參數敏感性分析Fig.10 Sensitivity analysis of strength surface parameters

Kong-Fang模型中損傷參數影響到壓縮和拉伸損傷累積的快慢，包括壓縮損傷參數d1，d2，α和拉伸損傷參數εfrac，損傷參數對高壓狀態方程的影響，如圖11所示。可以看出影響較小可以忽略。

圖11 損傷參數敏感性分析Fig.11 Sensitivity analysis of damage-related parameters

通過本章敏感性分析可得出，材料模型中的破壞面(即偏量部分)對高壓狀態方程的影響有限，因此在一維平面應變的高壓力狀態中偏量和球量可近似解耦，且在缺少實驗數據情況下細觀模型中骨料和水泥石的破壞面可不必精確確定，采用Kong-Fang模型中自動生成的參數即可滿足要求。

5.2 水灰比

由第3章可知，水泥石的組分配合比和孔隙率由水灰比決定，調整水灰比會同時改變水泥石的破壞面和狀態方程，本節考慮另一典型水灰比0.56，參考已有實驗數據[31]，該水灰比下的水泥石fc=35 MPa和E=30.7 GPa，其他強度面參數按照Kong-Fang模型自動生成得到。該水灰比下水泥石組分配合比可根據德國水泥協會的建議得到，進而基于本文提出的改進的混合法則得到其高壓狀態方程。

在相同的骨料級配下，基于細觀模型得到的兩種水灰比下的混凝土高壓狀態方程，如圖12所示。由圖12可知：水灰比對混凝土高壓狀態方程的影響非常明顯，隨著水灰比減小，水化作用產生的毛細孔數量大幅減少，水泥石中凝膠產物占比變大，混凝土組分間粘結力增強，混凝土材料的加載體積模量增大，即在相同體積應變下壓力更高。初始粒子速度800 m/s工況下t=17 μs時刻的壓力云圖，如圖13所示。由圖13可知：沖擊波在低水灰比下傳播更快，進一步說明了其加載體積模量更高。

圖12 不同水灰比混凝土高壓狀態方程Fig.12 EoS of concrete with different water-cement ratio

圖13 不同水灰比壓力云圖Fig.13 Pressure contour of different water-cement ratio

6 結 論

基于作者前期提出的混凝土兩相介質細觀模型和混凝土類材料的Kong-Fang材料模型，本文在細觀尺度上對混凝土高壓狀態方程的確定、基于體積分數的混合法則、沖擊波在兩相介質細觀尺度上的傳播規律以及強度面和水灰比敏感性等方面進行了詳細探討，主要研究工作和結論如下：

(1)提出了一種改進的基于體積分數的混合法則確定水泥石的高壓狀態方程，在壓實階段，仍采用原有體積分數混合法則；而在密實段，由于水泥凝膠晶體結構發生變化等因素，原有混合法則不再適用，基于已有實驗數據可得到了密實段的加載模量。

(2)在細觀尺度上分析了混凝土中沖擊波傳播規律，發現數值模擬和已有實驗中測得的壓力波動現象是由混凝土材料的非均質性造成的。

(3)參數敏感性分析結果表明，用于描述混凝土偏量的破壞面參數對高壓狀態方程影響較小，在一維平面應變的高壓力狀態中偏量和球量可近似解耦，因此在基于該方法確定混凝土高壓狀態方程時，破壞面可不必精確確定，極大簡化了模型參數標定的復雜度。而水灰比決定了水泥石的組分配合比和孔隙率，對高壓狀態方程影響較大，因此采用該方法時需明確混凝土的水灰比。

(4)通過細致了解混凝土中骨料和水泥石一維平面應變高壓力下的球量關系，可將本文基于細觀模型確定混凝土高壓狀態方程的方法推廣至任意配合比混凝土，數值模擬結果得到了已有平板撞擊實驗數據的較好驗證。