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

    砂巖核磁共振響應(yīng)模擬及受限擴(kuò)散

    2016-07-28 09:31:58郭江峰謝然紅鄒友龍
    地球物理學(xué)報(bào) 2016年7期
    關(guān)鍵詞:砂巖

    郭江峰, 謝然紅*, 鄒友龍

    1 油氣資源與探測國家重點(diǎn)實(shí)驗(yàn)室,中國石油大學(xué)(北京), 北京 102249 2 地球探測與信息技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室,中國石油大學(xué)(北京), 北京 102249

    ?

    砂巖核磁共振響應(yīng)模擬及受限擴(kuò)散

    郭江峰1,2, 謝然紅1,2*, 鄒友龍1,2

    1 油氣資源與探測國家重點(diǎn)實(shí)驗(yàn)室,中國石油大學(xué)(北京), 北京1022492 地球探測與信息技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室,中國石油大學(xué)(北京), 北京102249

    摘要本文運(yùn)用隨機(jī)游走方法模擬了砂巖儲(chǔ)層中流體的核磁共振(NMR)響應(yīng)及其受限擴(kuò)散現(xiàn)象.通過改變數(shù)字巖心的分辨率模擬生成不同孔隙尺寸的砂巖,研究了不同孔隙尺寸砂巖飽含水時(shí)流體擴(kuò)散系數(shù)隨擴(kuò)散時(shí)間的變化關(guān)系,同時(shí)模擬了砂巖飽和單相流體和兩相流體的NMR響應(yīng);研究了流體的受限擴(kuò)散系數(shù)與橫向弛豫時(shí)間T2的關(guān)系,分析了表面弛豫率和膠結(jié)指數(shù)對(duì)潤濕相流體受限擴(kuò)散系數(shù)線位置的影響,并將其用于解釋砂巖儲(chǔ)層的D-T2分布.結(jié)果表明:孔隙流體的擴(kuò)散系數(shù)會(huì)隨擴(kuò)散時(shí)間的增加而逐漸減小并趨于定值.隨著巖石孔隙尺寸的減小,受限擴(kuò)散現(xiàn)象越明顯,受限擴(kuò)散對(duì)巖石NMR響應(yīng)的影響也越大.潤濕相流體受限擴(kuò)散系數(shù)線的位置受巖石膠結(jié)指數(shù)和表面弛豫率的影響較大.由于潤濕相流體擴(kuò)散系數(shù)減小,導(dǎo)致D-T2分布中潤濕相流體信號(hào)偏離其自由擴(kuò)散系數(shù)線,需要利用流體的受限擴(kuò)散系數(shù)線準(zhǔn)確識(shí)別D-T2分布中的潤濕相流體.

    關(guān)鍵詞砂巖; NMR響應(yīng); 隨機(jī)游走方法; 受限擴(kuò)散; 擴(kuò)散系數(shù)線

    1引言

    隨著石油工業(yè)的發(fā)展,勘探開發(fā)對(duì)象從常規(guī)儲(chǔ)層逐漸轉(zhuǎn)變?yōu)榉浅R?guī)儲(chǔ)層,其中重要的一個(gè)特點(diǎn)是儲(chǔ)層巖石的孔隙尺寸逐漸變小,其核磁共振(NMR)響應(yīng)也會(huì)發(fā)生變化,研究不同孔隙尺寸巖石的NMR響應(yīng)具有重要意義.近年發(fā)展的數(shù)字巖心技術(shù)為巖石的NMR響應(yīng)研究提供了經(jīng)濟(jì)有效的解決途徑.巖石NMR響應(yīng)模擬可以采用隨機(jī)游走方法(Toumelin et al.,2003, 2004,2007;Jin et al.,2009;Talabi et al.,2009;蔡淑惠等,2009;肖立志等,2012;Alghamdi et al.,2013;成家杰等,2013;鄒友龍,2013;Tan et al.,2014),有限元方法(Hagsl?tt et al.,2003)和有限差分方法(Zientara and Freed,1980),而隨機(jī)游走方法具有更高的靈活性且比較容易實(shí)現(xiàn),最重要的是該方法能夠適用于復(fù)雜的巖石孔隙結(jié)構(gòu),所以本文選用了隨機(jī)游走方法模擬巖石的NMR響應(yīng).

    巖石孔隙中的流體由于受到孔隙壁的影響,其擴(kuò)散系數(shù)隨擴(kuò)散時(shí)間變化而變化(Valfouskaya et al.,2006;Luo et al.,2015;張宗富,2013),具體表現(xiàn)為隨著擴(kuò)散時(shí)間的增加,擴(kuò)散系數(shù)逐漸減小并趨于一個(gè)常數(shù)值.Minh等(2003)定性地描述了受限擴(kuò)散對(duì)D-T2分布中水信號(hào)的影響,受限擴(kuò)散使D-T2分布中水信號(hào)向擴(kuò)散系數(shù)減小的方向移動(dòng),此時(shí)水信號(hào)與油信號(hào)容易混淆,給D-T2分布的解釋帶來困難.Zielinski等(2010)研究了受限擴(kuò)散對(duì)D-T2分布計(jì)算碳酸鹽巖飽和度的影響,根據(jù)計(jì)算的流體受限擴(kuò)散系數(shù)線可以準(zhǔn)確識(shí)別流體及獲取流體的飽和度.Minh等(2012)將該方法擴(kuò)展應(yīng)用于受限擴(kuò)散更明顯的頁巖儲(chǔ)層D-T2分布中,可以準(zhǔn)確識(shí)別頁巖儲(chǔ)層中的流體.

    本文應(yīng)用隨機(jī)游走方法模擬砂巖的NMR響應(yīng),通過改變數(shù)字巖心分辨率來模擬不同孔隙尺寸的砂巖,并研究砂巖中流體擴(kuò)散系數(shù)隨擴(kuò)散時(shí)間的變化,探究受限擴(kuò)散和含水飽和度對(duì)砂巖NMR響應(yīng)的影響.研究了流體的擴(kuò)散系數(shù)與橫向弛豫時(shí)間T2的變化關(guān)系,分析了影響因素,并驗(yàn)證了隨機(jī)游走模擬得到的結(jié)果.本文的研究結(jié)果對(duì)非常規(guī)儲(chǔ)層流體D-T2分布解釋有重要的指導(dǎo)意義.

    2NMR響應(yīng)數(shù)值模擬基礎(chǔ)

    2.1理論基礎(chǔ)

    在飽和流體介質(zhì)中,磁化矢量的衰減遵循Bloch-Torrey方程(Torrey,1956)

    (1)

    邊界條件為

    (2)

    初始時(shí)刻的磁化矢量為

    (3)

    式中ρ為表面弛豫率,m(r,t)表示磁化矢量,γ為質(zhì)子的旋磁比,Bz為靜磁場,T2為橫向弛豫時(shí)間,M(0)是初始時(shí)刻總的磁化強(qiáng)度,VP是總的孔隙體積.

    孔隙流體的橫向弛豫時(shí)間T2由體弛豫時(shí)間、表面弛豫時(shí)間和擴(kuò)散弛豫時(shí)間組成,且流體的體弛豫、表面弛豫和擴(kuò)散弛豫都是同時(shí)發(fā)生的,也是彼此獨(dú)立的(Talabi,2008).由此可知,流體NMR的弛豫信號(hào)強(qiáng)度M(t)隨時(shí)間的變化關(guān)系為

    (4)

    式中MB(t)、MS(t)、MG(t)分別是t時(shí)刻的體弛豫信號(hào)強(qiáng)度、表面弛豫信號(hào)強(qiáng)度、擴(kuò)散弛豫信號(hào)強(qiáng)度,M0是初始時(shí)刻總的弛豫信號(hào)強(qiáng)度.

    隨機(jī)游走方法通過模擬大量粒子的布朗運(yùn)動(dòng)求解方程(1)來實(shí)現(xiàn)流體磁化強(qiáng)度衰減的模擬.由(4)式可知,模擬NMR信號(hào)強(qiáng)度的衰減可以通過分別模擬體弛豫信號(hào)強(qiáng)度、表面弛豫信號(hào)強(qiáng)度和擴(kuò)散弛豫信號(hào)強(qiáng)度的衰減來實(shí)現(xiàn),其對(duì)應(yīng)的方法步驟如圖1所示.

    圖1 NMR隨機(jī)游走方法流程圖Fig.1 Flowchart of NMR random-walk method

    (5)

    沒有被吸收的粒子將在界面上發(fā)生彈性碰撞,并保持相位和幅度的連續(xù)變化.從而整個(gè)數(shù)值模擬的過程中,表面弛豫信號(hào)的強(qiáng)度可以表示為

    (6)

    式中N(t)為t時(shí)刻沒有被吸收的粒子總數(shù),N(0)是初始時(shí)刻的粒子總數(shù).

    (7)

    式中Normal()為高斯隨機(jī)函數(shù).為了符合CPMG(Carr-Purcell-Meiboom-Gill)脈沖序列采集要求,當(dāng)t=(n+1/2)TE時(shí),相位發(fā)生反轉(zhuǎn),即φ(t)=-φ(t);當(dāng)t=nTE時(shí),采集CPMG自旋回波信號(hào).此時(shí)擴(kuò)散弛豫信號(hào)強(qiáng)度是所有自旋散相的余弦和,即

    (8)

    2.2方法驗(yàn)證

    選取球形孔隙驗(yàn)證隨機(jī)游走方法的有效性.首先構(gòu)造球形孔隙,如圖2所示,正方體內(nèi)部的球?yàn)榭紫读黧w部分,其余的為顆粒骨架.對(duì)于球形孔隙,磁化強(qiáng)度衰減可表達(dá)為

    圖2 球形孔隙模型Fig.2 Sperical pore model

    (9)式中M0為初始磁化強(qiáng)度,設(shè)為1.假定孔隙半徑r為2 μm,體弛豫時(shí)間T2B為3.0 s,表面弛豫率ρ為30 μm·s-1,磁場梯度G為0.3 T·m-1,TE為1 ms,流體擴(kuò)散系數(shù)D為2.5×10-9m2·s-1.由(9)式計(jì)算得到的理論磁化強(qiáng)度衰減曲線和由隨機(jī)游走方法模擬得到的磁化強(qiáng)度衰減曲線如圖3所示,可以看出結(jié)果一致性較好,驗(yàn)證了方法的有效性.

    圖3 球形孔隙中隨機(jī)游走方法模擬的 磁化強(qiáng)度衰減曲線與理論值的對(duì)比Fig.3 Comparison of magnetization decay curve in an sperical pore obtained from random-walk simulation with theoretical curve

    3砂巖NMR響應(yīng)數(shù)值模擬

    3.1單相流

    本文所用的砂巖為英國帝國理工學(xué)院已公開的砂巖S,孔隙度為16.9%,分辨率為9.1 μm,體素為3003.模擬時(shí)改變巖心的分辨率,數(shù)值分別為9.1 μm、0.91 μm和0.091 μm,這樣可以在不改變孔隙度的情形下,改變巖心樣品的孔隙尺寸.假設(shè)砂巖為水潤濕,且孔隙中完全含水,其數(shù)字巖心模型如圖4a所示,黑色部分表示顆粒骨架,紅色部分表示孔隙中的水;水的擴(kuò)散系數(shù)為2.5×10-9m2·s-1,體弛豫時(shí)間為3.0 s,砂巖的表面弛豫率ρ為30 μm·s-1.

    圖4 砂巖數(shù)字巖心模型 (a) 含水飽和度為100%的砂巖; (b) 含水飽和度75%和含油飽和度25%的砂巖; (c) 含水飽和度50%和含油飽和度50%的砂巖; (d) 含水飽和度25%和含油飽和度75%的砂巖.Fig.4 Digital core models of sandstone (a) Sandstone saturated with 100% water; (b) Sandstone saturated with 75% water and 25% oil; (c) Sandstone saturated with 50% water and 50% oil; (d) Sandstone saturated with 25% water and 75% oil.

    在0.3 T·m-1的梯度場中,用隨機(jī)游走方法模擬了不同分辨率下砂巖中水的擴(kuò)散系數(shù)隨擴(kuò)散時(shí)間的變化,如圖5所示.從圖中可以看出: (1) 隨著擴(kuò)散時(shí)間的增加,砂巖中水的擴(kuò)散系數(shù)逐漸減小,并趨于定值,該定值隨著分辨率的提高逐漸減小; (2) 經(jīng)過相同的擴(kuò)散時(shí)間后,隨著分辨率的提高(砂巖孔隙尺寸逐漸減小),孔隙流體擴(kuò)散系數(shù)逐漸減小,即受限擴(kuò)散現(xiàn)象越明顯.

    用隨機(jī)游走方法模擬了回波間隔TE=0.6 ms時(shí)不同分辨率情形下砂巖在0.3 T·m-1的梯度場下的CPMG自旋回波串,如圖6a所示,隨著分辨率的提高, 砂巖完全飽和水時(shí)回波串幅度衰減加快.利用奇異值分解法(Dunn and Latorraca,1999)分別對(duì)圖6a中的回波串反演得到各自的T2分布,如圖6b所示,隨著分辨率的提高,砂巖完全飽和水時(shí)T2分布的峰值對(duì)應(yīng)的T2值逐漸減小,說明砂巖孔隙尺寸越來越小.

    圖5 砂巖中水的擴(kuò)散系數(shù)隨擴(kuò)散時(shí)間的變化Fig.5 Changes of water diffusion coefficient with time in sandstone

    圖6 砂巖完全飽和水時(shí)的(a)回波串及(b) T2分布Fig.6 Echo trains (a) and T2 distributions (b) of sandstone fully saturated with water

    用隨機(jī)游走方法模擬回波間隔分別為0.6 ms、1.2 ms、2.4 ms、4.8 ms、9.6 ms、12 ms、20 ms和40 ms時(shí)不同分辨率情形下完全飽和水的砂巖在0.3 T·m-1的梯度場下的CPMG自旋回波串,對(duì)應(yīng)的回波個(gè)數(shù)NE=int(1200 ms/TE),運(yùn)用多回波串聯(lián)合反演方法(謝然紅等,2009)對(duì)8組回波串反演得到不同分辨率下砂巖的D-T2分布,如圖7所示,圖中紅色實(shí)線、藍(lán)色實(shí)線、綠色實(shí)線分別代表氣、水和油的自由擴(kuò)散系數(shù)線(下文D-T2分布中的線條含義也是如此).從圖中可以看出,隨著分辨率的提高,D-T2分布中水信號(hào)逐漸偏離水的自由擴(kuò)散系數(shù)線,且水信號(hào)對(duì)應(yīng)的T2值逐漸減小,由于當(dāng)砂巖孔隙尺寸逐漸減小時(shí),水在孔隙中隨機(jī)運(yùn)動(dòng)時(shí)受到巖石顆粒表面的影響越大,即受限擴(kuò)散現(xiàn)象越來越明顯,導(dǎo)致D-T2分布中水信號(hào)向擴(kuò)散系數(shù)減小方向移動(dòng).

    圖7 砂巖完全飽和水時(shí)的D-T2分布 (a) 分辨率為9.1 μm; (b) 分辨率為0.91 μm; (c) 分辨率為0.091 μm.Fig.7 D-T2 distributions of sandstone fully saturated with water (a) Resolution of 9.1 μm; (b) Resolution of 0.91 μm; (c) Resolution of 0.091 μm.

    3.2兩相流

    不同含水飽和度的砂巖數(shù)字巖心模型如圖4所示,黑色部分表示骨架,紅色部分表示孔隙中的水,藍(lán)色部分表示孔隙中的油.設(shè)水的自由擴(kuò)散系數(shù)為2.5×10-9m2·s-1,體弛豫時(shí)間為3.0 s;而油的自由擴(kuò)散系數(shù)為1×10-10m2·s-1,體弛豫時(shí)間為0.2 s;砂巖為水潤濕,其表面弛豫率ρ=30 μm·s-1,分辨率取0.91 μm;在0.3 T·m-1的梯度場中,設(shè)置回波間隔為0.6 ms、1.2 ms、2.4 ms、4.8 ms、9.6 ms、12 ms、20 ms和40 ms,對(duì)應(yīng)的回波個(gè)數(shù)NE=int(1200 ms/TE),分別對(duì)不同含水飽和度砂巖采用隨機(jī)游走方法模擬得到8組CPMG自旋回波串,運(yùn)用多回波串聯(lián)合反演方法得到其D-T2分布,如圖8所示.從圖中可以看出,隨著含水飽和度的減小,D-T2分布中水信號(hào)偏離水的自由擴(kuò)散系數(shù)線越來越大;而油信號(hào)在油的自由擴(kuò)散系數(shù)線上,其位置不受含油飽和度的影響.這是因?yàn)殡S著含水飽和度的減小,水所占據(jù)孔隙的表面積與體積比值越來越大,粒子運(yùn)動(dòng)過程中與巖石顆粒發(fā)生碰撞的概率越來越大,它們的運(yùn)動(dòng)方向更容易發(fā)生改變,導(dǎo)致受限擴(kuò)散現(xiàn)象越來越明顯;而油由于不接觸巖石顆粒表面,所以不存在受限擴(kuò)散現(xiàn)象.

    圖8 砂巖不同含水飽和度時(shí)D-T2分布(分辨率為0.91 μm) (a) 含水飽和度75%和含油飽和度25%的砂巖; (b) 含水飽和度50%和含油飽和度50%的砂巖; (c) 含水飽和度25%和含油飽和度75%的砂巖.Fig.8 D-T2 distributions of the sandstone with different water saturation (resolution of 0.91 μm) (a) Sandstone saturated with 75% water and 25% oil; (b) Sandstone saturated with 50% water and 50% oil; (c) Sandstone saturated with 25% water and 75% oil.

    4流體的受限擴(kuò)散系數(shù)線

    4.1流體的受限擴(kuò)散系數(shù)線特征

    附錄A推導(dǎo)出潤濕相流體的擴(kuò)散系數(shù)D隨橫向弛豫時(shí)間T2的變化關(guān)系.假設(shè)巖石飽含單相流體(氣、水或油)時(shí),且其潤濕性隨著飽含的單相流體變化而變化,根據(jù)(A9)式分別在其D-T2分布中繪出各自流體的受限擴(kuò)散系數(shù)線,如圖9所示,圖中紅色、藍(lán)色和綠色實(shí)線分別代表氣、水和油的自由擴(kuò)散系數(shù)線;紅色和藍(lán)色虛線分別代表氣和水的受限擴(kuò)散系數(shù)線;由于不同黏度的油對(duì)應(yīng)不同的自由擴(kuò)散系數(shù)和體弛豫時(shí)間,圖中從上到下的三條綠色虛線分別代表1 cp的輕質(zhì)油(D0,oil=5×10-10m2·s-1, T2B,oil=1.0 s)、10 cp的中等黏度油(D0,oil=5×10-11m2·s-1, T2B,oil=0.1 s)和100 cp的稠油(D0,oil=5×10-12m2·s-1, T2B,oil=0.01 s)的受限擴(kuò)散系數(shù)線.可以看出,潤濕相流體的受限擴(kuò)散系數(shù)線在D-T2分布中表現(xiàn)為:隨著T2值的減小,流體受限擴(kuò)散系數(shù)逐漸減小并趨于穩(wěn)定,且不同流體的受限擴(kuò)散系數(shù)線變化趨勢一致.主要原因如下:在只考慮表面弛豫,忽略體弛豫和擴(kuò)散弛豫時(shí),巖石的表面弛豫率是一定的,這樣T2值減小的原因只能是孔隙尺寸的減小.依據(jù)受限擴(kuò)散的機(jī)理可知,隨著孔隙尺寸的減小,自旋粒子擴(kuò)散過程中與巖石顆粒發(fā)生碰撞的概率就會(huì)增大,從而導(dǎo)致單位時(shí)間步長內(nèi)粒子隨機(jī)運(yùn)動(dòng)方向被改變的概率增大,巖石孔隙中流體受限擴(kuò)散現(xiàn)象會(huì)越明顯,表現(xiàn)為受限擴(kuò)散系數(shù)越來越小;當(dāng)孔隙尺寸減小到一定程度時(shí),在單位時(shí)間步長內(nèi)粒子隨機(jī)運(yùn)動(dòng)方向被改變的概率逐漸趨于1,此時(shí)巖石孔隙中流體的受限擴(kuò)散系數(shù)不再減小,而趨于一個(gè)定值.不同類型流體作為潤濕相在相同孔隙尺寸的孔隙中擴(kuò)散機(jī)理是相同的,所以它們的擴(kuò)散系數(shù)隨T2的變化趨勢也一致.

    圖9 流體的受限擴(kuò)散系數(shù)線分布特征Fig.9 Distribution features of fluid restricted diffusion coefficient lines

    4.2影響因素分析

    流體受限擴(kuò)散系數(shù)線的位置不是一成不變的,它與很多因素有關(guān),例如巖石膠結(jié)指數(shù)和表面弛豫率等.下面具體討論巖石膠結(jié)指數(shù)m和表面弛豫率ρ對(duì)水的受限擴(kuò)散系數(shù)線的影響.圖10a中藍(lán)色虛線和紅色虛線分別代表m=2、ρ=3 μm·s-1及m=3、ρ=3 μm·s-1時(shí)水的受限擴(kuò)散系數(shù)線,從圖中可以看出,當(dāng)T2較小時(shí),膠結(jié)指數(shù)m對(duì)孔隙中潤濕相流體(水)受限擴(kuò)散系數(shù)線位置的影響占主導(dǎo)地位,且膠結(jié)指數(shù)m越大,受限擴(kuò)散越顯著,流體受限擴(kuò)散系數(shù)線偏離自由擴(kuò)散系數(shù)線越遠(yuǎn).這是因?yàn)閹r石膠結(jié)指數(shù)m越大,孔隙結(jié)構(gòu)越復(fù)雜,孔隙網(wǎng)絡(luò)彎曲度越大,粒子擴(kuò)散過程中與巖石顆粒發(fā)生碰撞的概率就會(huì)越大,巖石孔隙中流體受限擴(kuò)散就會(huì)越明顯.圖10b中綠色虛線和藍(lán)色虛線分別代表m=2、ρ=30 μm·s-1及m=2、ρ= 3 μm·s-1時(shí)水的受限擴(kuò)散系數(shù)線,從圖中可以看出,當(dāng)膠結(jié)指數(shù)m一定時(shí),即巖石孔隙結(jié)構(gòu)一定,隨著T2的逐漸增加,巖石表面弛豫率ρ不同,孔隙中潤濕相流體(水)受限擴(kuò)散系數(shù)線的位置不同,表現(xiàn)為表面弛豫率ρ越小,水的受限擴(kuò)散系數(shù)線偏離自由擴(kuò)散系數(shù)線越遠(yuǎn).

    圖10 水的受限擴(kuò)散系數(shù)線影響因素分析 (a) 膠結(jié)指數(shù);(b) 表面弛豫率.Fig.10 Analysis on influencing factors for water restricted diffusion coefficient line (a) Cementation index;(b) Surface relaxivity.

    4.3 砂巖D-T2分布中水信號(hào)的識(shí)別

    由第3節(jié)親水砂巖的模擬結(jié)果可知,當(dāng)砂巖孔隙尺寸較小時(shí),水的受限擴(kuò)散使D-T2分布中水信號(hào)的位置向擴(kuò)散系數(shù)減小的方向移動(dòng),隨含水飽和度的減小,D-T2分布中水信號(hào)偏離水的自由擴(kuò)散系數(shù)線越來越大,此時(shí)依據(jù)水的自由擴(kuò)散系數(shù)線識(shí)別水信號(hào)變得困難,下面研究依據(jù)水的受限擴(kuò)散系數(shù)線來識(shí)別水信號(hào).親水砂巖有效表面弛豫率ρeff隨含水飽和度的改變而變化,通過(A7)式計(jì)算出不同含水飽和度時(shí)該砂巖的有效表面弛豫率,如表1所示,隨含水飽和度的減小,有效表面弛豫率ρeff逐漸減小.通過(A9)式計(jì)算分辨率為0.91 μm的砂巖不同含水飽和度時(shí)的受限擴(kuò)散系數(shù)線,如圖11所示,藍(lán)色虛線代表水的受限擴(kuò)散系數(shù)線,可以看出,隨著含水飽和度的減小,藍(lán)色虛線在D軸上的截距越來越小,即受限擴(kuò)散現(xiàn)象越來越明顯;D-T2分布中水信號(hào)落在藍(lán)色虛線上,說明依據(jù)水的受限擴(kuò)散系數(shù)線可以準(zhǔn)確識(shí)別水信號(hào).

    圖11 砂巖D-T2分布中水信號(hào)的識(shí)別 (a) 含水飽和度為100%的砂巖; (b) 含水飽和度為75%和含油飽和度25%的砂巖; (c) 含水飽和度為50%和 含油飽和度50%的砂巖; (d) 含水飽和度25%和含油飽和度75%的砂巖.Fig.11 Water signal identification in D-T2 distributions of sandstone (a) Sandstone saturated with 100% water; (b) Sandstone saturated with 75% water and 25% oil; (c) Sandstone saturated with 50% water and 50% oil; (d) Sandstone saturated with 25% water and 75% oil.

    含水飽和度SW100%75%50%25%有效表面弛豫率ρeff(μm·s-1)30.0022.1813.3110.52

    5結(jié)論

    (1) 通過改變巖心分辨率模擬生成不同孔隙尺寸的砂巖,隨著砂巖孔隙尺寸逐漸變小,砂巖中潤濕相流體受限擴(kuò)散現(xiàn)象越明顯,其擴(kuò)散系數(shù)隨擴(kuò)散時(shí)間的增加逐漸減小并趨于定值,該定值隨著孔隙尺寸的減小而減小.

    (2) 當(dāng)砂巖飽含單相流體時(shí),隨孔隙尺寸的減小,其T2分布峰值逐漸向短弛豫時(shí)間方向移動(dòng),其D-T2分布中流體信號(hào)向擴(kuò)散系數(shù)減小方向以及短弛豫時(shí)間方向移動(dòng);當(dāng)砂巖含有兩相流體時(shí),隨潤濕相流體飽和度的減小,潤濕相流體受限擴(kuò)散越來越明顯,D-T2分布中潤濕相流體信號(hào)逐漸偏離自由擴(kuò)散系數(shù)線,非潤濕相流體信號(hào)的位置不變.

    (3) 巖石膠結(jié)指數(shù)和表面弛豫率都影響D-T2分布中潤濕相流體受限擴(kuò)散系數(shù)線的位置.當(dāng)T2較小時(shí),膠結(jié)指數(shù)對(duì)潤濕相流體受限擴(kuò)散系數(shù)線位置的影響占主導(dǎo)地位,膠結(jié)指數(shù)越大,受限擴(kuò)散越明顯,即潤濕相流體受限擴(kuò)散系數(shù)線偏離自由擴(kuò)散系數(shù)線越遠(yuǎn);當(dāng)巖石膠結(jié)指數(shù)一定時(shí),隨著T2的逐漸增大,表面弛豫率對(duì)潤濕相流體受限擴(kuò)散系數(shù)線位置的影響表現(xiàn)為:表面弛豫率越小,流體受限擴(kuò)散系數(shù)線偏離自由擴(kuò)散系數(shù)線越遠(yuǎn).

    (4) 受限擴(kuò)散使D-T2分布中潤濕相流體信號(hào)偏離其自由擴(kuò)散流體線,需要利用潤濕相流體的受限擴(kuò)散系數(shù)線準(zhǔn)確識(shí)別D-T2分布中的潤濕相流體.

    附錄A流體受限擴(kuò)散系數(shù)線推導(dǎo)

    流體在無限介質(zhì)中擴(kuò)散為隨機(jī)的布朗運(yùn)動(dòng),其均方根位移由愛因斯坦方程(Einstein,1956)決定,為

    (A1)

    式中r(t)是任意t時(shí)刻粒子的位置,r(0)是初始時(shí)刻粒子的位置,D0為粒子的自由擴(kuò)散系數(shù).

    巖石孔隙中的流體由于受到孔隙壁的影響,流體中粒子的運(yùn)動(dòng)表現(xiàn)為受限擴(kuò)散,孔隙中的流體擴(kuò)散系數(shù)隨著擴(kuò)散時(shí)間變化而變化,其在三維空間中的擴(kuò)散系數(shù)表達(dá)式(Mitra et al.,1992)為

    (A2)當(dāng)孔隙中的流體在t時(shí)間內(nèi)擴(kuò)散距離小于孔隙半徑,即擴(kuò)散時(shí)間t較小時(shí),其擴(kuò)散系數(shù)表達(dá)式(Mitra et al.,1993)為

    (A3)

    式中S/V是孔隙的表面積與體積的比值,當(dāng)孔隙中存在多種流體時(shí),S和V分別是各自流體的表面積與體積.

    當(dāng)孔隙中的流體在t時(shí)間內(nèi)擴(kuò)散距離大于孔隙半徑,即擴(kuò)散時(shí)間t較大時(shí),其擴(kuò)散系數(shù)趨于孔隙的曲折度,即

    (A4)

    式中,φ是孔隙度,m是膠結(jié)指數(shù),τ是曲折度.

    根據(jù)Padé近似(Hürlimann et al.,1994)可將以上兩種情況綜合表示為

    (A5)

    巖石孔隙中的流體橫向弛豫過程受體弛豫、表面弛豫和擴(kuò)散弛豫三種機(jī)理的作用,孔隙流體的橫向弛豫時(shí)間T2可表示為

    (A6)

    式中ρeff指有效的表面弛豫率,在只含潤濕相流體時(shí),ρeff=ρ;當(dāng)飽含兩種流體時(shí),

    (A7)

    式中Sg,1是潤濕相流體與巖石顆粒表面的接觸面積,S1,2是潤濕相流體與非潤濕相流體的接觸面積.

    (A8)

    將式(A4)和式(A8)代入式(A5)中,并用T2取代t,可以得到潤濕相流體的擴(kuò)散系數(shù)隨橫向弛豫時(shí)間T2的變化關(guān)系為

    (A9)

    References

    Alghamdi T M, Arns C H, Eyvazzadeh R Y. 2013. Correlations between NMR-relaxation response and relative permeability from tomographic reservoir-rock images. SPE Reservoir Evaluation & Engineering, 16(4): 369-377.

    Cai S H, Chen Q L, Cai C B. 2009. Monte Carlo simulation of liquid self-diffusion with surface relaxation. Chinese Journal of Computational Physics (in Chinese), 26(1): 42-48.

    Cheng J J, Xiao L Z, Xu W, et al. 2013. Simulating NMR responses in porous media based on reconstructed digital rock image. Chinese Journal of Magnetic Resonance (in Chinese), 30(3): 336-344.

    Dunn K J, Latorraca G A. 1999. The inversion of NMR log data sets with different measurement errors. Journal of Magnetic Resonance, 140(1): 153-161.

    Einstein A. 1956. Investigations on the Theory of the Brownian Movement. New York: Dover Publications.

    Hagsl?tt H, J?nsson B, Nydén M, et al. 2003. Predictions of pulsed field gradient NMR echo-decays for molecules diffusing in various restrictive geometries. Simulations of diffusion propagators based on a finite element method. Journal of Magnetic Resonance, 161(2): 138-147.

    Hürlimann M D, Helmer K G, Latour L L, et al. 1994. Restricted diffusion in sedimentary rocks. Determination of surface-area-to-volume ratio and surface relaxivity. Journal of Magnetic Resonance, Series A, 111(2): 169-178.

    Jin G D, Torres-Verdín C, Toumelin E. 2009. Comparison of NMR simulations of porous media derived from analytical and voxelized representations. Journal of Magnetic Resonance, 200(2): 313-320.

    Luo Z X, Paulsen J, Vembusubramanian M, et al. 2015. Restricted diffusion effects on nuclear magnetic resonance DT 2 maps. Geophysics, 80(2): E41-E47.

    Minh C C, Heaton N, Ramamoorthy R, et al. 2003. Planning and interpreting NMR fluid-characterization logs. ∥SPE Annual Technical Conference and Exhibition. Denver, Colorado: Society of Petroleum Engineers. Minh C C, Crary S F, Zielinski L, et al. 2012. 2D-NMR applications in unconventional reservoirs.∥SPE Canadian Unconventional Resources Conference. Calgary, Alberta, Canada: Society of Petroleum Engineers.Mitra P P, Sen P N, Schwartz L M, et al. 1992. Diffusion propagator as a probe of the structure of porous media. Physical Review Letters, 68(24): 3555-3558.

    Mitra P P, Sen P N, Schwartz L M. 1993. Short-time behavior of the diffusion coefficient as a geometrical probe of porous media. Physical Review B, 47(14): 8565-8574.

    Talabi O A. 2008. Pore-scale simulation of NMR response in porous media [Ph. D. thesis]. Department of Earth Science and Engineering, Imperial College London.

    Talabi O, AlSayari S, Iglauer S, et al. 2009. Pore-scale simulation of NMR response. Journal of Petroleum Science & Engineering, 67(3-4): 168-178. Tan M J, Xu J J, Zou Y L, et al. 2014. Nuclear magnetic resonance(NMR) microscopic simulation based on random-walk: theory and parameters analysis. Journal of Central South University, 21(3): 1091-1097.

    Torrey H C. 1956. Bloch equations with diffusion terms. Physical Review, 104(3): 563-565.

    Toumelin E, Torres-Verdin C, Chen S. 2003. Modeling of multiple echo-time NMR measurements for complex pore geometries and multiphase saturations. SPE Reservoir Evaluation & Engineering, 6(4): 234-243. Toumelin E, Torres-Verdin C, Sun B Q, et al. 2004. A numerical assessment of modern borehole NMR interpretation techniques. ∥SPE Annual Technical Conference and Exhibition. Houston, Texas: Society of Petroleum Engineers.

    Toumelin E. 2006. Pore-scale petrophysical models for the simulation and combined interpretation of nuclear magnetic resonance and wide-band electromagnetic measurements of saturated rocks [Ph. D. thesis]. The University of Texas at Austin.

    Toumelin E, Torres-Verdín C, Sun B Q, et al. 2007. Random-walk technique for simulating NMR measurements and 2D NMR maps of porous media with relaxing and permeable boundaries. Journal of Magnetic Resonance, 188(1): 83-96.

    Valfouskaya A, Adler P M, Thovert J F, et al. 2006. Nuclear magnetic resonance diffusion with surface relaxation in porous media. Journal of Colloid and Interface Science, 295(1): 188-201.

    Xiao L Z, Xie R H, Liao G Z. 2012. The Principle and Methods of Nuclear Magnetic Resonance Well Logging of Complex Hydrocarbon Reservoir in China (in Chinese). Beijing: Science Press. Xie R H, Xiao L Z, Liu J J, et al. 2009. A method for multiple echo trains jointing inversion of NMR relaxation measurements. Chinese J. Geophys. (in Chinese), 52(11): 2913-2919, doi: 10.3969/j.issn.0001-5733.2009.11.027.

    Zhang Z F. 2013. Research on 3D Laplace NMR method and simulation of diffusion for water molecular in micropore [Master′s thesis] (in Chinese). Beijing: China University of Petroleum (Beijing).

    Zheng L H, Chiew Y C. 1989. Computer simulation of diffusion-controlled reactions in dispersions of spherical sinks. Journal of Chemical Physics, 90(1): 322-327.

    Zielinski L, Ramamoorthy R, Minh C C, et al. 2010. Restricted diffusion effects in saturation estimates from 2D diffusion-relaxation NMR maps. ∥SPE Annual Technical Conference and Exhibition. Florence, Italy: Society of Petroleum Engineers. Zientara G P, Freed J H. 1980. Spin-echoes for diffusion in bounded, heterogeneous media: A numerical study. Journal of Chemical Physics, 72(2): 1285-1292.

    Zou Y L. 2013. Sedimentary rock reconstruction and pore-scale NMR response simulation [Master′s thesis] (in Chinese). Beijing: China University of Geosciences (Beijing).

    附中文參考文獻(xiàn)

    蔡淑惠, 陳巧龍, 蔡聰波. 2009. 具有表面弛豫的液體自擴(kuò)散的Monte Carlo模擬. 計(jì)算物理, 26(1): 42-48.

    成家杰, 肖立志, 許巍等. 2013. 基于巖石重構(gòu)圖像的核磁共振響應(yīng)模擬. 波譜學(xué)雜志, 30(3): 336-344.

    肖立志, 謝然紅, 廖廣志. 2012. 中國復(fù)雜油氣藏核磁共振測井理論與方法. 北京: 科學(xué)出版社.

    謝然紅, 肖立志, 劉家軍等. 2009. 核磁共振多回波串聯(lián)合反演方法. 地球物理學(xué)報(bào), 52(11): 2913-2919, doi: 10.3969/j.issn.0001-5733.2009.11.027.

    張宗富. 2013. 三維核磁共振方法研究與水分子擴(kuò)散模擬[碩士論文]. 北京: 中國石油大學(xué)(北京).

    鄒友龍. 2013. 沉積巖重構(gòu)及其孔隙尺度NMR響應(yīng)模擬[碩士論文]. 北京: 中國地質(zhì)大學(xué)(北京).

    (本文編輯何燕)

    基金項(xiàng)目國家自然科學(xué)基金(41272163)及高等學(xué)校博士學(xué)科點(diǎn)專項(xiàng)科研基金(20130007110012)資助.

    作者簡介郭江峰,男,1991年生,博士研究生,主要從事巖石物理及核磁共振測井方法研究.E-mail:jiangfeng_guo@163.com E-mail:xieranhong@cup.edu.com

    *通訊作者謝然紅,女,1966年生,教授,主要從事巖石物理、核磁共振測井方法及測井儲(chǔ)層評(píng)價(jià)等方面的研究及教學(xué)工作.

    doi:10.6038/cjg20160733 中圖分類號(hào)P631

    收稿日期2015-06-24,2016-06-07收修定稿

    Simulation of NMR responses in sandstone and restricted diffusion

    GUO Jiang-Feng1,2, XIE Ran-Hong1,2*, ZOU You-Long1,2

    1StateKeyLaboratoryofPetroleumResourcesandProspecting,ChinaUniversityofPetroleum,Beijing102249,China2KeyLaboratoryofEarthProspectingandInformationTechnology,ChinaUniversityofPetroleum,Beijing102249,China

    AbstractNuclear magnetic resonance (NMR) responses and restricted diffusion phenomenon of the fluids in sandstone were simulated by the Random-Walk method. Sandstone of different pore sizes was generated by changing digital cores′ resolution, permitting to study how the fluid diffusion coefficient changes with time when the sandstone of different pore sizes were saturated with water, and to simulate the NMR responses of the either single or two phase fluids in the sandstone. Then the relationship between the fluid restricted diffusion coefficient and the T2 relaxation time was researched, and the influence on restricted diffusion phenomenon of wetting phase fluid posed by the surface relaxivity and cementation index of the sandstone was analyzed. Finally the analysis results were used to interpret the D-T2 distributions of the sandstone. It was concluded that the diffusion coefficient of fluid gradually decreases and tends to be constant with increase of diffusion time. The smaller pore size of the rock, the restricted diffusion phenomenon is more obvious, and the restricted diffusion has a greater influence on NMR responses. The position of the restricted diffusion coefficient line of wetting phase fluid is affected by the rock cementation index and surface relaxivity. The wetting phase fluid signal in the D-T2 distributions deviates from the unrestricted diffusion coefficient line due to the decreasing diffusion coefficient of the wetting phase fluid, hence the restricted diffusion coefficient line of fluid should be used to recognize the wetting phase fluid in D-T2 distributions.

    KeywordsSandstone; NMR responses; Random-Walk method; Restricted diffusion; Diffusion coefficient line

    郭江峰, 謝然紅, 鄒友龍. 2016. 砂巖核磁共振響應(yīng)模擬及受限擴(kuò)散. 地球物理學(xué)報(bào),59(7):2703-2712,doi:10.6038/cjg20160733.

    Guo J F, Xie R H, Zou Y L. 2016. Simulation of NMR responses in sandstone and restricted diffusion. Chinese J. Geophys. (in Chinese),59(7):2703-2712,doi:10.6038/cjg20160733.

    猜你喜歡
    砂巖
    三種不同粒徑砂巖的強(qiáng)度與破壞特征
    CSAMT法在柴北緣砂巖型鈾礦勘查砂體探測中的應(yīng)用
    火星上的漩渦層狀砂巖
    砂巖:黏結(jié)在一起的沙子
    基于DEM的地貌特征分析與類型劃分——以砒砂巖區(qū)為例
    低滲砂巖地層因素的應(yīng)力敏感研究
    裂縫性致密砂巖氣藏入井流體傷害規(guī)律
    賀蘭口砂巖吸水率的研究
    多裂縫系統(tǒng)致密砂巖氣藏的水鎖及應(yīng)力敏感疊加傷害研究:以川西蓬萊鎮(zhèn)組致密砂巖為例
    蘇里格氣田致密砂巖氣層識(shí)別難點(diǎn)及方法評(píng)述
    三级毛片av免费| 香蕉丝袜av| 啪啪无遮挡十八禁网站| 18美女黄网站色大片免费观看| 国产精品秋霞免费鲁丝片| 正在播放国产对白刺激| 免费在线观看视频国产中文字幕亚洲| 人人妻,人人澡人人爽秒播| 日本精品一区二区三区蜜桃| 欧美精品一区二区免费开放| 国产人伦9x9x在线观看| 老熟妇乱子伦视频在线观看| 波多野结衣av一区二区av| 久久久国产精品麻豆| 真人一进一出gif抽搐免费| 天堂影院成人在线观看| 亚洲av五月六月丁香网| 美女大奶头视频| 亚洲av电影在线进入| 亚洲一卡2卡3卡4卡5卡精品中文| 高潮久久久久久久久久久不卡| 黑人巨大精品欧美一区二区mp4| 99re在线观看精品视频| 欧美激情高清一区二区三区| 国产成人欧美在线观看| 国产成人精品久久二区二区免费| 99久久人妻综合| 国产一区二区在线av高清观看| 色精品久久人妻99蜜桃| 校园春色视频在线观看| 咕卡用的链子| 欧美在线黄色| 亚洲五月天丁香| 97超级碰碰碰精品色视频在线观看| 国产麻豆69| 波多野结衣一区麻豆| 视频区图区小说| 老司机福利观看| 久久伊人香网站| 一本综合久久免费| 午夜精品国产一区二区电影| 久久精品人人爽人人爽视色| 亚洲欧美激情综合另类| 精品高清国产在线一区| 女人被躁到高潮嗷嗷叫费观| 五月开心婷婷网| 麻豆国产av国片精品| 我的亚洲天堂| 久久久久久人人人人人| 99国产精品免费福利视频| 成人手机av| 最近最新中文字幕大全电影3 | www国产在线视频色| 欧美日韩国产mv在线观看视频| 自线自在国产av| 国产一区二区激情短视频| 久久亚洲真实| 国产精品一区二区免费欧美| 丝袜在线中文字幕| 十八禁网站免费在线| 久久久国产成人免费| 亚洲 欧美一区二区三区| 999久久久精品免费观看国产| 欧美激情 高清一区二区三区| 嫁个100分男人电影在线观看| 国产aⅴ精品一区二区三区波| 成人特级黄色片久久久久久久| 欧美激情 高清一区二区三区| 国产精华一区二区三区| 人人澡人人妻人| 热re99久久精品国产66热6| 美女高潮喷水抽搐中文字幕| 别揉我奶头~嗯~啊~动态视频| 日日爽夜夜爽网站| 999久久久国产精品视频| 久久香蕉国产精品| 午夜老司机福利片| 啦啦啦 在线观看视频| 69av精品久久久久久| 精品日产1卡2卡| 日韩精品中文字幕看吧| 一a级毛片在线观看| 男女之事视频高清在线观看| 国产激情久久老熟女| 水蜜桃什么品种好| 午夜福利欧美成人| 免费女性裸体啪啪无遮挡网站| 亚洲一区高清亚洲精品| 欧美 亚洲 国产 日韩一| 天天影视国产精品| 69av精品久久久久久| 欧美中文日本在线观看视频| 黑人欧美特级aaaaaa片| 亚洲一区高清亚洲精品| 男女下面插进去视频免费观看| 午夜精品在线福利| 在线播放国产精品三级| 一本综合久久免费| 一区二区三区激情视频| 丝袜在线中文字幕| 成人三级黄色视频| 精品人妻在线不人妻| 长腿黑丝高跟| 欧美老熟妇乱子伦牲交| 国产亚洲欧美在线一区二区| 18禁美女被吸乳视频| 99在线人妻在线中文字幕| 岛国视频午夜一区免费看| 岛国视频午夜一区免费看| 亚洲熟妇熟女久久| 在线天堂中文资源库| 满18在线观看网站| 精品熟女少妇八av免费久了| 成熟少妇高潮喷水视频| 国产精品1区2区在线观看.| 亚洲成人国产一区在线观看| 日本wwww免费看| 欧美日韩国产mv在线观看视频| 亚洲国产毛片av蜜桃av| 久久精品国产亚洲av高清一级| 欧美日韩亚洲国产一区二区在线观看| 成人特级黄色片久久久久久久| 亚洲精品美女久久av网站| 久久这里只有精品19| 99国产精品一区二区三区| 两个人免费观看高清视频| 久久性视频一级片| 久久 成人 亚洲| 91成人精品电影| 高清黄色对白视频在线免费看| 国产片内射在线| 国产无遮挡羞羞视频在线观看| 99在线人妻在线中文字幕| av片东京热男人的天堂| 国产精品久久久久成人av| 国产不卡一卡二| 午夜福利在线观看吧| 久久久久久久久免费视频了| 亚洲精品一区av在线观看| 亚洲av熟女| 最新美女视频免费是黄的| 身体一侧抽搐| 欧美激情高清一区二区三区| 欧美久久黑人一区二区| 男人舔女人的私密视频| 免费女性裸体啪啪无遮挡网站| 国产精品久久电影中文字幕| 精品一区二区三区视频在线观看免费 | 久久性视频一级片| 久久 成人 亚洲| 美女大奶头视频| 久久久久国内视频| 可以在线观看毛片的网站| 日韩免费av在线播放| 久久中文字幕一级| 成人三级做爰电影| 国产精品亚洲一级av第二区| 亚洲狠狠婷婷综合久久图片| 在线播放国产精品三级| 欧美日韩亚洲高清精品| 欧美激情高清一区二区三区| 级片在线观看| 神马国产精品三级电影在线观看 | 久久精品91蜜桃| av有码第一页| 欧美另类亚洲清纯唯美| 国产无遮挡羞羞视频在线观看| 亚洲人成77777在线视频| 欧美一级毛片孕妇| 亚洲三区欧美一区| 一个人免费在线观看的高清视频| 成人18禁在线播放| 精品一区二区三区视频在线观看免费 | 国产色视频综合| 丝袜美足系列| 久久天堂一区二区三区四区| 久久久久久亚洲精品国产蜜桃av| 久久久久久免费高清国产稀缺| 露出奶头的视频| 美女扒开内裤让男人捅视频| 亚洲 国产 在线| 日本a在线网址| 国产av一区在线观看免费| 成人av一区二区三区在线看| 高清黄色对白视频在线免费看| 国产成人av激情在线播放| 村上凉子中文字幕在线| 麻豆av在线久日| 国产激情欧美一区二区| 国产真人三级小视频在线观看| 午夜福利在线观看吧| 亚洲精品久久午夜乱码| 久久国产乱子伦精品免费另类| av超薄肉色丝袜交足视频| 50天的宝宝边吃奶边哭怎么回事| 亚洲少妇的诱惑av| 一个人观看的视频www高清免费观看 | 婷婷丁香在线五月| 亚洲欧美激情在线| 美女高潮到喷水免费观看| 性欧美人与动物交配| 亚洲专区国产一区二区| 久久久久久大精品| √禁漫天堂资源中文www| 女人被躁到高潮嗷嗷叫费观| 在线看a的网站| 久久午夜综合久久蜜桃| 欧美另类亚洲清纯唯美| 又黄又粗又硬又大视频| 精品国产美女av久久久久小说| 中文字幕人妻熟女乱码| 国产精品久久久久成人av| 天堂动漫精品| 级片在线观看| 国产一区在线观看成人免费| 精品国产一区二区久久| 国产精品98久久久久久宅男小说| 国产av精品麻豆| 国产亚洲欧美精品永久| 亚洲熟女毛片儿| 久久中文字幕人妻熟女| 久久午夜综合久久蜜桃| 久久人人精品亚洲av| 十八禁人妻一区二区| 国产伦一二天堂av在线观看| 美女高潮到喷水免费观看| www.精华液| 97碰自拍视频| 制服诱惑二区| 欧美精品亚洲一区二区| 成人免费观看视频高清| 国内久久婷婷六月综合欲色啪| 嫩草影视91久久| 视频在线观看一区二区三区| 欧美性长视频在线观看| 亚洲成av片中文字幕在线观看| 国产97色在线日韩免费| 丰满饥渴人妻一区二区三| 久久人妻福利社区极品人妻图片| 亚洲精品av麻豆狂野| avwww免费| 国产精品秋霞免费鲁丝片| 亚洲专区国产一区二区| 亚洲九九香蕉| 99国产精品99久久久久| 国产一卡二卡三卡精品| 狂野欧美激情性xxxx| 一进一出抽搐gif免费好疼 | 黑人巨大精品欧美一区二区蜜桃| 国产精华一区二区三区| www日本在线高清视频| 99在线人妻在线中文字幕| 日本三级黄在线观看| 男女午夜视频在线观看| 成人黄色视频免费在线看| 亚洲午夜精品一区,二区,三区| 国产免费现黄频在线看| 制服人妻中文乱码| 国产一区在线观看成人免费| 女人被狂操c到高潮| 少妇粗大呻吟视频| 国产精品 国内视频| 国产亚洲精品一区二区www| 88av欧美| 久久久久精品国产欧美久久久| 亚洲精品一区av在线观看| 国产色视频综合| 欧美日韩av久久| 叶爱在线成人免费视频播放| 中文字幕高清在线视频| 久久中文看片网| 久久久国产一区二区| 精品人妻1区二区| 麻豆国产av国片精品| 亚洲专区国产一区二区| 性色av乱码一区二区三区2| 一级毛片女人18水好多| 另类亚洲欧美激情| 性欧美人与动物交配| 18禁观看日本| 在线观看免费午夜福利视频| 欧美日韩黄片免| а√天堂www在线а√下载| 看片在线看免费视频| 99热只有精品国产| 欧洲精品卡2卡3卡4卡5卡区| 久久狼人影院| 精品日产1卡2卡| 国产男靠女视频免费网站| 麻豆国产av国片精品| 黑人巨大精品欧美一区二区mp4| 亚洲国产精品一区二区三区在线| 免费在线观看视频国产中文字幕亚洲| av在线天堂中文字幕 | 宅男免费午夜| 黑人猛操日本美女一级片| 亚洲国产精品999在线| 国产成人一区二区三区免费视频网站| 亚洲成人久久性| 桃红色精品国产亚洲av| av免费在线观看网站| 欧美日韩亚洲综合一区二区三区_| 一级毛片精品| www.999成人在线观看| 亚洲精品国产色婷婷电影| 国产一区在线观看成人免费| 老司机福利观看| 色在线成人网| 首页视频小说图片口味搜索| 国产三级在线视频| 18禁黄网站禁片午夜丰满| 999精品在线视频| 丝袜人妻中文字幕| 美女国产高潮福利片在线看| 免费人成视频x8x8入口观看| 日韩精品中文字幕看吧| 国产伦人伦偷精品视频| 在线国产一区二区在线| 国产亚洲欧美在线一区二区| 91精品国产国语对白视频| 91精品三级在线观看| 成人18禁高潮啪啪吃奶动态图| 日韩欧美一区视频在线观看| 成年女人毛片免费观看观看9| 久久久久亚洲av毛片大全| 一区在线观看完整版| 亚洲欧美激情在线| 18美女黄网站色大片免费观看| 激情在线观看视频在线高清| 亚洲aⅴ乱码一区二区在线播放 | 69精品国产乱码久久久| 亚洲中文字幕日韩| 视频在线观看一区二区三区| 国产亚洲av高清不卡| 最近最新免费中文字幕在线| 激情在线观看视频在线高清| 一级片免费观看大全| 在线十欧美十亚洲十日本专区| 国产有黄有色有爽视频| 9191精品国产免费久久| 午夜福利在线免费观看网站| 在线观看舔阴道视频| 午夜福利在线免费观看网站| 精品久久久精品久久久| 一个人免费在线观看的高清视频| 在线观看舔阴道视频| 亚洲狠狠婷婷综合久久图片| 精品卡一卡二卡四卡免费| 成人免费观看视频高清| 麻豆一二三区av精品| 长腿黑丝高跟| 黄色毛片三级朝国网站| 美女高潮到喷水免费观看| 国产精品电影一区二区三区| 亚洲精品粉嫩美女一区| 一级毛片精品| 国产一区二区三区在线臀色熟女 | 夜夜爽天天搞| av免费在线观看网站| 精品福利永久在线观看| 老司机在亚洲福利影院| 久久久久久久久久久久大奶| 九色亚洲精品在线播放| 久久久久久久精品吃奶| 老司机在亚洲福利影院| 亚洲国产欧美一区二区综合| 侵犯人妻中文字幕一二三四区| 久久亚洲真实| 亚洲欧美激情在线| 国产精华一区二区三区| 国产无遮挡羞羞视频在线观看| 两性夫妻黄色片| 久久久久久免费高清国产稀缺| 欧美精品一区二区免费开放| 12—13女人毛片做爰片一| 欧美色视频一区免费| 亚洲一区二区三区色噜噜 | 精品福利观看| 亚洲人成电影观看| 亚洲精品在线美女| 中文欧美无线码| 动漫黄色视频在线观看| 真人一进一出gif抽搐免费| 免费在线观看影片大全网站| 国产亚洲精品久久久久5区| 欧美乱色亚洲激情| 国产精品自产拍在线观看55亚洲| 超碰成人久久| 色精品久久人妻99蜜桃| 中文亚洲av片在线观看爽| 国产精品成人在线| 欧美亚洲日本最大视频资源| 热99re8久久精品国产| 日韩欧美国产一区二区入口| 亚洲成国产人片在线观看| 一级a爱片免费观看的视频| 成人三级做爰电影| 久久亚洲精品不卡| 久久午夜综合久久蜜桃| 国产精品一区二区三区四区久久 | 真人一进一出gif抽搐免费| 亚洲国产欧美日韩在线播放| 精品久久久久久电影网| 亚洲国产毛片av蜜桃av| 亚洲男人的天堂狠狠| 日韩欧美一区视频在线观看| 国产在线精品亚洲第一网站| 琪琪午夜伦伦电影理论片6080| 丝袜美足系列| 亚洲七黄色美女视频| 老司机在亚洲福利影院| 国产成人影院久久av| 不卡一级毛片| 极品人妻少妇av视频| 99久久人妻综合| 女人被狂操c到高潮| 久久狼人影院| 精品第一国产精品| 欧美激情极品国产一区二区三区| 精品久久久精品久久久| 国产成人精品无人区| 欧美精品亚洲一区二区| 国产精华一区二区三区| 午夜精品国产一区二区电影| 在线av久久热| 精品人妻在线不人妻| 丰满迷人的少妇在线观看| 长腿黑丝高跟| 好看av亚洲va欧美ⅴa在| 最近最新中文字幕大全免费视频| 一边摸一边做爽爽视频免费| 久99久视频精品免费| 日韩精品青青久久久久久| 操美女的视频在线观看| 精品国产一区二区三区四区第35| 久久影院123| 久久久水蜜桃国产精品网| 亚洲国产欧美网| 亚洲精品中文字幕一二三四区| 一二三四在线观看免费中文在| 男人的好看免费观看在线视频 | 精品午夜福利视频在线观看一区| 亚洲激情在线av| 脱女人内裤的视频| 俄罗斯特黄特色一大片| 亚洲国产欧美网| 日本黄色日本黄色录像| 一级毛片精品| 不卡一级毛片| 一级作爱视频免费观看| 午夜免费成人在线视频| 法律面前人人平等表现在哪些方面| 国产蜜桃级精品一区二区三区| 国产成人av激情在线播放| 久久久国产精品麻豆| 国产亚洲精品一区二区www| 国产精品久久久久久人妻精品电影| 亚洲欧美激情在线| 午夜老司机福利片| 久久狼人影院| 黄色怎么调成土黄色| 久久久精品欧美日韩精品| 亚洲一区中文字幕在线| 国产色视频综合| 黑人巨大精品欧美一区二区蜜桃| av网站在线播放免费| 他把我摸到了高潮在线观看| 久久国产亚洲av麻豆专区| 亚洲九九香蕉| 精品少妇一区二区三区视频日本电影| 国产精品美女特级片免费视频播放器 | 天天躁狠狠躁夜夜躁狠狠躁| 亚洲欧洲精品一区二区精品久久久| 国产欧美日韩精品亚洲av| 精品国产国语对白av| 久久99一区二区三区| 老熟妇乱子伦视频在线观看| 国产欧美日韩一区二区精品| 日韩成人在线观看一区二区三区| www.熟女人妻精品国产| av有码第一页| 曰老女人黄片| 亚洲色图 男人天堂 中文字幕| 国产av在哪里看| 国产亚洲欧美精品永久| 亚洲av成人av| 欧美性长视频在线观看| 欧美乱妇无乱码| 神马国产精品三级电影在线观看 | 亚洲精品美女久久av网站| 女人爽到高潮嗷嗷叫在线视频| 免费人成视频x8x8入口观看| 成人av一区二区三区在线看| www国产在线视频色| 1024香蕉在线观看| 色综合婷婷激情| 日韩av在线大香蕉| 午夜视频精品福利| 成熟少妇高潮喷水视频| 一二三四社区在线视频社区8| 国产伦一二天堂av在线观看| 99riav亚洲国产免费| 亚洲av成人不卡在线观看播放网| 夫妻午夜视频| videosex国产| 国产精品1区2区在线观看.| 精品熟女少妇八av免费久了| 国产不卡一卡二| 免费观看精品视频网站| 亚洲七黄色美女视频| 妹子高潮喷水视频| 国产单亲对白刺激| 亚洲久久久国产精品| 无限看片的www在线观看| 18禁裸乳无遮挡免费网站照片 | 91在线观看av| 国产成人影院久久av| 99国产极品粉嫩在线观看| 色尼玛亚洲综合影院| videosex国产| 变态另类成人亚洲欧美熟女 | 咕卡用的链子| 高清毛片免费观看视频网站 | av福利片在线| 日韩精品免费视频一区二区三区| xxxhd国产人妻xxx| 高清黄色对白视频在线免费看| 身体一侧抽搐| 自拍欧美九色日韩亚洲蝌蚪91| 国产精品一区二区三区四区久久 | 熟女少妇亚洲综合色aaa.| 天堂√8在线中文| 中文字幕精品免费在线观看视频| 欧美乱妇无乱码| 成人三级黄色视频| 免费观看精品视频网站| 亚洲av熟女| 一进一出抽搐gif免费好疼 | 国产一区二区在线av高清观看| 日韩精品中文字幕看吧| 在线观看日韩欧美| av欧美777| 精品高清国产在线一区| 国产黄色免费在线视频| 亚洲精品在线观看二区| 国产日韩一区二区三区精品不卡| 99精国产麻豆久久婷婷| 在线观看www视频免费| 亚洲成国产人片在线观看| 天堂动漫精品| 一边摸一边抽搐一进一出视频| 久9热在线精品视频| av在线播放免费不卡| 亚洲 欧美一区二区三区| 人成视频在线观看免费观看| 成人国产一区最新在线观看| 级片在线观看| 少妇裸体淫交视频免费看高清 | 丝袜人妻中文字幕| 亚洲美女黄片视频| 怎么达到女性高潮| 69av精品久久久久久| 日本免费一区二区三区高清不卡 | 在线观看66精品国产| 97超级碰碰碰精品色视频在线观看| 黄片小视频在线播放| 免费在线观看视频国产中文字幕亚洲| 欧美人与性动交α欧美精品济南到| 美女午夜性视频免费| 麻豆一二三区av精品| 亚洲一区二区三区欧美精品| 色在线成人网| 亚洲 欧美一区二区三区| 老鸭窝网址在线观看| 午夜免费成人在线视频| 久久人人爽av亚洲精品天堂| 男女下面插进去视频免费观看| 91字幕亚洲| 成人影院久久| 欧美日韩乱码在线| 丝袜在线中文字幕| 日韩欧美在线二视频| 国产av精品麻豆| 亚洲一区二区三区色噜噜 | 无遮挡黄片免费观看| 免费观看精品视频网站| 18禁美女被吸乳视频| 18禁观看日本| 国产aⅴ精品一区二区三区波| 天堂中文最新版在线下载| 国产精品综合久久久久久久免费 | 久久国产亚洲av麻豆专区| 1024视频免费在线观看| 日韩成人在线观看一区二区三区| 亚洲色图综合在线观看| 久久99一区二区三区| 国产深夜福利视频在线观看| 亚洲情色 制服丝袜| cao死你这个sao货| av中文乱码字幕在线| 日韩大尺度精品在线看网址 | 搡老乐熟女国产| 久久久精品欧美日韩精品| 身体一侧抽搐| www.熟女人妻精品国产| 黄频高清免费视频| 好男人电影高清在线观看| 两个人看的免费小视频| 欧美成人性av电影在线观看| 999久久久国产精品视频| 国产aⅴ精品一区二区三区波| 最近最新中文字幕大全免费视频| 久久欧美精品欧美久久欧美| 9色porny在线观看| 一个人观看的视频www高清免费观看 |