亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        船用燃機(jī)低壓轉(zhuǎn)子-軸承系統(tǒng)非線性動(dòng)力學(xué)分析

        2019-07-22 12:10:02李默劉永葆王強(qiáng)
        軸承 2019年11期
        關(guān)鍵詞:龐加萊游隙徑向

        李默,劉永葆,王強(qiáng),b

        (海軍工程大學(xué) a.動(dòng)力工程學(xué)院;b.兵器工程學(xué)院,武漢 430033)

        滾動(dòng)軸承因摩擦因數(shù)小,質(zhì)量輕,易于安裝和冷卻等優(yōu)點(diǎn),廣泛應(yīng)用于船用燃?xì)廨啓C(jī)和航空發(fā)動(dòng)機(jī)等轉(zhuǎn)子系統(tǒng),為轉(zhuǎn)子系統(tǒng)提供良好的穩(wěn)定支承。但軸承變剛度振動(dòng)、非線性Hertz恢復(fù)力、徑向游隙等非線性因素會(huì)引發(fā)轉(zhuǎn)子系統(tǒng)的非線性動(dòng)力學(xué)行為。

        國(guó)內(nèi)外對(duì)轉(zhuǎn)子-軸承系統(tǒng)的非線性動(dòng)力學(xué)特性做了大量研究。文獻(xiàn)[1]考慮軸承徑向游隙、變剛度頻率等多個(gè)非線性因素,分析了轉(zhuǎn)速對(duì)轉(zhuǎn)子-軸承系統(tǒng)的非線性動(dòng)力學(xué)特性的影響;文獻(xiàn)[2]建立了兩自由度轉(zhuǎn)子-軸承系統(tǒng),分析了系統(tǒng)阻尼變化對(duì)系統(tǒng)分岔、混沌的影響;文獻(xiàn)[3]建立了轉(zhuǎn)子-軸承系統(tǒng)動(dòng)力學(xué)模型,結(jié)合Poincaré 映射圖等進(jìn)行了轉(zhuǎn)子系統(tǒng)非線性特性分析;文獻(xiàn)[4]考慮軸承波紋度等系統(tǒng)的非線性因素,分析了波紋度頻率對(duì)系統(tǒng)振動(dòng)的影響;文獻(xiàn)[5]分析了不同轉(zhuǎn)速下轉(zhuǎn)子-軸承系統(tǒng)的分岔、混沌等非線性動(dòng)力學(xué)特性;文獻(xiàn)[6]考慮軸承的接觸非線性和游隙非線性建立了偏置轉(zhuǎn)子系統(tǒng)動(dòng)力學(xué)模型,分析了徑向游隙對(duì)偏置轉(zhuǎn)子系統(tǒng)的非線性動(dòng)力學(xué)影響;文獻(xiàn)[7]建立了三自由度對(duì)稱支承的轉(zhuǎn)子-軸承系統(tǒng)動(dòng)力學(xué)模型,分析了徑向游隙對(duì)系統(tǒng)振動(dòng)的影響;文獻(xiàn)[8]建立了單自由度剛性轉(zhuǎn)子非線性動(dòng)力學(xué)模型,分析了徑向游隙對(duì)系統(tǒng)非線性動(dòng)力學(xué)的影響;文獻(xiàn)[9]建立單自由度轉(zhuǎn)子-軸承系統(tǒng),分析軸承變剛度振動(dòng)和強(qiáng)迫振動(dòng)耦合振動(dòng)對(duì)轉(zhuǎn)子系統(tǒng)非線性動(dòng)力學(xué)特性的影響;文獻(xiàn)[10]建立單跨轉(zhuǎn)子系統(tǒng)動(dòng)力學(xué)模型,分析了轉(zhuǎn)速、擠壓油膜間隙對(duì)系統(tǒng)非線性動(dòng)力學(xué)特性的影響;文獻(xiàn)[11]建立轉(zhuǎn)子-軸承系統(tǒng)有限元模型,分析了結(jié)構(gòu)參數(shù)對(duì)連續(xù)質(zhì)量的轉(zhuǎn)子系統(tǒng)非線性動(dòng)力學(xué)特性的影響;文獻(xiàn)[12]建立了陀螺儀轉(zhuǎn)子-軸承系統(tǒng)非線性動(dòng)力學(xué)模型,分析了轉(zhuǎn)子-軸承系統(tǒng)的阻尼剛度對(duì)系統(tǒng)非線性動(dòng)力學(xué)特性的影響;文獻(xiàn)[13-15]建立了偏置轉(zhuǎn)子系統(tǒng)非線性動(dòng)力學(xué)模型,分析了偏置距離和徑向游隙對(duì)轉(zhuǎn)子系統(tǒng)非線性動(dòng)力學(xué)特性的影響;文獻(xiàn)[16-17]建立了四支承雙跨轉(zhuǎn)子系統(tǒng)動(dòng)力學(xué)模型,分析了軸承支承下的故障轉(zhuǎn)子非線性動(dòng)力學(xué)特性;文獻(xiàn)[18-19]建立了三支承雙跨轉(zhuǎn)子系統(tǒng)動(dòng)力學(xué)模型,分析了滑動(dòng)軸承支承下故障轉(zhuǎn)子系統(tǒng)的非線性動(dòng)力學(xué)特性。

        上述研究對(duì)轉(zhuǎn)子-軸承系統(tǒng)進(jìn)行了非線性動(dòng)力學(xué)分析,但僅針對(duì)簡(jiǎn)單的單跨單盤雙支承轉(zhuǎn)子-軸承系統(tǒng),結(jié)構(gòu)參數(shù)單一,未根據(jù)實(shí)際轉(zhuǎn)子-軸承結(jié)構(gòu)進(jìn)行建模。鑒于此,以某船用燃?xì)廨啓C(jī)轉(zhuǎn)子-軸承系統(tǒng)為例,建立1套深溝球軸承和2套圓柱滾子軸承的三支承雙跨形式轉(zhuǎn)子-軸承系統(tǒng),分析支承軸承轉(zhuǎn)速和徑向游隙對(duì)系統(tǒng)非線性特性的影響。

        1 滾動(dòng)軸承建模

        為便于對(duì)轉(zhuǎn)子-軸承系統(tǒng)進(jìn)行非線性動(dòng)力學(xué)分析,做如下假設(shè):1)軸與軸承內(nèi)圈剛性連接;2)軸承座與軸承外圈剛性連接;3)滾動(dòng)體在滾道內(nèi)自由滾動(dòng),不產(chǎn)生滑動(dòng)。滾動(dòng)軸承簡(jiǎn)化模型如圖1所示,圖中:O1為軸承不受載時(shí)的軸質(zhì)心;O2為軸承受載后的軸質(zhì)心;Ri為內(nèi)圈半徑;Re為外圈半徑;Gr為徑向游隙;j表示第j個(gè)滾動(dòng)體;θj為第j個(gè)滾動(dòng)體位置角。

        根據(jù)幾何關(guān)系可得

        (1)

        式中:Z為滾動(dòng)體數(shù)量;t為運(yùn)動(dòng)時(shí)間。

        圖1 滾動(dòng)軸承簡(jiǎn)化模型

        第j個(gè)滾動(dòng)體接觸變形量為

        χj=xcosθj+ysinθj-Gr,

        (2)

        式中:x,y為受載后軸質(zhì)心坐標(biāo)。

        接觸只能產(chǎn)生正向接觸壓力,當(dāng)χj<0時(shí),不會(huì)產(chǎn)生壓力,僅χj>0時(shí)才會(huì)有壓力,引入Heaviside函數(shù),即

        (3)

        根據(jù)Hertz接觸理論得到軸承恢復(fù)力為

        (4)

        式中:Fbx,Fby分別為軸承在x,y方向的恢復(fù)力;Kb為接觸剛度;Fj為第j個(gè)滾動(dòng)體所受載荷;n為載荷-變形指數(shù),對(duì)于深溝球軸承,n=3/2,對(duì)于圓柱滾子軸承,n=10/9。

        2 轉(zhuǎn)子-軸承系統(tǒng)建模

        根據(jù)某船用燃?xì)廨啓C(jī)低壓轉(zhuǎn)子-軸承系統(tǒng)的實(shí)際支承形式,將轉(zhuǎn)子系統(tǒng)簡(jiǎn)化為三支承雙跨的三盤轉(zhuǎn)子系統(tǒng),如圖2所示。圖中:O為坐標(biāo)原點(diǎn);Ob1,Ob2,Ob3分別為軸承1,2,3的形心;Op1,Op2,Op3分別為盤1,2,3的形心;Oc1,Oc2,Oc2分別盤1,2,3的質(zhì)心;e1,e2,e3分別為盤1,2,3的偏心距。A段模擬壓氣機(jī)段,由2個(gè)剛性盤模擬,分別位于壓氣機(jī)段的三等分處;B段模擬渦輪段,由一個(gè)剛性盤模擬,位于渦輪段中點(diǎn):整個(gè)轉(zhuǎn)子系統(tǒng)由3套軸承支承(軸承1為深溝球軸承,軸承2,3為圓柱滾子軸承)。模型考慮轉(zhuǎn)子在軸承和圓盤處的集中質(zhì)量,轉(zhuǎn)子與軸承之間采用無質(zhì)量彈性軸,忽略轉(zhuǎn)子扭轉(zhuǎn)振動(dòng)和軸向振動(dòng),僅研究橫向振動(dòng),每個(gè)集中質(zhì)量只考慮x,y這2個(gè)方向的自由度。

        圖2 轉(zhuǎn)子-軸承系統(tǒng)簡(jiǎn)化模型

        假設(shè)3個(gè)盤偏心距相同,由于主要研究軸承對(duì)系統(tǒng)非線性響應(yīng)的影響,假設(shè)壓氣機(jī)段和渦輪段轉(zhuǎn)子結(jié)構(gòu)、材料相同,根據(jù)平面梁結(jié)構(gòu)關(guān)系得到軸段的剛度關(guān)系為:k1=k3=k/2,k2=2k,k4=k5=k,k=2.5×107N/m。根據(jù)Lagrange運(yùn)動(dòng)方程和轉(zhuǎn)子動(dòng)力學(xué)理論得到系統(tǒng)非線性動(dòng)力學(xué)微分方程組為

        (5)

        M=diag(mb1,mp1,mp2,mb2,mp3,mb3),

        C=diag(cb1,cp1,cp2,cb2,cp3,cb3),

        K=

        Fbx=(Fbx,0,0,Fmx,0,Frx)T,

        Fby=(Fby,0,0,Fmy,0,Fry)T,

        X=(xb1,xp1,xp2,xb2,xp3,xb3)T,

        Y=(yb1,yp1,yp2,yb2,yp3,yb3)T,

        Fex=(0,mp1e1ω2cos(ωt),mp2e2ω2cos(ωt),

        0,mp3e3ω2cos(ωt),0)T,

        Fey=(0,mp1e1ω2sin(ωt),mp2e2ω2sin(ωt),

        0,mp3e3ω2sin(ωt),0)T,

        式中:M為系統(tǒng)剛度矩陣;C為系統(tǒng)阻尼矩陣;K為系統(tǒng)剛度矩陣;Fbx,F(xiàn)by為軸承恢復(fù)力矩陣;Fex,F(xiàn)ey為轉(zhuǎn)子不平衡力矩陣;g為重力加速度;mb1,mb2,mb3分別為轉(zhuǎn)子在軸承1,2,3處的集中質(zhì)量;mp1,mp2,mp3分別為轉(zhuǎn)子在盤1,2,3處的集中質(zhì)量;cb1,cb2,cb3分別為轉(zhuǎn)子在軸承1,2,3處的阻尼;cp1,cp2,cp3分別為轉(zhuǎn)子在盤1,2,3處的阻尼;Flx,F(xiàn)mx,F(xiàn)rx為軸承1,2,3恢復(fù)力的水平分量;Fly,F(xiàn)my,F(xiàn)ry為軸承1,2,3恢復(fù)力的垂直分量;ω為角速度;xb1,xb2,xb3分別為軸承1,2,3的水平位移;xp1,xp2,xp3分別為盤1,2,3的水平位移;yb1,yb2,yb3分別為軸承1,2,3的垂直位移;yp1,yp2,yp3分別為盤1, 2,3的垂直位移。

        引入量綱一化參數(shù)

        (6)

        式中:γ是與軸承徑向游隙相當(dāng)?shù)囊粋€(gè)長(zhǎng)度量,數(shù)值選擇根據(jù)數(shù)值計(jì)算方便。

        對(duì)(5)式進(jìn)行量綱一化處理可得

        (7)

        3 系統(tǒng)非線性動(dòng)力學(xué)分析

        軸承主要結(jié)構(gòu)參數(shù)見表1,轉(zhuǎn)子參數(shù)為:mb1=1.2 kg,mb2=1.2 kg,mb3=1.2 kg,mp1=2.09 kg,mp2=11.00 kg,mp3=6.00 kg,cb1=2 200 N·s/m,cb2=2 200 N·s/m,cb3=2 200 N·s/m,cp1=1 050 N·s/m,cp2=1 050 N·s/m,cp3=1 050 N·s/m,k=2.5×107N/m,e1=e2=e3=10 μm。微分方程組具有強(qiáng)非線性,故采用四階Runge-Kutta法求解。變剛度頻率fvc=frRi/(Ri+Re)(fr為軸旋轉(zhuǎn)頻率)。

        表1 軸承主要結(jié)構(gòu)參數(shù)

        3.1 轉(zhuǎn)速對(duì)系統(tǒng)非線性動(dòng)力學(xué)特性的影響

        隨轉(zhuǎn)速逐漸增大,旋轉(zhuǎn)機(jī)械運(yùn)動(dòng)狀態(tài)發(fā)生復(fù)雜變化。以轉(zhuǎn)速n為變量,分析盤1水平位移的分岔非線性特性,結(jié)果如圖3所示。由圖3可知:隨轉(zhuǎn)速增大,盤1經(jīng)歷了擬周期、單周期、擬周期、周期二、周期四、周期二和單周期運(yùn)動(dòng)過程,說明系統(tǒng)具有復(fù)雜的非線性特性。

        圖3 盤1水平位移隨轉(zhuǎn)速變化的分岔圖

        在n=1 432 r/min時(shí)盤1的頻譜圖如圖4所示(頻率比λ為特征頻率與旋轉(zhuǎn)頻率的比值),由圖可知:系統(tǒng)有3個(gè)激振頻率,1個(gè)轉(zhuǎn)速激勵(lì)頻率和2個(gè)軸承變剛度頻率。在fr,fvc1(軸承1變剛度頻率),fvc2(軸承2,3變剛度頻率)處出現(xiàn)峰值,且變剛度頻率幅值大于旋轉(zhuǎn)頻率幅值,系統(tǒng)振動(dòng)頻率以軸承變剛度頻率為主。分析結(jié)果與文獻(xiàn)[20]一致。

        圖4 在n=1 432 r/min時(shí)盤1的頻譜圖

        盤1水平位移隨轉(zhuǎn)速變化的三維譜圖如圖5所示,由圖可知:隨轉(zhuǎn)速增大,盤1變剛度頻率幅值逐漸減小,旋轉(zhuǎn)頻率幅值逐漸增大,不平衡激勵(lì)為盤1的主要激勵(lì)。當(dāng)轉(zhuǎn)速遠(yuǎn)小于臨界轉(zhuǎn)速(8 000 r/min)時(shí),盤1變剛度頻率和旋轉(zhuǎn)頻率的作用差別不大,軸承變剛度頻率幅值與旋轉(zhuǎn)頻率幅值接近,且頻率不可相約,故盤1表現(xiàn)為擬周期運(yùn)動(dòng)。

        圖5 盤1水平位移隨轉(zhuǎn)速變化的三維譜圖

        不同轉(zhuǎn)速下盤1的龐加萊圖和頻譜圖如圖6所示,由圖可知:1)n=5 254 r/min時(shí),龐加萊圖為一個(gè)封閉曲線,頻譜圖有多個(gè)不可相約的主要峰值;2)n=9 745 r/min時(shí),龐加萊圖為一個(gè)點(diǎn),頻譜圖上主要是旋轉(zhuǎn)頻率及其倍頻;3)n=14 426 r/min時(shí),龐加萊圖為一個(gè)不規(guī)則的封閉曲線,頻譜圖上出現(xiàn)較多的低倍頻。

        圖6 盤1的龐加萊圖(左)和頻譜圖(右)

        隨轉(zhuǎn)速增加,不平衡激勵(lì)及其分頻越來越明顯,盤1逐漸出現(xiàn)分岔運(yùn)動(dòng)。更高的轉(zhuǎn)速下,盤1的龐加萊圖如圖7所示,由圖可知:在轉(zhuǎn)速為18 535,20 636,22 356,24 936 r/min時(shí)龐加萊圖分別為2個(gè)單點(diǎn)、4個(gè)單點(diǎn)、2個(gè)單點(diǎn)和1個(gè)單點(diǎn),對(duì)應(yīng)盤1的周期二、周期四、周期二和單周期運(yùn)動(dòng)。

        圖7 高轉(zhuǎn)速下盤1的龐加萊圖

        3.2 徑向游隙對(duì)系統(tǒng)非線性動(dòng)力學(xué)特性的影響

        3.2.1 3套軸承同時(shí)變化

        轉(zhuǎn)子-軸承系統(tǒng)非線性因素主要來自軸承恢復(fù)力,體現(xiàn)為軸承載荷-變形指數(shù)和恢復(fù)力分段性。根據(jù)Hertz接觸理論,載荷-變形指數(shù)與軸承類型有關(guān),恢復(fù)力分段性由徑向游隙決定,與軸承結(jié)構(gòu)參數(shù)有關(guān)。在此分析徑向游隙同時(shí)變化對(duì)系統(tǒng)非線性動(dòng)力學(xué)特性的影響。

        3套軸承徑向游隙相同,隨徑向游隙增大盤1水平位移的分岔圖如圖8所示,由圖可知:隨徑向游隙增大,盤1經(jīng)歷了單周期、周期二、周期四和周期二運(yùn)動(dòng)。當(dāng)徑向游隙小于6.4 μm時(shí),盤1進(jìn)行單周期運(yùn)動(dòng),此時(shí)徑向游隙較小,軸承恢復(fù)力非線性比較弱,軸承恢復(fù)力對(duì)系統(tǒng)非線性運(yùn)動(dòng)影響較??;徑向游隙為6.4~11.1 μm時(shí),盤1進(jìn)行周期二運(yùn)動(dòng),徑向游隙增大使盤1發(fā)生倍周期分岔;徑向游隙為11.1~ 25.9 μm時(shí),盤1進(jìn)行周期四運(yùn)動(dòng);最后進(jìn)行周期二運(yùn)動(dòng)。

        圖8 隨徑向游隙增大盤1水平位移的分岔圖

        不同徑向游隙時(shí)盤1水平位移的相圖如圖9所示,由圖可知: 當(dāng)徑向游隙為5 μm,相圖為一個(gè)封閉曲線;當(dāng)徑向游隙為8 μm,相圖為2個(gè)封閉曲線;當(dāng)徑向游隙為20 μm時(shí),相圖為4個(gè)封閉曲線;當(dāng)徑向游隙為28 μm時(shí),相圖為2個(gè)封閉曲線。不同徑向游隙時(shí)盤1的龐加萊圖如圖10所示,由圖可知:當(dāng)徑向游隙為5 μm時(shí),龐加萊圖為孤立的點(diǎn);當(dāng)徑向游隙為8 μm時(shí),龐加萊圖為2個(gè)孤立的點(diǎn);當(dāng)徑向游隙為20 μm時(shí),龐加萊圖為4個(gè)孤立的點(diǎn);當(dāng)徑向游隙為28 μm,龐加萊圖又是2個(gè)孤立的點(diǎn);分析結(jié)果與盤1的水平位移分岔圖對(duì)應(yīng)。

        圖9 不同徑向游隙時(shí)盤1的相圖

        圖10 不同徑向游隙時(shí)盤1的龐加萊圖

        3.2.2 單套軸承徑向游隙變化

        在實(shí)際工作過程中,由于系統(tǒng)質(zhì)量沿軸向分布不均或由于潤(rùn)滑不足等原因,滾動(dòng)體與內(nèi)、外滾道接觸而磨損,會(huì)使軸承徑向游隙發(fā)生較大的變化,偏離了原始設(shè)計(jì)值。根據(jù)某船用燃?xì)廨啓C(jī)實(shí)際工況,渦輪處軸承工作在高溫環(huán)境中時(shí),軸承徑向游隙變化較大,在此分析右端軸承(軸承3)徑向游隙變化對(duì)系統(tǒng)非線性動(dòng)力學(xué)特性的影響。

        隨軸承3徑向游隙變化,盤1水平位移的分岔圖如圖11所示。由圖可知:當(dāng)徑向游隙小于2.2 μm時(shí),盤1進(jìn)行單周期運(yùn)動(dòng);當(dāng)徑向游隙超過2.2 μm時(shí),盤1通過倍周期分岔進(jìn)入周期二運(yùn)動(dòng)狀態(tài),且大范圍游隙內(nèi)盤1進(jìn)行周期二運(yùn)動(dòng)。

        圖11 盤1水平位移隨軸承3徑向游隙變化的分岔圖

        軸承3不同徑向游隙時(shí)盤1的相圖和龐加萊截面圖如圖12所示,由圖可知:1)徑向游隙為2 μm時(shí),盤1相圖只有一個(gè)封閉的圓形,龐加萊圖上只有一個(gè)孤立的單點(diǎn),充分說明了此時(shí)盤1處于單周期運(yùn)動(dòng)狀態(tài);2)徑向游隙為20 μm時(shí),盤1的相圖有2個(gè)互相交叉的封閉圓形,龐加萊圖有2個(gè)孤立的點(diǎn),說明此時(shí)盤1處于周期二運(yùn)動(dòng)狀態(tài)。

        圖12 軸承3不同徑向游隙時(shí)盤1的相圖和龐加萊圖

        4 結(jié)論

        以某船用燃?xì)廨啓C(jī)為研究對(duì)象,建立了三支承雙跨三盤滾動(dòng)軸承支承動(dòng)力學(xué)模型,分析了轉(zhuǎn)速、徑向游隙對(duì)轉(zhuǎn)子-軸承系統(tǒng)非線性動(dòng)力學(xué)特性的影響,得出以下結(jié)論:

        1)隨轉(zhuǎn)速增大,軸承變剛度影響變?nèi)?,系統(tǒng)主要表現(xiàn)為轉(zhuǎn)頻的倍頻和分頻運(yùn)動(dòng),出現(xiàn)擬周期、單周期、周期二、周期四等復(fù)雜非線性運(yùn)動(dòng)。

        2)3套軸承徑向游隙同時(shí)增大時(shí),系統(tǒng)通過倍周期分岔途徑運(yùn)動(dòng),且軸承徑向游隙越大,系統(tǒng)運(yùn)動(dòng)越復(fù)雜,當(dāng)徑向游隙超過一定范圍,由于游隙過大,滾動(dòng)體與內(nèi)外圈接觸不充分,系統(tǒng)表現(xiàn)為周期二運(yùn)動(dòng)。在軸承設(shè)計(jì)時(shí),應(yīng)根據(jù)轉(zhuǎn)子系統(tǒng)工作環(huán)境溫度、載荷等,控制軸承徑向游隙小于6.4 μm,保證轉(zhuǎn)子系統(tǒng)處于穩(wěn)定的單周期運(yùn)動(dòng)狀態(tài)。

        3)由于轉(zhuǎn)子-軸承系統(tǒng)右端軸承徑向游隙增大,系統(tǒng)通過倍周期分岔途徑由單周期過渡到周期二運(yùn)動(dòng),該非線性特性可作為右端軸承徑向游隙增大的判斷依據(jù),為軸承故障診斷提供理論依據(jù)。

        猜你喜歡
        龐加萊游隙徑向
        龐加萊偵察術(shù)
        中外文摘(2022年7期)2022-05-17 09:36:42
        龐加萊偵察術(shù)
        青年文摘(2021年20期)2021-12-11 18:45:12
        龐加萊偵查術(shù)
        淺探徑向連接體的圓周運(yùn)動(dòng)
        RN上一類Kirchhoff型方程徑向?qū)ΨQ正解的存在性
        軸承游隙對(duì)隔爆電機(jī)的影響分析
        水平軸設(shè)計(jì)中軸承游隙的計(jì)算與選擇
        基于PID+前饋的3MN徑向鍛造機(jī)控制系統(tǒng)的研究
        一類無窮下級(jí)整函數(shù)的Julia集的徑向分布
        推挽式直流變換器的龐加萊映射圖分析
        国产一区二区自拍刺激在线观看| aaaaaa级特色特黄的毛片| 国产a v无码专区亚洲av| 91超碰在线观看免费| 国产熟女自拍视频网站| 日韩av中文字幕波多野九色| 精品福利一区二区三区蜜桃| 久久只精品99品免费久23| 亚洲国产韩国欧美在线| 亚洲成a人片在线观看导航| 日本大胆人体亚裔一区二区| 精品国产一区二区三区香| 中国孕妇变态孕交xxxx| 无码不卡av东京热毛片| 亚洲av色福利天堂| 日本中文字幕av网址| 国产天堂av在线播放资源| 国产av一区二区三区无码野战| 亚洲日韩国产一区二区三区在线| 一区二区无码中出| 日本精品国产1区2区3区| 精品国产一区二区三区香| 蜜臀av在线播放一区二区三区| 99香蕉国产精品偷在线观看| 亚洲无码夜夜操| 亚洲一区视频中文字幕| 亚洲精一区二区三av| 国产成人无码av一区二区| 真人男女做爰无遮挡免费视频| 亚洲国产精品一区二区第一| 丝袜美腿福利视频在线| 亚洲av网一区二区三区| 精品国产乱码久久久久久口爆网站| 国产在线天堂av| 国产二区中文字幕在线观看| 国产成人精品免费视频大全软件| 特级做a爰片毛片免费看无码| 亚洲欧美日韩中文字幕网址| 国内国外日产一区二区| 韩国av一区二区三区不卡| 亚洲熟女乱色综合亚洲av|