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

基于EEMD與分形盒維數(shù)的氣固兩相流流型識別

2016-10-10 02:43:41胡紅利
西北大學學報(自然科學版) 2016年2期
關(guān)鍵詞:信號

張 肖,李 林,胡紅利

(1.西安交通大學 電氣學院, 陜西 西安 710049;2.西安市產(chǎn)品質(zhì)量監(jiān)督檢驗院線纜室,陜西 西安 710065)

?

·數(shù)理科學·

基于EEMD與分形盒維數(shù)的氣固兩相流流型識別

張肖1,2,李林1,胡紅利1

(1.西安交通大學 電氣學院, 陜西 西安710049;2.西安市產(chǎn)品質(zhì)量監(jiān)督檢驗院線纜室,陜西 西安710065)

針對氣固兩相流氣力輸送管道中測控裝置后常見的幾種過渡流型,如層流、環(huán)流與核心流,文中通過靜電傳感器獲取靜電波動信號,利用總體經(jīng)驗模式分解和分形盒維數(shù)相結(jié)合的方法進行流型識別。首先使用總體經(jīng)驗模式分解對不同流型的靜電信號進行自適應分解從而獲得位于不同頻段的本征模式函數(shù),然后利用分形技術(shù)對所分解的IMF求取盒維數(shù),并將其作為特征量。最后通過提取的特征量對BP神經(jīng)網(wǎng)絡(luò)進行訓練和測試從而識別流型。實驗結(jié)果表明,該方法能夠有效識別氣固兩相流流型,識別率達94%。

氣固兩相流;總體經(jīng)驗模式分解;盒維數(shù);流型識別

氣固兩相流動普遍存在于工業(yè)生產(chǎn)過程中,如煤粉燃燒、火箭推力燃燒、醫(yī)藥、食品行業(yè)中的粉末處理過程及氣力輸送等工業(yè)領(lǐng)域。流型是氣固兩相流系統(tǒng)基本的特征參數(shù),它的變化會引起兩相流系統(tǒng)的流動特性、傳質(zhì)和傳熱的變化,直接影響兩相流系統(tǒng)的經(jīng)濟性、穩(wěn)定性和安全性。因此,對兩相流流型在線、實時和準確地識別具有很高的科學研究意義和應用價值[1-2]。由于氣固兩相流的流動特性十分復雜,使得流型識別成為工業(yè)現(xiàn)場測控的一個難題。目前針對氣固兩相流流型識別主要有兩類方法:一類是利用傳感器采集數(shù)據(jù),對管道截面成像,直接獲得管道內(nèi)流型及介質(zhì)分布,如ECT[3-4]和ERT[5-6]等方法;另一類是通過對傳感器獲取的流動信息進行信號分析獲得特征量,針對不同流型進行模型的訓練與識別,其信號分析方法有希爾伯特-黃變換[7]和小波包變換[8-9]等。

在氣固兩相流氣力輸送過程中,由于顆粒與顆粒、顆粒與管道之間的碰撞摩擦,顆粒與氣體之間的相對滑移以及靜電感應都可能使顆粒產(chǎn)生自然荷電現(xiàn)象,即運動的顆粒有靜電噪聲信號產(chǎn)生[10-11]。其靜電噪聲信號包含了大量的顆粒流動信息,本文通過靜電傳感器采集靜電噪聲信號,并對信號進行適當?shù)男畔⑻幚?即可獲得氣固兩相流流型信息。由于氣固兩相流流態(tài)形式多種多樣,十分復雜,使得靜電信號具有非線性非平穩(wěn)特性。總體經(jīng)驗模式分解(ensemble empirical mode decomposition, EEMD)[12]是Huang在經(jīng)驗模式分解[13]之后提出的又一種非線性非平穩(wěn)信號的時頻分析方法??傮w經(jīng)驗模式分解能夠利用噪聲在不同時間尺度上的連續(xù)性對不同時間尺度進行彌補,消除模態(tài)混疊。分形是一門以不規(guī)則事物為研究對象、探索復雜性的科學,所以它很自然地被用來描述設(shè)備振動信號的不規(guī)則性和復雜性[14]。分形理論是一種認識事物由粗到細的一個過程,這與EEMD對信號進行分解獲得位于不同特征時間尺度下本征模式函數(shù)的過程是一致的。由此啟發(fā)本文將EEMD與分形理論相結(jié)合,使用EEMD對信號進行處理后得到不同頻帶的本征模式函數(shù)(intrinsic mode function, IMF)信號,其盒維數(shù)的大小反映出信號的復雜度和非平穩(wěn)性,從而達到流型識別的目的。

1 總體經(jīng)驗模式分解

總體經(jīng)驗模式分解能夠根據(jù)信號的局部特征尺度有效的把信號分解為若干個IMF信號,它是一種噪聲輔助分析的方法。由于噪聲在時間上隨機分布,所以它在不同的時間尺度上具有連續(xù)性,可以彌補一些時間尺度的缺失。同時利用白噪聲的統(tǒng)計特性,對多次加噪聲處理后的分解結(jié)果求平均值能降低加入噪聲的影響。這樣避免了經(jīng)驗模式分解(empirical mode decomposition, EMD)方法中由于IMF的不連續(xù)性而造成的模式混淆現(xiàn)象。

EMD算法可以看作一種迭代過程,把原始數(shù)據(jù)看作r0,從j=1開始每循環(huán)一次j加1,其步驟如下。

1)找到信號rj-1的所有極大值和極小值,利用3次樣條插值獲得上包絡(luò)線和下包絡(luò)線v(t);

2)計算rj-1與上下包絡(luò)線均值之差,把它作為hj,1,

(1)

3)把hj,1作為步驟1),2)中的rj-1重復步驟1),2)得到hj,2,hj,3,…,直到hj,k滿足IMF的條件為止,從而得到IMFj,定義為cj,同時得到剩余信號為

rj=rj-1-cj,

(2)

rj形成了下次迭代的起始信號,從而形成了一個迭代的過程,當剩余信號rj為單調(diào)函數(shù)或者足夠小時,迭代停止,得到若干個IMF信號和一個剩余信號,如公式(3)所示,

(3)

EEMD是對信號多次加入噪聲,對每次加噪聲的信號進行經(jīng)驗模式分解,最終對多次分解結(jié)果求平均的一個過程,i代表加入白噪聲的次數(shù),其步驟如下:

1)對原始信號x(t)加入一組白噪聲分量ni(t),得到信號Xi(t),

Xi(t)=x(t)+ni(t);

(4)

2)利用EMD算法把加入噪聲后的信號分解為IMF信號,

(5)

3)重復步驟1)和2),每次加入的白噪聲具有相同幅值,每重復一次對i加1,直到i等于提前設(shè)定的值N;

4)對分解的IMF分量總體求平均作為最后的結(jié)果

(6)

式(6)中,cj(t)表示最終平均后得到的第j個IMF分量, N為加入白噪聲的次數(shù),加入白噪聲的次數(shù)與幅值滿足式(7),

(7)

其中,εn為加入噪聲引起的分解誤差,ε為加入噪聲的幅值,N為加入白噪聲的次數(shù)。Huang指出對于大部分非線性非平穩(wěn)信號,加入噪聲的幅值為信號標準差的0.2倍,如果信號主要包含高頻分量加入噪聲幅值可減小,如果主要包含低頻分量加入噪聲幅值可適當增大。

為了驗證EEMD在克服模態(tài)混疊現(xiàn)象上的有效性,對包含間斷信號的正弦信號分別使用EEMD和EMD進行分解,信號可以表示為

x(t)=sin(20πt)+n(t)。

(8)

式中,n(t)是幅值為0.2的間斷信號。結(jié)果如圖1和圖2所示。圖1中,由于間斷信號的存在,改變了原來的極值點分布,使得該間斷信號與部分低頻正弦波分解到imf1中,造成了明顯的模態(tài)混疊現(xiàn)象,導致分解結(jié)果不能反映出信號中固有的振動分量。圖2中,由于噪聲提供了各個頻段的振動背景,imf1,imf2和imf3中的振動模式分量主要由噪聲提供。其中imf1主要包含噪聲中的高頻振動。同時從imf2的幅值可以看出,間斷信號在該模式中被很好地提取出來,這是因為imf2提供的時間尺度與間斷信號相近。imf3是時間尺度位于imf2和imf4之間的一個振動分量,從幅值可以看出間斷信號對該分量影響已經(jīng)很小。在imf4中由于噪聲在該尺度上的分量相對于信號本身的能量很小,噪聲在該分量上的振動對原始正弦信號的影響可以忽略,所以原始正弦信號得到了很好的保留。

圖1 仿真信號的EMD分解結(jié)果Fig.1 Decomposition results of simulation signal using EMD algorithm

圖2 仿真信號的EEMD分解結(jié)果Fig.2 Decomposition results of simulation signal using EEMD algorithm

2 分形盒維數(shù)及其計算

2.1分形盒維數(shù)

若離散信號x(m)?U(m=1,2,…,N0)是n維歐式空間上的閉集。將Rn劃分為盡可能細的網(wǎng)格,若NΔ是網(wǎng)格寬度為Δ的離散空間上覆蓋集的最少網(wǎng)格個數(shù),則盒維數(shù)(boxdimension, db)定義為[14]

(9)

2.2盒維數(shù)的計算

由于離散信號x(m)的最高分辨率為采樣間隔Δt,即最小的網(wǎng)格應大于Δt,因此不能按照定義中的極限進行求解。實際中采用近似計算方法,在分形對象的無標度區(qū)內(nèi),將一系列尺寸為方形的網(wǎng)格對信號進行覆蓋,得到各尺度下的有效覆蓋的網(wǎng)格計數(shù)NkΔ,通過最小二乘法得到lgkΔ-lgNkΔ的擬合直線,其斜率就是該形體的分形盒維數(shù)D[15]。

令NkΔ為網(wǎng)格寬度為kΔ的離散空間上集合的網(wǎng)格計數(shù),可由公式(8)和公式(9)求得

l=1,2,…,N0/k, k=1,2,…,K,K

(10)

則信號x(m)的網(wǎng)格計數(shù)NkΔ為

NkΔ=P(kΔ)/kΔ。

(11)

然后在lgkΔ-lgNkΔ圖中確定線性好的一段為信號無標度區(qū)。如果無標度區(qū)的起點和終點分別為k1,k2,則在此區(qū)域有如下線性關(guān)系

lgNkΔ=-dblgkΔ+bk1≤k≤k2。

(12)

最后利用最小二乘法可求得信號的盒維數(shù)為

k1≤k≤k2。

(13)

一維離散信號的盒維數(shù)是介于1和2之間的一個分數(shù),信號越復雜維數(shù)越大。所以盒維數(shù)可以作為描述信號復雜性和不規(guī)則性的無量綱指標。

3 靜電信號采集系統(tǒng)

本文搭建了如圖3所示的吸送式氣固兩相流氣力輸粉試驗平臺[3]。該試驗平臺主要包括:載氣系統(tǒng)、給粉系統(tǒng)、流型發(fā)生器、測試段及收集系統(tǒng)。

圖3 氣力輸送裝置示意圖Fig.3 Schematic diagram of a pneumatic conveying system

ZF98工業(yè)吸塵器作為氣源安裝在氣力輸送系統(tǒng)的末端,能夠在管道內(nèi)部形成負壓。固相物料由孔徑可調(diào)的漏斗給粉(實驗所用的固相媒質(zhì)為經(jīng)由細篩子篩選,并用清水洗去塵土,烘干的Φ0.02~0.5mm細沙),然后物料和空氣在混合器中混合,形成連續(xù)均勻的氣固兩相流。在測試段之前安裝流型發(fā)生器,用來產(chǎn)生工業(yè)現(xiàn)場中常見的3種流型:層流、環(huán)流和核心流。測試段管道由有機玻璃制成,長度為0.6m,內(nèi)徑為36mm,外徑為40mm。距離流型發(fā)生器5cm處玻璃管外壁安裝靜電傳感器。流體流經(jīng)測試段后,將被送入旋風除塵器,在旋風除塵器中物料與空氣分離,被分離出來的物料進入旋風除塵器下方的回收器。實驗所采用的ZD98工業(yè)吸塵器,其容積為100L,空氣流量108L/s,由3個電機驅(qū)動,根據(jù)實驗所需要的動力選擇開啟電機的數(shù)量,調(diào)節(jié)合適的閥門開度,可以在負壓管道中形成不同的濃度以及速度。實驗中所使用的靜電傳感器的物理模型如圖4所示。

圖4 圓環(huán)型靜電傳感器結(jié)構(gòu)Fig.4 The construction of the ring shaped electrostatic sensor

靜電傳感器主要由圓環(huán)型電極、絕緣管道和金屬屏蔽層3部分組成靜電傳感器,通過法蘭將靜電傳感器的兩端與其余管道連接組成。為了提高靜電傳感器的抗干擾能力,將金屬管道接地以便起到屏蔽和抗干擾的作用。傳感器的結(jié)構(gòu)參數(shù)由表1給出。

表1靜電傳感器的結(jié)構(gòu)參數(shù)

Tab.1The parameters of the electrostatic sensor

參數(shù)名參數(shù)參數(shù)名參數(shù)參數(shù)名參數(shù)管道壁厚2mm電極軸向長度10mm屏蔽罩直徑60mm絕緣管道材料玻璃電極材料黃銅片屏蔽罩長度30mm管道外徑40mm電極厚度0.1mm屏蔽罩厚度0.1mm

4 實驗分析

利用靜電采集系統(tǒng)對層流、環(huán)流以及核心流以2 000Hz的采樣速率各采集100組數(shù)據(jù)。首先利用EEMD對實驗中采集的數(shù)據(jù)進行分解,加入白噪聲的次數(shù)為100次,幅值為信號標準差的0.2倍,圖5列舉出其中一組結(jié)果。從圖中可以看出隨著階次的增高,IMF中包含的時間尺度變大,復雜度降低。

EEMD具有自適應劃分頻帶以及二分濾波器[16]的作用,所以在3種流型下,同階IMF的頻譜位于相近的頻段范圍。不同階次IMF盒維數(shù)代表著該頻段范圍內(nèi)振動分量的復雜性和不規(guī)則性。如果3種靜電信號在不同頻段所包含的振動分量具有不同的復雜度,則能夠在各階IMF的盒維數(shù)上反映出來。列舉6組層流、環(huán)流和核心流靜電信號分解結(jié)果的盒維數(shù),得到圖6所示的結(jié)果。

圖5 層流、環(huán)流、核心流信號的分解結(jié)果Fig.5 Decomposition results of stratified flow, annular flow and core flow

圖6 3種流型靜電信號的IMF盒維數(shù)Fig.6 Fractal box dimensions of IMFsfor three flow regimes′ electrostatic signals

通過圖6可以看出,在列舉的6組對比中,層流比環(huán)流與核心流的盒維數(shù)大;核心流的盒維數(shù)與環(huán)流相近,但部分略大于環(huán)流的盒維數(shù)。對所有分解結(jié)果的盒維數(shù)求均值得到表2,對比均值可以看出,在總體趨勢上,層流各階IMF的盒維數(shù)比環(huán)流和核心流相應階次IMF的盒維數(shù)大,核心流各階IMF的盒維數(shù)比環(huán)流的略大。因此我們將各階IMF的盒維數(shù)作為神經(jīng)網(wǎng)絡(luò)的輸入進行流型識別。

表23種流型IMF盒維數(shù)的均值

Tab.2The averages of Fractal box dimensions of IMFs for three flow regimes

IMF階數(shù)層流環(huán)流核心流11.61101.46241.462321.51661.42071.444731.48451.42781.449341.43271.36861.379951.35171.28291.301961.27391.21471.232371.19251.13751.145181.11801.07501.082091.06771.03761.0440101.03451.03041.0291

對3種流型各100組數(shù)據(jù)進行處理,提取出300組IMF盒維數(shù)向量。300組數(shù)據(jù)中,3種流型各取70組用于訓練,其余的各30組用于測試。為了方便判斷,將層流編碼為[1,0,0],環(huán)流編碼為[0,1,0],核心流編碼為[0,0,1]作為目標輸出。神經(jīng)網(wǎng)絡(luò)訓練時的輸入為信號10階IMF的盒維數(shù)。實驗中BP神經(jīng)網(wǎng)絡(luò)的隱層神經(jīng)元個數(shù)為6個,隱層傳遞函數(shù)選為具有任意階數(shù)導數(shù)的非線性函數(shù)logsig,輸出層傳遞函數(shù)為線性函數(shù)purelin,訓練函數(shù)為traindx。訓練參數(shù)定義如下:最大訓練步數(shù)100,訓練目標誤差為0.000 01。構(gòu)建好的網(wǎng)絡(luò)如圖7所示。

圖7 神經(jīng)網(wǎng)絡(luò)模型Fig.7 Model of neural network

將90組測試樣本輸入到訓練好的神經(jīng)網(wǎng)絡(luò)模型,定義識別結(jié)果中大于0.5的位碼取1,小于0.5的位碼取0,把神經(jīng)網(wǎng)絡(luò)的輸出與實際流型的編碼對比,將結(jié)果繪制在表3中。

表3神經(jīng)網(wǎng)絡(luò)識別結(jié)果

Tab.3Identifications results from neural network

測量流型識別結(jié)果層流環(huán)流核心流層流2901環(huán)流0282核心流0327

根據(jù)表3中BP神經(jīng)網(wǎng)絡(luò)輸出結(jié)果可知,環(huán)流與核心流相對更容易識別錯誤,層流特征比較明顯,其總體識別率在94%,達到了很好的識別效果,驗證了該方法的有效性。

5 結(jié) 論

本文利用EEMD與分形盒維數(shù)相結(jié)合的特征提取方法, 對氣固兩相流進行流型識別。 EEMD能夠有效地對信號進行自適應分解而且抑制了模態(tài)混疊,分形盒維數(shù)反映了結(jié)果在各個頻段的復雜度。由于不同流型在某一相同頻段中可能包含不同的振動分量,所以盒維數(shù)能夠作為流型識別特征量。試驗結(jié)果證明該方法能有效地用于氣固兩相流流型識別,且識別率達94%。

[1]白博峰,王忠,郭烈錦,等.多相流流型識別技術(shù)與系統(tǒng)開發(fā)[J].西安交通大學學報(自然科學版),2003,37(3):306-309.

[2]HU H L, ZHANG J,DONG J, et al. Identification of gas-solid two-phase flow regimes using Hilbert-Huangtransformand neural-network techniques[J]. Instrumentation Science and Technology, 2011,39(2):198-210.

[3]周云龍,高云鵬,衣得武,等. ECT系統(tǒng)關(guān)鍵問題研究進展[J]. 化工自動化及儀表,2011(5):503-509.

[4]胡紅利,張娟,陳夏. 用于兩相流測量的ECT圖像重構(gòu)技術(shù)研究[J]. 工業(yè)儀表與自動化裝置,2010(2):100-103.

[5]李霞,冀海峰,黃志堯,等. 氣液兩相流ERT系統(tǒng)圖像邊緣檢測方法研究[J]. 儀器儀表學報,2005,26(增):55-56.

[6]董峰,姜之旭,喬旭彤, 等. 基于ERT技術(shù)的垂直管道兩相流流型識別[J]. 儀器儀表學報,2004,25(4):457-461.

[7]周云龍,王強,孫斌, 等. 基于希爾伯特——黃變換與Elman神經(jīng)網(wǎng)絡(luò)的氣液兩相流流型識別方法[J]. 中國電機工程學報,2007,27(11):50-56.

[8]孫斌,周云龍,張玲, 等. 基于小波包分解和Kohonen神經(jīng)網(wǎng)絡(luò)的氣液兩相流流型識別方法[J]. 熱能動力工程,2005,20(1):48-51+106.

[9]周云龍,張學清,李峰. 基于HMM和小波包分解的氣液兩相流流型識別[J]. 自動化儀表,2009(8):42-46.

[10] YAN Y. Guide to the flow measurement of particulate solids in pipelines, part 1: Fundamentals and principles[J]. Powder Handling & Processing, 2001, 13(4): 343-352.

[11] MEHRANI P, BI H T, GRACE J R. Electrostatic behavior of different fines added to a Faraday cup fluidized bed[J]. Journal of Electrostatics, 2007, 65(1): 1-10.

[12] WU Z, HUANG N E. Ensemble empirical mode decomposition: A noise-assisted data analysis method[J]. Advances in Adaptive Data Analysis,2009,1(1): 1-41.

[13] HUANG N E, SHEN Z, LONG S R, et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[C]//Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences. The Royal Society, 1998: 903-995.

[14] 訾艷陽,胥永剛,何正嘉. 離散振動信號分形盒維數(shù)的改進算法和應用[J]. 機械科學與技術(shù),2001,20(3):373-375.

[15] 郝研,王太勇,萬劍, 等. 分形盒維數(shù)抗噪研究及其在故障診斷中的應用[J]. 儀器儀表學報,2011,32(3):540-545.

[16] FLANDRIN P, RILLING G, GONCALVES P. Empirical mode decomposition as a filter bank [J]. Signal Processing Letters, IEEE, 2004, 11(2): 112-114.

(編輯李靜)

Identification of gas-solid two-phase flow regimes using EEMD and fractal box dimension

ZHANG Xiao1,2, LI Lin1, HU Hong-li1

(1.College of Electrical Engineering, Xi′an Jiaotong University, Xi′an 710049, China;2.Xi′an Supervision & Inspection Institute of Product Quality Cable Compartment, Xi′an 710065, China)

For three common transitional flow regimes behind the detection and control devices in the pneumatic conveying pipeline of gas-solid two-phase flow, namely stratified flow, annular flow and core flow, the electrostatic fluctuation signals were obtained through electrostatic sensor. In this paper, Ensemble Empirical Mode Decomposition (EEMD) and fractal box dimension are combined to identify the flow regimes. Firstly, via EEMD algorithm, the electrostatic signals of different flow regimes were decomposed into a series of Intrinsic Mode Function (IMF) in different frequency channels correspondingly. Then fractal box dimensions of IMFs were calculated and treated as the feature parameters. At last, the BP neural network was trained and tested using the feature parameters to identify the flow regimes. The experiment shows this method of combing EEMD and fractal box dimension can identify the flow regimes effectively, and the identification rate reaches 94%.

gas-solid two-phase flow; ensemble empirical mode decomposition; fractal box dimension; flow regime identification

2015-10-20

國家自然科學基金資助項目(51177120);國家“863”計劃基金資助項目(2009AA04Z130);電力設(shè)備電氣絕緣國家重點實驗室主任基金資助項目(EIPE14132)

張肖,女,陜西藍田人,從事氣固兩相流相關(guān)研究。

O359.2

A

10.16152/j.cnki.xdxbzr.2016-02-005

猜你喜歡
信號
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
7個信號,警惕寶寶要感冒
媽媽寶寶(2019年10期)2019-10-26 02:45:34
孩子停止長個的信號
《鐵道通信信號》訂閱單
基于FPGA的多功能信號發(fā)生器的設(shè)計
電子制作(2018年11期)2018-08-04 03:25:42
基于Arduino的聯(lián)鎖信號控制接口研究
《鐵道通信信號》訂閱單
基于LabVIEW的力加載信號采集與PID控制
Kisspeptin/GPR54信號通路促使性早熟形成的作用觀察
主站蜘蛛池模板: 国产精品九九视频| 久爱午夜精品免费视频| 伊在人亚洲香蕉精品播放 | 成人免费视频一区二区三区| 91丝袜乱伦| 99久久免费精品特色大片| 亚洲黄色激情网站| 最新日本中文字幕| 欧美午夜性视频| 欧美精品1区2区| 欧美精品H在线播放| www.亚洲天堂| 久久精品aⅴ无码中文字幕 | 性视频久久| 中文国产成人精品久久| 国产一区二区免费播放| 亚洲精品免费网站| 日本国产在线| 91亚瑟视频| 麻豆精品在线视频| 国产高潮视频在线观看| www欧美在线观看| 精品99在线观看| 伊人成人在线视频| 亚洲精选无码久久久| 国产精品久久精品| 欧美在线网| 一级毛片免费观看不卡视频| 91丝袜在线观看| 久久国产亚洲欧美日韩精品| 色婷婷综合激情视频免费看| 一级福利视频| 国产乱人免费视频| 国产99精品视频| 国产香蕉在线视频| 欧美性猛交xxxx乱大交极品| 91精品国产自产在线老师啪l| 狠狠躁天天躁夜夜躁婷婷| 国产丝袜一区二区三区视频免下载| 不卡视频国产| 亚洲中文无码av永久伊人| 一区二区三区高清视频国产女人| 日韩123欧美字幕| 青青草一区| 亚洲国产日韩欧美在线| 视频一区亚洲| 亚洲国产精品国自产拍A| 精品丝袜美腿国产一区| 欧美亚洲综合免费精品高清在线观看| 国产导航在线| 国产精品久久久久久久久| 国产麻豆91网在线看| 国产后式a一视频| 久久窝窝国产精品午夜看片| 国产网友愉拍精品| 久久综合九色综合97网| 午夜精品国产自在| 亚洲成A人V欧美综合| 欧美日韩激情在线| 狠狠干综合| 亚洲av无码专区久久蜜芽| 国产第八页| 97国产精品视频自在拍| 美女内射视频WWW网站午夜 | 美女视频黄又黄又免费高清| 欧美一级专区免费大片| 欧洲av毛片| 国产精品亚洲专区一区| 国产99精品久久| 激情无码视频在线看| 欧美一区精品| 播五月综合| 最新日韩AV网址在线观看| 免费在线播放毛片| 午夜a级毛片| 国产白浆一区二区三区视频在线| 婷婷午夜影院| 欧美日韩成人在线观看| 2021精品国产自在现线看| 国产丝袜91| 国产亚洲精品97在线观看| 国产激情无码一区二区三区免费|