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

某供水工程TBM刀盤破巖過程動靜態響應特性

2016-06-27 10:00:12夏毅敏陳卓林賚貺田彥朝唐露吳才章楊妹
哈爾濱工程大學學報 2016年5期

夏毅敏,陳卓,林賚貺,田彥朝,唐露,吳才章,楊妹

(1.中南大學 高性能復雜制造國家重點實驗室,湖南 長沙 410083; 2.中南大學 機電工程學院,湖南 長沙 410083)

某供水工程TBM刀盤破巖過程動靜態響應特性

夏毅敏1,2,陳卓2,林賚貺2,田彥朝2,唐露2,吳才章2,楊妹2

(1.中南大學 高性能復雜制造國家重點實驗室,湖南 長沙 410083; 2.中南大學 機電工程學院,湖南 長沙 410083)

摘要:刀盤是硬巖掘進機(TBM)的關鍵部件,為了研究刀盤刀具與掘進面耦合作用下載荷分布與變化規律,建立了某供水工程TBM刀盤-巖石相互作用的三維分析模型,基于顯式動力學方法模擬TBM刀盤掘進過程,研究掌子面巖體的損傷失效狀態以及刀盤刀具的動態掘進載荷。研究得到了全盤刀具三向力載荷及刀盤總載荷(總推力、總扭矩),其中刀盤扭矩平均值1 560 kN·m,最大值5 880 kN·m,推力平均值3 444 kN,最大值7 296 kN;通過工程實例驗證,得到實測載荷值稍大于仿真值,總推力誤差為15.74%,總扭矩誤差為15.38%,驗證了仿真結果的有效性。在此基礎上,對刀盤的動靜態響應特性進行研究,得到刀盤在極限工況下最大應力100 MPa,最大變形0.695 mm,確定刀盤有足夠的剛度強度。

關鍵詞:TBM刀盤;數值模擬;動態掘進特性;工程驗證;靜態力學特性;靜態響應;動態響應

網絡出版地址:http://www.cnki.net/kcms/detail/23.1390.U.20160411.0845.014.html

隨著隧道以及地下工程快速發展,盾構以及TBM等大型掘進裝備得到了廣泛應用[1-3]。刀盤刀具是TBM掘進機的關鍵部件,與盾構相比,TBM刀盤結構更為復雜,工作環境更為惡劣,刀具損耗大、換刀頻率高,承受大扭矩、大推力以及沖擊載荷的作用,振動十分劇烈[4-6]。對滾刀與掘進界面的相互作用進行研究,確定刀盤刀具的載荷特性是TBM設計的基礎。

目前在盾構與TBM刀盤刀具載荷方面,國內外學者做了大量研究。Jung-Woo Cho等[7]利用AUTODYN-3D,采用Drucker-Prager材料模型模擬巖石,對盤形滾刀直線切割巖石進行了模擬。Alessio Nardi等[8]應用離散元軟件來分析巖石被侵入發生破壞時的力學響應特征。李剛等[9]人通過建立巖石的損傷本構模型進行仿真分析,得到了不同巖石的最優切割參數。蘇翠俠等[10]利用ABAQUS模擬盾構掘進過程,分析了刀盤掘進載荷影響因素,得到了刀盤系統界面載荷沿半徑方向非線性分布的規律。霍軍周等[11]通過改進TBM刀盤支撐筋結構,提升了刀盤的靜動態性能。

本文結合工程實際,針對某引水工程刀盤的工作特性,利用數值模擬方法,建立了刀盤-巖石相互作用的三維分析模型,以Drucker-Prager準則作為巖石的屈服準則,應用包含單元刪除功能的損傷失效準則模擬切屑的形成與分離,基于顯式積分算法實現TBM刀盤掘進全物理過程的直接數值模擬,研究刀盤刀具與掘進面耦合作用下載荷分布與變化規律,分析該刀盤破巖過程的動靜態力學響應特性,并利用工程進行了部分驗證。

1刀盤動態掘進數值模型的建立

1.1巖石材料模型及損傷失效準則

巖石材料特性將對盤形滾刀載荷產生巨大影響。采用擴展的Drucker-Prager模型,引入參數K,其控制方程如下所示:

(1)

(2)

式中:t為相關應力參數;r為偏應力第三應力不變量;K為單軸拉伸屈服應力與單軸壓縮屈服應力比值;d為凝聚力;β為材料摩擦角;q為Mises等效應力;p=-1/3trace(σ)為平均壓應力。

在圖1中,當材料的壓縮強度與拉伸強度相同,即K=1,t=q,擴展Drucker-Prager準則退化為傳統的Drucker-Prager準則。

圖1 偏平面內線性擴展Drucker-Prager屈服準則Fig.1 Typical yield surfaces of linear model in deviatoric plane

設ωs為描述塑性變形隨等效塑性應變遞增的狀態變量,當ωs為1時,達到初始破壞點B,而后材料開始生成裂紋,剛度逐漸衰減,直至材料失去承載能力:

(3)

圖2 巖石損傷應力-應變響應曲線Fig.2 Description of isotropic damage mechanism

1.2刀盤破巖過程的顯式積分算法

TBM刀盤掘進是一個有著連續的動態接觸關系、包含材料的破壞和失效的復雜動態過程。

利用顯式動力學對其進行模擬,在時間段開始時(t時刻),求解整個刀盤系統的動力平衡方程:

(4)

在當前時間段開始時(t時刻)刀盤加速度為:

(5)

刀盤加速度是由中心差分法的時間積分得到的,即假定加速度為常數以求得速度的變化,用這個速度的變化值加上前一個時間段中點的速度來確定當前時間段的中點速度:

(6)

刀盤速度沿時間積分的結果加上此時間段開始時的位移,即為時間段結束時的位移:

(7)

由于方程的求解是非耦合的,不必同時求解聯立方程,因此大大節省了刀盤破巖過程的求解時間和存儲空間。

1.3刀盤掘進過程數值計算模型

1.3.1刀盤切削巖石三維模型

刀盤結構較為復雜,在動力學計算過程中,為了提高計算速度,將刀盤模型適當簡化,并將其與建立的巖石模型裝配,如圖3所示。刀盤直徑7.93 m,其中心滾刀8把(1#~8#),正面滾刀32把(9#~40#),邊緣滾刀11把(41#~51#)。

圖3 刀盤-巖石三維實體模型Fig.3 Model of cutterhead-rock

1.3.2邊界條件及網格劃分

在實際施工過程中,TBM刀盤旋轉切削巖石同時受到后方液壓缸的推進力,以緩慢的速度向前推進,為模型施加的邊界條件如下:

1)初始狀態時盤形滾刀與巖石剛好接觸。

2)刀盤旋轉速度為8 r/min。

3)刀盤的推進速度為6 mm/r。

4)約束巖石模型外邊界的位移自由度,保持待開挖巖石表面為自由面。

5)利用非光滑接觸模擬滾刀巖石的相互作用。

6)設定仿真時間為18 s,刀盤勻速旋轉,切削巖石三圈。

TBM刀盤材料為Q345鋼,隧道地質條件為混合花崗巖地層,刀盤及巖石的材料參數如表1所示。模型采用8節點六面體單元對實體單元進行網格劃分。

表1 刀盤以及巖石材料參數表

2刀盤動態掘進數值模型的建立

2.1巖石材料模型及損傷失效準則

由刀盤結構可以看出,滾刀始終高于刀盤面板,

在掘進過程中滾刀首先與掌子面接觸并成為整個切削過程的主體。初始掘進階段,在滾刀推力作用下,掌子面巖體發生了彈性變形;隨著切深的增加,應力值相應增大,由于巖石為脆性材料,在極短的時間內,巖石發生塑性變形;當達到強度極限后(對應于圖2中的B點),滾刀正前方巖體開始產生局部損傷(如圖4所示)。

損傷發生后,巖體剛度開始衰減,隨著刀盤的旋轉掘進,滾刀與巖體相互作用增強,損傷程度加劇,當損傷變量D的值達到1時,(對應于圖2中的C點),相應巖體因完全失去承載能力而被剝離刪除,在掌子面上形成了與滾刀切削軌跡相匹配的一系列同心圓溝槽(如圖5所示)。

圖4 t=1.5 s時巖體損傷云圖Fig.4 Damage contour of rock when t=1.5 s

圖5 t=9 s時巖體損傷云圖Fig.5 Damage contour of rock when t=9 s

2.2刀具刀盤動態掘進載荷分布與變化規律

對TBM刀盤破碎巖石過程進行模擬,可得到每個刀號滾刀和刀盤的載荷變化情況。

圖6 25#滾刀三向力載荷隨時間的變化規律Fig.6 three-dimensional mechanical force of 25# disc cutter

圖6為25#滾刀三向力載荷隨時間的變化曲線,由圖6可知:在t從0~2 s內,滾刀與巖石接觸,載荷逐漸增大,隨著滾刀繼續前進,巖石發生塑性變形并被壓實形成密實核,最終巖石隨著裂紋的增加而發生破碎,載荷大幅下降,完成了一個滾刀的破巖周期,滾刀切削巖石的過程就是上述破巖過程的不斷重復發生,呈現出階躍破碎特性;25#滾刀垂直力、滾動力、側向力的平均值分別為91.7、14.8、1.78 kN,峰值分別為130.7、79.8、20.8 kN,垂直力和滾動力的峰值可達平均值的1.5~5倍,因此滾刀軸承的選擇應同時考慮滾刀載荷平均值和峰值;滾刀刀刃兩側都受到巖石反作用力的作用,因兩側載荷方向相互抵消,側向力的合力較小。

圖7 不同刀位滾刀三向力分布曲線Fig.7 Distribution of three-dimensional mechanical load  ondifferent cutter

(a) 刀盤扭矩曲線

(b) 刀盤推力曲線圖8 刀盤總載荷變化規律Fig.8 Variation of cutterhead total load

對于過渡滾刀和邊緣滾刀,其安裝半徑差別不大,主要是安裝傾角的不同引起三向力的差異。隨著安裝傾角的增大,邊緣滾刀貫入度減小而導致垂直力減小,最外側邊緣滾刀的垂直力只有正滾刀垂直力的1/3左右;滾動力波動相對較小;由于存在安裝傾角,故過渡滾刀及邊緣滾刀側向力比中心滾刀與正滾刀的側向力要大得多,側向力指向刀盤中心,隨著安裝傾角的增大而明顯增大。

TBM掘進過程中刀盤不同刀位滾刀三向力平均值如圖7所示。由圖7可知,刀盤掘進載荷分布具有如下規律:對于中心滾刀和正滾刀,主要影響因素為刀具安裝半徑。隨著刀具安裝半徑的增大,垂直力增大,刀盤外區正滾刀(30#~40#)垂直力為中心滾刀垂直力的1.2倍;滾動力變化規律與垂直力趨勢一致,但波動較為平緩;側向力隨安裝半徑增大而減小,由于安裝半徑越大,滾刀越近似于直線切割,兩側巖石對滾刀的作用也越來越均衡。

圖8為刀盤總載荷(總推力、總扭矩)隨時間變化規律。隨著刀具滾壓破巖,刀盤載荷開始逐漸波動。在掘進的初始階段,滾刀與巖石開始接觸,接觸面匹配不佳,刀盤掘進載荷呈現緩慢波動;伴隨巖體剝離,刀具與前方巖體形成較好的接觸關系,且操作參數維持穩定,刀盤進入穩定掘進階段。整個掘進過程中,刀盤扭矩平均值為1 560 kN·m,最大值為5 880 kN·m;推力平均值為3 444 kN,最大值為7 296 kN。

2.3掘進載荷工程驗證

某引水工程地層巖性主要為侏羅系小南溝組凝灰巖和新太古代混合花崗巖。弱風化和微風化小南溝組凝灰巖、微風化混合花崗巖的物理力學性質參見表中,弱風化混合花崗巖的飽和單軸抗壓強度建議值50~70 MPa。

表2巖石物理力學性質參數表

Table 2Physical and mechanical properties parameter of rock

風化程度及巖性密度/(g/cm3)抗壓強度/MPa吸水率彈性模量/GPa泊松比弱風化粉砂巖2.7156.570.7460.200.29微風化粉砂巖2.6963.180.4843.400.13微風化混合花崗巖2.6559.860.4062.500.16

圖9 工程實際掌子面形貌Fig.9 Actual morphology of project's heading face

掌子面形貌如圖9所示,對比圖5可知,巖體經刀盤切削后形貌與動力學數值模擬所得巖體損傷形貌有一定的相似性,宏觀驗證了模擬巖體失效的有效性。取該掘進機2014年某天施工55 min,累計掘進1環掘進過程中刀盤載荷數據,刀盤實測扭矩1 800 kN·m,實測推力3 986 kN。將模擬得到的結果與實測結果進行對比,可知:實測載荷值稍大于仿真載荷值,總推力誤差為15.74%,總扭矩的誤差為15.38%,驗證了動力學仿真結果的準確性。

3刀盤動靜態力學響應特性

根據動態掘進特性分析,設置安全系數為2,將動態模擬中所得全盤刀具載荷乘以安全系數,分別等效施加于對應滾刀刀箱上,以此耦合為刀盤極限工況,在此工況下對刀盤力學響應特性進行研究。

3.1刀盤靜態力學響應特性分析

3.1.1刀盤強度分析

如圖10,刀盤面板應力分布基本在40 MPa范圍以內,由刀盤中心區域以傘狀向周圍發散,在刀具密集區域呈現相對較大的應力,在刀具稀少區域呈現較小應力分布,應力較大或應力集中主要出現在圖中。其中,中心刀區由于安裝了較多的滾刀,引起受力較大,主中心刀區應力最大為30 MPa;正滾刀密集區匯集了4把19寸滾刀,導致刀箱間面板面積余量較小且受載嚴重,產生40 MPa的應力;人孔區與正滾刀相距近,導致兩者結合處產生較大的應力分布,并向周圍擴展,又由于人孔的圓形設計,導致人孔周圍出現一定量的應力集中,最大為30 MPa。

圖10 刀盤面板應力分布云圖Fig.10 Stress distribution contour of cutterhead’s panel

圖11 刀盤支撐筋板處應力分布云圖Fig.11 Stress distribution contour of cutterhead’s support plates

如圖11所示,刀盤支撐筋板處的應力基本分布在40~75 MPa區域,為整個刀盤應力最大的區域,這是因為該位置除承擔全盤刀具傳遞的垂直載荷,還承擔主要的受扭作用,尤其在支撐筋板與刀盤背板的連接處,應力最大,同時,全刀盤最大應力100 MPa也發生在此區域。該最大應力低于刀盤材料Q345鋼的屈服極限,安全系數達到3.45。

3.1.2刀盤剛度分析

刀盤變形如圖12所示,由圖中可以看出,刀具變形主要集中在圖中標注的3個區域:中心刀區,正滾刀密集區,邊刀區。其中,中心刀區變形最大,這是因為刀盤變形主要由該區域承受的垂直載荷所致,而中心區域布刀數量最多,承擔大量的垂直載荷,中心區最大變形量達到0.695 mm;邊刀區的刀盤變形次之,因為該區域刀具分布密度較高,變形量較大,最大達到0.406 mm;此外,圖中標注2位置處正刀區刀具分布密度相較于其他正刀區較大,因此該區域產生的變形量也較大,達到0.348 mm。其他區域的變形量相對于以上三塊區域較小,按照布刀密度的減小依次遞減。

刀盤最大變形量為0.695 mm,與刀盤直徑相比,僅為刀盤直徑的0.008 8%,驗證了刀盤滿足剛度要求。

圖12 刀盤變形分布云圖Fig.12 Deformation distribution contour of cutterhead

3.2刀盤模態特性分析

刀盤的模態動態分析中,提取刀盤30階模態,通過頻率提取分析,得到主運動方向(Z方向)上總有效質量占可運動質量的比例是92%,因此提取15階振型是必需且足夠的。

通過模態分析,獲得每階模態的刀盤固有頻率以及振型,統計刀盤15階模態的固有頻率如表3。

表3 各階模態固有頻率表

3.3刀盤共振可能性分析

將圖9、10的刀盤動態掘進載荷進行傅里葉變換,采樣頻率1 000 Hz,幅頻曲線如圖13所示。

在圖13中A、B、C、D所示為刀盤高幅值載荷對應點,均發生在12 Hz頻率范圍以內,其中A點為刀盤載荷最大幅值點,對應頻率僅為0.167 Hz,可見刀盤發生大幅值載荷頻率很低,利于刀盤壽命;同時,由圖可知,刀盤載荷在12 Hz頻率以后趨于逐漸減小,根據表3得知刀盤1階最小固有頻率為58.4 Hz,載荷頻率到達該頻率時振幅可忽略不計,刀盤發生共振的可能性極小。

(a) 推力幅頻曲線

(b) 扭矩幅頻曲線圖13 刀盤載荷幅頻曲線Fig.13 Amplitude-frequency curve of cutterhead load

4結論

本文在有限元環境下建立了TBM刀盤系統綜合分析的數值仿真模型:

1)建立了TBM刀盤三維破巖仿真模型,對刀盤破巖過程的載荷分布及其響應特性進行了研究;

2)分析獲得動態掘進過程中整盤刀具三向力以及刀盤總載荷,刀盤轉動一圈過程中,扭矩平均值1 560 kN·m,最大值5 880 kN·m,推力平均值3 444 kN,最大值7 296 kN。通過實際工程驗證,得到實測載荷值稍大于仿真載荷值,總推力誤差為15.74%,總扭矩誤差為15.38%,驗證了仿真結果的準確性;

3)通過刀盤靜態力學特性分析,得到極限工況下刀盤所受最大應力為100 MPa,低于Q345的屈服強度,安全系數達到3.45,刀盤最大變形量0.69 mm,僅為刀盤直徑的0.008 8%,刀盤具有足夠的強度和剛度;

4)通過刀盤瞬時模態動態分析,獲得了刀盤各階振型、固有頻率,通過與動態掘進載荷的頻域分析對比驗證刀盤發生共振的可能性極小。

參考文獻:

[1]宋克志, 王夢恕. 復雜巖石地層隧道掘進機操作特性分析[J]. 土木工程學報, 2012, 45(5): 176-181.

SONG Kezhi, WANG Mengshu. Analysis of the operation characteristics of tunnel boring machine in complex rocks[J]. China civil engineering journal, 2012, 45(5): 176-181.

[2]張照煌, 喬永立. 全斷面巖石掘進機盤形滾刀布置規律研究[J]. 工程力學, 2011, 28(5): 172-177.

ZHANG Zhaohuang, QIAO Yongli. Research on the layout of TBM disc cutter[J]. Engineering mechanics, 2011, 28(5): 172-177.

[3]王虹. 我國綜合機械化掘進技術發展40a[J]. 煤炭學報, 2010, 35(11): 1815-1820. WANG Hong. The 40 years developmental review of the fully mechanized mine roadway heading technology in China[J]. Journal of China coal society, 2010, 35(11): 1815-1820.

[4]向前進, 申智杰, 劉曉波. 大直徑泥水盾構刀盤刀具應用技術[J]. 現代隧道技術, 2011, 48(4): 157-161.

XIANG Qianjin, SHEN Zhijie, LIU Xiaobo. Applied techniques for cutterhead and cutters of a large-diameter slurry shield[J]. Modern tunnelling technology, 2011, 48(4): 157-161.

[5]王旭, 趙羽, 張寶剛, 等. TBM滾刀刀圈磨損機理研究[J]. 現代隧道技術, 2010, 47(5): 15-19.

WANG Xu, ZHAO Yu, ZHANG Baogang, et al. Research on the ring wear mechanism of TBM disc cutter[J]. Modern tunnelling technology, 2010, 47(5): 15-19.

[6]ALBER M. Stress dependency of the Cerchar abrasivity index (CAI) and its effects on wear of selected rock cutting tools[J]. Tunnelling and underground space technology, 2008, 23(4): 351-359.

[7]CHO J W, JEON S, YU S W, et al. Optimum spacing of TBM disc cutters: a numerical simulation using the three-dimensional dynamic fracturing method[J]. Tunnelling and underground space technology, 2010, 25(3): 230-244.

[8]NARDIN A, SCHREFLER B. Numerical simulation of rock behaviour through a discrete model[J]. International journal of solids and structures, 2004, 41(21): 5945-5965.

[9]李剛, 王宛山. 全斷面巖石掘進機滾刀間距優化研究[J]. 哈爾濱工程大學學報, 2013, 34(7): 913-917.

LI Gang, WANG Wanshan. Study on optimum spacing of TBM disc cutters[J]. Journal of Harbin engineering university, 2013, 34(7): 913-917.

[10]蘇翠俠, 王燕群, 蔡宗熙, 等. 盾構刀盤掘進載荷的數值模擬[J]. 天津大學學報, 2011, 44(6): 522-528.

SU Cuixia, WANG Yanqun, CAI Zongxi, et al. Numerical simulation of excavation load on cutterhead in shield tunneling machine[J]. Journal of Tianjin university, 2011, 44(6): 522-528.

[11]霍軍周, 楊靜, 孫偉, 等. TBM刀盤支撐筋結構設計及靜動態特性分析[J]. 哈爾濱工程大學學報, 2014, 35(7): 883-888.

HUO Junzhou, YANG Jing, SUN Wei, et al. Structure design and static/dynamic analysis of TBM cutterhead supporting ribs[J]. Journal of Harbin engineering university, 2014, 35(7): 883-888.

本文引用格式:

夏毅敏,陳卓,林賚貺,等. 某供水工程TBM刀盤破巖過程動靜態響應特性[J]. 哈爾濱工程大學學報, 2016, 37(5): 732-737.

XIA Yimin, CHEN Zhuo, LIN Laikuang, et al. Static and dynamic response characteristics of a TBM cutterhead's rock-breaking process: a case study of a diversion project[J]. Journal of Harbin Engineering University, 2016, 37(5): 732-737.

Static and dynamic response characteristics of a TBM cutterhead's rock-breaking process: a case study of a diversion project

XIA Yimin1,2, CHEN Zhuo2, LIN Laikuang2, TIAN Yanzhao2, TANG Lu2,WU Caizhang2,YANG Mei2

(1. State Key Laboratory of High Performance Complex Manufacturing, Central South University, Changsha 410083, China; 2. Collegeof Mechanical and Electrical Engineering, Central South University, Changsha 410083, China)

Abstract:The cutterhead is the most important part of a tunnel-boring machine (TBM). To determine the load distribution and load variations under the action of coupling of the cutterhead and the cutters, a three-dimensional model of interactions between a water diversion project’s TBM cutterhead and the rock was established. Based on an explicit dynamics algorithm, the entire excavation process of the TBM cutterhead was simulated and the state of the tunnel face and the dynamic tunneling load on the cutterhead were obtained. By analysis, the distribution and variation of the three-dimensional mechanical loads on each disc cutter and the total load on the cutterhead were obtained. The average torque of the cutterhead was 1 560 kN·m, and the maximum was 5 880 kN·m. The average thrust was 3 444 kN, and the maximum was 7 296 kN. Verification provided by an actual project showed that the real measured load was slightly larger than the simulation values. The total thrust error was 15.74%, and the total torque error was 15.38%; therefore, the effectiveness of the simulation results was verified. Based on this, research into the static and dynamic response characteristics of the cutterhead showed that the cutterhead maximum bearing stress was 100 MPa, and the maximum deformation was 0.695 mm; therefore, the cutterhead was judged to be of sufficient strength and stiffness.

Keywords:TBM cutterhead; numerical simulation; dynamic tunneling properties; project verification; static mechanical properties;static response; dynamic response

收稿日期:2015-01-30.

基金項目:國家自然科學基金資助項目(51475478);國家863發展計劃基金資助項目(2012AA041801);湖南省戰略新興產業重大科技攻關項目(2012GK4068).

作者簡介:夏毅敏(1967-),男,教授,博士生導師. 通信作者:夏毅敏,E-mail:xiaymj@mail.csu.edu.cn.

DOI:10.11990/jheu.201501048

中圖分類號:U455.31

文獻標志碼:A

文章編號:1006-7043(2016)05-0732-06

網絡出版時間:2016-04-11.

主站蜘蛛池模板: 亚洲成人精品久久| 无码免费的亚洲视频| 日韩无码一二三区| 米奇精品一区二区三区| 制服丝袜无码每日更新| 九色综合伊人久久富二代| 人妻中文久热无码丝袜| 永久成人无码激情视频免费| 亚洲天堂久久| 免费观看国产小粉嫩喷水| 在线观看国产精品日本不卡网| 四虎永久免费地址在线网站| 乱码国产乱码精品精在线播放| 亚洲综合色在线| 99在线免费播放| 欧美色综合久久| 亚洲另类色| 综合亚洲网| 精品国产成人三级在线观看| 青青青伊人色综合久久| 干中文字幕| 波多野结衣在线se| 免费欧美一级| 波多野结衣在线se| 中文字幕在线免费看| 国产成人无码AV在线播放动漫 | 国产电话自拍伊人| 99热精品久久| 亚洲电影天堂在线国语对白| 亚洲中文在线看视频一区| 国产成人1024精品下载| 草草影院国产第一页| 色悠久久久久久久综合网伊人| 精品色综合| 国产免费网址| AV熟女乱| 中文字幕永久视频| 日韩精品成人在线| 国产区网址| 伊大人香蕉久久网欧美| 亚洲 欧美 中文 AⅤ在线视频| 欧美人与动牲交a欧美精品| 久久综合一个色综合网| 国产真实乱子伦视频播放| 成人久久精品一区二区三区 | 成人年鲁鲁在线观看视频| 欧美一级大片在线观看| 久久久四虎成人永久免费网站| 五月天久久综合| 免费国产黄线在线观看| 国产欧美日韩视频一区二区三区| 91毛片网| 国产va在线观看免费| 最新国产精品第1页| 成人午夜视频免费看欧美| 青草精品视频| 欧洲日本亚洲中文字幕| 國產尤物AV尤物在線觀看| 美女高潮全身流白浆福利区| 黄色片中文字幕| 久久网欧美| 天天综合网色中文字幕| 午夜福利网址| 91国内在线观看| 青草娱乐极品免费视频| 亚洲九九视频| 九色视频线上播放| 最近最新中文字幕在线第一页| 国产精品观看视频免费完整版| 久久亚洲国产最新网站| 无码人妻免费| 国产美女无遮挡免费视频| 日韩在线网址| 在线观看国产小视频| 亚洲v日韩v欧美在线观看| 无码国产伊人| 伊人五月丁香综合AⅤ| 久久综合亚洲鲁鲁九月天| 毛片手机在线看| 天天做天天爱天天爽综合区| 亚洲午夜福利精品无码不卡| 丁香六月综合网|