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

基于微分算子的局部特征尺度分解方法及其應用

2015-10-29 06:45:48吳占濤程軍圣李寶慶鄭近德
中國機械工程 2015年17期
關鍵詞:模態信號方法

吳占濤 程軍圣 李寶慶 鄭近德

1. 湖南大學汽車車身先進設計制造國家重點實驗室,長沙,4100822.安徽工業大學,馬鞍山,243032

基于微分算子的局部特征尺度分解方法及其應用

吳占濤1程軍圣1李寶慶1鄭近德2

1. 湖南大學汽車車身先進設計制造國家重點實驗室,長沙,4100822.安徽工業大學,馬鞍山,243032

針對局部特征尺度分解(LCD)方法在目標信號所含頻率分量較接近時分解能力降低,易出現模態混淆現象,從而導致內稟尺度分量失去原有物理意義的問題,提出了基于微分算子的局部特征尺度分解(DOLCD)。DOLCD在對目標信號進行分解前,先將目標信號進行一階微分,則在頻率比一定的情況下,可提升分解能力及抑制模態混淆能力。研究了DOLCD方法的原理,通過仿真信號模型將DOLCD與LCD的分解能力進行對比分析,結果表明,DOLCD方法在提高分解能力,抑制模態混淆等方面具有一定的優越性,并將DOLCD方法應用于轉子不對中故障的診斷,結果表明該方法有效。

局部特征尺度分解;微分算子;故障診斷;內稟尺度分量;轉子

0 引言

經驗模態分解(empirical mode decomposition,EMD)是一種自適應的非平穩信號處理方法[1-2],EMD方法自提出后在很多領域得到了廣泛的應用[3]。然而,EMD方法存在包絡過沖和欠包絡、迭代計算量大、端點效應以及模態混淆等問題[4]。作為對EMD方法的改進,程軍圣等[5]提出了另一種自適應時頻分析方法——局部特征尺度分解(local characteristic-scale decomposition, LCD)。LCD能夠自適應地將一個復雜信號分解為若干個瞬時頻率具有物理意義的內稟尺度分量(intrinsic scale component,ISC)之和,從而得到目標信號完整的時頻分布。與EMD相比,LCD在減小擬合誤差、提高計算速度和抑制模態混淆等方面有一定優勢,已被應用于信號分析和機械故障診斷等領域,取得了較好的效果[6-7]。LCD方法理論基礎尚不完善,實際應用中難免出現模態混淆現象,即同一尺度或頻率的信號被分解到多個不同的ISC分量當中,或者一個ISC分量當中出現了尺度或頻率差異較大的信號,降低了LCD方法的分解能力,導致分解得到的ISC分量失去原有的物理意義。

為解決LCD的模態混淆問題,鄭近德等[8]基于向目標信號中成對地添加符號相反的白噪聲信號,對加噪信號依據頻率高低逐層進行集成平均分解的思路,提出了部分集成局部特征分解(partly ensemble local characteristic-scale decomposition,PELCD)方法。PELCD對LCD分解中的模態混淆有很好的抑制作用,但由于PELCD不能自適應地選擇添加白噪聲的幅值和數目,若參數選擇不合適,可能會產生偽分量。本文從另一種思路,解決LCD模態混淆問題。

改變目標信號中低頻分量與高頻分量的幅值比,可提升LCD方法的分解能力,進而有效抑制模態混淆的產生。借鑒劉海波等[9]在對EMD改進時提出的利用微分算子增強EMD分解能力的思路,在對目標信號進行LCD分解前,先將目標信號進行一階微分,則在頻率比一定的情況下,可降低低頻分量信號與高頻分量信號的幅值比,提升LCD方法的分解能力及抑制模態混淆能力?;诖?,筆者提出了基于微分算子的局部特征尺度分解(differential operator based local characteristic-scale decomposition, DOLCD)。

本文研究了DOLCD方法的原理,通過仿真信號將DOLCD與LCD的分解能力進行對比分析,結果表明,DOLCD方法分解能力范圍要明顯大于LCD方法的分解能力范圍,增強了抑制模態混淆能力。采用DOLCD方法對具有不對中故障的轉子振動位移信號進行了分析,結果表明DOLCD能夠有效地將高頻不對中成分與轉頻成分等進行分離,從而實現轉子故障診斷。

1 DOLCD方法

楊宇等[10]對LCD方法的分解能力從理論層面上進行了深入研究。研究結果表明LCD方法的分解能力主要由各分量信號的頻率比、幅值比決定,在頻率比一定的情況下,降低低頻分量信號與高頻分量信號幅值比,可提高LCD方法的分解能力,抑制模態混淆;而各分量信號的初始相位以及單個分量迭代次數對LCD方法的分解能力影響均不大。為提高LCD方法的分解能力,抑制模態混淆,本文提出了基于微分算子的局部特征尺度分解方法。對目標信號x(t)(t>0),DOLCD方法的主要分解步驟如下:

(3)對y1(t)進行積分,得到z1(t)。

(4)對z1(t)進行LCD分解,得到目標信號x(t)的第一個有效ISC分量I1(t)與剩余信號r1(t)。

(5)將分量I1(t)與目標信號x(t)分離,得到剩余信號s1(t),即

s1(t)=x(t)-I1(t)

(1)

再將s1(t)視為目標信號,重復上述步驟(1)~(5),直到sn(t)為單調或常函數,依次得到ISC分量I1(t),I2(t),…,In(t)和趨勢項r1(t),r2(t),…,rn(t)。令

R(t)=r1(t)+r2(t)+…+rn(t)

(2)

(6)原目標信號x(t)被分解為n個ISC分量和一個趨勢項R(t)之和,即

(3)

由于信號的初始相位對LCD的分解能力影響不大,為了討論方便、不失一般性,設輸入目標信號x(t)由兩種模態構成,即

x(t)=x1(t)+x2(t)=A1sin(2πf1t)+A2sin(2πf2t)

(4)

其中,f1>f2。低頻分量信號x2(t)與高頻分量信號x1(t)的頻率比F為

F=f2/f1

(5)

低頻分量信號x2(t)與高頻分量信號x1(t)的幅值比a為

a=A2/A1

(6)

對式(4)兩邊進行一階微分得

2πf2A2cos(2πf2t)

(7)

b=(f2/f1)(A2/A1)

(8)

由于f1>f2,則b

2 DOLCD分解能力

2.1分解能力評估模型與評價指標

為了評價DOLCD方法的分解能力,需確定一個分解能力的評價指標。參考文獻[10]定義式(9)為分解能力評價指標(evaluationindexofdecompositioncapacity),記為E:

(9)

式中,xi(t)為單分量信號的真實值;Ii(t)為分解得到的單分量信號;T為信號長度。

為了保證分解得到的分量全部準確,分解能力評價指標由下式確定:

E=max[Ei]

(10)

由于式(10)可能出現E>1的情況,為了便于研究,對E>1的情況,則令E=1,即E∈[0,1]。E值越大說明分解得到的Ii(t)越偏離真實分量xi(t),分解能力越弱;E值越小說明分解得到的分量Ii(t)越接近于真實分量xi(t),分解能力越強。文獻[10-11]研究表明,當E≤0.05時,可認為混合信號被準確分解。

2.2分解能力對比

為了不失一般性,仍采用式(4)所示的信號模型為分析對象,并設A1=4,f1=80,頻率比F∈[0.1,1],幅值比a∈[0.1,∞],信號采樣頻率fs=8192Hz。則在不同的幅值比a、頻率比F條件下得到DOLCD方法的分解能力評價指標E值如圖1所示。

圖1 幅值比a和頻率比F對E值的影響

由圖1可以看出,當a大致滿足a<1.1時,E值曲面形成一個坡面,該坡面在0.1≤F<0.7的部分較為平坦,此時E≤0.05;隨著F增大,坡面逐漸變陡,在F接近1的部分坡面又變得相對平坦,但E值接近1。當a>1.1,F≤0.37時,E值曲面也較平坦,此時E≤0.05。當a>1.1,F>0.37時,隨著F增大,會出現非常陡峭的曲面。

當E≤0.05時,表示DOLCD方法具有較好的分解效果,此時可認為兩個分量信號能夠被有效分離。當0.05

圖2 DOLCD方法分解能力

由圖2可以看出,E=1和E=0.05兩條等高線將整個a-F平面劃分為3個區域:Ⅰ區為E=1區域,表示DOLCD方法完全不能將兩個單分量分離;Ⅱ區為0.05

采用LCD方法對式(4)所示的信號模型進行分解,運行參數與圖2相同,其分解能力區域分布圖如圖3所示。參照圖2,繪制E值分別為0.05和1的兩條等高線,將a-F平面劃分為3個區域。與圖2對比可以發現,圖3中的區域Ⅲ明顯小于圖2中的區域Ⅲ。表明在相同的采樣頻率fs、最大頻率fmax和相同最大迭代次數等運行條件下,DOLCD方法分解能力范圍要明顯大于LCD方法的分解能力范圍,增強了抑制模態混淆的能力。

圖3 LCD方法分解能力

2.3分解能力的影響因素

需要說明的是,圖3所示的LCD分解能力圖與文獻[10]有所差異,這主要是因為信號模型的采樣頻率fs,以及fs與最大頻率fmax的關系對LCD分解效果的影響。圖2和圖3中,采樣信號模型的采樣頻率fs=8192Hz,高頻分量信號最大頻率fmax=80Hz。仍采用式(4)所示的信號模型,分別更改信號模型的采樣頻率fs和最大頻率fmax,分析DOLCD方法和LCD方法的分解能力。當fs=4096Hz,高頻分量信號最大頻率fmax=80Hz時,以及當fs=8192Hz,高頻分量信號最大頻率fmax=200 Hz時,DOLCD方法和LCD方法的分解能力分別如圖4a、圖4c和圖4b、圖4d所示。

對圖4分析可得,盡管圖4a和圖4c的區域Ⅲ與圖2中的區域Ⅲ,以及圖4b和圖4d的區域Ⅲ和圖3中的區域Ⅲ均有所不同,但圖4a和圖4c的區域Ⅲ均明顯大于圖4b和圖4d的區域Ⅲ。這也說明了DOLCD方法的分解能力要強于LCD方法。

3 仿真分析與應用實例

3.1仿真分析

為了說明DOLCD方法的優越性、不失一般性,仍采用式(4)所示的信號模型為分析對象,其中A1=2,A2=18,f1=16,f2=4,即x(t)=2sin(32πt)+18sin(8πt),t∈[0,1],其時域波形如圖5所示。

(a)fs=4096 Hz、fmax=80 Hz時,DOLCD方法的分解能力

(b)fs=4096 Hz、fmax=80 Hz時,LCD方法的分解能力

(c)fs=8192 Hz、fmax=200 Hz時,DOLCD方法的分解能力

(d)fs=8192 Hz、fmax=200 Hz時,LCD方法的分解能力圖4 兩組不同運行參數對應的兩種方法的分解能力

(a)仿真信號x1(t)

(b)仿真信號x2(t)

(c)仿真信號x(t)圖5 仿真信號x(t)及其各成分的時域波形

信號x(t)是一個大幅值的低頻信號與小幅值的頻率略高信號的疊加信號,其幅值比a=A2/A1=9,頻率比F=f2/f1=0.25。根據對圖3的分析,信號x(t)處于圖3的區域Ⅰ,LCD無法將信號x(t)的兩分量有效分離。對圖6所示信號x(t)的LCD分解結果分析可得,LCD分解得到第一個分量I1(t)與目標信號x(t)的幅值、頻率均區別不大;第二個分量I2(t)的幅值很小,從其時域圖上無法看出有物理意義的頻率成分。對其第一個分量I1(t)的幅值譜(圖7)進行分析,可明顯看出分量I1(t)包含了目標信號x(t)的兩種頻率成分??梢姴捎肔CD對信號x(t)進行分解時,第一個ISC分量I1(t)中包含有頻率差異較大的不同信號成分,出現了嚴重的模態混淆,LCD未能實現對仿真信號x(t)的有效分解,這與前面的分析是一致的。

(a)第一個ISC分量I1(t)

(b)第二個ISC分量I2(t)

(c)剩余項R2(t)圖6 仿真信號x(t)的LCD分解結果

圖7 LCD分解信號x(t)得到的分量I1(t)的幅值譜

對圖2進行分析,信號x(t)處于圖2中區域Ⅲ,采用DOLCD可有效實現信號x(t)的兩分量分離,信號x(t)的DOLCD的分解結果如圖8所示。

(a)第一個ISC分量I1(t)

(b)第二個ISC分量I2(t)

為了進一步比較兩種方法的分解效果,論文還對兩種方法得到的兩分量與對應真實分量的相關系數[11]進行了對比,相關系數指標值越大,表示分解的準確性越高。DOLCD得到的兩分量的相關系數指標值為0.9999和1,LCD得到的兩分量的相關系數指標值為0.1104和0.1713。這說明與LCD方法相比,DOLCD方法在分解精確度方面表現出一定的優越性。

3.2應用實例

為了進一步說明DOLCD方法的有效性與實用性,論文將其應用于不對中故障[12]的轉子振動位移實驗信號分析。實驗裝置如圖9所示,轉子轉速為3000 r/min,轉頻fr=50 Hz,實驗采樣頻率fs=2048 Hz,采樣時長為0.5 s。實驗數據的時域波形如圖10所示。其幅值譜如圖11所示。采用剛性聯軸器連接的轉子不對中時,強制連接所產生的力矩,使轉子發生彎曲變形;轉軸內阻以及轉軸表面與旋轉體內表面之間的摩擦而產生的相對滑動,使轉子產生自激旋轉周期性振動,其振動頻率為轉子轉頻fr的2倍,常伴頻率為轉頻fr的1倍頻,轉子振動位移信號會產生一個以轉頻為幅值調制頻率的高頻分量[13]。由圖10可以看出,高頻的分量信號被淹沒在強大的背景信號中。從圖11中主要看到的是與轉頻fr相關的基頻分量,看不出不對中故障信息。為提取高頻不對中故障信息,分別采用DOLCD和LCD對實驗數據進行分解,分解結果如圖12和圖13所示,限于篇幅只畫出前三個分量。

圖9 轉子實驗裝置示意圖

圖10 實驗數據的時域波形

圖11 實驗數據的幅值譜

(a)第一個ISC分量I1(t)

(b)第二個ISC分量I2(t)

(c)第三個ISC分量I3(t)圖12 實驗數據的DOLCD分解結果

(a)第一個ISC分量I1(t)

(b)第二個ISC分量I2(t)

(c)第三個ISC分量I3(t)圖13 實驗數據的LCD分解結果

由圖12可以看出,DOLCD方法對實驗數據分解得到的第一個分量I1(t)是高頻背景噪聲信號,第二個分量I2(t)具有明顯的調制特征,包含了主要故障信息,分量I3(t)是與轉頻有關的背景信號。從圖13中可以看出,LCD方法對實驗數據分解得到的第一個分量I1(t)也是高頻背景噪聲信號,第二個分量I2(t)也能看出調制特征,第三個分量I3(t)是與轉頻有關的背景信號;但分量I2(t)和I3(t)均出現了嚴重的波形失真和模態混淆(如圖中實線所圈部分)。采用希爾伯特變換(Hilbert transform, HT)分別估計兩種方法得到的第二個分量的包絡譜,結果分別如圖14和圖15所示。

圖14 DOLCD的分量I2(t)的包絡譜

圖15 LCD的分量I2(t)的包絡譜

由圖14所示的DOLCD分量I2(t)的包絡譜圖可以看到明顯的轉頻fr的2倍頻和轉頻fr的1倍頻,這與轉子發生不對中故障時的頻率特征相吻合[14]。圖15所示的LCD分量I2(t)的包絡譜圖中未出現明顯的轉頻fr的2倍頻,無法有效識別實驗數據所包含的故障類型。

對上述實驗數據分析說明,相對LCD方法,DOLCD能更為有效地將高頻不對中信號從強大的背景信號中提取出來,實現不對中故障信號、背景信號和噪聲信號的分離,從而實現轉子故障診斷,同時也有效抑制了模態混淆。

4 結論

(1)采用仿真信號模型將DOLCD與LCD分解能力進行對比分析,發現DOLCD方法分解能力范圍要明顯大于LCD方法的分解能力范圍,分解結果準確性更高,也有效地抑制了模態混淆。

(2)將DOLCD方法應用于轉子不對中故障信號的分析,分解得到了高頻調制特征分量和低頻分量,有效地實現了轉子不對中故障的診斷。

論文提出的DOLCD方法不需要事先選擇基函數,是一種自適應的信號分析方法,可有效地應用于機械設備故障診斷。本文僅是簡單的從信號模型兩組不同的采樣頻率和最大頻率對兩種方法分解能力的影響進行了對比分析。筆者后續將進一步對信號模型不同的運行參數對LCD分解能力的影響進行系統性的研究。

[1]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[J]. Proceedings of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences, 1998, 454(1971): 903-995.

[2]Huang N E, Wu Z.A Review on Hilbert Huang Transform: Method and Its Applications to Geophysical Studies[J]. Reviews of Geophysics, 2008, 47(2): 1-23.

[3]Lei Yaguo, Lin Jing, He Zhengjia, et al. A Review on Empirical Mode Decomposition in Fault Diagnosis of Rotating Machinery[J]. Mechanical Systems and Signal Processing, 2013, 35(1): 108-126.

[4]時培明, 蔣金水, 劉彬, 等. 基于邊界特征尺度匹配延拓的EMD改進方法及應用[J]. 中國機械工程, 2014,25(12): 1616-1623.

Shi Peiming, Jiang Jinshui, Liu Bin, et al. Improved Method of EMD and Its Applications Based on Boundary Characteristic Scale Matching Extension Method[J]. China Mechanical Engineering, 2014,25(12): 1616-1623.

[5]程軍圣, 鄭近德, 楊宇. 一種新的非平穩信號處理方法——局部特征尺度分解法[J]. 振動工程學報, 2012, 25(2): 215-220.

Cheng Junsheng, Zheng Jinde, Yang Yu. A Nonstationary Signal Analysis Approach——the Local Characteristic-scale Decomposition Method[J]. Journal of Vibration Engineering, 2012, 25(2): 215-220.

[6]Zheng Jinde, Cheng Junsheng, Yang Yu. Generalized Empirical Mode Decomposition and Its Applications to Rolling Element Bearing Fault Diagnosis[J]. Mechanical Systems and Signal Processing, 2013, 40(1): 136-153.

[7]鄭超, 郭奇, 郭麗杰, 等. 基于局部特征尺度分解的旋轉機械故障欠定盲源分離方法研究[J]. 燕山大學學報, 2014, 38(2): 168-174.

Zheng Chao, Guo Qi, Guo Lijie, et al. Underdetermined Blind Source Separation Method of Rotating Machinery Faults Based on Local Characteristic-scale Decomposition[J]. Journal of Yanshan University, 2014, 38(2):168-174.

[8]鄭近德, 程軍圣, 楊宇. 部分集成局部特征尺度分解:一種新的基于噪聲輔助數據分析方法[J]. 電子學報, 2013, 41(5): 1130-1135.

Zheng Jinde, Cheng Junsheng, Yang Yu. Partly Ensemble Local Characteristic-scale Decomposition: A New Noise Assisted Data Analysis Method[J]. Acta Electronica Sinica, 2013, 41(5): 1130-1135.

[9]劉海波, 玄志武. 利用微分算子增強EMD算法頻帶分解能力[J]. 振動與沖擊, 2013, 32(18): 133-135.

Liu Haibo, Xuan Zhiwu. Improving Frequency-band Separating Ability of EMD with a Differential Operator[J]. Journal of Vibration and Shock, 2013, 32(18): 133-135.

[10]楊宇, 曾鳴, 程軍圣. 局部特征尺度分解方法及其分解能力研究[J]. 振動工程學報, 2012, 25(5): 602-609.

Yang Yu, Zeng Ming, Cheng Junsheng. Research on Local Characteristic-scale Decomposition and Its Capacities[J]. Journal of Vibration Engineering, 2012, 25(5):602-609.

[11]Rilling G, Flandrin P. One or Two Frequencies? The Empirical Mode Decomposition Answers[J]. IEEE Transactions on Signal Processing, 2008, 56(1): 85-95.

[12]孫斌, 劉立遠, 雷偉. 基于正交局部保持映射的轉子故障診斷方法[J]. 中國機械工程, 2014, 25(16): 2219-2224.

Sun Bin, Liu Liyuan, Lei Wei. Rotor Fault Diagnosis Method on Orthogonal Locality Preserving Projection[J]. China Mechanical Engineering, 2014, 25(16): 2219-2224.

[13]盛兆順, 尹琦嶺. 設備狀態監測與故障診斷技術及應用[M]. 北京: 化學化工出版社, 2003.

[14]劉楊, 太興宇, 姚紅良, 等. 雙盤轉子軸承系統不對中-碰摩耦合故障分析[J]. 振動、測試與診斷, 2013, 33(5): 819-823.

Liu Yang, TaiXingyu, Yao Hongliang, et al. Study on Misalignment-rubbing Coupling Fault of Dual-disk Rotor-bearing System Including the Impact of Motor Coupling[J]. Journal of Vibration, Measurement & Diagnosis, 2013, 33(5): 819-823.

(編輯王旻玥)

Local Characteristic-scale Decomposition Method Based on Differential Operator and Its Applications

Wu Zhantao1Cheng Junsheng1Li Baoqing1Zheng Jinde2

1.State key Laboratory of Advanced Design and Manufacture for Vehicle Body, Hunan University,Changsha,410082 2.Anhui University of Technology,Ma’anshan,Anhui,243002

A novel non-stationary signal method differential operator based LCD(DOLCD) was proposed for improving the LCD method, of which, when the target signal contained components with frequencies close to each other, the decomposition ability was decreased, mode mixing phenomenon would easily occur and the physical meaning of intrinsic scale components decomposed were lost. A differential operator was used to differentiate the target signal in time domain, then the differentiated target signal was decomposed with DOLCD, then in the certain frequency ratio, the decomposition ability was enhanced, and the mode mixing was inhibited. The paper firstly studied the theory of DOLCD, then simulation experiments were used to compare the decomposition ability of DOLCD with LCD. The results indicate that DOLCD is more efficient in improving the decomposition ability, and inhibiting the mode mixing than LCD. Finally, the proposed method was applied to diagnose the rotor with misalignment faults successfully which indicated the effectiveness of DOLCD.

local characteristic-scale decomposition(LCD); differential operator; fault diagnosis; intrinsic scale component(ISC); rotor

2015-03-15

國家自然科學基金資助項目(51375152);湖南省科技計劃資助項目(2014WK3005)

TH165.3;TH911.7DOI:10.3969/j.issn.1004-132X.2015.17.003

吳占濤,男,1982年生。湖南大學機械與運載工程學院博士研究生。主要研究方向為動態信號處理及機械設備故障診斷。程軍圣(通信作者),男,1968年生。湖南大學機械與運載工程學院教授、博士、博士研究生導師。李寶慶,男,1984年生。湖南大學機械與運載工程學院博士研究生。鄭近德,男,1986年生。安徽工業大學機械工程學院講師、博士。

猜你喜歡
模態信號方法
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
基于FPGA的多功能信號發生器的設計
電子制作(2018年11期)2018-08-04 03:25:42
用對方法才能瘦
Coco薇(2016年2期)2016-03-22 02:42:52
基于LabVIEW的力加載信號采集與PID控制
國內多模態教學研究回顧與展望
四大方法 教你不再“坐以待病”!
Coco薇(2015年1期)2015-08-13 02:47:34
捕魚
基于HHT和Prony算法的電力系統低頻振蕩模態識別
由單個模態構造對稱簡支梁的抗彎剛度
計算物理(2014年2期)2014-03-11 17:01:39
主站蜘蛛池模板: 国产乱码精品一区二区三区中文| 国产精品香蕉| 乱系列中文字幕在线视频| 在线国产综合一区二区三区 | 色135综合网| 美女一级毛片无遮挡内谢| 国国产a国产片免费麻豆| 国产精品自在在线午夜| 露脸真实国语乱在线观看| 91在线丝袜| 国产高潮流白浆视频| 欧美另类第一页| 欧亚日韩Av| 午夜福利网址| 亚洲精品制服丝袜二区| 99精品国产自在现线观看| 天天综合网色| 久久永久免费人妻精品| 国产成人精品日本亚洲| 国产女人在线| 搞黄网站免费观看| www成人国产在线观看网站| 国产十八禁在线观看免费| 亚洲男人的天堂在线观看| 国产91av在线| 精品视频第一页| 亚洲无线国产观看| 亚洲色欲色欲www在线观看| 伊人五月丁香综合AⅤ| 亚洲人成网址| 欧美第一页在线| 91色国产在线| 亚洲欧洲日产国产无码AV| 99视频在线观看免费| 亚洲最黄视频| 91精品啪在线观看国产60岁 | www.youjizz.com久久| 无码电影在线观看| 99草精品视频| 国产中文在线亚洲精品官网| 天堂成人在线| 狠狠色成人综合首页| 人妻一区二区三区无码精品一区| 亚洲日本中文字幕乱码中文 | 日韩精品欧美国产在线| 亚洲人成人伊人成综合网无码| 国产一级毛片高清完整视频版| av在线人妻熟妇| 久久亚洲天堂| 欧美精品1区| 亚洲无码高清视频在线观看| 这里只有精品在线播放| 一级不卡毛片| 国产国产人在线成免费视频狼人色| 亚洲啪啪网| 欧美激情网址| 久爱午夜精品免费视频| 国产真实乱子伦视频播放| 嫩草在线视频| 国产免费黄| 国产va免费精品| 久久美女精品| 黄片在线永久| 乱系列中文字幕在线视频| 免费一级毛片不卡在线播放| 欧美成人午夜在线全部免费| 国产av色站网站| 日本a∨在线观看| 91视频国产高清| 无码中字出轨中文人妻中文中| 九九热这里只有国产精品| 免费看a毛片| 一级毛片免费播放视频| 亚洲国产综合精品一区| 女人av社区男人的天堂| 日韩中文精品亚洲第三区| 成人va亚洲va欧美天堂| 亚洲av无码专区久久蜜芽| 成人综合在线观看| 欧美高清国产| 亚洲一区二区成人| av天堂最新版在线|