亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        煤分子結構對煤層氣吸附與擴散行為的影響

        2016-07-07 12:13:10王寶俊章麗娜凌麗霞章日光太原理工大學煤科學與技術教育部和山西省重點實驗室山西太原03004太原理工大學化學化工學院山西太原03004
        化工學報 2016年6期
        關鍵詞:擴散吸附

        王寶俊,章麗娜,凌麗霞,章日光(太原理工大學煤科學與技術教育部和山西省重點實驗室,山西 太原 03004;太原理工大學化學化工學院,山西 太原 03004)

        ?

        煤分子結構對煤層氣吸附與擴散行為的影響

        王寶俊1,章麗娜1,凌麗霞2,章日光1
        (1太原理工大學煤科學與技術教育部和山西省重點實驗室,山西 太原 030024;2太原理工大學化學化工學院,山西 太原 030024)

        摘要:對煤的芳香單元延展度、芳香單元堆砌層數(shù)以及包括不同缺陷和含氧官能團類型的表面結構對煤層氣吸附與擴散的影響進行了研究。采用Monte Carlo模擬方法及分子動力學模擬方法分別得到了煤層氣的吸附量與擴散系數(shù),模擬溫度為303 K,壓力為10 MPa。研究結果表明,單位質量的煤對甲烷的吸附量隨著芳香單元堆砌層數(shù)的增加而降低,缺陷和含氧官能團的存在不利于甲烷的吸附。甲烷的擴散隨著芳香單元延展度的增加呈現(xiàn)出一個N形的復雜變化過程,單缺陷和羰基的存在有利于煤層氣的擴散,煤結構中大的裂隙更有利于煤層氣的擴散。最后,基于煤層氣的微觀影響因素和宏觀運移行為提出了煤顆粒徑向不均質的煤層氣擴散微觀模型。

        關鍵詞:煤分子結構;煤層甲烷;吸附;擴散;分子模擬

        2015-11-26收到初稿,2016-01-07收到修改稿。

        聯(lián)系人及第一作者:王寶?。?964—),男,博士,教授。

        Received date: 2015-11-26.

        Foundation item: supported by the Major Projects of National Science and Technology (2011ZX05040-005-002-001), the National Natural Science Foundation of China (21276171, 21576178, 21276003) and the Program for the Innovative Talents of Higher Learning Institutions of Shanxi.

        引 言

        煤與瓦斯突出現(xiàn)象是大量的游離態(tài)煤層氣在煤體強度較弱面瞬間涌出的現(xiàn)象,這與地應力、氣體壓力及煤體本身的物理機械性能有關[1]。構造煤是指煤層受構造應力作用,原生結構、構造受到強裂破壞而產(chǎn)生碎裂、揉皺、擦光面等構造變動特征的煤。而原生結構煤則是指未受構造應力作用的煤。因此構造煤具有不同于原生結構煤的特殊的物理性質,如強度低、滲透率低、吸附量高和擴散系數(shù)大,使得構造煤是煤與瓦斯突出的高發(fā)生煤體[2]。

        關于構造煤的吸附與擴散行為已有大量的文獻報道。一致認為,構造煤的吸附量明顯高于原生結構煤[3-4],且Pan等[5]通過研究不同變形程度煤中甲烷的吸附量得到不同變形程度的煤體表現(xiàn)出不同的吸附等溫線和吸附量。通過對比構造煤和原生結構煤在相同實驗條件下的完全熱解生烴潛力,發(fā)現(xiàn)構造煤的熱解生烴潛力低,由此推斷,構造應力的作用使得構造煤在形成過程中產(chǎn)生了大量的烴類[6-7],故認為構造煤較相應的原生結構煤具有超前演化的特征。

        與原生結構煤相比,構造煤層中煤層氣在初期解吸速率和瞬間解吸總量上都有顯著提高[8],其中,初期解吸速率是原生結構煤的1.36~2.84倍,軟煤的解吸初速率明顯大于硬煤[9];構造煤在前10 min內的煤層氣解吸量已經(jīng)達到1 h解吸量的60%[8],解吸總量隨著煤樣粒度的減小而增加。

        構造煤宏觀結構變化的最根本因素是其微觀結構發(fā)生了變化,量子化學計算、Monte Carlo模擬和分子動力學模擬是從微觀上解決煤層氣吸附與擴散問題的有效工具。降文萍等[10]通過用不同數(shù)量的苯環(huán)來代表不同變質程度的煤,計算了不同變質程度的煤與甲烷氣體的相互作用,結果表明煤對甲烷的吸附屬于物理吸附,煤表面與甲烷分子的吸附能隨著煤變質程度的升高而增加。Qiu等[11]用3種不同的簇結構代表煤模型研究了甲烷在不同吸附構型及不同吸附位上的吸附,研究結果顯示六元環(huán)的中心位是最穩(wěn)定的吸附位,甲烷以正三角錐形吸附在表面時的吸附能最大,結構最穩(wěn)定,這與陳昌國等[12]的研究結果一致。而Liu等[13]采用石墨超胞來研究煤與甲烷分子的相互作用時,發(fā)現(xiàn)甲烷主要以正三角錐吸附在C C橋位上。近年來,分子動力學模擬逐漸成為解決煤體復雜物理行為的重要途徑,Hu等[14]采用分子模擬的方法研究了CO2和CH4在Wiser煤上的吸附與擴散。

        構造煤的吸附與擴散行為是由其大分子結構的變化引起的[15-16],相比于原生結構煤,構造煤在分子結構上有明顯的特征:石墨化程度較高,基本結構單元有序性增強,芳香疇環(huán)系統(tǒng)有序增加[7,17]。在脆性變形煤中芳香單元堆砌度(Lc)以及芳香單元堆砌層數(shù)(N)先減小后增加,而在韌性變形煤中則呈現(xiàn)相反的趨勢[18]。同時,在構造應力作用下,由于構造煤相比于同一時期的原生結構煤具有超前演化的特征,使得煤中含氧官能團(羥基、羧基、羰基和醚)的主要存在形式也發(fā)生了變化。其中,醚鍵在構造煤中幾乎消失[19],而羧基和羰基卻基本沒有變化。因此,亟需探討煤分子結構對煤層氣運移的影響。據(jù)此,有必要構建一系列可以代表煤某一方面特征的煤分子模型,用于研究導致構造煤特殊性質甚至引起煤與瓦斯突出的煤微晶結構參數(shù)和表面結構。

        根據(jù)實驗獲得的煤分子微晶結構參數(shù),以芳香單元延展度(La)、芳香單元堆砌層數(shù)(N)及表面缺陷和含氧官能團類型為變量,構建一系列可以代表煤某一方面特征的煤分子模型,采用Monte Carlo模擬方法和分子動力學模擬方法分別研究了煤層氣的吸附與擴散,通過分析煤分子結構與煤層氣運移的關系,剖析影響煤與瓦斯突出的分子結構,得到了煤層氣擴散微觀模型。

        1 模型構建與模擬參數(shù)

        1.1模型構建

        應力降解和應力縮聚是構造應力影響煤化作用的兩個基本過程。應力降解過程使得煤基本結構單元中側鏈和官能團發(fā)生斷裂,導致煤分子的化學組成發(fā)生變化。應力縮聚促使煤基本結構單元重排、有序疇增大、定向生長,最終導致芳香疇環(huán)體系加大[17,20]。因此,在構造煤形成的過程中芳香單元延展度(La),芳香單元堆砌度(Lc),芳香單元延展度與芳香單元堆砌度的比值(La/Lc),芳香單元堆砌層數(shù)(N)以及芳香環(huán)的個數(shù)(n)都發(fā)生了相應變化。一些研究者通過實驗方法得到了多個煤分子的微晶結構參數(shù)[16,18,21-23],每項研究的目的不同,有些針對不同的煤種,有些針對不同的變形程度?;谶@些結果得到了各個微晶結構參數(shù)的范圍,列于表1,基本涵蓋了煤微晶結構參數(shù)的變化范圍。

        表1 煤微晶結構參數(shù)的范圍Table 1 Microcrystal structure parameter ranges of coal

        另外,構造煤中缺陷的存在也是造成其物理行為變化的關鍵因素[21],在煤的演化過程中,含氧官能團的主要存在形式也發(fā)生著變化。因此,有必要構建含有不同缺陷類型和含氧官能團類型的煤結構模型,以獲得造成構造煤煤層氣吸附量高、擴散快的主要因素。

        圖1 煤微晶結構模型Fig.1 Models of microcrystal structures of coal

        根據(jù)表1中的數(shù)據(jù),獲得了一系列以煤的某一方面特征為變量的煤分子模型。模型采用N-La的命名方式,如3-5×5表示芳香單元堆砌層數(shù)為3,芳香單元延展度為5×5的煤分子模型。其中,以芳香單元延展度La為變量的有5個模型;以芳香單元堆砌層數(shù)N為變量的有7個模型;以表面所含缺陷類型為變量的有2個模型,其中V1表示單缺陷、V2表示次生缺陷;以表面所含含氧官能團類型為變量的有4個模型。由于3-9×9的表面可以很好地排除邊界效應的影響,因此各類模型的構建以煤分子模型3-9×9為基礎,進而分別研究了不同缺陷和含氧官能團類型的影響。所構建的煤分子模型共5+ 7? 1+ 2 + 4 = 17個。部分煤分子模型的結構見圖1,各模型的微晶結構參數(shù)列于表2中。

        1.2吸附模擬

        采用Monte Carlo統(tǒng)計方法[24]對煤層氣中甲烷分子的吸附行為進行研究,采用的模擬軟件為Materials Studio 5.5。甲烷分子及各個模型的幾何構型首先通過Forcite模塊中的Geometry Optimization任務完成,使體系能量達到最小,同時表面電荷呈電中性,得到其穩(wěn)定結構。再啟動Sorption模塊下的Fix pressure任務模擬煤分子模型對甲烷氣體分子的吸附量,得到固定壓力下甲烷分子在煤分子模型中的飽和吸附量。

        表2 所構建煤模型的微晶結構參數(shù)Table 2 Microcrystal structure parameters of established coal models

        整個模擬過程采用COMPASS力場[25-26],方法為Metropolis,Quality采用Medium,電子勢和范德華勢分別采用Ewald[27]和Atom based方法進行統(tǒng)計處理。模擬過程中前10000步用于使體系達到平衡,后100000步用來統(tǒng)計所需的熱力學性質。

        1.3擴散模擬

        關于甲烷單組分氣體分子在不同煤分子模型表面的擴散行為,采用Materials Studio 5.5軟件中的Discover模塊模擬。模擬過程是在周期性條件下進行的,采用COMPASS力場計算相互作用,各分子的起始速度由Maxwell-Boltzmann分布隨機產(chǎn)生。分子動力學模擬選用NVT系綜(恒定粒子數(shù)、恒定體積、恒定壓力),壓力為10 MPa,溫度為303 K。模擬時采用的控溫方法為Andersen[28]。總模擬時間為2000 ps,時間步長為1.0 fs,每1000步記錄一次體系的軌跡信息(幀,frame),便于進行下一步分析。

        完成分子動力學模擬后,分析原子軌跡可以獲得均方位移(mean square displacement,MSD),進一步則可以計算出氣體分子的擴散系數(shù)(D)。由于受結構外界環(huán)境的影響,在愛因斯坦擴散之前常常發(fā)生異常擴散(非愛因斯坦擴散)[29],而這一部分不用于分析擴散行為。擴散系數(shù)主要描述的是擴散粒子的隨機走動,通過愛因斯坦擴散公式[30]求得

        式中,Ri(t)為第i個分子或粒子在t時刻的位置;Ri(0)為第i個分子或粒子在初始時刻的位置。

        因此,有必要在計算擴散系數(shù)之前驗證此擴散是否符合愛因斯坦擴散[31]。確認的方法為:作關于MSD與t的對數(shù)曲線,如果為直線則表示模擬過程已經(jīng)收斂。此外,對數(shù)曲線的斜率m的大小可以作為判斷愛因斯坦擴散發(fā)生的依據(jù),當m<1時表示發(fā)生的是異常擴散,當m= 1時則表示發(fā)生了愛因斯坦擴散[32]。

        2 結果與討論

        2.1煤層氣的吸附行為

        甲烷依靠弱相互作用吸附在煤表面[33]。采用10 MPa的壓力,303 K的溫度來模擬當前礦井條件[3],根據(jù)Peng-Robinson方程[34]得到甲烷在上述壓力和溫度下的逸度為8287.1 kPa。模型首先經(jīng)過優(yōu)化得到最穩(wěn)定的構型,然后通過固定壓力下的模擬得到給定壓力下模型所能吸附的甲烷的最大量。

        吸附模擬可以得到模型所能吸附氣體分子的最大量,即絕對吸附量,它包含了吸附在多孔材料孔壁上以吸附相形式存在的氣體分子,以及存在于多孔材料孔內以氣相形式存在的氣體分子。而實際中用于比較的是超額吸附量,即只包含吸附在多孔材料孔壁上以吸附相形式存在的氣體分子。通過式(2)[35]可以將絕對吸附量(Nabs)轉化為超額吸附量(Nexc)

        式中,Nexc為超額吸附量;mcell為每個超晶胞的質量;Nabs為絕對吸附量,即每個晶胞中吸附的氣體分子數(shù);Vp為分子模型的空隙體積。煤分子結構的空隙體積可以通過Materials Studio中的Connolly[36]任務獲得,在此過程中Connolly半徑設為1 ?[37](1 ?= 10?10m);甲烷氣體在10 MPa、303K下的密度為75798.5 g·m?3。

        2.1.1芳香單元延展度對煤層氣吸附量的影響通過Connolly任務和Fix pressure任務分別得到了不同芳香單元延展度下煤的空隙體積(Vp)和煤層氣的絕對吸附量(Nabs),結果列于表3中。

        表3 不同芳香單元延展度下煤分子的空隙體積及對甲烷的吸附量Table 3 Void volumes of coal in different Laand its adsorption amount of methane

        可以看出,隨著芳香單元延展度的增大,煤分子結構的空隙體積由3184.70 ?3增加到7731.10 ?3,絕對吸附量由26.67 u.c.?1增加到70.01 u.c.?1,超額吸附量幾乎沒有變化,在147.90~152.33 cm3·g?1之間浮動。這是由于絕對吸附量表示的是單個晶胞中吸附甲烷氣體分子的數(shù)量,其數(shù)值與煤分子表面的吸附位點密切相關,甲烷氣體分子單層吸附在煤分子表面時,芳香單元延展度越大,所提供的吸附位點則越多;而超額吸附量表示的是單位質量的煤吸附甲烷氣體的體積,其數(shù)值不僅與表面的吸附位點有關,還與煤的分子質量及空隙體積有關。因此延展度由5×5增加到9×9時,由于吸附位點的增加,使得煤分子對甲烷氣體分子的絕對吸附量增加,同時由于單位質量的煤吸附位點不變,造成煤對甲烷氣體分子的超額吸附量基本不變。因此可以得出,煤分子結構中的芳香單元延展度并不是影響甲烷氣體分子實際吸附量的因素。

        2.1.2芳香單元堆砌層數(shù)對煤層氣吸附量的影響同樣地,對不同芳香單元堆砌層數(shù)的煤分子模型也進行了定壓吸附甲烷氣體分子的模擬。不同芳香單元堆砌層數(shù)下,煤分子的空隙體積及對甲烷氣體分子的絕對吸附量和超額吸附量見表4。

        表4 不同芳香單元堆砌層數(shù)下煤分子的空隙體積及對甲烷的吸附量Table 4 Void volumes of coal in different N and its adsorption amount of methane

        從表可見,芳香單元堆砌層數(shù)由3變化到9時,煤分子的空隙體積略微增加,整體變化不大,在3184.70~3411.73 ?3之間浮動。煤分子模型堆砌層數(shù)不同時,受延展度影響的表面吸附位點數(shù)沒有發(fā)生變化,因此對甲烷氣體分子的絕對吸附量基本保持不變,浮動范圍為23.98~26.74 u.c.?1。但增加芳香單元堆砌層數(shù),會導致煤分子的分子質量增加,當模型表面吸附的甲烷氣體分子數(shù)基本相當時(絕對吸附量相當),分子質量的增加會使得單位質量的煤對甲烷氣體的吸附體積減少,即超額吸附量降低。如表4所示,在不同芳香單元堆砌層數(shù)的煤分子中,煤對甲烷氣體分子的超額吸附量由N =3時的152.33 cm3·g?1降低到N =9時的41.32 cm3·g?1,且在堆砌層數(shù)較低時降低趨勢更為明顯。

        2.1.3缺陷及含氧官能團類型對煤層氣吸附量的影響隨后,為了研究煤分子表面上存在不同缺陷及含氧官能團對甲烷氣體分子吸附量的影響,對含有不同缺陷和含氧官能團的煤分子進行了定壓吸附甲烷氣體分子的模擬。不同表面結構下煤分子的空隙體積及對甲烷氣體分子的絕對吸附量和超額吸附量見表5。

        表5 不同表面結構下煤分子的空隙體積及對甲烷的吸附量Table 5 Void volumes of coal in different surface structure and its adsorption amount of methane

        缺陷和含氧官能團的存在對煤分子空隙體積的影響不大,單缺陷、次生缺陷下煤分子的空隙體積分別為7777.89和7755.33 ?3,其值略高于煤分子模型3-9×9(7731.10 ?3),含氧官能團存在下煤分子的空隙體積在7661.76~7718.98 ?3范圍內,整體略低于煤分子模型3-9×9。

        對比絕對吸附量可以發(fā)現(xiàn),缺陷和含氧官能團的存在會使單個晶胞中吸附的甲烷分子數(shù)減少,這是由于表面對甲烷氣體分子的吸附位點減少了。同樣,煤分子對甲烷的超額吸附量也并沒有由于表面結構的變化而發(fā)生很大的變化,但仍可以看出,單位質量的煤對甲烷氣體的吸附體積有略微降低。其中,煤分子模型3-9×9中單位質量的煤對甲烷氣體的飽和吸附量為148.35 cm3·g?1,模型3-9×9-V1 和3-9×9-V2對甲烷分子的超額吸附量分別為141.46和144.70 cm3·g?1;模型3-9×9-OH、3-9× 9-COOH、3-9×9-OCH3和3-9×9-C O對甲烷分子的超額吸附量分別為141.51、139.66、138.97和140.13 cm3·g?1。

        因此,當有缺陷和官能團存在時會降低煤層對甲烷分子的吸附量,這與Hao等[38]的觀點一致,表明煤的石墨化程度越高,其對煤層氣的吸附量也就越大。而構造煤在經(jīng)過一系列的應力降解和應力縮聚作用后,會大大降低結構中官能團及側鏈的含量[17],提高石墨化程度,導致石墨化程度高于同一時期原生結構煤的構造煤往往含有較大的煤層氣吸附量,容易發(fā)生煤與瓦斯突出現(xiàn)象。

        2.1.4影響煤層氣吸附量的分子結構將煤分子結構中的芳香單元延展度、芳香單元堆砌層數(shù)及表面含有的缺陷和含氧官能團類型作為變量,構建了17個不同的煤分子模型。通過不同煤分子結構下甲烷氣體分子在表面的吸附情況,可以得到影響煤層氣宏觀吸附量的微觀結構。其中,芳香單元延展度通過改變煤表面的吸附位點來影響煤層氣的絕對吸附量;芳香單元堆砌層數(shù)通過改變煤的分子質量來影響煤層氣的超額吸附量??傮w來說,芳香單元延展度并不會影響單位質量的煤對煤層氣的吸附體積,而隨著芳香單元堆砌層數(shù)的增加及表面缺陷和官能團的存在,單位質量的煤對煤層氣的吸附體積呈逐漸降低的趨勢,且在較低層數(shù)下,煤分子中芳香單元堆砌層數(shù)的變化對吸附量的影響較大。

        然而,構造應力一方面使構造煤分子結構的聚合、縮合作用增強,導致芳香單元延展度和堆砌度增加[7,17],另一方面造成煤層氣吸附量增加[3-4],這與本節(jié)模擬結果不符。因此可以判斷:煤分子結構并不是影響煤層氣吸附量高低的主要因素,造成構造煤煤層氣高吸附量的主要因素應該是納米級孔結構的變化。

        2.2煤層氣的擴散行為

        擴散系數(shù)是基于統(tǒng)計原理獲得的,為了提高精度,模擬過程中必須保證有一定數(shù)量的樣本數(shù)才能避免明顯的漲落現(xiàn)象,這就需要將原來的單胞擴大。因此,為了確定分子動力學模擬的初始模型,首先要將飽和吸附甲烷分子后的煤分子模型的單胞進行不同程度的擴大。選取的煤分子模型分別是分子結構中芳香單元堆砌層數(shù)和芳香單元延展度最小的模型3-5×5,芳香單元堆砌層數(shù)最大的模型9-5×5及芳香單元延展度最大的模型3-9×9。3個結構的單胞分別擴大到1×1×1、2×1×1和2×2×1,得到的9個結構含有不同的樣本數(shù),分別對其進行了20次分子動力學模擬。表6列出了各個單胞中樣本數(shù)與擴散系數(shù)的合格數(shù)之間的關系。

        表6 不同樣本數(shù)下擴散系數(shù)的合格數(shù)Table 6 Qualified number of diffusion coefficient under different number of samples

        表中的數(shù)據(jù)顯示當樣本數(shù)接近90時,數(shù)據(jù)可重復性較為可靠。芳香單元延展度最小、單胞吸附甲烷分子數(shù)最少的模型x-5×5滿足樣本數(shù)大于90時,需要使單胞擴大到2×2×1。因此為了保證所有模型進行動力學模擬的樣本數(shù)大于90,選擇飽和吸附甲烷分子后、單胞為2×2×1的結構作為擴散模擬的初始構型。

        將達到飽和吸附后的煤結構進行分子動力學模擬,模擬條件為303 K、10 MPa,得到了以不同芳香單元延展度、不同芳香單元堆砌層數(shù)及表面不同缺陷和含氧官能團類型為變量的煤結構的擴散系數(shù),見表7。

        表7 不同煤結構模型中甲烷的擴散系數(shù)Table 7 Diffusion coefficients of methane in different models

        2.2.1不同芳香單元延展度下煤層氣的擴散系數(shù)以芳香單元延展度為變量,對結構進行多次分子動力學模擬,獲得擴散系數(shù),取平均值,即得到甲烷氣體分子在各結構中的平均擴散系數(shù)。從數(shù)據(jù)分析可得,芳香單元延展度的增大使得甲烷在煤分子中的擴散系數(shù)增大。在芳香單元延展度為5×5時,甲烷的擴散系數(shù)僅為1.19×10?10m2·s?1,芳香單元延展度增加到6×6時,甲烷的擴散系數(shù)顯著增加,達到1.93×10?10m2·s?1。甲烷的擴散系數(shù)隨著芳香單元延展度的持續(xù)增加而不斷增加,當芳香單元延展度為9×9時,擴散系數(shù)已經(jīng)達到了3.19×10?10m2·s?1,但增加量要比延展度較低時有所降低。構造煤在形成過程中受到構造應力的作用,使得芳香單元延展度增加,有序性增強,而芳香單元延展度的增加恰恰是造成煤層氣擴散系數(shù)增大的關鍵因素,由此可以解釋構造煤中煤層氣擴散速率高于相應原生結構煤的現(xiàn)象[8]。

        2.2.2不同芳香單元堆砌層數(shù)下煤層氣的擴散系數(shù)芳香單元堆砌層數(shù)是構造應力影響下煤分子結構的重要變量。從表7可以看出,芳香單元堆砌層數(shù)增加過程中,甲烷氣體分子在煤分子中的擴散系數(shù)發(fā)生明顯的變化。隨芳香單元堆砌層數(shù)的增加,甲烷氣體分子的擴散系數(shù)表現(xiàn)出先增大后減小再增大的特點,整體呈“N”形趨勢,如圖2所示。當芳香單元堆砌層數(shù)從3增加到5時,甲烷氣體分子在煤分子中的擴散系數(shù)由1.19×10?10m2·s?1增加到 2.77× 10?10m2·s?1;然而隨著芳香單元堆砌層數(shù)繼續(xù)增加,甲烷氣體分子的擴散系數(shù)從模型5-5×5中的2.77×10?10m2·s?1降低到模型7-5×5中的1.75× 10?10m2·s?1;當芳香單元堆砌層數(shù)從7增加到9時,甲烷氣體分子的擴散系數(shù)又呈增加的趨勢,從模型7-5×5中的1.75×10?10m2·s?1增加到模型9-5×5中的2.43×10?10m2·s?1。此外,當芳香單元堆砌層數(shù)較低時,其變化對甲烷氣體擴散行為的影響更加明顯。堆砌層數(shù)從3增加到5,甲烷氣體的擴散系數(shù)增加了1.58×10?10m2·s?1;堆砌層數(shù)從5增加到7,甲烷氣體的擴散系數(shù)降低了1.02×10?10m2·s?1;堆砌層數(shù)從7增加到9,甲烷氣體的擴散系數(shù)增加了0.68×10?10m2·s?1。

        圖2 堆砌層數(shù)變化下甲烷的擴散行為Fig.2 Diffusion behavior of methane with changes of number of basic structural unit

        本節(jié)只模擬了10 MPa、303 K下不同芳香單元堆砌層數(shù)下煤層氣的擴散行為,雖然得出煤層氣的擴散系數(shù)隨芳香單元堆砌層數(shù)的增加呈“N”形趨勢,但由于外界環(huán)境條件的局限性,并不能充分說明煤層氣的擴散系數(shù)隨芳香單元堆砌層數(shù)的增加先增加后減小再增加,而只能說明其不是一個單調的變化過程,而是一個受多種因素共同作用的復雜過程。

        2.2.3表面缺陷及含氧官能團對煤層氣擴散系數(shù)的影響煤分子表面單缺陷的存在對甲烷擴散的影響是最大的,擴散系數(shù)為4.39×10?10m2·s?1,而在近石墨化的煤分子模型3-9×9中,甲烷氣體分子的擴散系數(shù)為3.19×10?10m2·s?1,因此單缺陷的存在使得甲烷在煤分子中的擴散系數(shù)增加了1.20×10?10m2·s?1,有利于甲烷氣體的擴散。然而,煤分子表面次生缺陷位的存在卻降低了甲烷氣體分子的擴散系數(shù),僅為2.70×10?10m2·s?1,故次生缺陷的存在不利于甲烷氣體的擴散。

        分析煤分子表面不同含氧官能團存在時甲烷氣體分子的擴散系數(shù),可以發(fā)現(xiàn),僅當表面有羰基存在的情況下,甲烷的擴散系數(shù)要比在類石墨化的分子模型3-9×9中高,其擴散系數(shù)為3.37×10?10m2·s?1,表面羧基的存在很大程度上降低了甲烷的擴散系數(shù),其擴散系數(shù)為2.39×10?10m2·s?1,這與構造煤中含有較多的羰基,較少的羧基、醚鍵有關[19]。表面存在羥基和甲氧基時,甲烷氣體分子的擴散系數(shù)分別為3.06×10?10和3.12×10?10m2·s?1,與近石墨化的模型3-9×9相差不大。

        2.2.4不同空隙下煤層氣的擴散為了研究宏觀裂隙對煤層氣擴散行為的影響,將煤結構模型3-5×5的真空層設置為不同厚度,分別為13、18和23 ?。這個真空層厚度可以較好模擬芳香層片間隙的大小和宏觀裂隙的大小。模擬過程首先將25個甲烷分子分別吸附在3個不同真空層厚度的煤模型中,再將單胞擴大為2×2×1,保證擴散的樣本數(shù)符合要求,進行多次重復的動力學模擬。不同真空層厚度下煤層氣的擴散系數(shù)見表8。

        表8 不同真空層厚度下甲烷的擴散系數(shù)Table 8 Diffusion coefficients of methane in differentvacuum layer

        真空層厚度越大代表實際情況下煤層氣運移的宏觀裂隙越大。當真空層厚度為13 ?時,甲烷分子的擴散系數(shù)為1.10×10?10m2·s?1;當煤分子結構中的真空層厚度增加到18 ?時,即煤體宏觀裂隙增加,此時甲烷分子的擴散系數(shù)為4.72×10?10m2·s?1,相比于真空層厚度為13 ?時的擴散系數(shù)增加了3.62×10?10m2·s?1;當結構中的真空層厚度繼續(xù)增加到23 ?時,甲烷分子的擴散系數(shù)為7.64×10?10m2·s?1,相比于真空層厚度為18 ?時的擴散系數(shù)增加了2.92×10?10m2·s?1??梢钥闯稣婵諏雍穸嚷晕⒃黾?,就能引起煤層氣擴散系數(shù)的大幅度提高,因此煤結構中大的裂隙更有利于煤層氣的擴散。

        2.2.5影響煤層氣擴散行為的分子結構通過對不同結構擴散系數(shù)的模擬可得,煤層氣的擴散系數(shù)與芳香單元延展度的變化呈正比,由此可以解釋構造煤中煤層氣擴散速率高于相應原生結構煤的現(xiàn)象[8-9]。煤層氣的擴散系數(shù)隨著芳香單元堆砌層數(shù)的增加并不是一個單調的變化過程,而是一個受多種因素共同作用的復雜過程。煤分子表面含有不同缺陷和含氧官能團時,對煤層氣的擴散系數(shù)有不同程度的影響,其中,單缺陷的存在能夠增加煤層氣的擴散系數(shù),而次生缺陷的存在則降低煤層氣的擴散系數(shù)。比較4種含氧官能團對煤層氣擴散行為的影響發(fā)現(xiàn),只有當表面存在羰基的情況下,擴散系數(shù)要高于類石墨化的煤分子模型3-9×9,表面羧基的存在很大程度上降低了煤層氣的擴散系數(shù),而當表面存在羥基和甲氧基時,煤層氣的擴散系數(shù)與近石墨化的煤分子模型3-9×9相差不大。這與構造煤中含有較多的羰基,較少的羧基、醚鍵是符合的。真空層厚度越大,煤層氣的擴散系數(shù)越大。因此,煤結構中大的裂隙更有利于煤層氣的擴散。

        2.3煤層氣擴散微觀模型

        李志強等[39]指出煤粒中煤層氣擴散系數(shù)的大小取決于溫度、流體黏度與分子大小,并與擴散分子流動的平均速度平方呈正比。因此煤層氣的擴散系數(shù)并非常數(shù),而是隨時間延長呈動態(tài)衰減,基于這一關鍵特點,通過如下假設提出了煤顆粒徑向異性擴散系數(shù)模型。

        假設:①煤粒為同徑均質、異徑異性的徑向異性球體;②煤基質(擴散)孔隙系統(tǒng)由多尺度、大小不一的非均質多級孔隙組成,擴散系數(shù)假設為隨著時間延長而變化的變量D;③擴散過程采用Fick定律描述。

        在煤顆粒模型中,徑向結構是不均勻的(同徑同質),同時假定煤顆粒中賦存煤層氣分子的擴散在徑向上是先顆粒表層后顆粒內層(先外后內)。那么,在煤顆粒徑向向外方向的芳香層片延展度和間距增大(可能對應著孔隙的增大),含氧官能團和結構缺陷減少。這些因素均可共同導致沿徑向煤大分子對煤層氣分子束縛作用的減小,引起煤顆粒中煤層氣擴散系數(shù)隨徑向增加,從而導致煤層氣分子擴散系數(shù)表現(xiàn)出先大后小隨時間衰減的實驗結果。

        為了將煤層氣運移的微觀影響因素與宏觀表現(xiàn)相結合,在煤顆粒徑向異性擴散系數(shù)模型的基礎上,依據(jù)煤分子結構與外界環(huán)境對煤層氣吸附及運移的影響結果,提出了煤顆粒徑向不均質的煤層氣擴散微觀模型。如圖3所示。

        3 結 論

        依據(jù)實驗獲得的一系列煤的微晶結構參數(shù),構建了17個以煤的某一方面特征為變量的煤分子模型,通過Monte Carlo模擬方法與分子動力學模擬相結合,分別獲得各個煤分子結構對煤層氣的吸附量和擴散系數(shù),用于研究煤層氣的吸附與擴散。模擬溫度為303 K,壓力為10 MPa。模擬結果表明:芳香單元延展度并不會影響單位質量的煤對煤層氣的吸附體積;而隨著芳香單元堆砌層數(shù)的增加及表面缺陷和官能團的存在,單位質量的煤對煤層氣的吸附體積呈逐漸降低的趨勢,且在較低層數(shù)下,芳香單元堆砌層數(shù)的變化對吸附量的影響較大。但是,煤分子結構并不是影響煤層氣吸附量高低的主要因素,造成構造煤煤層氣高吸附量的主要因素應該是納米級孔結構的變化。

        圖3 煤顆粒徑向不均質的煤層氣擴散微觀模型Fig.3 Micro model of coalbed methane diffusion for coal particles of radial heterogeneous

        隨著芳香單元延展度的增加,煤層氣在表面的擴散系數(shù)增大,且在延展度較小時所引起的變化更加明顯;隨著芳香單元堆砌層數(shù)的增加,煤層氣在表面的擴散系數(shù)并不是一個單調的變化過程,而是一個受多種因素共同作用的復雜過程;單缺陷和含氧官能團C O的存在,提高了煤層氣的擴散速率;煤結構中大的裂隙更有利于煤層氣的擴散。

        在煤顆粒徑向異性擴散系數(shù)模型的基礎上,依據(jù)煤分子結構與外界環(huán)境對煤層氣吸附及運移的影響結果,提出了煤顆粒徑向不均質的煤層氣擴散微觀模型,進而將煤層氣運移的微觀影響因素與宏觀表現(xiàn)相結合。

        References

        [1]朱興珊, 徐鳳銀. 論構造應力場及其演化對煤和瓦斯突出的主控作用 [J]. 煤炭學報, 1994, 19 (3): 304-313. DOI: 10.13225/ j.cnki.jccs.1994.03.011. ZHU X S, XU F Y. The controlling effect of tectonic stress field and its evolution on coal and gas outburst [J]. J. China Coal Soc., 1994, 19 (3): 304-313. DOI: 10.13225/j.cnki.jccs.1994.03.011.

        [2]LEI D J, LI C W, ZHANG Z M, et al. Coal and gas outburst mechanism of the “Three Soft” coal seam in western Henan [J]. Min. Sci. Technol., 2010, 20 (5): 712-717. DOI: 10.1016/ S1674-5264(09)60268-4.

        [3]JU Y W, JIANG B, HOU Q L, et al. Behavior and mechanism of the adsorption desorption of tectonically deformed coals [J]. Chinese Sci. Bull., 2009, 54 (1): 88-94. DOI: 10.1007/s11434-008-0412-4.

        [4]王向浩, 王延斌, 高莎莎, 等. 構造煤與原生結構煤的孔隙結構及吸附性差異 [J]. 高校地質學報, 2012, 18 (3): 528-532. DOI: 10.16108/j.issn1006-7493.2012.03.003. WANG X H, WANG Y B, GAO S S, et al. Differences in pore structures and absorptivity between tectonically deformed and undeformed coals [J]. Geol. J. China U., 2012, 18 (3): 528-532. DOI: 10.16108/j.issn1006-7493.2012.03.003.

        [5]PAN J N, HOU Q L, JU Y W, et al. Coalbed methane sorption related to coal deformation structures at different temperatures and pressures [J]. Fuel, 2012, 102: 760-765. DOI: 10.1016/j.fuel.2012.07.023.

        [6]曹代勇, 李小明, 魏迎春, 等. 構造煤與原生結構煤的熱解成烴特征研究 [J]. 煤田地質與勘探, 2005, 33 (4): 39-41. CAO D Y, LI X M, WEI Y C, et al. Study on pyrogenation hydrocarbon generation character of deformed coal and undeformed coal [J]. Coal Geol. Explor., 2005, 33 (4): 39-41.

        [7]張玉貴, 張子敏, 曹運興. 構造煤結構與瓦斯突出 [J]. 煤炭學報, 2007, 32 (3): 281-284. DOI: 10.13225/j.cnki.jccs.2007.03.013. ZHANG Y G, ZHANG Z M, CAO Y X. Deformed-coal structure and control to coal-gas outburst [J]. J. China Coal Soc., 2007, 32 (3): 281-284. DOI: 10.13225/j.cnki.jccs.2007.03.013.

        [8]李云波, 張玉貴, 張子敏, 等. 構造煤瓦斯解吸初期特征實驗研究[J]. 煤炭學報, 2013, 38 (1): 15-20. DOI: 10.13225/ j.cnki.jccs.2013.01.017. LI Y B, ZHANG Y G, ZHANG Z M, et al. Experimental study on gas desorption of tectonic coal at initial stage [J]. J. China Coal Soc., 2013, 38 (1): 15-20. DOI: 10.13225/j.cnki.jccs.2013.01.017.

        [9]魏建平, 陳永超, 溫志輝. 構造煤瓦斯解吸規(guī)律研究 [J]. 煤礦安全, 2008, (8): 1-3. DOI: 10.13347/j.cnki.mkaq.2008.08.013. WEI J P, CHEN Y C, WEN Z H. Study of gas desorption laws of tectonically coal [J]. Safety in Coal Mines, 2008, (8): 1-3. DOI: 10.13347/j.cnki.mkaq.2008.08.013.

        [10]降文萍, 崔永君, 張群, 等. 不同變質程度煤表面與甲烷相互作用的量子化學研究結果 [J]. 煤炭學報, 2007, 32 (3): 292-295. DOI:10.13225/j.cnki.jccs.2007.03.016. JIANG W P, CUI Y J, ZHANG Q, et al. The quantum chemical study on different rank coals surface interacting with methane [J]. J. China Coal Soc., 2007, 32 (3): 292-295. DOI: 10.13225/j.cnki.jccs.2007.03.016.

        [11]QIU N X, XUE Y, GUO Y, et al. Adsorption of methane on carbon models of coal surface studied by the density functional theory including dispersion correction (DFT-D3) [J]. Comput. Theor. Chem., 2012, 992: 37-47. DOI: 10.1016/j.comptc.2012.04.024.

        [12]陳昌國, 魏錫文, 鮮學福. 用從頭計算法研究煤表面與甲烷分子的相互作用 [J]. 重慶大學學報(自然科學版), 2000, 23 (3): 77-79. CHEN C G, WEI X W, XIAN X F. Ab initio study on the interaction between CH4and the coal surface [J]. J. Chongqing U., 2000, 23 (3): 77-79.

        [13]LIU X Q, XUE Y, TIAN Z Y, et al. Adsorption of CH4on nitrogenand boron-containing carbon models of coal predicted by density-functional theory [J]. Appl. Surf. Sci., 2013, 285: 190-197. DOI:10.1016/j.apsusc.2013.08.035.

        [14]HU H X, LI X C, FANG Z M, et al. Small-molecule gas sorption and diffusion in coal: molecular simulation [J]. Energy, 2010, 35 (7): 2939-2944. DOI: 10.1016/j.energy.2010.03.028.

        [15]JU Y W, LI X S. New research progress on the ultrastructure of tectonically deformed coals [J]. Prog. Nat. Sci., 2009, 19 (11): 1455-1466. DOI: 10.1016/j.pnsc.2009.03.013.

        [16]JU Y W, YAN Z F, LI X S, et al. Structural characteristics and physical properties of tectonically deformed coals [J]. J. Geol. Res., 2012, 2012: 1-14. DOI: 10.1155/2012/852945.

        [17]CAO D Y, LI X M, ZHANG S R. Influence of tectonic stress on coalification: stress degradation mechanism and stress polycondensation mechanism [J]. Sci. China Ser. D, 2007, 50 (1): 43-54. DOI: 10.1007/s11430-007-2023-3.

        [18]LI X S, JU Y W, HOU Q L, et al. Response of macromolecular structure to deformation in tectonically deformed coal [J]. Acta Geol. Sin.-Engl., 2013, 87 (1): 82-90. DOI: 10.1111/1755-6724.12032.

        [19]JU Y W, JIANG B, HOU Q L, et al.13C-NMR spectra of tectonic coals and the effects of stress on structural components [J]. Sci. China Ser. D, 2005, 48 (9): 1418-1437. DOI: 10.1360/04yd0199.

        [20]WANG Z R, JIANG B, CHEN L X. The influence of intensity and properties of tectonic stress on the process of coalification [J]. J. Coal Sci. Eng., 2012, 18 (2): 158-162. DOI: 10.1007/s12404-012-0209-6.

        [21]LI X S, JU Y W, HOU Q L, et al. Spectra response from macromolecular structure evolution of tectonically deformed coal of different deformation mechanisms [J]. Sci. China Earth Sci., 2012, 55 (8): 1269-1279. DOI: 10.1007/s11430-012-4399-y.

        [22]姜波, 秦勇, 金法禮. 高溫高壓實驗變形煤XRD結構演化 [J]. 煤炭學報, 1998, 23 (2): 188-192. DOI: 10.13225/j.cnki.jccs.1998.02.016. JIANG B, QIN Y, JIN F L. XRD analysis of the structural evolution of deformed coal samples tested under high temperature and high confined pressure [J]. J. China Coal Soc., 1998, 23 (2): 188-192. DOI: 10.13225/j.cnki.jccs.1998.02.016.

        [23]屈爭輝, 姜波, 汪吉林. 構造煤結構演化及其應力-應變環(huán)境 [J].高校地質學報, 2012, 18 (3): 453-459. DOI: 10.16108/j.issn1006-7493.2012.03.017. QU Z H, JIANG B, WANG J L. Evolution of textures and stress-strain environments of tectonically-deformed coals [J]. Geol. J. China U., 2012, 18 (3): 453-459. DOI: 10.16108/j.issn1006-7493.2012.03.017.

        [24]CRACKNELL R F, NICHOLSON D, QUIRKE N. A grand canonical Monte-Carlo study of Lennard-Jones mixtures in slit shaped pores [J]. Mol. Phys., 1993, 80 (4): 885-897. DOI: 10.1080/00268979300102741.

        [25]SUN H. COMPASS: An ab initio force-field optimized for condensed-phase applications-overview with details on alkane and benzene compounds [J]. J. Phys. Chem. B, 1998, 102 (38): 7338-7364. DOI: 10.1021/jp980939v.

        [26]BUNTE S W, SUN H. Molecular modeling of energetic materials: the parameterization and validation of nitrate esters in the COMPASS forcefield [J]. J. Phys. Chem. B, 2000, 104 (11): 2477-2489. DOI: 10.1021/jp991786u.

        [27]LIU Y L, HU C J, ZHAO C C. Efficient parallel implementation of Ewald summation in molecular dynamics simulations on multi-core platforms [J]. Comput. Phys. Commun., 2011, 182 (5): 1111-1119. DOI: 10.1016/j.cpc.2011.01.007.

        [28]ANDERSEN H C. Molecular dynamics simulations at constant pressure and/or temperature [J]. J. Chem. Phys., 2008, 72 (4): 2384-2393. DOI: 10.1063/1.439486.

        [29]MüLLER-PLATHE F, ROGERS S C, VAN GUNSTEREN W F. Computational evidence for anomalous diffusion of small molecules in amorphous polymers [J]. Chem. Phys. Lett., 1992, 199 (3/4): 237-243. DOI: 10.1016/0009-2614(92)80112-O.

        [30]EINSTEIN A. Investigations on the Theory of the Brownian Movement [M]. New York: Dover Publication, Inc., 1956.

        [31]YANG J Z, LIU Q L, WANG H T. Analyzing adsorption and diffusion behaviors of ethanol/water through silicalite membranes by molecular simulation [J]. J. Membrane Sci., 2007, 291 (1/2): 1-9. DOI:10.1016/ j.memsci.2006.12.025.

        [32]LIU Q L, HUANG Y. Transport behavior of oxygen and nitrogen through organasilicon- containing polystyrenes by molecular simulation [J]. J. Phys. Chem. B, 2006, 110 (35): 17375-17382. DOI: 10.1021/jp063174x.

        [33]王繼仁, 鄧存寶, 鄧漢忠. 煤與瓦斯突出微觀機理研究 [J]. 煤炭學報, 2008, 33 (2): 131-135. DOI: 10.13225/j.cnki.jccs.2008.02.016. WANG J R, DENG C B, DENG H Z. Study on the microcosmic mechanism for coal-gas outburst [J]. J. China Coal Soc., 2008, 33 (2): 131-135. DOI: 10.13225/j.cnki.jccs.2008.02.016.

        [34]PENG D Y, ROBINSON D B. A new two-constant equation of state [J]. Ind. Eng. Chem. Fundamen., 1976, 15 (1): 59-64. DOI: 10.1021/ i160057a011.

        [35]JHON Y H, CHO M, JEON H R, et al. Simulations of methane adsorption and diffusion within alkoxy-functionalized IRMOFs exhibiting severely disordered crystal structures [J]. J. Phys. Chem. C, 2007, 111 (44): 16618-16625. DOI: 10.1021/jp0749470.

        [36]CONNOLLY M L. Analytical molecular surface calculation [J]. J. Appl. Cryst., 1983, 16: 548-558. DOI: 10.1107/S0021889883010985.

        [37]ZEITLER T R, ALLENDORF M D, GREATHOUSE J A. Grand canonical Monte-Carlo simulation of low-pressure methane adsorption in nanoporous framework materials for sensing applications [J]. J. Phys. Chem. C, 2012, 116 (5): 3492-3502. DOI: 10.1021/jp208596e.

        [38]HAO S X, WEN J, YU X P, et al. Effect of the surface oxygen groups on methane adsorption on coals [J]. Appl. Surf. Sci., 2013, 264: 433-442. DOI: 10.1016/j.apsusc.2012.10.040.

        [39]李志強, 王登科, 宋黨育. 新擴散模型下溫度對煤粒瓦斯動態(tài)擴散系數(shù)的影響 [J]. 煤炭學報, 2015, 40 (5): 1055-1064. DOI: 10.13225/j.cnki.jccs.2014.1218. LI Z X, WANG D K, SONG D Y. Influence of temperature on dynamic diffusion coefficient of CH4into coal particles by new diffusion model [J]. J. China Coal Soc., 2015, 40 (5): 1055-1064. DOI: 10.13225/j.cnki.jccs.2014.1218.

        Effects of coal molecular structure on adsorption and diffusion behaviors of coalbed methane

        WANG Baojun1, ZHANG Lina1, LING Lixia2, ZHANG Riguang1
        (1Key Laboratory of Coal Science and Technology of Ministry of Education and Shanxi Province, Taiyuan University of Technology, Taiyuan 030024, Shanxi, China;2College of Chemistry and Chemical Engineering, Taiyuan University of Technology, Taiyuan 030024, Shanxi, China)

        Abstract:The effects of microcrystal structures, such as the extension of basic structural unit (La) and the number of basic structural unit (N), as well as the surface structures including different types of defects and oxygen-containing functional groups, on the adsorption and diffusion of coalbed methane were investigated. The adsorption amounts were simulated by using the Monte Carlo (MC) method, and the diffusion coefficients of coalbed methane were simulated by using molecular dynamics (MD) simulation. All of the simulations were conducted in the condition of 303 K and 10 MPa. It can be concluded that the adsorption capacity of methane on unit mass of coal decreases with N increasing, as well as the existence of defects and oxygen-containing functional groups. The diffusion coefficient of methane increased with Laextending and shows a complex process influenced by increasing N. The existence of single-defect and C O was favorable of diffusion of coalbed methane. The coalbed methane was easy to diffuse in larger pore of coal. Finally, the micro model of coalbed methane diffusion for the coal particles of radial heterogeneous was proposed base on the micro-influenced factors and macro performance of the migration of coalbed methane.

        Key words:coal molecular structure; coalbed methane; adsorption; diffusion; molecular simulation

        中圖分類號:O 647.32

        文獻標志碼:A

        文章編號:0438—1157(2016)06—2548—10

        DOI:10.11949/j.issn.0438-1157.20151780

        基金項目:國家科技重大專項項目(2011ZX05040-005-002-001);國家自然科學基金項目(21276171,21576178,21276003);山西省青年學術帶頭人項目。

        Corresponding author:Prof. WANG Baojun, wangbaojun@tyut.edu.cn

        猜你喜歡
        擴散吸附
        一種吸附膜的制備及其吸附性能
        未來英才(2016年13期)2017-01-13 18:38:08
        Fe(Ⅲ)負載改性核桃殼對Cu2+吸附研究
        炭氣凝膠的制備與吸附氫氣性能
        核桃殼基吸附材料的制備及性能測試
        科技視界(2016年25期)2016-11-25 12:30:06
        混合無線多跳網(wǎng)絡安全性問題探討
        基于混沌理論的自適應參數(shù)圖像加密算法
        科技資訊(2016年18期)2016-11-15 07:45:11
        活化粉煤灰在煤焦油加工廢水預處理中的應用
        四環(huán)素類抗生素的環(huán)境行為研究進展
        科技視界(2016年11期)2016-05-23 08:07:47
        山東省產(chǎn)業(yè)集聚存在的問題與演進路徑探討
        山東省港口發(fā)展效果評價
        中國水運(2014年12期)2015-01-09 06:00:27
        一个色综合中文字幕人妻激情视频| 国产又黄又爽视频| 久久99久久99精品免观看不卡 | 无码不卡一区二区三区在线观看| 日本精品一区二区三区在线播放| 国产熟女露脸91麻豆| 天天噜日日噜狠狠噜免费| 亚洲国产欧美日韩欧美特级| 欧美黑人又粗又大久久久| 中文字幕乱码亚洲无线| 一级黄色一区二区三区| 秘书边打电话边被躁bd视频| 亚洲av无码成人yellow| 伊人色综合九久久天天蜜桃| 四季极品偷拍一区二区三区视频| 巨人精品福利官方导航| 免费人成毛片乱码| 91在线视频视频在线| 丁香婷婷激情视频在线播放| 久久久久波多野结衣高潮| 无码 制服 丝袜 国产 另类 | 日本a爱视频二区三区| 少妇人妻在线无码天堂视频网| 黑人巨大videos极度另类| 日本一区二区三区看片 | 亚洲精品国产精品系列| 国产亚洲精品综合一区| 亚洲欧美在线观看| 欧美z0zo人禽交欧美人禽交| 久久久9色精品国产一区二区三区| 人妻少妇艳情视频中文字幕| 欧美一性一乱一交一视频| 性导航app精品视频| 人妖系列在线免费观看| 四虎成人精品在永久免费| 中文字幕无码不卡一区二区三区 | 国产真实乱对白精彩| 亚洲色欲大片AAA无码| 国产麻豆极品高清另类| 久久亚洲日韩精品一区二区三区| 亚洲av第一成肉网|