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

基于格林函數(shù)的輪軌缺陷對輪軌垂向動態(tài)相互作用的研究

2017-04-27 01:47:58張寶安鄭靜周和超
關鍵詞:振動

張寶安, 鄭靜, 周和超

(1.中車青島四方機車車輛股份有限公司a.國家工程研究中心;b.技術中心, 山東青島266000;2.同濟大學鐵道與城市軌道交通研究院, 上海201804)

基于格林函數(shù)的輪軌缺陷對輪軌垂向動態(tài)相互作用的研究

張寶安1a, 鄭靜1b, 周和超2

(1.中車青島四方機車車輛股份有限公司a.國家工程研究中心;b.技術中心, 山東青島266000;2.同濟大學鐵道與城市軌道交通研究院, 上海201804)

在輪軌垂向耦合振動模型中,把鋼軌考慮成離散支撐的連續(xù)的鐵木辛科梁,用一系懸掛和二系懸掛把在鋼軌上滾動的車輪、1/4構(gòu)架質(zhì)量塊和1/8車體質(zhì)量塊連接起來作為簡化的車輛模型。為了高效精確地研究不同車輪和鋼軌缺陷對輪軌垂向動態(tài)相互作用的影響,對鋼軌在移動載荷作用下的運動方程進行拉普拉斯變換,用格林函數(shù)法求解其在頻域上的響應解,再通過反傅里葉變換求得鋼軌在時域內(nèi)的振動響應。考慮不同工況下車輪的扁疤、多邊形化缺陷和鋼軌的波磨缺陷,用非線性赫茲接觸理論和數(shù)值迭代方法計算輪軌垂向沖擊力。結(jié)果表明:應用格林函數(shù)法可以高效準確地求解輪軌垂向動態(tài)響應,為在高頻范圍內(nèi)研究輪軌作用提供了理論基礎。

鐵木辛科梁;格林函數(shù);輪軌垂向沖擊力;車輪扁疤;車輪多邊形化;鋼軌波磨

引言

車輪圓周非圓化現(xiàn)象廣泛存在于鐵路運輸現(xiàn)場。制動熱傷損和滾動接觸疲勞會引起扁疤、剝離等車輪滾動圓局部非圓化,車輪半徑沿整個圓周呈周期性變化是車輪多邊形化的車輪滾動圓全周非圓化的典型特征。車輪圓周的局部非圓化將會引起極大的輪軌高頻沖擊載荷,車輪圓周的全局非圓化將引起低頻輪軌接觸力的增大,導致輪軌和車輛簧下結(jié)構(gòu)的損傷。車輪圓周非圓化會對車輛-軌道系統(tǒng)造成很大的沖擊作用,影響車輛零件的強度,加劇輪對的損壞,還對軌道產(chǎn)生極大的影響,加劇輪軌間磨損,影響輪軌服役時間,對行車穩(wěn)定性、安全性有很大影響,并且產(chǎn)生的沖擊噪聲或增大的滾動噪聲又將影響乘客的舒適度。鋼軌波磨的存在導致車輛軌道結(jié)構(gòu)激烈的振動,不僅影響車輛和軌道結(jié)構(gòu)的使用壽命,而且嚴重的鋼軌波磨會導致重大脫軌事故的發(fā)生,增加了鐵路養(yǎng)護部門的維修工作量和維修費用。

因此,各國學者在輪軌垂向相互作用的數(shù)值模擬、輪軌缺陷形成機理方面展開了一些研究。Nielsen等分別用有限元方法[1]、復模態(tài)疊加法[2]和模態(tài)子結(jié)構(gòu)的方法[3]建立了輪軌振動模型,研究了不同波長的扁疤產(chǎn)生的最大輪軌沖擊力,確定了對應的輪軌接觸丟失的車輛運行速度。Wu和Thompson等[4-6]用數(shù)值解析的方法建立軌道的鐵木辛科梁模型,分析了在連續(xù)兩層彈性支撐和離散支撐上的軌道受到不同尺寸的車輪扁疤沖擊力和沖擊噪音。Baeza[7]分析比較了赫茲模型和三維非赫茲模型中車輪扁疤產(chǎn)生的動態(tài)響應。文獻[8]應用車輛-軌道耦合動力學理論分析了車輪扁疤激擾下輪軌沖擊噪聲的特性。翟婉明[9]將車輪扁疤等沖擊型激擾模型作為激擾源輸入到車輛-軌道耦合系統(tǒng)模型中。Johansson[10]研究了車輪多邊形化對輪軌垂向動態(tài)作用力的影響。Meinke[11]等人數(shù)值模擬了初始理想圓形車輪在激勵作用下逐漸向多邊形的發(fā)展過程。谷永磊等建立高速鐵道車輛-板式無砟軌道動力學數(shù)值分析模型,采用現(xiàn)場測試得到的高速鐵路鋼軌波磨數(shù)據(jù)作為系統(tǒng)激勵,研究不同深度的鋼軌波磨對高速鐵路輪軌相互作用、車輛運行穩(wěn)定性的影響[12],并研究了高速鐵路軌道結(jié)構(gòu)振動特性與鋼軌波磨間的內(nèi)在關系以及軌道振動特性對鋼軌波磨影響的機理[13]。陳光雄等基于摩擦系統(tǒng)自激振動理論[14-15]建立了由輪對和鋼軌組成的輪軌系統(tǒng)摩擦自激振動模型,使用該模型對輪軌系統(tǒng)進行了摩擦自激振動發(fā)生趨勢的仿真分析[16]。

本文為了探究輪軌垂向動態(tài)相互作用,把軌道考慮成離散支撐的連續(xù)的鐵木辛科梁,在鋼軌上移動車輪通過一系懸掛連接一質(zhì)量塊作為簡單的車輛模型,建立了輪軌垂向耦合振動模型。通常,相比頻域范圍,鋼軌的時域振動模型需要更多的計算工作量[17],所以對鋼軌在移動載荷作用下的運動方程進行拉普拉斯變換,用格林函數(shù)法求解其在頻域上的響應解,再通過反傅里葉變換求得鋼軌在時域內(nèi)的振動響應。用非線性赫茲接觸理論和數(shù)值迭代方法分析計算不同工況下車輪扁疤、多邊形化缺陷和鋼軌波磨缺陷引起的輪軌垂向動態(tài)沖擊力。

1輪軌垂向動態(tài)響應的求解

輪軌垂向耦合振動模型包括車輛子系統(tǒng)和軌道子系統(tǒng)。車輛子系統(tǒng)由兩個質(zhì)量塊(1/8車體和1/4構(gòu)架)、一系懸掛、二系懸掛和車輪組成的車輛簡化模型,軌道子系統(tǒng)由鋼軌、扣件、軌枕和道砟組成,輪軌垂向耦合振動的物理模型如圖1所示。車輛子系統(tǒng)以速度V沿著軌道縱向正方向前進,其中,Mc是1/8車體質(zhì)量,Mt是1/4構(gòu)架質(zhì)量,Mw是車輪質(zhì)量,并且分別只考慮它們各自的垂向位移自由度(Zc,Zt和Zw),K1和C1分別是一系懸掛系統(tǒng)的剛度和阻尼,K2和C2分別是二系懸掛系統(tǒng)的剛度和阻尼。在軌道子系統(tǒng)中,Kpv和Cpv分別為扣件的剛度和阻尼,Kbv和Cbv分別為道砟的剛度和阻尼,Kpα和Cpα分別為扣件的回轉(zhuǎn)剛度和回轉(zhuǎn)阻尼,si是第i個軌枕到參考點的距離,si+1是第i+1個軌枕到參考點的距離,兩者間距為軌距d,并且只考慮鋼軌垂向位移和橫截面旋轉(zhuǎn)位移兩個自由度,考慮軌枕的垂向位移和點頭位移兩個自由度。其他物理參數(shù)為:楊氏模量E,剪切模量G,密度ρ,鋼軌的單位長度質(zhì)量mr,鋼軌的橫截面積S,截面慣性矩I,剪切系數(shù)κ,軌枕質(zhì)量之半Ms,軌枕慣性質(zhì)量矩之半Js。假設車輛子系統(tǒng)和軌道子系統(tǒng)都是對稱的,所以只研究整個輪軌垂向耦合振動的一側(cè)。

車輛子系統(tǒng)的運動方程[18]:

(1)

其中,{z}=[zc(t)zt(t)zw(t)]T,{p}=[0 0P0-P(t)]T,P0是軸重之半,P(t)為輪軌垂向接觸力,M為質(zhì)量矩陣,C為阻尼矩陣,K為剛度矩陣。

圖1輪軌垂向耦合振動模型

軌道子系統(tǒng)的運動方程:

[-P(t)δ(x-Vt) 0]T

(2)

(3)

δ(x-ξ)[1 0]T

(4)

(5)

其中,Csv=Cpv+Cbv,Ksv=Kpv+Kbv。

由式(4)和式(5)可以得到:

(6)

(7)

其中

格林函數(shù)的Γ矩陣是以下等式的解:

(8)

其中

Γ=Γ(x,x′,ω)=diag(Γw(x,x′,ω),

Γθ(x,x′,ω)),Γw(x,x′,ω)=Γθ(x,x′,ω)=

(9)

2輪軌垂向沖擊力的計算

(10)

對于不同的輪軌缺陷,輪軌相對位移Δz會有不同的表達形式。

2.1車輪缺陷引起的輪軌相對位移

主要研究車輪扁疤和車輪多邊形化兩種車輪缺陷。車輪扁疤主要分為新扁疤和舊扁疤,車輪扁疤和車輪多邊形化的示意圖如圖2所示。車輪的新扁疤和舊扁疤的深度都為e,車輪新扁疤的寬度為l0,車輪舊扁疤的寬度為l,則車輪的新扁疤和舊扁疤分別產(chǎn)生的輪軌相對位移Δzn和Δzr為:

(11)

(12)

車輪多邊化導致的輪軌相對位移Δzp為:

(13)

其中,Δr0是車輪多邊形化的幅值,n為車輪多邊形的階數(shù),Rw是車輪的名義滾動圓半徑。

圖2車輪缺陷

2.2鋼軌缺陷引起的輪軌相對位移

只考慮鋼軌的波磨缺陷導致的輪軌相對位移的變化,把鋼軌波磨考慮成正弦函數(shù)形式r=r0sinkx,其中r0是鋼軌波磨的幅值,k為波數(shù),車輪在含有波磨的鋼軌上的不同位置,其對應產(chǎn)生的輪軌相對位移也不同,如圖3所示。

圖3車輪在波磨鋼軌上的不同位置

車輪運行在波磨鋼軌的不同位置會產(chǎn)生不同的輪軌相對位移Δzc,用分段函數(shù)表示Δzc,

(14)

3算例與分析

在數(shù)學軟件Maple里對上述各種缺陷產(chǎn)生的輪軌垂向沖擊力進行了編程計算,其中計算參數(shù)如下:E=210GPa,G=81GPa,κ=0.4,ρ=7850kg/m3,mr=60kg/m,S=7.69×10-3m2,I=30.55×10-6kgm2,Kpv=50×106N/m,Kpα=1×107Nm/rad,Cpv=2×105Ns/m,Cpα=1×104Nms/rad,Kbv=1×109N/m,Cbv=1×106Ns/m,Ms=129kg,Js=1.28kgm2,d=0.6m,s1=0.3m,V=60m/s,P0=70kN,Mc=4000kg,Mt=1800kg,Mc=1000kg,K1=1.176×106N/m,C1=1.96×104Ns/m,K2=0.2205×106N/m,C2=9.8×103Ns/m,Rw=0.46m。為了更加清楚地研究車輪缺陷和鋼軌缺陷對輪軌垂向動態(tài)相互作用的影響,不考慮軌道不平順譜對輪軌垂向沖擊力的影響。仿真結(jié)果如圖4~圖7所示。

圖4不同新扁疤深度引起的輪軌沖擊力

圖5不同舊扁疤深度引起的輪軌沖擊力

圖6不同車輪多邊形化幅值導致的輪軌沖擊力

圖7不同鋼軌波磨幅值導致的輪軌沖擊力

分別計算了車輪新扁疤和舊扁疤在0.5mm、1.0mm、1.5mm和2.0mm四種不同深度引起的輪軌垂向沖擊力,計算軌道長度為60m,時間步長為0.05ms。圖4為四種不同車輪新扁疤深度在第37號軌枕和第38號軌枕之間引起的輪軌垂向沖擊力,為了方便分析對比,圖5為四種不同車輪舊扁疤深度在第37號軌枕和第38號軌枕之間引起的輪軌垂向沖擊力,從圖4和圖5可以看出,隨著車輪扁疤深度的增加,輪軌垂向沖擊力最大值也隨之增加,車輪新扁疤的數(shù)學模型會引起車輪滾過新扁疤區(qū)域后的輪軌垂向沖擊力在輪軌靜載荷附近的振蕩,最后再收斂到輪軌靜載荷,而車輪舊扁疤的數(shù)學模型則不會引起輪軌垂向沖擊力在輪軌靜載荷附近的振蕩,從而說明車輪舊扁疤更能精確地模擬輪軌垂向沖擊力的變化。

對車輪多邊形化在0.05mm、0.10mm和0.15mm三種不同幅值下的輪軌垂向沖擊力進行了計算,時間步長為0.05ms。圖6為車輪多邊形化的不同幅值在第85號軌枕和第95號軌枕之間導致輪軌沖擊力的變化。從圖中可以看出,同一速度下的輪軌垂向沖擊力的周期一致,但是輪軌垂向沖擊力的幅值是與車輪多邊形化的幅值有關,隨著車輪多邊形化的幅值的增加,輪軌垂向沖擊力也會相應的增大,基本上呈正比例增加,并在輪軌靜載荷處上下波動。

對鋼軌波磨的不同幅值引起的輪軌沖擊力變化進行了計算,分別計算了鋼軌波磨幅值為20μm、40μm和60μm三種不同幅值下輪軌沖擊力的變化,鋼軌波磨的波長為0.05m,計算時間步長為0.1ms。圖7為鋼軌波磨的不同幅值在第37號軌枕和第38號軌枕之間引起的輪軌沖擊力變化。從圖中可以看出,由于鋼軌波磨的波長是一致的,所以鋼軌垂向沖擊力的變化周期是一致的,但不同的鋼軌波磨幅值會產(chǎn)生不同的輪軌沖擊力的大小,類似于不同幅值的車輪多邊形化引起的輪軌沖擊力,隨著鋼軌波磨幅值的增加,輪軌垂向沖擊力也隨之增大,其輪軌沖擊力的最大值基本上與三種鋼軌波磨幅值呈線性關系,并在輪軌靜載荷處上下波動。圖7與圖6的不同之處是,圖7的輪軌垂向沖擊力變化并不呈現(xiàn)光滑的諧波變化,這可能與鋼軌波磨的波長大小有關,另外,即使兩種缺陷的幅值相同,兩種缺陷產(chǎn)生的最大輪軌垂向力也不同,鋼軌波磨會產(chǎn)生相對比較大的輪軌垂向沖擊力,因此鋼軌的波磨相比車輪多邊形化會對輪軌垂向沖擊力產(chǎn)生更大的影響。

4結(jié)論

(1) 建立了輪軌垂向耦合振動模型,把鋼軌考慮成離散支撐的連續(xù)的鐵木辛科梁,用格林函數(shù)法求解其在頻域上的響應解,再通過反傅里葉變換求得鋼軌在時域內(nèi)的振動響應,可以從高頻結(jié)構(gòu)振動的角度研究輪軌垂向相互作用,克服了傳統(tǒng)動力學中只考慮輪軌振動的低頻成分,忽略中高頻振動影響的缺陷。

(2) 隨著車輪扁疤深度的增加,輪軌垂向沖擊力的最大值也隨之線性地增加,車輪新扁疤的數(shù)學模型會引起車輪滾過新扁疤區(qū)域后的輪軌垂向沖擊力在輪軌靜載荷附近的振蕩,最后再收斂到輪軌靜載荷,而車輪舊扁疤的數(shù)學模型則不會引起這種振蕩,從而說明車輪舊扁疤更能精確地模擬輪軌垂向沖擊力的變化。

(3) 相比車輪多邊形化,鋼軌的波磨則會對輪軌垂向沖擊力產(chǎn)生更大的影響,也就是說,即使兩種缺陷的幅值相同,各自導致的輪軌垂向力不同,鋼軌波磨缺陷會產(chǎn)生相對比較大的輪軌垂向沖擊力。

[1] NIELSEN J C O,IGELAND A.Vertical dynamic interaction between train and track-influence of wheel and track imperfections[J].Journal of Sound and Vibration,1995,187(5):825-839.

[2] NIELSEN J C O,OSCARSSON J.Simulation of dynamic train-track interaction with state-dependent track properties[J].Journal of Sound and Vibration,2004,275(3):515-532.

[3] BAEZA L,RODA A,NIELSEN J C O.Railway vehicle/track interaction analysis using a modal substructuring approach[J].Journal of Sound and Vibration,2006,293(1-2):112-124.

[4] WU T X,THOMPSON D J.A hybrid model for the noise generation due to railway wheel flats[J].Journal of Sound and Vibration,2002,251(1):115-139.

[5] WU T X,THOMPSON D J.On parametric excitation of the wheel/track system[J].Journal of Sound and Vibration,2004,278(4):725-747.

[6] WU T X,THOMPSON D J.The effects of track non-linearity on wheel/rail impact[J].Proceedings of the Institution of Mechanical Engineers,Part F:Journal of Rail and Rapid Transit,2004,218(1):1-15.

[7] BAEZA L,RODA A,CARBALLEIRA J,et al.Railway train-track dynamics for wheelflats with improved contact models[J].Nonlinear Dynamics,2006,45(3):385-397.

[8] 楊新文,翟婉明.車輪扁疤激起的輪軌沖擊噪聲機理分析[J].振動與沖擊,2009,28(8):46-49.

[9] 翟婉明.車輛-軌道耦合動力學[M].3版.北京:中國科學出版社,2007.

[10] JOHANSSON A,ANDERSSON C.Out of round railway wheels-a study of wheel polygonalization through simulation of three-dimensional wheel-rail interaction and wear[J].Vehicle System Dynamics,2005,43(8):539-559.

[11] MEINKE P,MEINKE S.Polygonalization of wheel treads caused by static and dynamic imbalances[J].Journal of Sound and Vibration,1999,227(5):979-986.

[12] 谷永磊,趙國堂,金學松,等.高速鐵路鋼軌波磨對車輛-軌道動態(tài)響應的影響[J].中國鐵道科學,2015,36(4):27-31.

[13] 谷永磊,趙國堂,王衡禹,等.軌道振動特性對高速鐵路鋼軌波磨的影響[J].中國鐵道科學,2016,37(4):42-47.

[14] CUI X L,CHEN G X,YANG H G,et al.Study on rail corrugation of a metro tangential track with Cologne-egg type fasteners[J].Vehicle System Dynamics,2016,54(3):353-369.

[15] HU W P,WANG P,CHEN G X,et al.Experimental study on corrugation of a sliding surface caused by frictional self-excited vibration[J].Tribology Transactions,2016,59(1):8-16.

[16] 陳光雄,崔曉璐,王科.高速列車車輪踏面非圓磨耗機理[J].西南交通大學學報,2016,51(2):244-250.

[17] PIERINGER A,KROPP W,NIELSEN J C O.A time domain model for wheel/rail interaction aiming to include non-linear contact stiffness and tangential friction[J].Notes on Numerical Fluid Mechanics and Multidisciplinary Design,2008,99:285-291.

[18] 張寶安,陸正剛,王恒亮.柔性輪對的輪軌接觸參數(shù)分析[J].同濟大學學報:自然科學版,2013,41(4):577-582.

Investigation on Wheel-Rail Vertical Dynamic Interaction Due to Wheel-rail Imperfections Based on Green’s Function Method

ZHANGBaoan1a,ZHENGJing1b,ZHOUHechao2

(1a.National Engineering Research Center for High-speed EMU; 1b. R&D Center, CRRC Sifang Co.,Ltd, Qingdao 266000, China; 2.Institute of Railway & Urban Rail Transit, Tongji University, Shanghai 201804, China)

The wheel-rail vertical coupled vibration model was established, in which the rail was considered as a continuous Timoshenko beam discretely supported. The lumped masses of wheel, a quarter of bogie frame and one eighth of carbody were connected using primary and secondary suspensions respectively, and this subsystem was taken as a simplified vehicle model. In order to efficiently and accurately investigate the effects of wheel and rail defects on vertical dynamic interaction between wheel and rail, Laplace transform was implemented for the motion equations of the rail under a moving load, and those motion equations were solved in the frequency domain by Green's function method, then inverse Fourier transform was carried out to obtain vibration response of rail in the time domain. Wheel-rail vertical impact forces due to wheel flat, wheel polygonalization and rail corrugation in various cases were calculated by numerical iterations according to the nonlinear Hertzian contact theory. Results show that the vertical wheel-rail dynamic response can be solved by the Green function method accurately and efficiently, and it provides a theoretical foundation for the research of wheel-rail interaction in the range of high frequency.

Timoshenko beam; Green’s function; vertical wheel-rail impact force; wheel flat; wheel polygonalization; rail corrugation

2017-01-06

中國鐵路總公司科技研究開發(fā)計劃項目(2016G008-F);上海市浦江人才計劃(16PJ1409500)

張寶安(1987-),男,山東日照人,工程師,博士,主要從事車輛動力學方面的研究,(E-mail)zhangbaoan@cqsf.com

1673-1549(2017)02-0026-06

10.11863/j.suse.2017.02.06

U271.91

A

猜你喜歡
振動
振動的思考
科學大眾(2023年17期)2023-10-26 07:39:14
某調(diào)相機振動異常診斷分析與處理
大電機技術(2022年5期)2022-11-17 08:12:48
振動與頻率
天天愛科學(2020年6期)2020-09-10 07:22:44
This “Singing Highway”plays music
具非線性中立項的廣義Emden-Fowler微分方程的振動性
中立型Emden-Fowler微分方程的振動性
基于ANSYS的高速艇艉軸架軸系振動響應分析
船海工程(2015年4期)2016-01-05 15:53:26
主回路泵致聲振動分析
UF6振動激發(fā)態(tài)分子的振動-振動馳豫
計算物理(2014年2期)2014-03-11 17:01:44
帶有強迫項的高階差分方程解的振動性
主站蜘蛛池模板: 好紧好深好大乳无码中文字幕| 日韩欧美中文字幕在线韩免费| 亚洲无码高清一区| 视频一区视频二区日韩专区 | 欧美性天天| 青草精品视频| 国产幂在线无码精品| 欧亚日韩Av| 国产不卡一级毛片视频| 十八禁美女裸体网站| 九色91在线视频| 九色综合视频网| 亚洲人成在线免费观看| 一区二区日韩国产精久久| 青青青视频免费一区二区| 欧美日韩v| 国产三级a| 天天躁夜夜躁狠狠躁躁88| 99re热精品视频中文字幕不卡| 成人在线观看不卡| 三级国产在线观看| 直接黄91麻豆网站| 国产午夜小视频| 99在线免费播放| 国产在线精彩视频二区| 国产激情无码一区二区APP| 国产成人亚洲无码淙合青草| 日韩无码视频网站| 亚洲av片在线免费观看| 91麻豆久久久| 国产综合在线观看视频| 亚洲国产成人无码AV在线影院L| 亚洲天堂首页| 尤物午夜福利视频| 久久国产精品77777| 91视频青青草| 香港一级毛片免费看| 亚洲色无码专线精品观看| 手机在线国产精品| 她的性爱视频| 亚洲精品成人片在线播放| 久久国产精品电影| 91破解版在线亚洲| 欧美激情伊人| 特级毛片免费视频| 91亚洲精选| 久久精品娱乐亚洲领先| 免费啪啪网址| vvvv98国产成人综合青青| 婷五月综合| 久久婷婷五月综合色一区二区| 国产h视频在线观看视频| 亚洲性日韩精品一区二区| 国产成人AV综合久久| 日韩中文欧美| 制服无码网站| 婷婷中文在线| 国产亚洲视频免费播放| 亚洲精品视频在线观看视频| 98超碰在线观看| 亚洲精品天堂在线观看| 日本午夜影院| 在线国产欧美| 久久精品人人做人人| 91免费在线看| 区国产精品搜索视频| 国产成人成人一区二区| 不卡色老大久久综合网| 人人妻人人澡人人爽欧美一区| 91av成人日本不卡三区| 美女扒开下面流白浆在线试听| 午夜久久影院| 欧美日韩国产综合视频在线观看| 天天摸夜夜操| 欧美黄网站免费观看| 亚洲熟女中文字幕男人总站| 国产91蝌蚪窝| 中文字幕啪啪| 巨熟乳波霸若妻中文观看免费| 国外欧美一区另类中文字幕| 国产女人在线观看| 免费无遮挡AV|