• <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)述
    videos熟女内射| 你懂的网址亚洲精品在线观看 | 免费观看的影片在线观看| 美女内射精品一级片tv| 精品久久久久久久末码| 亚洲伊人久久精品综合 | 欧美一区二区国产精品久久精品| 国产三级中文精品| 欧美激情国产日韩精品一区| 亚洲无线观看免费| 国产精品人妻久久久久久| 少妇熟女aⅴ在线视频| 中国美白少妇内射xxxbb| 亚洲天堂国产精品一区在线| 欧美变态另类bdsm刘玥| 欧美丝袜亚洲另类| 91精品国产九色| 日日干狠狠操夜夜爽| 神马国产精品三级电影在线观看| 夜夜看夜夜爽夜夜摸| 欧美精品一区二区大全| 亚洲欧美成人精品一区二区| 麻豆一二三区av精品| 中文字幕av在线有码专区| 日韩欧美在线乱码| 久久亚洲国产成人精品v| 国产午夜精品论理片| 亚洲国产日韩欧美精品在线观看| 国产精品乱码一区二三区的特点| www.av在线官网国产| 大香蕉久久网| 成人亚洲欧美一区二区av| 男女啪啪激烈高潮av片| 亚洲国产精品sss在线观看| 久久精品久久久久久噜噜老黄 | 美女xxoo啪啪120秒动态图| 国产黄色小视频在线观看| 欧美激情在线99| 亚洲综合精品二区| 精品人妻熟女av久视频| 中文乱码字字幕精品一区二区三区 | 亚洲精品自拍成人| 草草在线视频免费看| 99九九线精品视频在线观看视频| 听说在线观看完整版免费高清| 噜噜噜噜噜久久久久久91| 听说在线观看完整版免费高清| 少妇人妻精品综合一区二区| 欧美变态另类bdsm刘玥| 久久99精品国语久久久| 久久综合国产亚洲精品| 女人被狂操c到高潮| 中文资源天堂在线| 亚洲成人中文字幕在线播放| 卡戴珊不雅视频在线播放| 亚洲av不卡在线观看| 伊人久久精品亚洲午夜| 搡女人真爽免费视频火全软件| 久久6这里有精品| 嘟嘟电影网在线观看| 国产精品99久久久久久久久| 在线观看一区二区三区| 一边摸一边抽搐一进一小说| 成人一区二区视频在线观看| 国产极品天堂在线| 欧美性感艳星| 欧美另类亚洲清纯唯美| 国产精品久久久久久久电影| 国产成年人精品一区二区| 国产亚洲av嫩草精品影院| 老师上课跳d突然被开到最大视频| 在线观看66精品国产| 菩萨蛮人人尽说江南好唐韦庄 | 色5月婷婷丁香| 蜜臀久久99精品久久宅男| 久久久久久大精品| 免费看日本二区| 内射极品少妇av片p| 精品人妻视频免费看| 久久99热这里只频精品6学生 | 午夜日本视频在线| 91在线精品国自产拍蜜月| 欧美成人免费av一区二区三区| 少妇熟女aⅴ在线视频| 久久精品久久久久久噜噜老黄 | 亚洲精品aⅴ在线观看| 日韩三级伦理在线观看| 久久精品夜色国产| 日韩国内少妇激情av| 哪个播放器可以免费观看大片| 99国产精品一区二区蜜桃av| 观看免费一级毛片| 又粗又爽又猛毛片免费看| 综合色av麻豆| 老司机福利观看| 白带黄色成豆腐渣| 亚洲自偷自拍三级| 免费在线观看成人毛片| 伦理电影大哥的女人| 精品久久久久久电影网 | 91午夜精品亚洲一区二区三区| 男女国产视频网站| 老司机福利观看| 欧美最新免费一区二区三区| 国产亚洲91精品色在线| 久久精品久久久久久久性| 国产精品乱码一区二三区的特点| 免费不卡的大黄色大毛片视频在线观看 | 精品一区二区三区视频在线| 久久久久久久久久成人| 99热网站在线观看| 亚洲成人久久爱视频| 尾随美女入室| eeuss影院久久| 性插视频无遮挡在线免费观看| 精品人妻一区二区三区麻豆| 国内少妇人妻偷人精品xxx网站| 91久久精品国产一区二区三区| 成人鲁丝片一二三区免费| 熟妇人妻久久中文字幕3abv| 少妇人妻一区二区三区视频| 天天躁日日操中文字幕| 观看免费一级毛片| 亚洲精品aⅴ在线观看| 少妇的逼好多水| 国产片特级美女逼逼视频| 国产免费又黄又爽又色| 边亲边吃奶的免费视频| 淫秽高清视频在线观看| 99热全是精品| 白带黄色成豆腐渣| 我的女老师完整版在线观看| 精品人妻熟女av久视频| 欧美一级a爱片免费观看看| 中文亚洲av片在线观看爽| 久久久久久久久大av| av女优亚洲男人天堂| 一级毛片aaaaaa免费看小| 我的老师免费观看完整版| 国产亚洲91精品色在线| 啦啦啦啦在线视频资源| 99久久九九国产精品国产免费| 日本午夜av视频| 青青草视频在线视频观看| 亚洲精品456在线播放app| 99久久九九国产精品国产免费| 最近最新中文字幕免费大全7| 国产真实乱freesex| 国产成人一区二区在线| 亚洲伊人久久精品综合 | 99视频精品全部免费 在线| av免费在线看不卡| 色视频www国产| 一个人看的www免费观看视频| 精品午夜福利在线看| 麻豆精品久久久久久蜜桃| 国产精品国产高清国产av| 国产精品一及| 韩国高清视频一区二区三区| 国产av码专区亚洲av| 成年女人永久免费观看视频| 国产久久久一区二区三区| 最近中文字幕2019免费版| 日本wwww免费看| 亚洲在线观看片| 可以在线观看毛片的网站| 亚洲成人av在线免费| 女人久久www免费人成看片 | 久久99热这里只有精品18| or卡值多少钱| av女优亚洲男人天堂| 国产成人精品一,二区| 亚洲精华国产精华液的使用体验| 亚洲在久久综合| 国产午夜福利久久久久久| 亚洲欧美日韩无卡精品| 一边亲一边摸免费视频| 最后的刺客免费高清国语| 一级av片app| 赤兔流量卡办理| 亚洲欧美精品自产自拍| 91av网一区二区| 欧美一区二区国产精品久久精品| 18禁在线无遮挡免费观看视频| 欧美丝袜亚洲另类| 99久久精品一区二区三区| 人体艺术视频欧美日本| 色综合亚洲欧美另类图片| 少妇猛男粗大的猛烈进出视频 | 亚洲人成网站高清观看| 亚洲国产高清在线一区二区三| 国产国拍精品亚洲av在线观看| 精品一区二区三区视频在线| 亚洲av中文字字幕乱码综合| 国产乱人视频| 精品久久久久久久久av| 日本一二三区视频观看| 亚洲欧美精品综合久久99| 国产在线男女| 亚洲电影在线观看av| 搡女人真爽免费视频火全软件| 亚洲欧美日韩卡通动漫| 啦啦啦啦在线视频资源| 国产成人午夜福利电影在线观看| 中文字幕av成人在线电影| 高清日韩中文字幕在线| 日本免费a在线| 桃色一区二区三区在线观看| 在线播放无遮挡| 日韩人妻高清精品专区| 亚洲综合精品二区| 大香蕉久久网| 少妇的逼好多水| 午夜福利成人在线免费观看| 波多野结衣高清无吗| 毛片一级片免费看久久久久| 欧美日韩国产亚洲二区| 九色成人免费人妻av| 欧美极品一区二区三区四区| 国产极品天堂在线| 亚洲精品国产av成人精品| 欧美日本视频| 少妇的逼好多水| 亚洲久久久久久中文字幕| 中文字幕av成人在线电影| 18禁在线无遮挡免费观看视频| 国产又色又爽无遮挡免| 免费不卡的大黄色大毛片视频在线观看 | 免费看av在线观看网站| 久久精品人妻少妇| 欧美色视频一区免费| 国内精品一区二区在线观看| 久久久精品欧美日韩精品| 国产 一区 欧美 日韩| 久久99蜜桃精品久久| 在现免费观看毛片| av福利片在线观看| 女人十人毛片免费观看3o分钟| 啦啦啦韩国在线观看视频| 亚洲欧洲国产日韩| 超碰97精品在线观看| 久久精品91蜜桃| av在线播放精品| 精品久久久久久久人妻蜜臀av| www.色视频.com| 国产精品一二三区在线看| 色网站视频免费| 一级毛片久久久久久久久女| 久久久久久久久中文| 国产精品一及| 18禁动态无遮挡网站| ponron亚洲| 一区二区三区免费毛片| 欧美不卡视频在线免费观看| 天美传媒精品一区二区| 亚洲成av人片在线播放无| 国产欧美另类精品又又久久亚洲欧美| 亚洲精品乱码久久久v下载方式| 日本五十路高清| 精品欧美国产一区二区三| 少妇高潮的动态图| 嫩草影院新地址| 赤兔流量卡办理| 亚洲精品乱久久久久久| 亚洲精品乱码久久久v下载方式| 免费看a级黄色片| 男插女下体视频免费在线播放| 欧美激情国产日韩精品一区| 一级二级三级毛片免费看| 天堂网av新在线| 少妇丰满av| 在线免费观看不下载黄p国产| 免费观看在线日韩| 91久久精品国产一区二区成人| 在线播放无遮挡| 嫩草影院新地址| 久久精品熟女亚洲av麻豆精品 | 久久99热这里只有精品18| 干丝袜人妻中文字幕| 亚洲电影在线观看av| 一二三四中文在线观看免费高清| 亚洲人成网站在线观看播放| 精品久久久久久久久av| 中文字幕久久专区| 国产高清视频在线观看网站| 国产精品麻豆人妻色哟哟久久 | 18禁在线无遮挡免费观看视频| 中文字幕av在线有码专区| 久久久久久久久久黄片| 久久精品国产鲁丝片午夜精品| av播播在线观看一区| 乱系列少妇在线播放| 小说图片视频综合网站| 成人午夜高清在线视频| 搞女人的毛片| 观看美女的网站| 天堂av国产一区二区熟女人妻| 国产私拍福利视频在线观看| 成年av动漫网址| 亚洲自拍偷在线| 久久精品夜色国产| 老司机影院毛片| av又黄又爽大尺度在线免费看 | 婷婷色麻豆天堂久久 | 久久这里有精品视频免费| 蜜臀久久99精品久久宅男| 美女内射精品一级片tv| 国产激情偷乱视频一区二区| 精品熟女少妇av免费看| or卡值多少钱| 又爽又黄无遮挡网站| 久久久久久久久久久免费av| 我的女老师完整版在线观看| 国产精品一区www在线观看| 中国美白少妇内射xxxbb| 精品久久久久久电影网 | 亚洲电影在线观看av| 欧美性猛交黑人性爽| 欧美成人精品欧美一级黄| 人体艺术视频欧美日本| or卡值多少钱| 久久精品国产亚洲av天美| 国产免费一级a男人的天堂| 一级二级三级毛片免费看| 亚洲天堂国产精品一区在线| 亚洲真实伦在线观看| videossex国产| 天堂av国产一区二区熟女人妻| 热99在线观看视频| 成年av动漫网址| 97人妻精品一区二区三区麻豆| 麻豆精品久久久久久蜜桃| 久久久亚洲精品成人影院| 亚洲国产精品成人久久小说| 成人毛片a级毛片在线播放| 男人和女人高潮做爰伦理| 插逼视频在线观看| 午夜免费激情av| 一个人看的www免费观看视频| 欧美激情在线99| 国产精品乱码一区二三区的特点| 日日摸夜夜添夜夜爱| 久久久久国产网址| 一个人观看的视频www高清免费观看| 国产精品乱码一区二三区的特点| 99国产精品一区二区蜜桃av| 一个人看的www免费观看视频| 国产免费视频播放在线视频 | 婷婷色综合大香蕉| 亚洲人成网站在线播| 亚洲精品成人久久久久久| 国产一级毛片在线| av国产免费在线观看| 婷婷色av中文字幕| 亚洲欧美精品专区久久| 干丝袜人妻中文字幕| av播播在线观看一区| 久久精品久久久久久久性| 国产精品福利在线免费观看| 欧美成人精品欧美一级黄| 少妇人妻一区二区三区视频| 亚洲av电影不卡..在线观看| 免费播放大片免费观看视频在线观看 | 日韩一区二区视频免费看| 国产高清视频在线观看网站| 美女被艹到高潮喷水动态| 精品久久国产蜜桃| 国产在线男女| 深爱激情五月婷婷| 免费看光身美女| 免费观看性生交大片5| 两性午夜刺激爽爽歪歪视频在线观看| 91在线精品国自产拍蜜月| 日韩,欧美,国产一区二区三区 | av.在线天堂| 国产亚洲av嫩草精品影院| 天堂中文最新版在线下载 | 亚洲国产欧美在线一区| 国产探花在线观看一区二区| 亚洲国产欧美人成| 国产精品不卡视频一区二区| 日本av手机在线免费观看| 亚洲无线观看免费| 久久6这里有精品| 成年版毛片免费区| av国产久精品久网站免费入址| 三级男女做爰猛烈吃奶摸视频| 97热精品久久久久久| 欧美日韩综合久久久久久| 一边摸一边抽搐一进一小说| 大话2 男鬼变身卡| 少妇人妻精品综合一区二区| 91aial.com中文字幕在线观看| 一级av片app| 精品久久久久久久久亚洲| 亚洲国产成人一精品久久久| 亚洲伊人久久精品综合 | 精品国产三级普通话版| 99久久成人亚洲精品观看| 国产成人91sexporn| 级片在线观看| 国产精品一区二区三区四区久久| 日韩一区二区视频免费看| 国产亚洲最大av| 国产精品人妻久久久久久| 欧美一区二区亚洲| 伦理电影大哥的女人| 欧美xxxx黑人xx丫x性爽| 国产熟女欧美一区二区| 在线免费十八禁| 日韩精品青青久久久久久| 国产亚洲一区二区精品| 国内揄拍国产精品人妻在线| 老师上课跳d突然被开到最大视频| 亚洲av免费在线观看| 国产淫语在线视频| 99九九线精品视频在线观看视频| 日韩国内少妇激情av| 女人十人毛片免费观看3o分钟| 国产精品1区2区在线观看.| 一个人观看的视频www高清免费观看| 免费av观看视频| 国产精品一区二区三区四区久久| 午夜免费男女啪啪视频观看| 色网站视频免费| 亚洲国产成人一精品久久久| 午夜精品一区二区三区免费看| av播播在线观看一区| 非洲黑人性xxxx精品又粗又长| 精品不卡国产一区二区三区| 丝袜美腿在线中文| 日本熟妇午夜| 有码 亚洲区| 直男gayav资源| 精华霜和精华液先用哪个| 国产69精品久久久久777片| 日本五十路高清| 国产中年淑女户外野战色| 日本免费一区二区三区高清不卡| 精品欧美国产一区二区三| h日本视频在线播放| 国产精品精品国产色婷婷| 欧美高清成人免费视频www| or卡值多少钱| 看片在线看免费视频| 男女下面进入的视频免费午夜| 亚洲av不卡在线观看| 五月伊人婷婷丁香| 少妇熟女欧美另类| 波多野结衣高清无吗| 秋霞在线观看毛片| 纵有疾风起免费观看全集完整版 | 3wmmmm亚洲av在线观看| 韩国av在线不卡| 欧美精品一区二区大全| 男插女下体视频免费在线播放| 国产老妇女一区| 国产亚洲av嫩草精品影院| 日本-黄色视频高清免费观看| 日本一二三区视频观看| 成年女人看的毛片在线观看| 亚洲四区av| 我的老师免费观看完整版| 美女cb高潮喷水在线观看| 少妇熟女aⅴ在线视频| 一二三四中文在线观看免费高清| 亚洲人成网站高清观看| 久久婷婷人人爽人人干人人爱| 美女高潮的动态| 国产 一区精品| 亚洲欧美精品综合久久99| 亚洲精品成人久久久久久| 欧美另类亚洲清纯唯美| 少妇丰满av| 91精品国产九色| 神马国产精品三级电影在线观看| 成年免费大片在线观看| 汤姆久久久久久久影院中文字幕 | 一本一本综合久久| 国产极品精品免费视频能看的| 最后的刺客免费高清国语| 小说图片视频综合网站| 亚洲综合精品二区| 欧美区成人在线视频| 天美传媒精品一区二区| 久久久欧美国产精品| 亚洲内射少妇av| 国产色爽女视频免费观看| 久久精品夜夜夜夜夜久久蜜豆| 美女高潮的动态| 日产精品乱码卡一卡2卡三| 久久久国产成人免费| 亚洲第一区二区三区不卡| 国产成人a∨麻豆精品| 99久久九九国产精品国产免费| 日本爱情动作片www.在线观看| 久久久久免费精品人妻一区二区| 色哟哟·www| 一级毛片aaaaaa免费看小| 中文精品一卡2卡3卡4更新| 国产在视频线在精品| 国内揄拍国产精品人妻在线| 91精品国产九色| 91av网一区二区| 国产v大片淫在线免费观看| 国产av在哪里看| 午夜激情福利司机影院| 亚洲精品亚洲一区二区| 久99久视频精品免费| 欧美一区二区精品小视频在线| 国产真实乱freesex| 亚洲真实伦在线观看| 亚洲国产成人一精品久久久| 久久99热这里只有精品18| 国产在视频线在精品| 国产精品av视频在线免费观看| 国产亚洲av片在线观看秒播厂 | 亚洲aⅴ乱码一区二区在线播放| 国内精品美女久久久久久| 久久99热6这里只有精品| 久久久成人免费电影| 精品久久久久久久人妻蜜臀av| 日韩av在线大香蕉| 22中文网久久字幕| 内地一区二区视频在线| 麻豆国产97在线/欧美| 国产精品国产高清国产av| av女优亚洲男人天堂| 国产一区二区在线观看日韩| 国产精品电影一区二区三区| 两个人视频免费观看高清| 1000部很黄的大片| 美女高潮的动态| 最新中文字幕久久久久| av线在线观看网站| www日本黄色视频网| 国产在视频线在精品| 国产欧美另类精品又又久久亚洲欧美| 国产av一区在线观看免费| 国产成人免费观看mmmm| 久久久精品大字幕| 99久久中文字幕三级久久日本| 99热精品在线国产| 日韩强制内射视频| 成人欧美大片| 综合色av麻豆| 免费看光身美女| 午夜老司机福利剧场| av在线蜜桃| 韩国av在线不卡| 黄片wwwwww| 高清在线视频一区二区三区 | 毛片女人毛片| 国产一级毛片七仙女欲春2| 插阴视频在线观看视频| 插逼视频在线观看| 99久久精品国产国产毛片| 永久免费av网站大全| 三级毛片av免费| 欧美日韩在线观看h| 成人亚洲精品av一区二区| 日本wwww免费看| 国产色爽女视频免费观看| 九色成人免费人妻av| 亚洲成av人片在线播放无| 综合色丁香网| 久久精品国产亚洲av涩爱| 国产成人a∨麻豆精品| 高清午夜精品一区二区三区| 一本一本综合久久| 亚洲激情五月婷婷啪啪| 久久人人爽人人片av| 亚洲精品日韩在线中文字幕| 免费黄色在线免费观看| 男人舔奶头视频| 久久草成人影院| 亚洲伊人久久精品综合 | 91aial.com中文字幕在线观看| 亚洲在线自拍视频| 丝袜美腿在线中文| 波多野结衣高清无吗| 老女人水多毛片| 观看美女的网站| 成人毛片a级毛片在线播放| 色尼玛亚洲综合影院| 深爱激情五月婷婷| 1000部很黄的大片| 色尼玛亚洲综合影院| 一级毛片电影观看 | 国产高清有码在线观看视频| 亚洲精品aⅴ在线观看| 久久精品国产鲁丝片午夜精品| 我要看日韩黄色一级片| 99视频精品全部免费 在线| 亚洲精品456在线播放app| 国产一区二区亚洲精品在线观看| 99视频精品全部免费 在线| 久热久热在线精品观看| 成人特级av手机在线观看| 欧美区成人在线视频| 欧美成人精品欧美一级黄| 久久久久久大精品| 成年女人看的毛片在线观看| 亚洲人成网站在线播| 欧美性猛交╳xxx乱大交人| 国产免费一级a男人的天堂| 别揉我奶头 嗯啊视频| 青春草亚洲视频在线观看| 免费观看在线日韩| 蜜桃久久精品国产亚洲av| 国产精品日韩av在线免费观看| 中文亚洲av片在线观看爽| 亚洲av电影不卡..在线观看| 精品久久国产蜜桃| 小说图片视频综合网站|