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

基于人工化學反應優化的SVM及旋轉機械故障診斷

2015-10-28 09:58:22羅頌榮程軍圣Hunglinh
中國機械工程 2015年10期
關鍵詞:故障診斷分類優化

羅頌榮 程軍圣 Hunglinh A O

1.湖南文理學院,常德,415003  2.湖南大學,長沙,410082

基于人工化學反應優化的SVM及旋轉機械故障診斷

羅頌榮1,2程軍圣2Hunglinh A O2

1.湖南文理學院,常德,4150032.湖南大學,長沙,410082

針對支持向量機(SVM)的參數優化問題,結合人工化學反應優化算法的優點,提出了基于人工化學反應優化算法的支持向量機(ACROA_SVM)方法;然后利用標準數據驗證了ACROA_SVM方法的有效性和優越性;最后,結合局部均值分解信號分析和能量矩特征提取,將ACROA_SVM方法應用于旋轉機械故障診斷中。分析結果表明,ACROA_SVM方法不但具有較高的故障診斷精度和較好的泛化能力,而且時間消耗短,故障診斷效率高,有利于實現在線智能故障診斷。

支持向量機;人工化學反應優化算法;旋轉機械;故障診斷

0 引言

機械故障診斷實質上是模式識別的問題。目前使用最廣泛的模式識別方法有聚類分析法、人工神經網絡(artificial neural network, ANN)法和支持向量機(support vector machine, SVM)法[1-3]。聚類分析法缺乏通用性,同時計算量較大。ANN法具有較強的自組織、自學習能力和非線性模式分類能力,但需要大量典型故障樣本,而在機械故障診斷的工程實踐中,往往典型的故障樣本稀缺;同時神經網絡具有過學習,結構和類型的選擇過分依賴于先驗知識等局限性,這些局限都會嚴重影響識別精度[4]。SVM是Vapnik[5]基于統計學習理論提出的機器學習方法,它通過合適的核函數將非線性分類空間的樣本映射到高維線性分類空間中,特別適合處理非線性、高維度、小樣本模式識別問題,較好地克服了人工神經網絡的缺陷,但SVM的核函數及其參數會嚴重影響SVM分類器性能。針對旋轉機械故障診斷中SVM核參數的優化問題,國內外學者研究了不少基于優化方法[6-8]的SVM模型,并應用于機械故障診斷中,取得了良好的效果[9-11]。但這些優化方法計算復雜且費時,不利于機械故障在線實時診斷。

人工化學反應優化算法(artificial chemical reaction optimization algorithm, ACROA)是2011年由Alatas Bilal提出的一種新的自適應的優化算法。該算法具有魯棒性好、輸入參數少、計算量小、計算時間短等優點[12-13]。本文針對故障診斷的實時性問題,結合人工化學反應優化算法的優點,提出了一種基于人工化學反應優化算法的支持向量機(support vector machine based on artificial chemical reaction optimization algorithm, ACROA_SVM)。隨后,利用標準數據,將ACROA_SVM和遺傳算法支持向量機(GA_SVM)、粒子群優化支持向量機(PSO_SVM)進行了性能對比分析,分析結果表明,ACROA_SVM在計算時間、分類精度上均優于GA_SVM方法和PSO_SVM方法。最后作為應用實例,將提出的AROCA-SVM方法應用于旋轉機械系統的典型元件——軸承的故障診斷中,結果表明,與GA_SVM、PSO_SVM相比,ACROA_SVM方法的分類準確度更高,耗費時間更短,更能滿足機械故障在線診斷的實時性要求。

1 基于化學反應的優化算法

基于化學反應的優化算法通過模擬一個容器內化學反應過程來解決全局和局部搜索優化問題。假設在固定容積的容器中有N種不同種類的化學反應物Ri(i=1,2,…,N),這些化學反應物可以通過M種方式進行化學反應。將Ri以一定的方式進行編碼,不同的編碼對應不同的反應規則,反應規則代表著化學反應物之間的相互作用。通常的編碼方式有實數編碼、二進制編碼、字符串編碼等。首先,設定初始反應物,然后,基于反應物的濃度和勢能來選擇反應物,并通過化學反應消耗或者產生反應物,并根據目標函數來自適應地選擇和更新反應物種群。當不再有新的反應物產生時,化學系統的熵達到最大時,或者達到最大迭代次數,滿足算法終止條件時,運算終止。AROCA可概括成以下步驟:

(1)優化問題和算法參數的初始化。一個優化問題可以概括為

minf(x)

s.t.xj∈[xjl,xju]

其中,f(x)為目標函數,變量x=(x1,x2,…,xj,…,xN),xjl為xj的下限,xju為xj的上限。

(2)設定初始反應物和計算焓值或熵值。在這一階段,采用統一種群方法將初始反應物均勻地初始化到可行搜索區間。首先,設定兩種反應物R`0=(u1,u2,…,un),R`1=(l1,l2,…,ln),n為反應物的長度,此時,分割因子k=1。然后,將分割因子不斷增加1,若分割因子k=2,則產生額外的兩種反應物R`2和R`3:

R`2=(r*u1,r*u2,…,r*un/2,r*ln/2+1,r*ln/2+2,…,r*ln)

R`3=(r*l1,r*l2,…,r*ln/2,r*un/2+1,r*un/2+2,…,r*un)

其中,r為隨機數,且0≤r≤1。若分割因子k=3,則由R`0和R`1產生額外的6(23-2)種反應物R`4、R`5、R`6、R`7、R`8、R`9,依此類推。且

R`4=(r*u1,r*u2,…,r*u2n/3,r*l2n/3+1,

r*l2n/3+2,…,r*ln)

R`5=[r*u1,r*u2,…,r*un/3,r*ln/3+1,…,r*l2n/3,

r*u2n/3+2,…,r*un]

R`6=[r*u1,r*u2,…,r*un/3,r*ln/3+1,…,r*ln]

R`7=[r*l1,r*l2,…,r*ln/3,r*un/3+1,…,r*un]

R`8=[r*l1,r*l2,…,r*ln/3,r*un/3+1,…,

r*u2n/3,r*l2n/3+1,…,r*ln]

R`9=[r*l1,r*l2,…,r*l2n/3,r*u2n/3+1,…,r*un]

假設初始反應物的種類數為R,種群的大小為P,只要滿足R

(1)

置換反應產生的新的反應物為

(2)

(3)

λtd+1=2.3(λtd)2sin(πλtd)

(4)

若R`1是氧化還原反應產生的新的反應物,則

(5)

(6)

(4)反應物更新。在這一步中,通過化學平衡測試來更新反應物。各物質的質量分數不再發生改變的狀態稱為平衡狀態。在平衡狀態下,能使化學系統的焓值減小或熵增大的生成物被選擇作為新反應物,同時排除不良的反應物,從而進行反應物更新。

(5)算法終止條件判斷。符合條件,算法終止,否則轉入步驟(3)。

滿足終止條件時,算法終止。通常的終止條件,如最大迭代次數,或者焓最小、熵最大。否則轉入第(3)步和第(4)步。ACROA流程見圖1。

2 基于化學反應優化算法的支持向量機

SVM是泛化能力較強的機器學習算法。基礎的SVM通過求解二次最優規劃問題尋找一個最優分類超平面,保證每類樣本到最優分類超平面的距離最大,將兩類線性可分的訓練樣本最大限度地分開。當訓練樣本線性不可分時,SVM算法引入懲罰系數或者松弛因子放寬約束,同時引入核函數,將原始樣本從原始非線性空間Rn映射到高維線性特征空間F,從而借助線性分類來解決非線性分類問題。設Strain={(x1,y1),(x2,y2),…,(xN,yN)}為訓練樣本,N為樣本數,xi∈Rn為樣本點,yi∈{1,-1}為樣本分類標簽。則最優分類超平面的判別函數為

圖1 ACROA流程圖

(7)

式中,x為待分樣本;b為偏差;αi為拉格朗日乘積算子;K(x,xi)為核函數。

拉格朗日乘積算子αi通過求解二次最優規劃問題得到:

(8)

(9)

其中,C為懲罰系數或者松弛因子,C的確定需要折中考慮分類邊界最大化和誤差最小化。常用的核函數K(x,xi)有多項式函數、線性函數、高斯徑向基函數、Sigmoid函數等。高斯徑向基函數具有良好的非線性映射能力[14],本文選用高斯徑向基函數作為SVM核函數:

(10)

SVM參數C和σ對SVM的分類性能影響較大,需要進行優化[15]。本文將ACROA應用于SVM參數C和σ的優化求解,提出了ACROA_SVM方法,以SVM分類精度作為適應度函數來優化參數C和σ。ACROA_SVM方法的流程見圖2。首先根據具體問題確定以分類精度為適應度函數、設定參數C和σ的初始值(即初始反應物);然后利用ACROA優化參數C和σ,隨后利用數據子集訓練SVM模型,同時利用測試子集測試SVM模型,計算分類精度;最后以分類精度為依據判斷D和σ是否滿足要求,若滿足要求,則輸出最佳參數C和σ,若不滿足要求,則繼續采用ACROA方法優化,直至滿足要求,循環終止。

圖2 ACROA_SVM方法的流程圖

3 標準數據分析

為了驗證本文提出的ACROA_SVM方法的有效性,從UCI網站上下載標準數據進行分析。實驗下載使用了兩種數據集,分別為iris數據集和thyroid數據集。iris數據集包含3類樣本,分別為Setosa、Versicolor和Virginica,每類樣本各50組數據,共有150組數據,每組數據4個屬性值。thyroid數據集包含3類樣本,共有215組數據,分別為nornal數據150組、hyperthyroid數據35組和hypothyroid數據30組。實驗中抽取三分之二樣本作為訓練樣本,剩下的三分之一樣本作為測試樣本,分別采用ACROA_SVM、GA_SVM、PSO_SVM對以上兩種數據集進行對比分析。由于SVM本質是一種二分類方法,本文采用“一對多”分類方式實現多分類,即針對3類數據設計了兩個分類器——SVM1和SVM2,SVM多分類擴展流程見圖3。 對比分析結果見表1。從分析結果可見,AROCA_SVM不但和GA_SVM、PSO_SVM一樣具有較好的分類能力,而且ACROA_SVM的優化運算速度快,運行時間大大縮短。

圖3 SVM多分類流程圖

方法iris數據集thyroid數據集訓練樣本數目測試樣本數目分類精度(%)時間(s)訓練樣本數目測試樣本數目分類精度(%)時間(s)AROCA_SVM111436100751625310011AROCA_SVM27624100381424310016GA_SVM1114361001491625394.3102GA_SVM2762410021214243100124PSO_SVM1114361001651625388.7524PSO_SVM276241009714243100378

4 ACROA_SVM在旋轉機械故障診斷中的應用

軸承是旋轉機械中的主要元件之一,據統計,旋轉機械30%的故障問題是因軸承引起的。因此,本次實驗采用美國CaseWesternReserveUniversity電氣工程實驗室軸承數據來驗證ACROA_SVM方法能否有效地應用于旋轉機械故障診斷。實驗中,采用1.5kW的三相感應電機驅動風機運轉。電機輸入軸支撐軸承型號為6205-2RSSKF。通過電火花加工技術設置凹坑來模擬單點故障。將帶磁性座底的振動加速度傳感器安放在電機驅動端軸承座上。輸入軸轉速為1750r/min,信號采樣頻率為12kHz,采集7種狀態振動加速度信號各30組數據(共210個樣本)。每個樣本數據長度為2048。實驗中軸承的7種狀態為:正常、輕微外圈故障、輕微內圈故障、輕微滾動體故障、嚴重外圈故障、嚴重內圈故障、嚴重滾動體故障。輕微故障和嚴重故障的凹坑直徑分別為0.18mm、0.53mm。

4.1LMD信號分析與特征提取

圖4為軸承處于不同運行狀態時原始信號時域波形。由圖4可見軸承故障振動加速度信號為非平穩信號,內圈故障和外圈故障振動加速度信號具有明顯的周期沖擊,且具有復雜的調頻-調幅特征,但滾動體故障振動加速度信號沖擊與調制特征幾乎被大量的噪聲淹沒,很難與正常狀態相區別。因此需要采用時頻信號分析方法進行處理。常用的時頻分析方法有小波分析、經驗模態分解(empirical mode decomposition, EMD)、局部均值分解(local mean decomposition, LMD)方法等。但小波變換和EMD方法都有一定的局限性。局部均值分解( local mean decomposition, LMD)是一種新的數據驅動的自適應時頻分析方法,相比EMD方法,具有迭代次數少、端點效應不明顯、得到的虛假分量少等優點;該方法將多分量的調幅調頻信號x(t)分解為若干個頻率由高到低的乘積函數(product function, PF)分量cp,即x(t)可表示為[16]

(11)

式中,uk(t)為信號的平均趨勢。

圖4 軸承不同狀態原始振動加速度信號時域波形

實驗中對每個樣本進行LMD分解,分別得到頻率由高到低的數個PF分量,并利用相關分析方法篩選出真實的PF分量[17]。圖5為內圈故障振動加速度信號的前4個高頻PF分量的時域波形,可見各PF分量的頻率由高到低;而且第一個PF分量的沖擊與調制特征最明顯,即高頻PF分量蘊含了豐富的故障信息,與軸承故障機理理論相符。因此選取前4個高頻PF分量做進一步分析。存在不同類型局部故障的軸承振動信號的能量分布存在明顯的差別,在某些頻帶內能量將增加,而另一些頻帶內PF能量將減少。然而,能量僅能反映信號的強度,不能反映信號的非平穩特性。為了準確獲取信號的本質特征,提取PF能量矩作成特征向量, PF能量矩定義為

Mi=∫t|ci(t)|2dt

(12)

圖5 內圈故障振動加速度信號前4個PF分量時域波形

PF能量矩是PF能量在時間軸的積分,PF能量矩同時關注了信號能量在頻率軸和時間軸上的分布,相對于能量,更能有效地提取信號的本質特征。PF能量矩特征向量計算步驟如下:

(1)采用LMD方法將軸承原始振動加速度信號進行分解,得到n個PF分量ci(i=1,2,…,n)。并以相關系數法篩選出真實PF分量cj(j=1,2,…,m,m≤n)。

(2)選擇包含故障特征的PF分量。對于離散時間序列,按式(12)計算Mj,設第j個分量為1×N維的向量,Δt為采樣時間間隔,N為采樣點數,k為采樣點,則每個真實PF分量能量矩為

(13)

(3)將PF分量能量矩Mj組成特征向量:

M=(M1,M2,…,Mm)

4.2ACROA-SVM故障診斷

在特征提取之后,采用ACROA_SVM進行模式分類識別。本次實驗中對7類軸承狀態的所有原始數據首先進行LMD分解;然后提取PF的能量矩特征,得到7類特征向量樣本集,其中每類樣本集包含30個樣本。分別從每類的樣本集中隨機抽取三分之一樣本作為訓練樣本,訓練ACROA_SVM,得到符合要求的分類器;最后將剩下的三分之二樣本作為測試樣本來測試評價ACROA_SVM的分類性能。故障診斷模型見圖6。

圖6 故障診斷模型

SVM方法本質是一種二分類方法,為此采用“一對多”的方式擴展實現SVM多分類;本次實驗中有7類,設計6個分類器SVM1、SVM2、…、SVM6,軸承故障多分類器設計流程如圖7所示。分析結果見表2~表4。對比分析可見,采用基于LMD分解的PF的能量矩作特征向量,ACROA_SVM,GA_SVM、PSO_SVM均能達到很好的模式識別效果,但ACROA_SVM的時間消耗最短,最有利于實現機械故障在線診斷。

圖7 軸承故障診斷流程圖

分類器訓練樣本數測試樣本數參數C參數σ分類精度(%)時間消耗(s)ACROA_SVM17014011115.71005.8ACROA_SVM2601201.890.340938.5ACROA_SVM35010012375.91007.5ACROA_SVM4408010614.71005.6ACROA_SVM5306062.815.91002.6ACROA_SVM620401140.2301002.2

表3 GA_SVM方法的故障診斷結果

表4 PSO_SVM方法分析結果

5 結論

(1)SVM方法特別適合于樣本稀缺的機械故障診斷,但SVM方法的性能受核參數的影響較大。針對SVM的參數優化問題,首次引入新的優化算法——ACROA,提出了ACROA_SVM方法。

(2)采用UCI標準數據,與現有的GA_SVM、PSO_SVM進行對比分析,結果表明,ACROA_SVM不但具有良好的分類識別效果,而且相比現有GA_SVM、PSO_SVM方法,運算更簡單、更快速。

(3)結合LMD信號分析與PF的能量矩特征提取,將提出的ACROA_SVM方法應用于旋轉機械故障診斷,結果表明了方法的有效性和優越性,有利于實現智能在線故障診斷。

(4)ACROA_SVM可以推廣至其他智能故障診斷領域,將ACROA方法應用多分類SVM模型,值得下一步研究。

[1]陳安華, 周博, 周會福, 等 基于PCA與蟻群算法的機械故障聚類診斷方法[J].中國機械工程, 2013, 24(24): 3333-3337.

Chen Anhua, Zhou Bo, Zhang Huifu, et al. Clustering Method of Mechanical Fault Diagnosis Based on PCA and Ant Colony Algorithm[J].China Machanical Engineering, 2013, 24(24): 3333-3337.

[2]Yuan Shengfa, Chu Fulei. Fault Diagnosis Based on Support Vector Machines with Parameter Optimisation by Artificial Immunisation Algorithm[J]. Mechanical Systems and Signal Processing, 2007, 21(3): 1318-1330.

[3]Wang Huaqing, Chen Peng. Intelligent Diagnosis Method for Rolling Element Bearing Faults Using Possibility Theory and Neural Network[J]. Computers & Industrial Engineering, 2011, 60(4):511-518.

[4]Wang C, Kang Yuan, Shen Pingchen, et al. Applications of Fault Diagnosis in Rotating Machinery by Using Time Series Analysis With Neural Network[J]. Expert Systems with Applications, 2010, 37(2): 1696-1702.

[5]Vapnik V, Chapelle O. Bounds on Error Expectation for Support Vector Machines[J]. Neural Computation, 2000, 12(9): 2013-2036.

[6]Shi X H,Liang Y C, Lee H P, et al. An Improved GA and a Novel PSO-GA Based Hybrid Algorithm[J]. Information Processing Letters, 2005, 93(5): 255-261.

[7]Shelokar P S, Siarry P, Jayaraman V K, et al. Particle Swarm and Ant Colony Algorithms Hybridized for Improved Continuous Optimization[J]. Applied Mathematics and Computation, 2007, 188(1): 129-142.

[8]Xiao Jing, Li Liangping, A Hybrid Ant Colony Optimization for Continuous Domains[J]. Expert Systems with Applications, 2011, 38(9): 11072-11077.

[9]An Senjian, Liu Wanquan, Venkatesh S. Fast Cross-validation Algorithms for least Squares Support Vector Machine and Kernel Ridge Regression[J]. Pattern Recognition, 2007, 40(8): 2154-2162.

[10]Diosan L, Rogozan A, Pecuchet J P. Improving Classification Performance of Support Vector Machine by Genetically Optimising Kernel Shape and Hyper-Parameters[J]. Applied Intelligence, 2010, 36(2): 280-294.

[11]張弦, 王宏力. 基于粒子群優化的最小二乘支持向量機在時間序列預測中的應用[J]. 中國機械工程, 2011, 22(21): 2572-2576.

Zhang Xian, Wang Hongli. LSSVM Based on PSO and Its Application to Time Series Prediction[J].China Machanical Engineering, 2011, 22(21): 2572-2576.

[12]Alatas B. A Novel Chemistry Based Metaheuristic Optimization Method for Mining of Classification Rules[J]. Expert Systems with Applications, 2012, 39: 11080-11088.

[13]Alatas B. ACROA: Artificial Chemical Reaction Optimization Algorithm for Global Optimization[J]. Expert Systems with Applications, 2011, 32: 13170-13180.

[14]Lai Xingwu, Tse P W, Zhang Guicai, et al. Classification of Gear Faults Using Cumulants and the Radial Basis Function Network[J]. Mechanical Systems and Signal Processing, 2004, 18: 381-389.

[15]Zhang Xiaoli, Chen Xuefeng, He Zhengjia. An ACO-based Algorithm for Parameter Optimization of Support Vector Machines[J]. Expert Systems with Applications, 2010, 37: 6618-6628.

[16]Jonathan S S. The Local Mean Decomposition and Its Application to EEG Perception Data[J]. Journal of the Royal Society Interface, 2005, 2(5): 443-454.

[17]Peng Zhike, Peter W, Chu Fulei. A Comparison Study of Improved Hilbert-Huang Transform and

Wavelet Transform:Application to Fault Diagnosis for Rolling Bearing[J]. Mechanical Systems and Signal Processing, 2005, 19(1): 974-988.

(編輯王艷麗)

SVM Based on ACROA and Its Applications to Rotating Machinery Fault Diagnosis

Luo Songrong1,2Cheng Junsheng2Hunglinh A O2

1.Hunan University of Arts and Science,Changde,Hunan,415003 2.Hunan University,Changsha,410082

Firstly,in view of SVM parameters optimization problem, combination to the advantage of ACROA, a new classification model, called ACROA_SVM was presented herein. Furthermore, the effectiveness and superiority of the ACROA_SVM model was identified via benchmark datasets, which was downed from the sit web of UCI. Lastly, combination to local mean decomposition and energy moment feature extraction, ACROA_SVM was served as approach of pattern recognition to identify rotating machinery fault types. The experimental results show ACROA_SVM method has higher precision, better generalization ability of fault diagnosis, and less time consumption, higher efficiency of fault diagnosis, which is conducive to realize online intelligent fault diagnosis.

support vector machine(SVM); artificial chemical reaction optimization algorithm(ACROA); rotating machinery; fault diagnosis

2014-01-17

國家自然科學基金資助項目(51175158, 51075131);湖南省教育廳科研項目(14C0789);湖南省“十二五”重點建設學科項目(湘教發2011[76])

TH165.3;TN911.7DOI:10.3969/j.issn.1004-132X.2015.10.006

羅頌榮,女,1973年生。湖南文理學院機械工程學院副教授,湖南大學機械與運載工程學院博士研究生。主要研究方向為機械設備狀態監控與故障診斷、動態信號處理與分析。發表論文10余篇。程軍圣,男,1968年生。湖南大學機械與運載工程學院教授、博士研究生導師。Hunglinh AO,男,1974年生。湖南大學機械與運載工程學院博士研究生。

猜你喜歡
故障診斷分類優化
超限高層建筑結構設計與優化思考
房地產導刊(2022年5期)2022-06-01 06:20:14
民用建筑防煙排煙設計優化探討
關于優化消防安全告知承諾的一些思考
分類算一算
一道優化題的幾何解法
分類討論求坐標
數據分析中的分類討論
教你一招:數的分類
因果圖定性分析法及其在故障診斷中的應用
基于LCD和排列熵的滾動軸承故障診斷
主站蜘蛛池模板: 亚洲婷婷六月| 黄色网址手机国内免费在线观看| 日韩国产黄色网站| 狼友视频一区二区三区| 中文字幕日韩久久综合影院| 91无码人妻精品一区二区蜜桃| 亚洲免费播放| 呦视频在线一区二区三区| 成人国产一区二区三区| 91青青草视频| 国产剧情一区二区| 久久96热在精品国产高清| 456亚洲人成高清在线| 欧美啪啪网| 久久青草精品一区二区三区| 亚洲一区二区无码视频| 综合亚洲网| 国产男女免费视频| 国产亚洲视频免费播放| 久久精品亚洲热综合一区二区| 在线精品亚洲一区二区古装| 91在线激情在线观看| 久久大香香蕉国产免费网站| 日韩大乳视频中文字幕 | 精品一区二区三区视频免费观看| 日本一区二区不卡视频| 国产第一页免费浮力影院| 国产欧美亚洲精品第3页在线| 国产高颜值露脸在线观看| 亚洲成肉网| 亚洲综合亚洲国产尤物| 国产人成乱码视频免费观看| 久久九九热视频| 久久精品欧美一区二区| 免费看的一级毛片| 91精品国产91欠久久久久| 亚洲免费毛片| 国产十八禁在线观看免费| 狠狠干欧美| 精品国产自| 中文字幕久久亚洲一区| 孕妇高潮太爽了在线观看免费| 成人a免费α片在线视频网站| 国产综合无码一区二区色蜜蜜| 一级爱做片免费观看久久| 青青青草国产| 久久99精品久久久久久不卡| 无码精油按摩潮喷在线播放 | 香蕉精品在线| 免费无码网站| 亚洲精品777| 久久无码高潮喷水| 激情六月丁香婷婷四房播| 国产另类乱子伦精品免费女| 免费大黄网站在线观看| a级毛片免费看| 五月综合色婷婷| 人妻21p大胆| 看国产毛片| 国产黑丝视频在线观看| 乱色熟女综合一区二区| 亚洲国产成人精品无码区性色| 99久久国产综合精品女同| 91丝袜美腿高跟国产极品老师| 亚洲最大福利网站| 国产一级无码不卡视频| 国产在线精品人成导航| 又大又硬又爽免费视频| 色老头综合网| 中文字幕1区2区| 国产精品第5页| 全色黄大色大片免费久久老太| 亚洲码一区二区三区| 九色综合伊人久久富二代| 无码aaa视频| www.youjizz.com久久| 精品国产一区二区三区在线观看| 好紧太爽了视频免费无码| 玖玖精品在线| 精品人妻无码区在线视频| 日本91视频| a在线观看免费|