張 彬，曾凡桂，王德璋，康官先，張曉雨，康天合

(1.太原理工大學(xué) 原位改性采礦教育部重點(diǎn)實驗室，山西 太原 030024; 2.太原理工大學(xué) 礦業(yè)工程學(xué)院，山西 太原 030024; 3.山西晉城無煙煤礦業(yè)集團(tuán)有限責(zé)任公司，山西 晉城 048006; 4.太原理工大學(xué) 煤與煤系氣地質(zhì)山西省重點(diǎn)實驗室，山西 太原 030024; 5.太原理工大學(xué) 安全與應(yīng)急管理學(xué)院，山西 太原 030024)

沁水盆地南部的煤層氣開發(fā)已實現(xiàn)了小型商業(yè)化，但仍存在產(chǎn)量低和衰減快的問題[1-2]，往往需要壓裂改造提高抽采率[3]，而煤層本身的力學(xué)性質(zhì)是決定壓裂效果的關(guān)鍵因素[4]。陳立超和王生維[5]利用測井?dāng)?shù)據(jù)反演了鄭莊區(qū)塊煤巖的彈性模量、剪切模量、體積模量和泊松比等參數(shù)。馮晴等[6]實驗得出彈性模量和泊松比隨埋深的變化規(guī)律。實際上，沁水盆地?zé)o煙煤儲層是普遍含氣的[7]，含甲烷條件下的煤體力學(xué)性質(zhì)更符合生產(chǎn)實際。辛程鵬等[8]對新景礦含瓦斯突出煤進(jìn)行了常規(guī)三軸和分段變速加載力學(xué)試驗研究，結(jié)果表明:隨瓦斯壓力降低，煤體強(qiáng)度和彈性模量均增大;煤體在分段變速加載路徑下的強(qiáng)度普遍增大，峰值軸向應(yīng)變、峰值環(huán)向應(yīng)變絕對值和峰值體積應(yīng)變絕對值也普遍增大，失穩(wěn)破壞瞬間應(yīng)力跌落和能量釋放更加劇烈。孟筠青等[9]采用分子模擬的方法對趙莊煤大分子模型進(jìn)行力學(xué)性質(zhì)研究，并分析了徑向分布曲線和鍵角分布曲線，發(fā)現(xiàn)煤大分子中的環(huán)狀結(jié)構(gòu)在抵抗破壞過程中起到的重要作用，提出了研究納米尺度下煤的力學(xué)性質(zhì)的新思路，是揭示氣體與固體相互作用機(jī)制的有力手段[10]。筆者以沁水盆地15號煤為研究對象，對其結(jié)構(gòu)進(jìn)行表征并構(gòu)建其大分子結(jié)構(gòu)模型。在模型的基礎(chǔ)上，采用巨正則系綜蒙特卡洛(GCMC)對甲烷在無煙煤中的吸附特性進(jìn)行研究，獲得吸附量、吸附熱和吸附構(gòu)型;采用分子動力學(xué)的方法對吸附甲烷后的構(gòu)型進(jìn)行力學(xué)性質(zhì)計算，從分子角度揭示含甲烷無煙煤的力學(xué)性質(zhì)，為沁水盆地煤層氣抽采提供一些微觀的理論參數(shù)。

1 研究方法和研究方案

分子模擬研究含甲烷無煙煤力學(xué)性質(zhì)的方法為:① 無煙煤化學(xué)結(jié)構(gòu)表征;② 基于表征數(shù)據(jù)，建立無煙煤大分子結(jié)構(gòu)模型;③ 模擬無煙煤大分子結(jié)構(gòu)模型吸附甲烷分子;④ 計算分析含甲烷無煙煤構(gòu)型的力學(xué)性質(zhì)。

1.1 無煙煤化學(xué)結(jié)構(gòu)表征

目前，常采用工業(yè)分析、元素分析、核磁共振碳譜(13C-NMR)和X射線光電子能譜(XPS)等手段來表征煤體結(jié)構(gòu)，如馬汝嘉等[11]表征了陜西鳳縣高煤級煤的化學(xué)結(jié)構(gòu)。筆者也采用相似的方法對沁水盆地南部無煙煤進(jìn)行表征。

表1為無煙煤工業(yè)分析和元素分析的測試結(jié)果。測試依據(jù)為《煤的工業(yè)分析方法》(GB/T212—2008)和《煤的元素分析方法》(GB/T31391—2015)進(jìn)行無煙煤樣品的工業(yè)分析(以空氣干燥基為標(biāo)準(zhǔn))和元素分析(以干燥無灰基為標(biāo)準(zhǔn))。

表1 無煙煤工業(yè)分析和元素分析Table 1 Proximate and ultimate analysis of anthracite sample %

圖1 無煙煤樣的13C-NMR頻譜分峰擬合Fig.1 13C-NMR spectral peak fitting of anthracite sample

表2 13C-NMR 化學(xué)位移歸屬Table 2 13C-NMR chemical shift attribution

表3 無煙煤13C-NMR結(jié)構(gòu)參數(shù)計算結(jié)果Table 3 Calculation results of 13C-NMR structural parameters of anthracite %

圖2分別給出了無煙煤XPS總能譜以及C1s,N1s,S2p能譜范圍的分峰擬合曲線。測試儀器為ESCALA B250型的X射線光電子能譜儀。利用Avantage軟件分別對C,N，S元素的能譜進(jìn)行分峰擬合確定O,N，S元素不同存在形式的相對含量。表4給出了O,N，S元素不同化合態(tài)的分布情況。

1.2 無煙煤大分子結(jié)構(gòu)模型的構(gòu)建

根據(jù)上述表征結(jié)果，得出沁水盆地15號無煙煤的分子式為C228H77N3O13S3，共由324個原子組成，其大分子模型的結(jié)構(gòu)參數(shù)和芳香結(jié)構(gòu)見表5，6。構(gòu)建的沁水盆地?zé)o煙煤結(jié)構(gòu)平面模型，如圖3所示。將5個無煙煤的結(jié)構(gòu)平面模型沿三維周期性組合，構(gòu)建無煙煤大分子結(jié)構(gòu)模型，如圖4所示。無煙煤大分子結(jié)構(gòu)模型的三維尺寸為2.74 nm ×2.74 nm ×2.74 nm。

圖2 無煙煤XPS總能譜以及C1s,N1s,S2pFig.2 XPS spectrum and fitting curves of C1s,N1s,S2p of anthracite

表4 無煙煤雜原子分布特征Table 4 Distribution characteristics of hybrid atoms in anthracite

表5 無煙煤大分子模型的結(jié)構(gòu)參數(shù)Table 5 Parameter of molecular model unit of anthracite

表6 無煙煤大分子模型的芳香結(jié)構(gòu)
Table 6 Aromatic structure of macromolecular model of anthracite

圖3 沁水盆地15號無煙煤的結(jié)構(gòu)平面模型Fig.3 Structural plane model of No.15 Anthracite in Qinshui Basin

圖4 沁水盆地?zé)o煙煤的大分子結(jié)構(gòu)模型Fig.4 Macromolecular structure model of No.15 anthracite in Qinshui Basin

1.3 無煙煤吸附甲烷模擬方法和方案

目前，許多研究表明蒙特卡洛(GCMC)分子模擬方法是研究氣體在多孔介質(zhì)中性質(zhì)的有效手段[14,15]。因此，筆者通過Materials Studio軟件中的Sorption模塊，采用蒙特卡洛(GCMC)法研究甲烷在無煙煤中的吸附特性[16]。模擬力場為COMPASS，方法為Metropolis，Quality為Customized，電子勢和范德華勢分別采用Ewald和Atom based方法進(jìn)行統(tǒng)計處理，模擬1×107步，其中前5×106步用于使體系達(dá)到平衡，后5×106步用于計算吸附量和吸附熱等熱力學(xué)參數(shù)。

考慮目前煤層氣資源量以埋深2 000 m為臨界值，其原巖應(yīng)力近似于30 MPa。因此，模擬計算時，設(shè)置壓力為1,5,10,20,30 MPa，設(shè)置溫度為50 ℃。在Materials Studio軟件中需要輸入逸度而非壓力，通過Peng-Robinsion方程實現(xiàn)壓力和逸度的轉(zhuǎn)換[17]。

1.4 含甲烷無煙煤力學(xué)性質(zhì)模擬方法

對固體材料而言，通常利用楊氏模量、體積模量、剪切模量和泊松比等參數(shù)來表征其力學(xué)性質(zhì)。為更好的研究無煙煤的含氣力學(xué)性質(zhì)，基于上述蒙特卡洛法得到1,5,10,20,30 MPa下甲烷在無煙煤中的吸附構(gòu)型，采用分子動力學(xué)的方法，在Materials Studio軟件中Forcite模塊的Mechanical Properties程序中以恒定應(yīng)力的方式對構(gòu)型進(jìn)行單軸拉伸和純剪切變形計算，其中最大應(yīng)變?yōu)?.003，壓強(qiáng)為1個大氣壓(0.1 MPa)，迭代次數(shù)為10 000。在體系始終保持穩(wěn)定的前提下，材料的力學(xué)參數(shù)[18]可表示為

(1)

(2)

G=μ

(3)

(4)

式中，Y為彈性模量，GPa;λ和μ分別為Lame常數(shù)，GPa;K為體積模量，GPa;G為剪切模量，GPa;ν為泊松比。

2 甲烷在無煙煤中的吸附特性

2.1 模型的正確性驗證

為驗證模型的正確性，普遍將模擬數(shù)據(jù)和實驗數(shù)據(jù)進(jìn)行對比。圖5為無煙煤平面模型譜和固體核磁譜的對比?？梢钥闯?，無煙煤平面模型的核磁譜與上述實驗的固體核磁譜具有較好的一致性，對模型進(jìn)行反演計算得出無煙煤的元素含量，與之前實驗煤樣進(jìn)行元素分析時，含量的結(jié)果非常接近，與理論測定無煙煤的元素含量范圍[19]相符，詳見表7。采用分子動力學(xué)模擬得到的無煙煤模型密度為1.526 g/cm3，符合前人通過自然伽馬測井測定的密度(1.39～1.76 g/cm3)[20]，也滿足前人理論計算得到的無煙煤密度(1.3～1.7 g/cm3)[21]。根據(jù)理想氣體定律，用TALU和MYERS[22]的方法模擬計算了無煙煤的孔體積:

圖5 無煙煤平面模型譜和固體核磁譜的對比Fig.5 Comparison of plane model spectrum and solid state nuclear magnetic resonance spectrum of anthracite

其中，R為氣體常數(shù);Nm為單位質(zhì)量吸附劑(Mm)的吸附探針氦氣的分子數(shù);T為溫度;P為壓力。模擬得到的無煙煤大分子結(jié)構(gòu)模型的孔體積為0.018 6 cm3/g，與戚靈靈等[23]采用液氮吸附測得寺河礦煤的孔體積(0.016 4 cm3/g)相差不大。因此，通過上述分析說明本文所建立的無煙煤平面模型和大分子結(jié)構(gòu)模型是合理的，能夠保證后續(xù)計算的正確性。

表7 模擬計算、實驗測定和理論范圍的無煙煤元素質(zhì)量分?jǐn)?shù)、密度和孔體積Table 7 Element content,density and pore volume of anthracite in simulated calculation,experimental measurement and theoretical range

2.2 吸附量和吸附構(gòu)型

圖6為甲烷在無煙煤中的吸附量隨壓力的變化曲線?？梢钥闯?，甲烷在無煙煤中吸附量呈先增大后趨于穩(wěn)定的變化規(guī)律，符合Langmuir模型[24]。壓力PL為1,5,10,20和30 MPa時，甲烷在無煙煤中的吸附量VL分別為9.64,15.64,17.68,19.67,21.03個/晶胞，呈先增大后趨于穩(wěn)定的變化規(guī)律，符合Langmuir模型。通過Langmuir模型對數(shù)據(jù)進(jìn)行擬合，得出甲烷在無煙煤中的飽和吸附量為22.4個/晶胞，Langmuir壓力為1.12 MPa。圖7為20 MPa時甲烷在無煙煤中的吸附構(gòu)型，其他壓力條件也具有相似的構(gòu)型，但是吸附甲烷數(shù)量不同。

圖6 甲烷在無煙煤中的吸附量隨壓力的變化曲線Fig.6 Curve of methane adsorption in coal with pressure

圖7 20 MPa時甲烷在無煙煤中的吸附構(gòu)型Fig.7 Adsorption configuration of methane in anthracite at 20 MPa

2.3 吸附位

徑向分布函數(shù)(RDFs)是描述原子密度與原子間距離的函數(shù)，可以用來反映煤等多孔介質(zhì)材料中的不同原子與甲烷等氣體分子的相互作用強(qiáng)度，判斷其吸附位[25]表達(dá)式為

(5)

式中，gij(r)為在分子調(diào)圍距離r處，分子j的局部密度;dN為從r到r+dr的分子j的數(shù)量;ρj為分子j的密度。

圖8為壓力20 MPa時，無煙煤模型中的碳、氫、氧、氮和硫原子與甲烷分子的徑向分布函數(shù)曲線，可以看出第1個RDFs波峰出現(xiàn)在距離為0.11 nm處，為CH4—C的波峰，峰值為5.78。其次是在距離為0.238 nm和0.43 nm處分別由CH4—H和CH4—N產(chǎn)生的波峰，CH4—O和CH4—S引起的波峰距離較遠(yuǎn)，波峰相對較低。因此，芳香碳、吡啶型氮和吡咯型氮以及羧基是甲烷分子主要吸附位。甲烷分子更容易吸附在由芳香環(huán)和脂肪結(jié)構(gòu)組成的孔隙中，且更容易被孔隙中靠近芳香結(jié)構(gòu)的一側(cè)吸附。

圖8 無煙煤與甲烷分子的徑向分布函數(shù)曲線Fig.8 Radial distribution function of anthracite and methane molecule

2.4 等量吸附熱

在甲烷吸附過程中，等量吸附熱能夠反映無煙煤與甲烷間能量釋放的變化關(guān)系，可通過Clausius-Clapeyron方程計算求得[26]。圖9為甲烷在無煙煤中的等量吸附熱隨吸附壓力的變化規(guī)律?？梢钥闯?，甲烷在無煙煤中的等溫吸附熱qst隨吸附壓力P的升高呈對數(shù)規(guī)律下降，擬合得到的對數(shù)表達(dá)式為qst=-0.71× lnP+25.05，表明甲烷在無煙煤中吸附過程中釋放的熱量隨著壓力的升高逐漸減小，并最終達(dá)到平衡，這是因為甲烷在低壓力時率先占據(jù)無煙煤表面的高能吸附位。在吸附壓力范圍內(nèi)，甲烷在煤中的吸附熱為22.65～25.00 kJ/mol，遠(yuǎn)小于42 kJ/mol，屬于物理吸附。

圖9 甲烷在煤中的等量吸附熱隨壓力的變化曲線Fig.9 Curve of isosteric heat of methane adsorption in coal with pressure

3 無煙煤的含氣力學(xué)性質(zhì)

3.1 力學(xué)性質(zhì)

表8為不同含氣量無煙煤的力學(xué)參數(shù)。圖10為無煙煤的力學(xué)參數(shù)隨吸附量Q的變化規(guī)律。由圖10可以看出，體積模量K、彈性模量E和剪切模量G具有一致的變化規(guī)律，均隨著甲烷吸附量的增大呈對數(shù)規(guī)律降低，而泊松比ν隨吸附量的增大呈線性規(guī)律增大，擬合方程分別為

K=ln(881.04-36.65Q)

(6)

E=ln(234.15-8.07Q)

(7)

G=ln(7.71-0.14Q)

(8)

ν=ln(0.31+0.002Q)

(9)

王晉等[27]研究了不同含氣飽和度下煤巖力學(xué)性質(zhì)的變化趨勢，對其數(shù)據(jù)進(jìn)行擬合(圖11)，發(fā)現(xiàn)煤的彈性模量與含氣飽和度呈對數(shù)衰減規(guī)律，與本文的研究結(jié)論相一致。當(dāng)甲烷含氣量為21.03個/晶胞時，無煙煤的體積模量、彈性模量和剪切模量分別為4.12,4.29和1.45 GPa，同比于不含甲烷的無煙煤，體積模量、楊氏模量和剪切模量分別降低了38.5%,24.4%和27.1%。這表明吸附甲烷能夠顯著降低無煙煤的力學(xué)強(qiáng)度，且吸附量越大影響程度越高。

圖10 無煙煤的力學(xué)參數(shù)隨吸附量的變化曲線Fig.10 Curves of the mechanical parameters of coal with the adsorption capacity

圖11 煤的彈性模量隨含氣飽和度的變化曲線Fig.11 Variation curve of elastic modulus of coal with gas saturation

3.2 結(jié)果討論與機(jī)理分析

高保彬等[28]采用自主研發(fā)的含瓦斯煤巖受載變形破壞試驗裝置，通過實驗室試驗研究了不同瓦斯壓力條件下無煙煤的力學(xué)性質(zhì)，設(shè)置瓦斯壓力為0～2 MPa，發(fā)現(xiàn)煤的力學(xué)強(qiáng)度隨瓦斯壓力的增大而降低，與本文的結(jié)論相一致。但由于其所設(shè)置壓力范圍偏小，以此結(jié)果進(jìn)行無煙煤力學(xué)性質(zhì)評價有一定的局限性。本試驗延伸了瓦斯壓力范圍(20 MPa)，定量分析了無煙煤力學(xué)性質(zhì)和吸附量之間的關(guān)系，得出了無煙煤的體積模量、楊氏模量和剪切模量隨著吸附量的變化規(guī)律。

已有研究表明，多孔介質(zhì)材料吸附甲烷后會發(fā)生膨脹變形，其計算公式[16]為

a=(Vi-V0)/V0

(10)

其中，a為膨脹率;Vi，V0分別為煤的原始體積和吸附甲烷后的體積。圖12為無煙煤的體積和膨脹率隨含氣量的變化規(guī)律。由圖12可知，無煙煤的體積和膨脹率隨吸附量的增加呈指數(shù)規(guī)律增大。與未吸附甲烷的無煙煤相比，當(dāng)甲烷含氣量為21.03個/晶胞時，無煙煤的膨脹率達(dá)到了9.012%。根據(jù)周動等[29]的研究，發(fā)現(xiàn)無煙煤的膨脹變形量隨煤結(jié)構(gòu)尺度的減小而顯著增大。當(dāng)煤結(jié)構(gòu)尺度減小時，膨脹變形量的降幅增大，這是由于在較小統(tǒng)計尺度下煤結(jié)構(gòu)的非均質(zhì)性強(qiáng)，結(jié)構(gòu)間力學(xué)強(qiáng)度與膨脹變形能力差異較大，膨脹變形量均較大;而在較大統(tǒng)計尺度下煤結(jié)構(gòu)的非均質(zhì)性弱，結(jié)構(gòu)間力學(xué)強(qiáng)度與膨脹變形能力差異較小，膨脹較小。煤體的膨脹變形主要受到吸附壓力與結(jié)構(gòu)尺度的雙重影響。ZHANG等[30]采用分子模擬的方法研究了含水煙煤吸附甲烷后的膨脹變形，發(fā)現(xiàn)當(dāng)含水率為3.0%和壓力為1 MPa時，煙煤的膨脹率達(dá)到3.1%。本文的研究是不考慮含水狀態(tài)下的含甲烷力學(xué)性質(zhì)，甲烷的吸附位沒有被水占據(jù)，能夠吸附更多的甲烷，產(chǎn)生更大的膨脹變形。因此，本文在納米尺度條件下，無煙煤吸附甲烷的膨脹率是可以接受的。

表9為不同含氣量無煙煤模型的能量對比。可以看出，隨著含氣量的增大，煤體會發(fā)生吸附膨脹，導(dǎo)致無煙煤模型的總能量降低，說明無煙煤基質(zhì)之間的相互作用力降低，進(jìn)而導(dǎo)致無煙煤的強(qiáng)度降低，抵抗變形的能力減弱。表9結(jié)果表明，雖然成鍵能與非成鍵能相比降低幅度小，但鍵能、鍵角能、扭轉(zhuǎn)能、反轉(zhuǎn)能等均有不同程度的降低，這說明甲烷的吸附將導(dǎo)致鍵能的變化，使其鍵強(qiáng)降低。分子間作用力是物質(zhì)穩(wěn)定的重要物理因素，表中數(shù)據(jù)表明范德華能在吸附甲烷的過程中能量降低值最大，對無煙煤總能量的降低起主導(dǎo)作用，表明范德華能是保持煤體結(jié)構(gòu)和力學(xué)性質(zhì)穩(wěn)定的主導(dǎo)因素。當(dāng)含氣量分別為9.64，15.64，17.68，19.67和21.03個/晶胞時，無煙煤的范德華能分別降低了25.02，27.97，28.99，30.86和32.60 kJ/mol。因此，含甲烷無煙煤體積模量、彈性模量、剪切模量等降低、泊松比增大的內(nèi)在原因是吸附甲烷會降低煤分子的鍵能和分子間相互作用力。這種內(nèi)在關(guān)系值得進(jìn)一步的理論與實驗研究。

圖12 無煙煤的體積和膨脹率隨吸附量的變化曲線Fig.12 Variation of volume and expansion rate of anthracite with gas content

表9 不同含氣量無煙煤模型的能量對比Table 9 Energy comparison of anthracite models with different gas contents kJ/mol

4 結(jié) 論

(1)無煙煤樣大分子結(jié)構(gòu)中主要由芳烴類碳組成，平均縮聚芳香結(jié)構(gòu)由3～4個苯環(huán)的鏈型連接或環(huán)型連接。

(2)甲烷在無煙煤中吸附量呈先增大后趨于穩(wěn)定的變化規(guī)律，符合Langmuir模型，甲烷在無煙煤中的飽和吸附量為22.4個/晶胞，Langmuir壓力為1.12 MPa。甲烷分子更容易吸附在由芳香環(huán)和脂肪結(jié)構(gòu)組成的孔隙中，且更容易被孔隙中靠近芳香結(jié)構(gòu)的一側(cè)吸附。

(3)無煙煤的體積模量、彈性模量和剪切模量均隨著甲烷吸附量的增大呈對數(shù)規(guī)律降低，而泊松比隨吸附量的增大呈線性增大。

(4)隨著含氣量的增大，煤體會發(fā)生吸附膨脹，無煙煤模型的總能量降低，表明無煙煤基質(zhì)之間的相互作用力降低，進(jìn)而導(dǎo)致無煙煤的強(qiáng)度降低，抵抗變形的能力減弱。范德華能是保持煤體結(jié)構(gòu)和力學(xué)性質(zhì)穩(wěn)定的主導(dǎo)因素。

(5)分子模擬是研究煤等多孔介質(zhì)材料在含氣狀態(tài)下力學(xué)性質(zhì)的一種有效手段。