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

直筒-錐段型鋼結構冷卻塔平均風荷載及靜風響應分析

2017-04-21 01:09:38柯世堂杜凌云劉東華馬兆榮
振動與沖擊 2017年7期
關鍵詞:模型

柯世堂, 杜凌云, 劉東華, 馬兆榮

(1. 南京航空航天大學 土木工程系,南京 210016;2. 廣東省電力設計研究院有限公司 土木工程技術部,廣州 210663)

直筒-錐段型鋼結構冷卻塔平均風荷載及靜風響應分析

柯世堂1, 杜凌云1, 劉東華2, 馬兆榮2

(1. 南京航空航天大學 土木工程系,南京 210016;2. 廣東省電力設計研究院有限公司 土木工程技術部,廣州 210663)

作為一種新穎的典型風敏感結構,直筒-錐段型鋼結構冷卻塔的動力特性和風致受力性能亟待研究。以國內擬建的某超大型鋼結構冷卻塔(189 m)為例,基于有限元方法分別建立主筒、主筒+加強桁架、主筒+加強桁架+附屬桁架(鉸接)、主筒+加強桁架+附屬桁架(固接)四種鋼結構冷卻塔模型,并對比分析其動力特性及傳力路徑;然后基于計算流體動力學(CFD)技術進行直筒-錐段型冷卻塔表面平均風荷載數值模擬,有分別加載規范和數值模擬風壓對四種模型進行風致響應分析,對比研究增設加強桁架、附屬桁架及與主筒和地面不同連接方式對直筒-錐段型鋼結構冷卻塔動力特性和靜風響應的影響。主要研究結論可為我國此類超大型鋼結構冷卻塔的結構選型和抗風設計提供依據。

直筒-錐段型鋼結構冷卻塔;計算流體動力學數值模擬;動力特性;靜風響應;參數分析

隨著國內經濟的提升和電機總量的不斷擴大,冷卻塔的建設日益朝著高大化、多樣化的趨勢發展[1-2],作為一種新穎的冷卻塔結構形式,鋼結構冷卻塔施工快、可回收且不受冬季溫度限制,在國外已有十余例成功應用于電廠間接空冷建設[3],而在國內尚處于萌芽階段,僅有新疆某直筒型鋼結構冷卻塔在建。傳統的冷卻塔[4]多以鋼筋混凝土為基本材料,塔筒以雙曲薄壁結構形式為主,通常簡化為板殼單元進行受力計算;鋼結構冷卻塔通過不同截面形式的鋼管組成復雜的桁架體系,其構件受力形式與梁單元更為接近。故與傳統鋼筋混凝土冷卻塔相比,鋼結構冷卻塔總體上阻尼更小、頻率更低、且振型及傳力路徑復雜,對風荷載更加敏感,其風致動力放大效應和受力性能問題也更為突出。

針對大型冷卻塔的抗風問題,國內外學者對大型混凝土雙曲冷卻塔的風荷載[5-7]、干擾效應[8-9]、風致響應[10-11]、整體和局部穩定性[12]等方面進行了大量的研究,相關成果很好地指導了混凝土雙曲冷卻塔的抗風設計,然而鮮有涉及直筒-錐段型鋼結構冷卻塔的風荷載和受力性能研究。

鑒于此,本文以國內擬建的某超大直筒-錐段型鋼結構冷卻塔(189 m)為例,分別建立主筒、主筒+加強桁架、主筒+加強桁架+附屬桁架(鉸接)、主筒+加強桁架+附屬桁架(固接)四種結構形式的鋼結構冷卻塔有限元模型,然后采用計算流體動力學(Computation Fluid Dynamics, CFD)。在此基礎上,對比研究增設加強桁架、附屬桁架及與地面不同連接方式對直筒-錐段型鋼結構冷卻塔動力特性和靜風響應的影響,主要結論可為此類超大型鋼結構冷卻塔的抗風設計提供參考依據。

1 動力特性分析

1.1 工程簡介及有限元建模

該大型鋼結構冷卻塔整體塔高189 m,由主筒、加強桁架和附屬桁架組成,主要構件信息如表1所示。為便于表達,將主筒、主筒+加強桁架、主筒+加強桁架+附屬桁架(鉸接)、主筒+加強桁架+附屬桁架(固接)四種結構形式對應的計算模型分別簡稱為模型一、模型二、模型三和模型四。采用大型通用軟件ANSYS進行四種模型的有限元建模,分別如圖1所示,所有桿件均采用BEAM188單元,主筒底部斜桿與地面固接,模型二中加強桁架外圍節點與主筒對應節點共用,模型三中附屬桁架上部和下部分別與主筒和地面鉸接,模型四中附屬桁架上部和下部分別與主筒和地面固接。

表1 大型鋼結構冷卻塔構件信息列表

圖1 鋼結構冷卻塔有限元模型示意圖Fig.1 The finite element model of cooling tower for steel structure

1.2 動力特性分析

考慮到表面蒙皮的重量會對后續模態和靜風響應分析造成一定影響,通過增大密度的方法考慮其附加質量,將鋼的密度乘以1.15的系數后作為鋼結構冷卻塔的等效密度[13]。

圖2和圖3分別給出了四組鋼結構冷卻塔模型基階和第50階典型振型圖。由圖可見,模型一振型復雜且變形較大,并呈現出明顯的三維特征,環向和豎向諧波數隨振型階數的增大而增多;增設加強桁架使得模型環向和豎向變形約束明顯;由于附屬桁架每榀之間的聯系橫桿剛度較弱,使得模型三和模型四的主要變形明顯從主筒轉移到附屬桁架,其表現為附屬桁架的環向扭轉,且固接相對鉸接變形略小。

圖2 四種不同模型基階振型示意圖Fig.2 The first step vibration chart of four models

圖3 四種不同模型第50階振型示意圖Fig.3 The 50th step vibration chart of four models

圖4給出了四種鋼結構冷卻塔模型前200階頻率隨階數的變化曲線。由圖4可知:①增設加強桁架的模型二自振頻率最大,在20階之前頻率增長迅速,之后逐漸呈現線性規律;②模型一的前10階自振頻率最小,自振頻率幾乎呈線性規律分布;③附屬桁架的增設減小了模型20階之后的自振頻率,模型三和模型四的前10階自振頻率極為相近,固接相對鉸接的第10階~第35階頻率較大,第35階之后的頻率均小于鉸接模型;④模型三在第10階~第35階存在一個頻率平緩段,隨階數增大基本保持不變,而模型四的頻率平緩段則相對較長,從第10階~第65階基本保持不變。

圖4 四種不同模型前200階頻率圖Fig.4 The first 200 steps of frequency of four models

2 平均風荷載數值模擬

2.1 參數設置

由于模型一和二為標準圓形截面,故直接采用建筑結構荷載規范[14]中的風壓推薦曲線,而模型三和模型四需采用CFD技術進行平均風荷載數值模擬。按照實際尺寸建模,以30%的透風率考慮百葉窗工作開啟狀態,為了兼顧計算效率和精度,將整個計算域劃分為外圍區域和局部加密區域。外圍區域形狀規整,采用高質量的結構化網格進行劃分,局部加密區域內含冷卻塔模型,采用非結構化網格進行劃分,總網格數量約640萬。圖5給出了計算域及模型網格劃分。數值模擬計算參數設置如表2所示。

表2 數值模擬計算參數

圖5 計算域及模型網格劃分示意圖Fig.5 The sketch of computational field mesh

2.2 結果分析

通過數值模擬計算得到直筒-錐段型鋼結構冷卻塔模型表面平均風荷載分布,將其分為下部錐段風壓和上部直筒段風壓。圖6和圖7分別給出了兩部分風壓沿環向和子午向分布曲線并與火工規范[15]加肋雙曲線混凝土冷卻塔和建筑荷載規范圓截面構筑物的數據進行對比??梢园l現:①CFD數值模擬得出的冷卻塔下部錐段平均與雙曲冷卻塔曲線較為吻合,側風面區別壓力系數值差別較大,錐段底部受氣流影響風壓略偏離于規范數據;②上部直筒段的平均風壓系數沿環向分布與圓截面構筑物的壓力系數曲線十分接近,但側風區和背風區數值有明顯區別,且端部繞流特性使得塔頂壓力系數偏離規范較遠;③對比結果同時驗證了本文數值方法用于鋼結構冷卻塔風荷載模擬的可靠性。

圖8給出了用于下文風荷載輸入的CFD數值模擬上部主筒及下部錐段數值模擬風壓系數平均值與兩種規范的對比曲線。由圖8可知:①兩種規范數據相比,圓截面構筑物背風區風壓平臺較小,負壓極值出現滯后且數值較大;②下部錐段的負壓極值點和分離點對應角度與火工規范中雙曲冷卻塔一致,但在背風區負壓數值要明顯小于規范值;③上部直筒段風壓系數分布曲線與荷載規范中圓截面構筑物基本吻合。表3具體列出了下部錐段和上部直筒段計算風壓分別與火工規范和建筑荷載規范的差值,冷卻塔塔頂、中部和底部流場平均速度云圖見圖9。

圖6 下部錐段平均風壓系數與規范對比曲線Fig.6 Wind pressure coefficient for conic section below

圖7 上部直筒段平均風壓系數與規范對比曲線Fig.7 Wind pressure coefficient for cylinder section upper

圖8 數值模擬平均風壓系數與規范曲線對比圖Fig.8 Average wind pressure coefficient for the model

環向角度/(°)下部錐段上部直筒段環向角度/(°)下部錐段上部直筒段00.00270.15171800.09240.0316120.06770.07281920.09240.0314240.13820.00612040.09980.0074360.16120.12162160.10710.0462480.18430.23702280.15440.1525600.21890.13922400.20170.2587720.25350.04142520.20200.0944840.02580.05562640.20230.0698960.20190.06982760.02560.05561080.20180.09442880.25350.04141200.20170.25873000.21910.13921320.15440.15253120.18480.23701440.10710.04623240.16170.12161560.09980.00723360.13860.00611680.09240.03183480.06790.0728

(a) z=37 m

(b) z=107 m

(c) z=187 m圖9 冷卻塔壓力系數及流場速度云圖Fig.9 The nephograms of pressure coefficient and speed for cooling tower and flow field

3 靜風響應分析

3.1 計算參數選取

該擬建冷卻塔位于B類地貌,基本風壓取為0.45 kPa。風荷載計算時風振系數取1.9,群塔放大系數取1.1。模型一和模型二采用建筑荷載規范中圓截面構筑物風壓系數,模型三和模型四采用數值模擬風壓系數進行計算?;陟o風響應計算結果,提取主筒、加強桁架和附屬桁架外部所有節點的徑向位移及典型構件應力響應,研究不同結構及附屬桁架與地面連接方式對鋼結構冷卻塔靜風響應的影響。

3.2 位移響應

3.2.1 主筒位移響應

圖10給出了四種模型的主筒徑向位移沿環向和子午向變化云圖,表4提煉出了四種模型主筒徑向位移最大和最小值。對比發現:①四種模型的主筒徑向位移均關于0°迎風角對稱,最大負位移在0°角附近,最大正位移出現在±70°~100°;②模型一的主筒變形最為顯著,極值29.49 cm和-33.59 cm位于塔頂處;模型二的主筒變形急劇減小,極值4.44 cm和-3.14 cm位于40~60 m高度之間;模型三和模型四的徑向位移分布模式相似,極值均位于第一層桁架所在處,且模型四的徑向位移整體上小于模型三,附屬桁架的設置使主筒大部分區域的徑向位移顯著減小。

表4 四種模型主筒節點徑向位移極值

Tab. 4 Extremum radial displacement of nodes for the main cylinder cm

徑向位移模型編號模型一模型二模型三模型四極大正值29.494.444.383.53極大負值-33.59-3.14-6.18-5.00

圖10 四種模型主筒徑向位移沿環向和子午向變化云圖Fig.10 The radial displacement along the ring and meridian nephogram for the main cylinder

3.2.2 加強桁架位移響應

圖11給出了模型二~模型四的三層加強桁架的徑向位移變化曲線,表5給出了加強桁架徑向位移極值大小及方位對比。分析可知:①模型二第一層加強桁架的迎風面與背風面徑向位移與模型三和四相比較為接近,但側風面徑向位移相對較大,且隨著高度的增大,三種模型加強桁架的區別越來越明顯;②第三層和第五層加強桁架節點徑向位移模型三最大,模型四次之,模型二最小,端部效應使得第五層加強桁架的背風面位移由負轉正,區別于其他幾層;③三種模型位移極值均出現在第二層加強桁架的迎風面和側風面上,模

型二加強桁架極值位移最小,模型四加強桁架受力優于模型三。

表5 模型二~模型四加強桁架節點徑向位移極值

Tab. 5 Extremum radial displacement of nodes for the stiffening truss cm

徑向位移模型編號模型一模型二模型三極大正值1.984.383.54極大負值-2.08-5.98-4.82出現位置二層290°/0°二層70°/0°二層290°/0°

圖11 模型二~模型四加強桁架徑向位移沿環向和子午向變化云圖Fig.11 The radial displacement along the ring and meridian nephogram for the stiffening truss

3.2.3 附屬桁架位移響應

圖12給出了模型三和四附屬桁架徑向位移云圖。對比分析可知,附屬桁架與主筒和地面的連接方式對大型直筒錐段型冷卻塔附屬桁架徑向位移分布趨勢的影響微弱,但鉸接方式下附屬桁架徑向位移值大于固接方式。

圖12 模型三和模型四附屬桁架徑向位移沿環向和子午向變化云圖Fig.12 The nephogram of radial displacement for the accessory truss

3.3 應力響應

四種冷卻塔模型的各桿件梁單元軸向應力和Y向彎曲應力極值大小及方位對比見表6。由表6可知:①模型一和模型二桿件單元軸向應力極值均位于主筒49.5 m高度處的迎風面和側風面橫桿上,而模型一彎曲應力極大和極小值分別位于冷卻塔頂部迎風面和側風面的橫桿上,設置加強桁架后,彎曲應力極值急劇減小并轉移至模型中下部;②設置附屬桁架后,應力極值從主筒橫桿轉移至附屬桁架頂部的外部斜桿上,附屬桁架與主筒和地面的連接方式對應力分布方位產生很大的影響,綜合比較可知鉸接方式下應力極值略大于固接。

表6 四種模型梁單元軸向應力和Y向彎曲應力對比

4 結 論

基于CFD和有限元技術,對超大直筒-錐段型鋼結構冷卻塔表面平均風荷載進行了數值模擬,并對比研究了增設加強桁架、附屬桁架及與主筒和地面不同連接方式對直筒-錐段型鋼結構冷卻塔動力特性和靜風響應的影響,得到主要結論如下:

(1) 主筒結構振型的三維特征顯著且變形最大,增設加強桁架后主筒徑向變形約束明顯,其自振頻率最大且隨階數近似線性分布,增設附屬桁架后主要變形明顯轉至附屬桁架,固接下結構整體變形相對鉸接下略小,隨階數增大均出現明顯的頻率平緩段。

(2) CFD數值模擬得到的冷卻塔下部錐段平均風壓分布曲線的負壓極值點和分離點對應角度與火工規范中雙曲冷卻塔一致,但在背風區負壓數值要明顯小于規范值;上部直筒段風壓系數分布曲線與荷載規范中圓截面構筑物基本吻合。

(3) 四種模型的主筒徑向位移均關于0°迎風角對稱,最大正、負位移均位于側風面和迎風面。模型一主筒變形最為顯著,極值出現在塔頂;模型二主筒位移明顯減小,極值下移至塔體中下部,且加強桁架極值位移最?。辉O置附屬桁架后主筒位移急劇減小,加強桁架位移相對模型二較大,最大徑向位移轉移至附屬桁架且分布趨勢相似,但模型三位移值略大于模型四。

(4) 四種模型的單元軸向應力和彎曲應力極值均位于迎風面和背風面,但模型一和二的軸向應力極值出現在主筒中下部橫桿上,模型一彎曲應力極值位于塔頂,設置加強桁架后下移至模型中下部,模型三和模型四應力極值轉移至附屬桁架頂部的外部斜桿上且前者略大于后者,與主筒和地面的連接方式對應力極值的分布方位產生較大的影響。

本文成果的知識產權歸廣東省電力設計研究院有限公司所有。

[ 1 ] 張軍鋒, 陳淮, 李靜斌. 基于反應譜方法的雙曲冷卻塔地震響應特征分析[J]. 振動與沖擊, 2015, 34(10):128-134. ZHANG Junfeng, CHEN Huai, LI Jingbin. Seimic response features of hyperboloidal cooling towers based on response spectral analysis[J]. Journal of Vibration and Shock, 2015, 34(10):128-134.

[ 2 ] KE S T , GE Y J, ZHAO L, et al. Wind-induced vibration characteristics and parametric analysis of large hyperbolic cooling towers with different feature sizes[J]. Structural Engineering and Mechanics, 2015, 54(5): 891-908.

[ 3 ] IZADI M, BARGI K. Natural draft steel hyperbolic cooling towers: optimization and performance evaluatim [J]. Structural Design of Tall ﹠ Special Buildings, 2013,23(9):713-720.

[ 4 ] 柯世堂, 侯憲安, 趙林,等. 特大型冷卻塔風荷載和風振響應參數分析:自激力效應[J]. 土木工程學報, 2012, 45(12): 45-53. KE Shitang, HOU Xian’an, ZHAO Lin, et al. Parameter analysis of wind loads and wind induced responses for super-large cooling towers: self-excited force effects[J]. China Civil Engineering Journal, 2012, 45(12): 45-53.

[ 5 ] 劉若斐, 沈國輝, 孫炳楠. 大型冷卻塔風荷載的數值模擬研究[J]. 工程力學, 2006, 23(增刊1): 177-183. LIU Ruofei, SHEN Guohui, SUN Bingnan. Numerical simulation study of wend load on large hyperbolic cooling towers[J]. Engineering Mechanics, 2006, 23(Sup1): 177-183.

[ 6 ] 鄒云峰, 牛華偉, 陳政清. 超大型冷卻塔單塔外表面風荷載三維效應及其設計取值[J]. 振動與沖擊, 2013, 32(24): 76-82. ZOU Yunfeng, NIU Huawei, CEHN Zhengqing. Three-dimensional effect and design values of outer surface wind loading for a single super-large cooling tower[J].Journal of Vibration and Shock, 2013, 32(24): 76-82.

[ 7 ] KE S T, GE Y J, ZHAO L, et al. A new methodology for analysis of equivalent static wind loads on super-large cooling towers[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2012, 111(3): 30-39.

[ 8 ] 鄒云峰, 牛華偉, 陳政清. 基于完全氣動彈性模型的冷卻塔干擾效應風洞試驗研究[J]. 湖南大學學報(自然科學版), 2013, 40(12):2-7. ZOU Yunfeng, NIU Huawei, CEHN Zhengqing. Wind tunnel test on wind-induced interference effect of cooling towers based on full areo-elastic model[J]. Journal of Hunan University (Natural Sciences), 2013, 40(12):2-7.

[ 9 ] 沈國輝, 余關鵬, 孫炳楠,等. 大型冷卻塔雙塔干擾的風洞試驗研究[J]. 振動與沖擊, 2011, 30(3):109-114. SHEN Guohui, YU Guanpeng, SUN Bingnan, et al. Study on interference effects of two large hyperbolic cooling towers by using wind tunnel test[J]. Journal of Vibration and Shock, 2011, 30(3):109-114.

[10] 沈國輝, 余關鵬, 孫炳楠,等. 大型冷卻塔風致響應的干擾效應[J]. 浙江大學學報(工學版), 2012, 46(1):33-45. SHEN Guohui, YU Guanpeng, SUN Bingnan, et al. Interference effect of wind-induced response on large hyperbolic cooling tower[J]. Journal of Zhejiang University(Engineering Science), 2012, 46(1):33-45.

[11] 鄒云峰, 李壽英, 牛華偉. 雙曲冷卻塔等效靜力風荷載規范適應性研究[J]. 振動與沖擊, 2013, 32(11): 100-105. ZOU Yunfeng, LI Shouying, NIU Huawei. Adaptability study on a China’s code for equivalent static wind load of hyperbolic cooling tower[J].Journal of Vibration and Shock, 2013, 32(11): 100-105.

[12] 張軍鋒, 葛耀君, 趙林. 基于風洞試驗的雙曲冷卻塔靜風整體穩定研究[J]. 工程力學, 2012, 29(5):187-196. ZHANG Junfeng, GE Yaojun, ZHAO Lin. Study on global aerostatic stability of hyperboloidal cooling towrs based on the wind tunnel tests[J]. Engineering Mechanics, 2012, 29(5):187-196.

[13] 柯世堂. 鋼結構冷卻塔數值風洞模擬試驗及風振系數研究報告[R]. 南京航空航天大學風工程試驗研究中心, 2015.

[14] 建筑結構荷載規范:GB 50009—2012 [S]. 北京: 中國建筑工業出版社, 2012.

[15] 火力發電廠水工設計規范: DL/T 5339—2006 [S]. 北京: 中國電力出版社, 2006.

Analysis of average wind pressure and responses under static wind for a cylinder-conic section steel cooling tower

KE Shitang1, DU Lingyun1, LIU Donghua2, MA Zhaorong2

(1. Department of Civil Engineering, Nanjing University of Aeronautics and astronautics, Nanjing 210016,China;2. Technology Department of Civil Engineering, Guangdong Electric Power Design Institute, Guangzhou 210663, China)

As a kind of new typical wind-sensitive structure, dynamic characteristics and wind-induced mechanical performance of a cylinder-conic section steel cooling tower are urgent to be studied. With a proposed super large steel cooling tower (189 m) in our country as an example, firstly, four steel cooling tower models induding a main cylinder, a main cylinder+ stiffening truss, a main cylinder+ stiffening truss+ accessory truss(hinged), and a main cylinder+ stiffening truss+ accessory truss(fixed) were established based on the finite element method. Their dynamic characteristics and force transmission paths were comparatively analyzed.Secondly, the numerical simulation of wind load on the surface of the cylinder-conic section steel cooling tower was conducted by using the computational fluid dynamics (CFD) method. Finally, with the specified and simulated wind pressure loads, respectively, wind-induced responses of four models were analyzed. The influences of adding stiffening truss and accessory truss and different ways of connection among accessory truss, the main cylinder and ground on dynamic characteristics and responses under static wind of the tower were comparatively studied. The studying conclusions provided a scientific basis for the structure form selection and wind-resistance design of this kind of super large steel cooling towers in our country.

cylinder-conic section steel cooling tower; CFD numerical simulation; dynamic characteristics; responses under static wind; parametric analysis

國家自然科學基金(51208254);中央高?;究蒲袠I務費專項資金資助(NS2016019);中國博士后科學基金(2013M530255; 1202006B);江蘇省優秀青年基金(BK20160083)

2015-10-28 修改稿收到日期: 2016-02-16

柯世堂 男,博士,副教授,1982年生 E-mail:verylyzh@163.com

TU279.7+41

A

10.13465/j.cnki.jvs.2017.07.023

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 自拍偷拍欧美| 免费一级毛片在线播放傲雪网| 亚洲综合在线最大成人| 又爽又大又黄a级毛片在线视频| 日韩国产欧美精品在线| 波多野结衣一区二区三区四区视频| 国产亚洲一区二区三区在线| 99re精彩视频| 精品99在线观看| 国产精品专区第1页| 国产精品毛片一区| 99久久精品国产麻豆婷婷| 日韩av高清无码一区二区三区| 国产亚洲精品自在久久不卡| 欧美三级不卡在线观看视频| 久久中文电影| 无码丝袜人妻| 尤物国产在线| 在线免费观看a视频| a国产精品| 午夜毛片免费看| 亚洲无码高清免费视频亚洲| 五月综合色婷婷| 精品免费在线视频| 久操中文在线| 伊人精品成人久久综合| 国产无吗一区二区三区在线欢| 日本亚洲最大的色成网站www| 狠狠ⅴ日韩v欧美v天堂| 久久五月天国产自| 国产性生交xxxxx免费| 欧美日韩免费| 国产精女同一区二区三区久| 亚洲日本韩在线观看| 无码国内精品人妻少妇蜜桃视频| 亚洲国产成人无码AV在线影院L| 亚洲AⅤ综合在线欧美一区| 国产成人无码Av在线播放无广告| 欧美色图第一页| 伊人久久久久久久| 午夜在线不卡| 国产网站免费| 亚洲人成成无码网WWW| 久草中文网| 欧美中文一区| 91欧美亚洲国产五月天| 亚洲精品高清视频| 久久人人妻人人爽人人卡片av| 亚洲福利网址| 亚洲色图欧美视频| 亚洲永久视频| 国产香蕉一区二区在线网站| 国产精品一区在线麻豆| aⅴ免费在线观看| 亚洲成在线观看 | 国产黑丝视频在线观看| 国产女人综合久久精品视| 亚洲美女一区| 欧美天堂在线| 亚洲无码高清视频在线观看 | 91久久大香线蕉| 国产99在线观看| 欧美午夜视频| 91久久天天躁狠狠躁夜夜| 亚洲中文字幕国产av| 久久夜色精品国产嚕嚕亚洲av| 成人免费黄色小视频| 国产精品成人啪精品视频| a级毛片在线免费观看| 国产永久在线视频| 亚洲区欧美区| 日韩视频免费| 99热这里只有精品免费国产| 欧美一级在线| 88国产经典欧美一区二区三区| 亚洲福利视频网址| 色婷婷在线影院| 国产麻豆另类AV| 国产成人亚洲毛片| 青青热久麻豆精品视频在线观看| 十八禁美女裸体网站| 国产在线一二三区|