• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    煤儲層粗糙割理中煤層氣運(yùn)移機(jī)理數(shù)值分析

    2014-06-07 05:55:20祝一搏鄭軍領(lǐng)董佳斌
    煤炭學(xué)報 2014年9期
    關(guān)鍵詞:運(yùn)移格子維數(shù)

    金 毅,祝一搏,吳 影,鄭軍領(lǐng),董佳斌,李 翔

    (1.河南理工大學(xué)資源環(huán)境學(xué)院,河南焦作 454003;2.中原經(jīng)濟(jì)區(qū)煤層(頁巖)氣河南省協(xié)同創(chuàng)新中心,河南焦作 454003)

    煤儲層粗糙割理中煤層氣運(yùn)移機(jī)理數(shù)值分析

    金 毅1,2,祝一搏1,吳 影1,鄭軍領(lǐng)1,董佳斌1,李 翔1

    (1.河南理工大學(xué)資源環(huán)境學(xué)院,河南焦作 454003;2.中原經(jīng)濟(jì)區(qū)煤層(頁巖)氣河南省協(xié)同創(chuàng)新中心,河南焦作 454003)

    有效描述粗糙割理中煤層氣運(yùn)移機(jī)理是實現(xiàn)其產(chǎn)能準(zhǔn)確評估的基本前提。結(jié)合裂隙的立方定律以及廣義Kozeny-Carman孔-滲方程,于理論層面推導(dǎo)了考慮內(nèi)、外摩擦所致的裂-滲關(guān)系模型。其中內(nèi)摩擦效應(yīng)以煤層氣運(yùn)移的水文彎曲度來描述,外摩擦則通過引入端面曲折率來定義。借助分形理論實現(xiàn)了割理端面粗糙幾何的定量描述,并采用格子Boltzmann方法模擬了煤層氣運(yùn)移過程,分析了端面分形維、絕對粗糙度以及相對粗糙度對煤層氣輸運(yùn)特征的影響。在此基礎(chǔ)上,對比分析了新裂-滲方程解析值同數(shù)值模擬滲透率之間的關(guān)系。結(jié)果表明,考慮了內(nèi)、外摩擦效應(yīng)的新方程能有效描述微觀粗糙割理中煤層氣的運(yùn)移規(guī)律,并且物理意義明確。

    煤儲層;粗糙割理;煤層氣;裂隙立方定律;Kozeny-Carman方程;格子Boltzmann方法

    煤是一種具有雙重孔隙結(jié)構(gòu)的多孔介質(zhì),其中基質(zhì)孔隙是煤層氣的主要儲集場所,而裂隙及割理網(wǎng)絡(luò)主宰著煤儲層輸運(yùn)屬性。煤層氣產(chǎn)能評估的一個重要物理參數(shù)就是煤儲層介質(zhì)的滲透率,因此定量描述割理中煤層氣運(yùn)移規(guī)律就成為了預(yù)測其產(chǎn)能的基礎(chǔ)。研究表明[1-4]:光滑平行端面所組成的裂隙空間中,等溫、層流條件下流體絕對滲透率是裂隙開度b的函數(shù):

    其中,kabs為絕對滲透率;H為裂隙橫截面最小外包的高度;φ為孔隙度;β為形狀因子。光滑裂隙中β=3, b=H且φ=1,式(1)所示的立方定律可簡化為kabs= b2/12的形式。

    然而,自然裂隙端面幾何形貌異常粗糙,其起伏高度滿足分形統(tǒng)計特征,研究表明這種粗糙端面對裂隙輸運(yùn)屬性有十分重要的影響,特別是在微觀裂隙[4-10]中。而煤中割理形成于不同尺度、強(qiáng)度的地質(zhì)作用及內(nèi)部過程,其端面幾何同樣表現(xiàn)出分形粗糙特征,因此采用經(jīng)典的立方定律(式(1))評估煤儲層滲透率必然導(dǎo)致理論預(yù)測同實際測試結(jié)果之間的巨大偏差。實際應(yīng)用中,為了降低實測值同理論預(yù)測結(jié)果之間的差異,一種所謂的粗糙度因子(fr)被引入到立方定律中[2]。

    式中,kr為修正后的滲透率;〈b〉為粗糙裂隙的等效開度;fr為一無量綱參數(shù),且滿足fr≥1的關(guān)系,如果裂隙端面光滑平行,則有fr=1,此時式(2)等價于式(1)。

    細(xì)觀式(2)可知,引入fr的根本目的是為了削弱粗糙裂隙中的等效開度〈b〉。如果〈b〉?σ時(其中σ為端面絕對粗糙度,常以端面高度均方根rmh表示),因宏觀流總體上的直線運(yùn)動特征,故削弱效應(yīng)是可以接受的;但當(dāng)〈b〉≈σ或〈b〉?σ時,因裂隙空間的幾何限制,流體運(yùn)移過程除了受阻于外部摩擦(邊界摩擦)外,曲折流效應(yīng)會大大提升層流線間的接觸長度,故流體間的內(nèi)摩擦作用被增強(qiáng)。

    因此,探明微觀割理中粗糙端面對流體運(yùn)移的多重阻礙作用,包括開度的削弱、流-固間的外摩擦增量以及彎曲流所致的內(nèi)摩擦增量,并定量描述微觀幾何同其輸運(yùn)屬性之間的關(guān)系就成為了煤層氣產(chǎn)能評估的基礎(chǔ)。

    由于微裂隙流的重要性,在過去的數(shù)十年中,不同的學(xué)科對這種裂-滲關(guān)系展開了大量研究,包括實驗測試、理論分析以及數(shù)值模擬等[3,11-19]。實驗研究受測試技術(shù)、研究尺度和實驗環(huán)境的影響明顯;而現(xiàn)有的理論方程均于裂隙的簡單假設(shè)條件下推導(dǎo)而得(式(1)),對于復(fù)雜的裂隙而言往往失去效用。

    正因如此,如今更多的研究轉(zhuǎn)向數(shù)值模擬方向。因為在數(shù)值分析中,任何相關(guān)的影響因素都可被輕易選取或忽略,以便對耦合環(huán)境下各主要因素進(jìn)行獨(dú)立的分析并定量描述。同時數(shù)值模擬能有效克服實驗測試中的種種限制,如今已成為發(fā)掘微觀裂隙中流體運(yùn)移機(jī)理的重要手段[20-21]。其中,新近發(fā)展起來的格子Boltzmann方法(LBM)引起了學(xué)術(shù)界的高度重視。相對于其他方法而言,LBM優(yōu)勢明顯,主要表現(xiàn)為復(fù)雜邊界的處理能力以及“自下而上”由質(zhì)到本的演化策略[22-28]。越來越多的證據(jù)表明[12,29-31],LBM是當(dāng)今發(fā)掘復(fù)雜流控制機(jī)理的有效手段。同時,國內(nèi)學(xué)者利用該方法,從不同的側(cè)面模擬,分析了煤儲層中煤層氣的運(yùn)移規(guī)律[32-37],但均基于宏觀或介觀尺度上,微觀尺度下的研究較少[38-39],而考慮煤層氣曲折流效應(yīng)的微裂隙流方面的研究鮮見報道。

    為了全面理解煤儲層微割理中端面幾何分形粗糙特征對煤層氣運(yùn)移行為的影響,筆者首先從理論層面推導(dǎo)了粗糙裂隙空間中裂-滲關(guān)系,隨后采用Syn-Frac軟件[40]構(gòu)建了不同粗糙特征的割理復(fù)合拓?fù)淠P?并在孔隙尺度下采用LBM模擬了煤層氣的運(yùn)移過程。基于以上實驗結(jié)果,系統(tǒng)分析了微觀物性參數(shù)對流體運(yùn)移的控制作用并驗證了新裂-滲方程的有效性。

    1 理論和方法

    1.1 割理端面幾何特征及其分形描述

    傳統(tǒng)歐氏幾何主要用于描述規(guī)則結(jié)構(gòu),其維數(shù)為整數(shù)。而分形幾何則以廣義的非整數(shù)維(分形維)概念實現(xiàn)了對復(fù)雜多變的非規(guī)則結(jié)構(gòu)的有效描述[42]。已有研究表明,自然裂隙端面幾何多表現(xiàn)出多尺度、自相似的粗糙特性[21,41-45]。而割理端面的微觀復(fù)雜幾何可由3個參數(shù)于統(tǒng)計層面等效描述[40,46],包括端面分形維數(shù)Df、端面絕對粗糙度σ以及隨機(jī)相集合{ψ},即

    其中,x為平面坐標(biāo);F為裂隙端面構(gòu)造函數(shù);Z(x)為x處的端面高度。其中分形維數(shù)控制端面的高、低頻成分分布,即Df越大,高頻成分越多,沿割理展布方向Z(x)的變化越快(圖1),微觀結(jié)構(gòu)更復(fù)雜;σ是一個獨(dú)立于頻率成分的尺度因子或粗糙度因子,σ越大,則垂直于割理展布方向上的起伏越大,端面越粗糙(圖2)。而引入隨機(jī)相集合{φ}的目的是確保裂隙端面的復(fù)雜無序特征。

    圖1 分形維對割理端面粗糙影響(σ=1.0 mm)Fig.1 Effects of fractal dimension on the roughness of cleat surface with same scale factor(σ=1.0 mm)

    圖2 分形維確定情況下,裂隙端面高度均方根對其粗糙度的影響Fig.2 Effect of rmh on the roughness of surface profiles with same fractal dimensions

    1.2 粗糙割理中煤層氣運(yùn)移的LBM模擬

    1.2.1 LBM算法及邊界處理

    LBM是一種介觀尺度上的流體運(yùn)移模擬方法,在其格子系統(tǒng)中,空間D被離散成等邊長規(guī)則格子(通常為簡單立方體),時間t離散為間隔h的序列,速度空間分解為向量集{ci},且確保cih能連接相鄰格點[47]。如常用的D3Q19格子模型采用19個速度向量來離散格子的速度空間[48],除了部分靜止粒子外,其他速度向量必須確保其中6個粒子能到達(dá)最近鄰域,12粒子到達(dá)次近鄰域。記具有速度ci的粒子質(zhì)量密度分布為fi(x,t),LBM算法通過操作fi(x,t)確保局部流體密度ρ(x,t)滿足

    類似地,動量密度M(x,t)可表示為

    因此可得格子的平均流速為

    綜上所述,描述以上算法的格子Boltzmann方程(LBE)可表達(dá)為

    LBE中,通過碰撞算子Ωi修改不同位置粒子密度分布({fi}代表不同位置粒子分布的集合)以確保質(zhì)量和動量的守恒(式(8)),碰撞后的粒子質(zhì)量分布隨后流向鄰近的位置。

    雖然LBE滿足了流體流動模擬過程中的一些基本要求,即質(zhì)量、能量守恒以及局部性。但因離散速度的限制,式(8)無法保證伽利略不變性以及各向同性動量的傳輸。而基于單一松弛時間的簡化碰撞函數(shù)[25,48]——Bhatnagar-Gross-Krook模型(LBGK)有效地克服了以上的問題,其碰撞算子為式中,tlbm和(x,t)分別為無量綱松弛時間和平衡分布函數(shù)。

    選用D3Q19型的三維LBGK模型,其速度向量如圖3所示。其中,D3Q19型的三維LBGK模型的平衡態(tài)分布函數(shù)為

    其中,cs為聲速;ωi為權(quán)重系數(shù)。出于對稱性原因,權(quán)重系數(shù)應(yīng)滿足條件ωi>0和∑iωi=1,以上限定確保了ωi僅是依賴于速度ci的絕對值,而與方向無關(guān)。同時,權(quán)重的設(shè)定必須使平衡分布函數(shù)滿足如下屬性:

    圖3 D3Q19格子模型的離散速度Fig.3 Velocity directions of D3Q19 lattice model

    基于以上約束條件,可得D3Q19格子模型的權(quán)重系數(shù)分別為1/3(i=0),1/18(i=1~6)和1/36(i= 7~18)。

    同時,因割理端面無序及異常復(fù)雜的微觀幾何所致的流-固邊界上幾乎不存在連續(xù)流動[23]。在LBM模擬中,這種情況常采用反彈模型作為物理邊界條件,為了使問題簡化同時又不失一般性,使用完全反彈的策略來近似割理表面煤層氣的邊界響應(yīng)[29]。

    1.2.2 物理模型同格子系統(tǒng)間的轉(zhuǎn)換

    在模擬過程中,割理復(fù)合拓?fù)淇臻g離散為格子系統(tǒng),以0表示孔隙相,1表示固相。為了降低計算誤差,割理水平展布方向的離散分辨率設(shè)為1 024個格子單位,而實際物理尺寸確定為1 μm。二值轉(zhuǎn)換如圖4所示。

    圖4 二維裂隙空間的二值化表達(dá)Fig.4 Binary representation of two-dimensional fracture space

    本文數(shù)值研究中,煤層氣被假想為單組分的理想氣體,LBGK模型為流體運(yùn)移演化控制體系。為保證數(shù)值模擬的穩(wěn)定性,無量綱松弛時間τlbm取值1,該條件下完全反彈邊界模型所產(chǎn)生的誤差極小[49]。

    同時,為確保割理中煤層氣的層流特征及達(dá)西定律的有效性,壓力梯度設(shè)置為ΔP≤10-5Pa/m,模擬滲透率可由式(13)獲得。

    翻轉(zhuǎn)課堂是“課前”和“課堂”兩個學(xué)習(xí)環(huán)節(jié)緊密配合的教學(xué)組織形式,體現(xiàn)了“學(xué)生主體,教師主導(dǎo)”的教育理念,課堂上主要以相應(yīng)的學(xué)習(xí)活動為學(xué)習(xí)發(fā)生的主要場所。楊開城教授在《以學(xué)習(xí)活動為中心的教學(xué)設(shè)計理論:教學(xué)設(shè)計理論新探索》一書中提到,“以學(xué)習(xí)活動為中心的教學(xué)設(shè)計理論將學(xué)習(xí)活動作為組成教學(xué)系統(tǒng)的基本單元”。何克抗教授“雙主理念”的教學(xué)設(shè)計模式主要針對的是課堂上45分鐘的教學(xué)組織形式。鑒于此,我們需要結(jié)合“以學(xué)習(xí)活動為中心”的教學(xué)設(shè)計理論,在“雙主理念”教學(xué)設(shè)計模式的基礎(chǔ)上進(jìn)行更加具體地調(diào)整、修改和發(fā)展而形成真正適合翻轉(zhuǎn)課堂的教學(xué)設(shè)計模式(如圖1)。

    其中,klbm為格子系統(tǒng)下無量綱滲透率;ΔP=-dP/dL為壓力梯度,L為宏觀流方向長度;υ為流體黏度,因tlbm=1,故有υ=0.166 7;U為格子平均速度。為確保格子系統(tǒng)模擬的無量綱滲透率klbm同割理實際物理滲透率kns之間的等效性,需采用式(14)[49]進(jìn)行轉(zhuǎn)換。

    其中,Lphysical為割理的實際物理尺寸,有量綱參數(shù); Llbm為離散后對應(yīng)的格子數(shù)量。本文割理水平尺寸為1 μm,離散后格子數(shù)量為1 024個,因此一個格子邊長對應(yīng)的實際尺寸為1/1 024 μm,這不但基本滿足了對粗糙裂隙幾何形貌的有效描述,同時也保證了數(shù)值模擬的精度。

    2 理論及數(shù)值分析

    2.1 耦合曲折效應(yīng)的裂-滲關(guān)系

    如前所述,在微觀裂隙中端面粗糙幾何除了會削弱等效開度的大小外,還會引入內(nèi)摩擦增量,此時內(nèi)摩擦因子大于1,記其為無量綱參數(shù)fτ,因此可得此時滲透率為

    與此同時,廣義的多孔介質(zhì)孔-滲關(guān)系,即Kozeny-Carman(K-C)方程,也包括對裂-滲關(guān)系的描述。其基本形式為其中,kkc為多孔介質(zhì)滲透率;C0為待定系數(shù);τ為水文彎曲度;S為比表面積,在裂隙系統(tǒng)中的定義為S=(As1+As2)/(HA0),其中A0為裂隙端面水平投影面積,As1和As2為復(fù)合端面的面積。因端面相同的巖石力學(xué)性質(zhì)及相似的經(jīng)歷,故滿足〈As1〉=〈As2〉,記為As,其中〈…〉表示期望或均值。

    結(jié)合裂隙空間的物性特征,式(16)可改寫為如下的形式:

    如前所述,外摩擦來源于裂隙端面的粗糙以及由此引入的端面曲折率τs;而曲折流效應(yīng)導(dǎo)致了流體之間相對接觸面的增大。結(jié)合式(18),可將外摩擦因子fr和內(nèi)摩擦因子fτ分別表示為

    進(jìn)而得粗糙裂隙中裂-滲方程為

    如果裂隙由光滑且水平的端面組成,則有τ=1, τs=1,φ=1,此時式(21)退化為的形式;對于端面幾何粗糙特征相同的裂隙空間,隨著開度〈b〉的不斷增加,φ→1,τ→1,此時fτ→1,fr→τs,此時表現(xiàn)出粗糙端面對等效開度的削弱作用。

    2.2 數(shù)值模擬結(jié)果與分析

    利用 Glover[40]分形粗糙裂隙端面描述方法(式(3))模擬了不同物性特征的煤巖割理復(fù)合拓?fù)淠P?隨后,在孔隙尺度下利用LBM模擬了煤層氣的運(yùn)移過程。

    圖5展示了分形維數(shù)Df、絕對粗糙度σ對裂隙端面幾何的影響。對比分析可知隨著Df的增加,裂隙延展方向上高度起伏頻率明顯加劇,如圖5(a),(b)所示,雖然其具有相同的σ值,但因端面分形維不同,圖5(b)的粗糙性明顯增大,對于微觀裂隙而言,必然導(dǎo)致端面曲折率τs的提升,流阻增大;而σ主要控制裂隙端面高度起伏的強(qiáng)度,即σ越大,端面起伏高度越大,此時如果限定裂隙最小外包H,隨著σ的增加,必然導(dǎo)致裂隙的平均開度〈b〉及孔隙度φ的減小,而τs和τ均會增加。

    圖5 三維割理端面的分形模擬Fig.5 Fractal representations of 3D cleat surfaces

    圖6為LBM模擬的粗糙割理中煤層氣運(yùn)移特征其中的流線為流體運(yùn)移的軌跡,灰體為固相基質(zhì)。流體軌跡直觀說明了粗糙割理中煤層氣的曲折流效應(yīng)。

    為了能更好地理解分形維數(shù)Df,σ以及σ/〈b〉對煤層氣運(yùn)移彎曲度的影響,基于LBM模擬的流場數(shù)據(jù)以及裂隙的微觀物性參數(shù)重點分析了σ不變,不同相對粗糙度σ/〈b〉條件下水文彎曲度τ同分形維數(shù)Df的關(guān)系,以及保持分形維數(shù)Df不變條件下,相對粗糙度對流體水文彎曲度增量的影響,其中水文彎曲度采用文獻(xiàn)[50]局部彎曲度平均法計算而得。

    圖7為相對粗糙度σ/〈b〉分別為0.400,0.533以及0.800時,水文彎曲度隨割理端面分形維的變化趨勢。在相對粗糙度不變的情況下,水文彎曲度隨端面分形維的增大而線性增加,并均表現(xiàn)出分段的特征:當(dāng)Df≤2.5時,Df對水文彎曲度τ的影響較弱,但當(dāng)Df>2.5時,這種影響明顯加劇。與此同時,相對粗糙度σ/〈b〉越小,水文彎曲度τ越低,該趨勢表明:在σ固定的情況下,隨著割理開度〈b〉的增加,煤層氣運(yùn)移逐漸趨于直線流動,或者是開度〈b〉不變時,裂隙的流阻會隨著割理端面粗糙度σ的降低而減小。

    圖6 粗糙割理中煤層氣運(yùn)移的曲折流效應(yīng)Fig.6 Tortuous effect of CBM flow through cleat with rough surfaces

    圖7 水文彎曲度隨分形維數(shù)變化的分段特征Fig.7 Dual linear relationships between hydraulic tortuosity and surface fractal dimension for fractures

    為探明相對粗糙度對水文彎曲度增量的影響,割理端面分形特征被剔除,即Df不變時相對粗糙度σ/〈b〉同彎曲度增量τ-1的關(guān)系,如圖8所示。結(jié)果表明:當(dāng)σ/〈b〉→0時,τ-1→0,即開度〈b〉足夠大時,割理粗糙端面對流體曲折流效應(yīng)可以忽略不計;隨著σ/〈b〉的增加,端面粗糙特征對流體運(yùn)移軌跡的控制作用逐漸增強(qiáng)。

    而圖8中的分段特征主要是由于開度〈b〉接近或小于絕對粗糙度σ時(σ/〈b〉=1.5附近),裂隙微觀幾何對流體的限制作用增強(qiáng)、流體的動力學(xué)特征被削弱。此時水文彎曲度會逐漸接近但永遠(yuǎn)低于幾何彎曲度,對于確定的割理系統(tǒng)而言,幾何彎曲度固定,這也是當(dāng)σ/〈b〉>1.5時,水文彎曲度隨σ/〈b〉的線性變化趨緩的根本原因。

    圖8 水文彎曲度增量同相對粗糙度之間的關(guān)系Fig.8 Relationship between relative surface roughness and the increment of hydraulic tortuosity

    依據(jù)現(xiàn)有的滲透率預(yù)測模型,影響割理滲透率的另一個因素是流-固間的外摩擦效應(yīng),如式(21)所示。結(jié)合粗糙割理端面的分形假設(shè),分析了絕對粗糙度σ、端面分形維數(shù)Df對端面曲折率τs的影響,如圖9所示。其中σ對τs是一種線性影響,而分形維數(shù)Df的影響雖然也表現(xiàn)出線性特征,但是影響程度在Df≈2.5時出現(xiàn)突變。

    為驗證式(21)的有效性,對比分析了LBM模擬滲透率kns同式(2)和(21)解析滲透率之間的關(guān)系。它們之間的主要差別就是式(21)中有效地包含了內(nèi)、外摩擦增量,分別以水文彎曲度τ以及割理端面曲折率τs來定量描述。

    采用式(21)解析估算時,形狀因子β取值3,同時由于待定系數(shù)C0主要受控于割理端面的分形特征,因此針對不同的Df分別分析了式(2)和(21)的解析滲透率同模擬滲透率kns之間的關(guān)系,結(jié)果如圖10所示。

    圖10分別為割理端面分形維數(shù)Df=2.3,2.5和2.7的情況。當(dāng)Df較小時,式(2)同LBM模擬滲透率kns之間的相關(guān)性較好,但其理論預(yù)測值均大于模擬值;隨著分形維數(shù)的增加,其預(yù)測結(jié)果同模擬滲透率之間不再線性相關(guān),表現(xiàn)出隨相對粗糙度的增加而相關(guān)性變差的特點,這種趨勢表明,當(dāng)端面絕對粗糙度〈σ〉固定時,隨著開度〈b〉的減小,式(2)的預(yù)測值比實際模擬滲透率大且呈增加的趨勢(圖10(c))。而圖中式(21)理論預(yù)測值同實際模擬的kns之間不但存在很好的線性關(guān)系,并且這種關(guān)系不因端面分形維數(shù)Df而改變。雖然在圖10中,式(21)中的C0不同,如Df=2.3,2.5,2.7時分別為0.67,0.38和0.20,其主要原因是分形維對形狀因子的影響。

    圖9 割理端面絕對粗糙度和割理端面分形維數(shù)同端面曲折率的關(guān)系Fig.9 Relationships between absolute roughness,surface fractal dimension and surface tortuosity of cleat surfaces

    圖10 LBM模擬滲透率同理論預(yù)測模型式(2),式(21)解析值之間的關(guān)系Fig.10 Relationship between the permeability by LBM and those estimated by Eq.(2)and Eq.(21)

    3 結(jié)論與討論

    (1)構(gòu)建了粗糙割理中考慮內(nèi)、外摩擦增量的煤層氣裂-滲方程。

    (2)確定了內(nèi)摩擦效應(yīng)同水文彎曲度、外摩擦效應(yīng)同割理端面曲折率之間的定量關(guān)系。

    (3)對比不同端面幾何對割理輸運(yùn)屬性的影響表明:水文彎曲度同割理端面幾何分形維及相對粗糙度之間均存在分段線性關(guān)系;端面曲折率同其絕對粗糙度之間滿足線性關(guān)系,但同分形維之間同樣表現(xiàn)出分段線性相關(guān)特征。

    本次研究只考慮了單割理中煤層氣運(yùn)移規(guī)律,實際應(yīng)用中需要考慮煤儲層中割理的空間分布密度、端面幾何粗糙屬性以及割理開度等因素。在此基礎(chǔ)上,結(jié)合式(21)所表征的裂-滲關(guān)系定能為煤層氣產(chǎn)能評估提供更加準(zhǔn)確的評估及預(yù)測。同時在后續(xù)研究中,要有效構(gòu)建割理開度同割理端面幾何之間的尺度不變關(guān)系,如水文彎曲度、割理端面曲折率同開度之間的冪律關(guān)系,為煤層氣產(chǎn)能評估及煤儲層輸運(yùn)屬性的科學(xué)評價提供理論支持及明確參數(shù)解釋。

    [1] Eker E,Akin S.Lattice Boltzmann simulation of fluid flow in synthetic fractures[J].Transport in Porous Media,2006,65(3):363-384.

    [2] Witherspoon P A,Wang J S Y,Iwai K,et al.Validity of cubic law for fluid flow in a deformable rock fracture[J].Water Resources Research,1980,16(6):1016-1024.

    [3] Jacob B.Dynamics of fluids in porous media[Z].New York:Elsevier Science,1972.

    [4] Snow D T.Anisotropie permeability of fractured media[J].Water Resources Research,1969,5(6):1273-1289.

    [5] Zou M Q,Yu B M,F Y J,et al.A Monte Carlo method for simulating fractal surfaces[J].Physica A:Statistical Mechanics and its Applications,2007,386(1):176-186.

    [6] Schmittbuhl J,Gentier S,Roux S.Field measurements of the roughness of fault surfaces[J].Geophysical Research Letters,1993, 20(8):639-641.

    [7] Wang J S Y,Cox B L.Fractal surfaces:Measurement and applications in the earth sciences[J].Fractals,1993,1(1):87-115.

    [8] Vickers B C,Neuman S P,Sully M J,et al.Reconstruction and geostatistical analysis of multiscale fracture apertures in a large block of welded tuff[J].Geophysical Research Letters,1992, 19(10):1029-1032.

    [9] Brown S R.Fluid flow through rock joints:The effect of surface roughness[J].Journal of Geophysical Research,1987,92(B2): 1337-1347.

    [10] Brown S R,Scholz C H.Broad bandwidth study of the topography of natural rock surfaces[J].Journal of Geophysical Research,1985, 90(B14):12512-12575.

    [11] Chen L,Kang Q,Robinson B A,et al.Pore-scale modeling of multiphase reactive transport with phase transitions and dissolution-precipitation processes in closed systems[J].Physical Review E, 2013,87(4):43306.

    [12] Cai J C,Hu X Y,Standnes D C,et al.An analytical model for spontaneous imbibition in fractal porous media including gravity[J].Colloids and Surfaces A:Physicochemical and Engineering Aspects,2012,414:228-233.

    [13] Mitra A,Harpalani S,Liu S.Laboratory measurement and modeling of coal permeability with continued methane production:Part 1-Laboratory results[J].Fuel,2012,94:110-116.

    [14] Pan Z J,Connell L D.Modelling permeability for coal reservoirs:A review of analytical models and testing data[J].Internatial Journal of Coal Geology,2012,92:1-44.

    [15] Liu J S,Chen Z W,Elsworth D,et al.Interactions of multiple processes during CBM extraction:A critical review[J].International Journal of Coal Geology,2011,87(3-4):175-189.

    [16] Cai Y D,Liu D M,Yao Y B,et al.Geological controls on prediction of coalbed methane of No.3 coal seam in Southern Qinshui Basin, North China[J].International Journal of Coal Geology,2011,88(2 -3):101-112.

    [17] Xu H,Sang S,Fang L,et al.Production characteristics and the control factors of surface wells for relieved methane drainage in the huainan mining area[J].Acta Geologica Sinica-English Edition, 2011,85(4):932-941.

    [18] Connell L D,Lu M,Pan Z J.An analytical coal permeability model for tri-axial strain and stress conditions[J].International Journal of coal Geology,2010,84(2):103-114.

    [19] Fu X H,Qin Y,Wang G G X,et al.Evaluation of coal structure and permeability with the aid of geophysical logging technology[J].Fuel,2009,88(11):2278-2285.

    [20] Croce G,Paola D A,Nonino C.Three-dimensional roughness effect on microchannel heat transfer and pressure drop[J].International Journal of Heat and Mass Transfer,2007,50(25-26):5249-5259.

    [21] Croce G,Paola D A.Numerical analysis of roughness effect on microtube heat transfer[J].Superlattices and Microstructures,2004, 35(3-6):601-616.

    [22] Nourgaliev R R,Dinh T N,Theofanous T G,et al.The lattice Boltzmann equation method:theoretical interpretation,numerics and implications[J].International Journal of Multiphase Flow,2003, 29(1):117-169.

    [23] Succi S.The lattice Boltzmann equation for fluid dynamics and beyond[M].Oxford:Clarendon Press,2001.

    [24] Kandhai D,Vidal D J E,Hoekstra A G,et al.Lattice-Boltzmann and finite element simulations of fluid flow in a SMRX Static Mixer Reactor[J].International Journal for Numerical Methods in Fluids,1999,31(6):1019-1033.

    [25] Chen S,Doolen G D.Lattice Boltzmann method for fluid flows[J].Annual Review of Fluid Mechanics,1998,30(1):329-364.

    [26] Koponen A,Kataja M,Timonen J.Permeability and effective porosity of porous media[J].Physical Review E,1997,56(3):3319-3325.

    [27] Ladd A J C.Numerical simulations of particulate suspensions via a discretized Boltzmann equation.Part 2.Numerical results[J].Journal of Fluid Mechanics,1994,271(1):311-339.

    [28] Ladd A J C.Numerical simulations of particulate suspensions via a discretized Boltzmann equation.Part 1.Theoretical foundation[J].Journal of Fluid Mechanics,1994,271(1):285-309.

    [29] 郭照立,鄭楚光.格子Boltzmann方法的原理及應(yīng)用[M].北京:科學(xué)出版社,2009.

    [30] Qian J,Li Q,Yu K,et al.A novel method to determine effective thermal conductivity of porous materials[J].Science in China Series E:Technological Sciences,2004,47(6):716-724.

    [31] Cai J C,Yu B M,Zou M Q,et al.Fractal analysis of invasion depth of extraneous fluids in porous media[J].Chemical Engineering Science,2010,65(18):5178-5186.

    [32] 金 毅.煤微觀結(jié)構(gòu)描述及其輸運(yùn)屬性數(shù)值模擬分析[D].北京:北京大學(xué),2011.

    Jin Yi.Quantitative description of the microstructure of coal reservoir and numerical simulation of its transport property[D].Beijing: Peking University,2011.

    [33] 陸秋琴.地下煤礦瓦斯運(yùn)移數(shù)值模擬及積聚危險性評價研究[D].西安:西安建筑科技大學(xué),2010.

    Lu Qiuqin.Research on numerical simulation of gas migration and risk evaluation of gas aggregation in underground coal mine[D].Xi’an:Xi’an University of Architecture and Technology,2010.

    [34] 邢玉飛.基于LBM地下煤層瓦斯運(yùn)移仿真實現(xiàn)[D].西安:西安建筑科技大學(xué),2009.

    Xing Yufei.Simulation of gas flow in underground coal seam based on LBM[D].Xi’an:Xi’an University of Architecture and Technology,2009.

    [35] 滕桂榮,譚云亮,高 明,等.基于LBM方法的裂隙煤體內(nèi)瓦斯抽放的模擬分析[J].煤炭學(xué)報,2008,33(8):914-919.

    Teng Guirong,Tan Yunliang,Gao Ming,et al.Simulation of gas drainage in fissured coal based on Lattice Boltzmann method[J].Journal of China Coal Society,2008,33(8):914-919.

    [36] 朱益華,陶 果,方 偉,等.低滲氣藏中氣體滲流Klinkenberg效應(yīng)研究進(jìn)展[J].地球物理學(xué)進(jìn)展,2007,22(5):1591-1596.

    Zhu Yihua,Tao Guo,Fang Wei,et al.Research progress of the Klinkenberg effcet in tight gas reservoir[J].Progress in Geophysics,2007,22(5):1591-1596.

    [37] 滕桂榮,譚云亮,高 明.基于lattice Boltzmann方法對裂隙煤體中瓦斯運(yùn)移規(guī)律的模擬研究[J].巖石力學(xué)與工程學(xué)報, 2007,26(S1):3503-3508.

    Teng Guirong,Tan Yunliang,Gao Ming.Simulation of gas seepage in fissured coal based on Lattice Boltzmann method[J].Chinese Journal of Rock Mechanics and Engineering,2007,26(S1):3503-3508.

    [38] 金 毅,宋慧波,潘結(jié)南,等.煤微觀結(jié)構(gòu)三維表征及其孔-滲時空演化模式數(shù)值分析[J].巖石力學(xué)與工程學(xué)報,2013, 32(S1):2632-2641.

    Jin Yi,Song Huibo,Pan Jienan,et al.Three-dimensional representation of coal’s microstructure and numerical analysis of its porepermeability spatial-temporal evolution mode[J].Chinese Journal of Rock Mechanics and Engineering,2013,32(S1):2632-2641.

    [39] Jin Y,Song H B,Hu B,et al.Lattice Boltzmann simulation of fluid flow through coal reservoir’s fractal pore structure[J].Science China Earth Sciences,2013,56(9):1519-1530.

    [40] Glover P W J,Matsuki K,Hikima R,et al.Synthetic rough fractures in rocks[J].Journal of Geophysical Research,1998,103(B5): 9609-9620.

    [41] Mandelbrot B B.The fractal geometry of nature[M].New York: Macmillan,1983.

    [42] Zhang C B,Chen Y P,Deng Z L,et al.Role of rough surface topography on gas slip flow in microchannels[J].Physical Review E, 2012,86(1):16319.

    [43] Chen Y P,Zhang C B,Shi M H,et al.Role of surface roughness characterized by fractal geometry on laminar flow in microchannels [J].Physical Review E,2009,80(2):26301.

    [44] Auradou H,Drazer G,Boschan A,et al.Flow channeling in a single fracture induced by shear displacement[J].Geothermics,2006, 35(5-6):576-588.

    [45] Madadi M,Sahimi M.Lattice Boltzmann simulation of fluid flow in fracture networks with rough,self-affine surfaces[J].Physical Review E,2003,67(2):26309.

    [46] Brown S R.Simple mathematical model of a rough fracture[J].Journal of Geophysical Research,1995,100(B4):5941-5952.

    [47] Dunweg B,Schiller U D,Ladd A J C.Statistical mechanics of the fluctuating lattice Boltzmann equation[J].Physical Review E, 2007,76(3):36701-36704.

    [48] Qian Y H,d'Humieres D,Lallemand P.Lattice BGK models for Navier-Stokes equation[J].Europhysics Letters,1992,17(6): 479.

    [49] Sukop M C,Thorne D T.Lattice Boltzmann modeling:An introduction for geoscientists and engineers[M].Berlin:Springer,2006.

    [50] Jin Y,Dong J B,Li X,et al.Kinematical measurement of hydraulic tortuosity of fluid flow in porous media[J].International Journal of Modern Physics C,2015,26(2):1550017.

    Numerical investigation of migration mechanism for coal-bed methane flow through cleats with rough surfaces in coal reservoir

    JIN Yi1,2,ZHU Yi-bo1,WU Ying1,ZHENG Jun-ling1,DONG Jia-bin1,LI Xiang1

    (1.School of Resource and Environment,Henan Polytechnic University,Jiaozuo 454003,China;2.Collaborative Innovation Center of Coalbed Methane and Shale Gas for Central Plains Economic Region,Henan Province,Jiaozuo 454003,China)

    To estimate the productivity of coal-bed methane(CBM),the premise is the effective description of CBM migration mechanism.The authors have developed a new fracture-permeability relative model by coupling the general Kozeny-Carman equation and the classical cubic law of fracture flow,which takes into account the effects of external and internal friction.In the model developed,the hydraulic tortuosity is responsible for the internal friction,and the surface tortuosity for external friction effect,which is defined as the area ratio between cleat rough surface and the projected one on the macro-flow plane.Consequently,the authors modelled the cleats composed of rough surfaces characterized by the theory of fractal geometry,and simulated the CBM flow through them by lattice Boltzmann method.In addition,the authors investigated the effects,such as surface fractal dimension,absolute surface roughness,as well as the relative surface roughness,on the CBM flow.Based on the results obtained,the validity of the newly derived fracturepermeability relationship is verified by comparative analyses.The results reveal that the newly developed model accounting for the effects of external and internal friction,can effectively describe the migration law of CBM through rough cleats.

    coal reservoir;rough cleat;coal-bed methane;cubic law of fracture flow;Kozeny-Carman equation;lattice Boltzmann method

    P618.11

    A

    0253-9993(2014)09-1826-09

    2014-05-21 責(zé)任編輯:常 琛

    國家自然科學(xué)基金資助項目(41102093,41472128);山西省煤層氣聯(lián)合基金資助項目(2012012002)

    金 毅(1979—),男,湖北鄂州人,副教授,博士。Tel:0391-8316286,E-mail:jinyi2005@hpu.edu.cn

    金 毅,祝一搏,吳 影,等.煤儲層粗糙割理中煤層氣運(yùn)移機(jī)理數(shù)值分析[J].煤炭學(xué)報,2014,39(9):1826-1834.

    10.13225/ j.cnki.jccs.2014.8003

    Jin Yi,Zhu Yibo,Wu Ying,et al.Numerical investigation of migration mechanism for coal-bed methane flow through cleats with rough surfaces in coal reservoir[J].Journal of China Coal Society,2014,39(9):1826-1834.doi:10.13225/j.cnki.jccs.2014.8003

    猜你喜歡
    運(yùn)移格子維數(shù)
    β-變換中一致丟番圖逼近問題的維數(shù)理論
    曲流河復(fù)合點壩砂體構(gòu)型表征及流體運(yùn)移機(jī)理
    一類齊次Moran集的上盒維數(shù)
    東營凹陷北帶中淺層油氣運(yùn)移通道組合類型及成藏作用
    數(shù)格子
    填出格子里的數(shù)
    格子間
    女友(2017年6期)2017-07-13 11:17:10
    關(guān)于齊次Moran集的packing維數(shù)結(jié)果
    格子龍
    涉及相變問題Julia集的Hausdorff維數(shù)
    美女高潮的动态| 建设人人有责人人尽责人人享有的 | 一个人观看的视频www高清免费观看| 男人和女人高潮做爰伦理| 亚洲精品456在线播放app| 亚洲av成人精品一区久久| 伦理电影大哥的女人| 少妇人妻一区二区三区视频| 国产精品爽爽va在线观看网站| 亚洲精品影视一区二区三区av| 亚洲精品中文字幕在线视频 | 高清日韩中文字幕在线| 天天躁日日操中文字幕| 狠狠精品人妻久久久久久综合| 亚洲国产日韩欧美精品在线观看| 1000部很黄的大片| 午夜激情久久久久久久| 亚洲精品日本国产第一区| 高清在线视频一区二区三区| h日本视频在线播放| 亚洲精品,欧美精品| 色综合站精品国产| 乱系列少妇在线播放| 国产免费又黄又爽又色| 亚洲欧美一区二区三区黑人 | 成人欧美大片| ponron亚洲| xxx大片免费视频| 久久韩国三级中文字幕| 内地一区二区视频在线| 伦理电影大哥的女人| 3wmmmm亚洲av在线观看| 两个人视频免费观看高清| 麻豆av噜噜一区二区三区| 18+在线观看网站| 中文天堂在线官网| 丰满人妻一区二区三区视频av| 嫩草影院入口| av线在线观看网站| 国产精品一区二区在线观看99 | a级毛色黄片| 波多野结衣巨乳人妻| 国内少妇人妻偷人精品xxx网站| 国产精品人妻久久久影院| 亚洲综合精品二区| 亚洲av在线观看美女高潮| 亚洲色图av天堂| 黄片wwwwww| 午夜亚洲福利在线播放| 99久国产av精品| 亚洲在久久综合| 国产黄频视频在线观看| 超碰av人人做人人爽久久| 久久久久性生活片| 久久久久免费精品人妻一区二区| 日韩欧美精品v在线| 国产综合精华液| 国产亚洲一区二区精品| 十八禁网站网址无遮挡 | 日韩 亚洲 欧美在线| 丝瓜视频免费看黄片| 国产欧美日韩精品一区二区| 老司机影院毛片| 国产黄频视频在线观看| 真实男女啪啪啪动态图| 在线播放无遮挡| 欧美日韩一区二区视频在线观看视频在线 | 水蜜桃什么品种好| 中文字幕制服av| 精品人妻一区二区三区麻豆| 亚洲三级黄色毛片| 精品国内亚洲2022精品成人| 国产在线男女| 天天一区二区日本电影三级| 国产精品一区二区在线观看99 | 色综合亚洲欧美另类图片| 国产黄a三级三级三级人| 久久久久久国产a免费观看| 日韩成人av中文字幕在线观看| 欧美成人一区二区免费高清观看| 久久久久精品性色| 亚洲精品中文字幕在线视频 | 国产黄片美女视频| 草草在线视频免费看| 嫩草影院入口| 国产高清三级在线| 三级国产精品欧美在线观看| 男女边吃奶边做爰视频| 青春草国产在线视频| 日韩成人伦理影院| 日本色播在线视频| av卡一久久| 99久国产av精品| 国产精品.久久久| 国产一区二区三区av在线| 校园人妻丝袜中文字幕| 丝瓜视频免费看黄片| 人体艺术视频欧美日本| 2018国产大陆天天弄谢| 国产中年淑女户外野战色| 国产成人午夜福利电影在线观看| 亚洲av中文av极速乱| 99九九线精品视频在线观看视频| 国产精品蜜桃在线观看| av福利片在线观看| 七月丁香在线播放| 亚洲精品日本国产第一区| 80岁老熟妇乱子伦牲交| 国产真实伦视频高清在线观看| 成人午夜高清在线视频| 乱码一卡2卡4卡精品| 日本免费在线观看一区| 亚州av有码| 精品欧美国产一区二区三| 联通29元200g的流量卡| 国产精品福利在线免费观看| 在线观看免费高清a一片| 极品少妇高潮喷水抽搐| av网站免费在线观看视频 | 国产一区有黄有色的免费视频 | 一区二区三区免费毛片| 免费少妇av软件| 22中文网久久字幕| 久久久精品欧美日韩精品| 亚洲成人一二三区av| 日韩国内少妇激情av| 又大又黄又爽视频免费| 美女cb高潮喷水在线观看| 男女那种视频在线观看| 日韩av在线免费看完整版不卡| 欧美一区二区亚洲| 亚洲不卡免费看| 欧美日韩一区二区视频在线观看视频在线 | 欧美日韩精品成人综合77777| 又大又黄又爽视频免费| 久久久久网色| 久久久欧美国产精品| 国产午夜精品久久久久久一区二区三区| 久久精品久久精品一区二区三区| 国产色爽女视频免费观看| 国产永久视频网站| 一个人免费在线观看电影| 欧美成人午夜免费资源| 麻豆久久精品国产亚洲av| 狂野欧美激情性xxxx在线观看| 免费看日本二区| 亚洲成人av在线免费| 中文精品一卡2卡3卡4更新| 狂野欧美白嫩少妇大欣赏| 搞女人的毛片| 最近最新中文字幕免费大全7| 男人爽女人下面视频在线观看| 免费看不卡的av| 哪个播放器可以免费观看大片| 毛片一级片免费看久久久久| av在线老鸭窝| 只有这里有精品99| av国产久精品久网站免费入址| 联通29元200g的流量卡| 美女主播在线视频| 日本欧美国产在线视频| 久久久久免费精品人妻一区二区| 高清日韩中文字幕在线| 少妇熟女欧美另类| 网址你懂的国产日韩在线| 777米奇影视久久| 国产精品一二三区在线看| 少妇高潮的动态图| 高清av免费在线| 日本黄色片子视频| 国产av码专区亚洲av| 男女国产视频网站| 亚洲av在线观看美女高潮| 国产女主播在线喷水免费视频网站 | 亚洲av成人精品一区久久| freevideosex欧美| 亚洲婷婷狠狠爱综合网| 80岁老熟妇乱子伦牲交| 日韩制服骚丝袜av| 国产美女午夜福利| 看十八女毛片水多多多| 夫妻午夜视频| 久久精品国产亚洲av天美| 国产女主播在线喷水免费视频网站 | 精品久久久噜噜| 久久99精品国语久久久| 最后的刺客免费高清国语| 亚洲av中文字字幕乱码综合| 亚洲精品成人久久久久久| 高清欧美精品videossex| 国产亚洲精品av在线| 国产精品av视频在线免费观看| 777米奇影视久久| 少妇人妻精品综合一区二区| av福利片在线观看| 舔av片在线| 国产在线一区二区三区精| 22中文网久久字幕| 好男人视频免费观看在线| 一级毛片aaaaaa免费看小| 五月伊人婷婷丁香| 国产精品嫩草影院av在线观看| 91狼人影院| 亚洲综合色惰| 国产视频内射| 亚洲第一区二区三区不卡| 成人漫画全彩无遮挡| 高清在线视频一区二区三区| 神马国产精品三级电影在线观看| 97热精品久久久久久| 少妇的逼水好多| 又大又黄又爽视频免费| 亚洲av免费高清在线观看| 国产在线男女| 日韩大片免费观看网站| 色视频www国产| 日本免费在线观看一区| 91在线精品国自产拍蜜月| 欧美一级a爱片免费观看看| 91av网一区二区| 国产成人a∨麻豆精品| 久久久久久久久久久丰满| 成人午夜高清在线视频| 国产老妇伦熟女老妇高清| 高清毛片免费看| 精品久久久久久久人妻蜜臀av| 在线a可以看的网站| 精品人妻视频免费看| 性插视频无遮挡在线免费观看| 日韩一区二区视频免费看| 久久精品夜夜夜夜夜久久蜜豆| 乱系列少妇在线播放| 亚洲成色77777| 噜噜噜噜噜久久久久久91| 日日啪夜夜爽| 亚洲欧美精品专区久久| 精品一区二区三区视频在线| 99久久九九国产精品国产免费| 亚洲欧洲日产国产| 精品久久久久久久久久久久久| 亚洲精品日韩在线中文字幕| 国产亚洲av嫩草精品影院| 欧美+日韩+精品| av又黄又爽大尺度在线免费看| 国产精品伦人一区二区| 中文天堂在线官网| 久久久久九九精品影院| 免费大片黄手机在线观看| 日韩av免费高清视频| 丝瓜视频免费看黄片| 日产精品乱码卡一卡2卡三| 午夜福利视频1000在线观看| 女人久久www免费人成看片| 亚洲内射少妇av| 亚洲久久久久久中文字幕| 精品人妻偷拍中文字幕| 日本一二三区视频观看| 亚洲色图av天堂| 听说在线观看完整版免费高清| 丝袜喷水一区| 亚洲美女视频黄频| 久热久热在线精品观看| 汤姆久久久久久久影院中文字幕 | 夜夜爽夜夜爽视频| 精品久久久久久电影网| 极品教师在线视频| 久久99精品国语久久久| 国产黄色视频一区二区在线观看| 熟女电影av网| 午夜福利视频精品| 国产成人freesex在线| 亚洲欧洲日产国产| 国产色婷婷99| 色吧在线观看| 九九爱精品视频在线观看| 午夜福利在线在线| 国产美女午夜福利| 1000部很黄的大片| 一级毛片我不卡| 十八禁网站网址无遮挡 | 国产精品一区二区三区四区免费观看| 一级毛片久久久久久久久女| 国产淫语在线视频| 色吧在线观看| 久久精品国产鲁丝片午夜精品| 亚洲成人精品中文字幕电影| 非洲黑人性xxxx精品又粗又长| 日本免费在线观看一区| 成人无遮挡网站| 久久精品国产亚洲av涩爱| 久久久精品94久久精品| 日韩一区二区视频免费看| 国产 亚洲一区二区三区 | 日本黄色片子视频| 男女边摸边吃奶| 日本一二三区视频观看| 精品国产一区二区三区久久久樱花 | 亚洲色图av天堂| 超碰av人人做人人爽久久| a级一级毛片免费在线观看| 在线免费观看的www视频| 欧美成人午夜免费资源| 精华霜和精华液先用哪个| 三级毛片av免费| 国产免费福利视频在线观看| 一区二区三区乱码不卡18| 色综合色国产| 丝袜美腿在线中文| 看免费成人av毛片| 日韩成人av中文字幕在线观看| 午夜免费男女啪啪视频观看| 亚洲欧洲日产国产| 午夜免费观看性视频| 国产欧美另类精品又又久久亚洲欧美| 国产综合懂色| 久久久久久久久大av| 欧美日韩亚洲高清精品| 黄片wwwwww| 免费看av在线观看网站| 中文欧美无线码| 国产麻豆成人av免费视频| 国产免费一级a男人的天堂| 精品久久久精品久久久| 成人亚洲精品一区在线观看 | 97人妻精品一区二区三区麻豆| 麻豆乱淫一区二区| 婷婷色麻豆天堂久久| 精品一区二区三区视频在线| 边亲边吃奶的免费视频| 国产精品久久久久久精品电影小说 | 精品久久久久久久久久久久久| 黄片wwwwww| 一级av片app| 亚洲精品久久久久久婷婷小说| 国产精品美女特级片免费视频播放器| 男的添女的下面高潮视频| 日韩亚洲欧美综合| 亚洲av日韩在线播放| 毛片女人毛片| 天堂中文最新版在线下载 | 精品人妻熟女av久视频| 男的添女的下面高潮视频| 国产精品熟女久久久久浪| 熟女电影av网| 亚洲av免费高清在线观看| 夜夜看夜夜爽夜夜摸| 最后的刺客免费高清国语| 波野结衣二区三区在线| 一本久久精品| 亚洲欧美一区二区三区黑人 | 一级毛片aaaaaa免费看小| 国产精品国产三级国产av玫瑰| 又爽又黄a免费视频| 在线观看一区二区三区| 亚洲人成网站在线播| av免费观看日本| 啦啦啦啦在线视频资源| av免费观看日本| 嘟嘟电影网在线观看| 中文欧美无线码| 三级国产精品片| 日韩精品青青久久久久久| 欧美xxxx性猛交bbbb| 午夜福利高清视频| 51国产日韩欧美| 国产永久视频网站| 少妇熟女欧美另类| 91aial.com中文字幕在线观看| 建设人人有责人人尽责人人享有的 | 亚洲最大成人av| 亚洲不卡免费看| 日韩av在线免费看完整版不卡| 亚洲精品影视一区二区三区av| 又爽又黄a免费视频| 日韩一本色道免费dvd| av在线播放精品| 亚洲av成人精品一区久久| 亚洲成人中文字幕在线播放| 日韩亚洲欧美综合| 免费大片黄手机在线观看| 国产av不卡久久| 欧美 日韩 精品 国产| 久久精品人妻少妇| 秋霞在线观看毛片| 日韩欧美精品v在线| 一级av片app| 日本与韩国留学比较| 久久精品久久久久久久性| 2021天堂中文幕一二区在线观| 国产视频内射| 日韩中字成人| 蜜桃久久精品国产亚洲av| 国产精品av视频在线免费观看| 国产在线一区二区三区精| 美女国产视频在线观看| 国产成人a区在线观看| 中文精品一卡2卡3卡4更新| 九九爱精品视频在线观看| 国产极品天堂在线| 最近中文字幕高清免费大全6| av线在线观看网站| 你懂的网址亚洲精品在线观看| 国产 一区精品| 婷婷色综合www| 亚洲电影在线观看av| 久久国内精品自在自线图片| 纵有疾风起免费观看全集完整版 | 中国美白少妇内射xxxbb| 久久精品综合一区二区三区| 免费av观看视频| 我要看日韩黄色一级片| 国产探花极品一区二区| 少妇猛男粗大的猛烈进出视频 | 久久久成人免费电影| 高清在线视频一区二区三区| 国产成人精品婷婷| 日韩电影二区| av在线天堂中文字幕| 国产成人freesex在线| 国产成人a∨麻豆精品| 人体艺术视频欧美日本| 在线免费十八禁| 午夜日本视频在线| 如何舔出高潮| 一级毛片 在线播放| 精品久久久噜噜| 久久久久久国产a免费观看| 久久亚洲国产成人精品v| 日本爱情动作片www.在线观看| 国产精品.久久久| 午夜老司机福利剧场| av.在线天堂| 国产精品一区www在线观看| 国产av码专区亚洲av| 亚洲精品影视一区二区三区av| 婷婷色av中文字幕| 色尼玛亚洲综合影院| 婷婷色综合大香蕉| 高清欧美精品videossex| 成人一区二区视频在线观看| 国产精品一及| 国产在视频线在精品| 美女主播在线视频| 99热6这里只有精品| 国产精品久久久久久久电影| av一本久久久久| 精品久久久噜噜| 美女黄网站色视频| 国产精品伦人一区二区| freevideosex欧美| 老司机影院成人| 色网站视频免费| 深爱激情五月婷婷| 青春草亚洲视频在线观看| 国产一级毛片在线| 人妻夜夜爽99麻豆av| 亚洲欧美清纯卡通| www.av在线官网国产| 有码 亚洲区| 又黄又爽又刺激的免费视频.| 91av网一区二区| 久久精品国产亚洲网站| 亚洲电影在线观看av| 毛片一级片免费看久久久久| 熟妇人妻久久中文字幕3abv| 国产黄a三级三级三级人| 男女那种视频在线观看| 国产精品av视频在线免费观看| 国产一级毛片在线| 人人妻人人澡欧美一区二区| 精品人妻熟女av久视频| 国产黄片美女视频| 观看美女的网站| 久久亚洲国产成人精品v| 国产精品日韩av在线免费观看| 成人综合一区亚洲| 色吧在线观看| 精品久久久久久久久久久久久| 男女下面进入的视频免费午夜| 成年免费大片在线观看| 男人舔奶头视频| 精品一区二区三区视频在线| 精品国产一区二区三区久久久樱花 | 日本猛色少妇xxxxx猛交久久| 亚洲国产成人一精品久久久| 国产真实伦视频高清在线观看| 国产视频内射| 2022亚洲国产成人精品| 五月天丁香电影| 插逼视频在线观看| 69人妻影院| 日韩一区二区视频免费看| 国产一区亚洲一区在线观看| av黄色大香蕉| 国产精品久久久久久精品电影| 99热这里只有精品一区| 精品久久久久久久久亚洲| 晚上一个人看的免费电影| 黄片无遮挡物在线观看| 亚洲精品成人av观看孕妇| 成人毛片60女人毛片免费| 99久久精品国产国产毛片| 亚洲国产高清在线一区二区三| 久久久久九九精品影院| 亚洲人成网站在线观看播放| 一边亲一边摸免费视频| 精品一区二区三卡| av又黄又爽大尺度在线免费看| 丝袜喷水一区| 97超视频在线观看视频| 精品久久久久久久久久久久久| 欧美bdsm另类| 午夜激情福利司机影院| 国产一区亚洲一区在线观看| 国产乱人偷精品视频| 亚洲伊人久久精品综合| 建设人人有责人人尽责人人享有的 | 亚州av有码| 免费少妇av软件| 国产成人免费观看mmmm| 亚洲精品乱码久久久v下载方式| 少妇被粗大猛烈的视频| 韩国高清视频一区二区三区| 日韩av不卡免费在线播放| 乱系列少妇在线播放| 永久免费av网站大全| 精品欧美国产一区二区三| 日本黄大片高清| 在现免费观看毛片| 三级国产精品欧美在线观看| 最近手机中文字幕大全| 国产精品一区二区三区四区免费观看| 两个人视频免费观看高清| 免费看不卡的av| 亚洲最大成人av| 高清视频免费观看一区二区 | 人妻一区二区av| 欧美日韩视频高清一区二区三区二| 国产亚洲精品av在线| 18禁在线播放成人免费| 欧美日本视频| 日本爱情动作片www.在线观看| 久久鲁丝午夜福利片| 卡戴珊不雅视频在线播放| 久久鲁丝午夜福利片| 欧美zozozo另类| 看十八女毛片水多多多| 欧美xxxx黑人xx丫x性爽| 久久久欧美国产精品| 少妇的逼好多水| 波多野结衣巨乳人妻| 深爱激情五月婷婷| 亚洲最大成人中文| 69av精品久久久久久| 高清午夜精品一区二区三区| 毛片一级片免费看久久久久| 特级一级黄色大片| 观看美女的网站| 人妻一区二区av| 色综合站精品国产| 99热这里只有是精品50| 色综合亚洲欧美另类图片| 狠狠精品人妻久久久久久综合| 五月伊人婷婷丁香| av线在线观看网站| 亚洲怡红院男人天堂| 久久久久久久大尺度免费视频| 国产午夜精品久久久久久一区二区三区| 国产美女午夜福利| 午夜福利高清视频| 水蜜桃什么品种好| 神马国产精品三级电影在线观看| 简卡轻食公司| 午夜福利网站1000一区二区三区| 乱系列少妇在线播放| 日韩一本色道免费dvd| 国产黄色视频一区二区在线观看| 午夜福利高清视频| 国产麻豆成人av免费视频| 亚洲激情五月婷婷啪啪| 免费黄色在线免费观看| 天天躁夜夜躁狠狠久久av| 午夜激情福利司机影院| 婷婷色麻豆天堂久久| 美女高潮的动态| 男人舔女人下体高潮全视频| 日本一本二区三区精品| 精品久久国产蜜桃| 91狼人影院| 久久97久久精品| 国产av码专区亚洲av| 成人av在线播放网站| 午夜福利在线在线| 男人狂女人下面高潮的视频| 精品国产三级普通话版| av黄色大香蕉| 国产成人福利小说| 亚洲,欧美,日韩| av卡一久久| 国产精品一区二区性色av| 久久精品国产亚洲网站| 国产一区二区在线观看日韩| 亚洲图色成人| 中文精品一卡2卡3卡4更新| 91精品一卡2卡3卡4卡| 欧美日韩视频高清一区二区三区二| 免费大片黄手机在线观看| 国模一区二区三区四区视频| 亚洲欧美一区二区三区黑人 | 性色avwww在线观看| 纵有疾风起免费观看全集完整版 | 亚洲三级黄色毛片| 亚洲最大成人av| 大香蕉久久网| 久久午夜福利片|