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

群性波列作用下系泊船舶的水動力響應數值模擬

2015-04-25 10:11:48馬小劍孫昭晨張志明楊國平1江蘇科技大學船舶與海洋工程學院江蘇鎮江100大連理工大學海岸和近海工程國家重點實驗室遼寧大連11604中交水運規劃設計院有限公司北京100007
船舶力學 2015年7期

馬小劍,孫昭晨,劉 思,張志明,楊國平,周 豐(1江蘇科技大學 船舶與海洋工程學院,江蘇 鎮江 100;大連理工大學 海岸和近海工程國家重點實驗室,遼寧 大連11604;中交水運規劃設計院有限公司,北京100007)

群性波列作用下系泊船舶的水動力響應數值模擬

馬小劍1,2,孫昭晨2,劉 思2,張志明3,楊國平3,周 豐3
(1江蘇科技大學 船舶與海洋工程學院,江蘇 鎮江 212003;2大連理工大學 海岸和近海工程國家重點實驗室,遼寧 大連116024;3中交水運規劃設計院有限公司,北京100007)

該文根據隨機波群的模擬方法,模擬得到具有同一有效波高,有效周期但不同群性的隨機波列,研究了不同群高參數GFH和群長參數GLF的波列作用下系泊船的運動響應及纜繩張力。數值計算結果表明,系泊船橫移運動量和纜繩張力呈現隨波浪群高的增大而增大的趨勢,升沉運動量也略有增加,而波浪群長的變化對系泊船的運動響應及纜繩張力基本沒有影響。

波群;系泊船;群高因子;群長因子

0 引 言

系泊船舶在波浪作用下的運動響應直接關系到碼頭結構、船舶以及貨物裝卸作業的安全,對其準確確定有著非常重要的意義。目前系泊浮體的相關設計常采用代表波法,即選某一特征波研究其對建筑物、系泊浮體等的作用力及響應,而在系泊船的運動響應及纜繩張力研究時僅僅通過模擬特征波要素和波譜來模擬不規則波,沒有考慮波浪群性對透空式碼頭前系泊船的運動量響應的影響。很多學者針對隨機波列作用下系泊船舶的運動響應做了很多工作[1-6]。

已有研究表明,系泊浮體結構所受波浪力及其運動響應,在連續的超過設計波高波列作用時,相比單個設計波高波浪將增加。海浪的群性是波浪場的一重要特征,Nolte和Hsu[7]將隨機波列中一群連續的高波稱為波群,其對海上系泊浮體結構、防波堤的動力反應有直接影響。Johnson等[8]對波群對海上防波堤的破壞作用作了研究。Murray等[9]研究了波群對漂浮方箱駁船慢漂運動的影響。Balaji[10]對波群作用下浮標運動響應作了物理模型實驗研究,實驗結果表明波浪群性對浮標運動響應影響較大。Xu 等[11]將波包譜的概念用于波群模擬中,以給定波譜和波包譜為靶譜的方法模擬波群,認為波群的特征可以由頻譜和波譜和波包譜兩者決定。

本文采用模擬波浪頻譜和波包譜的方法,模擬出具有同一有效波高、有效周期但不同群性的隨機波列,以探討不同群高參數GFH和群長參數GLF的波列作用下系泊船的運動響應和纜繩張力。

1 數值模型

1.1 坐標系

為了描述系泊船舶的運動,建立如圖1所示兩個坐標系:靜坐標系xyz和動坐標系x′y′z′,兩個坐標的原點均取在系泊船的重心位置,兩坐標系在初始狀態時重合,系泊船舶的線位移運動為隨船坐標系x′y′z′相對固定坐標系xyz的移動,定義為△=(△X,△Y,△Z),船舶的角度位移為船舶坐標軸x′y′z′繞相應固定坐標系xyz軸的轉角,定義為θ= (θx,θy,θz)。

圖1 坐標系Fig.1 Coordinate system

1.2 運動方程

船舶在碼頭前系泊時受到護弦和纜繩的約束,以上兩者的作用力是非線性的,而且在船舶的運動過程中會出現纜繩松弛的現象,因此,系泊船舶的運動也是非線性的,此時的運動響應問題不能用頻域方法求解,而需要在時域內求解。時域方法有間接時域法和直接時域法兩種,考慮到計算量,本文采用了間接時域法,即利用頻域下的水動力系數通過Fourier變換在時域內求解系泊船舶的非線性運動方程。

其中:x(t)為船體的位移向量(縱移,橫移,升沉,橫搖,縱搖,回轉),M為船體質量矩陣,a∞為船體的附加質量矩陣,C為船體的靜水恢復力系數矩陣,k(t- τ)為船體遲滯函數矩陣,船體附加質量、阻尼系數和頻域波浪力采用三維源匯分布法計算[12-13],)為船體所受外力,f(t)=波浪力+纜力+護弦力。下面依次介紹這三種外力的計算方法。

1.3 波浪力的計算

波浪力的時間歷程由對時域內波浪力的脈沖響應函數與波面做卷積求得:

其中:

其中:Hi(ω)為單位振幅規則波作用下物體的一階波浪激振力,hi(ω)為時域內一階脈沖響應函數,Fi(t)為一階時域波浪力,η1(t)為隨機波面,本文考慮波面群性大小,按照文獻[11]的方法模擬得到給定群性的波浪,該方法簡述如下:

目前生成隨機波群信號時一般假定波譜為窄譜,首先由Jonswap譜公式:

其中:α為能量尺度參量,γ為譜峰升高因子,ωm為譜峰頻率,σ為峰形參數,ω為頻率,g為重力加速度,s(ω)為頻譜。生成靶譜根據波浪疊加法產生一隨機波列波面(4):

其中:η(t)為波面,t為時間,ε為初相位,△ωi為頻率間隔。然后通過對由Jonswap譜得到的隨機波面進行希爾特變換得到η′(t),從而依據下式:

得到波包線的相位函數φ(t),根據給定的群高因子GFH(波包線的方差和均值之比)和群長因子GLF(波譜峰頻和波包譜峰頻之比):

其中:σA和A′(t)分別為波包線的方差和均值,GFH和GLF分為群高因子和群長因子,GFH是俞聿修[14-16]建議描述群高的參數,GFH趨于零時波列為規則波,fp和fPA為波譜的峰頻和波包譜峰頻,所以波包線的周期相比波周期越大,則波群越長。本文采用文獻[17]中建議的經驗波包譜作為靶譜,即(9)式計算波包譜SA(f):

當GFH>0.7,fP/fPA=10~28,當GFH<0.7,fP/fPA=5~15。同理通過(4)式可以產生一隨機過程,波包線的均值A′由波浪頻譜得到,因此波包線為最后由波包線A1(t)和相位函數φ(t)依據(11)式可構造出波群波列η1(t):

得到波群波面后對其逐步進行修正,使得波浪靶譜和波包靶譜與模擬頻譜和波包譜吻合良好,如圖2 (a)和圖2(b)所示。

根據上述方法及步驟,模擬出具有不同群高及群長特征的波列。圖2(a)、(b)給出了其中一組模擬波列的計算頻譜、計算波包譜與相應的靶譜的對比圖,從圖中可知靶譜和模擬譜吻合良好。圖3(a)、(b)給出了GLF為10.0,而群高因子GFH分別為0.39,0.89情況下,所模擬出來的波群波列。圖4(a)、(b)給出了GFH為0.71,而群長因子GLF分別為9.82,20.8情況下,所模擬出來的波群波列。從圖上可以看出,按照上述方法可以模擬出所要求的不同群高和不同群長特征的波列。群高主要表示波群內外波振幅的差別,其統計特征量可反映波包的起伏幅度,從圖3可以看出,群高參數GFH越小,波浪包絡越平穩,GFH=0時,即為規則波情況。同時,群長主要表示群內大波個數,其統計特征量可反映波浪包絡的變化快慢,或者說是連續大波的長度,從圖4可以看出,通過給定不同的群長參數GLF,群長大的波浪包絡線內大波更多。

圖2 波浪模擬譜和靶譜的比較Fig.2 Comparison between the simulation spectrum and the target spectrum

圖3 (b)GFH=0.89,GLF=10.0波群波面Fig.3(b)Wave grouping surface GFH=0.89,GLF=10.0

4(a)GFH=0.71,GLF=9.82波群波面Fig.4(a)Wave grouping surface GFH=0.71,GLF=9.82

圖4 (b)GFH=0.71,GLF=20.80波群波面Fig.4(b)Wave grouping surface GFH=0.71,GLF=20.80

1.4 纜繩對系泊船的作用力及力矩

根據系泊系統的位移—張力關系即導纜孔至系纜柱的空間距離的變化,依據Wilson公式求得纜力和力矩:

其中:CP為纜繩彈性系數,d為纜繩直徑,ε為纜繩變形量。

對于纜繩系統,t時刻纜繩系在碼頭系船柱上的端點坐標O(x,y,z),纜繩系在船體上的端點坐標),纜繩原始長度:

t+1時刻隨著船體運動到下一個位置,纜繩系船柱上的端點坐標仍為O(x,y,z),在船體上的端點運動到下一位置,坐標變為O2),纜繩長度:

則t+1時刻纜繩變形長度:

3.環境營造:讓學員置身于一個良好的具身學習情境。具身認知理論認為,認知不僅源自大腦和身體本身,還包含與所處特定情境中身體與該情境互動的結果。萊夫(lave)也明確地指出:“情境意味著,在特殊性和普遍性的許多層面上,一個特定的社會實踐和活動系統中社會過程的其他方面具有多重的交互聯系。”[10]身體一旦離開了情境,就無所謂交互,更不可能產生新的認知構建。所以,培訓是人腦、身體與特定情境交互的知識體驗過程,學員的身體只有沉浸在良好的情境中,調動大腦的參與,才能更好地獲得知識和能力的建構、內化和遷移。而培訓中的情境應該包括:心理情境、物理情境和虛擬情境。

△s=s-s1

ε=△s/s

1.5 護舷對船舶的撞擊力和力矩

碼頭上通常設有護舷,用以吸收船舶靠泊時的能量,保護船舶和碼頭結構。設第j個護舷的坐標為,護舷對船舶的作用力和關于重心的力矩為

式中:Ej為護舷剛度;為第j個護舷對系泊船的力矩。

2 模型驗證及數值計算結果與分析

2.1 模型驗證

為了驗證數值計算模型,取一漂浮半球,半徑為a,自由漂浮于H=3a的水中,球的重力與其浮力相等,重心與球心重合。圖5(a)給出了規則波浪作用下漂浮半球水平和垂向一階力的計算結果,并與滕斌[18]的計算結果作了比較。圖中f1,f3分別為波浪對漂浮半球一階水平力和垂向力。從圖中可以看出,本文計算漂浮半球的波浪力結果與文獻結果吻合很好。

圖5 (a)漂浮半球一階作用力Fig.5(a)The first order wave force of a floating hemisphere

圖5 (b)漂浮半球縱移運動時歷Fig.5(b)The surge of a floating hemisphere

圖5 (c)漂浮半球升沉運動時歷Fig.5(c)The heave of a floating hemisphere

圖5 (d)漂浮半球縱搖運動時歷Fig.5(d)The pitch of a floating hemisphere

2.2 數值模型

數值算例中船舶、碼頭系纜布置及環境荷載方向如圖6所示。布纜方式采用2:2:2:2:2:2方式,即首纜2根(圖6中1號纜),首橫纜2根(圖6中2號纜),首倒2根(圖6中3號纜),尾倒2根(圖中6中4號纜),艉橫纜2根(圖6中5號纜),尾纜2根(圖6中6號纜)。碼頭頂高程9.5 m,碼頭底高程-21.0 m,計高水位1.84 m,入射波浪為橫浪。

圖6 碼頭纜繩布置圖Fig.6 Mooring arrangement of the model experiment

為了考察系泊船在不同群性波浪作用下的運動響應,按照1.3節所述數值模擬了波浪要素為H1/3=1.2m,T1/3=9.3s,不同群高和群長的隨機波浪。文中列出了部分波列作用下的數值模擬結果,其中,波包譜通過對模擬波面希爾特變換后計算得到。

表1 船舶主要尺寸及系纜尺寸Tab.1 Dimensions of the 80kdwt and mooring

圖7 (a)GFH=0.93與GFH=0.39,頻譜對比圖(GLF=13.5)Fig.7(a)Wave spectrum(GFH=0.93,GFH=0.39)

圖7 (b)GFH=0.93與GFH=0.39波包譜對比圖(GLF=13.5)Fig.7(b)Wave envelope spectrum(GFH=0.93,GFH=0.39)

圖8 不同群高波列作用橫移的計算頻譜對比圖(GFH=0.93與GFH=0.39,GLF=13.5)Fig.8 Comparison of the wave spectrum of sway with variety group height factor(GFH=0.93與GFH=0.39,GLF=13.5)

圖9 (b)GLF=20.88與GLF=9.82波包譜對比圖 (GFH=0.80)Fig.9(b)Wave envelope spectrum(GLF=20.88, GLF=9.82)

圖10 不同群長波列作用下橫移的計算頻譜對比圖 (GLF=9.82與GLF=20.88,GFH=0.80)Fig.10 Comparison of the wave spectrum of sway with variety group length factor(GLF=9.82與GLF=20.88,GFH=0.80)

圖7 給出了GLF為13.5時,GFH分別為0.93和0.39的波群波面的計算頻譜、計算波包譜。從圖上可以明顯看出,波面的頻譜基本接近,而對于群高大的波群波面,其波包譜能量相對較大。為了進一步探討不同群高波面作用下對系泊船橫移和升沉運動的頻域特征的影響,圖8給出了兩組波群波面作用下系泊船的橫移運動量的計算頻譜對比圖。從圖中還可以看出,在峰頻基本一致條件下,群高參數大的波列作用下的橫移計算頻譜峰值大。

圖9給出了GFH為0.8時,GLF分別為9.82和20.88的波群波面的計算頻譜,計算波包譜。可以看出,對于群長大的波群波面,其波包譜能量相對較大。同樣地,圖10為GFH為0.8時,GLF分別為9.82和20.88的波群波面作用下系泊船的橫移運動量的計算頻譜。從圖10中可以看出群長小波列作用下的橫移計算頻譜譜峰值大,峰頻基本一致。

2.3 數值計算結果分析

圖11至圖14分別給出了系泊船的纜繩平均張力和護舷的撞擊力在不同群高及不同群長波列作用下的變化情況。從圖11和圖12可知,總體上來講,系泊船的平均纜繩張力隨著群高因子的增大而增大,其中橫纜相比增幅更明顯,同時,隨著群長因子的增大而略微減小,且減小幅度不大,即群長因子對平均纜繩張力影響較不明顯。從圖13和14可見,撞擊力變化很小,即群高因子和群長因子對護舷撞擊力影響較小。

圖11 纜繩平均張力與GFH的關系Fig.11 Relationship between the mean line force and GFH

圖12 纜繩平均張力與GLF的關系Fig.12 Relationship between mean line force and GLF

圖13 護舷撞擊力與GFH的關系Fig.13 Relationship between impact force of fender and GFH

圖14 護舷撞擊力與GLF的關系Fig.14 Relationship between impact force of fender and GLF

圖15 至圖18分別給出了系泊船的橫移和升沉運動量在不同群高即不同群長波列作用下的變化情況。圖中系泊船的橫移和升沉運動量計算頻譜的零階距以無因次形式給出,其中M0ξ2,M0ξ3代表橫移及升沉計算頻譜的零階距,M0η波面的頻譜零階距。從圖15和圖16可以看出,總體上,系泊船橫移運動量隨著群高因子增大而增大,隨著群長因子的增大而略微減小。與此同時,從圖17和圖18可以看出,系泊船升沉運動量隨著群高因子增大而略微增大,隨著群長因子的增大而略微減小。相比橫移運動量而言,系泊船的升沉運動量受波浪群性的影響小。因為橫移運動量的慣性相比升沉運動量的慣性小,在同等外力作用下運動量更大。

圖15 橫移與GFH的關系Fig.15 Relationship between sway and GFH

圖16 橫移與GLF的關系Fig.16 Relationship between sway and GLF

圖17 升沉與GFH的關系Fig.17 Relationship between heave and GFH

圖18 升沉與GLF的關系Fig.18 Relationship between heave and GLF

圖19 有效波高與GFH的關系Fig.19 Relationship between GFH and significant wave heights

圖20 最大波高與GFH的關系Fig.20 Relationship between GFH and max heights

圖21 波群對橫移運動量的作用圖Fig.21 Effect of grouping wave on sway

圖22 波群對升沉運動量的作用Fig.22 Effect of grouping wave on heave

為了探討不同群性的波浪引起系泊船運動響應變化的原因,下面我們就模擬波面的群性特征與其統計特征波高的關系作了分析。圖19和圖20給出了不同群高參數的波面其有效波高和最大波高值。這里按上跨零點法計算統計波高。從圖上可以看出,隨波面群高的增大,有效波高和最大波高都逐漸增大,且最大波高值受群高影響更大些。可見,波浪統計特征與群高有關,這是由于按相同波浪要素模擬得到的不同群高的波面,其波浪形狀變化很明顯,群高因子較小時,波面趨于平緩,故有效波高和最大波高相對較小,群高較大時,有效波高和最大波高相對較大。波浪群高特征的不同引起的波面統計特征的差異是導致系泊船的運動響應受群高參數變化影響的原因所在。

為了進一步說明考慮群性的波浪對系泊船運動響應的影響,圖21、22給出波浪波群與未考慮波浪群性的差異。本文對兩種波浪作用下不同波高的無量綱的運動響應進行了對比,波浪群性對系泊船的運動響應有較大的影響,相同波要素下比未考慮波浪群性波浪對系泊船運動量增幅明顯。

3 結 論

本文通過模擬具有相同的有效波高、有效周期,但不同群高和群長的隨機波浪,研究了系泊船在不同群性的波浪作用下的運動響應。數值計算結果表明,波浪群性對系泊船的運動影響較大,其中群高因子GFH較大的隨機波列作用下系泊船的運動響應也較大,群長因子GLF較大的隨機波列作用時系泊船的運動響應反而略有減小。此外,纜繩張力隨著波列的群高因子的增大而增大,隨著群長因子的增大而略微減小,且減小幅度不大,即群長因子對平均纜繩張力影響較不明顯,同時,波浪的群性對護舷撞擊力影響較小。

[1]López M,Iglesias G.Long wave effects on a vessel at berth[J].Applied Ocean Research,2014,47:63-72.

[2]Sakakibara S,Kubo M.Characteristics of low-frequency motions of ships moored inside ports and harbors on the basis of field observations[J].Marine Structures,2008,21:196-223.

[3]Sasa K,Incecik A.Numerical simulation of anchored ship motions due to wave and wind forces for enhanced safety in offshore harbor refuge[J].Ocean Engineering,2012,44:68-78.

[4]Cummins W E.The impulse response function and ship motions[J].Schiffstechnik,1962,47:101-109.

[5]Bingham H B.A hybrid Boussinesq-panel method for predicting the motion of a moored ship[J].Coastal.Engineering, 2000,40(1):21-38.

[6]Kwak M,Moon Y,Pyun C.Computer simulation of moored ship motions induced by harbour resonance in Pohang new harbor[C].In:Lynett P,Smith J M,Proceedings of the 33rd International Conference on Coastal Engineering-ICCE. 2012:2156-1028.

[7]Nolte K G,Hsu F H.Statistics of ocean wave groups[C].Fourth Annual Offshore Technology Conf paper,OTC1688,1972: 637-644.

[8]Johnson R R,Mansard E P D,Ploeg J.Effects of wave grouping on breakwater stability[C]//Proceedings of Sixteenth Coastal Engineering Conference.Hamburg,Germany,1978,3:2228-2243.

[9]Murray J J,Wilkie B P,Muggeeridge D B.An experimental and analytical study fo the effects of wave grouping on the slowly varying drift oscillation of a floating rectangular barge[J].Ocean Engineering,1987,14(4):255-274.

[10]Balaji R,Sannasiraj S A,Sundar V.Detection of wave groups from the motion behavior of a discus buoy[J].Journal of Hydro-environoment Research,2008:195-205.

[11]Xu D L,Hou W,Zhao M,Wu J.The statistical simulation of wave groups[J].Applied Ocean Research,1993,15:217-226.

[12]Van Oortmerssen G.The motions of a moored ship in waves[C].Proceedings of the 2nd Int Conf of Behavior of Offshore Structures,1979:339-356.

[13]鄒志利.港口內靠碼頭系泊船運動的計算[[J].海洋工程,1995,13(3):25-35. Zou Z L.A mathematical model simulation of a ship moored against a quay in harbor[J].Ocean Engineering,1995,13(3): 25-35.

[14]俞聿修.隨機波浪及其工程應用[M].大連:大連理工大學出版社,2003.

[15]俞聿修,桂滿海.海浪的群性及其主要特征參數[J].海洋工程,1998,16(3):9-21.Yu Yuxiu,Gui Manhai.Numerical simulation and physical simulation of sea wave groups[J].Journal of Dalian University of Technology,1998,38(1):86-91.

[16]俞聿修,桂滿海.波群的數值模擬和物理模擬[J].大連理工大學學報,1998,38(1):86-91. Yu Yuxiu,Gui Manhai.Numerical simulation and physical simulation of sea wave groups[J].Journal of Dalian University of Technology,1998,38(1):86-91.

[17]劉 思,柳淑學,俞聿修.改進的海浪波包譜[J].水道港口,2010,31(4):229-235. Liu Si,Liu Shuxue,Yu Yuxiu.An improved wave envelope spectrum[J].Journal of Waterway and Harbor,2010,31(4): 229-235.

[18]滕 斌.波浪對三維浮體的二階作用[J].水動力學研究與進展,1995,10(3):316-327. Teng Bin.Second order wave action on floating bodies[J].Journal of Hydrodynamics,1995,10(3):316-327.

A numerical study on the hydrodynamic response of moored ship under the action of wave grouping

MA Xiao-jian1,2,SUN Zhao-chen2,LIU Si2,ZHANG Zhi-ming3,YANG Guo-ping3,ZHOU Feng3
(1.School of Naval Architecture and Ocean Engineering,Jiangsu University of Science and Technology,Zhenjiang 212003, China;2.State Key Laboratory of Coastal and Offshore Engineering,Dalian University of Technology,Dalian 116024, China;3.China Communications Planning&Design Institute for Water Transportation,Beijing 100007,China)

Waves with the same wave parameters,such as significant wave height and period,but with different wave groupiness are simulated,the responses of moored ship under the action of wave groups with different group height factor GFH and group length factor GLF are discussed.The numerical results show that the sway,heave motions and mooring line force increase with an increase in GFH,and the influence of GLF on the motions of a moored ship and mooring line force is nearly neglectable.

wave grouping;moored ship;GFH;GLF

P75

A

10.3969/j.issn.1007-7294.2015.07.006

1007-7294(2015)07-0797-13

2014-12-23

國家自然科學基金(NSF51309124);江蘇科技大學船舶與海洋工程學院青年基金項目;國家創新研究群體科學基金(5092100)資助

馬小劍(1982-),男,博士,講師,E-mail:immaxiaojian@126.com;

孫昭晨(1960-),男,教授,E-mail:sunzc@dlut.edu.cn。

主站蜘蛛池模板: 国产又爽又黄无遮挡免费观看| 日韩第八页| 99热这里只有精品在线播放| 欧美精品H在线播放| 亚洲综合极品香蕉久久网| 丁香婷婷激情网| 国产永久无码观看在线| 成年人午夜免费视频| 成年午夜精品久久精品| 亚洲中文制服丝袜欧美精品| 国产成本人片免费a∨短片| 这里只有精品免费视频| 中国黄色一级视频| 欧美精品影院| 免费Aⅴ片在线观看蜜芽Tⅴ| 亚洲精品免费网站| 青青草国产精品久久久久| 区国产精品搜索视频| av在线无码浏览| 国模沟沟一区二区三区| 欧美.成人.综合在线| 久久这里只有精品2| 成人在线不卡| 国产av一码二码三码无码| 国产第二十一页| 18禁色诱爆乳网站| 国产一区二区三区精品久久呦| 91啪在线| 精品久久久无码专区中文字幕| 91啪在线| 欧美日韩中文国产| 国产一二三区在线| 99无码中文字幕视频| 人妻免费无码不卡视频| 国产亚洲第一页| 一级毛片网| 国产精品9| 国产又黄又硬又粗| 欧美专区日韩专区| 国产午夜小视频| 成人免费午夜视频| 美女黄网十八禁免费看| 青草视频免费在线观看| 一本色道久久88亚洲综合| 欲色天天综合网| 亚洲人成人无码www| 国产成人福利在线视老湿机| 一区二区三区四区在线| 亚洲美女一区二区三区| 亚洲91在线精品| 成人在线综合| 中文字幕无码av专区久久| 一区二区日韩国产精久久| 免费观看男人免费桶女人视频| 国产一区二区三区日韩精品| 日韩免费毛片视频| 久久公开视频| 超薄丝袜足j国产在线视频| 国产成人精品高清不卡在线 | 欧美国产日本高清不卡| 亚洲天堂.com| 亚洲全网成人资源在线观看| 亚洲视频一区| 伊人色天堂| 国产老女人精品免费视频| 91视频99| 欧美中文字幕在线二区| 91亚洲免费| 九九热视频精品在线| 国产肉感大码AV无码| 中文字幕自拍偷拍| 亚洲激情99| 91最新精品视频发布页| 亚洲精品片911| 国产新AV天堂| 国产免费一级精品视频| 天堂av综合网| vvvv98国产成人综合青青| 9啪在线视频| a网站在线观看| 精品无码一区二区三区电影| 91福利国产成人精品导航|