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

基于分形理論的渦輪葉片特征提取

2015-04-28 03:34:52馮馳胡楊王兆豐
應(yīng)用科技 2015年4期
關(guān)鍵詞:特征故障

馮馳,胡楊,王兆豐

1.哈爾濱工程大學(xué)信息與通信工程學(xué)院,黑龍江哈爾濱150001 2.西安航空發(fā)動機(jī)(集團(tuán))有限公司,陜西西安710021

基于分形理論的渦輪葉片特征提取

馮馳1,胡楊1,王兆豐2

1.哈爾濱工程大學(xué)信息與通信工程學(xué)院,黑龍江哈爾濱150001 2.西安航空發(fā)動機(jī)(集團(tuán))有限公司,陜西西安710021

采用了葉片溫度這一葉片質(zhì)量的重要指標(biāo),對數(shù)據(jù)預(yù)處理并進(jìn)行特征提取,為渦輪葉片建立起特征模型。基于分形理論提取葉片溫度信號的3種分形維數(shù)特征,結(jié)合K-means聚類分析和ReliefF算法計(jì)算各特征值的權(quán)重,從而建立起渦輪葉片的溫度特征模型,實(shí)現(xiàn)對渦輪葉片故障的早期預(yù)警。統(tǒng)計(jì)結(jié)果表明,該特征模型能夠較好地反映出處于故障狀態(tài)的渦輪葉片的狀態(tài)。

渦輪葉片;特征提取;分形理論; K-means; ReliefF

網(wǎng)絡(luò)出版地址: http://www.cnki.net/kcms/detail/23.1191.U.20150727.0857.001.html

燃?xì)廨啓C(jī)作為國家高技術(shù)水平和科技實(shí)力的重要標(biāo)志之一被廣泛應(yīng)用于航空國防和工業(yè)發(fā)電領(lǐng)域。在它的工作過程中,發(fā)動機(jī)噴出高速噴氣流,同時帶動壓氣機(jī)和渦輪繼續(xù)旋轉(zhuǎn),渦輪葉片長期工作在高溫高壓的惡劣環(huán)境下,其自身溫度過熱使壽命受到極大影響,因此對渦輪葉片工作狀態(tài)的監(jiān)測具有十分重要的意義。分形理論最基本特點(diǎn)是用分形的數(shù)學(xué)工具來描述研究客觀事物。它跳出了一維的線、二維的面、三維的立體乃至四維時空的傳統(tǒng)束縛,更加趨近復(fù)雜系統(tǒng)的真實(shí)屬性與狀態(tài)的描述,更加符合客觀事物的多樣性與復(fù)雜性[1]。綜上所述,文中將分形理論應(yīng)用于渦輪葉片模型的建立,作為葉片故障診斷的理論依據(jù)。文中對渦輪葉片旋轉(zhuǎn)一個周期的溫度信號進(jìn)行預(yù)處理,以得到單個葉片的溫度信號,提取3種分形維數(shù)組成的渦輪葉片特征向量,同時引入K-means聚類分析方法對葉片進(jìn)行分類,根據(jù)分類結(jié)果應(yīng)用ReliefF算法計(jì)算各分形維數(shù)所占權(quán)重,繼而形成每個葉片的特征模型。渦輪葉片輻射測溫系統(tǒng)的輸出信號為電壓信號,電壓經(jīng)擬合后變換為溫度,為方便起見,文中均直接采用電壓表示葉片溫度,并未進(jìn)行電壓溫度轉(zhuǎn)換。

1 單葉片溫度數(shù)據(jù)分割

文中所采集的每組數(shù)據(jù)包括多個渦輪旋轉(zhuǎn)周期,短時間內(nèi)(100 ms),每個葉片的溫度分布幾乎不變,因此可將多個周期的數(shù)據(jù)對齊取平均以減少噪聲帶來的影響。然而因數(shù)據(jù)的每個渦輪周期之間不可避免地存在相位偏差,為提高對齊精度,需先對數(shù)據(jù)進(jìn)行10倍插值。文中采用的三次樣條插值法具有良好的收斂性與穩(wěn)定性,又有二階光滑性。使用三次樣條法插值前后溫度數(shù)據(jù)局部對比如圖1所示。從圖中可以看出,插值后數(shù)據(jù)沒有出現(xiàn)畸變,而是變得更加平滑。

圖1 插值前后波形對比

以每周期溫度數(shù)據(jù)的最大值為基準(zhǔn)分割出17個周期,計(jì)算余下每周期與第一周期數(shù)據(jù)的相關(guān)性,并按照相關(guān)性最強(qiáng)點(diǎn)調(diào)整對齊。圖2為調(diào)整前后對比圖。

圖2 對齊前后波形對比

對齊之后取平均,對這一周期的平均溫度信號按極大值分割法(以每個葉片波形的最大的極大值點(diǎn)為基準(zhǔn)逐一分割)分出每個渦輪葉片的波形,如圖3所示。

圖3 葉片波形分割圖

文中算法流程如圖4所示。

圖4 特征提取算法流程

2 特征值計(jì)算

2.1分形維數(shù)

分形維數(shù)[2]是分形意義上由標(biāo)度關(guān)系得出的一個定量數(shù)值,可以定量地描述分形集的不規(guī)則程度和復(fù)雜程度。

2.1.1分形的標(biāo)度不變性

使用分形理論的前提是需分析的信號具有分形性質(zhì)。事實(shí)上,滿足一定條件的動力學(xué)系統(tǒng)都會產(chǎn)生分形;反之,判別任意形體是否分形,沒有必要從這個形體的產(chǎn)生機(jī)制上入手,只要判斷它是否具有標(biāo)度不變性就可以了。

所謂標(biāo)度不變性是指無論測量尺度如何改變,所測量對象的特性(如形態(tài)特性、復(fù)雜程度、不規(guī)則性、統(tǒng)計(jì)特性等)均不發(fā)生變化。除了嚴(yán)格的數(shù)學(xué)模型以外,如Koch曲線,對于實(shí)際的分形集來說,這種標(biāo)度不變性只在一定的范圍內(nèi)適用。通常把標(biāo)度不變性適用的空間稱之為該分形體的無標(biāo)度區(qū)間。

設(shè)平面R2內(nèi)有圖形F,在平面內(nèi)作間距為δ的方格網(wǎng),則F與方格網(wǎng)相交的方格數(shù)Nδ(F )稱為圖形F在標(biāo)度(分辨率)δ下的盒數(shù)。如果存在一個標(biāo)度范圍(δ1,δ2),在此范圍內(nèi)log Nδ(F )與log δ保持大約恒定的斜率,則這個范圍就稱為無標(biāo)度區(qū)。在無標(biāo)度區(qū)內(nèi),可以將該形體視為分形。即雙對數(shù)關(guān)系曲線中存在一段直線(或近似直線),則該段就是無標(biāo)度區(qū)。在該段直線的標(biāo)度范圍(δ1,δ2)內(nèi),可以認(rèn)為波形F是分形,且這段直線的平均斜率就是分形的盒維數(shù),即

2.1.2分形特性驗(yàn)證實(shí)驗(yàn)

基于上述理論,文中將標(biāo)度定位于(1,80)實(shí)驗(yàn)驗(yàn)證渦輪葉片溫度曲線是否具有分形特性。從實(shí)驗(yàn)所得的log Nδ(F )-log δ雙對數(shù)曲線來看不易準(zhǔn)確判斷定位是否為直線段,為此作出斜率與標(biāo)度的關(guān)系曲線,如圖5所示。

圖5 logNδ(F )-logδ雙對數(shù)曲線及曲線斜率

圖中可看出,δ在(1,15)的范圍內(nèi),斜率值基本保持不變,故可認(rèn)為無標(biāo)度區(qū)為(1,15),此處無需精確無標(biāo)度區(qū)范圍,只需證明其存在無標(biāo)度區(qū)即可。綜上,渦輪葉片的溫度信號具有分形特性,事實(shí)上,根據(jù)分形的定義,文中提取的葉片溫度信號的分形特征所反映出的是葉片在分形算法設(shè)定區(qū)域內(nèi)的波形分布特征。

2.1.3分形特征值計(jì)算

在設(shè)備的故障診斷方面,單重分形維數(shù)(例如盒維數(shù))只能從單一測度出發(fā)描述設(shè)備狀態(tài)信號的分形特征,并不能全面地反映其完整特性。因此文中選用多重分形,計(jì)算渦輪葉片溫度信號的分形特征,取權(quán)重因子q=0、1、2,所對應(yīng)的分形維數(shù)分別為容量維數(shù)(盒維數(shù))、信息維數(shù)、關(guān)聯(lián)維數(shù)[3]。部分葉片的3種分形維數(shù)計(jì)算結(jié)果如表1所示。

表1 第1~5號渦輪葉片的3種分形維數(shù)

3 分配特征權(quán)重

當(dāng)一個特征能夠清晰刻畫形體的時候,類似形體的特征應(yīng)是相近的,而非類似的形體特征相距較遠(yuǎn)。因此,具有相似特性的葉片的特征向量之間距離較近,反之較遠(yuǎn)。而好的特征應(yīng)該使屬于同一類的葉片樣本接近,使不同類的葉片樣本之間遠(yuǎn)離。鑒于這種思想,文中使用最常用的K-means聚類算法將葉片分為N類,根據(jù)結(jié)果使用ReliefF算法計(jì)算各分形維數(shù)所占權(quán)值。

3.1改進(jìn)的K-means聚類算法

聚類算法是給予數(shù)據(jù)自然上的相似劃分,要求得到的結(jié)果是每個類別內(nèi)部數(shù)據(jù)盡可能的相似而類別之間要盡可能大的存在差異。K-means算法是常用的基于劃分的聚類算法,首先隨機(jī)選擇k個對象作為初始的k個簇的質(zhì)心,然后將其余對象根據(jù)其與各個簇的質(zhì)心的距離分配到最近的簇;最后重新計(jì)算各個簇的質(zhì)心,不斷重復(fù)此過程,直到目標(biāo)函數(shù)最小為止[4]。

3.1.1改進(jìn)方法

K-means算法存在以下缺陷[5]:

1)聚類個數(shù)k需要預(yù)先給定,而k值選定難估計(jì);

2)對初始聚類中心的選取具有依賴性,算法常常因此陷入局部最優(yōu)解。

基于上述缺陷,文中采用了改進(jìn)的K-means算法:

1)最佳聚類數(shù)的選擇。根據(jù)DB Index準(zhǔn)則,DBk值越小,說明聚類的效果越好,一般情況下,最佳聚類個數(shù)不會超過個,因而迭代算法可以在2~之間進(jìn)行,在文中選擇2~10,選擇DB指標(biāo)最小的聚類數(shù)為最佳聚類數(shù)。

2)初始聚類中心的選擇。普通K-means聚類方法使用隨機(jī)選擇的聚類中心,這使得算法極易陷入局部最優(yōu)解。文中通過選舉機(jī)制產(chǎn)生初始聚類中心。具體做法是首先計(jì)算最大最小樣本點(diǎn)之間的距離d。指定常整數(shù)p,k≤p≤N。在N個模式中,找到p個相距較遠(yuǎn)的樣本點(diǎn)。再對p個樣本點(diǎn)計(jì)算以本身為球心、以d(指定的常數(shù))為半徑的球內(nèi)所包含的樣本點(diǎn)個數(shù)。選舉個數(shù)較多的前k個樣本點(diǎn)作為初始聚類中心。

3.1.2聚類效果演示

經(jīng)計(jì)算,根據(jù)DB Index準(zhǔn)則,計(jì)算得最佳聚類數(shù)為9,圖6所示的聚類初始中心點(diǎn)與更新中心點(diǎn)對比圖,可見改進(jìn)的K-means聚類方法沒有使聚類結(jié)果陷入局部最優(yōu)解。

圖6 聚類初始中心點(diǎn)與更新中心點(diǎn)圖

最后的聚類結(jié)果分類如圖7所示,表2記錄了各分類所包含的渦輪葉片號,基于分類結(jié)果,文中采用ReliefF算法計(jì)算各分形維數(shù)對分類所做的貢獻(xiàn)值,以貢獻(xiàn)值大小作為該種分形維數(shù)的權(quán)重。

圖7 聚類結(jié)果分類

表2 各類別所包含的葉片代號(共86個葉片)

3.2基于聚類結(jié)果的ReliefF算法

ReliefF算法從訓(xùn)練集D中隨機(jī)選擇一個樣本R,然后從和R同類的樣本中尋找最近鄰的k個樣本Hii=1,2,…,k (),稱為Near Hits,從和R不同類的樣本中尋找最近鄰的k個樣本Mi(C ) i=1,2,…,k (),稱為Near Misses,然后根據(jù)式(1)更新某個特征A的權(quán)重W(A)。

式中: m為重復(fù)次數(shù),diff(A,R,S )指的是樣本R和S在特征A上的差,p(C )指樣本落在集合C中的概率[6-8]。特征的權(quán)重越大,表示該特征的分類能力越強(qiáng),反之,表示該特征分類能力越弱。經(jīng)計(jì)算,葉片的3種分形維數(shù)所對應(yīng)的權(quán)值如表3所示。

表3 3種分形維數(shù)所對應(yīng)的權(quán)值

將葉片各特征乘以其所對應(yīng)的權(quán)值后,再按改進(jìn)的K-means算法計(jì)算分類后的DB指標(biāo),指標(biāo)從0.286 8降到0.155 7,可見ReliefF算法能有效地提高聚類性能,并且能定量分析各維特征對分類的貢獻(xiàn)程度,使葉片溫度信號特征更具有代表性。

3.3模擬故障信號及其識別實(shí)驗(yàn)

建立葉片特征模型的最終目的是為了能夠在葉片發(fā)生故障的初期進(jìn)行預(yù)警,提示操作人員停機(jī)檢查。文中主要研究葉片故障中2種情況:一類是葉片氣道大面積堵塞導(dǎo)致整體溫度過高。整體溫度過高又分為3種情況:前后氣道均堵塞,中部正常導(dǎo)致的前緣溫度過高,尾緣溫度也過高;葉片前緣氣道堵塞導(dǎo)致前緣溫度過高,使得尾緣氣流通過量增大導(dǎo)致尾緣溫度過低;前緣溫度過低,尾緣溫度過高。另一類是葉片局部毛細(xì)氣道堵塞或者隔熱涂層小面積脫落導(dǎo)致局部小范圍溫度過高。隔熱涂層的發(fā)射率大概在0.7左右,而金屬發(fā)射率在0.9左右,因此若葉片的隔熱涂層脫落會引起葉片局部的溫度偏高,通過斯特藩—玻爾茲曼定律,葉片工作溫度在650~750℃的時候,隔熱涂層脫落會使得葉片溫度上升40℃左右,轉(zhuǎn)換成電壓即上升0.14 V左右,為模擬涂層脫落故障的漸變過程,文中僅將電壓調(diào)整為上升0.8 V。整體溫度偏高的故障較為少見,常見故障以局部故障為主[9-10]。

文中使用數(shù)學(xué)軟件MATLAB R2014a版本進(jìn)行模擬計(jì)算,MATLAB具有高效的數(shù)值計(jì)算及符號計(jì)算功能,且提供了大量方便實(shí)用的處理工具箱。隨機(jī)選取正常的第11號葉片溫度數(shù)據(jù),模擬6個模擬故障信號,其中1號和2號故障數(shù)據(jù)是通過將正常數(shù)據(jù)乘以Kaiser窗獲得的,其中一個Kaiser窗的衰減參數(shù)設(shè)為0.3,另一個設(shè)為0.4,來模擬葉片前緣、尾緣溫度都過高的情況; 3、4號故障數(shù)據(jù)是通過將正常葉片數(shù)據(jù)分別通過幅度為0.01和-0.01的正弦窗獲得的,用來模擬葉片前緣溫度高、尾緣溫度低和前緣溫度低、尾緣溫度高的情況; 5、6號故障數(shù)據(jù)是將正常數(shù)據(jù)加上一個以第50數(shù)據(jù)點(diǎn)和第100數(shù)據(jù)點(diǎn)為中心,寬度為50個數(shù)據(jù)點(diǎn),高度為0.8的高斯窗,用來模擬第50數(shù)據(jù)點(diǎn)和第100數(shù)據(jù)點(diǎn)發(fā)生局部故障的情況。1、3、5、6號故障數(shù)據(jù)模擬前后對比如圖8所示,圖中實(shí)線為正常數(shù)據(jù)波形,而虛線則為故障模擬數(shù)據(jù)波形。

表4所示的是11號葉片溫度信號在未進(jìn)行求和取平均時在前10個周期的波形所對應(yīng)的乘權(quán)值后的分形特征值及因這10個周期11號葉片波形變化導(dǎo)致的分形特征值變化的范圍。由表4可以看出,3種分形維數(shù)的變化范圍均不超過2%,故將σ=-2%+ 2% []設(shè)為置信區(qū)間,用來判斷葉片是否處于正常狀態(tài),此置信區(qū)間是根據(jù)文中所用數(shù)據(jù)設(shè)置的故障識別門限,因文中的數(shù)據(jù)量有限,置信區(qū)間與實(shí)際可能存在出入,具體應(yīng)用時可根據(jù)大量數(shù)據(jù)樣本對故障識別門限進(jìn)行修改。

表4 11號葉片10個周期的分形特征及其波動范圍

與11號葉片類似,文中使用同樣的方式模擬出其他85個葉片對應(yīng)的6種故障數(shù)據(jù),即每種故障均有86個模擬樣本。將各故障類型的分形特征與其原始正常波形的分形特征相比較,統(tǒng)計(jì)故障分形特征落在置信區(qū)間之外不同范圍的葉片個數(shù),并計(jì)算識別率組成表5。設(shè)需判別的某模擬葉片的某分形特征為Dxi,其所對應(yīng)的正常葉片數(shù)據(jù)xi同一分形特征為dxi,則識別率η計(jì)算公式如下:

表5 各故障類型分形特征與正常相比偏離超過2%的葉片個數(shù)統(tǒng)計(jì)表

由表可以看出,模擬的各類故障波形與正常波形相比,均能在分形特征上體現(xiàn)出顯著差異,且在關(guān)聯(lián)維數(shù)上體現(xiàn)得尤其明顯,信息維數(shù)其次,容量維數(shù)最次,且當(dāng)增大故障變化的幅度時,識別率能夠達(dá)到更高。因此,分形特征能夠幫助監(jiān)測人員識別故障葉片以及預(yù)測發(fā)生何種故障。對于不同種故障,3種分形特征會體現(xiàn)出不同的識別情況,因模擬的故障類型有限,不能對各分形維數(shù)最適于識別何種故障類型進(jìn)行一一討論。值得一提的是,若使用將原正常波形平移或壓縮伸展的方式來模擬故障信號的話,提取出的故障信號的分形特征與正常信號相比,將不會發(fā)生變化,這是因?yàn)橛?jì)算分形維數(shù)時,算法會先將波形伸展壓縮在一個等邊長的正方形區(qū)域里,這導(dǎo)致無論將原始波形做任何伸展平移,經(jīng)計(jì)算后都會產(chǎn)生同樣的結(jié)果,事實(shí)上,葉片發(fā)生故障時不會僅僅是正常波形的單純平移伸展,所以這種情況可以忽略。綜上,文中所研究的渦輪葉片分形特征能顯著地反映出葉片故障情況,可以較好地應(yīng)用于葉片的故障識別中。但故障識別門限以及故障識別類型需大量的具有代表性的正常葉片波形和故障葉片波形樣本,因研究條件有限故文中不做深入討論。

4 結(jié)束語

針對渦輪葉片溫度信號進(jìn)行了插值、對齊、分割等數(shù)據(jù)預(yù)處理,通過分形理論提取了溫度數(shù)據(jù)的3種分形特征,并引用經(jīng)改善過的K-means算法結(jié)果和應(yīng)用ReliefF算法,使得葉片特征按重要性進(jìn)行權(quán)值分配,最后進(jìn)行了常見葉片故障波形仿真。經(jīng)統(tǒng)計(jì)結(jié)果顯示,文中提取的特征值能夠有效地反映出葉片故障情況,可以應(yīng)用于葉片的故障識別中。但因本課題研究條件有限,不能深入探究故障識別的具體門限值,接下來的研究重點(diǎn)將放在基于大量具有代表性葉片波形樣本的故障識別方向上。

[1]趙玉春.基于混沌分形與模糊聚類的滾動軸承故障診斷[D].秦皇島:燕山大學(xué),2011: 36-42.

[2]郝研.分形維數(shù)特性分析及故障診斷分形方法研究[D].天津:天津大學(xué),2012: 35-51.

[3]李兆飛.振動故障分形特征提取及診斷方法研究[D].重慶:重慶大學(xué),2013: 36-48.

[4]王千,王成,馮振元,等.K-means聚類算法研究綜述[J].電子設(shè)計(jì)工程,2012(7) : 21-24.

[5]袁方,周志勇,宋鑫.初始聚類中心優(yōu)化的K-means算法[J].計(jì)算機(jī)工程,2007,33(3) : 65-66.

[6]張麗新,王家,趙雁南,等.基于Relief的組合式特征選擇[J].復(fù)旦學(xué)報(bào),2004,43(5) : 893-897.

[7]張勇.基于ReliefF算法的模糊聚類新算法[J].應(yīng)用技術(shù),2009(1) : 43-46.

[8]李曉嵐.基于Relief特征選擇算法的研究與應(yīng)用[D].大連:大連理工大學(xué),2013: 11-30.

[9]劉英乾.渦輪葉片故障診斷與模擬研究[D].哈爾濱:哈爾濱工程大學(xué),2013: 37-39.

[10]劉大響.航空發(fā)動機(jī)葉片故障及預(yù)防研討會論文集[C]//北京:航空航天工業(yè)出版社,2005: 91-97.

Feature extraction of turbine blades based on the fractal theory

FENG Chi1,HU Yang1,WANG Zhaofeng2

1.College of Information and Communication Engineering,Harbin Engineering University,Harbin 150001,China 2.Xi’an Aero-engine (Group) Ltd.,Xi’an 710021,China

The data is preprocessed and features are extracted to establish a feature model of turbine blades by blade temperature,which is an important indicator of the turbine blade's quality.Based on the fractal theory,the three fractal dimension features of blade temperature signals are extracted,the weight of each feature is calculated by combining with K-means clustering analysis and ReliefF algorithm,and thereby establish the temperature feature model of turbine blades,achieving early warning of the turbine blade failure.Statistical results show that the feature model can effectively reflect the state of turbine blades with failure.

turbine blades; feature extraction; fractal theory; K-means; ReliefF

TK473

A

1009-671X(2015) 04-064-06

10.3969/j.issn.1009-671X.201411006

2014-11-14.網(wǎng)絡(luò)出版日期: 2015-07-27.

黑龍江省自然科學(xué)基金資助項(xiàng)目(F201413).

馮馳(1961-),男,教授;胡楊(1990-),女,碩士研究生.

胡楊,E-mail: huyang900218@163.com.

猜你喜歡
特征故障
抓住特征巧觀察
新型冠狀病毒及其流行病學(xué)特征認(rèn)識
故障一點(diǎn)通
如何表達(dá)“特征”
不忠誠的四個特征
抓住特征巧觀察
奔馳R320車ABS、ESP故障燈異常點(diǎn)亮
故障一點(diǎn)通
故障一點(diǎn)通
故障一點(diǎn)通
主站蜘蛛池模板: 中文无码毛片又爽又刺激| 呦女精品网站| 青青草国产精品久久久久| 欧美日韩一区二区三区在线视频| 欧美第一页在线| 亚洲码一区二区三区| 91区国产福利在线观看午夜| 色综合五月婷婷| 成年人福利视频| 国产精品妖精视频| 国产流白浆视频| 69精品在线观看| 波多野结衣的av一区二区三区| 欧美一区二区三区不卡免费| 亚洲第一天堂无码专区| 久久精品aⅴ无码中文字幕| 日本不卡视频在线| 在线精品视频成人网| 亚洲男人天堂网址| 456亚洲人成高清在线| 在线观看免费人成视频色快速| 四虎永久在线视频| 国产日韩久久久久无码精品| 国产靠逼视频| 99无码中文字幕视频| 青青青国产视频| 在线亚洲精品自拍| 免费国产福利| 久久精品女人天堂aaa| 久久精品亚洲热综合一区二区| 国产黄色片在线看| 少妇精品网站| 日本人真淫视频一区二区三区| 亚洲性视频网站| 国产精品美女免费视频大全| 永久免费精品视频| 国产精品浪潮Av| 国产成人精品高清在线| a级毛片在线免费| 中文无码影院| 2021亚洲精品不卡a| 色婷婷色丁香| 欧美视频在线播放观看免费福利资源 | 久久精品欧美一区二区| 成人在线综合| 欧美一级专区免费大片| 亚洲色图综合在线| 在线观看亚洲国产| 国产欧美日韩另类| 欧美劲爆第一页| 四虎永久免费网站| 亚洲中文字幕无码爆乳| 久久99国产精品成人欧美| 最新日韩AV网址在线观看| 成人免费一级片| 国产乱码精品一区二区三区中文| 在线观看国产精品一区| 欧洲极品无码一区二区三区| 亚洲欧美日韩中文字幕在线| 免费在线视频a| 极品尤物av美乳在线观看| 国产高清在线精品一区二区三区| 免费又黄又爽又猛大片午夜| 911亚洲精品| 波多野结衣一区二区三区四区视频 | 国产美女精品在线| 欧美亚洲网| 欧美日韩北条麻妃一区二区| 欧美有码在线| 免费看的一级毛片| 中文字幕无码电影| 久久青草热| 毛片免费观看视频| 中国黄色一级视频| 最新国产网站| 在线免费无码视频| 免费国产福利| 欧美全免费aaaaaa特黄在线| 另类专区亚洲| 国产精品.com| 欧美午夜一区| 国产激情无码一区二区APP|