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

遮擋山體對輸電塔線體系風振疲勞的影響

2013-09-10 11:01:52張春濤李正良范文亮汪之松
振動與沖擊 2013年10期
關鍵詞:風速

張春濤,李正良,2,范文亮,2,汪之松,2,孫 毅

山地地形中,由于遮擋山體的存在致使山地風場的平均風速、湍流度、空間相關性、脈動風等更加復雜[1-2]。結構風荷載計算和風致響應分析時均不能采用一般平地風場的計算方法。因此,對周圍地形中存在遮擋山體干擾的高聳結構進行抗風設計和研究時,必須考慮具體山地情況。同時,特高壓大跨越輸電塔大多建設在地勢較高的山坡或江河兩邊較高的山頂上,受風場影響巨大,設計中對抗風亦提出了更高的要求[3-4]。

風荷載作用下,一方面輸電塔與輸電線的風致振動相互影響、共同作用,形成復雜的動力耦合體系;另一方面特高壓輸電塔線體系同時兼有高聳結構和大跨度結構的共同特點。因此,輸電塔線體系對風荷載作用十分敏感,容易發生振動疲勞破壞和動態坍塌。于是,近年來輸電塔線體系的風振疲勞得到了獲得了廣泛關注。Havard等[5]對多個服役15年左右的輸電塔進行觀察發現風荷載致使橫擔上出現了疲勞裂紋,并進一步分析得到風荷載和導線運動是疲勞裂紋產生的主要原因;Repetto等[6-9]對風荷載作用下高細結構的動力響應、風振疲勞等進行了一系列的研究,得到了計算疲勞壽命的建議方法,但動力分析時未對風向風速分布進行全面考慮;王之宏[10]、屠海明等[11]、鄧洪洲等[12]分別通過等效應力法、雨流法和等效窄帶法對格構式桅桿結構進行了風振疲勞研究;謝強等[13]通過氣彈性風洞試驗研究了1000kV特高壓輸電塔線體系在均勻流場和紊流場的風致響應,并首次采用光纖光柵測得了塔體及導線等的動應變;李正良等[14]、汪之松等[15-16]、肖正直[17]基于風洞試驗研究了特高壓輸電塔線體系的時域和頻域風振響應,并利用Miner損傷理論分析了輸電塔的風振疲勞特性,預測了風荷載作用下該塔線體系的剩余壽命。

綜上述及,雖然上述研究者對輸電塔線體系的風振響應、風振疲勞特性等做出了卓有成效的研究,為后續研究和實際工程應用均奠定了良好的基礎。但是不難發現上述研究均是在平地風場中對輸電塔線體系的相關性能進行了研究,而對山地風場作用下輸電塔線體系風振響應和風振疲勞等特性的研究尚屬起步階段,資料相對匱乏。然而,實際工程中特高壓輸電塔常建在地勢較高的山坡或山頂,致使整個服役期均處于受山地地形影響的風場中。為此,本文通過風洞試驗研究了在多個山體中,遮擋山體的高度、坡度和間距對受擾山體風場特性的影響規律,同時,將風向風速的統計概率與時域和頻域疲勞累積損傷理論相結合,對特高壓輸電塔線體系進行了風振疲勞分析,研究了其風振疲勞壽命隨遮擋山體坡度、高度和間距的變化規律。

1 風洞試驗

1.1 試驗概況

本次試驗在中國空氣動力研究與發展中心低速所的FL-11風洞中進行,風洞截面尺寸為1.4 m×1.4 m,試驗段長為6.3 m。試驗采取余弦平方型山體模,如圖1所示;模型幾何縮尺比為1∶1 000,風速縮尺比約1∶2.5;試驗中利用DSM3400電子掃描閥系統采集試驗數據,采樣頻率為156 Hz,每個測點采用時間步數為4 096步,排管高0.7 m,沿其高度平均分布20根探針,如圖2所示。

圖1 試驗模型示意圖Fig.1 Schematic diagram of specimens

圖2 試驗測點位置Fig.2 Distribution of test measuring points

試驗中,由不同幾何尺寸的山體模型(如表1所示)模擬遮擋山體(施擾山體)對受擾山體風場的影響。本次試驗從遮擋山體的5種坡度、3種高度和與受擾山體的5種遮擋間距共三類地形分別研究了受擾山體平均風和脈動風加速比ΔS的變化情況,即在平地和山地同一高度z處,山地平均風速U(z)與平地平均風速U0(z)之差與平地平均風速的比值,表達式為:

表1 山體模型及試驗工況Tab.1 Test conditions of hilly terrain wind filed

1.2 試驗結果分析

1.2.1 遮擋山體坡度

試驗中,分別以高度相同坡度不同的DX2~DX6五個模型作為遮擋山體,以模型DX2和DX6為受擾山體。同時,試驗中干擾和受擾山體坡腳相連,即各工況中遮擋間距始終為0。山體模型DX2和DX6在不同坡度的遮擋山體的干擾作用下山頂和背風面山腳風場特性的變化規律如圖3~圖5所示。

不難發現,受擾山頂近地面47 mm到121 mm高度范圍內,平均風速加速比ΔS曲線隨遮擋山體坡度增加而遞增,但隨高度增加ΔS增勢減緩;隨受擾山體坡度增加,ΔS出現最大值所對應的遮擋山體坡度增大,但ΔS的最大值卻逐漸減小。圖4中受擾山體DX2背風面山腳的ΔS在121 mm高度以下為負值,具有明顯的減速效應,而在121 mm到269 mm高度范圍內產生了較小的加速效應;隨坡度增加,121 mm高度以下的近地面ΔS增大,121 mm高度以上的ΔS逐漸減小,而且在遮擋山體坡度為0.75時ΔS均為0。DX6山體背風面山腳的近地面平均風效應逐漸減弱,在遮擋山體坡度為0.375時,ΔS亦均為0。當坡度繼續增大時,加速比變為正數后,并不斷增加。在坡度大于0.75后,ΔS已基本保持不變,而且121 mm以上測點的ΔS隨遮擋山體的坡度增大稍有減小。由圖5中各曲線可知,當受擾山體坡度較大時,背風面脈動風遮擋效應隨遮擋山體坡度增大而顯著增加,脈動風加速比隨遮擋山體坡度減小而增大;受擾山體坡度較小時,雖然遮擋山體坡度較大時遮擋效果亦比較明顯,但脈動風加速比卻隨遮擋山體坡度減小而增大,隨遮擋山體坡度增大而減小。

圖3 山頂平均風加速比隨遮擋山體坡度的變化Fig.3 Effect of slope of occluding hill on speed-up ratio at crest

圖4 背風面山腳平均風加速比隨遮擋山體坡度的變化Fig.4 Effect of slope of occluding hill on speed-up ratio at hill foot in lee side

圖5 背風面脈動風加速比隨遮擋山體坡度的變化Fig.5 Effect of slope of occluding hill on fluctuate speed-up ratio at hill base in lee side

1.2.2 遮擋山體高度

考察遮擋山體坡度對受擾山體風場的影響時,以坡度相同,高度不同的山體模型DX1、DX4、DX7(阻塞度為3.57%)為遮擋山體,并仍以DX2和DX6為受擾山體,各工況中遮擋距離亦為0。

由圖6~圖8可以看出,隨遮擋山體高度增加,受擾山體山頂及背風面山腳僅近地面10 mm高度的ΔS稍有變化,而其余高度的ΔS幾乎相等;同時,脈動風最大加速比亦不隨山體高度發生變化。由此可見,遮擋山體高度對受擾山體的風場影響小。

圖6 山頂平均風加速比隨遮擋山體高度的變化Fig.6 Effect of height of occluding hill on speed-up ratio at crest

圖7 背風面山腳加速比隨遮擋山體高度的變化Fig.7 Effect of height of occluding hill on speed-up ratio at hill foot in lee side

圖8 背風面脈動風加速比隨遮擋山體高度的變化Fig.8 Effect of height of occluding hill on fluctuate speed-up ratio at hill base in lee side

1.2.3 遮擋山體間距

試驗中,對遮擋距離為 0、200 mm、400 mm、600 mm、800 mm的DX2和DX6兩類山體的相同模型分別進行了試驗。

由圖9~圖10可知,不同遮擋間距的兩山體坡度較大時,下風向模型受干擾嚴重,遮擋間距為0時山頂ΔS達到最大值,但隨間距增加ΔS迅速減小,間距大于400 mm和近地高度大于150 mm時ΔS均不再隨間距發生變化;坡度較小時,遮擋效應也相應減小,間距大于600 mm和高度大于150 mm時ΔS才不再發生變化。坡度較大而且間距小時,受擾山體背風面山腳近地面的ΔS均為正值,并且在121 mm高度處最大,但隨遮擋間距增大,近地面4個高度的ΔS迅速減小,并由加速變為減速;高度大于158 mm和間距大于800 mm時,ΔS基本趨于穩定,遮擋效應基本消失。圖11則反映出山體坡度較大且遮擋間距較小時,背風面山腳的脈動風速最大加速比受到一定遮擋干擾,隨遮擋間距增大,最大加速比基本重合。

圖9 山頂平均風加速比隨遮擋距離的變化Fig.9 Effect of occlusion distance on speed-up ratio at crest

圖10 背風面山腳平均風加速比隨遮擋山體間距的變化Fig.10 Effect of height of occluding hill on speed-up ratio at hill foot in lee side

圖11 背風面脈動風加速比隨遮擋山體距離的變化Fig.11 Effect of occlusion distance on fluctuate speed-up ratio at hill base in lee side

2 風向風速的聯合統計概率

重慶市氣象局提供了重慶市1980年1月1日至2010年6月30日的日最大風速值(10 min平均時距)和相應風向記錄(共16個方位,正北為0度方向,順時針方向為正)的實測數據。由此數據,一方面可整理計算出16個風向所對應的風速概率直方圖,并結合對數正態分布、Gumbel分布、Gamma分布、兩參數 Weibull分布等典型概率密度函數可統計出各風向所有風速的分布概率。因篇幅所限,此處僅給出風向2(11.25°~33.76°)、風向6(101.25°~123.76°)、風向10(191.25°~213.76°)以及風向14(281.25°~303.76°)四個方向的結果,如圖12所示。另一方面,由方位記錄法根據風速數據可統計出風向頻度值,如表2所示。最后,結合風向頻度及風速概率即可得到整個計算期內風向風速的聯合統計概率,部分風向風速聯合分布統計概率如表3所示。

圖12 不同風向內風速概率直方圖和典型概率密度函數的比較Fig.12 Comparison with typical distributions for wind speed in wind direction

3 輸電塔線體系風振疲勞計算

3.1 風荷載處理

結構上某點X(x,y,z)的瞬時風速U(X,t)是平均風速(X)與脈動風速U(X,t)={u(X,t),v(X,t),w(X,t)}之和。顯然,風荷載是由平均風引起的靜風力、脈動風引起的抖振力和氣動耦合產生的自激力三部分組成。其中,靜風荷載由模型風洞實驗獲得的靜力三分力系數確定;抖振力按準定常氣動公式計算,不考慮氣動導納的修正;自激力在輸電塔線體系的風振響應亦不予考慮。并且脈動風速一般看作一個零均值的平穩高斯過程,根據準定常假設認為來流在結構表面重分布的影響可以忽略不計[17]。為此,基于準定常假設和不考慮氣動導納的修正時,X(x,y,z)點所受的脈動風荷載僅與該點來流的速度有關,即

表2 風向頻度函數表(%)Tab.2 Frequency function of wind direction(%)

表3 部分風向風速的聯合分布的統計概率(%)Tab.3 The joint statistic probability value of wind speech for each wind direction(%)

式中:ρ、CD(M)、A(M)分別表示空氣密度,結構的截面的阻力系數和結構的迎風面積。

3.2 風振疲勞時域計算方法

首先根據重慶市氣象局提供的重慶地區1980年1月1日至2010年6月30日的風速數據,按荷載規范采用指數型風剖面[18],根據文獻[17]的建議方法即由平均風速和地貌指數生成平地風速功率譜;然后,由前述風洞試驗數據對風速功率譜進行修正,得到相應工況下山地風速功率譜;最后,利用公式(3)和公式(4)將各風速功率譜轉化為氣動力并對輸電塔線體系進行時程分析,得到各風向中不同風速下輸電塔線結構體系中各桿件的風致響應后,通過ANSYS進行對比分析即可尋找出最不利受力桿件(關鍵桿件)的應力時程,由此可計算出各風向中各平均風速模擬出的風場在整個計算期對構件造成的損傷,然后再結合該風速和其所在風向的聯合概率即可求得該風向中該風速對應功率譜對結構在整個計算期內造成的實際損傷,最后再將各風向中所有風速對構件的損傷相疊加即為整個計算期內結構的實際風振疲勞損傷。

具體計算過程為,由上述過程可得各風向中不同風速下輸電塔線結構體系中最不利受力桿件(關鍵桿件)的應力時程[σi,j](風向i=1,2…,16;風速j=1,2,….,26)。但是計算過程中桿件應力均值隨風速變化而變化,因此,疲勞分析時應利用Goodman模型對其進行修正可得到等效應力幅矩陣:

顯然,將式(5)代入下式

可得到各風向中不同風速下輸電塔的疲勞壽命矩陣,即:

同時,經統計分析可得到整個計算期內式(5)中各等效應力幅對應的輸電塔循環次數:

根據Miner累積損傷理論,即:

可得到各風向和風速下桿件的損傷矩陣

然后,表3中根據實測數據得到的風向風速聯合分布統計概率可表示為:

于是,結合式(10)和式(11)可以得到整個計算周期T內輸電塔的損傷:

將上式右邊矩陣中的每一項相加即可得到T時間內總的損傷為:

再由Miner累積損傷理論,即可得輸電塔的最終疲勞壽命:

3.3 風振疲勞頻域計算方法

將風向風速聯合分布的統計概率Pi,j引入理想窄帶疲勞損傷計算方法中[15],對其進行改進。

利用3.2節中模擬的風場對輸電塔線體系進行時程分析,由此可計算出各風向上不同風速中輸電塔線結構體系中最不利受力桿件(關鍵桿件)的窄帶隨機應力響應s(t),應力幅值由s到s+δs發生的循環次數在整個隨機循環中的概率為:

式中:f(s)為應力幅值s所對應的循環次數的概率密度函數,理想窄帶隨機應力響應過程中,因應力峰值s符合Rayleigh分布,則:

若整個計算期T內應力循環次數為NT,那么,由s到s+δs應力循環的次數為:

并且應力幅值s對應的疲勞壽命為N(s),由此可得,該應力作用下結構的累積損傷為:

式中:C和k均為材料參數。

顯然,整個計算期內i風向上j風速的所有應力幅對結構的累積損傷為:

式中:Pi,j為表3中的風向風速聯合分布的統計概率。

將式(13)代入式(16),i風向上的j風速所有應力幅對結構的累積損傷可進一步表示為:

式中,Г(·)為Gamma函數。

不難發現,結合表3中各風向風速的聯合分布的統計概率,在整個計算期T內,窄帶應力響應過程所形成的損傷可表示為:

根據式(21)即可計算出結構在整個計算期T內的損傷后,再采用式(14)亦可得到結構的疲勞壽命。

3.4 計算結果分析

文中基于工程實例及風洞試驗[1-3],取特高壓大跨越輸電塔的總高為108 m,桿件設計強度為345 MPa,構件材料疲勞參數取C為0.41×1012及β為3,服役期為100年,有限元模型如圖13所示。

圖13 輸電塔分段示意圖Fig.13 The schematic diagram of transmission tower

通過時程分析得到如圖13所示部位四根關鍵桿件兩端的應力最大,然后利用3.2和3.3節給出的風振疲勞計算方法計算出各工況下山頂和背風面山腳塔線體系四根最不利受力桿件的平均疲勞壽命,從而預測出各工況下輸電塔線體系的疲勞壽命,如表4所示。

由表4中數據可知,同一高度不同坡度的遮擋山體對輸電塔疲勞壽命存在明顯影響,當坡度較小或較大時受擾山體風場干擾嚴重,輸電塔損傷較大,疲勞壽命較短,且背風面山腳疲勞壽命最短;相同坡度不同高度的遮擋山體對受擾山體風場無明顯干擾,輸電塔疲勞壽命未發生較大變化;遮擋山體間距較小時,下風向受擾山體風場受干擾嚴重,輸電塔疲勞壽命較短,但隨著遮擋間距受擾山體風場逐漸趨于穩定,疲勞壽命隨之延長。同時,對比雨流法(時域)和理想窄帶(頻域)兩種方法的計算結果,不難發現頻域法得到的疲勞損傷較大,疲勞壽命相對保守,與文獻[15]得到的規律相吻合。

綜上述及,高度相同坡度不同的遮擋山體以及不同間距的遮擋山體均會對受擾山體的風場造成一定程度的干擾作用,致使輸電塔損傷最大,疲勞壽命變短。為此,輸電塔選址時應盡量避免風場干擾嚴重的不利位置。

表4 不同工況下輸電塔疲勞壽命(單位:年)Tab.4 Fatigue life of transmission tower under different hilly terrain wind field(year)

4 結論

本文基于風洞試驗對通過13種工況測試了遮擋山體坡度、高度以及間距對受擾山體山頂及背風面山腳風場的影響情況。同時,將風向風速聯合分布的統計概率與疲勞損傷理論相結合,得到了考慮風向風速的疲勞損傷計算方法。在此基礎上,根據試驗數據擬合出了各工況的風場,并對輸電塔進行了風致響應計算和疲勞損傷分析,得到了輸電塔疲勞壽命隨遮擋山體坡度呈先增后減的復雜變化趨勢,且山體背風面山腳疲勞壽命最短;坡度相同而高度不同的遮擋山體未對輸電塔疲勞壽命造成明顯變化;遮擋山體間距較小時輸電塔疲勞壽命最短,但隨遮擋間距增加受擾山體風場逐漸趨于穩定,疲勞壽命隨之延長。

[1]孫 毅.山地風場中超高層建筑風致振動研究[D].重慶:重慶大學,2010.

[2]李正良,孫 毅,黃漢杰,等.山地風場中超高層建筑風荷載風荷載幅值特性試驗研究[J].建筑結構學報,2010,31(6):171-178.LIZheng-liang, SUN Yi, HUANG Han-jie, etal.Experimental research on amplitude characteristics of wind loads of super tall buildings in hilly terrain field[J].Journal of Building Structures,2010,31(6):171-178.

[3]李正良,魏奇科,孫 毅.山地地形對輸電塔風振響應的影響[J].電網技術,2010,34(11):214-219.

LI Zheng-liang,WEI Ji-ke,SUN Yi.Influence of mountain topography on response to wind-induced vibration of transmission tower[J].Power System Technology,2010,34(11):214-219.

[4]李正良,魏奇科,黃漢杰,等.山地超高層建筑風致響應研究[J].振動與沖擊,2011,30(5):43-48.

LI Zheng-liang,WEI Qi-ke,HUANG Han-jie,et al.Windinduced response of super tall building in hilly terrain[J].Journal of Vibration and Shock,2011,30(5):43-48.

[5]Havard D G,Perry O C.Lattice tower member fatigue and its controlusing a noveldamping scheme [J]. Power Engineering Society Summer Meeting,IEEE,Seattle,WA,USA,200,4:2560-2565.

[6]Repetto M P,Solari G.Dynamic alongwind fatigue of slender vertical structure [J].Engineering Structure,2001,23:1622-1633.

[7]Repetto M P,Solari G.Dynamic crosswind fatigue of slender vertical structures[J].Wind and Structures,2002,5(6):527-542.

[8]Repetto M P,Solari G.Directional wind-induced fatigue of slender vertical structure [J]. Journal of Structure Engineering,2004,130(7):1032-1040.

[9]Repetto M P,SolariG.Bimodalalongwind fatigueof structures[J].Journal of Structure Engineering,2006,132(6):899-908.

[10]王之宏.桅桿結構的風振疲勞分析[J].特種結構,1994,11(3):3-8.

WANG Zhi-hong.Analysis of wind-induced fatigue in guyed steel masts[J].Special Structures,1994,11(3):3-8.

[11]屠海明,鄧洪洲.基于頻域的桅桿結構風振疲勞分析[J].特種結構,1999,15(4):34-36.

TU Hai-ming,DENG Hong-zhou.The study of guyed mast fatigue based on stress spectrum under wind vibration[J].Special Structures,1999,15(4):34-36.

[12]鄧洪洲,溫應龍,何鵬飛.格構式桅桿順風向風振疲勞可靠性分析[J].特種結構,2004,21(2):31-35.

DENG Hong-zhou,WEN Ying-long,HE Peng-fei.Analysis of fatigue reliability of the lattice guyed-mast along wind vibration[J].Special Structures,2004,21(2):31-35.

[13]謝 強,嚴承涌.1000 kV特高壓交流同塔雙回輸電塔線耦聯體系風洞試驗[J].高電壓技術,2010,36(4):900-906.

XIE Qiang,YAN Cheng-yong.Wind tunnel test on 1000kV UHV AC double circuit transmission tower-conductor coupling system[J].High Voltage Engineering,2010,36(4):900-906.

[14]李正良,肖正直,韓 楓,等.1000 kV漢江大跨越特高壓輸電塔線體系氣動彈性模型的設計與風洞試驗[J].電網技術,2008,32(12):1-5.

LIZheng-liang,XIAO Zheng-zhi,HAN Feng, etal.Aeroelastic model design and wind tunnel tests of 1000 kV Hanjiang long span transmission line system [J].Power System Technology,2008,32(12):1-5.

[15]汪之松.特高壓輸電塔風振響應及風振疲勞性能研究[D].重慶:重慶大學,2009.

[16]汪之松,李正良,肖正直,等.1 000 kV雙回路特高壓輸電塔等效靜風荷載研究[J].電網技術,2009,33(14):6-12.

WANG Zhi-song,LI Zheng-liang,XIAO Zheng-zhi,et al.Equivalent along-wind direction static wind load of transmission tower for 1000 kV double circuit transmission lines on same tower[J],Power System Technology,2009,33(14):6-12.

[17]肖正直.特高壓輸電塔風振響應及等效風荷載研究[D].重慶:重慶大學,2009.

[18]GB 50009-2001建筑結構荷載規范[S].北京:中國建筑工業出版社,2006.

猜你喜歡
風速
邯鄲市近46年風向風速特征分析
基于Kmeans-VMD-LSTM的短期風速預測
基于最優TS評分和頻率匹配的江蘇近海風速訂正
海洋通報(2020年5期)2021-01-14 09:26:54
基于時間相關性的風速威布爾分布優化方法
陜西黃土高原地區日極大風速的統計推算方法
陜西氣象(2020年2期)2020-06-08 00:54:38
基于GARCH的短時風速預測方法
快速評估風電場50年一遇最大風速的算法
風能(2016年11期)2016-03-04 05:24:00
考慮風切和塔影效應的風力機風速模型
電測與儀表(2015年8期)2015-04-09 11:50:06
GE在中國發布2.3-116低風速智能風機
考慮風速分布與日非平穩性的風速數據預處理方法研究
主站蜘蛛池模板: 国产AV毛片| 精品人妻无码中字系列| 青青草a国产免费观看| 色妞www精品视频一级下载| 精品一區二區久久久久久久網站| 国产又色又爽又黄| 免费看的一级毛片| 国产又粗又爽视频| 国产91小视频在线观看| 免费国产好深啊好涨好硬视频| 91精品国产自产91精品资源| 久久精品无码中文字幕| 精品三级在线| 天天综合网站| 中文字幕久久亚洲一区| 丁香六月激情婷婷| 天堂va亚洲va欧美va国产 | 国产精品成人免费综合| 日韩高清欧美| 99精品欧美一区| 国产在线一区视频| lhav亚洲精品| 国产日本欧美在线观看| 成人精品区| 中文字幕在线看| 国产综合色在线视频播放线视| 久久精品视频一| 国产成人福利在线| 国产美女精品人人做人人爽| 男女猛烈无遮挡午夜视频| 日韩AV无码免费一二三区| 国产00高中生在线播放| 97精品国产高清久久久久蜜芽| 国产成人精品在线| 欧美国产在线一区| 欧美亚洲日韩中文| 国产精品色婷婷在线观看| 71pao成人国产永久免费视频| 伊人成人在线| 狼友视频一区二区三区| 又粗又硬又大又爽免费视频播放| 国产精品香蕉在线| 伊人天堂网| 久久久久国产精品嫩草影院| 久久这里只有精品66| 欧美福利在线观看| 无码aⅴ精品一区二区三区| 666精品国产精品亚洲| 色AV色 综合网站| 久久久久久国产精品mv| 欧美成人一级| 日本一区中文字幕最新在线| 毛片在线区| 美女国产在线| 婷婷色婷婷| 女高中生自慰污污网站| 99久久精品免费观看国产| 午夜色综合| 亚洲AV无码久久精品色欲| 伊人久久大香线蕉aⅴ色| 高清免费毛片| 三级国产在线观看| 中文字幕免费播放| 91精品国产综合久久不国产大片| 精品伊人久久久香线蕉| 国产一区二区精品高清在线观看| 国产综合在线观看视频| 亚洲精选无码久久久| 视频一本大道香蕉久在线播放| 亚洲综合专区| 欧美视频二区| 国产精品偷伦在线观看| 青青草原国产| 性色在线视频精品| 国产在线一区二区视频| 992tv国产人成在线观看| 四虎国产成人免费观看| 九九线精品视频在线观看| 亚洲精品不卡午夜精品| 久久午夜夜伦鲁鲁片不卡| 亚洲乱亚洲乱妇24p| 亚洲成人播放|