999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

考慮儲層改造體積頁巖氣藏復合模型

2015-08-01 02:38:16中國石油大學華東地質資源與地質工程博士后流動站山東青島266580中國石油大學華東石油工程學院山東青島266580
東北石油大學學報 2015年2期
關鍵詞:區域模型

(1.中國石油大學(華東)地質資源與地質工程博士后流動站,山東青島 266580; 2.中國石油大學(華東)石油工程學院,山東青島 266580)

在時間上,對方程式(10)的求解采用隱式向后差分格式,并利用上一個時間步k時刻的基質擬壓力值,得到時間步k+1時刻裂縫系統的擬壓力值,再代入基質的控制方程式(8),得到k+1時刻基質的擬壓力值;因此,首先計算裂縫系統在k+1時刻的擬壓力值為

(1.中國石油大學(華東)地質資源與地質工程博士后流動站,山東青島 266580; 2.中國石油大學(華東)石油工程學院,山東青島 266580)

針對頁巖氣藏中水平井結合體積壓裂開采、吸附氣和游離氣共存的方式,建立考慮儲層改造體積的頁巖氣藏復合模型,定義新的參數表征基質中吸附解吸氣量與游離氣彈性釋放量的比值,且將儲層分為人工主裂縫區域、儲層改造區域和未改造區域,其中人工主裂縫基于離散裂縫模型降維處理,儲層改造區域為雙孔雙滲模型,未改造區域為單孔隙介質模型;模型采用有限元方法進行求解,與雙重介質解析解對比驗證算法的正確性.結果表明:頁巖氣藏水平井體積壓裂復合模型主要存在主裂縫周圍線性流、過渡區域擬穩態、竄流階段、未改造區域的擬徑向流動和到達邊界后的擬穩態等5個主要流動階段,且考慮吸附解吸后,定產量生產所需壓差小,壓力波傳播到邊界時間長,壓力導數曲線凹槽更加明顯,定井底流壓生產時壓裂水平井產量更大,穩產時間更長;儲層改造體積越大,到達區域擬穩態流越晚,可判定儲層改造體積;Langmuir吸附體積越大,壓力波傳播越慢,所需壓差越小,壓力導數曲線凹槽越深,頁巖氣藏穩產時間越長,產量越大,但產量的增幅越來越小.

頁巖氣藏;體積壓裂;儲層改造區域;復合模型;有限元法

0 引言

頁巖氣藏資源量豐富、潛力巨大已成為研究的熱點[1-2],頁巖內吸附氣與游離氣共存[3-5],不同于常規的氣藏(游離氣),也不同于煤層氣(主要以吸附狀態賦存于煤基質孔隙)[6].因此,在頁巖氣藏中,需要同時考慮游離氣的彈性釋放和吸附氣的解吸脫附;隨著水平井鉆井、同步壓裂、微地震裂縫診斷等技術的進步與應用[7-8],頁巖氣藏的開采增長迅速.2006年,Mayerhofer M等[9-10]采用微地震技術研究Banett頁巖中壓裂裂縫變化時,首次提出儲層改造體積,并指出增大儲層改造體積的技術思路與常規的雙翼對稱裂縫不同,體積壓裂方式在形成一條或多條主裂縫的同時,對天然裂縫、巖石層理進行溝通,從而實現儲層在長、寬、高三維方向的全面改造[11].目前,頁巖氣藏的研究大部分基于單孔隙或多重孔隙介質模型[12-14],不能準確刻畫頁巖氣藏儲層改造體積的形成.為了更好地研究頁巖氣藏體積壓裂的滲流特征,針對頁巖氣藏構建體積壓裂復合模型,筆者把模型分為3個不同的主要滲流區域,并對模型進行求解,分析考慮儲層改造體積壓裂水平井流動形態及各因素對流動形態的影響,為頁巖氣藏開采和產能預測提供有效的技術支持.

1 體積壓裂復合

1.1 物理模型

頁巖氣藏儲層孔隙度滲透率較低,需采用水平井鉆井增大含油氣層與井筒的接觸面積,并利用多級水力壓裂技術,在主裂縫周圍開啟天然裂縫形成儲層改造區域,因此,對于頁巖氣藏考慮吸附解吸過程,建立復合油藏分段壓裂水平井模型(SRV為主裂縫周圍區域,見圖1).假定三維盒狀封閉油氣藏內有一口水平井,儲層在x、y、z方向的長度分別為Xe、Ye、h,水平井長度為L,對水平井進行人工壓裂改造,主裂縫條數為Nf.

頁巖氣藏復合模型將儲層分為人工主裂縫、儲層改造區和未改造區3個滲流區域,且各個區域滿足的基本假設條件:(1)人工主裂縫,有限導流且服從達西滲流規律.(2)儲層改造體積,次生裂縫被開啟且相互連接,區域由頁巖基質和微裂縫介質構成,滿足雙重介質模型.(3)頁巖基巖表面吸附大量的甲烷氣體,且氣體的吸附量與壓力滿足Langmuir等溫吸附公式.(4)頁巖基質系統中,氣體以游離相和吸附相形式共存,壓力下降時氣體排出包括游離氣膨脹和吸附氣解吸兩部分.(5)微裂縫系統內氣體主要以游離相的形式出現,為主要運移通道,遵循達西滲流規律.(6)未改造區域.該區域頁巖基質滲透率和孔隙度與改造區域內頁巖基質系統內的相同.(7)水平井裸眼或射孔完井,當定總產量生產時,分析井底壓力動態;當定井底流壓時,計算水平井的不穩定產能.

1.2 數學模型

圖1 復合油藏分段壓裂水平井物理模型Fig.1 Physical model diagram of fractured horizontal well in composite reservoir

其中,ΨmD、ΨpD和ΨsD分別為頁巖基質系統、主裂縫系統和SRV區域次生裂縫系統的擬壓力值;Ψ0為氣藏的初始擬壓力值;λ為竄流系數;η為基巖與微裂縫總的彈性儲能;ω為彈性儲能比;Tsc為標準狀況下溫度;psc為標準狀況下壓力;q為定產量生產時氣體流量;h為氣藏厚度;s、m分別為次生裂縫系統和頁巖基質系統的孔隙度;Ks、Km分別為次生裂縫系統和頁巖基質系統的滲透率;μ為氣體的黏度;t為時間;Cg為氣體的壓縮系數;α為竄流的形狀因子;tD為無因次時間;xD、yD、zD分別為在x、y、z方向的無因次距離.

1.2.1 人工主裂縫

假定氣體在主裂縫內服從達西滲流規律,主裂縫與基巖交界處的壓力處處相等,則主裂縫的無因次數學模型可表示為

式中:Kp為人工主裂縫系統的滲透率;qpD=qp/qt為人工主裂縫流入水平井筒的無因次流量值;p為人工主裂縫的孔隙度;Cp、Cm、Cs分別為人工主裂縫、基巖和微裂縫系統的壓縮系數;MD、M′D分別為射孔處和地層任意一點處的無因次位置;Ωp為人工裂縫邊界;hD為無因次厚度.

1.2.2 儲層改造區域

定義一個新的參數β,表示單位時間內由于基質壓力降低,基質骨架吸附解吸氣量與游離氣彈性能釋放量的比值,即

式中:Qads為基質單位骨架的甲烷吸附量;ρm、ρs分別為基質系統和人工裂縫系統氣體的密度;Vstd為氣體在標準狀況下的摩爾體積;R為理想氣體常數;T為氣藏溫度;VL為Langmuir體積;pL為Langmuir壓力;Z為氣體壓縮因子;pm為頁巖基質系統壓力.在計算過程中參數采用上一步的基質壓力,因此,對于某一具體時間步參數β是已知的.

假定儲層改造區域滿足雙孔雙滲模型,流體在次生裂縫和基巖系統中滿足達西滲流規律,介質系統之間為擬穩態竄流,則該區域的無因次數學模型為

1.2.3 儲層未改造區域

未改造區域內的基巖與SRV區域雙重介質內的基巖系統相同,滿足達西滲流規律,基巖系統內邊界壓力與雙重介質系統內次生裂縫系統的壓力處處相等,外邊界為封閉邊界,則未改造區域無因次數學模型為

2 有限元求解及驗證

2.1 有限元求解

對盒狀頁巖氣藏水平井體積壓裂復合模型,采用非結構化網格進行剖分,利用有限元方法對模型進行求解,其中主裂縫基于離散裂縫模型采用二維三角形單元,采用四面體單元剖分儲層改造區域和未改造區域,分別對3個主要滲流區域進行單元特性分析.

(1)人工主裂縫基于離散裂縫模型進行顯式處理[15-16],將三維的板狀裂縫化為二維的裂縫面,采用二維三角形單元進行剖分,得到二維主裂縫面單元特性矩陣為

(2)儲層改造區域微裂縫和基質系統的四面體單元特性矩陣為

同理,未改造區域內的基質系統單元特性矩陣為

將氣藏整個區域內的單元特性矩陣和列陣進行組合得到氣藏的整體矩陣和列陣.首先根據結點局部序號的對應關系,將單元特性矩陣里的元素移置到與整體序號對應的位置上,并在沒有元素的位置設置為0.記結點總數為Np,則基質和微裂縫的結點壓力值可表示為Ψm=[Ψm,1,Ψm,2,…Ψm,Np]T,Ψs=[Ψs,1,Ψs,2,…Ψs,Np]T,其中只有儲層改造區域才存在微裂縫系統,故未改造區域處雙重介質系統微裂縫的壓力值Ψs,i=0.為了方便起見,當結點屬于儲層改造區域時方程組記為

當結點屬于儲層未改造區域時方程組化為

將人工主裂縫單元矩陣和列陣與雙重介質微裂縫系統進行耦合,得到儲層改造區域內微裂縫系統的壓力計算方程組為

其中,裂縫系統整體矩陣和列陣的計算化為人工主裂縫與微裂縫的組合,即

在時間上,對方程式(10)的求解采用隱式向后差分格式,并利用上一個時間步k時刻的基質擬壓力值,得到時間步k+1時刻裂縫系統的擬壓力值,再代入基質的控制方程式(8),得到k+1時刻基質的擬壓力值;因此,首先計算裂縫系統在k+1時刻的擬壓力值為

再計算基質系統不同區域的k+1時刻的擬壓力值為

2.2 正確性驗證

由于一般的解析方法無法直接得到體積壓裂復合模型壓裂水平井的解,為了驗證有限元算法的正確性,將模型退化為雙重介質壓裂水平井模型,并與經典Zerzar模型[17]中壓裂水平井解析解進行對比.盒狀封閉氣藏水平井和人工裂縫的無因次參數為:氣藏大小XeD=12,YeD=12,hD=0.1,水平井長度LD=1,5條裂縫在X方向無因次坐標(-0.4,-0.2,0,0.2,0.4),雙重介質參數竄流系數λ=1,彈性儲能比ω=0.1,無因次裂縫半長xfD=0.05.假設只有人工裂縫向水平井井筒供液,對比文中模型計算的有限元數值解與Zerzar模型中壓裂水平井解析解,得到壓裂水平井流動形態見圖2.由圖2可知:兩模型的壓力及壓力導數基本相同,表現出裂縫線性流、裂縫徑向流、地層線性流、竄流階段及擬徑向流動階段,說明算法正確.

圖2 有限元數值解與Zerzar解析解對比Fig.2 Comparing Zerzar classic analytical solution with our model by finite element method

3 水平井動態分析

為了討論頁巖氣藏復合模型下壓裂水平井的壓力及產能動態,盒狀封閉氣藏水平井和人工裂縫基本參數,油藏大小XeD=12,YeD=12,hD=0.1;SRV區域,水平井長度LD=1,5條裂縫在X方向無因次坐標(-0.4,-0.2,0,0.2,0.4),無因次裂縫半長xfD=0.15,假設人工裂縫全部穿透儲層,主裂縫滲透率與次生裂縫滲透率的比值Kp/Ks=100,主裂縫開度αp=10-3,SRV區域雙重介質竄流系數λ=0.5,彈性儲能比ω=0.1,基巖系統孔隙度m=0.050,次生裂縫系統孔隙度f=0.005,氣藏初始壓力p0=10MPa,定產氣量q=0.1m3/s,次生裂縫系統滲透率Kf=10-15m2,基巖系統滲透率Km=10-18m2,則基巖與次生裂縫滲透率的比值Km/Ks=0.001.

頁巖氣藏基本參數:甲烷的摩爾質量Mg=0.016kg/mol,甲烷的摩爾體積Vstd=0.022 37m3/mol,密度ρs=2 600kg/m3,理想氣體常數R=8.314J/(K·mol),油藏溫度T=323K,標準狀況下溫度Tsc=273 K,標準狀況下壓力psc=1.013 25×105Pa.

3.1 吸附解吸特征

為了研究吸附解吸特性對頁巖氣藏動態的影響,分別計算不考慮吸附解吸特性(ρ=0)和考慮吸附解吸特性時壓裂水平井動態特征,頁巖氣藏的基本參數不變,取頁巖對甲烷氣體的Langmuir吸附體積VL=5×10-4m3/kg,Langmuir吸附壓力pL=12MPa,分別計算定產量生產時水平井的壓力動態及定井底流壓生產時水平井的產量動態(見圖3).

頁巖氣藏復合模型下壓裂水平井主要分為5個流動階段(見圖4):(1)SRV區域主裂縫周圍線性流.由于SRV區域微裂縫張開,滲透率較大,在流動早期,主要為SRV區域微裂縫流體流向人工主裂縫,表現為人工裂縫周圍的線性流動,壓力導數曲線表現為斜率1/2的直線.(2)過渡區域擬穩態.隨著SRV區域壓力的降低,壓力波往外傳播,到達改造與未改造區域的邊界時,由于未改造區域的滲透率小,按原來的壓降速度無法滿足定產量生產要求,因此,有類似于封閉邊界的形態,井底壓差快速上升,出現區域的擬穩態流動,壓力和壓力導數表現為斜率1的直線段.(3)竄流階段.在復合油藏出現過渡區域擬穩態的同時,SRV區域的基質系統也開始向微裂縫系統供液,表現為基質系統向裂縫系統的竄流.(4)未改造區域的擬徑向流動.若油藏區域足夠大時,壓力波傳播到離SRV較遠處,在試井曲線上表現為以儲層改造區域為中心的擬徑向流動,壓力導數曲線表現為平行于橫坐標的直線.(5)擬穩態流動階段.由于外邊界封閉,壓力波傳播到邊界,流動達到擬穩態流動階段,壓力降迅速增加,壓力導數為斜率1的直線段.

圖3 考慮和不考慮吸附解吸時復合氣藏壓力和產能曲線Fig.3 Pressure and rate curves of wether consider the adsorption and desorption or not in composite shale gas reservoir

圖4 體積壓裂復合模型水平井主要流動形態Fig.4 Main flow stages of fractured horizontal well in the composite model

由圖3可知:當在頁巖氣藏中考慮吸附解吸的影響時,對流動早期的影響不明顯;當SRV區域的流體發生竄流時,由于吸附解吸氣體補足,使得雙重介質模型中次生裂縫系統所占的彈性儲能比例減小,壓力導數曲線下凹更加明顯,地層過渡到徑向流動的時間越晚,壓力波傳播到外部封閉邊界的時間也越延遲;在產能曲線上可見考慮吸附解吸后,在竄流階段之后氣藏中吸附氣體開始解吸,使得壓裂水平井氣井產能的穩產時間更長.因此,頁巖氣藏產能評價過程中必須考慮吸附解吸的影響.

3.2 SRV大小

為了說明儲層改造體積對頁巖氣藏壓裂水平井流動形態的影響,在其他參數不變的情況下,分別取無因次裂縫半長xfD=0.05,0.10,0.15進行計算,得到頁巖氣藏體積壓裂復合模型在不同儲層改造體積下,壓裂水平井壓力和壓力導數曲線見圖5.由圖5可知:在其他參數一定的情況下,儲層改造區域SRV大小影響壓力波到達過渡區域擬穩態流動的時間,儲層改造區域越大,壓力波傳播到邊界所需要的時間越長,進入區域擬穩態的時間越早,因此可根據早期頁巖氣藏體積壓裂后壓力導數曲線上出現斜率為1的直線段時間判定SRV的大小.

3.3 Langmuir吸附體積

頁巖對甲烷的吸附量滿足Langmuir吸附等溫公式.在Langmuir壓力不變的情況下,隨著Langmuir吸附體積VL的增大,頁巖對甲烷的吸附能力增強.為了研究頁巖Langmuir吸附體積對壓裂水平井壓力及產量的影響,在Langmuir壓力pL=5MPa時,取3個不同Langmuir體積VL,即VL=(2,6,10)× 10-4m3/kg,分別計算定產量生產時壓力及壓力導數曲線及定井底流壓生產時的產量曲線(見圖6).由圖6可知:隨著Langmuir吸附體積增大,SRV內次生裂縫系統和頁巖基質系統發生竄流時,擬穩態竄流的凹槽更深,到達未改造區域過渡到擬徑向流的時間越晚;同時,壓力波傳播到邊界的時間也越晚.隨著Langmuir吸附體積的增大,發生竄流以后的產能越大,產量延遲的時間也越長;隨著Langmuir吸附體積的線性增加,壓力波影響和產能穩產時間的增幅逐漸減小.

圖5 不同SRV大小時復合油藏壓裂水平井壓力動態分布Fig.5 Dynamic pressure curves of composite fractured horizontal well model with different SRV sizes

圖6 不同Langmuir體積下頁巖氣藏復合壓力和產能曲線Fig.6 Pressure and rate response affected by Langmuir absorption volume in composite shale gas reservoir

4 結論

(1)頁巖氣藏水平井體積壓裂復合模型將儲層分為人工主裂縫、儲層改造和未改造3個區域,采用有限元方法求解,得到模型的存在SRV主裂縫周圍線性流、過渡區域擬穩態、竄流階段、未改造區域的擬徑向流動和到達邊界后的擬穩態5個流動階段.

(2)吸附解吸過程對頁巖氣藏分段壓裂水平井壓力和產能的影響表現在基質與天然裂縫發生竄流后,且考慮吸附解吸過程后,水平井定產量生產所需的壓降值比不考慮時要小,壓力波到達邊界的時間更長,壓力導數曲線凹槽更加明顯;同時,考慮吸附解吸過程后,頁巖氣藏壓裂水平井產能越大,穩產時間更長.

(3)儲層改造區域越大,則壓力波到達過渡區域擬穩態時間越晚,在測試過程中可根據壓力波到達區域擬穩態的時間,判斷壓力波到達儲層改造區域邊界,計算儲層改造體積的大小.

(4)Langmuir吸附體積越大,頁巖氣藏發生竄流后定產量生產所需壓差越小,壓力波傳播到邊界越晚,壓力導數曲線下凹更加明顯;考慮吸附解吸后定井底流壓生產時壓裂水平井產能越大,穩產時間也越長;隨著Langmuir吸附體積線性增大,頁巖氣藏壓裂水平井壓力影響幅度逐漸越小,產能增加幅度也越小.

):

[1] U S Energy Information Administration.World shale gas resources:An initial assessment of 14regions outside the united states[M].Washington:US Department of Energy,2011.

[2] 陳尚斌,朱炎銘,王紅巖,等.中國頁巖氣研究現狀與發展趨勢[J].石油學報,2010,31(4):689-694.

Chen Shangbin,Zhu Yanming,Wang Hongyan,et al.Status and trends of shale gas in China[J].Acta Petrolei Sinica,2010,31(4):689-694.

[3] Hill D G.Nelson CR gas productive fractured shales:An overview and update[J].Gas Tips,2000,6(2):4-18.

[4] John P V.Geomechanical studies of the barnett shale[D].Texas:Stanford University,2011.

[5] 孫海,姚軍,孫致學,等.頁巖氣數值模擬技術進展及展望[J].油氣地質與采收率,2012,19(1):46-49.

Sun Hai,Yao Jun,Sun Zhixue,et al.Recent development and prospect on numerical simulation of shale gas reservoirs[J].Petroleum Geology and Recovery Efficiency,2012,19(1):46-49.

[6] 陳霞,劉洪林,王紅巖,等.沁水盆地含水煤層氣藏的氣體滲流特性[J].石油學報,2011,32(3):500-503.

Chen Xia,Liu Honglin,Wang Hongyan,et al.Gas-seepage characteristics of coalbed reservoirs with different water saturation in the Qinshui basin[J].Acta Petrolei Sinica,2011,32(3):500-503.

[7] 鄒才能,朱如凱,吳松濤,等.常規與非常規油氣聚集類型、特征、機理及展望——以中國致密油和致密氣為例[J].石油學報,2012,33(2):173-187.

Zou Caineng,Zhu Rukai,Wu Songtao,et al.Types,characteristics,genesis and prospects of conventional and unconventional hydrocarbon accumulations:Taking tight oil and tight gas in China as an instance[J].Acta Petrolei Sinica,2012,33(2):173-187.

[8] 張抗.頁巖油氣發展的中國之鑒——致密油氣和煤層氣[J].天然氣工業,2013,33(4):18-25.

Zhang Kang.Lessons for shale oil &gas development from that of tight oil &gas and coalbed methane gas in China[J].Natural Gas Industry,2013,33(4):18-25.

[9] Mayerhofer M,Lolon E,Youngblood J,et al.Integration of microseismic fracture mapping results with numerical fracture network production modeling in the Barnett shale[C].SPE Annual Technical Conference and Exhibition,2006.

[10] Mayerhofer M,Lolon E,Warpinski N,et al.What is Stimulated rock volume[C].SPE Shale Gas Production Conference,2008.

[11] Waters G,Dean B,Downie R,et al.Simultaneous hydraulic fracturing of adjacent horizontal wells in the woodford shale[C].SPE Hydraulic Fracturing Technology Conference,2009.

[12] 劉曉旭,楊學鋒,陳遠林,等.頁巖氣分段壓裂水平井滲流機理及試井分析[J].天然氣工業,2013,33(12):1-5.

Liu Xiaoxu,Yang Xuefeng,Chen Yuanlin,et al.Seepage behavior and well testing in horizontal shale gas wells under multi-stage fracking[J].Natural Gas Industry,2013,33(12):1-5.

[13] Zhao Yulong,Zhang Liehui,Zhao Jinzhou,et al."Triple porosity"modeling of transient well test rate decline analysis for multi-fractured horizontal well in shale gas reservoirs[J].Journal of Petroleum Science and Engineering,2013(110):253-262.

[14] Guo Jingjing,Zhang Liehui,Wang Haitao,et al.Presssure transient analysis for multi-stage fractured horizontal wells in shale gas reservoirs[J].Transport in Porous Media,2012,93(3):635-653.

[15] 姚軍,王子勝,張允,等.天然裂縫性油藏的離散裂縫網絡數值模擬方法[J].石油學報,2010,31(2):284-288.

Yao Jun,Wang Zisheng,Zhang Yun,et al.Numerical simulation method of discrete fracture network for naturally fractured reservoirs[J].Acta Petrolei Sinica,2010,31(2):284-288.

[16] Yao Jun,Huang Zhaoqin,Li Yajun,et al.Discrete fracture-vug network model for modeling fluid flow in fractured vuggy porous media[C].SPE 130287,2010.

[17] Zerzar A,Bettam Y.Interpretation of multiple hydraulically fractured horizontal wells in closed systems[C].SPE 84888,2003.

考慮儲層改造體積頁巖氣藏復合模型

樊冬艷1,2,姚 軍2,金 強1,孫 海2,曾 慧2

DOI 10.3969/j.issn.2095-4107.2015.02.010

TE332

A

2095 4107(2015)02 0077 08

2014 13 30;編輯:關開澄

國家自然科學基金重點項目(51234007);山東省自然科學基金項目(ZR2014EL016,ZR2014EEP018);中國博士后科學基金項目(2014M551989)

樊冬艷(1985-),女,博士,主要從事水平井及壓裂水平井滲流理論及產能方面的研究.

book=84,ebook=87

猜你喜歡
區域模型
一半模型
永久基本農田集中區域“禁廢”
今日農業(2021年9期)2021-11-26 07:41:24
分割區域
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
關于四色猜想
分區域
FLUKA幾何模型到CAD幾何模型轉換方法初步研究
基于嚴重區域的多PCC點暫降頻次估計
電測與儀表(2015年5期)2015-04-09 11:30:52
主站蜘蛛池模板: 好吊色妇女免费视频免费| 成人毛片免费观看| 无码日韩精品91超碰| 国产二级毛片| 欧美在线网| 欧洲亚洲一区| 手机在线免费不卡一区二| 亚亚洲乱码一二三四区| 亚洲一区毛片| 国产精品国产主播在线观看| 亚洲91在线精品| 国产91特黄特色A级毛片| 狠狠亚洲婷婷综合色香| 91黄视频在线观看| 99热这里只有精品国产99| 免费毛片在线| 日日拍夜夜操| 少妇高潮惨叫久久久久久| 国产91久久久久久| 国内精品久久久久鸭| 欧美一级在线看| 亚洲国产精品国自产拍A| 99久久亚洲精品影院| 国产午夜无码片在线观看网站| 狠狠亚洲五月天| 婷婷六月色| 9999在线视频| 91视频免费观看网站| 国产成人精品高清在线| 91精品国产自产在线观看| 色成人综合| 国产成人三级在线观看视频| 97亚洲色综久久精品| 国产一区二区三区精品久久呦| 国产免费福利网站| 2021国产精品自产拍在线观看| 亚洲成A人V欧美综合天堂| 国产剧情国内精品原创| 青青操国产| 久久男人资源站| 97超爽成人免费视频在线播放| 色综合久久久久8天国| 国产一区二区免费播放| 国产高清国内精品福利| 国产在线精彩视频二区| 日韩高清中文字幕| 五月丁香伊人啪啪手机免费观看| 欧美a在线| 国产91在线|日本| 五月天在线网站| 香港一级毛片免费看| 青青热久免费精品视频6| 波多野结衣国产精品| 手机在线免费不卡一区二| 亚洲精品片911| 精品偷拍一区二区| 小说区 亚洲 自拍 另类| 99热最新在线| 美女视频黄又黄又免费高清| 在线免费观看AV| 久久久精品国产SM调教网站| 毛片在线区| 久久久久国产一级毛片高清板| 自慰高潮喷白浆在线观看| 久久综合九色综合97婷婷| 91久久性奴调教国产免费| 久青草免费在线视频| 亚洲 日韩 激情 无码 中出| 91久久大香线蕉| 成人国产精品一级毛片天堂 | 嫩草国产在线| 国产激情无码一区二区APP | 亚洲国产高清精品线久久| 婷婷六月激情综合一区| 国产精品手机在线播放| 国产xx在线观看| 凹凸精品免费精品视频| 天天激情综合| 久久亚洲国产一区二区| 国内黄色精品| 狠狠干综合| 国产91特黄特色A级毛片|