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

考慮磨耗的鋼軌疲勞裂紋萌生壽命預測仿真

2016-05-15 07:14:10王少鋒姜俊楠
鐵道學報 2016年7期
關鍵詞:裂紋模型

周 宇, 張 杰, 王少鋒, 姜俊楠, 于 淼

(1.同濟大學 道路與交通工程教育部重點實驗室,上海 201804;2.華東交通大學 軌道交通學院,江西 南昌 330013)

滾動接觸疲勞(Rolling Contact Fagitue,RCF)裂紋和磨耗是影響鋼軌壽命的主要傷損,兩者相互影響,共同控制鋼軌壽命[1]。大多數針對疲勞裂紋[2-4]和磨耗[5-6]的研究都是分別開展的,沒有討論彼此之間的關系。已有的討論兩者相互關系的文獻是從經驗和現場觀測[7-8]或者實驗室輪輪雙盤滾動試驗[9-10]來研究,缺少理論探索。在非鐵路領域,MADGE等[11-12]采用Archard磨耗模型和臨界平面法裂紋萌生預測模型分析了磨耗對Ti-6Al-4V合金疲勞裂紋萌生和擴展的影響。LEEN等[13]通過三維有限元模型研究了典型的航空發動機花鍵聯軸器的疲勞和接觸磨耗。ZHANG等[14]基于有限元方法對人工髖關節假肢之間的磨耗-疲勞裂紋相互影響進行了分析。而鋼軌在輪軌反復作用下,其表面材料的磨耗和疲勞總是同時存在的,特別是硬質鋼軌具有磨耗小、裂紋萌生早、擴展快等特點,已經嚴重影響使用壽命,有必要建立理論模型對鋼軌的疲勞裂紋和磨耗的相互關系進行研究,探索兩者的平衡關系和減緩措施。

本文根據Archard磨耗理論和臨界平面理論分別建立鋼軌磨耗模型和疲勞裂紋萌生預測模型,通過磨耗計算、型面變化和平滑、規定磨耗量的型面迭代、裂紋萌生預測和疲勞損傷累積,將車輪通過次數引起的疲勞裂紋萌生和磨耗發展的過程結合起來,建立考慮磨耗的鋼軌疲勞裂紋萌生壽命預測方法,預測磨耗影響下的疲勞裂紋萌生特征。

1 考慮磨耗的鋼軌疲勞裂紋萌生壽命預測方法

1.1 磨耗及型面預測方法

Archard磨耗模型認為材料的磨耗體積與硬度成反比,與法向力、滑動距離成正比[15],即

Vm/D=K·N/H

( 1 )

式中:Vm為材料磨耗體積;D為滑動距離;N為輪軌法向力;H為材料的硬度;K為磨耗系數,其值由滑動距離和法向壓力決定,本文根據文獻[15-16]中的磨耗系數取平均值。

當考慮接觸斑面積及接觸斑的黏著-滑動區分布時,式( 1 )中的N由接觸應力代替,則可以計算出接觸斑內滑動區任意點的磨耗量,即磨耗深度。

接觸點與接觸斑關系見圖1。當一個車輪滾過鋼軌表面任意O點時,該點受到不同時刻輪軌接觸斑A1、A2、A3…的分別碾壓作用。假設該車輪在各個時刻產生的輪軌接觸斑為穩態形式,即其輪軌正壓力、蠕滑率、蠕滑力、接觸面積、黏著-滑動區分布等參數不變。為計算一個車輪滾過O點的累積磨耗量,將可能的接觸斑區域劃分成m×n個單元格,如圖2所示,O點的累積磨耗量就是所有接觸斑通過該點所在單元格引起的磨耗量的累加。例如,假設該點位于接觸斑中心,其累積磨耗量為經過該點所有縱向單元格磨耗量的累加,如圖2中通過O點的縱向黑色條帶。

因此,接觸斑內任意點的磨耗為

( 2 )

式中:Δzys為接觸斑第y個縱向條帶的累積磨耗量;Δzs(x,y)為單元格(x,y)處的磨耗量;m為接觸斑縱向單元格數量;n為接觸斑橫向單元格數量。

當一個轉向架的前后2個車輪均作用在一股鋼軌上時,如2個車輪作用在曲線外軌上,相應的輪軌接觸斑位置、尺寸、黏著-滑動區分布及鋼軌型面磨耗量見圖3。其中,接觸斑和磨耗深度圖的坐標系原點在軌頭中心線上,橫向為鋼軌橫斷面方向,縱向為列車運行方向。進一步,當1節車廂2個轉向架4個輪對作用時,任意橫斷面鋼軌磨耗量為4個車輪引起的累積磨耗量之和。

由圖3可知,1個車輪作用下的鋼軌磨耗量最大不到1.8×10-6mm。為提高計算效率,定義當1節車輛反復作用下1股鋼軌上的4個車輪引起的累積磨耗量之和再乘以車輪通過次數得到總磨耗量。當鋼軌最大總磨耗量(磨耗深度)小于0.04 mm時,型面不發生變化;當最大總磨耗量達到0.04 mm時,型面發生變化。因此,仿真中鋼軌型面連續變化采用每次最大磨耗量為0.04 mm的型面分段迭代實現。這時將鋼軌型面各點垂直方向減掉對應的磨耗量,形成磨耗后的鋼軌型面。

特別地,當曲線半徑小于600 m時,外輪-外軌特別是導向輪與外軌常發生2點接觸(轉向架自由內接通過曲線時)。這時,軌肩-車輪踏面接觸點按上述方法計算磨耗量,而輪緣-軌側的接觸點由于為全滑動狀態,主要引起外軌側磨和型面改變,所以根據輪緣-軌側接觸斑面積、接觸斑法向力(為輪軌橫向力在接觸斑內的分布),按庫倫摩擦理論計算出接觸斑切向力分布,并進一步計算輪緣-軌側接觸點處引起的磨耗量,作為軌側的磨耗量和磨耗控制點。

由于磨耗模型認為磨耗只發生在分布于接觸斑兩側和后部的滑動區,鋼軌型面在2個位置會引起較大磨耗,使得磨耗量直接疊加到前一個鋼軌型面上造成型面不規則,因此,采用三次插值樣條曲線法將磨耗后的型面進行平滑處理[17]。

以輪軌1點接觸為例,型面平滑控制參數見圖4。

圖4中的控制參數具體為:

(1) 控制點。起點Ps為軌頂中心;終點Pe為軌頂面垂直向下16 mm處;磨耗計算得到的2個磨耗量峰值點P1和P2。

(2) 下降量d1和d2分別為P1和P2的磨耗深度,即標準75 kg/m鋼軌型面對應位置上的垂直磨耗。

(3) 下降量所在位置L1和L2分別為P1和P2的橫坐標。

(4) 樣條曲線控制點密度。在圖4所示坐標系下,橫坐標0~34 mm范圍內控制點間隔1 mm;橫坐標34~36 mm范圍內控制點間隔0.1 mm。

由于鋼軌磨耗型面沿x軸方向的變化是無突變漸變過程,定義磨耗型面上各個控制點的下降方式為

( 3 )

相鄰2個控制點(點i和點i+1)間的曲線方程為

y=F(x)=yi+Ci,1(xi+1-xi)+

Ci,2(xi+1-xi)2+Ci,3(xi+1-xi)3

( 4 )

式中:x為控制點的橫坐標;Ci,1、Ci,2、Ci,3均為方程的系數。

令式( 4 )的二階導數F″(xi)=Si。由于F(x)的一階導數連續,即F′(xi-0)=F′(xi+0)。當相鄰2個點之間的距離相等即hi=xi-xi-1=hi+1=xi+1-xi=h時,可得

Si-1+4Si+Si+1=6(yi+1+yi-1-2yi)/h2

( 5 )

因此,式( 4 )中的系數為

( 6 )

此外,圖4中起點和終點處的一階導數應該連續以保證磨耗部分與鋼軌型面剩余部分曲線相切。當x=xi和x=x0時,從式( 5 )可得出其邊界條件為

( 7 )

式中:μi=hi/(hi+hi+1)=(xi-xi-1)/(xi+1-xi-1),λi=hi+1/(hi+hi+1)=(xi+1-xi)/(xi+1-xi-1)=1-μi,當2個相鄰點橫坐標距離相等時,μi=λi=1/2,i=1,2,…,n-1。

至此,磨耗型面的平滑可由式(4)~式(7)得到。2點接觸時,增加軌側工作邊的側磨點作為第3個磨耗控制點。

1.2 裂紋萌生壽命預測方法

對于每一個磨耗型面,其對應的疲勞裂紋預測方法從臨界平面法和能量密度法得到,其中疲勞參量為

( 8 )

式中:〈 〉為MacCauley括號,〈σmax〉=0.5(|σmax|+σmax);σmax為裂紋面上的最大正應力;Δε為裂紋面上正應變幅值;Δτ和Δγ分別為裂紋面上剪應力幅值和剪應變幅值;J為材料參數。

疲勞裂紋萌生壽命預測時根據式( 8 )中正應力部分和剪切應力部分的關系,選取不同的裂紋萌生壽命預測公式,即

FPijmax=

( 9 )

1.3 材料疲勞破壞累積方法

假設第i個型面被第i+1個型面替換前,共有ni次車輪通過,則第i個型面上第j點的無量綱的疲勞損傷為

Dij=ni/Nfij

(10)

假如型面上的第j點沒有在磨耗過程中被磨掉,根據Miner疲勞法則,當∑Dj=D1j+D2j+…+Dij=DCR=1時,認為在第i個型面的第j點上裂紋萌生,裂紋萌生壽命為

(11)

式中:m為裂紋萌生時型面因磨耗而被替換的總次數。

1.4 考慮磨耗的鋼軌疲勞裂紋萌生壽命預測流程

綜上所述,考慮磨耗的鋼軌疲勞裂紋萌生壽命預測步驟如下:

Step1以鋼軌初始型面(如標準型面)為鋼軌上道后的原始狀態(i=1),采用多體動力學軟件建立車輛-軌道動力學模型,計算相應的輪軌力、接觸斑位置和面積;采用Kalker輪軌蠕滑模型,計算輪軌法向力、蠕滑力及其在輪軌接觸斑內的分布;當外輪輪緣-外軌軌側出現第2點接觸時,用庫倫摩擦理論計算接觸斑面積內法向力與切向力分布;當型面因磨耗而發生替換后(i=i+1),同樣根據上述方法計算對應型面的輪軌接觸斑和輪軌力。

Step2采用式( 1 )、式( 2 )計算型面各點的磨耗量,當型面最大磨耗量沒有超過設定的磨耗量(本文是0.04 mm),型面不替換,同時累加車輪通過次數;當型面最大磨耗量超過設定的磨耗量,將型面各點磨耗量疊加到初始型面(或第i個型面)對應位置,采用式( 3 )~式( 7 )平滑型面,替換初始型面(或第i個型面),得到第i+1個型面(計數為i=i+1),同時得到該磨耗階段的累積車輪通過次數,從而實現型面磨耗和分段迭代。

Step3結合鋼軌初始型面(或第i個型面),采用有限元方法建立鋼軌全局模型和局部模型,全局模型用于計算局部模型的約束條件;局部模型用于施加接觸斑內應力分布和約束條件,計算出軌頭的應力應變分布。

Step5采用式(10)、式(11)計算疲勞階段損傷Dij和疲勞累積損傷∑Dj。

Step6若疲勞累積損傷∑Dj小于DCR,則各點均沒有萌生裂紋,繼續按Step1~Step5條計算軌頭各點的疲勞累積損傷;若疲勞累積損傷∑Dj等于或大于DCR,則認為在該點萌生裂紋,對應的各個磨耗階段的累積車輪通過次數之和就是裂紋萌生壽命。

考慮磨耗的疲勞裂紋萌生預測流程見圖5。

2 仿真結果

仿真計算時選用我國重載鐵路小半徑曲線常用的75 kg/m、U75V熱處理鋼軌。鋼軌的基本力學參數取自文獻[2]。鋼軌為曲線外軌,曲線半徑800 m,均衡超高,列車運行速度60 km/h。鋼軌表面摩擦系數取0.3。

2.1 型面磨耗預測結果檢驗

對現場曲線外軌型面進行了測量,測量時間分別為新軌上道和通過總重約10×107kN時,期間車輪通過次數約6×105次。將2次測量的型面數據插值平滑處理,得到通過總重為1×107kN(車輪通過次數約6×104次)時的鋼軌磨耗型面,作為實測磨耗型面。同時,采用1.1節的磨耗及型面預測方法預測車輪通過約6×104次時的鋼軌型面。預測磨耗型面與實測型面的對比及其一階導數變化趨勢見圖6。

由圖6(a)可知,預測的磨耗型面與實測型面相差很小。對這2個型面進行離散分析,相同橫坐標處的縱坐標偏差平均值為0.009 mm,各點偏差的方差為7.42×10-5。進一步分析軌距邊一側型面的一階導數,對比兩者的曲線變化趨勢,如圖6(b)所示,兩者一階導數偏差的平均值為0.002,一階導數偏差的方差為4.44×10-6,說明預測的磨耗型面與實測磨耗型面基本接近。

2.2 考慮磨耗的裂紋萌生位置及壽命

在考慮磨耗的情況下,預測得到鋼軌型面經過7次迭代,即7次型面磨耗,裂紋出現萌生。考慮磨耗的裂紋萌生壽命計算結果見表1。

表1 考慮磨耗的裂紋萌生壽命計算結果

由表1可知,隨著磨耗的增加以及型面的迭代,接觸斑內應力呈增加趨勢,且初期緩慢增加,后期快速增加,最終導致疲勞累積達到臨界值(DCR=1)。

考慮磨耗情況下鋼軌軌頭疲勞累積損傷量分布見圖7。可知,疲勞損傷主要分布在軌頂面以下12 mm、距離軌頂中心10~25 mm處。其中最大疲勞累積損傷量所處位置就是裂紋萌生的位置。考慮磨耗和不考慮磨耗2種方法預測出的外軌裂紋萌生壽命及位置見表2。

表2 考慮磨耗和不考慮磨耗情況下外軌疲勞裂紋萌生壽命及位置

工況裂紋萌生壽命/(105次)萌生位置/mm鋼軌截面橫向鋼軌截面垂向不考慮型面磨耗5.3918.432.95考慮型面磨耗2.1716.872.42現場觀測1.63~3.54(約3.0~6.5×107kN)[18]——文獻[19]2.60~6.58—鋼軌表面

由表2可知,考慮型面磨耗時曲線外軌疲勞裂紋萌生壽命小于不考慮磨耗時的預測結果,與現場觀測的結果更為接近。這時裂紋萌生于鋼軌亞表面,距離軌頭表面深度約2~3 mm,考慮磨耗時預測的萌生位置略靠近軌頂中心。

3 結果分析

根據表2中的裂紋萌生位置找出對應在鋼軌上的位置,見圖8。其中,A點為考慮鋼軌型面磨耗情況下的裂紋萌生位置,P點為不考慮型面磨耗時的裂紋萌生位置。為了分析考慮磨耗時的P點疲勞損傷,假定在型面磨耗情況下該位置點為P′。

單次車輪引起的疲勞損傷分布見圖9,其中橫坐標表示因磨耗達到設定的磨耗量而替換型面的次數,也反映裂紋萌生前的不同磨耗階段。由圖9可知:

(1) 在標準型面時(型面第1個磨耗階段),P(P′)點的單次疲勞損傷大于A點的單次疲勞損傷,說明最初輪軌作用對P(P′)點的損傷較大,印證了當不考慮磨耗時,P點的單次疲勞損傷以及累積疲勞損傷都較大,裂紋最終萌生于P點而不是A點。

(2) 在第1~4個型面磨耗階段(第4個磨耗階段對應的車輪通過次數約為0.91×105~1.22×105,換算成通過總重約為1.67~2.24×107kN),A點的的單次疲勞損傷都較小,但其單次疲勞損傷呈緩慢增加趨勢,P′點也具有相似的規律。說明在這個時期,輪軌幾何關系匹配比較好,單次疲勞損傷較小且發展緩慢,裂紋萌生相對較慢。

(3) 第5個型面磨耗階段(車輪通過次數約為1.22×105~1.54×105,換算成通過總重約為2.24~2.83×107kN),A點的單次疲勞損傷比P′點大49%。說明這個磨耗階段,由于軌距角和軌肩的型面被磨損降低,輪軌接觸點略向軌頂中心移動。同時,輪軌接觸斑長短半軸之比略有增大,導致接觸應力增加,使得A點疲勞損傷加快增長。

(4) 隨著型面的進一步磨耗(第6、7個磨耗階段),輪軌作用關系越來越差,A點的單次疲勞損傷迅速增加,導致A點比P′點更早達到疲勞臨界損傷值,因此,考慮磨耗時A點對應的車輪通過次數較不考慮磨耗時P點的車輪通過次數小約60%,A點先達到疲勞累積損傷臨界值。

將每個磨耗階段的A點和P′點單次疲勞損傷乘以對應的車輪通過次數,可以得到2點的疲勞累積損傷,再結合不考慮磨耗時P點的疲勞累積損傷(以單次損傷為斜率隨車輪通過次數累積而呈線性增加趨勢),得到疲勞累積損傷與車輪通過次數的關系,見圖10。

由圖10可知:在第1個磨耗階段(車輪通過次數約為0.3×105以下),P′(P)點的疲勞累積損傷為0.055,A點的疲勞累積損傷為0.043,說明若沒有磨耗,P′(P)點將以這個累積損傷達到疲勞破壞,對應車輪通過次數為5.39×105次;但是由于輪軌接觸會因為型面的磨耗而發生變化,進一步導致軌頭各點應力應變的變化,因此,A點和P′點的疲勞累積損傷隨著磨耗逐漸增加,且A點增加得更快;到了第5個磨耗階段(車輪通過次數約為1.22×105~1.54×105次),A點的疲勞累積損傷為0.27,大于P′點的疲勞累積損傷0.22,并持續快速增大;最后,當A點疲勞累積損傷達到1時,P′點的疲勞累積損傷為0.62,A點為裂紋萌生點,而不是不考慮磨耗時的P′(P)點。裂紋萌生在A點時對應的車輪通過次數為2.17×105次。

綜上所述,在考慮磨耗的裂紋萌生壽命預測方法中,由于考慮了型面磨耗引起的輪軌接觸位置和鋼軌應力應變變化,找到了鋼軌內部更容易發生疲勞累積直至破壞的點,比不考慮磨耗的裂紋萌生壽命預測方法得到的裂紋萌生位置更合理、裂紋萌生壽命更接近現場觀測。

4 結論

(1) 將磨耗模型和裂紋萌生預測模型相結合,根據不同的磨耗階段分別進行鋼軌應力應變計算和疲勞損傷累積,建立了考慮磨耗的鋼軌疲勞裂紋萌生壽命預測方法,裂紋萌生預測更接近現場觀測結果,并可分析磨耗和疲勞裂紋萌生的相互關系。

(2) 在疲勞裂紋萌生之前,磨耗使得輪軌經常接觸區域的型面發生變化,造成接觸斑逐漸向軌頂中心方向偏移,使得處在接觸斑偏移方向的材料點的疲勞累積損傷加快,最終更早地萌生裂紋;而遠離接觸斑影響范圍的、在不考慮磨耗時認為是裂紋萌生點的材料點則不會萌生裂紋。

(3) 預測磨耗情況下的U75V熱處理鋼軌在800 m半徑曲線外軌時的裂紋萌生壽命約為4.5×107kN,對應車輪通過次數2.17×105次。當輪軌摩擦系數為0.3時,裂紋萌生在鋼軌亞表面,距離軌頂面2.42 mm。

(4) U75V熱處理鋼軌在上道初期的一定階段內輪軌關系較穩定,磨耗和疲勞裂紋發展緩慢。該階段大約是鋼軌上道至車輪通過次數約為1.22×105次,對應通過總重約為2.24×107kN。

參考文獻:

[1] ZHOU Yu,WANG Shaofeng,WANG Tianyi,et al.Field and Laboratory Investigation of the Relationship Between Rail Head Check and Wear in a Heavy-haul Railway[J].Wear,2014,315(1/2):68-77.

[2] 王少鋒,周宇,許玉德,等.基于蠕滑理論的鋼軌臨界平面疲勞參量仿真[J].鐵道學報,2014,36(4):65-70.

WANG Shaofeng,ZHOU Yu,XU Yude,et al.Simulation of Fatigue Parameters of Rail Critical Plane Based on Wheel-rail Creep Theory [J].Journal of the China Railway Society,2014,36(4):65-70.

[3] RINGSBERG J W.Life Prediction of Rolling Contact Fatigue Crack Initiation[J].International Journal of Fatigue,2001,23(7):575-586.

[4] JIANG Y,SEHITOGLU H.A Model for Rolling Contact Failure[J].Wear,1999,224(1):38-49.

[5] ARCHARD J F.Wear Control Handbook:Wear Theory and Mechanisms [M].New York:ASME,1980:161-78.

[6] LI Z L,ZHAO X,DOLLEVOET R,et al.Differential Wear and Plastic Deformation as Causes of Squat at Track Local Stiffness Change Combined with Other Track Short Defects[J].Vehicle System Dynamics,2008,46(1):237-246.

[7] KALOUSEK J,MAGEL E.Achieving a Balance:The Magic Wear Rate[J].Railway Track & Structures,1997,93(5):50-52.

[8] MAGEL E,RONEY M,KALOUSEK J,et al.The Blending of Theory and Practice in Modern Rail Grinding[J].Fatigue and Fracture of Engineering Materials and Structures,2003,26(10):921-929.

[9] FLETCHER D I,BEYNON J H.Equilibrium of Crack Growth and Wear Rates During Unlubricated Rolling-sliding Contact of Pearlitic Rail Steel[J].Proceedings of the Institution of Mechanical Engineers,Part F:Journal of Rail and Rapid Transit,2000,214(2):93-105.

[10] DONZELLA G,FACCOLI M,MAZZA,et al.Progressive Damage Assessment in the Near-surface Layer of Railway Wheel-rail Couple Under Cyclic Contact[J].Wear,2011,271(1/2):408-416.

[11] MADGE J J,LEEN S B,MCCOLL I R,et al.Contact-evolution Based Prediction of Fretting Fatigue Life:Effect of Slip Amplitude[J].Wear,2007,262(9/10):1 159-1 170.

[12] MADGE J J,LEEN S B,SHIPWAY P H.A Combined Wear and Crack Nucleation-propagation Methodology for Fretting Fatigue Prediction[J].International Journal of Fatigue,2008,30(9):1 509-1 528.

[13] LEEN S B,HYDE T H,RATSIMBA C,et al.An Investigation of the Fatigue and Fretting Performance of a Representative Aero-engine Spline Coupling[J].Journal of Strain Analysis for Engineering Design,2002,37(6):565-583.

[14] ZHANG T,HARRISON N M,MCDONNELL P F,et al.A Finite Element Methodology for Wear-fatigue Analysis for Modular Hip Implants[J].Tribology International,2013,65(3):113-127.

[15] LI Z L,KALKER J J.Simulation of Severe Wheel-Rail Wear[C]// Proceedings of the 6th International Conference on Computer Aided Design,Manufacture and Operation in the Railway and Other Advanced Mass Transit Systems.UK:WIT Press,1998:393-402.

[16] JENDEL T,BERG M.Contact Mechanics:Prediction of Wheel Wear for Rail Vehicles-Methodology and Verification[M].Berlin:Springer Netherlands,2002:229-236.

[17] 王天一.重載鐵路小半徑曲線鋼軌特殊型面設計方法研究[D].上海:同濟大學,2014:35-40.

[18] WANG J X,XU Y D,LIAN S L,et al.Probabilistic Prediction Model for Initiation of RCF Cracks in Heavy-haul Railway[J].International Journal of Fatigue,2011,33(2):212-216.

[19] 王少鋒.重載鐵路鋼軌疲勞裂紋萌生及路徑演變規律研究[D].上海:同濟大學,2014:131-136.

猜你喜歡
裂紋模型
一半模型
裂紋長度對焊接接頭裂紋擴展驅動力的影響
一種基于微帶天線的金屬表面裂紋的檢測
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
Epidermal growth factor receptor rs17337023 polymorphism in hypertensive gestational diabetic women: A pilot study
微裂紋區對主裂紋擴展的影響
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉換方法初步研究
預裂紋混凝土拉壓疲勞荷載下裂紋擴展速率
主站蜘蛛池模板: 国产成人高清在线精品| 国产在线专区| 免费xxxxx在线观看网站| 四虎永久免费网站| 成人在线观看一区| 国产成人永久免费视频| 色综合久久久久8天国| 国产日本一区二区三区| 亚洲欧美日韩中文字幕在线| 亚洲免费成人网| 亚洲国产看片基地久久1024| 免费国产黄线在线观看| 精品国产自在现线看久久| 国产噜噜噜视频在线观看 | 欧美亚洲欧美区| 日韩欧美中文字幕在线精品| 欧美亚洲国产一区| 亚洲无码高清一区| 玖玖精品在线| 玖玖精品视频在线观看| 免费午夜无码18禁无码影院| 人妻中文久热无码丝袜| 91久久青青草原精品国产| 国产又色又爽又黄| 在线日韩日本国产亚洲| 国产精品视频第一专区| 欧美在线观看不卡| 秋霞国产在线| 久久久久免费看成人影片| 伊人无码视屏| 999国产精品永久免费视频精品久久 | 国产成人夜色91| 亚洲成人精品在线| 久久综合婷婷| 欧美性爱精品一区二区三区| 久久夜夜视频| 婷婷六月综合网| 国产乱子伦手机在线| 91久久偷偷做嫩草影院| 久久99精品久久久久久不卡| 极品av一区二区| 欧美日韩在线第一页| 亚洲一区二区三区国产精华液| 免费一级成人毛片| 凹凸国产分类在线观看| 2021天堂在线亚洲精品专区| 成人国产精品网站在线看| 亚洲午夜国产片在线观看| 国产91色在线| A级毛片无码久久精品免费| 亚洲国产成人精品青青草原| 99久久精品久久久久久婷婷| 久久久久久久久久国产精品| 青青操视频在线| 无码av免费不卡在线观看| 国产免费高清无需播放器| 欧洲亚洲欧美国产日本高清| 久久黄色小视频| 国产精品99在线观看| 亚洲综合精品第一页| 日韩av资源在线| 久久人人妻人人爽人人卡片av| 免费在线a视频| 456亚洲人成高清在线| 亚洲精品无码av中文字幕| 欧美激情二区三区| 国产综合网站| 激情综合图区| 日韩毛片在线播放| 波多野结衣视频网站| 国产av无码日韩av无码网站| 亚洲成aⅴ人在线观看| 天天做天天爱夜夜爽毛片毛片| 午夜国产理论| 九色视频一区| 国模极品一区二区三区| 国产手机在线小视频免费观看 | 久久综合成人| 国产青青草视频| 成人在线天堂| 凹凸精品免费精品视频| 色婷婷综合在线|