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

中介軸承復(fù)合故障動力學(xué)建模與振動特征分析

2022-12-01 10:25:32艾辛平劉麗麗張鳳玲
振動與沖擊 2022年22期
關(guān)鍵詞:故障模型

田 晶,艾辛平,劉麗麗,張鳳玲,王 志

(沈陽航空航天大學(xué) 遼寧省航空推進系統(tǒng)先進測試技術(shù)重點實驗室,沈陽 110136)

中介軸承是雙/三轉(zhuǎn)子航空發(fā)動機轉(zhuǎn)子系統(tǒng)的關(guān)鍵部件。中介軸承工作在高溫、高壓、高轉(zhuǎn)速和重載的惡劣工況下,極易發(fā)生機械故障。中介軸承一旦出現(xiàn)故障,將會起整臺發(fā)動機發(fā)生事故。建立準確的中介軸承動力學(xué)模型,分析中介軸承故障特征,對發(fā)展中介軸承故障診斷具有重要意義。

Walters[1]首次開展動力學(xué)建模嘗試,隨后國內(nèi)外專家學(xué)者針對滾動軸承動力學(xué)建模開展了大量研究。Gupta[2-3]建立了對軸承各部件運動及受力狀態(tài)進行描述的動力學(xué)模型。Patel等[4-5]提出了一種深溝球軸承6自由度動力學(xué)模型,采用脈沖序列和延遲脈沖序列組合的形式對恒定載荷下單點故障和多點故障進行仿真分析。Patil等[6]基于非線性接觸理論將滾道和滾子之間的接觸考慮成一個非線性接觸彈簧,提出了一種用于分析局部缺陷的動力學(xué)模型。張建軍等[7]建立的軸承故障6自由度模型,考慮了由軸承間隙引起的變剛度振動,使用單元諧振器模擬故障軸承元件的高頻固有振動。剡昌鋒等[8-9]基于非線性Hertz理論,考慮了彈流潤滑的影響,建立了深溝球軸承單點故障2自由度模型,深入分析了彈流潤滑對滾動體與滾道之間接觸剛度的影響。董振振[10]同樣以Hertz接觸理論為基礎(chǔ),采用脈沖序列對滾動軸承存在的單點故障和復(fù)合故障情況進行分析,建立了考慮撞擊力的動力學(xué)模型。王凱[11]對滾動體經(jīng)過滾道的全過程進行了有限元分析,在此基礎(chǔ)之上建立了相應(yīng)的4自由度動力學(xué)模型。Cao等[12-13]將Gupta的軸承模型與剛性滾動體模型相結(jié)合,建立了球軸承轉(zhuǎn)子系統(tǒng)的動力學(xué)模型。該模型考慮了轉(zhuǎn)子-軸承-軸承座系統(tǒng)的轉(zhuǎn)子彈性變形,研究了在配合間隙影響下軸承系統(tǒng)的振動特性和穩(wěn)定性。結(jié)果表明配合間隙會引起局部沖擊現(xiàn)象。Chen等[14]提出了一種基于Hertz接觸理論的振動信號模型,探究了振動信號的變化過程。在滾動軸承動力學(xué)建模方法基礎(chǔ)上,結(jié)合中介軸承的工作特點,學(xué)者們也提出了多種中介軸承動力學(xué)建模方法。徐可君等[15]建立了一種中介軸承內(nèi)圈存在單一故障的動力學(xué)模型,該模型考慮軸承的徑向游隙以及內(nèi)、外圈的徑向載荷。曹宏瑞等[16-17]以Gupta復(fù)雜動力學(xué)模型為基礎(chǔ),分別研究了高速運轉(zhuǎn)下普通軸承和中介軸承的時變剛度受靜載荷和動載荷的影響規(guī)律,還針對中介軸承的磨損故障和局部缺陷故障建立了相應(yīng)的動力學(xué)模型,并對這兩種故障所帶來的振動響應(yīng)進行了分析。研究表明,基于非線性Hertz接觸理論建立的滾動軸承動力學(xué)模型能夠有效模擬滾動軸承的故障特征,同時相較于其他方式,采用考慮時變位移激勵的建模方式可以更加準確的描述帶有局部缺陷軸承的實際工作狀態(tài)。

本文對中介軸承振動故障特性進行分析,結(jié)合現(xiàn)有滾動軸承動力學(xué)建模方法,提出了一種基于非線性Hertz接觸理論的4自由度動力學(xué)模型,實現(xiàn)了中介軸承內(nèi)圈多點故障、外圈多點故障和內(nèi)外圈多點復(fù)合故障的振動特征分析,模型中采用時變位移激勵描述滾子的多點局部缺陷的沖擊,并通過試驗驗證了所建模型的有效性與準確性。

1 故障動力學(xué)模型建立

1.1 非線性Hertz接觸模型

航空發(fā)動雙轉(zhuǎn)子通常采用內(nèi)、外圈同向或反向旋轉(zhuǎn)方式工作,中介軸承位于高、低壓轉(zhuǎn)子之間,內(nèi)、外圈隨高、低壓轉(zhuǎn)子旋轉(zhuǎn),比普通軸承具有更為復(fù)雜的動力學(xué)特性。本文以內(nèi)、外圈反向旋轉(zhuǎn)的情況為例,在非線性Hertz接觸理論的基礎(chǔ)上,采用時變位移激勵對軸承內(nèi)圈多點故障、外圈多點故障、內(nèi)外圈復(fù)合故障進行模擬,建立4自由度多點故障動力學(xué)模型。模型中假設(shè)滾動體不發(fā)生打滑,內(nèi)圈缺陷和外圈缺陷的形式相同。中介軸承徑向位移圖如圖1所示。

圖1 中介軸承徑向位移圖Fig.1 Radial displacement diagram of inter-shaft bearing

本文采用非線性Hertz接觸力-變形關(guān)系來計算滾子和接觸點之間的彈性趨近量,公式如式(1)所示[18]

Q=Kδn

(1)

式中:K為Hertz接觸剛度;δ為徑向變形量。對于圓柱滾子軸承n=10/9,K=8.06×104l8/9,其中l(wèi)為滾子長度。

因為中介軸承的內(nèi)、外圈均存在振動,所以滾子在內(nèi)外圈之間總的變形量δ為

δ=(xo-xi)cosθbi+(yo-yi)sinθbi-Cr-Ht

(2)

式中:xo和yo分別為外圈中心沿X軸和Y軸方向上的位移;xi和yi為內(nèi)圈中心沿X軸和Y軸方向上的位移;θbi為第i個滾子中心與X軸的夾角;Cr為滾子和滾道之間的初始徑向間隙;Ht為滾子經(jīng)過缺陷時產(chǎn)生的時變位移激勵。

(3)

(4)

式中:θb0為第1個滾子在初始時刻相對于X軸的夾角;Z為中介軸承滾子個數(shù);ωc,ωi和ωo分別為中介軸的保持架、內(nèi)圈和外圈的角速度;db為滾子直徑;Dm為中介軸承節(jié)徑;α為中介軸承接觸角。

由式(1)可得到第i個滾子的非線性Hertz接觸力為

Qi=K[(xo-xi)cosθbi+
(yo-yi)sinθbi-Cr-Ht]10/9

(5)

因此,中介軸承總的非線性Hertz接觸力在X方向和Y方向上的分量分別為

(6)

1.2 故障動力學(xué)模型

基于非線性Hertz接觸理論,建立4自由度中介軸承復(fù)合故障的動力學(xué)模型如下

(7)

式中:Mi為內(nèi)圈質(zhì)量;Mo為外圈質(zhì)量;Fi內(nèi)圈所受偏心載荷;Fo為外圈所受偏心載荷;W為軸承所受徑向載荷;Ki為 內(nèi)圈與滾子接觸剛度;Ko為外圈與滾子接剛度;Ci為軸承徑內(nèi)圈向間隙;Co為軸承徑外圈向間隙;λ為滾子與滾道是否接觸的開關(guān)量,表達式為

(8)

2 軸承故障時變位移激勵

在滾子經(jīng)過缺陷的整個過程中,位移激勵先增大后減小,此過程與半正弦函數(shù)相似,所以采用分段函數(shù)來描述時變位移激勵。圖2所示為滾子經(jīng)過缺陷的過程圖,圖2中滾子經(jīng)過缺陷時產(chǎn)生的最大附加位移H如式(9)所示。

(9)

式中:db為滾子直徑;b為缺陷寬度。時變位移激勵是指滾子在公轉(zhuǎn)過程中,滾子中心在徑向方向上位移隨時間的變化量。

圖2 滾子經(jīng)過缺陷過程圖Fig.2 Process of roller through the defect

2.1 內(nèi)圈多點故障時變位移激勵

本文研究內(nèi)圈同時存在兩點缺陷的情況。中介軸承內(nèi)、外圈反向旋轉(zhuǎn)時滾子和內(nèi)圈缺陷的位置關(guān)系,如圖3所示。圖3中:ωi,ωo和ωc分別為內(nèi)圈、外圈和保持架的角速度;σ01為故障的缺陷角度;σ02為另一個故障的缺陷角度,規(guī)定逆時針為正方向。

圖3 內(nèi)圈多點故障滾子和缺陷位置關(guān)系圖Fig.3 Location relation diagram of roller and defect in inner ring with multi-point faults

假設(shè)ωo>ωi,根據(jù)保持架角速度ωc與內(nèi)圈和外圈角速度的大小關(guān)系,可以判斷其保持架角速度ωc的值處于內(nèi)圈角速度和外圈角速度之間,方向與外圈方向一致。經(jīng)過時間t,兩種故障的位置角度分別為θσ1和θσ2。時變位移激勵Hmi可由式(10)表示為

(10)

(11)

式中:mod為求余;H為滾子經(jīng)過缺陷時產(chǎn)生的最大附加位移;θb0為初始時刻第1個滾子的位置角度;θσ01和θσ02分別為初始時刻兩種故障右側(cè)的位置角度。

2.2 外圈多點故障時變位移激勵

當(dāng)中介軸承外圈同時存在兩點故障,內(nèi)、外圈反向旋轉(zhuǎn)時,滾子經(jīng)過內(nèi)圈缺陷產(chǎn)生的時變位移激勵Hmo如式(12)所示

(12)

(13)

θφ1=ωot+θφ01

(14)

θφ2=ωot+θφ02

(15)

式中:mod為求余;H為最大附加位移;θb0為第1個滾子在初始時刻的位置角;θφ01和θφ02分別為在初始時刻兩故障右側(cè)相對X軸的角度。

圖4 外圈多點故障滾子和缺陷位置關(guān)系圖Fig.4 Location relation diagram of roller and defect in outer ring with multi-point fault

2.3 內(nèi)外圈復(fù)合故障時變位移激勵

當(dāng)中介軸承內(nèi)圈存在單點故障,內(nèi)、外圈反向旋轉(zhuǎn)時,滾子經(jīng)過內(nèi)圈缺陷產(chǎn)生的時變位移激勵Hti如式(16)、式(17)所示

(16)

(17)

當(dāng)中介軸承外圈存在單點故障,內(nèi)、外 圈反向旋轉(zhuǎn)時,滾子經(jīng)過外圈缺陷產(chǎn)生的時變位移激勵Hto如式(18)、式(19)所示

(18)

(19)

基于單點故障的時變位移激勵,當(dāng)設(shè)置內(nèi)、外圈各一個缺陷時,得出滾子在內(nèi)外圈之間總的變形量δ為

δ=(xo-xi)cosθbi+
(yo-yi)sinθbi-Cr-Hti+o

(20)

式中,Hti+o為中介軸承內(nèi)、外圈復(fù)合故障的時變位移激勵。Hti和Hto以及Hti+o之間的關(guān)系為Hti+o=Hti+Hto。將其代入式(2)中可對復(fù)合故障進行描述。

3 試驗驗證

3.1 試驗設(shè)備與軸承參數(shù)

中介軸承支承結(jié)構(gòu)主要應(yīng)用在雙轉(zhuǎn)子航空發(fā)動機中。中介軸承的外圈連接于高壓轉(zhuǎn)子,內(nèi)圈與低壓轉(zhuǎn)子相連,既可以同向旋轉(zhuǎn)又可以反向旋轉(zhuǎn)。根據(jù)中介軸承的實際工作狀態(tài)搭建試驗臺,如圖5所示。

圖5 中介軸承試驗臺Fig.5 Inter-shaft bearing test rig

試驗采用西門子公司的LMS SCADAS Recorder SCR202系統(tǒng)采集振動信號,該采集系統(tǒng)具有2個轉(zhuǎn)速通道,8個振動信號采集通道,最大采樣頻率為102 400 Hz,24位的采樣精度,如圖6所示。

圖6 LMS數(shù)據(jù)采集儀Fig.6 LMS data acquisition instrument

試驗采用6個加速度傳感器,型號均為333B30,傳感器的采樣范圍是0~6 400 Hz,本文試驗設(shè)置的采樣頻率為4 096 Hz,頻率分辨率為0.1 Hz。3個傳感器安裝在靠近低壓電機的軸承座上,靈敏度分別為98.6 mV/g,99 mV/g和97.6 mV/g,如圖7(a)所示;其余3個安裝于中介軸承的軸承座上,靈敏度分別為99.3 mV/g,98.3 mV/g和100.4 mV/g,如圖7(b)所示,用于檢測軸承3個方向的振動。

3.2 軸承結(jié)構(gòu)參數(shù)

試驗采用NSK公司NU202EM型號的圓柱滾子軸承作為試驗軸承。軸承參數(shù)如表1所示。使用電花火線切割的加工方式,在正常的軸承內(nèi)、外圈上,加工出縱向貫穿整個滾道的矩形故障,故障的寬度和深度均為1 mm。

圖7 加速度傳感器布置示意圖Fig.7 Acceleration sensor layout schematic

表1 NU202EM軸承參數(shù)Tab.1 NU202EM bearing parameters

3.3 振動特征分析

3.3.1 中介軸承內(nèi)圈多點故障

本文設(shè)置內(nèi)圈轉(zhuǎn)速為200 r/min,外圈轉(zhuǎn)速為1 000 r/min,根據(jù)文獻[19]的中介軸承故障特征頻率公式計算出此轉(zhuǎn)速下中介軸承的故障特征頻率理論值,如表2所示,其中反向旋轉(zhuǎn)指內(nèi)、外圈反向旋轉(zhuǎn)的情況。

表2 內(nèi)圈多點故障特征頻率理論值Tab.2 Theoretical value of characteristic frequency of inner ring multipoint fault

基于本文所建立的模型計算出內(nèi)圈多點故障的時域波形和包絡(luò)譜[20],如圖8所示。從圖8中可以看出,仿真模擬的故障特征頻率Fi為132.1 Hz,與表2中理論計算故障特征頻率值132 Hz基本一致。圖9表示為試驗測得的時域波形和包絡(luò)譜圖。對比仿真分析與試驗結(jié)果可以發(fā)現(xiàn),在時域波形中可以看出中介軸承內(nèi)圈存在多點故障時,信號具有周期性的脈沖特征。在包絡(luò)譜中可以看到明顯的轉(zhuǎn)頻及其倍頻、故障特征頻及其倍頻。由于滾子與缺陷的接觸處存在周期性變化,故障信號會出現(xiàn)以故障頻率為中心,轉(zhuǎn)頻及其倍頻為邊帶的轉(zhuǎn)頻調(diào)制幅值頻率現(xiàn)象。圖9中試驗測得的故障特征頻率為131 Hz,而理論模型仿真分析的故障特征頻率為132.1 Hz,對比發(fā)現(xiàn)試驗數(shù)據(jù)與理論數(shù)值誤差小于1%。且轉(zhuǎn)頻、倍頻及邊帶頻均與模擬的包絡(luò)譜相同,因而驗證了該模型的有效性與準確性。

圖8 反向旋轉(zhuǎn)時內(nèi)圈多點故障時域波形和包絡(luò)譜仿真結(jié)果Fig.8 Simulation results of time domain waveform and envelope spectrum of inner ring with multi-point fault in reverse rotation

圖9 反向旋轉(zhuǎn)時內(nèi)圈多點故障時域波形和包絡(luò)譜試驗結(jié)果Fig.9 Experimental results of time domain waveform and envelope spectrum of inner ring with multi-point fault in reverse rotation

圖10為滾子與缺陷的位置關(guān)系圖,圖10為提取的局部時域信號細節(jié)圖。圖10中的脈沖疊加現(xiàn)象是由于內(nèi)圈含有2個故障,同一個滾子在相隔的時間內(nèi)和2個故障接觸形成脈沖。假設(shè)圖11中t1時刻的脈沖是由圖10中第1個滾子和第1個故障產(chǎn)生的,則t2時刻的脈沖是由第2個滾子和第2個故障形成的,通過內(nèi)圈和保持架的轉(zhuǎn)速關(guān)系,可以計算出這2個脈沖的時間間隔為0.005 79 s。同理,t3與t1時刻的脈沖均是由第1個故障產(chǎn)生,2個脈沖的時間差值為1/132=0.007 58 s。而通過圖11可計算t2-t1=0.005 7 s,t3與t1的時間間隔為0.007 6 s,兩組數(shù)值基本一致,證明了模型的準確性。

圖10 反向旋轉(zhuǎn)內(nèi)圈多點故障運動分析圖Fig.10 Motion analysis diagram of inner ring with multi-point fault reverse rotation

圖11 仿真模擬反向旋轉(zhuǎn)內(nèi)圈多點故障時域波形細節(jié)圖Fig.11 Simulation results of detailed time domain waveform of inner ring with multi-point fault in reverse rotating

3.3.2 外圈多點故障試驗驗證

在分析外圈多點故障時,設(shè)置內(nèi)圈轉(zhuǎn)速為800 r/min,外圈轉(zhuǎn)速為200 r/min,此轉(zhuǎn)速下中介軸承的故障特征頻率理論值,如表3所示。

表3 外圈多點故障特征頻率理論值Tab.3 Theoretical value of characteristic frequency of outer ring multipoint fault

外圈多點故障的仿真模擬的時域波形和包絡(luò)譜,如圖12所示。仿真模擬的故障特征頻率Fo為73.2 Hz,表3中的經(jīng)驗公式計算故障特征頻率理論值為73.3 Hz,兩者計算誤差為0.1 Hz。從時域波形中可以看出中介軸承外圈存在多點故障時具有周期性的脈沖特征。同一個故障與不同滾子接觸形成的脈沖為故障特征頻率的倒數(shù)。

圖12 反向旋轉(zhuǎn)時外圈多點故障時域波形和包絡(luò)譜仿真結(jié)果Fig.12 Simulation results of time domain waveform and envelope spectrum of outer ring with multi-point fault in reverse rotation

圖13為故障軸承模擬試驗獲得時域波形和包絡(luò)譜圖。圖13中存在周期性的脈沖信號是同樣由于滾子受到內(nèi)、外圈轉(zhuǎn)動及保持架轉(zhuǎn)動頻率的調(diào)制包絡(luò)譜中存在的故障頻率以及故障頻率的邊頻同樣證明該現(xiàn)象的發(fā)生。在包絡(luò)譜中的故障特征頻率試驗值Fo為73.1 Hz,與理論模型仿真模擬的故障特征頻率73.2 Hz同樣僅差0.1 Hz,試驗結(jié)果與理論的故障特征頻率值基本一致。證明了數(shù)值模型的準確性。

圖13 反向旋轉(zhuǎn)時外圈多點故障時域波形和包絡(luò)譜試驗結(jié)果Fig.13 Experimental results of time domain waveform and envelope spectrum of outer ring with multi-point fault in reverse rotation

圖14為滾子與缺陷的位置關(guān)系圖。圖15為提取的局部時域細節(jié)圖。與內(nèi)圈多點故障中滾子的運動相似,根據(jù)相對位置關(guān)系可計算t2與t1的時間間隔為0.003 22 s,t3與t1的時間間隔為故障特征頻率的倒數(shù),即1/73.3=0.013 64 s。而通過圖15可提取得t2-t1=0.003 3 s。同理,t3與t1的時間間隔為0.013 7 s,與分析值保持一致。

圖14 反向旋轉(zhuǎn)外圈多點故障運動分析圖Fig.14 Motion analysis diagram of outer ring with multi-point fault reverse rotation

圖15 仿真模擬反向旋轉(zhuǎn)外圈多點故障時域波形細節(jié)圖Fig.15 Simulation results of detailed time domain waveform of outer ring with multi-point fault in reverse rotating

3.3.3 內(nèi)外圈復(fù)合故障試驗驗證

設(shè)置中介軸承內(nèi)圈轉(zhuǎn)速為300 r/min,外圈轉(zhuǎn)速為1 400 r/min。根據(jù)田晶等的研究計算得到中介軸承內(nèi)、外圈故障頻率,如表4所示。

表4 內(nèi)、外圈復(fù)合故障特征頻率理論值Tab.4 Theoretical value of characteristic frequency of compound fault inner and outer rings

基于本文所建立的模型計算出內(nèi)、外圈復(fù)合故障的仿真模擬時域波形和包絡(luò)譜如圖16所示。在圖16中可以觀察到內(nèi)圈、外圈的故障特征頻率。內(nèi)圈故障特征頻率Fi為187.1 Hz,外圈故障特征頻率Fo為124.5 Hz,二者與表4中經(jīng)驗公式計算的故障特征頻率187.0 Hz和124.7 Hz的誤差值小于1%。圖17為故障軸承由試驗得到的時域波形和包絡(luò)譜圖。時域波形中可以明顯的觀察到存在兩種沖擊信號,包絡(luò)譜中同樣包含內(nèi)圈和外圈各自的轉(zhuǎn)頻、故障特征頻率和邊帶頻及其二倍頻。這說明,內(nèi)、外圈均有故障時,內(nèi)、外圈故障相互獨立,互不干擾,內(nèi)外圈復(fù)合故障信號為內(nèi)、外圈單點故障信號的疊加,且受轉(zhuǎn)頻調(diào)制。在包絡(luò)譜中試驗測得的內(nèi)圈故障特征頻率試驗值為188.8 Hz,外圈故障特征頻率試驗值為125.1 Hz,對比理論模型仿真得到的故障特征頻率值187.1 Hz和124.5 Hz,誤差同樣小于1%。因而驗證了該模型的有效性與準確性。

圖16 反向旋轉(zhuǎn)時內(nèi)、外圈復(fù)合故障時域波形和包絡(luò)譜仿真結(jié)果Fig.16 Simulation results of time domain waveform and envelope spectrum of compound fault inner and outer rings in reverse rotation

圖17 反向旋轉(zhuǎn)時內(nèi)、外圈復(fù)合故障時域波形和包絡(luò)譜試驗結(jié)果Fig.17 Experimental results of time domain waveform and envelope spectrum of compound fault inner and outer rings in reverse rotation

4 結(jié) 論

(1) 基于本文建立的中介軸承內(nèi)圈多點故障、外圈多點故障和內(nèi)外圈復(fù)合故障模型數(shù)值模擬獲得的故障特征頻率值與試驗值的誤差均小于1%,獲取的時域波形及包絡(luò)譜中的時頻信息分布規(guī)律與試驗一致,驗證了本文所建立模型的準確性。

(2) 同一個滾子在旋轉(zhuǎn)的整個過程中會與不同的故障相撞形成脈沖,兩次脈沖在時間上會有間隔,時域圖中表現(xiàn)出脈沖疊加的現(xiàn)象。由模型計算得到的同一故障、不同滾子形成的脈沖的時間間隔數(shù)值與理論分析值基本一致,且由同一故障與不同滾子相接觸形成的相鄰脈沖,其時間間隔為故障特征頻率的倒數(shù)。

(3) 中介軸承內(nèi)圈、外圈存在復(fù)合故障時,在時域波形和包絡(luò)譜中可以看到復(fù)合故障信號中包括內(nèi)圈故障信號及外圈故障信號,且每種信號皆受與其對應(yīng)的轉(zhuǎn)頻的調(diào)制。

猜你喜歡
故障模型
一半模型
重要模型『一線三等角』
重尾非線性自回歸模型自加權(quán)M-估計的漸近分布
故障一點通
3D打印中的模型分割與打包
奔馳R320車ABS、ESP故障燈異常點亮
FLUKA幾何模型到CAD幾何模型轉(zhuǎn)換方法初步研究
故障一點通
故障一點通
故障一點通
主站蜘蛛池模板: 日韩欧美中文字幕在线精品| 亚洲人成影视在线观看| 亚洲Aⅴ无码专区在线观看q| 国产一级在线播放| 日韩福利在线视频| 天天综合色网| 久久国产高潮流白浆免费观看| 凹凸国产熟女精品视频| 国产精品露脸视频| 欧美a√在线| 色欲色欲久久综合网| 日本欧美精品| 国产美女91视频| 97视频免费看| 亚洲综合极品香蕉久久网| 久久中文无码精品| 国产精品美女网站| 国产精品va免费视频| 亚洲精品成人7777在线观看| a毛片免费观看| 亚洲三级影院| 国产亚洲视频免费播放| 99人妻碰碰碰久久久久禁片| 欧美一区福利| 男人天堂亚洲天堂| 久久特级毛片| 色综合网址| 国产原创自拍不卡第一页| 国产a v无码专区亚洲av| 欧美人人干| 欧美日韩亚洲综合在线观看| 992tv国产人成在线观看| 欧美日韩第三页| 91福利免费| 国产美女丝袜高潮| 天天色天天操综合网| 免费毛片全部不收费的| 综合天天色| 99久久精品国产自免费| 美女毛片在线| 午夜福利在线观看入口| 国产欧美在线| 久久婷婷国产综合尤物精品| 欧美日韩高清在线| 国产免费a级片| 欧美日韩成人| 国产成人一区| 亚洲天堂福利视频| 亚洲午夜福利精品无码| 亚洲欧美在线综合一区二区三区| 国产日韩欧美视频| 麻豆精品在线视频| 亚洲成人精品在线| 久久亚洲综合伊人| 激情爆乳一区二区| 欧美中文字幕在线播放| 日韩大乳视频中文字幕| 色老头综合网| 久久久久无码国产精品不卡| 91久久天天躁狠狠躁夜夜| 欧美色图第一页| 久久久精品无码一二三区| 极品尤物av美乳在线观看| 国产真实乱了在线播放| 无码高潮喷水专区久久| 国产网站免费| 91最新精品视频发布页| 亚洲一区二区日韩欧美gif| 亚洲综合中文字幕国产精品欧美| 亚洲成人高清在线观看| 91精品综合| 免费观看无遮挡www的小视频| 特级毛片免费视频| 国产成人久视频免费| 亚洲成人播放| 亚洲精品无码不卡在线播放| 福利在线不卡| 国产精品网曝门免费视频| 亚洲成人77777| 91蜜芽尤物福利在线观看| 国产一区亚洲一区| 国产青榴视频在线观看网站|