亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        基于信息融合和廣義循環(huán)互相關(guān)熵的電機(jī)軸承故障診斷

        2022-02-16 01:29:24郝如江
        振動(dòng)與沖擊 2022年2期
        關(guān)鍵詞:振動(dòng)故障信號(hào)

        李 輝, 郝如江

        (1.天津職業(yè)技術(shù)師范大學(xué) 機(jī)械工程學(xué)院,天津 300222; 2.石家莊鐵道大學(xué) 機(jī)械工程學(xué)院,石家莊 050043)

        滾動(dòng)軸承廣泛應(yīng)用于各種類型電機(jī)中,其健康狀態(tài)對(duì)保證電機(jī)的正常運(yùn)轉(zhuǎn)具有重要作用[1]。在電機(jī)軸承故障檢測(cè)與診斷中,常常采用基于振動(dòng)信號(hào)處理的方法[2]。但在采集得到的振動(dòng)信號(hào)中,不僅包含大量的高斯和非高斯噪聲,而且還呈現(xiàn)一定的非線性[3]。為得到可靠的診斷效果,往往需要對(duì)傳感器采集的信號(hào)進(jìn)行降噪處理,如卡爾曼濾波、共振稀疏分解[4]、最大相關(guān)峭度解卷積[5]和經(jīng)驗(yàn)?zāi)B(tài)分解[6]等,而這些方法在處理非線性、非高斯噪聲信號(hào)時(shí),會(huì)造成性能衰退甚至失效[7]。基于信號(hào)高階統(tǒng)計(jì)量的方法[8-9],雖然對(duì)高斯噪聲具有免疫能力,但也難以有效抑制強(qiáng)脈沖噪聲,而且高階統(tǒng)計(jì)量計(jì)算方法復(fù)雜且計(jì)算量大,限制了高階統(tǒng)計(jì)量的實(shí)際應(yīng)用和推廣。近年來,在通信領(lǐng)域,基于信號(hào)分?jǐn)?shù)低階統(tǒng)計(jì)量[10-12]( fractional lower-order statistics, FLOS)和相關(guān)熵[13-15](Correntropy)的方法得到了廣泛應(yīng)用,在非高斯、非線性信號(hào)處理方面取得了較好的效果。相關(guān)熵是處理脈沖噪聲的有效方法,已在雷達(dá)和通信信號(hào)檢測(cè)、信號(hào)濾波、波達(dá)方向估計(jì)和時(shí)延估計(jì)等方面得到應(yīng)用和驗(yàn)證,取得了良好效果[16-18],其性能優(yōu)于FLOS。相關(guān)熵采用高斯核函數(shù),高斯核函數(shù)只有一個(gè)核長(zhǎng)調(diào)整參數(shù),只能調(diào)整高斯概率密度曲線的方差,而不能調(diào)整其“鐘形”形態(tài),因而自適應(yīng)性較差。為改善相關(guān)熵的性能,2016年,Chen等[19]以廣義高斯分布函數(shù)作為核函數(shù),提出了廣義相關(guān)熵的概念,廣義高斯分布函數(shù)具有形狀和核長(zhǎng)兩個(gè)參數(shù),因而具有更好的適應(yīng)性,已在通信信號(hào)自適應(yīng)濾波方面取得了良好效果[20-21],但還未在機(jī)電設(shè)備故障診斷領(lǐng)域得到應(yīng)用。本文將多傳感器信息融合[22]、廣義相關(guān)熵與循環(huán)平穩(wěn)信號(hào)處理方法相結(jié)合,提出了基于信息融合和廣義循環(huán)互相關(guān)熵的軸承故障診斷方法,并用仿真信號(hào)和電機(jī)軸承內(nèi)圈、外圈局部裂紋故障試驗(yàn)信號(hào)對(duì)提出方法進(jìn)行了驗(yàn)證,信號(hào)仿真和試驗(yàn)結(jié)果表明:本文提出的方法具有很強(qiáng)的降噪能力,為電機(jī)軸承故障檢測(cè)和診斷提供了一種有效方法。

        1 廣義循環(huán)互相關(guān)熵

        1.1 廣義高斯分布

        廣義高斯分布(generalized gaussian distribution, GGD)函數(shù)定義為[23]

        (1)

        式中,Γ(·)是Γ函數(shù),即

        (2)

        (3)

        式中,σ為方差。

        GGD共有兩個(gè)參數(shù):p是形狀指數(shù)(p>0),控制著概率密度函數(shù)的“形狀”,即衰減的速度;σ為方差(核長(zhǎng),σ>0)。

        1.2 廣義互相關(guān)熵

        對(duì)于任意兩個(gè)實(shí)隨機(jī)變量x和y(x,y∈R),它們的廣義互相關(guān)熵可定義為

        Vp,σ(x,y)=E[κp,σ(x,y)]

        (4)

        式中:E[·]為期望算子;κp,σ(·)為由式(1)定義的廣義高斯核函數(shù),即

        (5)

        式中,‖·‖為范數(shù)算子。

        (6)

        實(shí)隨機(jī)過程x(t)和y(t)(t∈R)的時(shí)變廣義互相關(guān)熵(cross correntropy)可表示為

        (7)

        式中,τ為時(shí)間滯后量。

        1.3 廣義循環(huán)互相關(guān)熵和譜密度

        (8)

        (9)

        在式(9)中,當(dāng)α=0時(shí),廣義循環(huán)互相關(guān)熵譜密度退化為傳統(tǒng)的功率譜密度。

        2 基于信息融合的廣義循環(huán)互相關(guān)熵軸承診斷方法

        利用兩個(gè)傳感器,采集兩路振動(dòng)信號(hào)x(t)和y(t),并對(duì)兩路信號(hào)進(jìn)行解卷積,以消除傳遞路徑的影響,然后計(jì)算兩路解卷積后振動(dòng)信號(hào)的廣義循環(huán)互相關(guān)熵譜密度,最后根據(jù)互相關(guān)熵譜識(shí)別軸承故障。其方法步驟如下:

        步驟1采集兩路振動(dòng)信號(hào)x和y;

        步驟2對(duì)振動(dòng)信號(hào)x和y進(jìn)行解卷積預(yù)處理;

        3 信號(hào)仿真

        滾動(dòng)軸承點(diǎn)蝕故障振動(dòng)響應(yīng)信號(hào)模型可以表示為[24-25]

        (10)

        式中:Ai為瞬態(tài)脈沖幅值;T為脈沖周期;τi為第i次沖擊相對(duì)于脈沖周期T的微小波動(dòng);n(t)為平穩(wěn)隨機(jī)噪聲。

        對(duì)于式(10)所示的滾動(dòng)軸承點(diǎn)蝕故障振動(dòng)響應(yīng)信號(hào)模型,瞬態(tài)沖擊s(t)是以系統(tǒng)的固有頻率fb為頻率的振蕩衰減信號(hào)。假定模擬軸承外圈故障振動(dòng)響應(yīng)信號(hào),此時(shí)Ai為常數(shù),式(10)可具體描述為

        (11)

        x(t)=x1(t)+n1(t)+n2(t)

        (12)

        式中:Ai為瞬態(tài)脈沖幅值;Ci為阻尼衰減因子;ti為沖擊持續(xù)的時(shí)間;θbi為初始相位;fbi為電動(dòng)機(jī)系統(tǒng)的共振頻率;n1(t)為零均值高斯噪聲;n2(t)為脈沖噪聲。

        函數(shù)Θ(t-ti)用來指定沖擊發(fā)生的時(shí)間,可用式(13)定義

        (13)

        設(shè)軸承外圈故障特征頻率為fo=110 Hz,電動(dòng)機(jī)系統(tǒng)的固有振動(dòng)頻率fb=1 000 Hz,瞬態(tài)沖擊振幅Ai=5,采樣頻率fs=6 000 Hz。

        信號(hào)仿真的目的主要有2個(gè)方面,首先通過式(11)表示的未加任何噪聲的軸承外圈仿真信號(hào),驗(yàn)證廣義循環(huán)互相關(guān)熵的頻譜特征,為軸承故障診斷奠定基礎(chǔ);其次通過式(12)的仿真信號(hào),驗(yàn)證廣義循環(huán)互相關(guān)熵的降噪性能。

        3.1 廣義循環(huán)相關(guān)熵頻譜特征仿真

        圖1(a)是式(11)表示的仿真信號(hào)的時(shí)域波形,圖1(b)是其FFT,從圖1(b)中可以看出,以系統(tǒng)的固有頻率1 000 Hz為中心,兩邊分布著以軸承外圈故障特征頻率fo=110 Hz為間隔的邊頻帶。

        圖1 軸承外圈故障仿真信號(hào)Fig.1 Simulative signal for bearing with outer race fault

        圖2、圖3分別為軸承外圈故障仿真信號(hào)廣義循環(huán)相關(guān)熵函數(shù)的二維和三維圖,從圖2和圖3中可以清楚地觀察到信號(hào)的循環(huán)頻率信息分布在循環(huán)頻率域的低、高兩個(gè)不同的頻率段,低頻段的循環(huán)頻率對(duì)應(yīng)軸承外圈故障頻率fo及其倍頻;高頻段的循環(huán)頻率是以二倍的載波頻率(α=2 000 Hz)為中心、以軸承外圈故障特征頻率fo為間隔的邊頻帶。因此,根據(jù)廣義循環(huán)相關(guān)熵函數(shù)可以準(zhǔn)確識(shí)別軸承故障。

        圖2 廣義循環(huán)相關(guān)熵(Contour圖)Fig.2 Generalized cyclic correntropy of simulative signal (Contour plot)

        圖3 廣義循環(huán)相關(guān)熵(三維圖)Fig.3 Generalized cyclic correntropy of simulative signal (3D plot)

        圖4和圖5分別為仿真信號(hào)廣義循環(huán)相關(guān)熵譜密度函數(shù)的二維和三維圖,在圖4中存在著許多獨(dú)立的譜峰,這些譜峰主要分布在f=0.5α和f=-0.5α兩條頻譜線及其平行線上;在頻譜中心(0,0)和2倍共振頻率處(±2 000,±1 000),這些平行線構(gòu)成菱形,菱形在水平方向的對(duì)角線長(zhǎng)度等于2fo,在垂直方向的對(duì)角線長(zhǎng)度等于fo,這種頻譜結(jié)構(gòu),清晰表達(dá)了軸承外圈故障頻譜特征,因此,根據(jù)這些頻譜特征可準(zhǔn)確識(shí)別軸承故障。

        圖4 廣義循環(huán)相關(guān)熵譜密度(Contour圖)Fig.4 Generalized cyclic correntropy spectral density of simulative signal (Contour plot)

        圖5 廣義循環(huán)相關(guān)熵譜密度(三維圖)Fig.5 Generalized cyclic correntropy spectral density of simulative signal (3D plot)

        3.2 廣義循環(huán)相關(guān)熵降噪性能仿真

        首先在軸承外圈故障仿真信號(hào)x1(t)中加入零均值高斯噪聲n1(t),信噪比為SNR=-3 dB,之后再隨機(jī)加入幾個(gè)幅值不同的脈沖信號(hào),以模擬脈沖噪聲。圖6(a)為式(12)表示的仿真信號(hào)時(shí)域波形,從圖6(a)可以看出,仿真信號(hào)完全被噪聲淹沒,因此從圖6(a)已完全看不出信號(hào)幅值的變化規(guī)律。圖6(b)為其FFT,在圖6(b)中也不能有效識(shí)別信號(hào)的頻率成分。

        圖6 軸承外圈故障仿真信號(hào)Fig.6 Simulative signal for bearing with outer race fault

        含噪聲仿真信號(hào)的廣義循環(huán)相關(guān)熵譜密度(p=3,σ=0.6),如圖7(為清晰看出頻譜結(jié)構(gòu),僅畫出了頻譜中心部分)和圖8所示。在圖7和圖8中可以清晰看出:在循環(huán)頻率σ和譜頻率f構(gòu)成的雙頻平面內(nèi),含噪聲仿真信號(hào)的能量主要分布在f=0.5α和f=-0.5α等兩條頻譜線上,譜峰幅值較大;由于噪聲的存在,在f=0.5α±1 000和f=-0.5α±1 000等四條譜線上也存在譜峰,但幅值較小,由這四條譜線形成較大的菱形,大菱形在水平方向的對(duì)角線長(zhǎng)度等于2fb,在垂直方向的對(duì)角線長(zhǎng)度等于fb,這種頻譜準(zhǔn)確表達(dá)了系統(tǒng)的共振頻率fb。通過上述分析表明:廣義循環(huán)相關(guān)熵譜具有很強(qiáng)的能量聚集性,能準(zhǔn)確表示強(qiáng)噪聲中信號(hào)的頻率成分,能有效抑制高斯噪聲和非高斯脈沖噪聲。

        圖7 廣義循環(huán)相關(guān)熵譜密度(Contour圖)Fig.7 Generalized cyclic correntropy spectral density of simulative signal (Contour plot)

        圖8 廣義循環(huán)相關(guān)熵譜密度(三維圖)Fig.8 Generalized cyclic correntropy spectral density of simulative signal (3D plot)

        為凸顯廣義循環(huán)相關(guān)熵的降噪能力,將廣義循環(huán)相關(guān)熵譜密度與傳統(tǒng)基于二階統(tǒng)計(jì)量的譜相關(guān)密度進(jìn)行對(duì)比,圖9和圖10給出了含噪聲仿真信號(hào)x(t)的譜相關(guān)密度圖[26],可以看出信號(hào)的能量分散在整個(gè)雙頻平面內(nèi),譜相關(guān)密度的能量聚集性很差,軸承外圈特征頻率和系統(tǒng)共振頻率已完全被噪聲掩蓋,難以有效識(shí)別。

        圖9 仿真信號(hào)譜相關(guān)密度(Contour圖)Fig.9 Spectral correlation density of simulative signal (Contour plot)

        圖10 仿真信號(hào)譜相關(guān)密度(三維圖)Fig.10 Spectral correlation density of simulative signal (3D plot)

        通過上述仿真信號(hào)可以得出以下結(jié)論:廣義循環(huán)相關(guān)熵具有很強(qiáng)的從強(qiáng)高斯噪聲和非高斯噪聲中提取信號(hào)特征的能力,基于廣義高斯核函數(shù)的廣義循環(huán)相關(guān)熵為高斯、非高斯噪聲的處理提供了一種嶄新的解決方法,其性能優(yōu)于傳統(tǒng)基于二階統(tǒng)計(jì)量的譜相關(guān)分析。

        4 電機(jī)軸承故障診斷

        采用美國(guó)凱斯西儲(chǔ)大學(xué)(Case Western Reserve University)軸承數(shù)據(jù)中心(Bearing Data Center)的網(wǎng)站上公布的試驗(yàn)數(shù)據(jù)[27]:試驗(yàn)軸承型號(hào)為:深溝球軸承6205-2RS JEM SKF,采樣頻率fs=12 000 Hz,電機(jī)負(fù)載為空載,電機(jī)轉(zhuǎn)速1 797 r/min(fr=29.95 Hz)。兩路信號(hào)為同一精度的兩個(gè)傳感器對(duì)同一目標(biāo)的不同部位進(jìn)行測(cè)量拾取的振動(dòng)信號(hào)。電機(jī)軸承內(nèi)圈故障振動(dòng)信號(hào)數(shù)據(jù)記錄號(hào)為105DE 和125DE,軸承外圈故障振動(dòng)信號(hào)數(shù)據(jù)記錄號(hào)為130DE和144DE,構(gòu)成兩路傳感器信號(hào),將兩路傳感器采集的振動(dòng)信號(hào)(以下分別用x和y表示)進(jìn)行信息融合,計(jì)算其廣義循環(huán)互相關(guān)熵譜。計(jì)算得到滾動(dòng)軸承內(nèi)圈、外圈故障特征頻率為[28]

        fi=162.185 Hz,fo=107.365 Hz

        4.1 軸承內(nèi)圈故障診斷

        圖11為滾動(dòng)軸承內(nèi)圈故障振動(dòng)信號(hào)x及其FFT。圖12為滾動(dòng)軸承內(nèi)圈故障振動(dòng)信號(hào)y及其FFT。在圖11(a)和圖12(a)中存在明顯的幅值調(diào)制現(xiàn)象,但根據(jù)時(shí)域波形及其FFT還不能識(shí)別軸承故障。

        圖11 軸承內(nèi)圈故障振動(dòng)信號(hào)x及其FFTFig.11 Vibration signal x and its FFT with inner race fault

        圖12 軸承內(nèi)圈故障振動(dòng)信號(hào)y及其FFTFig.12 Vibration signal y and its FFT with inner race fault

        為了驗(yàn)證廣義循環(huán)互相關(guān)熵在軸承故障診斷中的有效性,圖13給出了軸承內(nèi)圈故障振動(dòng)信號(hào)的廣義循環(huán)互相關(guān)熵譜的輪廓圖(p=4,σ=1,為清晰顯示頻譜結(jié)構(gòu),僅畫出了頻譜中心部分)。在圖13和圖14中可以清楚地看到:軸承內(nèi)圈故障特征fi頻率及其倍頻,主要沿頻譜直線f=±0.5α分布,逐漸由頻譜中心向外擴(kuò)展張成整個(gè)頻譜平面;以軸承內(nèi)圈故障特征頻率fi的譜峰構(gòu)成菱形,菱形對(duì)角線在水平方向上的長(zhǎng)度是2fi,在垂直方向上的長(zhǎng)度是fi,這種頻譜結(jié)構(gòu),清晰表達(dá)了軸承內(nèi)圈故障特征頻率fi。圖14給出了軸承內(nèi)圈故障振動(dòng)信號(hào)的廣義循環(huán)互相關(guān)熵譜的三維圖,從圖14中,能更清晰地看到軸承內(nèi)圈故障特征頻率fi的頻譜特征及其分布規(guī)律。

        圖13 軸承內(nèi)圈故障振動(dòng)信號(hào)廣義循環(huán)互相關(guān)熵譜Fig.13 Generalized cyclic cross correntropy spectral density of vibration signal with inner race fault (Contour plot)

        圖14 軸承內(nèi)圈故障振動(dòng)信號(hào)廣義循環(huán)互相關(guān)熵譜(三維圖)Fig.14 Generalized cyclic cross correntropy spectral density of vibration signal with inner race fault (3D plot)

        將廣義循環(huán)互相關(guān)熵譜密度與傳統(tǒng)基于二階統(tǒng)計(jì)量的譜相關(guān)密度進(jìn)行對(duì)比,圖15和圖16給出了軸承內(nèi)圈故障振動(dòng)信號(hào)的譜相關(guān)密度圖,可以看出傳統(tǒng)譜相關(guān)密度的能量聚集性很差,分辨率較低,不能有效識(shí)別軸承內(nèi)圈故障特征頻率fi。

        圖15 軸承內(nèi)圈故障振動(dòng)信號(hào)譜相關(guān)密度(Contour圖)Fig.15 Spectral correlation density of vibration signal with inner race fault (Contour plot)

        圖16 軸承內(nèi)圈故障振動(dòng)信號(hào)譜相關(guān)密度(三維圖)Fig.16 Spectral correlation density of vibration signal with inner race fault (3D plot)

        4.2 軸承外圈故障診斷

        圖17為滾動(dòng)軸承外圈故障振動(dòng)信號(hào)x及其FFT。圖18為滾動(dòng)軸承外圈故障振動(dòng)信號(hào)y及其FFT。圖19為其廣義循環(huán)互相關(guān)熵譜的輪廓圖、圖20為其廣義循環(huán)互相關(guān)熵譜的三維圖(p=4,σ=2)。從圖19和圖20可以清楚地看到:軸承外圈故障特征頻率fo及其倍頻,主要沿頻譜直線f=±0.5α分布,逐漸由頻譜中心向外擴(kuò)展張成整個(gè)頻譜平面;軸承外圈故障特征頻率fo及其倍頻的譜峰構(gòu)成菱形,菱形對(duì)角線在水平方向上的長(zhǎng)度是2fo,在垂直方向上的長(zhǎng)度是fo,這種頻譜結(jié)構(gòu),清晰表達(dá)了軸承外圈故障特征頻率fo的分布特征。

        圖17 軸承外圈故障振動(dòng)信號(hào)x及其FFTFig.17 Vibration signal x and its FFT with outer race fault

        圖18 軸承外圈故障振動(dòng)信號(hào)y及其FFTFig.18 Vibration signal y and its FFT with outer race fault

        圖19 軸承外圈故障振動(dòng)信號(hào)廣義循環(huán)互相關(guān)熵譜Fig.19 Generalized cyclic correntropy spectral density of vibration signal with outer race fault (Contour plot)

        圖20 軸承外圈故障振動(dòng)信號(hào)廣義循環(huán)互相關(guān)熵譜Fig.20 Generalized cyclic cross correntropy spectral density of vibration signal with outer race fault (3D plot)

        圖21和圖22給出了軸承外圈故障振動(dòng)信號(hào)的譜相關(guān)密度圖,從圖21和圖22可以看出傳統(tǒng)譜相關(guān)密度也不能有效識(shí)別軸承外圈故障特征頻率。

        圖21 軸承外圈故障振動(dòng)信號(hào)譜相關(guān)密度(Contour圖)Fig.21 Spectral correlation density of vibration signal with outer race fault (Contour plot)

        圖22 軸承外圈故障振動(dòng)信號(hào)譜相關(guān)密度(三維圖)Fig.22 Spectral correlation density of vibration signal with outer race fault (3D plot)

        通過上述分析可以看出:基于兩路傳感器振動(dòng)信號(hào)的廣義循環(huán)互相關(guān)熵具有信號(hào)解調(diào)功能,在由循環(huán)頻率α和譜頻率f構(gòu)成的雙頻平面內(nèi),能夠很清晰地刻畫軸承內(nèi)圈、外圈的故障特征,提高了軸承故障診斷的準(zhǔn)確性和可靠性,為從噪聲環(huán)境中提取軸承故障特征的有效方法,其性能優(yōu)于傳統(tǒng)的基于二階統(tǒng)計(jì)量的譜相關(guān)方法。

        5 結(jié) 論

        為提高電機(jī)軸承故障診斷的準(zhǔn)確性和可靠性,基于循環(huán)平穩(wěn)信號(hào)處理、相關(guān)熵理論和廣義高斯分布,提出了一種基于多傳感器信息融合的廣義循環(huán)互相關(guān)熵軸承故障診斷方法。先從廣義高斯分布入手,推導(dǎo)了廣義相關(guān)熵、廣義循環(huán)互相關(guān)熵和廣義循環(huán)互相關(guān)熵譜密度的計(jì)算公式;并利用軸承外圈故障仿真信號(hào)分析了軸承故障的廣義循環(huán)互相關(guān)熵譜的頻譜結(jié)構(gòu)及其特征,并利用該信號(hào)對(duì)廣義循環(huán)互相關(guān)熵譜的降噪性能進(jìn)行了驗(yàn)證,表明廣義循環(huán)互相關(guān)熵能夠有效處理高斯和非高斯噪聲,具有從強(qiáng)噪聲環(huán)境中提取信號(hào)頻譜特征的能力。最后將廣義循環(huán)互相關(guān)熵應(yīng)用于電機(jī)滾動(dòng)軸承內(nèi)圈、外圈故障診斷,試驗(yàn)結(jié)果表明:基于多傳感器信息融合的廣義循環(huán)互相關(guān)熵譜能有效提取軸承內(nèi)圈、外圈故障特征,其性能優(yōu)于傳統(tǒng)基于二階統(tǒng)計(jì)量的譜相關(guān)方法。

        猜你喜歡
        振動(dòng)故障信號(hào)
        振動(dòng)的思考
        信號(hào)
        鴨綠江(2021年35期)2021-04-19 12:24:18
        完形填空二則
        振動(dòng)與頻率
        故障一點(diǎn)通
        基于FPGA的多功能信號(hào)發(fā)生器的設(shè)計(jì)
        電子制作(2018年11期)2018-08-04 03:25:42
        中立型Emden-Fowler微分方程的振動(dòng)性
        奔馳R320車ABS、ESP故障燈異常點(diǎn)亮
        基于LabVIEW的力加載信號(hào)采集與PID控制
        故障一點(diǎn)通
        日本免费一区精品推荐| 亚洲精品午睡沙发系列| 乌克兰少妇xxxx做受6| 国产av无码专区亚洲aⅴ| 高清国产日韩欧美| 在线不卡中文字幕福利| 精品女厕偷拍视频一区二区| 久久久久夜夜夜精品国产| 99re久久精品国产| 色优网久久国产精品| 国产av无毛无遮挡网站| 久久久久亚洲av片无码| 欧美精品久久久久久久自慰| 国产男女乱婬真视频免费| 中文字幕精品一区二区的区别| 亚洲乱码日产精品一二三| 84pao强力打造免费视频34| 秀人网嫩模李梓熙大尺度| av免费不卡一区二区| 天天做天天爱天天综合网2021| 爱我久久国产精品| 精品一区二区亚洲一二三区| 国产自拍视频在线观看网站| 无码人妻一区二区三区在线视频| 激情人妻在线视频| 国产一区精品二区三区四区| 国产成人av一区二区三区在线观看 | 日本一区二区精品高清 | 欧美老妇人与禽交| 在线一区二区三区视频观看| 亚洲精品一区三区三区在线| 欧美 国产 综合 欧美 视频| jjzz日本护士| 成人在线观看视频免费播放| 欧美国产亚洲日韩在线二区| 天天操夜夜操| 91亚洲精品久久久蜜桃| 亚洲字幕中文综合久久| 无码精品a∨在线观看| 国产成年无码aⅴ片在线观看| 亚洲av日韩专区在线观看|