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

基于多穩態隨機共振的軸承微弱故障信號檢測*

2017-01-09 05:48:50陸寶春張登峰
振動、測試與診斷 2016年6期
關鍵詞:故障信號模型

馮 毅, 陸寶春, 張登峰

(南京理工大學機械工程學院 南京,210094)

?

基于多穩態隨機共振的軸承微弱故障信號檢測*

馮 毅, 陸寶春, 張登峰

(南京理工大學機械工程學院 南京,210094)

針對雙穩態隨機共振模型無法有效處理調制信號的缺點,提出了一種以包絡信號為輸入信號的自適應多穩態級聯隨機共振(adaptive multi-stable cascaded stochastic resonance,簡稱AMCSR)信號強化方法。首先,對振動信號進行包絡解調,依據包絡信號分布特點,選用與信號分布相匹配的多穩態隨機共振模型;然后,以故障特征頻率的頻譜幅值為指標,采用蟻群算法自適應地優化隨機共振模型參數;最后,以噪聲為強化源和驅動信號,通過級聯隨機共振方法對包絡信號中的故障特征頻率進行逐級強化,獲得故障特征成分的強化信號。對實測軸承振動信號的驗證結果表明,該方法能夠增強故障特征頻率成分,有效地提取被其他頻率成分淹沒的微弱故障信號。

滾動軸承; 包絡信號; 多穩態; 隨機共振

引 言

隨機共振(stochastic resonance,簡稱SR)是一種有效的微弱信號處理方法,它可以簡單理解為:通過非線性系統作用,噪聲能量向微弱信號轉移,從而加強了原本微弱的信號。隨機共振的產生必須滿足3個要素,即具有勢能壘的非線性系統、微弱輸入信號以及噪聲源,當三要素達到最佳匹配關系時,隨機共振效應最強[1-2]。軸承出現早期故障時,由于故障信號過于微弱,傳統方法對信號的處理結果并不理想,采用隨機共振方法可有效強化微弱信號。文獻[1]通過二次采樣變換,將僅在小參數范圍成立的經典隨機共振理論擴展到大參數范圍,為隨機共振的工程應用提供了理論依據。文獻[2]提出了級聯隨機共振方法,有效地提高了隨機共振對微弱信號的增強效果。文獻[3]采用蟻群算法優化雙穩態隨機共振參數,成功地提取了行星齒輪的故障信號。文獻[4-7]采用雙穩態隨機共振方法,提高軸承信號的信噪比,再通過經驗模態分解處理降噪后的信號,取得了一定效果。文獻[8]采用倒頻譜白化方法預處理輸入信號,提高了雙穩態模型對故障特征的強化效果。文獻[9]將Woods-Saxon勢函數引入隨機共振,優化了雙穩態隨機共振模型。文獻[10]提出了可增加噪聲利用率的多穩態模型,其結構特點可處理非對稱分布的包絡信號。

以上研究表明,雙穩態隨機共振方法在機械故障診斷中具有顯著的效果,但雙穩態隨機共振同時存在一定缺陷:a.故障軸承的振動信號具有調制特點,而故障信號特征頻率成分主要存在于包絡信號中,雙穩態模型將滿足正態分布的原始信號作為隨機共振系統的輸入,故障特征成分的強化效果并不明顯;b.故障信號為復合信號,以信噪比為優化指標時,無法準確反映復合信號中某一特定故障特征頻率的強化效果。針對問題a引入可處理非對稱分布信號的多穩態模型[10],以含有故障特征頻率的包絡信號作為輸入信號,從而克服雙穩態隨機共振無法處理調制信號的問題;針對問題b,以故障特征頻率的頻譜幅值為優化指標,通過蟻群算法進行模型參數優化,可實現故障特征頻率的選擇性強化。

基于以上分析,提出一種以故障頻率頻譜幅值為自適應優化指標,以包絡信號為輸入的自適應多穩態級聯隨機共振方法,可解決雙穩態隨機共振的兩種缺陷。試驗對比驗證了AMCSR方法的有效性。

1 雙穩態模型及其缺陷

1.1 雙穩態隨機共振模型

雙穩態模型[1-2,11]具有中間勢壘和兩個對稱勢阱,雙穩態勢函數表示為

(1)

對應的Langevin方程為

(2)

其中:s(t)為輸入信號;n(t)為均值為0的高斯白噪聲,其均方根為σ。

方程描述了布朗粒子同時受到輸入信號和噪聲驅動時,在雙勢阱中的過阻尼運動[1]。雙穩態隨機共振模型如圖1所示。

圖1 雙穩態隨機共振模型Fig.1 Stochastic resonance model potential function

1.2 雙穩態隨機共振的缺陷

采用文獻[11]中以信噪比為優化指標所獲得的最優雙穩態模型,通過輸入不同類型的信號,驗證雙穩態模型隨機共振的缺陷。

最優模型參數取值a=b=1,最優共振頻率(在已知參數的隨機共振系統中,唯一能得到最大程度強化的頻率)信號s(t)=0.5cos(2π×0.1t),即信號頻率f=0.1 Hz,噪聲均方根σ=5,采樣頻率fs=20 Hz,采樣點數N=2 000。

不同類型輸入信號的隨機共振強化效果如下:

1) 在單一信號成分條件下,圖2(a)中被噪聲完全淹沒的輸入信號在經過隨機共振系統處理后,輸出信號的信噪比大為提高,頻譜中0.1Hz頻率的幅值非常顯著,證明隨機共振系統實現了噪聲能量向目標信號能量的轉換,從而加強了目標信號的能量,如圖2(b)所示;

圖2 復合信號隨機共振Fig.2 Stochastic resonance of composite signal

2) 當輸入信號為s(t)=0.5cos(2π×0.1t)+0.6sin(2π×0.15t)+0.4sin(2π×0.05t),是以0.1 Hz為頻率中心的復合疊加信號,隨機共振仍能有效加強信號中的3種頻率分量,如圖2(c)所示;

3) 當輸入信號為s(t)=0.5cos(2π×0.1t)+0.6sin(2π×0.5t),信號成分0.5 Hz已遠離0.1 Hz頻帶中心,0.1 Hz分量仍得到有效強化,如圖2(d)所示;

4) 輸入信號為s(t)=|0.5sin(2π×0.05t)|×cos(2π×0.5t)時,較高頻率的0.5 Hz信號為主要成分,0.1 Hz調幅信號為包絡信號,隨機共振系統無法有效加強0.1 Hz頻率成分,如圖2(e)所示。

仿真結果表明雙穩態隨機共振存在以下缺陷:

1) 圖2(c)與圖2(b)結果對比表明,以信噪比為優化指標時,隨機共振系統會強化最優共振頻率附近的所有頻率成分,無法針對性地強化0.1 Hz的最優共振頻率;

2) 圖2(e)與圖2(d)結果對比表明,0.1Hz為包絡信號時,雙穩態隨機共振系統無法有效強化具有調制特性的頻率成分。

2 自適應多穩態級聯隨機共振

滾動軸承局部損傷對振動信號具有調制作用,而包絡信號可有效地反映此類故障特征[12]。包絡信號是非零均值的非對稱信號,若用雙穩態模型處理,在噪聲影響下,上包絡信號會越過中間勢壘進入負半區勢阱[1-2],則會引入無意義的頻率成分。

針對雙穩態隨機共振方法存在的兩種缺陷,筆者提出了“自適應多穩態級聯隨機共振”方法,可對包絡信號進行隨機共振處理,實現故障特征頻率的針對性強化。

2.1 多穩態隨機共振模型

文獻[10]提出了可增加噪聲利用率的多穩態模型,其結構特點可以處理非對稱信號。多穩態模型具有3個勢阱和2個勢壘的結構特點,如圖3所示,其勢函數可表示為

(3)

多穩態模型中的參數a,b均影響兩側勢阱深度,且作用效果相反,參數c影響兩勢壘高度,通過調節系統參數a,c以及采樣頻率fs,獲得最佳共振模型。

圖3 多穩態模型勢函數Fig.3 Multi-stable model potential function

2.2 包絡信號多穩態隨機共振

AMCSR方法具體處理步驟如下:

1) 采用Hilbert變換對振動信號進行解調,并提取軸承振動信號的上下包絡信號;

2) 依據包絡信號和包絡譜幅值以及故障特征頻率,按照大參數隨機共振尺度變換方法[1]確定噪聲強度以及初始共振參數范圍;

3) 以待檢測的故障特征頻率頻譜幅值為優化指標,采用蟻群算法進行參數優化得到最優共振參數;

4) 以噪聲為驅動信號,對上下包絡信號進行多穩態級聯隨機共振處理,獲得強化信號。

其中,以噪聲為驅動信號,即噪聲信號在任意時刻t的瞬時幅值為正(負),則以上(下)包絡信號作為t時刻的系統輸入信號,實現包絡信號在多穩態模型中的隨機共振。

(4)

其中:envu(t)為上包絡信號;envd(t)為下包絡信號。

級聯隨機共振[2]是以上一級隨機共振的輸出作為下一級隨機共振輸入的方式,將多個相同的隨機共振模型串聯,實現對微弱信號的逐級強化。多次仿真試驗結果表明,采用兩級隨機共振級聯,可以達到較理想的輸出強化效果。包絡信號多穩態隨機共振算法流程如圖4所示。

圖4 算法流程Fig.4 Algorithm flow

2.3 自適應參數優化

多穩態隨機共振系統對信號的處理效果受到多個參數的共同作用影響,若單一考慮某個參數則無法獲得復雜系統的最優解。蟻群算法源于螞蟻群體在覓食過程中,通過對路徑上信息素濃度的更新,選擇最優路徑的現象,蟻群算法可有效解決多參數的尋優問題。

“小參數”(信號頻率低于1Hz)隨機共振并不適用于實際滾動軸承信號,因而需要對多穩態模型共振參數a,b,c,fs進行尺度變換,以適應“大參數”條件(信號頻率遠大于1Hz)的隨機共振需求[1]。以文獻[10]中的“小參數”條件下多穩態模型中的參數值為基礎,采用文獻[1]中的尺度變換方法得到“大參數”條件下的a,b,c值,并以此為各參數優化區間的中間值,通過多次試驗,確定參數的優化區間(參數b固定不變)。同樣,以實際采樣頻率作為fs的優化區間中間值,確定優化區間。

采用蟻群算法對共振模型參數進行自適應優化,獲得最優參數。以滾動體故障為例,步驟如下:

1) 初始化ai→cj→fsk參數路徑信息素,將所有路徑的信息素置為相同初值τ;

2) 依據式(5)計算滾動體故障特征頻率fball;

3) m只螞蟻以隨機ai為起點,按照參數路徑,依據式(6)計算出的跳轉概率進行下一步路徑選擇,直至m只螞蟻均完成路徑中的3種參數選擇;

7) 依據更新后的路徑信息素,重復步驟3~5,直至達到迭代次數上限K,獲得全局最優參數。

上述步驟中,滾動體故障特征頻率為

(5)

其中:ω為軸的轉頻;D為軸承節徑;d為滾動體直徑;α為接觸角。

(6)

其中:n1,n2,n3分別為參數a,c,fs在尋優范圍內的取值點數。

滾動體故障特征頻率的頻譜幅值大小A由離散傅里葉變換得到

(7)

其中:n為故障特征頻率對應的頻域離散值序號;x(kΔt)為信號采樣值;Δt為采樣間隔;k為時域離散值的序號N為序列點數。

信息素的揮發與增強表示為

(8)

其中:k為當前迭代次數;ρ為信息素揮發因子。

3 示例驗證

3.1 自適應多穩態隨機共振算法驗證

算法驗證采用的數據來自美國Case Western Reserve University軸承數據中心。軸承類型為6205-ZRSJMESKF深溝球軸承,信號采樣頻率為12 kHz,采樣時間為1 s,驅動軸轉速為1 797 r/min。軸承故障都預先人為加工,加工的故障部位分別為軸承內圈、外圈以及滾動體,深度尺寸均為0.177 8 mm。通過經驗公式計算得到外圈、內圈和滾動體故障特征頻率分別為107.3,162.2和141.5 Hz(實測頻率存在誤差)。

外圈和內圈故障的包絡譜圖5(a),(b)中含有明顯的故障特征頻率,可判斷內圈、外圈存在故障。由于滾動軸承結構特點,滾動體在運行中存在公轉、自轉以及打滑,工作狀態相比外圈和內圈更為復雜,滾動體的特征信號非常微弱,往往被其他頻率成分淹沒,致使滾動體故障難以確定。

圖5(c)中難以觀察到明顯的故障特征,且與內圈、外圈故障特征具有明顯區別,可判斷該軸承不存在內圈、外圈故障,但可能存在滾動體故障。依據經驗公式,對待測信號中的滾動體特征頻率進行隨機共振增強檢測。

圖5 軸承信號包絡譜Fig.5 Bearing signal envelope spectrum

對待測軸承信號進行Hilbert變換得到上、下包絡信號,如圖6(a)所示。依據包絡信號幅值波動范圍,取噪聲參數σ=3.5。由于白噪聲頻譜為均勻分布,含噪聲的包絡信號頻譜中,所有頻率成分得到一定增強,但故障特征頻率仍不明顯,如圖6(b)所示。

圖6 待測包絡信號和含噪聲包絡譜Fig.6 Rolling envelope signal and spectrum

依據大參數隨機共振尺度變換方法[1],設置參數b=1 000固定不變,其余參數的優化區間為a∈[500,2 000],c∈[50,800],fs∈[4 000,16 000],取值步長分別為sa=2,sc=1,sf=5。

依據多次試驗輸出所得故障特征頻率幅值的平均值,為保證在算法趨于收斂時信息素更新的增強幅度逐漸減緩,設置螞蟻數m=100,全局路徑信息素初值τ=0.2,信息素揮發因子ρ=0.05,設置迭代次數為K=2 000。

算法收斂時得到全局最優參數,將最優二次采樣頻率fs固定,分析參數a,c對幅值A的影響。參數a,c同時改變時,輸出結果會出現局部極值,若單一考慮某一參數則無法獲得最優結果。尋優結果表明,在參數a∈[1 300,1 400],c∈[400,500]區域內存在全局最優解,蟻群算法對多參數優化具有良好效果,如圖7所示。

圖7 蟻群算法尋優結果Fig.7 Ant colony optimization algorithm results

經過蟻群算法優化后的隨機共振參數為a=1 364,c=483,fs=7 835。將最優參數賦值于隨機共振系統,采用最優模型以及最優二次采樣頻率對輸入信號進行處理,再按原有采樣頻率12 kHz繪制輸出信號頻譜。

依據最優共振模型參數,對含有噪聲信號的上、下包絡信號進行兩級多穩態隨機共振處理。由圖8多穩態級聯隨機共振輸出可知,第2級隨機共振輸出信號主要集中在以加速度幅值0,+2,-2 m/s2為中心的鄰域范圍內,對應多穩態模型的三勢阱結構,呈現出明顯的多穩態特征。

故障特征頻率的頻譜幅值從隨機共振增強之前的0.051 99經過第1級隨機共振增強至0.179 6,如圖8(a)所示。同時,頻譜中故障特征頻率的邊頻帶成分也受到一定程度增強,但故障特征頻率的增強不夠顯著。經過第2級隨機共振,故障特征頻率的頻譜幅值大幅度加強至0.228,其頻譜幅值明顯高于邊頻帶成分,如圖8(b)所示。經過兩級多穩態隨機共振,輸出信號中的噪聲成分減少,噪聲能量逐漸向故障特征頻率信號轉移。

對于滾動體故障特征頻率的倍頻282.8 Hz,采用同樣方法獲得相應的最優共振參數,并對其進行強化,輸出結果如圖9所示。

圖9 倍頻隨機共振輸出Fig.9 Fault multiplier stochastic resonance output

以上試驗數據驗證表明,AMCSR方法有效地增強了待測信號中的滾動體故障特征頻率,與正常軸承相比,增強后的待測軸承信號中含有明顯的滾動體故障特征頻率成分,可判斷待測軸承存在滾動體故障。

3.2 雙穩態隨機共振算法對比

文獻[5]以未經過包絡解調的待測軸承信號經過帶通濾波處理后作為雙穩態系統輸入,以信噪比為優化指標進行共振模型參數優化。圖10頻譜中故障特征頻率的幅值為0.000 121 6,與其他頻率成分相差102量級,能量極其微弱,且遠遠小于圖5(c)包絡譜中故障特征頻率的頻譜幅值0.003 34。

圖10 帶通濾波輸入信號和頻譜Fig.10 Band-pass filtering input signal and spectrum

圖11 參數與幅值關系Fig.11 Amplitude with parameters

將多穩態模型下的最優變尺度采樣頻率8 kHz作為雙穩態模型尋優的二次采樣頻率,繪制參數a∈[10,1 500],b∈[10,1 000]與故障特征頻率頻譜幅值的關系譜,如圖11所示。

由于輸入信號未經過解調,導致強化后故障特征頻率的頻譜幅值仍然遠低于其他頻率成分。對于具有調制特性的故障軸承信號,雙穩態隨機共振的處理效果并不理想,如圖12所示。

圖12 雙穩態模型隨機共振輸出和頻譜Fig.12 Stochastic resonance output of bistable model

4 結 論

1) 仿真信號分析與文獻算法驗證表明,以信噪比為優化指標的雙穩態隨機共振系統,會對共振頻帶中心內的頻率成分無選擇地加強,但對寬帶信號和調制信號的加強效果并不理想。

2) 實例驗證表明,上下包絡信號多穩態級聯隨機共振方法對存在調制特性的故障軸承信號具有良好的強化效果。

3) 以故障特征頻率頻譜幅值為優化指標,自適應地優化多穩態模型參數,能針對性地增強軸承滾動體故障特征頻率成分。

4) 隨機共振方法存在一定不足,即需要已知故障特征頻率的先驗知識,因此隨機共振適用于微弱信號的增強,而難以反映信號的細節特征。

[1] 冷永剛.雙穩調參高頻共振機理[J].物理學報,2011,60(2):1-7.

Leng Yonggang. Mechanism of high frequency resonance of parameter-adjusted bistable system[J].Acta Physica Sinica, 2011,60(2): 1-7.(in Chinese)

[2] 冷永剛,王太勇,郭焱.級聯雙穩系統的隨機共振特性[J].物理學報,2005,54(3):1118-1125.

Leng Yonggang, Wang Taiyong,Guo Yan. Stochastic resonance behaviors of bistable systems connected in series[J].Acta Physica Sinica,2005, 54(3): 1118-1125.(in Chinese)

[3] Lei Yaguo, Han Dong, Lin Jing.Planetary gearbox fault diagnosis using an adaptive stochastic resonance method[J].Mechanical Systems and Signal Processing, 2013, 38(1): 113-124.

[4] 張海如,王國富,張法全,等.改進的隨機共振和EMD混合模型用于轉子早期故障檢測[J].電機與控制學報,2014,18(2):83-89.

Zhang Hairu,Wang Guofu,Zhang Faquan,et a1. Conjoint model combining improved stochastic resonance and empirical mode decomposition for rotor incipient faults detection[J].Electric Machines and Control, 2014, 18(2): 83-89.(in Chinese)

[5] 韓東穎,丁雪娟,時培明.基于自適應變尺度頻移帶通隨機共振降噪的EMD多頻微弱信號檢測[J].機械工程學報,2013,49(8):10-18.

Han Dongying,Ding Xuejuan,Shi Peiming. Multi-frequency weak signal detection based on emd after denoising by adaptive rescaling frequency-shifted band-pass stochastic resonance[J].Journal of Mechanical Engineering, 2013, 49(8):10-18.(in Chinese)

[6] Li Jimeng, Chen Xuefeng, He Zhengjia. Adaptive stochastic resonance method for impact signal detection based on sliding window [J]. Mechanical Systems and Signal Processing, 2013, 36(2): 240-255.

[7] 雷亞國,韓冬,林京,等.自適應隨機共振新方法及其在故障診斷中的應用[J].機械工程學報,2012,48(7):62-67.

Lei Yaguo, Han Dong, Lin Jing,et a1.New adaptive stochastic resonance method and its application to fault diagnosis[J].Journal of Mechanical Engineering, 2012, 48(7): 62-67.(in Chinese)

[8] 張曉飛,胡蔦慶,胡雷,等.基于倒譜預白化和隨機共振的軸承故障增強檢測[J].機械工程學報,2012,48(23):83-89.

Zhang Xiaofei, Hu Niaoqing, Hu Lei, et a1.Enhanced detection of bearing faults based on signal cepstrum pre-whitening and stochastic resonance[J].Journal of Mechanical Engineering,2012, 48(23): 83-89.(in Chinese)

[9] Lu Siliang, He Qingbo, Kong Fanrong. Stochastic resonance with woods-saxon potential for rolling element bearing fault diagnosis[J].Mechanical Systems and Signal Processing, 2014, 45(2): 488-503.

[10]Li Jimeng , Chen Xuefeng , He Zhengjia. Multi-stable stochastic resonance and its application research on mechanical fault diagnosis[J].Journal of Sound and Vibration,2013, 332(22): 5999-6015.

[11]胡蔦慶.隨機共振微弱特征信號檢測理論與方法[M].北京:國防工業出版社,2012:31-59.

[12]馬川,李宏坤,趙利華.運用小波包峭度包絡的滾動軸承故障診斷[J].振動、測試與診斷,2011,31(6):720-723.Ma Chuan, Li Hongkun, Zhao Lihua. Rolling bearing fault diagnosis using wavelet packet-kurtosis-envelope[J].Journal of Vibration, Measurement and Diagnosis,2011,31(6):720-723.(in Chinese)

10.16450/j.cnki.issn.1004-6801.2016.06.021

*國家自然科學基金資助項目(51275245;61374133);江蘇省“六大人才高峰”計劃資助項目(2011-ZBZZ-011)

2014-10-16;

2014-10-29

THl65.3; TN911

馮毅,男,1989年10月生,博士生。主要研究方向為軸承故障診斷及機械系統信號處理。 E-mail:574102101@qq.com

猜你喜歡
故障信號模型
一半模型
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
重要模型『一線三等角』
完形填空二則
重尾非線性自回歸模型自加權M-估計的漸近分布
故障一點通
基于FPGA的多功能信號發生器的設計
電子制作(2018年11期)2018-08-04 03:25:42
3D打印中的模型分割與打包
奔馳R320車ABS、ESP故障燈異常點亮
基于LabVIEW的力加載信號采集與PID控制
主站蜘蛛池模板: AV老司机AV天堂| 在线免费看片a| 国产成人福利在线视老湿机| 国产拍在线| 亚洲娇小与黑人巨大交| 2021国产精品自拍| 久久精品这里只有精99品| 中文字幕一区二区人妻电影| 久久免费成人| 欧美一区二区三区不卡免费| 99精品免费在线| 性视频一区| 国产拍揄自揄精品视频网站| 91精品视频播放| 免费看美女自慰的网站| 国产精品亚洲天堂| 精品国产Av电影无码久久久| 国产91成人| 亚洲综合日韩精品| 免费观看男人免费桶女人视频| 国产www网站| 国产精品3p视频| 国产欧美精品一区二区| 四虎永久在线精品国产免费| 一本大道无码日韩精品影视| 国产成人综合欧美精品久久| 亚洲国产日韩欧美在线| 高清免费毛片| 亚洲欧美成人在线视频| 久久99这里精品8国产| 波多野结衣一级毛片| 亚洲最大情网站在线观看| 国产在线观看成人91| 欧美一级黄片一区2区| 美女潮喷出白浆在线观看视频| 亚洲五月激情网| 香港一级毛片免费看| 亚洲欧美一级一级a| 成人国产免费| 秋霞午夜国产精品成人片| 蜜臀av性久久久久蜜臀aⅴ麻豆 | 欧美色综合网站| 伊伊人成亚洲综合人网7777| 中文字幕永久视频| 亚洲高清在线播放| 亚洲AV成人一区二区三区AV| 综合人妻久久一区二区精品| 一本大道AV人久久综合| 国产爽爽视频| 国产成人精品一区二区三在线观看| 中文精品久久久久国产网址 | 亚洲精品在线影院| 亚洲欧美激情另类| 青青草国产精品久久久久| 亚洲欧美另类日本| 激情在线网| 欧美视频在线第一页| 青青草欧美| 国产中文在线亚洲精品官网| 久久精品视频亚洲| 久久综合成人| 2018日日摸夜夜添狠狠躁| 免费高清毛片| 亚洲福利网址| 凹凸国产熟女精品视频| 凹凸精品免费精品视频| 精品人妻无码区在线视频| 日韩免费视频播播| 99精品视频播放| 亚洲欧美不卡视频| 九色视频在线免费观看| 激情综合婷婷丁香五月尤物| 亚洲日韩欧美在线观看| 四虎永久免费地址| 国产日韩欧美中文| 亚洲第一极品精品无码| 亚洲欧美日韩成人高清在线一区| 国内精品小视频在线| 黄色三级毛片网站| 久久a毛片| 欧美成人看片一区二区三区 | 国产精品yjizz视频网一二区|