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

考慮接觸剛度的含間隙鉸接副動態磨損分析

2016-12-24 02:26:35朱愛斌何勝利鄒超陳渭
西安交通大學學報 2016年5期
關鍵詞:模型

朱愛斌,何勝利,鄒超,陳渭

(西安交通大學現代設計及轉子軸承系統教育部重點實驗室,710049,西安)

?

考慮接觸剛度的含間隙鉸接副動態磨損分析

朱愛斌,何勝利,鄒超,陳渭

(西安交通大學現代設計及轉子軸承系統教育部重點實驗室,710049,西安)

針對大多數含間隙鉸接副的磨損計算都非常復雜,精度和效率不能兼得,且很少能夠與含間隙鉸接副系統的動力學分析動態結合起來考慮表面接觸剛度對磨損影響的問題,提出了一種新的含間隙鉸接副磨損分析方法。基于無質量桿-彈簧阻尼模型建立考慮接觸剛度影響的含間隙鉸接副系統動力學分析模型,利用非對稱Winkler彈性基礎模型計算接觸壓力分布,釆用Archard磨損理論計算接觸表面磨損量,對接觸表面輪廓實時更新得出含間隙鉸接副的動態磨損量。分析結果表明,在表面接觸剛度較小時含間隙鉸接副的動態磨損嚴重,且隨系統轉速的變化呈現出不同的變化趨勢,反映了不同接觸剛度下含間隙鉸接副的動態磨損趨勢。該方法計算精度和計算效率較高且充分考慮了表面接觸剛度對含間隙鉸接副動態磨損的影響,對含間隙鉸接副系統的設計和動力學分析及磨損預測具有一定的指導意義。

接觸剛度;彈簧阻尼;非對稱Winkler彈性模型;動態磨損

工程中機械系統摩擦學與動力學行為是密切聯系且相互耦合的,即磨損會造成機構的幾何形狀及尺寸變化,引起機器動力學性能變化,如振動加劇這種性能的變化又反過來加重磨損[1]。由于制造誤差、磨損、材料變形和缺陷的存在,使得間隙不可避免地存在于各種運動副中。如人工關節中的間隙鉸接副[2],工業機器人中機械臂柱鉸、滑移鉸,現代航天器上的太陽能電池板、大型天線等空間結構上的鉸接結構等。1958年美國第一顆衛星“探險者一號”由于外伸天線的柔性未考慮,其剛柔耦合的動力學作用使得衛星翻轉而導致發射失敗[3]。工程機械中因鉸接副間隙過大、接觸碰撞嚴重、銷軸軸套的磨損導致工作裝置出現振動過大、噪聲異常等故障[4]。

含間隙鉸接副的磨損是一個復雜的動態過程,早期的研究多集中于摩擦副的接觸建模及磨損計算模型等方面。Kougut針對粗糙接觸表面的接觸壓力求解建立了基于有限元的改進彈塑性模型[5],P?dra分析構建了考慮接觸表面形貌的有限元磨損計算模型[6]。傳統的含間隙鉸接副接觸表面壓力分析多基于經典的Hertz接觸模型[7]。鉸接副一般為高共形、小間隙接觸,而Hertz模型是基于球面接觸提出的非共形接觸理論。Liu等已證實在間隙較小時Hertz計算結果與真實值偏差很大[8]。Lankarani等提出了既考慮接觸的非線性,又包含能量損失的較全面的接觸模型[9],但該模型求解過程復雜。

王國慶等結合間隙運動副碰撞接觸的非線性彈簧阻尼模型,提出了邊界潤滑條件下的間隙副接觸力模型[10],但沒有考慮到在高共形率接觸條件下,由Hertz模型發展而來的接觸力計算公式誤差較大的問題。林碧華基于無質量桿模型的含間隙鉸接副系統動力學分析,對含間隙鉸接副接觸表面的動態磨損量進行了計算[11],但忽略了接觸表面彈性的影響。宿月文等采用Winkler彈性基礎模型計算了接觸壓力分布[12],但傳統Winkler彈性接觸模型是將間隙鉸接副作為理想圓處理,其所得接觸壓力分布結果的精度相對比較差。

本文基于無質量桿-彈簧阻尼(Mass-less-Link and Spring-Damper,MLSD)模型[13],在考慮了接觸表面不同接觸剛度影響的基礎上,建立了含間隙鉸接副系統的動力學分析模型;利用非對稱Winkler彈性基礎模型計算了接觸表面的壓力分布,解決了高共形率的接觸壓力計算問題;采用Archard理論進行磨損計算,通過對接觸表面輪廓的實時更新,實現了不同接觸剛度下含間隙鉸接副的動態磨損分析。

1 含間隙鉸接副的動態磨損

1.1 考慮不同接觸剛度系數影響的動力學分析計算

含間隙鉸接副動態磨損算法的流程圖如圖1所示。基于MLSD的含間隙曲柄滑塊機構動力學分析模型如圖2所示。

圖1 不同接觸剛度影響下含間隙鉸接副動態磨損算法流程

圖2 含間隙的曲柄滑塊機構MLSD分析模型

基于MLSD模型,將滑塊與連桿間的間隙用長度為ξ、方向為γ的無質量桿CD表示,CD的徑向剛度系數為k,徑向阻尼系數為μ。為保證運動中銷軸與軸套始終接觸,在結構右端附加一個彈簧系統,彈簧彈性系數為ksp。假定曲柄做勻速轉動,角速度為ω,初始間隙為e。

機構位置約束方程為

(1)

機構動力學約束方程為

(2)

對間隙桿CD有

(3)

利用攝動法[14]的求解思路對φ2、γ、x3、FN做如下修正

(4)

式中:a為修正因子,是一個無窮小量。由于間隙很小,通常在10-4~10-5m范圍間,因此可以認為e=a為定值。

將式(4)修正值分別代入式(1)和式(2),通過求解無間隙標稱機構得到相應的含間隙機構的動力學響應數據,并對理想化數據進行修正處理后求得

(5)

式中:f1(m1,m2,m3,l1,l2,φ1,ω)和f2(m1,m2,m3,l1,l2,φ1,ω)是關于已知量的函數表達式。

根據線性彈簧阻尼接觸力計算公式,可以得到

(6)

(7)

再次將γ*、ξ的值帶入式(1)和式(2),求得表面接觸力FN和方位角γ。

1.2 接觸壓力分布計算

Winkler彈性基礎模型基于離散化的思想處理接觸問題,且不局限于高共形率的接觸要求,因此用該模型的算法結果必然比基于Hertz接觸理論的結果具有更高的精度[12]。該模型把接觸壓力用彈簧單元的接觸變形來表示,其計算效率要高于Hertz接觸模型。Winkler接觸模型的接觸壓力與變形量的關系為

(8)

式中:pi為接觸壓力;δi為法向彈性變形;Li為彈性層厚度;Ew為彈性層的復合彈性模量。

在接觸表面上所有彈簧合外力大小即為外載荷大小,即

(9)

式中:βi為第i個彈簧與FN的夾角;Ai為第i個彈簧處的接觸面積;n為接觸區域彈簧總數。

非對稱Winkler彈性基礎模型的接觸示意如圖3所示,其中接觸區域為劣弧C1C2,接觸區域的彈簧壓力總和Ftotal=FN。

圖3 非對稱Winkler彈性基礎模型接觸示意圖

本文采用文獻[15]提出的非對稱Winker彈性基礎模型進行壓力分布計算,基本求解思路如下。

(1)假設壓縮中心變形量為δi。

(2)鑒于磨損引起軸套內壁的非圓特性,圓弧函數已不再適用,所以利用銷軸外輪廓以及軸套內壁幾何結構的離散化數據,分別計算載荷方向左右兩側的受壓彈簧數量及變形量。

(3)根據受壓彈簧數及變形量,求解接觸表面的接觸壓力分布。

(4)計算受壓彈簧合力Ftotal,并與FN進行對比,如果二者平衡,則輸出接觸壓力分布大小,否則返回思路(1),修正初始變形量δi,繼續迭代,直到二者平衡。

1.3 磨損計算

Archard磨損計算式[16]為

(10)

式中:K為相對磨損系數;H為材料布氏硬度;F為接觸表面法向載荷;s為接觸表面相對滑移距離;V為磨損體積。

工程中考慮更多的是材料的法向侵入深度,因此對上式兩邊取接觸面積的微分得

(11)

式中:h為材料接觸表面的法向磨損深度;p為接觸表面的接觸壓力。

計算中通常將材料的磨損過程離散化為若干個過程,將每個過程磨損量連續累加,得出材料接觸表面的法向侵入深度。取式(11)對時間的微分得

dhi=Kdpidsi

(12)

(13)

式中:Kd為接觸材料的線性磨損系數,一般通過銷盤磨損實驗得出;dhi為每個離散化過程的單步磨損量。

在材料的磨損過程中,單步的磨損量非常小,其對接觸界面輪廓的影響不大,因此引入一個磨損加速因子M,把連續若干步的磨損深度計算視為單步磨損深度的若干倍,即

hi+1=hx-M+1+Mhx=hx-M+1+MKdpidsi

(14)

2 數據模擬分析

基于無質量桿-彈簧阻尼模型建立考慮不同接觸剛度影響的含間隙鉸接副系統動力學分析模型,利用非對稱Winkler彈性基礎模型計算接觸壓力分布,釆用Archard磨損理論計算接觸表面磨損量,對接觸表面輪廓實時更新得出含間隙鉸接副的動態磨損量,由此提出了一種考慮表面接觸剛度影響的含間隙鉸接副動態磨損分析方法。曲柄滑塊機構動力學分析模型見圖2,基本參數如表1、表2所示。

表1 曲柄滑塊機構基本參數

表2 銷軸與軸套基本參數

在文獻[11]的含間隙鉸接副磨損計算工作中,基于間隙桿模型的Winkler模型下的磨損計算相比有限元分析模型和實驗數據具有較好的精度。在相同的磨損周期計算過程中,基于間隙桿模型的Winkler模型的計算方法用時幾乎為有限元方法的1/100,計算效率較高。考慮接觸表面不同接觸剛度的本文算法與文獻[11]算法在不同狀態下的動態磨損結果對比如圖4所示。

(a)120 r/min

(b)240 r/min圖4 不同接觸剛度、10 000轉時本文與文獻[11]算法結果對比

結合圖4a和表3可以看出,當k=1×107N/m時,本文算法結果與文獻[11]基本接近,磨損曲線具有相同的變化趨勢和一致的磨損范圍。當k=1×104N/m時,本文算法的磨損峰值與文獻[11]有較大差異,本文磨損峰值小于文獻[11],同時磨損區域有所增大,磨損質量相對減小。結合圖4b和表4可以看到,當k=1×107N/m時,隨著轉速的提高,軸套表面的磨損量逐漸增大。在磨損的起始階段,文獻[11]算法的磨損結果偏大,而在最大磨損發生之后的區域里,本文算法略大。原因有兩方面:一是本文算法動力學分析計算基于了MLSD模型,同時考慮了運動副表面的基礎剛度影響,使得銷軸在最初與軸套接觸階段,隨著侵入深度的不斷增大,軸套表面發生彈性變形而產生彈性滯回,相對延長了作用時間,接觸表面的接觸力相對減小,磨損量相對減小,而在由最大侵入深度逐漸退出磨損區域時,彈性勢能釋放,接觸表面間隙成為影響接觸力大小的主導因素,彈性形變引起的實際接觸間隙增大導致接觸力增大,磨損量相對增大;二是本文考慮了磨損發生后軸套內表面的非圓特性,并且采用了非對稱Winkler彈性基礎模型來計算接觸壓力分布。由表4可以看到,本文算法的磨損峰值、磨損區域及磨損量遠大于文獻[11],接觸表面的彈性對磨損計算結果的影響較大,表明考慮表面接觸剛度的算法具有較廣的適用范圍。

表3 120 r/min下磨損計算結果對比

表4 240 r/min下磨損計算結果對比

由表3、表4還可以看到,考慮了表面接觸剛度,計算效率更高,這種時效性將在工程周期磨損計算中更好地體現出來。隨著系統轉速的提高,接觸表面磨損深度增大,磨損區域變寬,磨損質量增加。

轉速為120 r/min時,減小接觸表面剛度系數,表面磨損峰值減小,磨損量減小,磨損區域增大;轉速為240 r/min時,減小接觸表面剛度系數,表面磨損峰值增大,磨損量增加,磨損區域增大。由此可見,磨損區域的擴展是由于表面彈性影響所致,較大彈性的接觸表面,磨損區域較大,這種現象隨轉速的提高而越為明顯,與實際動態磨損現象相符。不同轉速下改變接觸表面剛度系數,磨損峰值和磨損量變化的原因如下:在相同外載荷下,同一剛度系數的接觸表面的彈性形變量(微增間隙)相同;轉速較低時,相對作用時間較長,尤其是達到磨損峰值的作用時間相對延長,滯回作用較明顯,接觸反力較小;微增間隙對接觸力的影響小于表面滯回作用,表面滯回作用的影響占據主導地位,表面接觸力相對減小,磨損峰值和磨損量減小;持續提高轉速,當轉速超過某一臨界值時,較高轉速下較短的接觸時間,尤其是達到磨損峰值的接觸時間相對縮短,接觸反力增大,表面滯回作用不明顯,微增間隙對磨損的影響超過了表面滯回的作用,微增間隙的影響占據了主導地位,接觸力相對增大,磨損峰值和磨損量增加。

不同轉速下接觸表面剛度系數對接觸表面的影響如圖5、表5和表6所示。可以看出:轉速為120 r/min、k≥1×106N/m時,表面磨損變化趨于穩定;轉速為240 r/min、k≥1×106N/m時,表面磨損變化趨于穩定,磨損峰值、磨損量和磨損區域不再發生變化。因此,在一定的轉速和不變的外載荷下,相同的磨損周期內持續增加接觸表面剛度系數,接觸表面的動態磨損量逐漸趨于穩定。文獻[11]適用于表面接觸剛度系數較大時的情形,而本文算法能夠統籌考慮不同接觸表面剛度系數的影響,且具有較高的精度。

(a)120 r/min

(b)240 r/min圖5 不同接觸表面剛度系數下磨損結果對比

k/N·m-1磨損深度/mm磨損寬度/(°)磨損質量/g1050.06001640.03491060.06041650.03511070.06041650.0351

表6 240 r/min不同剛度系數下磨損結果

k=1×107N/m時,本文算法在不同磨損周期和不同轉速下關節軸承表面磨損量的動態對比如圖6所示。由圖6可以看出,在相同轉速下,隨著磨損周期的延長,關節軸承表面的磨損區域發生了動態擴展。在該磨損區域,隨磨損周期的不斷拓展,銷軸與軸套表面振蕩周期延長,從而加速了軸承的失效,嚴重影響了關節軸承的使用壽命。

(a)240 r/min

(b)10 000轉圖6 不同狀態下本文算法結果

由表7磨損數據得,接觸表面的磨損質量和磨損深度隨磨損周期的延長而等比增長。分析表8數據可以看出,在相同的磨損周期中,轉速越高,磨損得越深,磨損區域越寬,磨損質量越大。對于高速運轉而言,關節軸承磨損更快,壽命周期更短。因此,含間隙鉸接副的磨損指標隨轉速和磨損周期的增加而不斷增大。

表7 240 r/min下不同仿真周期的磨損計算結果

表8 10 000轉時不同轉速的磨損計算結果

3 結 論

含間隙鉸接副的動態磨損計算是含間隙多體動力學和摩擦學耦合的交互影響的復雜工程問題,其蘊含了摩擦學所固有的宏、微觀多尺度特性,跨越了瞬態接觸和緩慢磨損的多時間尺度。

本文從系統耦合集成的角度融合各領域、多尺度的知識構建模型,基于無質量桿-線性彈簧阻尼模型,在考慮鉸接副表面接觸剛度影響的基礎上,對含間隙鉸接副系統進行了動力學分析,由此簡化了含間隙系統強非線性動力學微分方程組的求解過程,提高了算法的效率。同時,利用非對稱Winkler彈性基礎模型、采用離散化的思想,實現了高共形率接觸表面的非對稱壓力的分布計算,釆用Archard磨損理論計算含間隙運動副表面的磨損量、利用磨損量對間隙大小實時修正,得出不同接觸剛度下的含間隙鉸接副系統的動態磨損量變化趨勢。通過對比不同接觸剛度下接觸表面動態磨損趨勢的研究分析,表明本文算法可以反映不同接觸剛度下含間隙鉸接副表面的動態磨損趨勢。該研究方法在有效提高含間隙鉸接副動態磨損計算精度的同時,保證了計算的時效性,真實地反映了含間隙鉸接副接觸表面的動態磨損狀況,對工程機械設計和動力學分析及磨損預測具有重要理論意義和應用價值。

[1] 宿月文, 陳渭, 朱愛斌, 等. 鉸接副磨損與系統動力學行為耦合的數值分析 [J]. 摩擦學學報, 2009, 29(1): 50-54. SU Yuewen, CHEN Wei, ZHU Aibin, et al. Numerical analysis of wear coupling with dynamics behavior of machine [J]. Tribology, 2009, 29(1): 50-54.

[2] 王成燾. 天然與人工關節中的摩擦學問題 [J]. 醫用生物力學, 2009, 24(5): 317-324. WANG Chentao. Tribology problems in natural and artificial joint [J]. Journal of Medical Biomechanics, 2009, 24(5): 317-324.

[3] 馬興瑞, 王本利. 航天器動力學: 若干問題進展及應用 [M]. 北京: 科學出版社, 2001.

[4] 何清華, 張大慶, 郝鵬, 等. 液壓挖掘機工作裝置模型及控制的試驗研究 [J]. 中南大學學報: 自然科學版, 2006, 37(3): 542-546. HE Qinghua, ZHANG Daqing HAO Peng, et al. Model and experimental research on control of hydraulic excavator’s manipulator [J]. Journal of Central South University: Science and Technology, 2006, 37(3): 542-546.

[5] KOGUT L, ETSION I. A finite element based elastic-plastic model for the contact of rough surfaces [J]. Tribology Transactions, 2003, 46(3): 383-390.

[7] JOHNSON K L. Contact mechanics [M]. Cambridge, UK: Cambridge University Press, 1985.

[8] LIU Caishan, ZHANG Ke, YANG Rei. The FEM analysis and approximate model for cylindrical joints with clearances [J]. Mechanism & Machine Theory, 2007, 42(2): 183-197.

[9] LANKARANI H M, NIKRAVESH P E. A contact force model with hysteresis damping for impact analysis of multibody systems [J]. Journal of Mechanical Design, 1990, 112(3): 369-376.

[10]王國慶, 劉宏昭, 何長安. 含間隙曲柄滑塊機構運動副動態磨損研究 [J]. 機械強度, 2006, 28(6): 849-852. WANG Guoqing, LIU Hongzhao, HE Chang’an. Dynamic wear research of clearance joint in crank-rod mechanism [J]. Journal of Mechanical Strength, 2006, 28(6): 849-852.

[11]林碧華. 含間隙鉸接副關節軸承磨損壽命預測方法與實驗研究 [D]. 西安: 西安交通大學, 2013.

[12]宿月文, 陳渭, 郭彩霞. 應用Winkler彈性基礎模型的間隙鉸接副磨損預測 [J]. 摩擦學學報, 2012, 32(4): 320-324. SU Yuewen, CHEN Wei, GUO Caixia. Wear prediction of clearance revolute joints using Winkler surface model [J]. Tribology, 2012, 32(4): 320-324.

[13]SENEVIRATNE L D, EARLES S W E, FENNER D N. Analysis of a four-bar mechanism with a radially compliant clearance joint [J]. Proceedings of the Institution of Mechanical Engineers: Part C Journal of Mechanical Engineering Science, 1996, 210(3): 215-223.

[14]張建玲, 李哲, 白師賢. 求解含間隙機構運行微分方程的攝動法 [J]. 北京工業大學學報, 1997, 23(2): 1-7. ZHANG Jianling, LI Zhe, BAI Shixian. A perturbation method of solving the differential equations of the linkages with clearances [J]. Journal of Beijing Polytechnic University, 1997, 23(2): 1-7.

[15]李培. 基于動力學與摩擦學耦合的含間隙鉸接副多體系統研究 [D]. 西安: 西安交通大學, 2015.

[16]ARCHARD J F. Wear theory and mechanism [M]. New York, USA: ASME, 1980.

[本刊相關文獻鏈接]

朱愛斌,何大勇,鄒超,等.刀具磨損圖像視差圖的非標定方法.2016,50(3):8-15.[doi:10.7652/xjtuxb201603002]

劉雨薇,張衛正,葉福浩.單個微凸體接觸過程中的閃點溫度研究.2016,50(3):90-93.[doi:10.7652/xjtuxb201603014]

蒲偉,王家序,楊榮松,等.重載下準雙曲面齒輪傳動界面潤滑機理分析.2015,49(11):55-61.[doi:10.7652/xjtuxb2015 11010]

康樂,胡欲立,張克涵.水下磁諧振式無線電能傳輸系統的分析與設計.2015,49(10):41-47.[doi:10.7652/xjtuxb201510 007]

馬立峰,馬自勇,馬立東,等.利用塑性5εt與中性層偏移的二輥矯直輥型設計方法.2015,49(7):71-81.[doi:10.7652/xjtuxb201507013]

王勇,邢建東,馬勝強,等.Fe-B合金在鋅液中沖刷與腐蝕的交互作用.2015,49(7):92-97.[doi:10.7652/xjtuxb201507 016]

朱愛斌,胡浩強,何大勇,等.采用頻域融合方法的砂輪刀具磨損三維重構技術.2015,49(5):82-86.[doi:10.7652/xjtuxb201505013]

吳建華,雷源,王剛,等.往復式冰箱壓縮機曲軸動態特性與軸承潤滑計算分析.2015,49(2):55-61.[doi:10.7652/xjtuxb201502010]

鮑崇高,宋奕僑,侯書增,等.磨料表面微氧化對cBN磨具磨削性能的影響.2015,49(2):124-129.[doi:10.7652/xjtuxb 201502021]

張艷,鮑崇高,侯書增.金屬模具表面稀土催化共晶滲硼組織及其性能.2014,48(10):96-100.[doi:10.7652/xjtuxb201410 015]

(編輯 苗凌)

Dynamic Wear Analysis of Clearance Joint Considering Contact Stiffness

ZHU Aibin,HE Shengli,ZOU Chao,CHEN Wei

(Key Laboratory of Education Ministry for Modern Design & Rotor-Bearing System, Xi’an Jiaotong University, Xi’an 710049, China)

Most of wear calculations for clearance joint are very complicated and the accuracy and efficiency can not be satisfied simultaneously, and the dynamic wear analysis rarely combines with the dynamics analysis of clearance joint system considering the effect of different surface contact stiffness. A dynamic wear analysis for joint clearance considering different contact stiffness is proposed. The mass-less-link and spring-damping model is adopted to conduct dynamics analysis of clearance joint system considering different contact stiffness. The nonsymmetric Winkler surface model is taken to evaluate the contact pressure distribution, and the Archard wear theory is employed to calculate the worn mass of contact surface. The surface profile renews at real time with the wear calculation result, then the dynamic wear prediction of clearance joint is acquired. The dynamic wear trend for different contact stiffness is discussed. The proposed strategy fully considers the effects of contact stiffness on the dynamic wear of clearance joint, and reflects the dynamic wear trend of clearance joint surface for different contact stiffness with higher accuracy and efficiency.

contact stiffness; spring damping; nonsymmetric Winkler surface model; dynamic wear

10.7652/xjtuxb201605002

2015-11-26。 作者簡介:朱愛斌(1975—),男,副教授。 基金項目:國家自然科學基金資助項目(51175409)。

時間:2016-03-01

http:∥www.cnki.net/kcms/detail/61.1069.T.20160301.1018.016.html

TH117.1

A

0253-987X(2016)05-0012-07

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 久久久久青草大香线综合精品| 影音先锋亚洲无码| 久热中文字幕在线| 国产精品免费露脸视频| 狠狠亚洲五月天| 亚洲无码视频喷水| 女人18毛片水真多国产| 无码专区国产精品第一页| 久久亚洲国产视频| 久久中文字幕av不卡一区二区| h网站在线播放| 毛片在线看网站| 国产一区二区网站| 精品少妇人妻av无码久久| 国产三级精品三级在线观看| 亚洲人成网站色7799在线播放| 欧美视频在线播放观看免费福利资源| 欧美成人在线免费| 久久一本日韩精品中文字幕屁孩| 欧美黄网站免费观看| 欧美成人亚洲综合精品欧美激情| 97视频在线精品国自产拍| 尤物特级无码毛片免费| 午夜福利视频一区| 老色鬼欧美精品| 国产最新无码专区在线| 手机成人午夜在线视频| 日本免费福利视频| 亚洲人网站| 日韩欧美网址| 91丝袜美腿高跟国产极品老师| 免费看美女毛片| 成人午夜网址| 国产大片黄在线观看| 一级片一区| 成人免费午间影院在线观看| 夜夜操国产| 91无码网站| 亚洲中文字幕手机在线第一页| 亚洲综合天堂网| 国产jizzjizz视频| 狠狠综合久久| 日韩无码视频播放| 91在线激情在线观看| 欧美国产中文| 国产本道久久一区二区三区| 91精品最新国内在线播放| 国产午夜人做人免费视频中文| 少妇露出福利视频| 国产中文在线亚洲精品官网| 成人午夜免费观看| 国产精品无码AV片在线观看播放| 精品少妇人妻av无码久久| 免费中文字幕在在线不卡| 97在线免费视频| 午夜小视频在线| 亚洲色成人www在线观看| 精品偷拍一区二区| 黄色国产在线| 久久精品免费国产大片| 在线观看精品国产入口| 2021国产在线视频| 久久人妻系列无码一区| 国产乱人伦偷精品视频AAA| 亚洲人成色在线观看| 九色91在线视频| 成人无码一区二区三区视频在线观看| 激情综合五月网| 国产成+人+综合+亚洲欧美| 97综合久久| 91亚瑟视频| 日韩123欧美字幕| 高h视频在线| 五月婷婷综合网| 自拍偷拍欧美日韩| 国产超薄肉色丝袜网站| 高清不卡毛片| 亚洲国产日韩欧美在线| 亚洲乱码在线视频| 国产福利免费视频| 干中文字幕| 四虎精品国产永久在线观看|