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

管道中連續振動信號衰減模型與傳遞特性研究

2016-01-15 03:43:03劉均,袁峰
振動與沖擊 2015年19期

管道中連續振動信號衰減模型與傳遞特性研究

劉均1,2,袁峰1

(1. 哈爾濱工業大學電氣工程及自動化學院,哈爾濱150001; 2. 東北石油大學電氣信息工程學院,大慶163318)

摘要:借鑒電力線輸電理論中的阻抗和傳遞系數來描述連續壓力波動在管道中的傳遞過程,建立了基于傳遞矩陣的連續波衰減模型,利用特征阻抗來描述管道本體特征對波動傳遞的影響,并以此模型為基礎,分析了連續波在管道中傳遞時壓力和流量之間的關系、以及管道終端阻抗和傾角對波動傳遞的影響。通過對模型的仿真計算,得出了在管道中連續波信號傳遞時,管道內壓力波呈駐波分布的結論,同時分析了在不同終端阻抗下波動的分布情況;通過對管道內不同頻率信號的衰減情況分析,繪制了管道內信號傳遞的幅頻特性曲線,得出了隨著頻率的升高,管道內信號呈波動性衰減的結論。最后通過地面驗證實驗,說明了不同頻率信號在相同管道內的傳遞情況符合仿真分析,為隨鉆測量過程中通信頻率的選取提供了理論支持。

關鍵詞:幅頻特性;隨鉆測量;傳遞系數;水力阻抗;壓力波動

中圖分類號:TE242.9

文獻標志碼:A

DOI:10.13465/j.cnki.jvs.2015.19.013

Abstract:Here, impedance and transmission coefficient of the power line transmission theory were used to describe continuous pressure vibrations in transfer process of a pipeline, the continuous wave attenuation model was built based on the transfer matrix method. The characteristic impedance was adopted to describe the influence of the pipeline’s property on the wave propagation. According to this wave attenuation model, the relationship between pressure and flow was analyzed when the continuous waves moved along the pipeline, and the influence of terminal impedance and dip angle of the pipeline on the wave transmission was studied. Through simulation of the model, it was shown that the pressure waves are the standing wave distribution when the continuous wave signals transmit in the pipeline; at the same time, the transmission volatility under different terminal impedances is analyzed; the amplitude-frequency characteristic curve of the signal transmission in the pipe is mapped by analyzing the attenuation situations of different frequency signals in the pipeline, so the signal in the pipe is in the attenuation volatility with increase in frequency. The tests on the ground showed that the situations about signals of different frequencies transmission in the same pipe agreed well with those using the simulation analysis. The results provided a theoretical support for selecting the frequency of communication in the process of measurement while drilling.

Signal attenuation model and transfer characteristics for continuous vibration signals inside a pipeline

LIUJun1,2,YUANFeng2(1. School of Electrical Engineering & Automation, Harbin Institute of Technology, Harbin 150001, China;2. School of Electrical Engineering & Information, Northeast Petroleum University, Daqing Heilongjiang 163318, China)

Key words:amplitude-frequency characteristic; measurement while drilling; transfer coefficient; hydraulic impedance; pressure vibration

管道中水擊導致的壓力劇烈變化經常帶來劇烈振動和管道的損壞,但是在某些情況下,可控的水擊現象也可以用于數據的傳輸。在鉆井過程中,可以認為鉆柱是一根垂直或傾斜向下的中空管道,鉆井液通過鉆柱注入地下,并從鉆柱外面的環空返回地面,形成了鉆井液的循環過程。為了實時獲取井下信息,需要在地面和井下建立信息傳遞通道,如果在井下利用機械裝置,有目的產生壓力波動,并賦予波動一定的含義,在地面對波動信號進行測量并解碼,就可以將井下測量信息傳遞到地面,這種技術稱為無線(Measurement While Drilling,MWD),它是目前比較先進的井下測量技術[1-4]。目前國內外對于脈沖信號的傳遞與編碼技術已經比較成熟,研制并投入使用了多款以鉆井液壓力脈沖作為數據傳輸方式的儀器,但是使用脈沖信號傳遞數據時,由于信號基于時間編碼,導致數據傳遞速率較低[5-8]。為了進一步提高數據傳輸效率,國外一些鉆井設備公司(Schlumberger, Halliburton和Baker Hughes)已經開始研制基于連續波信號傳遞方式的MWD設備[3-4],連續波載波方式可以極大的提高數據傳輸速率。目前國內對于連續波信號傳遞的傳遞機理與編碼方式也開始了研究[9-18],從文獻[11]的摩阻模型和文獻[13]的傳輸特性模型可以得出管道中信號頻率越高衰減越嚴重的結論。但是波動信號在管道中傳遞時,其波動情況與很多因素有關,這些影響因素與連續波動信號的傳遞之間的關系還需要進一步研究,需要尋找更好的數學描述方法。本文從水擊理論出發建立了大傾角管道內連續波信號傳遞的矩陣計算模型,分析了平均摩擦力情況下管道內的壓力波分布和管道的幅頻特性,最后通過地面實驗,驗證了連續波動信號在管道中衰減特性,為MWD過程中信號頻率的選取提供了理論支持。

1連續波動信號的數學模型

目前對管道內單個壓力波動的數學描述是基于著名的水擊方程組,方程組由兩個方程構成,一個為運動方程,另一個為連續方程:

(1)

式中,V是流速,H是水力壓頭,β為管道與水平面的夾角,c是波速,D為管道水力直徑,t是時間,x是中間位置與波動起點的距離,f是與摩阻力相關的沿程阻力系數[6]。此方程組目前沒有解析解,在很多文獻中利用各種方法求出了它的數值解[6,8,12]。這些數值解說明了單個波動現象發生后管道中的波動情況,但是當波動信號持續產生后,只有數值解明顯是不夠的。

當管道中的壓力很大,但是流量Q(Q=AV,A為管道截面積)相對較小時,由于波速遠大于流速,所以有水頭在時間上的變化?H/?t遠大于其沿管路的變化V?H/?x,同時由于流速的變化V?V/?x也相對較小,因此如果忽略式(1)中V?H/?x和V?V/?x的影響,并將流速V用流量Q替代,可以將式(1)改寫為:

(2)

(3)

(4)

式中,R=fQ0/(gDA2),如果認為在波動傳遞過程中f不變,那么R代表了管道對波動的平均阻力,是流體穩定流動時管道的阻力。可以將(4)稱為管道內連續波動方程組。

2波動方程組的解

(5)

同時,從方程組(4)可以得到:

將它們代入式(5)可以得到:

(6)

從而得到壓頭波動在管道中關于時間和位置的微分方程。從式(6)可以看出,在管道中的任意位置存在著流量和水頭的連續波動。從式(6)還可以知道它的解為:

H′=(c1eγ1x +c2e-γ2x )ejωt

(7)

可以將式(7)寫成H′=h(x)ejωt的形式,h(x)代表了x處的壓頭波動振幅,式(7)描述的是管道中水頭的波動情況,式中c1和c2為積分常數,同時可以計算出γ1和γ2的表達式為:

可以看出當ω不變時,γ1和γ2是一個常復數,它們是管道材質與尺寸、流體特性以及管道傾角的集中體現,將γ1和γ2稱為特征常數,從式(7)可知,管道中的波動可以看作是向上游傳遞和向下游傳遞波動的疊加,它們的特征常數分別是γ1和γ2。特征常數的實部代表傳遞過程中幅值的衰減,虛部代表管道對相位的影響,gsinβ代表了重力對波動的影響,對于直井這種垂直向下的管道,有sinβ=1,對于斜井這種傾斜管道,則有0

(8)

將式(8)代入方程組(4)的第2式可以得到:

從而可以計算出Q′:

(9)

式(9)說明在管道中任意位置也存在著流量的波動。也可以將式(9)改寫為Q′=q(x)ejωt,q(x)可以看作是x位置處流量波動的振幅。

在管道中,由于是水頭的波動導致了流量的波動,同樣的,在電路中,電壓的波動也會導致電流的波動,根據二者的相似性,可以仿照電路中阻抗的定義來定義管道的阻抗為水頭與流量之比[19-20],定義管道阻抗為:

(10)

與輸電線路相似,也可以定義無限長管道的阻抗為特征阻抗,可以寫出管道的特征阻抗Zc的表達式:

(11)

從式(11)可以看出,特征阻抗Zc對于固定頻率的信號來說是不變的,代表了某個頻率下管道對波動信號的固有影響。

知道了c1,c2,可以寫出h(x)和q(x)的表達式為:

(12)

(13)

式(12)和式(13)式說明在管道中,只要知道了始端的壓頭和流量,就可以算出任意位置的壓頭和流量的振幅。

3連續壓力波動的傳遞矩陣

由于本文研究的是疊加在穩定水頭和流量上的變化量,根據流體力學知識,可以認為去除穩定不變量后的水頭振幅h與壓力波動振幅p是一致的,為了工程應用方便,用壓力振幅p來代替水頭振幅h。根據式(12)和式(13)可以定義出管道中連續壓力波的傳遞矩陣:

(14)

式中:

圖1 單一管道四端口等效模型 Fig.1 Four-port equivalent model of a single pipeline

如果不同材質與直徑的管道相連,則需要分別計算不同管道的傳遞矩陣、傳遞常數和特征阻抗。當一條管道和另一條管道連接時,相當于波動從一條管道的末端傳遞到另一條管道的始端,在接點處由于壓力和流量是連續的(qo=qi,po=pi),因此可以將連接點的傳遞矩陣表示為:

(15)

由于式(15)是單位矩陣,因此對于多條管道的串聯來說,整個管道的傳遞矩陣可以表示為:

M*=ML1ML2…MLn

(16)

第一條管道的輸入是整條管道的輸入,最后一條管道的輸出是整條管道的輸出,這些串聯管道之間的關系可以用圖2所示的多個串聯四端口模型表示。

圖2 多條管道串聯時的四端口等效模型 Fig.2 Four-port equivalent model of pipeline series

傳遞矩陣方法在分析機械振動傳遞過程中得到廣泛的應用,可以很方便的分析振動的各種特性[21-23]。當連續波動在串聯管道中傳遞時,根據圖2所示的模型和式(16),可以按照矩陣乘法將多條管道串聯時的傳遞矩陣簡化為等效傳遞矩陣M*,從而將傳遞過程表示為:

(17)

在MWD鉆井過程中,如果認為井下振動發生器是信號源,則通過式(17)就可以計算壓力波信號在鉆柱中的傳遞情況以及到達地面的振幅。當波動到達地面后,如果認為地面管道處于水平位置,可以取傳遞系數中的sinβ=0,說明波動傳遞不受重力的影響。

由于式(17)可以用于計算多條管道串聯后,輸入與輸出端的流量與壓力之間的關系,因此可以將(17)式稱為連續波動信號在管道中的傳遞模型。

4管道中連續波動信號幅值衰減分析

在MWD過程中,首先將需要傳遞的數據通過相位調制的方法調制到某一頻率正弦波動信號上,然后在地面測量管道中的壓力信號并對信號進行解調。在此過程中,最關心的是信號幅值的衰減,信號的衰減程度決定了數據傳遞的距離和信號處理的難度,如果定義發送端的壓力振幅為pi,接收端的壓力振幅為po,那么用于描述信號衰減程度的參數是po/pi,這個值越大說明收到的信號越強,信號傳遞效果越好。如果定義接收端的阻抗為Zl=po/qo,設管道長度為l,根據(14)可以推導出:

(18)

從式(18)可以看出,接收端壓力波動的幅值由pi和Zl決定。如果認為波動在傳遞中頻率不變,所受摩擦力不變,并且忽略流體中含氣量的變化,同時認為管道是剛性的,就可以運用前面推導的數學模型對連續波動信號在管道中的傳播過程進行理論分析與仿真。

4.1波動幅值沿管道的分布情況

根據(14)式和(18)式,如果知道波動過程中信號起始端的壓力和流量的幅值以及末端的終端阻抗,就可以計算管道中間位置x處壓力振幅與始端壓力振幅的關系:

(19)

根據前面A(x),B(x),C(x),D(x)的表達式,當x=l時,由于A(l)D(l)-B(l)C(l)=e(γ1-γ2)l ,可以看出(19)式與(18)式是一致的。假設管道直徑為80mm,長度為500m,動力粘度為20mPa·s,分別設定不同的終端阻抗,可以根據式(19)計算出沿管道的壓力波動幅值分布情況,當信號頻率為10Hz時其分布見圖3。在圖3中,橫坐標x/l是相對位置,x表示管道中間某個點到始端的距離,由于0≤x≤l,因此橫坐標的范圍為0~1,縱坐標是x/l位置處的壓力幅值px與始端x=0處的壓力幅值pi之比,在始端波動為周期性波動狀態下,管道中任意位置的壓力波動幅值是確定的,在不同位置波動幅值是不一樣的,幅值的變化呈駐波狀態。

圖3中三條曲線是三種不同終端阻抗下的壓力分布情況,當Zl=Zc時,管道中的壓力波動振幅呈現線性衰減,衰減斜率主要與管道中流體粘度和摩擦力相關。當在Zl>Zc時,整個管道中振幅波動較大,當Zl

圖3 10Hz時波動幅值沿管道的分布情況 Fig.3 The distribution of vibration amplitude at 10Hz alone the pipeline

從圖3還可以看出,在某些位置波動幅值接近0,說明在此位置波動非常小,這個現象說明,隨著井深的增加,有可能在某個深度收不到井下信號。

4.2管道中信號傳遞的幅頻分析

對于式(18),如果固定除頻率f以外的參數,讓頻率從0.1Hz變到50Hz,就可以計算出整個管道的50Hz以內的幅頻特性。

圖4 不同長度管道的幅值頻率特性曲線 Fig.4 The amplitude frequency curve in pipeline of different lengths

圖4為三種長度的管道在0~50Hz頻率振動情況下的幅頻特性圖,從圖中可以看出,在這三種長度的管道中,當振動頻率很低時,末端振幅與始端振幅比值接近1,但是隨著頻率的增加,不同長度管道對信號的衰減表現為不同的特性,在長度為50m的短管道中,末端幅值在整個頻率范圍內出現多個波峰,在7Hz,24Hz和42Hz達到峰值,在這三個頻率點附近,末端的振動幅值都大于始端的振動幅值,在長度為1000m的長管道中,波動衰減主要集中在低于10Hz的頻段,信號頻率超過10Hz時,末端信號振動幅值就已經非常微弱了。這種波動衰減現象也驗證了參考文獻[11,13]實驗中觀察到的波動現象和衰減現象。

從模型的幅值頻率特性可以看出,當輸入信號強度一定時,隨著信號頻率的增加,末端信號的強度呈現波動衰減現象,因此,可以將(17)式描述的數學模型稱為連續波動信號在管道中傳遞的衰減模型。

4.3流體運動粘度對信號傳遞的影響

流體振動信號在沿管道傳播過程中,可以看作是流體質點沿管道方向做往復振動,因此流體運動粘度對信號的傳遞也會有很大的影響。這個影響在宏觀上表現為流體內摩擦力,在水力學上可以用沿程阻力系數f來表示,f與運動粘度的關系可以由水力學中的達西公式計算,從而可以依據公式R=fQ0/(gDA2)計算出流阻。當固定管道中除粘度以外的參數時,可以繪制出不同粘度下管道的幅值頻率特性圖,圖5是長度為50m,管徑80mm的管道在三種粘度情況下的幅值頻率特性。

圖5 不同粘度流體在管道中的幅值頻率特性曲線 Fig.5 The amplitude frequency curve in pipeline at different viscosity

從圖5中可以看出,流體運動粘度會直接影響管道中波動幅值的大小,粘度越大,波動的峰值越小,而且粘度大的流體相對于粘度小的流體更早達到波動峰值。

5幅頻特性驗證實驗

為了驗證管道對連續波動的衰減特性,筆者利用地面實驗管道驗證了連續波動的幅頻特性,實驗裝置由恒壓水箱,管道,旋轉閥門和下游末端流量調節閥組成。管道上游為恒壓水箱,管道長度為67m,直徑27mm,管道下游安裝旋轉閥門作為波動發生器,旋轉閥門另一端接流量調節閥作為終端阻抗。在旋轉閥門附近安裝壓力傳感器,測量壓力記為pi,管道上游靠近水箱位置也安裝壓力傳感器,測量結果記為po,通過電機控制閥門轉動速度,就可以在管道中產生不同頻率的擾動波。

圖6 實驗數據與計算結果 Fig.6 Comparison of experimental data and calculated results

調整流量調節閥和水箱壓力,在管道中保持穩定流量,然后控制旋轉閥門在管道中產生從1Hz 到50Hz的擾動壓力波,波動頻率每次變化1Hz,分別記錄始端pi的波動Δpi和終端po的波動Δpo,就可以計算出衰減比值Δpo/Δpi,將50個測量結果繪制成曲線和計算結果相比較,比較結果見圖6,可以看出實驗數據的曲線和計算曲線是一致的,說明計算模型能夠很好的描述管道內連續波振動信號的傳遞與衰減。

在實驗過程中,隨著旋轉閥轉動頻率的增加,管道中流量會下降,為了保證管道流量的穩定,需要提高水箱的壓力來保證Δpi的穩定。

6結論

管道中的連續壓力波動在其沿管道傳遞過程中要受到多種因素的影響,這些因素主要包括管道的物理特性和尺寸,流體的粘度、壓力與流量,這些影響因素在本文研究的模型中都可以用傳遞系數和阻抗來描述,從而將復雜的問題簡單化,提供了分析波動信號的新手段。

本文推導的傳遞矩陣,不但適用于分析鉆柱這種大傾角管道,也適合分析地面普通管道中連續波動信號的傳遞。

傳遞矩陣分析方法可以將復雜的波動信號在管道中的傳遞過程用簡單的矩陣方式表達出來,適用于多種材質和尺寸的串聯管道中波動信號的分析,分析結果符合目前對傳遞現象的觀察與實驗。

如果要利用管道中的連續波動傳遞數據,則需要根據管道長度,流體特性選擇合適的載波頻率,否則很可能出現接收信號困難的情況。

參考文獻

[1]Montaron B A, Hache J M D, Voisin B. Improvements in MWD telemetry: The right data at the right time [R]. SPE,25356, 1993:337-346.

[2]石元會,劉志申,葛華,等. 國內隨鉆測量技術引進及現場應用[J].國外測井技術. 2009(01):9-13.

SHI Yuan-hui, LIU Zhi-shen,GE Hua, et al. Introduction and field application of Measurement While Drilling(MWD) technology inland[J]. World Well Logging Technology,2009(01):9-13.

[3]張紹槐,張潔. 21世紀中國鉆井技術發展與創新[J].石油學報.2001(06):63-68.

ZHANG Shao-huai, ZHANG Jie. The development and creation of drilling technology in China during the 21th century[J].Acta Petrolei Sinica. 2001(06):63-68.

[4]張春華,劉廣華. 隨鉆測量系統技術發展現狀及建議[J].鉆采工業,2010,33(1):31-35.

ZHANG Chun-hua, LIU Guang-hua. State of the art and development trend of MWD system[J]. Oil Drilling & Production Technology, 2010,33(1):31-35.

[5]蔡文軍,王平,祝遠征,等.機械式無線隨鉆測斜儀設計方案及關鍵技術[J].石油學報, 2006(27):103-106.

CAI Wen-jun, WANG Ping, ZHU Yuan-zheng, et al. Design scheme and key techniques for mechanical wireless inclinometer [J]. Acta Petrolei Sinica,2006(27):103-106.

[6]劉修善,蘇義腦.鉆井液脈沖信號的傳輸特性分析[J]. 石油鉆采工藝, 2000, 22(4) : 8-10.

LIU Xiu-shan, SU Yi-nao. Investigation on the transmission behaviors of drilling fluid pulse signal[J]. Oil Drilling & Production Technology, 2000,22(4):8-10.

[7]Klotz C, Kaniappan A, Thorsen A K. A new mud pulse telemetry systems reduce risks when drilling complex extended reach well[C]// SPE 2008, 15203.

[8]何樹山,劉修善.鉆井液正脈沖信號的衰減分析[J].鉆采工藝,2001, 24(6):1- 12.

HE Shu-shan,LIU Xiu-shan. Analysis of signal attenuation for positive drilling fluid pulse[J].Drilling & Production Technology, 2001, 24(6):1-12.

[9]沈躍,朱軍,蘇義腦,等.鉆井液壓力正交相移鍵控信號沿定向井筒的傳輸特性[J].石油學報,2011,32(3):340-345.

SHEN Yue, ZHU Jun, SU Yi-nao, et al.Transmission characteristics of the drilling fluid pressure quadrature phase shift keying signal along a directional wellbore[J].Acta Petrolei Sinica, 2011,32(3):340-345.

[10]Klotz C,Bond P,Wasserman I,et al.A new mud pulse telemetry system for enhanced MWD/LWD application[C]//2008, SPE 112683.

[11]王翔,王瑞和,紀國棟.井筒內鉆井液連續脈沖信號傳輸頻率相關摩阻模型[J].石油學報,2009,30(3):445-449.

WANG Xiang, WANG Rui-he, JI Guo-dong. Frequency dependent friction model for consecutive pulse signal of drilling fluid transmitting in borehole[J]. Acta Petrolei Sinica, 2009,30(3):445-449.

[12]石在虹, 劉修善. 井筒中鉆井信息的傳輸動態分析[J].天然氣工業, 2002,22(5):68-71.

SHI Zai-hong,LIU Xiu-shan. An analysis of drilling information transmission behavior in wellbore[J].Natural Gas Industry, 2002,22(5) :68-71.

[13]邊海龍,蘇義腦,李林,等.連續波隨鉆測量信號井下傳輸特性分析[J].儀器儀表學報,2011,32(5):983-988.

BIAN Hai-long, SU Yi-nao, LI Lin,et al. Downhole information transmission characteristic analysis of Measurement While Drilling(MWD) continuous wave signal[J].Chinese Journal of Scientific Instrument, 2011,32(5):983-988.

[14]沈躍,崔詩利,張令坦,等.鉆井液連續壓力波信號的延遲差動檢測及信號重構[J].石油學報,2013;34(2):353-358.

SHEN Yue,CUI Shi-li,ZHANG Ling-tan,et al. Delay differential detection and signal reconstruction of continuous pressure-wave signals of drilling fluid[J]. Acta Petrolei Sinica, 2013,34(2):353-358.

[15]劉瑞文,管志川,李春山.鉆柱振動信號的在線監測及應用[J].振動與沖擊,2013,32(1):60-68.

LIU Rui-wen,GUAN Zhi-chuan,LI Chun-shan.Drilling string vibration online monitoring and its application[J].Journal of Vibration and Shock,2013,32(1):60-68.

[16]Wang Jian-xun, Liu Hui-jin. A time-frequency mixed method for on-line monitoring of harmonics and interharmonics[C]//Proceedings of International Conference on Advanced Power System Automation and Protection,2011:28-233.

[17]Ronghui L, Erbin Y, Xiu Y. Analysis of transient harmonics in power systems based on wavelet packet transform[C]//Proceedings of 2010 The 3rd International Conference on Power Electronics and Intelligent Transportation System,2011,1:345-348.

[18]Jameson A. Time dependent calculations using multigrid, with applications to unsteady flow past airfoils and wings[J].AIAA Journal,1991,6: 1591-1596

[19]Boucher R F,Kitsiors E E.Simulation of fluid network dynamics by transmission line modeling, Proceedings of the institution of Mechanical Engineers[J]. Part C:Joumal of Mechanical Engineering Science,1986, 200:21-29.

[20]羅志昌.流體網絡理論[M].北京:機械工業出版社,1988:14-23.

[21]楊海根,芮筱亭,劉怡昕,等. 多體系統傳遞矩陣法分布式并行計算研究[J]. 振動工程學報,2014,27(1):9-15.

YANG Hai-gen,RUI Xiao-ting,LIU Yi-xin, et al. Study on distributed parallel computing of transfer matrix method for multibody systems[J].Journal of Vibration Engineering,2014,27(1):9-15.

[22]胡培民. 傳遞矩陣法在高頻振動分析中的應用[J].振動與沖擊,1996, 15(4):50-52.

HU Pei-ming .Analysis of high frequency vibration by transfer matrix method[J]. Journal of Vibration and Shock, 1996,15(4):50-52.

[23]李小燕,匡波,徐濟,等.網絡方法在管路流體動態仿真計算中的應用[J].核動力工程,2000,6(21):264-268.

LI Xiao-yan, KUANG,Bo, XU Ji,et al. Application of network method in calculation of dynamics of pipe-transmitted fluids simulation[J].Nuclear Power Engineering, 2000,6(21):264-268.

主站蜘蛛池模板: 一本大道东京热无码av| 毛片免费高清免费| AV无码一区二区三区四区| 亚洲不卡av中文在线| 亚洲综合久久成人AV| 国产精品视频观看裸模 | 国产又色又刺激高潮免费看| 中文字幕日韩欧美| 中文字幕亚洲综久久2021| 91久久青青草原精品国产| 又黄又湿又爽的视频| 成人国产一区二区三区| 亚洲天堂视频在线观看免费| 香蕉国产精品视频| 国产区福利小视频在线观看尤物| 亚洲国产系列| 亚洲成a人片| 午夜精品久久久久久久99热下载| 毛片久久久| 91在线精品麻豆欧美在线| 国产综合无码一区二区色蜜蜜| 国产第一页第二页| 毛片免费视频| 日韩在线视频网站| 自拍偷拍欧美| 波多野结衣在线se| a级毛片免费在线观看| 亚洲中文在线视频| 欧美色视频在线| 欧美专区在线观看| 亚洲午夜福利在线| 99九九成人免费视频精品| 精品国产一二三区| 中字无码av在线电影| 亚洲色图欧美一区| 久久免费成人| 国产精品美女自慰喷水| 少妇精品在线| 日本伊人色综合网| 国产杨幂丝袜av在线播放| 亚洲人人视频| 91青青草视频| 曰韩人妻一区二区三区| 91在线视频福利| 2021国产精品自产拍在线观看| 亚洲AⅤ波多系列中文字幕| 免费观看亚洲人成网站| 伊人色婷婷| 久久国产av麻豆| 国产色伊人| 久久a毛片| 欧美激情首页| 漂亮人妻被中出中文字幕久久| 在线免费亚洲无码视频| 69视频国产| 激情無極限的亚洲一区免费| 91欧美在线| 亚洲人成网站18禁动漫无码| 在线观看无码a∨| 国产精品高清国产三级囯产AV| 婷婷色一区二区三区| 婷婷色丁香综合激情| 真实国产乱子伦高清| 97在线观看视频免费| 色亚洲成人| 亚洲第一视频区| 日韩福利视频导航| 无码免费的亚洲视频| 日韩免费成人| 久久先锋资源| 国产乱人激情H在线观看| 精品一区二区久久久久网站| 欧美精品成人一区二区在线观看| 2021亚洲精品不卡a| 久久不卡国产精品无码| 国产精品三级专区| 久久人妻系列无码一区| 久久综合干| 国产日韩精品欧美一区灰| 国产国拍精品视频免费看| 亚洲无码高清免费视频亚洲| 国产精品无码制服丝袜|