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

基于EMD與多重分形去趨勢法的軸承智能診斷方法

2015-10-13 11:56:30賈峰武兵熊曉燕熊詩波
中南大學學報(自然科學版) 2015年2期
關鍵詞:趨勢振動故障

賈峰,武兵, 2,熊曉燕, 2,熊詩波, 2

?

基于EMD與多重分形去趨勢法的軸承智能診斷方法

賈峰1,武兵1, 2,熊曉燕1, 2,熊詩波1, 2

(1. 太原理工大學 機械電子工程研究所,山西 太原,030024;2. 太原理工大學 新型傳感器與智能控制教育部重點實驗室,山西 太原,030024)

引入經驗模態(tài)分解(EMD)方法去除故障信號的趨勢項,提出EMD-MFDFA(multifractal detrended fluctuation analysis)的多重分形分析方法,并通過仿真分析EMD方法去趨勢效果的有效性。然后將采用EMD-MFDFA方法提取的電機滾動軸承振動信號的多重分形特征向量作為訓練集,利用混合遺傳算法搜索全局最優(yōu)的能力優(yōu)化支持向量機(SVM)模型參數,建立電機滾動軸承的智能診斷模型。最后,通過對電機滾動軸承不同狀態(tài)的振動信號進行分析。研究結果表明:EMD-MFDFA方法能很好地揭示滾動軸承的振動信號多重分形特性,對滾動軸承正常狀態(tài)、單一故障與多故障耦合等狀態(tài)具有很強的辨識能力;所建立的智能診斷模型可以有效地診斷滾動軸承不同的故障狀態(tài),能夠作為滾動軸承故障在線監(jiān)測的有效工具。

多重分形;去趨勢波分析;經驗模態(tài)分解;遺傳算法;支持向量機

電機滾動軸承運行狀態(tài)往往直接影響整個系統(tǒng)的性能,因而研究其滾動軸承的故障狀態(tài)監(jiān)測診斷技術具有重要價值。電機滾動軸承的振動信號是包含著多種頻率成分的不規(guī)則信號,但是采集到的電機滾動軸承的振動信號的時域波形和變化趨勢會存在明顯的相似性即分形性質[1]。利用多重分形方法能精細地刻畫分形信號的局部結構,而且多重分形譜能夠很好地描述信號的結構特征和局部尺度行為[2],因而,可以將多重分形理論引入電機滾動軸承的故障診斷中。隨著多重分形理論的發(fā)展,Kantelhardt等[3]于2002年提出一種基于去趨勢波動分析(DFA)方法的多重分形分析方法——多重分形去趨勢法(MFDFA)。該方法利用DFA消除序列趨勢項的能力解決了標準多重分形分析處理非平穩(wěn)序列的欠缺問題,對非平穩(wěn)序列的處理更有效[4]。而且對平穩(wěn)時間序列,用去趨勢多重分形方法得到的結果與用標準方法得到的結果相同;而對于非平穩(wěn)的序列,其計算效果與其他方法相當[5]。因此,多重分形去趨勢法是一種更加實用的多重分形分析方法,已被應用于多個研究領域,取得了較好的分析效果。文獻[6]將MFDFA引入機械故障診斷領域,給出了DFA與MFDFA方法在分析機械故障信號上的優(yōu)缺點,表明MFDFA能夠更有效地揭示隱藏在多標度非平穩(wěn)信號中的動力學行為。文獻[7]利用MFDFA分析了滾動軸承的故障信號,利用實驗數據表明多重分形譜可以有效表征軸承故障信號的內在動力學機制,適合于軸承的故障特征提取。文獻[8]在利用MFDFA提取軸承振動信號的多重分形特征之后,使用支持向量機(SVM)建立診斷模型實現了故障診斷,但所建立的SVM模型只是單一特征的診斷模型,通用性較差。通常MFDFA方法假設信號的局部趨勢是多項式的形式,在消除序列趨勢項都是使用多項式擬合的方法,然而,機械振動信號的局部波動強烈且其局部趨勢都不是多項式的形式,所以,使用多項式進行擬合,很難準確地擬合信號的趨勢項。而經驗模態(tài)分解(EMD)是一種自適應的信號分析方法,依據數據自身的時間尺度特征來進行信號分解,能使復雜信號分解為有限個本征模函數(IMF)及1個余項,所分解出來的各IMF分量包含了原信號不同時間尺度的局部特征信號[9],余項即為信號的趨勢項。EMD方法理論上可以分解任何類型的信號,相對于多項式擬合,具有非常明顯的優(yōu)勢。而且使用MFDFA方法提取的信號多重分形特征有明確的物理意義,可以結合機器學習算法建立智能故障診斷模型,實現對軸承故障的在線診斷。支持向量機(support vector machine, SVM)作為近些年興起的機器學習方法,與其他方法相比,具有結構簡單、學習速度快、全局最優(yōu)、泛化性好等優(yōu) 點[10]。本文作者針對多項式不能合理擬合機械故障振動信號的趨勢項這一問題,引入EMD方法分解故障信號的趨勢項,提出了EMD-MFDFA的多重分形分析方法。然后利用EMD-MFDFA方法提取的電機滾動軸承多重分形特征,并結合混合遺傳算法優(yōu)化的SVM模型,建立電機滾動軸承的智能診斷模型。最后通過軸承的仿真信號及實驗信號,證明所提出的EMD-MFDFA方法及智能診斷模型的有效性。

1 理論基礎

1.1 MFDFA算法

對于時間序列x(=1,2,…,),設x序列是緊支集,MFDFA過程如下[4]。

1) 構造信號對于均值的累積離差()為

2) 將信號()以尺度劃分為互不重疊的N段數據:

3) 利用最小二乘法擬合每段數據,設y()為第段數據的階擬合多項式,為擬合多項式的階數,擬合的階數反映了數據去趨勢的效果,階數越高,去趨勢的效果越好。然后,求取原段數據與y()的均方誤差2(,):

式中:=1,2,…,;=1,2,…,N

4) 定義時間序列x的第階波動函數為

當=2時,MFDFA就退化為DFA。式(4)的物理意義為不同的描述不同程度波動對F()的影響。

5) 當時間序列x存在自相似特征,且時,波動函數F()與尺度之間存在穩(wěn)定的冪律關系為

式中:()為廣義Hurst指數,其值可在波動函數F()與尺度的對數關系圖上估算出來。

在描述多重分形的奇異譜時,()與基于配分函數方法中的標度指數()的關系為

再由統(tǒng)計物理中的legendre變換可以得到奇異指數和多重分形譜():

1.2 EMD方法

EMD方法是由美國NASA的黃鍔博士提出的一種信號分析方法,它依據數據自身的時間尺度特征來進行信號分解,無需預先設定任何基函數,在處理非平穩(wěn)及非線性數據上,具有非常明顯的優(yōu)勢[9]。因而將其引入MFDFA中,去除信號積累離差的趨勢項。其具體方法如下。

1) 找出信號()的局部極值點。

2) 把所有的局部極大值用三次樣條光滑地連接起來得到上包絡線+(),同樣地可以得到下包絡線?()。則可得到平均包絡線:

3) 求出差值函數z():

檢查z()是否滿足IMF條件,若不是,則把z()當作新的待處理量,重復做以上步驟,若z()滿足條件,則z()就是第1個IMF,記作1()。

4) 最后信號()被分解為個基本模式分量l()和一個余項r(),即

余項r()即為信號的趨勢項。

1.3 基于EMD的MFDFA算法

由于機械振動信號是由三角函數的疊加與相互調制得到的,且具有非平穩(wěn)的特性。使用多項式擬合這樣的信號,很容易會導致數據過擬合與欠擬合的現象發(fā)生,而EMD方法在理論上可以應用于任何類型的信號的分解,在處理非平穩(wěn)及非線性數據上,具有非常明顯的優(yōu)勢。所以,本文改進了MFDFA算法即EMD-MFDFA算法:在MFDFA算法第(3)步使用EMD方法去除信號累積離差量()的趨勢。

r()是第段累積離差的EMD分解后所得的余項,則均方誤差2(,)為

1.4 信號仿真

使用仿真的軸承故障振動信號進行仿真實驗。仿真信號的表達式為

式中:=900 Hz為固有頻率,=0.1為阻尼系數,rand[?1,1]表示在[?1,1]內的隨機數。采樣頻率為 20 kHz,采樣后數據1的數量為5 120點,如圖1所示。

圖1 軸承故障仿真信號

根據式(1)可得信號1的累積離差1,如圖2所示。以圖2可以看到:累積離差1局部波動很大,若使用多項式擬合1,則很容易出現過擬合及欠擬合的現象。

圖2 仿真信號的累積離差

對1進行MFDFA及EMD-MFDFA分析,由式(5)可知,在同一階波動下,lg2[F()]與lg2的關系為線性關系,即符合同一標度不變性原理。圖3所示為F()與的對數圖上,在各階波動下,lg2[F()]與擬合直線的殘差的均方根值(RMS)。從圖3可以看出:使用EMD方法去趨勢后,F()∝s(q)的冪律關系更穩(wěn)定,從而優(yōu)化了所提取的多重分形的譜特征。

1—MFDFA;2—EMD-MFDFA

2 EMD-MFDFA與SVM結合的智能診斷方法

2.1 EMD-MFDFA的特征提取方法

對信號進行EMD-MFDFA分析后,奇異指數反映了分形序列在某局部概率測試分布上的不均勻程度,奇異譜()是奇異指數的分維函數。由與()可得到描述多重分析譜的一些重要參數,這些參數能精細地表達多重分析時間序列動力學行為[6, 11]:1) 多重分形譜寬,Δ表示了信號多重分形性的強弱,可以反映振動信號波動的劇烈程度;2) 最大最小概率子集維數差別,Δ反映了概率測度所構成的子集中最大、最小概率處元素數量的比例,利用Δ可以以度量振動信號中波動平穩(wěn)與波動劇烈峰值數目的比例;3) 最小奇異指數min,min越小說明信號局部奇異性越強,即信號局部變化越劇烈,因而min在一定程度上也反映了振動信號的局部變化劇烈程度;4) 多重分形譜max()對應的奇異值0描述了信號的不規(guī)則性。

根據以上分析,在得到信號多重分形譜后,定義如下特征向量以提取多重分形特征:

2.2 建立智能診斷模型

支持向量機是建立在統(tǒng)計學習理論基礎上的一種新算法。自SVM出現以來,在故障診斷方面的應用越來越廣泛。SVM首先通過非線性變換將原始的樣本空間映射到高維的特征空間,然后在這個新空間中求取最優(yōu)線性分類面,進而對數據進行分類。而這種非線性變換通過定義適當的內積函數實現[10]。

在得到信號多重分形特征向量后,將其作為SVM的輸入,使用高斯徑向基作為SVM的核函數,從而實現軸承故障的智能診斷。由于懲罰因子與高斯核寬度是建立SVM智能診斷模型的重要參數,其取值決定著模型的診斷精度,本文使用遺傳算法對參數進行優(yōu)化[12]。盡管標準遺傳算法有很多優(yōu)點,但算法易陷入局部最優(yōu),而混合遺傳算法可以跳出局部最優(yōu)解,盡可能找到全局最優(yōu)解[13]。因而本文引入混合遺傳算法來確定與的最優(yōu)值,所使用的混合遺傳算法的流程圖見圖4。

圖4 混合遺傳算法流程圖

由上述分析,根據EMD-MFDFA與SVM原理,可將計算所得的作為SVM的輸入,建立智能診斷模型,診斷電機滾動軸承的故障。基于EMD-MFDFA的智能診斷模型的流程圖如圖5所示。

圖5 基于EMD-MFDFA的智能診斷模型的流程圖

3 實驗結果與分析

滾動軸承的故障振動信號是在Y160M2-8型電機軸承故障實驗臺上采集的。實驗是模仿實際工況進行的,電機軸承內圈與電機軸過盈配合,與軸無相對運動;外圈與電機端蓋配合,起支撐作用。由于軸承在電機內部,不能直接測量,但軸承振動直接傳遞到端蓋,這里將傳感器安裝在電機端蓋上,實現間接測量。滾動軸承型號為6309,加載器電流為0.5 A,采樣頻率為10 kHz。實驗數據包含正常、滾動體故障、內圈故障、滾動體與內圈的組合故障這4種故障狀態(tài)。由于垂直地面方向的振動信號包含信息最豐富,因而進行分析時選擇了傳感器的軸方向所采集到的信號,采集時間為21 s。圖6所示為電機軸承故障實驗臺示意圖,圖中傳感器為壓電式加速度傳感器,信號采集裝置為DEWETRON數據采集儀(型號為DEWE-2010)。所采集的4種故障振動信號如圖7所示。

圖6 電機軸承故障模擬實驗臺示意圖

(a) 正常狀態(tài);(b) 滾動體故障;(c) 內圈故障;(d) 滾動體與內圈的組合故障

選擇每類故障狀態(tài)振動信號的前5 120個數據進行EMD-MFDFA分析,參數最小值為64,最大值為1 024,取值間隔為48;參數最小值為?5,最大值5,取值間隔為0.1[14]。EMD方法的停機準則見參考文獻[15]。分析結果見圖8和圖9。

1—組合故障;2—滾動體故障;3—正常狀態(tài);4—內圈故障

1—組合故障;2—滾動體故障;3—正常狀態(tài);4—內圈故障

圖8所示為4種軸承振動信號的與廣義Hurst指數關系。從圖8可以看到:4種信號的廣義Hurst指數都是關于的曲線[14],說明這4種軸承狀態(tài)的振動信號都存在著多重分形的特性,從而導致振動信號出現多標度分形的特性。圖9所示為4種信號的多重分形奇異譜。從圖9可以看到:4種信號奇異譜具有不同的形狀、位置和譜寬,正是因為4種信號具有不同的多重分形特性,導致其內在的動力學機制不同,產生了譜差異。根據奇異譜的特征,可以提取信號的多重分形特性,從而進行有效的故障診斷。

由以上分析,將每種振動信號分為40段,每段5 120個數據點,分別使用EMD-MFDFA方法提取特征。然后與SVM原理結合,將所提取的多重分形特征向量作為SVM的輸入,進行電機軸承的智能診斷。

由于4種振動信號,每種信號有40個樣本,可以提取160組多重分形特征向量。將每種信號15組特征向量采用不同的方法建立不同SVM診斷模型,如表1所示。表1中:1表示每種信號有1組特征向量作為訓練集,2表示取每種信號有2組特征向量作為測試集。

表1 SVM模型信息

核函數選擇高斯徑向基核函數,然后對不同的SVM模型訓練。圖10所示為SVM-1模型使用混合遺傳算法尋找最佳參數時適應度函數值的進化曲線。

1—SVM-1;2—SVM-2;3—SVM-3

3個SVM的模型參數和的結果如表2所示。使用訓練好的SVM模型,對測試集分類,結果如表3所示。從表3可見:SVM-1可以有效區(qū)分正常軸承、單一故障軸承與多故障軸承,證明了EMD-MFDFA方法的有效性。從SVM-1與SVM-2的對比可見:使用EMD-MFDFA提取特征方法所建立的智能診斷模型的診斷精度更高,也說明去趨勢時使用EMD方法比使多項式法更適合于機械振動信號。從SVM-1與SVM-3的對比可見:在使用EMD-MFDFA方法建立智能診斷模型時,混合遺傳算法相對于標準遺傳算法能跳出局部極值點,收斂到全局最優(yōu),使SVM-1模型更適宜于在線智能診斷。

表2 SVM的參數值

表3 測試集的具體分類情況

注:/表示訓練樣本共有個,其中個訓練樣本分類正確。

4 結論

1) 利用EMD方法消除滾動軸承時間序列趨勢項,相比于多項式去趨勢更符合滾動軸承振動信號非平穩(wěn)的特性,去趨勢結果更準確。

2) 使用EMD-MFDFA方法提取滾動軸承不同故障狀態(tài)下振動信號的多重分形譜特征,與MFDFA方法相比,更能揭示滾動軸承多重分形特性,揭示不同狀態(tài)信號不同的內在動力學機制。該方法可以有效區(qū)分滾動軸承正常狀態(tài)、單一故障狀態(tài)與多故障狀態(tài),適合作為滾動軸承的特性提取方法。

3) 多重分形特征向量與SVM結合建立的智能診斷模型,可以有效診斷滾動軸承不同的故障狀態(tài),使用混合遺傳算法優(yōu)化的模型可以提高模型的診斷精度,且該智能診斷模型計算方便,可以作為滾動軸承故障狀態(tài)在線監(jiān)測的工具。

[1] YANG Junyan, ZHANG Youyun, ZHU Yongsheng, et al. Intelligent fault diagnosis of rolling element bearing based on SVMs and fractal dimension[J]. Mechanical Systems and Signal Processing, 2007, 21(5): 2012?2024.

[2] Lopes R, Betrouni N. Fractal and multifractal analysis: A review[J]. Medical Image Analysis, 2009, 13(4): 634?649.

[3] Kantelhardt J W, Zschiegner S A, Koscielny-Bunde E, et al. Multifractal detrended fluctuation analysis of nonstationary time series[J]. Physica A: Statistical Mechanics and its Applications, 2002, 316(1/4): 87?114.

[4] Talebinejad M, Chan A D C, Miri A. Fatigue estimation using a novel multi-fractal detrended fluctuation analysis-based approach[J]. Journal of Electromyography and Kinesiology, 2010, 20(3): 433?439.

[5] Norouzzadeh P, Rahmani B. A multifractal detrended fluctuation description of Iranian rial-US dollar exchange rate[J]. Physica A: Statistical Mechanics and its Applications, 2006, 367: 328?336.

[6] 林近山, 陳前. 基于多重分形去趨勢波動分析的齒輪箱故障特征提取方法[J]. 振動與沖擊, 2013, 32(2): 97?101.
LIN Jinshan, CHEN Qian. Fault feature extraction of gearboxes based on multifractal detrended fluctuation analysis[J]. Journal of Vibration and Shock, 2013, 32(2): 97?101.

[7] LIN Jinshan, CHEN Qian. Fault diagnosis of rolling bearings based on multifractal detrended fluctuation analysis and Mahalanobis distance criterion[J]. Mechanical Systems and Signal Processing, 2013, 38(2): 515?533.

[8] 李兆飛, 柴毅, 李華鋒. 多重分形去趨勢波動分析的振動信號故障診斷[J]. 華中科技大學學報(自然科學版), 2012, 40(12): 5?9, 17.
LI Zhaofei, CHAI Yi, LI Huafeng. Diagnosing faults in vibration signals by multifractal detrended fluctuation analysis[J]. Journal of Huazhong University of Science and Technology (Nature Science), 2012, 40(12): 5?9, 17.

[9] 楊世錫, 胡勁松, 吳昭同, 等. 旋轉機械振動信號基于EMD的希爾伯特變換和小波變換時頻分析比較[J]. 中國電機工程學報, 2003, 23(6): 102?107.
YANG Shixi, HU Jingsong, WU Zhaotong, et al. The comparison of vibration signals’ time-frequency analysis between EMD-based HT and WT method in rotating machinery[J]. Proceedings of the CSEE, 2003, 23(6): 102?107.

[10] Sapankevych N I, Sankar R. Time series prediction using support vector machines: A survey[J]. IEEE Computational Intelligence Magazine, 2009, 4(2): 25?38.

[11] 趙玲. 旋轉機械系統(tǒng)故障特征提取中的分形方法研究[D]. 重慶: 重慶大學機械工程學院, 2010.
ZHAO Ling. Research on fractal method of rotary mechanical fault feature extraction[D]. Chongqing: Chongqing University. College of Mechanical Engineering, 2010.

[12] Thomas M. Machine learning[M]. New York: McGraw Hill Higher Education, 1997: 249?273.

[13] WU Chihhung, KEN Yun, HUANG Tao. Patent classification system using a new hybrid genetic algorithm support vector machine[J]. Applied Soft Computing, 2010, 10(4): 1164?1177.

[14] Ihlen E A F. Introduction to multifractal detrended fluctuation analysis in Matlab[J]. Frontiers in Physiology, 2012, 141(3): 1?18.

[15] Rilling G, Flandrin P, Goncalvés P. On empirical mode decomposition and its algorithms[C]//Proc of IEEE EURASIP Workshop on Nonlinear Signal and Image Processing. [S.1.]: IEEE Press, 2003: 9?11.

Intelligent diagnosis of bearing based on EMD and multifractal detrended fluctuation analysis

JIA Feng1, WU Bing1, 2, XIONG Xiaoyan1, 2, XIONG Shibo1, 2

(1. Research Institute of Mechatronics Engineering, Taiyuan University of Technology, Taiyuan 030024, China;2. Key Laboratory of Advanced Transducers and Intelligent Control Systems of Ministry of Education,Taiyuan University of Technology, Taiyuan 030024, China)

The empirical mode decomposition (EMD) method was introduced to remove the fault signal local trends. A modified multifractal analysis method, i.e EMD-MFDFA (multifractal detrended fluctuation analysis), was proposed and the performance of the EMD-MFDFA method was proved by simulation. Then taking the multifractal feature vectors of motor bearing vibration signals extracted by the EMD-MFDFA method as the training sets, an intelligent diagnosis model was established to diagnose motor bearing early faults according to the support vector machine (SVM) theory, in which the parameters of SVM were optimized by using the hybrid genetic algorithm. Finally, the different states of the motor bearing vibration signals were analyzed. The results show that the EMD-MFDFA method can reveal the multifractal characteristics of the motor bearing signals and identify the bearing normal state, a single fault state and multiple faults state accurately. And the established model is able to diagnose bearing different fault effectively and suitable for on-line monitoring for motor bearing faults.

multifractal; detrended fluctuation analysis; empirical mode decomposition; genetic algorithm; support vector machine

TH165.3 TH17

A

1672?7207(2015)02?0491?07

2014?03?10;

2014?06?22

國家自然科學基金資助項目(51035007);山西省自然科學基金資助項目(2012011046-10)(Project (51035007) supported by the National Natural Science Foundation of China; Project (2012011046-10) supported by the Natural Science Foundation of Shanxi Province, China)

武兵,博士,副教授,從事機電系統(tǒng)故障診斷、動態(tài)測試研究;E-mail:wubing@tyut.edu.cn

10.11817/j.issn.1672-7207.2015.02.017

(編輯 趙俊)

猜你喜歡
趨勢振動故障
振動的思考
科學大眾(2023年17期)2023-10-26 07:39:14
趨勢
第一財經(2021年6期)2021-06-10 13:19:08
振動與頻率
天天愛科學(2020年6期)2020-09-10 07:22:44
故障一點通
中立型Emden-Fowler微分方程的振動性
初秋唇妝趨勢
Coco薇(2017年9期)2017-09-07 21:23:49
奔馳R320車ABS、ESP故障燈異常點亮
SPINEXPO?2017春夏流行趨勢
故障一點通
江淮車故障3例
主站蜘蛛池模板: 亚洲高清日韩heyzo| 国模沟沟一区二区三区| 国产另类乱子伦精品免费女| 国产成人一区二区| 91最新精品视频发布页| 久久精品嫩草研究院| 亚洲开心婷婷中文字幕| 午夜福利网址| 91外围女在线观看| 欧美不卡视频在线观看| 99re经典视频在线| 秘书高跟黑色丝袜国产91在线 | 久久精品娱乐亚洲领先| 亚洲精品无码不卡在线播放| 色窝窝免费一区二区三区| 亚洲欧洲国产成人综合不卡| 美女国内精品自产拍在线播放| 无码久看视频| 国产精品第三页在线看| 亚洲国产清纯| 国产大片黄在线观看| 日韩天堂在线观看| 国产丝袜第一页| 亚洲黄网在线| 国产69精品久久| 国产喷水视频| 亚洲精品波多野结衣| 国产在线精彩视频论坛| 好吊色妇女免费视频免费| 2048国产精品原创综合在线| 国产精品久久精品| 热久久这里是精品6免费观看| 亚洲婷婷在线视频| 国产主播在线一区| 狠狠亚洲五月天| 国产精品青青| 亚洲美女一区二区三区| 色婷婷视频在线| 亚洲va欧美ⅴa国产va影院| 亚洲成人在线网| 日本三级欧美三级| 国产人成乱码视频免费观看| 国产免费久久精品44| 欧美精品不卡| 国内99精品激情视频精品| 伊人中文网| 亚洲国产AV无码综合原创| 国产又色又爽又黄| 国产高清在线丝袜精品一区| 国产一区二区精品福利| 亚洲动漫h| 国产美女91视频| 亚洲国产清纯| 日韩欧美国产精品| 国产精品无码制服丝袜| 午夜国产精品视频| 精品夜恋影院亚洲欧洲| 国产小视频a在线观看| 青青草国产免费国产| 亚洲视频无码| 国产成熟女人性满足视频| 不卡午夜视频| 成人综合网址| 国产人成在线观看| 成人在线观看不卡| 激情無極限的亚洲一区免费| 丰满人妻久久中文字幕| 青青青伊人色综合久久| 久久久久国产一区二区| 99久久精彩视频| 色天天综合久久久久综合片| 就去吻亚洲精品国产欧美| 国产成人亚洲精品色欲AV| 亚洲午夜久久久精品电影院| 青青青国产精品国产精品美女| 亚洲精品老司机| 亚洲天堂区| 99久视频| 四虎永久免费在线| 伊人久久大香线蕉成人综合网| 欧美日韩激情在线| 亚洲天堂视频网站|