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

基于設計基準條件的臺風邊界層脈動風速譜建模方法研究

2016-01-15 02:09:33李利孝肖儀清周超英宋麗莉
振動與沖擊 2015年11期

李利孝,肖儀清,周超英,宋麗莉

(1.哈爾濱工業大學深圳研究生院,廣東 深圳 518055; 2.公共氣象服務中心中國氣象局,北京 100081)

第一作者李利孝男,博士,1984年生

基于設計基準條件的臺風邊界層脈動風速譜建模方法研究

李利孝1,肖儀清1,周超英1,宋麗莉2

(1.哈爾濱工業大學深圳研究生院,廣東深圳518055; 2.公共氣象服務中心中國氣象局,北京100081)

摘要:利用在博賀峙仔島觀測站實測的臺風黑格比數據,在莫寧-奧布霍夫相似理論框架內提出了一種數據驅動的考慮設計基準條件的臺風邊界層脈動風速譜建模方法。首先基于均勻各向同性湍流能譜需滿足的基本模型準則,將六參數脈動風速譜統一模型簡化為四參數譜模型;之后基于高頻區和低頻區脈動風速譜特征建立了四個譜參數的理論表達式,并利用實測的臺風黑格比數據將譜參數理論表達式進行參數化處理,進而提出基于設計基準條件的脈動風速譜建模方法;最后利用臺風黑格比過程中任意選取的四個1h時長樣本對上述脈動風速譜建模方法進行了驗證,表明其與實測樣本的風速譜具有較好的一致性。

關鍵詞:臺風;湍流;脈動風速譜;設計基準條件;結構抗風設計

基金項目:國家自然科學基金(51308168和51278161);中國博士后科學基金項目(2013M531045)

收稿日期:2014-01-08修改稿收到日期:2014-05-30

通信作者肖儀清男,博士,教授,1973年生

中圖分類號:TU973+.32

文獻標志碼:A

DOI:10.13465/j.cnki.jvs.2015.11.003

Abstract:A data-driven model was proposed to model wind velocity spectrum of typhoon boundary layer considering conditions of design reference. This model was established based on Moning-Obukhov similarity theory (MOST) and the field measured data of typhoon Hagupit on Zhizai Island, Bohe. Firstly, the six-parameter general model was simplified into a four-parameter spectral model by matching the model criteria of energy spectrum in homogeneous turbulent flow. Secondly, the theoretical expressions of the four spectral parameters were determined based on the features of wind velocity spectrum in lower and higher frequency ranges. Then the data-driven model considering conditions of design reference was established through parameterizing the four spectral parameters with the field measured data of typhoon Hagupit. Finally, the validation of the data-driven model was examined with arbitrarily selected four 1-hour length data sets in the process of typhoon Hagupit. It was shown that the wind velocity spectrum modeled using the proposed method agrees well with that obtained using the actual measured samples.

Modeling method for fluctuating wind velocity spectrum in typhoon boundary layer based on conditions of design reference

LILi-xiao1,XIAOYi-qing1,ZHOUChao-ying1,SONGLi-li2(1. Shenzhen Graduate School, Harbin Institute of Technology, Shenzhen 518055, China;2. Public Meteorological Service Center, China Meteorological Administration, Beijing 100081, China)

Key words:typhoon; turbulence; wind velocity spectrum; conditions of design reference; wind-resistance design

大量超高層建筑和超大跨度橋梁結構在我國東南沿海地區興建,每年遭受若干次臺風的襲擊,使得抵抗臺風作用成為這類結構在服役期安全性和舒適性的巨大挑戰。準確刻畫臺風風場結構特征和湍流特性是估計臺風風荷載和抗臺風設計的首要保證。風速剖面刻畫了平均風速隨高度的變化規律,是估算平均風荷載的重要參數,而風速譜則表征了湍動能在不同頻率區段上的分布特征,是刻畫脈動風荷載的關鍵參數[1-2]。

目前國內外大部分風荷載規范所采用的經驗風速譜模型[3-5]是利用非氣旋風場實測數據基于均勻各向同性湍流的能譜特征推導所得,然而在臺風場中,對流湍流和下沉作用將風場上部形成的大尺度渦漩向下輸運并影響近地面層的流場特性,使得臺風風場與常態風風場在結構特征和湍流特性均存在一定的差異[6-9]。針對臺風場脈動風速譜特征,國內外學者開展了大量的實測研究,部分學者的觀測發現臺風場脈動風速譜包含較多的低頻區能量,含能區特征尺度大于常態大風風場的特征尺度[10-11];而另一部分學者則獲得與之相反的結論[12-13]。鑒于此,部分學者采用最小二乘法對實測臺風場數據進行擬合以獲得適用于臺風影響區的脈動風速譜模型[14-15]。然而,采用這種擬合方法所建立的譜模型的參數變異性較大,對實測數據依賴性較高;同時所建立的譜模型不能直觀反映風場環境、粗糙度類別、高度以及風速等結構設計基本條件的影響,也不能反映大氣層結穩定度等氣象學要素對譜特性的影響[16]。

因此,本文將基于在博賀峙仔島觀測站實測的臺風黑格比數據,提出一種數據驅動的基于設計基準條件的臺風邊界層脈動風速譜建模方法,并利用臺風黑格比過程中任意選取的四個1h時長的樣本對上述建模方法的準確度進行驗證,從而為臺風影響區建筑結構的抗風設計提供借鑒和指導。

1脈動風速譜理論模型建立

均勻各向同性湍流的湍能譜可劃分如下三個區域:含能區、慣性子區和能量耗散區(見圖1)。含能區主要通過大尺度漩渦的脈動從平均流中獲取能量,其包含了流場絕大部分的湍動能;能量耗散區則主要是通過小尺度漩渦的相互摩擦碰撞作用進行能量的耗散;而在慣性子區,假定漩渦主要進行能量的傳遞作用,既不從平均流中獲取能量,也不進行能量的耗散[17]。所以在慣性子區,歸一化的脈動風速譜服從各向同性湍流的-2/3次律變化規律[18],

(1)

式中:Su(n)為脈動風速譜譜;u*為摩擦速度;φε為無量綱化風切力的莫寧奧布霍夫系數函數;Au為常數,取為0.27;n為自然頻率(Hz);f為無量綱化折算頻率,f=nΛ/U,其中Λ可為離地高度,積分尺度或一固定值。

基于慣性子區的脈動風速譜一般形式(式(1)),大量大氣中性層結下的經驗風速譜模型在此基礎被提出并廣泛應用于不同國家的抗風設計規范中,這些經驗譜模型可統一表達為六參數脈動風速譜廣義模型,即

(2)

式中:A、B、C、α、β和γ為六個待定參數;R為湍流比。

圖1 均勻各向同性湍流湍能譜示意圖Fig.1 Schematic of energy spectrum of homogeneous and isotropic turbulence

1.1模型準則

基于莫寧-奧布霍夫相似理論和均勻各向同性湍流的基本特征,脈動風速譜統一模型(式2)需滿足以下幾個基本準則:

(1)在慣性子區脈動風速譜需服從Kolmogorov的-2/3次律能量傳遞理論,即

αβ-γ=2/3

(3)

γ=1;A=4Cβ

(4)

(3)當頻率n趨近于0時,脈動風速譜Su(n)的導數亦趨近于0,由此可得

α≥1

(5)

將式(3)~式(5)代入式(2),則六參數脈動風速譜廣義模型可簡化為四參數譜模型,即

(6)

式中:A、B、C和α為四個待定參數。

1.2譜參數理論模型

無量綱化的順風向脈動風速譜nSu(n)/u*2在低頻區服從折算頻率f的+1次律變化規律,在慣性子區則服從折算頻率f的-2/3次律變化規律。因此無量綱化的脈動風速譜則可通過如下四個譜參數來確定:無量綱化脈動風速譜在譜能量達到最大值時對應的折算頻率fm和譜能量Gm(Gm=max(nSu(n)/σ2)),無量綱化風速譜在慣性子區系數Au,歸一化脈動風速譜的湍流比R。

上述四個譜參數則可由四參數譜模型(式(6))的四個待定參數來表達:

(1)無量綱化脈動風速譜在譜能量達到最大值時對應的折算頻率fm可通過對式(6)右邊部分求導得,

(7)

(2)歸一化風速譜在雙對數坐標系中的最大譜密度Gm則可表示為

(8)

(3)在慣性子區,脈動風速譜系數Au可由下式計算

(9)

(4)由順風向脈動風速譜所包含的能量等于順風向脈動風速的方差可得

(10)

式中

s≠0, -1, -2…

(11)

由上述分析可知,當描述風速譜形的四個譜參數確定之后,即可通過聯解式(7)~式(10)進而反算求得四參數譜模型(式(6))中的四個待定參數A、B、C和α,從而可建立一個相應的脈動風速譜的表達式。

2脈動風速譜參數模型建立

2.1數據來源

本文所采用的數據為2008年第14號臺風黑格比在博賀峙仔島觀測站實測得臺風資料。強臺風黑格比是自1996年以來登陸廣東的最強臺風,于2008-09-24T06∶45在廣東省西部電白縣陳村附近登陸,登陸時中心最大風力15級(48m/s)。強臺風黑格比數據是由廣東省氣象局在博賀峙仔島100 m觀測塔進行現場實測的所得,觀測站場地條件以及儀器設置等可見參考文獻[19-20]。

圖2所示為峙仔島觀測站60 m高度的超聲風速儀實測的強臺風黑格比的風速風向時程,由圖可見強臺風黑格比風眼經過了觀測站,在中心經過前后的風向轉角約為191°,臺風中心經過前的最大10min平均風速為45.88 m/s,臺風中心經過后的最大10min平均風速為40.11 m/s。實測的臺風黑格比144個樣本經過了數據質量控制和樣本選取準則處理后[20],最終有47個滿足平穩性要求的樣本用于本文的分析。

圖2 臺風黑格比10分鐘平均風速風向時程Fig.2 Time history of 10 min mean wind speeds and directions in typhoon Hagupit

2.2譜參數經驗模型建立

由1.2節的分析可知,只要確立了四個譜參數,即可通過式(7)~式(10)來確定參數A、B、C和α,本節將基于上節所選取的47個數據樣本,采用最小二乘法擬合建立各譜參數的經驗模型。由于大氣穩定度對風速譜在低頻區的能量分布較大程度的影響,在本節的分析中將同時考慮大氣層結穩定度對譜參數的影響。

在莫寧-奧布霍夫相似理論框架內,大氣層結穩定度可通過無量綱化的大氣穩定度系數來表達,即

(12)

在式(12)中,當z/L>0時為穩定層結,z/L<0時為不穩定層結,z/L=0時為中性層結。由于z/L=0的樣本基本上沒有,所以在本文中假定-0.1

基于莫寧-奧布霍夫相似理論分析時認為湍流比R在近面地層不隨高度發生變化,但隨著場地粗糙度的改變而改變,同時其還受到大氣層結穩定度的影響。在本文的分析中,先選取中性層結樣本研究湍流比隨粗糙度長度的變化規律,之后再分析湍流比隨大氣層結穩定的變化關系。

圖3(a)為實測的中性層結樣本中湍流比隨粗糙度長度的變化規律,為了比較臺風風場中風眼經過前后區域的湍流比差異,將所選樣本又劃分為中心經過前和中心經過后兩個類別,分別采用最小二乘擬合得湍流比隨粗糙度長度的關系如下兩式所示

TR(z0)=2.74-0.17ln(z0)

(13)

TR(z0)=0.65-0.28ln(z0)

(14)

圖3(b)為大氣穩定度參數對湍流比的影響,在分析大氣穩定度對湍流比的影響時,先利用式(13)扣除粗糙度長度對湍流比的貢獻,之所以選擇式(13)是因為在臺風中心經過后所選取到的近中性層結樣本基本上均位于外環流區域,不能代表臺風風場特性。疊加粗糙度長度對湍流比的影響和大氣穩定度對湍流比的影響之后,湍流比可表達為

R(z0,z/L)=2.74-0.17ln(z0)+0.93(z/L)1/3

(15)

圖3 湍流比TR隨粗糙度長度z0和大氣穩定度參數z/L的變化關系Fig.3 Variation of turbulence ratio with roughness length z0and atmospheric stability z/L

類似于湍流比,折算頻率fm也采用了相應的分析方法進行分析,圖4(a)所示為譜密度峰值Gm對應的自然頻率nm隨平均風速U的變化關系。在分析二者之間的關系時,只選用了近中性層結的樣本,利用最小二乘回歸得二者之間的關系為

nm=0.024 4exp(0.006 6U)

(16)

圖4(b)所示為大氣層結穩定度對頻率nm的影響。類似于分析大氣穩定度對湍流比的影響,在這里亦先利用式(16)扣除了平均風速對頻率nm的影響,綜合平均風速以及大氣穩定度對頻率nm的影響,可得如下表達式,

nm=0.024 4exp(0.006 6U)+

(17)

與之相對應的折算頻率fm則可表示為

fm=U-1z(0.024 4exp(0.006 6U)+0.020 8z/L)

(18)

式中:U為平均風速,0≤z/L≤0.4。

圖4 頻率nm隨平均風速U和大氣穩定度參數z/L的變化關系Fig.4 Variation of natural frequency nmwith mean wind speed U and atmospheric stability z/L

實測的歸一化風速譜的譜密度最大值Gm沒有表現出明顯的隨粗糙度長度或大氣穩定度的變化規律,所以在這里只分析了Gm隨平均風速U的變化關系,見圖5。利用最小二乘法擬合得Gm與平均風速U的變化關系如下式所示

Gm=0.748 5/U0.233 4

(19)

圖5 譜密度最大值Gm隨平均風速U的變化關系Fig.5 Variation of the maximum of normalized wind velocity spectra, Gm, with mean wind speed U

圖6為慣性子區脈動風速譜參數Au隨平均風速U的變化關系。由圖6可見,臺風風場順風向風速譜譜參數的實測值較為離散且大于理論推導的值0.26。譜參數的實測值是基于式(1)利用最小二乘法擬合慣性子區的實測風速譜而得到的,譜參數的實測值與平均風速的關系可用下式表示

Au=0.011 3U1.177 0+0.177 8,U≥17.2 m/s

(20)

圖6 慣性子區譜參數Au隨平均風速U的變化關系Fig.6 Variation of Auin inertial subrange with mean speed U

利用上述各譜參數的經驗模型(式(15)、式(18)、式(19)和式(20)),基于設計基準條件(設計基準風速、場地類別和大氣穩定度條件)可確定各譜參數的具體取值,進而利用1.2節所建立的理論模型(式(7)~式(10))反算得四參數風速譜模型(式(6))的四個待定參數,最后可確定唯一的風速譜表達式。

3模型有效性驗證

3.1數據選取

本文2.2節的譜參數模型是建立在平穩樣本的基礎上,為了使上述脈動風速譜建模方法不失一般性,在這里任意選取了強臺風黑格比風場不同部位(前外環流區FOV、前眼壁區FEW、后眼壁區BEW和后外環流區BOV)的四個時長為1h的樣本進行分析,所選的四個1h樣本的詳細湍流統計特性列于表1中。

3.2模型有效性驗證

基于表1所列的各樣本的平均風速及其湍流統計參數,分別采用式(15)、式(18)、式(19)和式(20)計算了各樣本的湍流比TR、折算頻率fm、歸一化風速譜譜密度最大值Gm和慣性子區譜參數Au(見表2)。

表1 四個1h時長樣本的湍流統計特性

表2 四個1h時長樣本的譜參數

基于表2所計算得譜參數,采用“1.2”節所建議的理論譜模型建立了強臺風黑格比不同部位的脈動風速譜表達式,如下各式所示:

(21)

(22)

(23)

(24)

圖7 ((a)、(b)、(c)和(d))分別給出了強臺風黑格比四個區域的實測風速譜、理論風速譜(式(21)~式(24))、von Karman譜和Kaimal譜進行對比。由圖7可見,采用本文建立的方法確定的風速譜(圖例中稱作理論風速譜)與實測風速譜吻合的非常好,尤其是在眼壁強風區(FEW和BEW)。在外環流區域(FOV和BOV),理論風速譜與實測風速譜有一定的差異,這可歸結為如下三個原因:第一個是由于樣本的非平穩性導致的,四個1h時長樣本是任意選取的,有些樣本不滿足平穩性要求;第二個可歸結為大氣穩定度的影響,本文在確定風速譜模型時選用了1h總樣本的穩定度參數,但是在1h內的6個10min時長的子樣本則有的是穩定層結,有的是不穩定層結;第三個原因可歸結為建立譜參數經驗模型的數據量相對較少,準確性需進一步提高。

圖7 實測風速譜與理論譜以及經驗風速譜對比Fig.7 Comparison of field measured spectra, theoretical model spectra and empirical spectra

4結論

通過對大氣邊界層順風向脈動風速譜統一表達式的分析,在莫寧-奧布霍夫相似理論的框架內和局部均勻各向同性湍流理論的基礎上提出了一種基于設計基準條件的臺風風場邊界層風速譜建模方法,并驗證了該方法的準確性,得出如下主要結論:

(1)通過對順風向脈動風速譜需滿足的基本條件分析,從理論上將六參數風速譜模型簡化為四參數風速譜模型,并提出了基于四個譜參數的建模方法;

(2)利用實測的強臺風黑格比的數據,建立了譜參數的數學模型,將其參數化為平均風速、粗糙度長度、大氣穩定度參數的函數;

(3)利用任意選取的四個一小時時長的臺風風場數據,驗證了該建模方法所確立的脈動風速譜模型的準確性,發現該理論模型與實測風速譜吻合較好,優于各經驗風速譜模型與實測臺風風速譜的一致性。

本文所提出臺風脈動風速譜建模方法是一種數據驅動模型,對譜參數經驗模型的準確性依賴度較高,所以有必要繼續擴大數據量以更為準確的模擬各譜參數經驗模型,更為精確的確立不同設計基準條件的風速譜模型。另外本文所建立的風速譜建模理論亦適用于豎向脈動風速譜的建模分析和非氣旋風場的風速譜建模。

參考文獻

[1]Li L, Kareem A, Xiao Y, et al. Wind profile and spectra in typhoon-prone regions in south China [C]. Miami, Florida, America: American Society of Civil Engineers, 2012.

[2]李利孝,肖儀清,宋麗莉,等. 基于風觀測塔和風廓線雷達實測的強臺風黑格比風剖面研究[J]. 工程力學,2012, 29(9): 284-293.

LI Li-xiao, XIAO Yi-qing, SONG Li-li, et al. Study on wind profile of typhoon Hagupit using wind observed tower and wind profile radar measurements[J]. Engineering Mechanics,2012, 29(9): 284-293.

[3]Kaimal J C, Wyngaard J C, Izumi Y, et al. Spectral characteristics of surface-layer turbulence[J]. Quarterly Journal of the Royal Meteorological Society,1972, 98(417): 563-589.

[4]Davenport A G. The spectrum of horizontal gustiness near the ground in high winds[J]. Quarterly Journal of the Royal Meteorological Society,1961, 87(372): 194-211.

[5]Karman T V. Progress in the statistical theory of turbulence[J]. Proc Natl Acod Sci U S A.,1948,34(11):530-539.

[6]Bradbury W, Deaves D M. The dependence of gust factor probabilities on convective activity: Analysis conducted for Eurotunnel[J]. Meteorological Applications,1994, 1(2): 159-164.

[7]Sparks P R, Huang Z. Gust factors and surface-to-gradient wind-speed ratios in tropical cyclones[J]. Journal of Wind Engineering and Industrial Aerodynamics,2001, 89(11/12): 1047-1058.

[8]Powell M D, Houston S H. Hurricane andrew’s landfall in south florida. part Ⅱ: surface wind fields and potential real-time applications[J]. Weather and Forecasting,1996, 11(3): 329-349.

[9]Li L, Xiao Y, Song L. Observed sub-hectometer-scale boundary layer rolls in surface layer in landfalling typhoons[C]//Robinson College, Cambridge, UK,2013.

[10]Schroeder J L, Smith D A. Hurricane Bonnie wind flow characteristics as determined from wemite[J]. Journal of Wind Engineering and Industrial Aerodynamics. 2003, 91(6): 767-789.

[11]Yu B, Gan Chowdhury A, Masters F. Hurricane wind power spectra, cospectra, and integral length scales[J]. Boundary-Layer Meteorology,2008, 129(3): 411-430.

[12]Li L, Xiao Y, Kareem A, et al. Modeling typhoon wind power spectra near sea surface based on measurements in the south China sea[J]. Journal of Wind Engineering and Industrial Aerodynamics,2012, 104-106: 565-576.

[13]Zhang J A. Spectral characteristics of turbulence in the hurricane boundary layer over the ocean between the outer rain bands[J]. Quarterly Journal of the Royal Meteorological Society,2010, 136(649): 918-926.

[14]胡曉紅,葛耀君,等. 上海“派比安”臺風實測結果的二維脈動風譜擬合[J]. 結構工程師, 2002(2): 41-47.

HU Xiao-bong, GE Yao-jun, PANG Jia-bin. The fitting of 2D fluctuating wind spectrum for the actual measurement data of “pha” typhoon in Shanghai[J].Structural Engineers,2002(2): 41-47.

[15]黃鵬,戴銀桃,王旭,等. 上海沿海地區近地風脈動風速譜及相干性研究[J]. 工程力學, 2014(4): 126-133.

HUANG Peng, DAI Yin-tao, WANG Xu,et al. Field measurement on power spectra and coherence characteristics of near-ground wind in shanghai coastal area[J].Engineering Mechanics, 2014(4): 126-133.

[16]李利孝. 熱帶氣旋邊界層風場剖面和工程湍流特性研究[D]. 哈爾濱:哈爾濱工業大學, 2012.

[17]Kaimal J C, Finnigan J J. Atmospheric boundary layer flows: their structure and measurement[M]. Oxford, New York: Oxford University Press, 1994.

[18]Simiu E, Scanlan R H. Wind effects on structure: fundamentals and applications to design[M]. Now York: John Wiley & Sons, 1996.

[19]Song L, Li Q S, Chen W, et al. Wind characteristics of a strong typhoon in marine surface boundary layer[J]. Wind and Structures,2012, 15(1): 1-15.

[20]肖儀清,李利孝,宋麗莉,等. 基于近海海面觀測的臺風黑格比風特性研究[J]. 空氣動力學學報, 2012, 30(3): 380-387.

XIAO Yi-qing, LI Li-xiao, SONG Li-li,et al.Study on wind characteristics of typhoon Hagupit based on offshore sea surface measurements[J].Acta Aerodynamica Sinica, 2012, 30(3): 380-387.

主站蜘蛛池模板: 国产精品女熟高潮视频| 91成人试看福利体验区| 在线国产你懂的| 40岁成熟女人牲交片免费| 国产免费人成视频网| 欧美天堂在线| 欧美日韩动态图| 色播五月婷婷| 在线观看免费AV网| 少妇精品网站| av在线5g无码天天| 国产在线97| 黄色免费在线网址| 99视频精品全国免费品| 成年人免费国产视频| 久久精品视频一| 国产精品hd在线播放| 国产精品专区第一页在线观看| 欧美成人精品在线| 国产玖玖视频| 97久久精品人人| 色悠久久综合| 亚洲天堂视频在线免费观看| 呦视频在线一区二区三区| 中文毛片无遮挡播放免费| 亚洲成a人在线观看| 国内嫩模私拍精品视频| 潮喷在线无码白浆| 在线a网站| 高潮毛片免费观看| 中国毛片网| 亚洲资源站av无码网址| 久久精品aⅴ无码中文字幕 | 国产一区在线观看无码| 午夜啪啪网| 午夜日韩久久影院| 精品综合久久久久久97超人该| 国产乱子伦手机在线| 国产在线欧美| 国产福利观看| 色综合手机在线| hezyo加勒比一区二区三区| 91年精品国产福利线观看久久| 美女被操91视频| 热久久这里是精品6免费观看| 国产人成午夜免费看| 国产在线无码一区二区三区| 国产综合欧美| 欧美国产综合视频| 人妻精品全国免费视频| 日韩一级二级三级| 99er精品视频| 1024你懂的国产精品| 亚洲开心婷婷中文字幕| 99视频在线免费| 中文字幕资源站| 国产精品久线在线观看| 91精品国产91久久久久久三级| 久久精品66| 91在线视频福利| 亚洲资源站av无码网址| 幺女国产一级毛片| 日本免费福利视频| 成年人国产视频| 在线观看亚洲天堂| 亚洲天堂.com| 成人精品午夜福利在线播放| 99人妻碰碰碰久久久久禁片| 成年午夜精品久久精品| 亚洲无码精品在线播放| 女人18毛片久久| 亚洲欧洲日韩综合| 中国黄色一级视频| 强乱中文字幕在线播放不卡| 最新国产在线| 精品无码人妻一区二区| 不卡的在线视频免费观看| 九九久久精品免费观看| 亚洲精选无码久久久| 中文字幕亚洲无线码一区女同| 国产18在线| 国产精品自在在线午夜|