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

變工況滾動軸承故障診斷方法綜述

2022-09-21 05:37:28胡春生李國利成芳娟
計算機工程與應用 2022年18期
關鍵詞:故障診斷振動特征

胡春生,李國利,趙 勇,成芳娟

寧夏大學 機械工程學院,銀川750021

在現代化裝備制造業中,旋轉機械設備正朝向重載化、自動化、智能化的趨勢不斷發展,其結構更加精密,功能更加復雜,是現代裝備不可或缺的重要組成部分[1-3]。滾動軸承常被用于支撐旋轉體、減小摩擦以及保證回轉精度等,是旋轉機械傳動裝置最重要的部件之一。由于軸承長時間高速旋轉,常常伴隨著磨損、高溫、重載等嚴苛的工作條件,因此故障頻發。研究顯示,軸承發生故障是很多重大安全事故的誘因,一旦軸承發生故障,輕則停機停產,重則造成人員傷亡,給生產生活帶來巨大災難[4]。因此,對軸承運行狀態進行實時監測以及故障診斷具有重大安全意義。經過國內外科研人員近40年的研究,軸承故障檢測與診斷技術發展了大量的有效方法,如溫度、油屑分析、振動、聲發射、靜電、超聲等[5]。由于振動數據便于分析并容易獲得,因此絕大多數的研究都是通過分析振動信號來進行故障診斷。

當前傳統恒定工況的軸承故障診斷方法主要分為兩種,一種是以信號分析為主的傳統故障診斷方式,如EMD[6]、集合經驗模態分解(EEMD)[7]和奇異值分解(SVD)[8]等,該方法依靠人工提取特征值,因此該方法十分依賴專家經驗;另一種為以數據驅動為主的機器學習方法,如RNN(包括LSTM和GRU)[9-10]、Logistic回歸[11]、支持向量機(SVM)[12]、人工神經網絡(ANN)[13]和CNN[14]等,該方法無需人工參與,可自動提取信號特征,且準確率大大提高,因此被廣泛應用于故障診斷領域。

然而在實際工程中,軸承運行環境是變轉速、變負載的。例如,汽車的車速和載重變化、風力發電機組風速風向的變化、生產線上產品需求的變化等因素都會導致其軸承運行工況的變化,變化的工況導致振動信號出現幅值變化、脈沖間隔變化、采樣相位不恒定等問題,從而造成故障特征的動態變化,僅靠單一手段進行特征提取的傳統方法無法解決變工況條件下數據分布的差異問題,給變工況條件下的故障診斷帶來巨大挑戰。當前解決變工況問題的方法主要有三種,分別是信號解調和信號處理、深度學習和遷移學習方法。本文回顧近五年變工況軸承故障診斷方法進展,簡要地介紹了各類算法的原理,總結其適用場景及優勢與不足,并對未來的研究方向做出展望,為后續研究指明方向。

1 信號解調與信號處理方法

當軸承處于健康的運行狀態時,產生的振動信號具有平穩性與隨機性,但是一旦軸承發生故障,將出現周期性的沖擊特性與一定的非平穩性[15]。由于基于時域的均方根(RMS)、偏度或峰度等統計量會因速度變化而變化,因此以這些統計量作為特征的故障診斷方法是無效的,在頻域上會引起頻域模糊現象。解決以上問題的方法主要有階次跟蹤以及時頻分析等方法。

1.1 階次跟蹤法

軸承振動信號一般通過在驅動端、風扇端安裝加速度傳感器,按照恒定的時間間隔測得,由于軸承轉動速度不斷發生變化,導致采集到的信號不再具有周期性與平穩性。階次跟蹤法是解決變速問題中應用廣泛的故障診斷方法之一。與按照固定時間間隔采樣方式不同,最初的階次跟蹤通過硬件對振動信號進行等角度重采樣,因此也被稱為硬件階次跟蹤(HOT)[16],該方法采集速度塊,可在線實現。然而,當軸承轉速變化較快時,需不斷動態調整采樣頻率及濾波器的截止頻率,導致出現一定程度的采樣延遲和誤差,而且,該方法需要可調模擬濾波器和模數轉化器等硬件支持,不僅增加了系統的復雜性及成本,同時限制了硬件階次跟蹤法的應用。

為解決硬件階次跟蹤的局限性,發展了計算階次跟蹤法[17]。該方法以固定時間間隔獲取振動信號并實時獲得轉速信息,然后根據轉速計脈沖計算出的旋轉相位通過數值插值算法對振動信號重采樣,從而將時域振動信號轉化為角域振動信號。通過信號分析等手段,提取故障特征,用于滾動軸承故障診斷。針對變速引起的頻域模糊現象,張亢等[18]使用階次跟蹤法(COT)將時域振動信號重采樣為角域振動信號后,對角域信號進行局部均值分解(LMD)得到多個乘積函數分量,通過分析各個乘積函數分量的頻譜判斷軸承的故障類型和損傷程度。SPRO、SPRI 和SPRR 等頻域特征僅在平穩運行條件下有效表征軸承缺陷,但是在變速條件下將會失效,Farhat等[19]基于計算階次分析技術,更新了SPRO、SPRI和SPRR 三個頻域特征,進而提取故障特征。楊武成[20]利用計算階次分析法對振動信號角域重采樣后,進一步利用Hilbert 包絡解調把低頻故障信號從高頻載波信號中解調出來,進而提取故障特征。尹學慧[21]認為Hilbert變換在軸承轉速波動較大的情況下并不適用,通過對實驗數據進行Envelope 函數包絡、Hilbert 變換包絡、平方包絡等包絡計算發現,基于Envelope函數的包絡能更有效地計算出變工況下故障信號的調制階次。基于此,通過階次跟蹤技術與Envelope函數包絡相結合,研究了軸承的故障機理和不同轉速下軸承的階次值,從而在變速條件下實現故障診斷。盡管計算階次跟蹤法擺脫了模擬濾波器等硬件的局限性,節約了成本以及系統復雜性,但是精度嚴重依賴于插值,并且需要轉速計獲取速度信息,導致該方法的應用受到了限制。

在軸承轉速信息不可用時,發展了無轉速計的階次跟蹤法[22],其中通過時頻分析獲得轉軸瞬時頻率的方法應用最為廣泛。Wang等[23]將振動信號經過短時傅里葉變換(STFT),使用基于幅度和的譜峰搜索算法從TFR中提取瞬時故障特征頻率(IFCF),然后基于階次跟蹤思想,對瞬時故障特征頻率(IFCF)進行重采樣,從而將非平穩時域信號轉換為平穩故障相位角度(FPA)域信號,最后將FPA 域信號轉換為故障特征階(FCO)域信號并識別故障類型。Huang等[24]認為基于瞬時故障特征頻率(IFCF)在TFR具有最大幅值的假設在某些條件下并不成立,因此提出了一種可靠的基于快速路徑優化的多重時頻曲線提取算法,該算法直接從TFR中提取包括IFCF、諧波和ISRF在內的T-F曲線,利用平均曲線曲線比來描述提取曲線之間的關系,通過比較平均曲線比與故障特征系數來實現故障診斷。Zhu等[25]認為角域隱藏著周期特征,為此將短時傅里葉變換與譜相關分析結合,提出了一種基于Teager-Kaiser能量算子(TKEO)和角域快速譜相關(Fast-SC)的故障特征提取方法。Zhao等[26]提出了一種無鍵相階次跟蹤法,首先使用無轉速階次跟蹤(TLOT)方法從時域信號中提取轉速信息,利用包絡階次譜(EOS)在階次域中恢復軸承特征頻率。實驗結果表明,通過結合TLOT 和EOS 的優點,該方法能夠在變速條件下有效、準確地識別出不同的軸承故障。Urbanek等[27]提出基于相位解調原理和聯合時頻分析的兩步法精確提取轉速曲線并完成模擬振動信號和實測振動信號的分析。Feng 等[28]開發了一種基于譜峭度的聯合時變幅頻解調頻譜以及改變故障特征頻率提取的時頻集中方法。

無轉速計的階次跟蹤法完全擺脫了硬件的限制,在一定程度上成功解決了速度信息不可用的問題,但是大部分研究都是基于短時傅里葉變換來估計轉速信息,然而,STFT是一種線性的固定分辨率的時頻分析方法,同時受Heisenberg-Gabor 不等式的約束,而變工況條件下軸承振動信號具有強時變、強非線性等特點,因此短時傅里葉變換無法從非穩定信號中精確估計瞬時轉頻,致使基于無轉速階次跟蹤法的故障診斷準確率低于基于硬件階次跟蹤法和計算階次跟蹤法。

1.2 時頻分析

盡管階次跟蹤法在一定程度上解決了變轉速問題,但是容易擾亂原始信號的脈沖響應,導致包絡畸變等問題,因此不能很好保證采樣后信號的穩態特征,國內外學者通過時頻分析來解決上述問題。

1.2.1 小波變換

由于快速傅里葉變換的窗口大小固定,無法針對頻率自適應。小波變換使用面積固定而形狀可變的窗函數通過多辨率分析能很好地平衡時間分辨率和頻率分辨率,其本質為通過小波基函數對信號進行濾波和加權,小波基的選擇是影響小波分析效果的關鍵。

(1)連續小波變換(CWT)

連續小波變換的數學原理為將振動信號與系列復共軛小波的卷積運算來實現,其數學表達式為:

其中,ψ*表示小波基函數ψ的復共軛;h表示尺度參數或膨脹參數;j表示平移參數。

王曉龍等[29]為解決軸承失效早期,特征信號微弱,并且受傳遞路徑衰減及環境噪聲影響等問題,使用不同尺度的小波系數重構信號,從而得到不同尺寸的信號分量,依據峭度對信號分量進行合并,并利用相關分析進行特征篩選,保留峭度值最大的分量。楊蕊等[30]認為峭度忽略了信號的周期特性,因此使用相關峭度來分析信號在各個周期間的相關程度,從而確定最優小波尺度。陳科百[31]認為原始信號沖擊成分能量決定連續小波分解系數的大小,信號與小波基越相似,分解出的系數值就越大,因此將信號進行連續小波變換后做小波尺度-能量譜圖,并對小波尺度-能量譜圖中的特殊尺度下的小波系數進行快速傅里葉變換,取得小波系數的頻譜,最后通過分析小波系數的頻譜從而提取特征頻率。

小波分析是通過不同尺度伸展的柔性窗口和原信號進行對比獲得小波系數,因而小波變換系數包含各種尺度下信號的信息,解決了傅里葉變換固定窗口的問題,因此可以處理突變和非穩定信號。然而在計算小波系數時,小波在每個時間尺度域上都是平滑移動的,因此連續小波變換具有較高的計算成本。

(2)離散小波變換(DWT)

對連續小波進行二階離散化,將公式(1)中的伸縮因子“h”替換為“2j”,平移參數“j”變為“2jk”,就可以得到離散小波變換:

離散小波變換可以將信號分解為高頻通道和低頻通道,分解后的高頻信號還可以繼續被分為高高頻和低高頻,因此離散小波變換可以有效減少冗余特征。

Li等[32]將振動信號經過離散小波變換后,計算各尺度信號的方差,最后根據對數方差的斜率估計基于小波的多尺度斜率特征,實現了軸承缺陷的識別以及對齒輪磨損的診斷。

小波分析的最大特點就是所使用的小波基函數是通過小波函數經過平移、伸縮變換得到的,具有多分辨分析的功能[33]。通過這種多分辨率分析,在高頻信號中獲得一個好的時間分辨率和較差的頻率分辨率,低頻信號中獲得較好的頻率分辨率和較高的時間分辨率。尤其在故障診斷領域,故障信號是瞬態的、非平穩的,小波分析的局部化能力能最大限度去分析信號在時域和頻域上的特征,從而準確刻畫出信號在任意時刻的瞬時頻率,這是傅里葉分析無法完成的。小波分析中不同的小波基函數,使其具有良好的多分辨率分析能力,能使時域分辨率和頻域分辨率達到最佳平衡狀態。而且,小波分析的時頻分辨率可以根據信號自適應調整。對于信號中比較陡峭的成分,可在該區域使用高分辨率的小波來增加頻率分辨能力,對于信號中平滑區域,小波分辨率可適當降低,從而達到時域頻域分辨率的一個平衡狀態[34]。小波1層分解后信號分為低頻信號、高頻信號兩個部分,但是在2 層小波分解時,小波變換只對低頻信號進行分解,高頻信號不做變化,針對以上不足,國內外學者發展了小波包變換。

(3)小波包變換(WPT)

小波包變換是一種經典的多時間分辨率時頻分析算法,信號經過小波包變換后得到低頻帶的近似系數和高頻帶的細節系數,彌補了小波分析只對低頻成分進行分解的缺陷。信號經過小波包分解后,通過對分解后和濾波后的分量的動態特性進行統計學計算來實現故障的分類。信息熵、排列熵、模糊熵等動力學指標是使用的統計學方法。張雄等[35]通過計算樣本小波包各子帶的散布熵值來構建特征矩陣,通過主成分分析進行特征選擇,較少冗余特征,最后采用Meanshift概率密度估計聚類中心位置。Rodriguez 等[36]認為擴展信息熵可以有效提高故障診斷準確率,因此提出一種平穩小波包結合色散熵和置換熵的故障診斷方法。

小波包分解的實質是對信號高頻部分進行小波分解,因此從本質上來說兩者只是在應用方面側重點不同。軸承振動信號在轉速比較高的情況下被載波調至高頻,實驗采集獲得的信號中高頻成分比較多,因此小波包分解比小波分解更適用于處理非平穩的振動信號[37-39]。波包分解對信號進行的分解更徹底,相比于小波分解更全面,對信號高頻部分繼續分解,使高頻分辨率增加。

1.2.2 經驗模態分解及其衍生算法

小波分析是通過不同尺度伸展的柔性窗口和原信號進行對比獲得小波系數,因而小波變換系數包含各種尺度下信號的信息,解決了傅里葉變換固定窗口的問題,因此可以處理突變和非穩定信號[40]。然而在計算小波系數時,小波在每個時間尺度域上都是平滑移動的,因此連續小波變換具有較高的計算成本。在處理非平穩信號時,只要小波基選取得當,小波包分解就會具有強大的多分辨率分析能力,但是如何選取小波基成為限制小波包分解應用的重要因素,經驗模態分解成功解決了該問題。

(1)經驗模態分解(EMD)

經驗模態分解將任意一個頻率成分復雜的信號看做一系列本征模態分量的疊加,信號的頻率信息、時域信息相互耦合,同時各本征模態函數相互獨立,因此對非平穩信號進行經驗模態分解實質上是一個解耦合的過程[41]。假設軸承故障信號為S(t),其定義如下:

其中,s(t)為振動信號,v(t)為噪聲。對用于含有噪聲的故障信號進行經驗模態分解步驟如下:

步驟1 確定故障信號S(t)的所有局部極值,應用樣條插值法得到故障信號的上包絡u(t)和下包絡l(t),通過下式計算上、下包絡的均值m1(t):

步驟2 故障信號S(t)與均值包絡相減,即可得到第一個分量h1(t):

步驟3 若分解后的分量h1(t)局部極值點與過零點數量差不大于1,其極大值與極小值可構成包絡線,且h1(t)的上包絡線與下包絡線均值為0,則h1(t)為本征模函數(IMF),若不符合,則重復步驟1與步驟2,符合上述條件后的h1(t)即為S(t)的第一個本征模函數(IMF)分量,記為c1(t)。分離后的信號如下式計算:

步驟4 將分離后的信號重新定義新的信號,并重復上述三個步驟,即可得到一系列本征模函數(IMF)分量ci(t)和殘余分量。即原始信號S(t)為各本征模函數分量ci(t)與殘差分量之和。

其中,r(t)被稱作余項。

楊建華等[42]應用級聯自適應分段線性隨機共振系統對振動信號進行降噪,然后對降噪后的信號進行經驗模態分解,從而提取故障特征。任學平等[43]將經驗模態分解與AR譜相結合來提取故障特征,該方法將信號進行經驗模態分解獲得固有模態分量后,通過自回歸模型進行特征提取,并由AR模型的參數和殘差的方差組成了故障特征向量矩陣,最后進行分類。

經驗模態分解無需借助小波基函數,只針對信號本身進行分解,消除了小波基的選取對信號分解的影響,并且分解后的本征模態分量都代表著原始信號的頻率成分,殘余分量rn(t) 代表了信號變化的平均趨勢。然而,在對信號進行經驗模態分解時,需確定信號的極值點,通過極值點擬合上、下包絡線,并求得平均包絡線,進而得到本征模態分量,由于噪聲的存在,導致不同類別信號的極值點很容易發生變化,從而引起混疊效應,此外經驗模態分解還存在邊界效應、對噪聲敏感等問題。

(2)集成經驗模態分解(EEMD)

為了解決經驗模態分解引起的混疊效應等缺陷,谷豪等[44]提出一種集合經驗模態分解方法。該方法通過在分解的過程中加入具有零均值統計特性的高斯白噪聲,以多次分解平均的方式來抵消信號的外來噪聲,從而避免異常事件模式與固有模式的混淆,改善了經驗模態分解的IMF分量的模態混淆現象,既保證了分解尺度分布的均勻性,又降低了異常事件的影響,凸顯了振動信號的故障特征。王瀟桐[45]通過對振動信號進行集成經驗模態分解和模糊熵來進行特征提取,并以拼接的方式與時域特征進行特征融合,最后得到一個10 維故障特征作為支持向量機的輸入,進而進行故障診斷。

時頻分析方法通過提取軸承信號中的時頻特征來進行故障診斷,因此信噪比越高,特征提取越容易。然而在故障發生早期,故障信號十分微弱,傳統的時頻分析方法很難從中提取特征,因此,國內外學者采用基于非線性系統的信號處理方法來解決上述問題。

1.3 基于非線性系統的信號處理

1.3.1 隨機振動

1981 年,意大利學者Benzi 等[46]首次提出隨機共振理論,用于研究古額爾德氣象冰川變化,經過幾十年的發展,現已廣泛應用于機械故障診斷、物理以及數學等領域。與傳統提出噪聲方法不同的是,隨機共振盡量保留原始信號,使噪聲與原始信號產生共振,從而將噪聲能量轉移到微弱的信號中,從而完成微弱故障特征提取。傳統的隨機共振僅能對微弱信號(f<<1 Hz)進行特征提取,然而,當前旋轉機械軸承發生故障時,其故障頻率大約在幾十到上千赫茲,某些高速列車的故障頻率甚至達到幾萬赫茲,因此該方法失效。為解決上述問題,發展了幾種微弱高頻特征檢測方法,如頻率轉移變尺度隨機共振[47]、二次采樣隨機共振[48]、多尺度噪聲調節隨機共振[49]以及多尺度雙穩陣列[50]等。影響隨機共振的因素有三個,分別為周期性驅動力、噪聲和隨機共振系統,但是在實際應用中,噪聲和信號本身往往無法控制,因此只能通過優化隨機共振系統來產生最佳的隨機共振。當前隨機共振系統按穩態類型主要分為單穩態、雙穩態與多穩態三種。

基于單穩態隨機共振系統的非線性勢函數,Zhang等[51]基于絕對勢和指數勢,構造了廣義指數型單井勢函數,然后數值分析了Levy 噪聲驅動下相應的指數型單阱自旋共振系統的特性。針對采用粒子群算法調整單穩態隨機共振系統的參數時,全局收斂時間較長問題,Liu等[52]提出了一種基于單穩態隨機共振的粒子群優化算法,將線性慣性加權函數改為余弦函數的非線性函數,以縮短全局收斂時間。

基于雙穩態隨機共振系統的非線性勢函數,田晶等[53]提出了基于容忍遺傳算法(TAGA)的自適應雙穩態隨機共振(BSR)的中介軸承故障診斷方法,該方法在傳統自適應遺傳算法中引入容忍度思想,建立一種容忍遺傳算法,采用容忍遺傳算法對雙穩態隨機共振系統的結構參數a、b進行優化,建立自適應雙穩態隨機共振系統對故障信號進行處理。尹進田等[54]提出一種基于狀態轉移算法的參數同步優化隨機共振新方法用于解決傳統隨機共振存在參數固定或僅對單個參數進行優化的不足問題。該方法采用移頻變尺度方法將大參數信號變換成小參數信號,通過狀態轉移算法對隨機共振系統的參數進行同步優化,以最大信噪比為優化目標,可同時獲取系統最優參數,從而有效削弱信號中的噪聲和增強微弱特征,從而實現早期故障準確診斷。以上方法在對雙穩隨機共振的系統參數進行調節時,會使勢函數的勢壘高度和勢阱間距同時發生變化。為了能直觀了解勢函數特征的變化對雙穩隨機共振的影響,劉進軍等[55]提出一種基于勢函數特征參數調節的隨機共振方法。該方法通過變量代換法實現勢函數特征參數(勢壘高度參數和勢阱間距參數)的解耦,因此有利于勢函數的調節。以上研究均基于對稱勢誘導隨機共振的研究,但是在非對稱勢誘導隨機共振中,阱深與阱寬非對稱性同時變化,基于此,譙自健等[56]研究阱深和阱寬單一非對稱雙阱勢誘導隨機共振的微弱征兆特征增強機理,提出了非對稱勢誘導隨機共振的機械重復瞬態特征提取方法。

基于多穩態隨機共振系統的非線性勢函數,馬強等[57]通過對傳統的雙穩態隨機共振進行研究改進,提出了一種三穩態隨機共振,并對三穩態隨機共振過阻尼布朗粒子的躍遷現象進行研究,分析不同參數下布朗粒子的函數變化,并采用變時間尺度方法對三穩態隨機共振進行優化,在經傅里葉變換進行時頻分析,由此提取出故障頻率。張剛等[58]將混合雙穩態隨機共振系統與Woods-Saxon模型相結合,提出混合三穩態隨機共振系統,然后利用自適應算法尋找混合型三穩隨機共振系統最優參數組合,對α噪聲背景下的仿真信號進行檢測,證明其優越性。

上述研究均基于單個非線性系統,由非線性系統控制理論可知,多個系統協同作用有可能對系統的輸出產生正向作用。王術光等[59]以級聯隨機共振為研究對象,采用人工魚群算法對雙穩態結構參數進行優化,進而設計隨機共振系統,并利用該系統對滾動軸承故障信號的數據進行處理。Li等[60]提出雙共振方法,在控制系統和被控系統中都輸入周期信號和噪聲,使得控制系統和被控系統均能產生隨機共振,因二者之間相互作用,增強了被控系統的共振效應,實現微弱信號的增強輸出。李偉等[61]將兩個單一周期勢系統線性耦合成一個新的系統,并采用粒子群優化算法自動匹配控制系統的相關參數,通過控制系統的輸出反饋調節被控系統的隨機共振,從而達到增強信號的目的。張剛等[62]通過將Duffing振子和Van der Pol 振子進行線性耦合,并將微分項作為耦合系統的反饋,從而增強耦合強度,形成一種Duffing與Van der Pol強耦合系統。為解決變換尺度的變化范圍的選取缺乏固定標準、參數自適應效率低、檢測到的目標信號不夠明顯等問題。張勇亮等[63]提出一種基于并聯自適應隨機共振的微弱信號檢測方法。該方法綜合考慮大參數信號發生隨機共振時變換尺度與采樣頻率之間的關系確定基于采樣頻率的變換尺度的最大變化范圍,通過尺度分段搜索、采用改進粒子群算法、共振系統并聯等方法綜合提高參數自適應的效率、突出檢測到的目標信號,實現強噪聲背景下微弱信號的快速、有效檢測。

1.3.2 混沌理論

1963年,美國氣象學家洛倫茲在研究天氣預報模型時,發現用計算機求解一個確定的微分方程時,其解是非周期的且具有隨機性,首次在一個確定性系統中從數值上得到混沌解,1975年,學者李天巖和約克首次將這種確定性方程產生隨機現象命名為混沌[64]。近年來,發展了幾種不同的混沌振子,在軸承故障診斷領域,Duffing振子是最常見的振子。

針對微弱信號特征提取困難問題,任學平等[65]改進了Duffing 振子的非線性項,成功將其應用于軸承故障診斷領域。李嶺陽等[66]將混沌振子檢測法應用于滾動軸承外圈、內圈和滾動體故障信號的檢測中,通過輸出相圖的變化來判斷故障信號是否存在,有效地實現了對滾動軸承故障信號的檢測。劉彬等[67]分析了微弱周期信號相位角對檢測系統的影響,提出采用多相位混沌振子陣列來消除微弱周期信號相位角對檢測系統的影響。針對現有混沌振子難以檢測頻率未知微弱信號這一難點,李國正等[68]提出利用Duffing、振子輸出值的方差峰值結合遺傳算法檢測淹沒在強噪聲背景中頻率未知微弱信號的一種新方法,該方法采用具有相位偏移的Duffing振子陣列覆蓋全相位,并結合遺傳算法,優化求解不同頻率輸入信號下系統輸出值方差的極值,以此得到待測信號頻率的方法。以上方法只能同時檢測某一頻率的待測信號,黃繼堯等[69]提出了一種可同時檢測多個未知頻率微弱信號的Duffing振子優化模型。該模型通過傅里葉與反傅里葉變換、頻率截取和尺度變換后,得到相圖在經過判決后輸出結果。

1.4 本章小結

基于信號處理方法的關鍵在于人工提取特征的優劣,且嚴重依賴于插值以及信號分析處理的方法,通過對振動信號進行一系列的處理、變換和分析,將原始振動信號轉化為時域信號、頻域信號以及時頻信號,再將這些人工提取的特征作為輸入,使用隱形馬爾可夫模型、支持向量機等機器學習方法進行故障分類。當前絕大部分的研究都是通過數值計算等手段,將高維數據特征轉化為低維特征,從而減小計算量,并作為模型輸入進行分類,因此該方法的分類性能嚴重依賴于人工特征值提取的好壞,且需要大量專家經驗,導致完全基于信號分析處理的方法應用局限性較大,并逐漸被取代。隨著深度學習的發展,國內外學者將深度學習應用到滾動軸承故障診斷領域。表1 總結了基于階次跟蹤和時頻分析方法的優勢、不足以及適用場景。

表1 基于階次跟蹤和時頻分析方法Table 1 Order tracking and time-frequency based analysis methods

2 深度學習

神經網絡省略了人工提取特征的繁瑣步驟,無需專家經驗,使用原始振動信號作為模型輸入,即可完成故障診斷,是智能提取特征值的代表方法。但是在處理變工況軸承故障診斷問題時,由于傳統神經網絡的網絡結構較為簡單,網絡層數較少,因此能夠處理的非線性運算水平較低,僅能提取淺層特征,為解決上述問題,提出了深度學習概念。與傳統神經網絡不同的是,深度學習通常擁有更為復雜的結構以及更深的網絡層數,底層特征逐層非線性映射到高維空間,特征表達隨著網絡層數的增加逐漸抽象化,最終形成層次化、抽象化的非線性特征[70]。并且網絡層數越深,非線性擬合能力越強,能夠處理的高維數據更加復雜。隨著深度學習的高速發展,目前應用于變工況軸承故障診斷的神經網絡主要有卷積神經網絡(CNN)、深度置信網絡(DBN)和自動編碼器。

2.1 卷積神經網絡

卷積神經網絡于1994 首先被提出,是一種經典的分層神經網絡[71]。因其權值共享等特點,致使其結構簡單,特征提取能力強,因此在軸承故障診斷領域有非常廣泛的應用。由于卷積神經網絡具有分層的特點,且每層網絡提取的輸入為二維數據,其分層結構可以智能學習振動信號的故障特征,是一種高效的端到端的學習系統,Guo等[72]提出一種自適應學習率的分層卷積神經網絡,第一層網絡用于識別故障的類型,第二層網絡用于識別故障的損傷程度,實驗結果表明其變化的學習率可以有效提升模型的學習能力。Xia等[73]認為通過提取單一傳感器數據的特征并不能有效減小不同工況下數據集的分布差異,因此提出一種多傳感器的故障診斷方法。通過對從多個傳感器測量的振動信號來提取不同工況下不變的故障特征,從而提升分類精度。張西寧等[74]認為傳統的卷積神經網絡在提取故障特征時,并不能只提取不同工況下的公共特征,因此提出一種小尺度卷積核以跳動的方式進行降采樣的方法,試圖通過改進傳統卷積神經網絡的模型結構來解決變工況問題。該方法使用滑動步長為2,激活函數為“relu”的卷積核代替傳統的池化層。該方法不僅可以使輸入數據的尺度降低一半,還可以在池化階段進行特征提取,通過訓練卷積核權重來選擇有效的特征。趙小強等[75]使用三個3×3的卷積核以堆疊和串聯的方式加入到殘差塊中,用于代替池化層,并將空洞卷積代替傳統卷積層來擴大感受野,殘差塊以跳躍的方式連接不同層,以此來進行特征融合。Wu等[76]認為卷積核尺寸是影響卷積神經網絡特征提取的關鍵,單一尺寸的卷積核無法同時平衡時間分辨率和頻率分辨率,因此提出一種多尺度卷積神經網絡來提取公共特征。以上方法都是通過人工調節參數的方法進行模型訓練,Wang 等[77]將蟻群算法與卷積神經網絡相結合,通過蟻群算法來優化卷積神經網絡的超參數,從而避免了人工調參的弊端。Ince等[78]提出一個實時監控系統用于約束卷積神經網絡中的權重,該方法不需要任何形式的預處理、特征轉換以及特征提取,省略人工調參和特征提取的繁瑣步驟,直接對原始振動信號進行評估。以上方法通過對原始振動信號進行特征提取進而進行故障診斷,而卷積神經網絡在處理圖像方面具有絕對的優勢,因此潘成龍等[79]通過小波變換,將一維振動信號轉化為二維時頻圖,以此進行故障診斷。陳里里等[80]將一維信號通過連續小波變換轉化為二維時頻圖后,將三種故障特征圖進行特征融合后輸入到卷積神經網絡中進行故障分類。張訓杰等[81]將卷積神經網絡與雙向門控循環單元相結合,由卷積神經網絡提取二維圖像的空間特征,由雙向門控循環單元提取時間特征,進而進行故障分類。

2.2 深度置信網絡

深度置信網絡兼顧生成模型和判別模型的優點,由可視層和隱藏層構成,是一種具有雙層結構深層混合網絡模型[82]。在軸承故障診斷領域,深度置信網絡依靠無監督學習自動提取特征,將原始信號映射到特征空間,并依靠監督學習將特征空間轉化為判別空間,從而進行故障類型判別,完成故障分類任務。

Tao等[83]最早將深度置信網絡應用到滾動軸承故障診斷領域,通過構建編碼器實現故障特征的提取,并將提取的低維故障特征作為模型輸入,最后通過最小化輸入輸出的能量實現故障類別的判別。單外平等[84]認為盡管深度置信網絡在提取故障特征時對輸入數據進行了降維,但是計算量仍然十分巨大,因此在自動提取故障特征后,在訓練判別器的過程中考慮了時間復雜度等因素,優先調節偏導數較大的參數,不僅提升了分類準確率,還大大降低了訓練耗時。以上方法從原始振動信號中提取故障特征,在一定程度上實現了故障類型判別,但是當工況發生變化時,該方法準確率大大降低。王圣杰等[85]認為單一的振動信號不足以表征變工況條件下的非平穩信號,因此提出一種多傳感器特征融合方法。該方法采用多個傳感器在不同位置同時測得軸承振動信號,再對測得的信號進行特征提取。胡永濤[86]認為當工況發生變化時,深度置信網絡自動提取的特征并不能解決訓練集與測試集之間的分布差異,因此提出一種多特征融合模型。該方法首先對原始振動信號進行FCM聚類分析,從而判定該樣本是否存在故障,并生成部分樣本的故障標簽,其次,分別提取振動信號的提取DTCWT 樣本熵MMEMD 能量熵和VMD 樣本熵,通過主成分分析方法,衡量三種特征的重要性,并給定一個權重,進而進行特征融合,得到一個同時包含三種特征的綜合故障特征,最后將綜合特征作為深度置信網絡模型的輸入,進行故障類型判別,實驗結果驗證了該方法的有效性。Chen 等[87]綜合多傳感器融合和多特征融合的優勢,提出一種多特征融合及多傳感器特征融合的DBN故障識別方法,準確率得到了進一步提升。Gan等[88]提出一種分層深度置信網絡方法,該方法首先提取原始振動信號的小波包能量作為模型輸入,第一層深度置信網絡用于識別故障出現位置,同時第一層的輸出作為第二層模型的輸入,用于判定損傷程度,綜合兩層網絡模型的結果,完成十分類任務。為了進一步提升模型分類準確率,深度置信網絡與其他優化算法相結合的方法不斷被提出。Liu等[89]將隨機梯度下降算法應用到深度置信網絡的訓練過程中,同時,為避免某些模型超參數的設置影響模型分類性能,采用粒子群算法優化模型超參數,進一步提升模型準確率。熊景鳴等[90]將深度置信網絡與支持向量機相結合,使用深度置信網絡提取故障特征后送入支持向量機中進行分類,取得了良好的效果。Shao 等[91]認為盡管深度置信網絡避免了人工提取特征等繁瑣的步驟,但是這種方法提取的特征維度過高,導致計算量很大,因此將主成分分析應用到特征提取過程中,對高維特征進行特征選擇后,降低了特征維度,在提升準確率的同時減小了訓練耗時。

2.3 基于自編碼器的軸承故障診斷

自動編碼器是一個輸入輸出相同、具有自我學習能力的神經網絡,由編碼層和解碼層組成。其中編碼層由輸入層和隱藏層構成,解碼層由隱藏層和重構層組成。網絡模型的編碼過程就是將輸入的高維特征經過激活函數的特征映射后,轉化為隱藏層的低維特征,解碼過程就是將隱藏層的低維特征經過激活函數重構為輸出目標。自編碼器的訓練過程,就是通過構建不同的優化函數,實現目標輸出與輸入之間的損失值達到最小。Jia等[92]使用一個五層的自動編碼器優先進行了故障的診斷。Yang等[93]使用將高斯核應用到深度自動編碼器中,實現了故障分類。童靳于等[94]采用最大相關熵代替傳統的均方誤差作為損失函數,同時,為減小重構誤差,加入了稀疏懲罰項和非負的約束因子收縮懲罰項,最后,通過灰狼優化算法對模型中的產參數進行優化,避免了人工調參等繁瑣復雜的任務,從而實現變工況軸承故障診斷。

2.4 本章小結

深度學習方法側重于挖掘數據層的規律,利用自身強大的特征提取能力,挖掘不同工況間相同的特征來解決變工況問題。該方法通過大量的數據訓練模型權重,智能挖掘故障特征,避免了人工提取特征等繁瑣的步驟,是實現變工況軸承故障診斷的有效方法。然而深度學習內部計算機制的可解釋性很差,缺乏理論支撐,導致無法清楚詮釋數據蘊含的領域知識[95],并且當輸入特征維度過高時,深度學習模型的模型計算成本巨大,導致需要花費大量的訓練時間成本,并出現特征同質化等問題,因此限制了深度學習方法的應用,遷移學習方法成功解決了上述問題。表2 總結了各類深度學習方法的特點以及優缺點。

表2 深度學習方法Table 2 Deep learning methods

3 遷移學習

傳統深度學習方法要求訓練集和測試集必須獨立同分布,當工況發生變化時,訓練集與測試集會出現明顯的域位移現象,因此限制了深度學習的發展。遷移學習可以在不同但相關聯的兩個領域內挖掘出公共特征。遷移學習已經成功應用于各個領域,在軸承故障診斷領域被廣泛應用的方法主要有三種,即領域自適應、參數遷移以及特征遷移。

3.1 領域自適應

3.1.1 基于數據分布的領域自適應

基于數據分布的領域自適應假設數據的特征空間與類別空間均一致,而源域及目標域間的條件分布及邊緣分布存在差異。通過減小域間的分布差異來解決變工況問題。根據數據分布,可分為邊緣分布自適應、條件分布自適應和聯合分布自適應三種。

邊緣分布自適應是通過減小源域及目標域間的邊緣概率分布距離來解決變工況問題,該方法最早由Sinno等[96]提出,并且通過MMD(maximum mean discrepancy)來計算域間分布距離。Chen等[97]基于遷移成分分析算法,提出一種變工況故障特征提取方法,在一定程度上解決了變域問題。Lu 等[98]提出一種具有領域自適應能力的深度網絡模型,該網絡模型采用自編碼器來提取故障特征,使用MMD計算源域及目標域間的分布差異并作為優化目標。以上方法僅考慮了對齊全局域,而忽略了類內的相似性,忽略源域及目標域之間的關系會導致丟失每個類別的細粒度信息,條件分布自適應通過對抗的思想來對齊子領域分布,從而解決變工況問題。吳靜然等[99]采用一維卷積神經網絡進行特征提取,并且通過最小化局部最大平均差異和分類器損失函數,進行相關子域的分布對齊。在江南軸承數據集上驗證了該方法的有效性。為了充分利用卷積神經網絡強大的圖像特征提取能力,董紹江等[100]采用連續小波變換,將一維軸承振動信號轉化為二維特征時頻圖,采用完成預訓練的Resnet-50網絡結構來提取源域及目標域共同的特征,使用局部最大均值差異來度量子領域自適應,通過計算目標域偽標簽以匹配條件分布距離來進行子領域自適應,從而減小域間的分布差異。Shen 等[101]認為邊緣分布和條件分布對遷移學習的貢獻不同,邊緣分布對齊全局域,而條件分布對齊子領域,兩種分布的相對重要性很難動態定量評估,因此單一使用某種分布無法滿足要求,針對以上缺陷,提出一種動態聯合分布對齊網絡。該網絡提出動態調整因子來定量比較MMD 和WCMMD 的相對重要性,從而同時對其條件分布和邊緣分布,同時,為加權條件最大差異,財通軟偽標簽來代替模型輸出的標簽來改變條件分布的機算,從而有效減小條件分布的距離。以上方法在一定程度上解決了變工況條件下深層特征的域位移現象,但是楊春柳[102]認為網絡較淺層仍然存在域位移現象,因此提出一種多層域自適應滾動軸承故障診斷方法,該方法將原始振動信號最為輸入,使用一維卷積神經網絡進行特征提取,并提出了多層域自適應和權重正則化項約束CNN 參數,從而減小因變工況而導致的特征分布差異。楊冰如等[103]采用預訓練好的ResNet18 進行特征提取,并將不同殘差塊提取的深、淺層特征計算MK-MMD 距離,以匹配邊緣分布差異,再將不同殘差塊提取的特征輸入到分類器中輸出偽標簽,從而縮小條件分布差異。楊春柳[104]認為當樣本數量過大時,使用MMD計算分布差異的成本巨大,而Wasserstein 距離幾乎在任何位置都是可微的,因此選用Wasserstein距離來度量域差異,以此來減小計算成本。

3.1.2 基于數據選擇的領域自適應

假設域間的類別空間一致,則認為域間的特征空間也一致。通過深度學習以及機器學習的手段,在不同域間提取相同的特征,以相同的特征作為模型的輸入,從而實現變工況軸承故障診斷的方法稱為基于特征選擇的領域自適應。

基于特征選擇的領域自適應的重點在于選擇不同域間相同的特征,Long 等[105]認為在遷移學習中源域數據量大小并不是影響模型準確率的唯一因素,因此對源域進行樣本選擇,同時結合邊緣條件自適應,從而達到變工況軸承故障診斷。Lu等[106]通過共享源域以及目標域特征的方法,實現無監督領域自適應,Tang等[107]提出一種EGR特征選擇方法用于度量特征與類別間的非線性關系,從而解決變域問題,但是該方法為半監督學習,在遷移學習的過程中目標域仍然需要小部分帶有標簽的數據。但是帶標簽的跨域數據集很難獲取,為解決該問題,郭亮等[108]提出一種跨數據集的變工況軸承故障診斷方法。該方法源域數據集為實驗室條件下獲取的帶標簽的數據集,采用帶領域適配正則約束項的一維卷積神經網絡進行特征知識的深度遷移適配,進而訓練模型,最后將在實際工作環境中測得的數據作為目標域進行測試模型,實驗結果驗證了該方法的有效性。Yu等[109]提出一種保留局部流行結構的傳遞分量分析方法,該方法通過故障敏感性以及特征相關性進行特征選擇,進而進行故障分類,通過兩個不同測試平臺收集的振動數據用于驗證該方法,實驗結果證明了該方法的可行性。

3.1.3 基于數據變換的領域自適應

基于特征變換的領域自適應的前提假設為源域及目標域在變換后的子空間中具有相似的分布,根據其變換的形式分為基于統計特征變換的統計特征對齊方法,以及基于流形變換的流形學習方法。

薛輝[110]將子空間分布對齊方法與概率分布自適應相結合,用于減小源域及目標域的分布差異,后又提出一種二階特征對齊的方法,通過學習二階特征變換矩陣,從而使源域以及目標域的歐氏距離最小。基于統計特征變換對齊方法可以在一定程度上解決域位移現象,但是基于統計特征的方法過于依賴專業經驗,Zhang等[111]將故障特征嵌入到流行子空間進行特征變換后,再進行故障分類,從而減小域位移現象,劉海寧等[112]采用一維卷積神經網絡進行特征提取后,將源域及目標域的邊緣分布和條件分布進行對齊,從而進行變工況軸承故障診斷。以上方法雖然在一定程度上解決了特征發散問題,但是其只在子空間進行流形特征學習,并未在特征變換后的子空間進行分布對齊,并且進行子空間學習后,特征發散問題并沒有完全消除,針對以上問題,王肖雨等[113]提出一種基于自適應噪聲完整經驗模態分解(CEEMDAN)與流形嵌入分布對齊的滾動軸承遷移故障診斷方法。該方法首先對原始振動信號進行自適應噪聲的完整經驗模態分解并進行特征選擇,取其峭度最大的前6個干內稟模態分量作為故障特征,將故障特征嵌入流行空間進行特征變換,并將變換后的流行特征進行動態分布對齊,實驗結果表明該方法可以有效解決特征發散問題。然而,童靳于等[114]認為以上方法僅僅通過最小化源域、目標域的分布差異來提升準確率,但是同一故障類別下樣本間的分布距離同樣影響遷移學習的效果,因此提出一種基于類內散度正則化的域自適應故障診斷方法,該方法將源域以及目標域嵌入流行子空間并動態對齊流行特征后,將類內散度引入到模型中,在減小源域以及目標域的分布差異后增大同一類別的聚集程度,使得模型分類準確率得到進一步提升。

3.2 參數遷移

參數遷移側重于對離散模型權重的調整,即使用源域對模型進行預訓練,再使用小部分帶有標簽的目標域數據集對模型參數進行微調,從而達到變工況條件下軸承的故障診斷。

Zhang 等[115]使用源域對模型進行訓練后,凍結除softmax 外的層,并通過小部分帶標簽的目標域數據微調softmax層參數,從而達到變工況故障診斷問題,Shao等[116]采用與訓練的CNN 網絡模型代替隨機初始化參數,僅改變最后全連接層節點的個數及權重來適應目標域的分類任務。Zhao等[117]認為不同尺度的卷積核可以提取不同的特征,在以往研究中特征提取的卷積核恒定,因此提出一種多尺度卷積神經網絡遷移學習框架,使用MS1、MS2 以及MS3 三個尺度卷積核進行特征提取,再權重遷移的過程中,三個卷積和權重保持不變,只改變全連接層的權重,從而實現變工況軸承故障診斷。

由于權重遷移在源域對模型進行預訓練后,需要部分有標記的目標域對模型權重進行微調,因此屬于半監督遷移學習,當無法獲得帶有標簽的目標域數據時,該方法將失效。

3.3 特征遷移

基于特征的遷移學習的主要任務在于尋找連接源域以及目標域的遷移函數來進行特征映射,通過學習特征映射來提取源域和目標域中可遷移的特征[118]。特征函數的作用是提取源域和目標域中不因轉速負載等運行環境發生變化而變化的公共特征,從而實現變工況軸承故障診斷,因此,特征函數的選擇是特征遷移的關鍵。

Wang等[119]提出一個可以將高維振動信號映射到低維潛在空間的主特征函數,在保持數據特征的同時有效減小了域位移,Han 等[120]結合神經網絡與聯合分布算法,提出一種具有搜索和標記功能的識別框架,該框架識別源域后,可以自適應目標域的條件分布,從而實現變工況故障診斷。Wen 等[121]提出一種三層稀疏自動編碼器,并以最大均值差異作為損失值來提取源域和目標域公共的特征。

3.4 本章小結

遷移學習將變工況問題轉化為域位移問題,通過特征提取、特征變換以及特征選擇等手段,成功解決了變工況導致的數據分布差異問題[122],在變工況軸承故障診斷領域取得了巨大的成功。然而,當前遷移學習仍然需要大量數據去訓練模型,當工況差異較大并且含有噪聲時,遷移學習的性能會顯著下降。表3總結了遷移學習方法的特點以及優缺點

表3 遷移學習方法Table 3 Transfer learning methods

4 結論與展望

4.1 結論

隨著社會的發展,旋轉機械越來越智能化、信息化,對旋轉機械的運行狀態在線監測具有重要意義,在此背景下軸承故障診斷技術得到了長足的發展。目前變工況故障診斷領域主要發展了信號分析處理、深度學習和遷移學習三種方法,本文綜述這些方法近五年的進展,主要得出以下四條結論:

(1)基于信號分析的階次跟蹤法將等時采樣轉化為等角度采樣,以此解決變速問題,因此計算量很大,且嚴重依賴于插值精度。并且該方法僅能解決變轉速問題,而無法解決變負載問題。

(2)基于信號分析的時頻分析方法通過將始于信號轉化為時頻域信號,人工提取特征值,以此挖掘不同工況中不變的特征。特征值的提取直接影響著故障診斷效果,因此十分依賴專家經驗,導致逐漸被取代。

(3)深度學習省略手動提取特征等繁瑣的工作,依靠自身強大的學習能力自動提取故障特征,實現變工況軸承故障診斷。但是深度學習方法需要大量帶標簽的數據訓練模型,數據量越大,模型分類性能越好。當不同工況間的差異較大或者工況過于復雜時,模型性能急劇下降。

(4)遷移學習方法不再受訓練集和測試集具有相同分布的約束,能夠在不同但是相關的兩個領域間挖掘出公共特征,能在一定程度上加快訓練過程并提高分類準確率,同時減少為測試集打標簽等繁瑣步驟,有效解決變工況問題。但是在實際工業現場環境中,由于外界的干擾等因素,故障診斷識別率均有不同程度的下降,實際工業場景下復雜工況以及實時數據在線診斷仍然還有很多問題需要研究和解決。

4.2 展望

盡管深度學習和遷移學習盡管在變工況軸承故障診斷領域取得了矚目的進展,但是仍然有需要改進或優化的空間,后續工作應從以下幾方面著手:

(1)多傳感器特征融合

信號處理方法均通過分析單個傳感器測得軸承振動信號來提取特征,因此存在諸多問題。采用多傳感器在多測點采集軸承運行數據可獲得更加豐富的信息,例如信號的物理屬性、空間屬性以及時間屬性,從而增加故障診斷的準確率與魯棒性,因此,基于多傳感器特征融合技術是未來故障診斷領域的發展趨勢。

(2)不均衡數據集

以上研究均基于均衡數據集的故障分類,但是在實際生產生活中,在絕大多數時間內軸承均處于正常運行狀態,一旦軸承發生故障,機械設備將停機檢修,導致可以獲取大量正常數據,而有故障的數據十分稀缺。基于數據驅動的深度學習方法可以取得良好性能的前提是需要大量帶標簽的數據訓練模型,直接將類別不均衡的樣本使用深度學習模型訓練時,會導致模型學習到的特征更加偏向于正常狀態的特征,但很難學習到故障類別的特征,導致出現嚴重的過擬合現象,因此不均衡樣本分類問題是未來的研究內容之一。

(3)多故障分類

目前幾乎所有的研究都是針對軸承的單故障分類,然而在實際工程中,軸承發生故障時可能出現多故障共存的情況,因此,軸承多故障分類是一個值得研究的領域。

(4)數字孿生技術

當前軸承故障均為實驗室環境下人工制造故障,因此不具備實際工程應用價值,而數字孿生技術在能源互聯網行業、電力行業等均取得了初步成果。因此通過數字孿生技術采集真實場景下軸承運行數據,進而進行軸承運行狀態監測及故障診斷是未來的發展趨勢。

猜你喜歡
故障診斷振動特征
振動的思考
科學大眾(2023年17期)2023-10-26 07:39:14
振動與頻率
天天愛科學(2020年6期)2020-09-10 07:22:44
如何表達“特征”
不忠誠的四個特征
當代陜西(2019年10期)2019-06-03 10:12:04
中立型Emden-Fowler微分方程的振動性
抓住特征巧觀察
因果圖定性分析法及其在故障診斷中的應用
UF6振動激發態分子的振動-振動馳豫
計算物理(2014年2期)2014-03-11 17:01:44
基于LCD和排列熵的滾動軸承故障診斷
基于WPD-HHT的滾動軸承故障診斷
機械與電子(2014年1期)2014-02-28 02:07:31
主站蜘蛛池模板: 国产乱子伦视频在线播放| 欧美色综合网站| 中文字幕欧美日韩| 亚洲高清资源| 中文纯内无码H| 欧美一级夜夜爽www| 全部免费毛片免费播放| 国产毛片基地| 在线免费无码视频| 毛片手机在线看| 99视频在线观看免费| 日韩黄色在线| 草逼视频国产| 最新国产在线| 国产综合色在线视频播放线视| 亚洲精品无码AV电影在线播放| 午夜视频在线观看免费网站| 本亚洲精品网站| 欧美国产在线一区| 喷潮白浆直流在线播放| 亚洲无码视频图片| 中文字幕永久在线观看| 亚洲国产成人麻豆精品| 一本一本大道香蕉久在线播放| 婷婷综合亚洲| 国产一区在线视频观看| 国产91成人| 午夜国产精品视频| 99热这里只有精品在线播放| 久久久波多野结衣av一区二区| 欧美97色| 国产精品久久久久久久久久98| 国产在线观看人成激情视频| 亚洲 欧美 中文 AⅤ在线视频| www.亚洲一区| 亚洲成年人片| 中国特黄美女一级视频| 亚洲欧美天堂网| 国产91av在线| 99久视频| 国产精品成人免费视频99| 黄色a一级视频| 亚洲码在线中文在线观看| 国产精品99一区不卡| 九九久久99精品| 亚洲国产理论片在线播放| 青青草原国产| 午夜在线不卡| 欧美亚洲日韩中文| 26uuu国产精品视频| 亚洲美女久久| 欧美精品在线视频观看| 国产一级妓女av网站| 免费看黄片一区二区三区| 91精品亚洲| 久久99蜜桃精品久久久久小说| 亚洲天堂在线视频| 久久动漫精品| 久久亚洲天堂| 91午夜福利在线观看| 国产18在线播放| 欧美亚洲国产精品第一页| 欧美中文字幕在线视频| 69综合网| 无码内射中文字幕岛国片| 久久窝窝国产精品午夜看片| 国产成年无码AⅤ片在线| 性欧美久久| 成人在线欧美| 婷婷亚洲最大| 国产精品美女自慰喷水| 久热中文字幕在线观看| 91精品久久久久久无码人妻| 国产精品毛片一区| 国产新AV天堂| 国产精品流白浆在线观看| 国产美女91视频| 呦视频在线一区二区三区| 欧美激情视频一区| 国模私拍一区二区| 久久久久中文字幕精品视频| 日本不卡在线播放|