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

海底多段懸跨管道渦激振動特性分析

2017-08-27 05:36:11徐萬海謝武德高喜峰馬燁璇
船舶力學 2017年8期
關鍵詞:模態振動

徐萬海,謝武德,高喜峰,馬燁璇

(天津大學水利工程仿真與安全國家重點實驗室,天津300072)

海底多段懸跨管道渦激振動特性分析

徐萬海,謝武德,高喜峰,馬燁璇

(天津大學水利工程仿真與安全國家重點實驗室,天津300072)

深海海底多跨管道相比于單跨,其渦激振動行為卻更為復雜。文章考慮了兩端跨肩和中間支撐處的管—土作用邊界條件,提出了多跨管道渦激振動預報模型,標定了模型參數,并重點分析了跨肩邊界條件、中間支撐條件和懸跨長度對管道渦激振動特性的影響。結果表明:較大的兩端邊界處扭轉彈簧彈性系數、較大的中間支撐處扭轉和拉伸彈簧彈性系數和較短的懸跨長度,會限制管道結構高階模態的激發,同時還發現了振動能量會在相鄰管跨之間傳遞的現象。

海底管道;流—固耦合;渦激振動;多跨

0 引言

海底管道是深海油氣開采輸運不可或缺的關鍵設備。由于深海海底地形復雜,地勢大量起伏,加之海流沖刷作用,直接鋪設于海底的管道頻繁出現懸跨,由此引發的懸跨管段渦激振動是導致管道結構疲勞損傷的重要因素[1]。

長期以來,人們更多地關注單跨管道渦激振動問題[2-4],缺乏對多跨管道渦激振動的深入研究。與單跨管道相比,多跨管道的邊界條件更加復雜、振動常呈現多模態、危險性更高等特點。近期,挪威的MARINTEK機構研究了多跨單模態和多跨多模態的渦激振動特性,采用特殊彈簧裝置近似兩端跨肩處的土體邊界對管道作用[5-6],部分實驗成果反映在最新推薦的海底管道設計規范DNV-RP-F105[7]上。李明明[8]數值模擬了跨肩長度、土體性質等因素對懸跨管段振動固有頻率的影響,并對深海多跨管道VIV疲勞進行了評估。Soni[9]通過特殊的彈簧裝置實驗模擬相鄰管跨之間的土體支撐,觀測了多跨管道鄰近跨間的相互作用。Koushan[10]在深入研究多段懸跨海底管道渦激振動后指出:準確地分析自由管跨段VIV,必須合理簡化跨肩處的土體邊界和多跨管道的相鄰跨間影響機制等。

現階段,學術界和工程界對海底多跨懸空管道渦激振動問題的認識水平還無法滿足深海工程實際需求。本文首次提出多跨管道渦激振動的預報模型,并討論多跨管道兩端跨肩邊界、中間支撐條件和相鄰管跨段長度等因素對渦激振動特性的影響機制。

圖1 多跨管道示意圖Fig.1 Multi-span pipelines

1 數學模型

1.1 管道結構方程

本文研究垂直來流方向上的橫流向渦激振動,忽略管道的幾何非線性和水動力對跨肩處管段的影響,如圖1所示的多段懸跨管道能夠簡化為歐拉-伯努利梁,管道的左端取為直角坐標系原點,結構的控制方程如下:

式中:m=ms+ma,ms為單位長度結構質量,ma(=CaπρD2/4)為附加質量,Ca是附加質量系數,D為管道的外徑,ρ為海水密度,y為橫流向的結構位移,t為時間,C=Cs+Cf為管道的阻尼(包括結構阻尼Cs和水動力阻尼Cf),T為管道軸向力,x為管道軸向的空間坐標,EI為結構的彎曲剛度,FL(=CLρDU2/2)是漩渦脫落引發的渦升力,CL為瞬時的升力系數,U為海流速度。

1.2 尾流振子方程

海流流經懸跨管段時,會誘發交替的漩渦脫落,工程上常采用尾流振子模型分析尾流與結構的相互作用,運用Van der pol方程描述漩渦脫落的尾跡特性。Facchinetti[11]提出了加速度耦合尾流振子模型,可很好地描述海底管道渦激振動這種流-固耦合問題,模型方程如下:

式中:變量q表示瞬時的升力系數CL與固定圓柱升力系數CL0之比,q=2CL/CL0,ωs(=2πStU/D)為漩渦脫落頻率,其中St是Strouhal常數。ε和A為尾流振子的模型參數,根據漩渦脫落的情況而選取不同的數值,它反映了尾流振子模型和管道渦激振動的耦合關系。徐萬海等[12]對尾流振子的模型參數進行了修定,彌補了Facchinetti模型僅能定性不能定量預報渦激振動的不足。當約化速度Vr(=2πU/(ωnD)),ωn為管道結構的圓頻率)處于不同范圍時,參數A的取值為[12]:

根據管道的最大響應幅值與實驗值相等的原則,確定模型參數ε,如下:

式中:SG為約化質量阻尼參數:

γ為遲滯因子,表達式為:

ζ為結構阻尼因子,CD為阻力系數,綜合利用(3)-(6)式,最終確定ε和A值。本文采用上述尾流振子模型描述漩渦脫落與管道結構之間的流-固耦合作用。

1.3 兩端跨肩和中間支撐處的邊界條件

土壤的性質、管土相互作用長度、管線掩埋深度和管線掩埋寬度對海底懸跨管道的渦激振動都有很大的影響。根據海底管線的實際情況,DNV-RP-F105[7]將端部跨肩處管土相互作用土壤分成了多個等級。工程界為了簡化分析過程,常將懸跨管道兩端約束簡化為簡支邊界或者固支邊界。實際上,兩端跨肩處的約束條件一般介于簡支狀態和固支狀態之間。因此,跨肩邊界可等效為彈性支承,以扭轉彈簧和拉伸彈簧約束模擬工程實際中的管土相互作用,假設海床水平,邊界條件方程為[13]:

其中:kt1和ks1分別為懸跨管道左端扭轉和拉伸彈簧彈性系數,kt2和ks2為右端扭轉和拉伸彈簧彈性系數,L=L1+L2,為懸跨管道的整體長度。

中間支撐部分近似為彈性支撐,ktm和ksm為支撐處的扭轉和拉伸彈簧彈性系數,并且在中間支撐處約束滿足懸跨L1和L2在連接處的位移、轉角、剪力和彎矩相等[14],如下所示:

其中:下角標“-”表示中間支撐處的左側,“+”表示中間支撐處的右側。

2 數值求解方法

本文以兩跨管道為例,設管道的振型函數為:

其中:φ1(x)是管道第一跨的振型函數,φ2(x)是管道第二跨的振型函數,即:

cij(i=1,2;j=1,2…4)是待求系數。

根據懸跨管道兩端邊界條件,以及中間支撐條件,可獲得管道的振型函數。運用模態展開法,分析懸跨管道的渦激振動特性,令:

將(18)式代入(1)和(2)式,根據振型正交性,對方程進行化簡:

3 模型正確性驗證

采用四階龍格-庫塔法對方程組(19)和(20)進行數值求解。海底懸跨管道渦激振動一般較難出現高階模態,因此僅取前四階模態(n=4)分析,即可得到合理結果。為了驗證多跨管道渦激振動尾流振子模型的合理性,以及數值算法的正確性,首先以MARINTEK的Koushan實驗(如圖2所示)[10]為算例,進行分析驗證。實驗中兩段懸跨管道長度相等,為4.126 m,管道外徑為0.035 11 m,模型的其它參數:軸向張力為50 N,單位長度質量為1.307 kg/m,彎曲剛度EI=203 Nm2。圖3對比了管道結構的最大響應幅值計算結果與實驗結果,DNV規范結果,以及渦激振動計算軟件VIVANA的模擬結果,通過分析可以發現:較低的流速條件下,多跨管道一階模態被激發,隨著流速增大,二階、三階模態相繼被激發,本文的預報模型計算結果比DNV規范略小,比VIVANA軟件結果稍大。計算過程中第四階模態一直未被激發,也充分證明了取前四階模態進行分析的合理性。模型預報結果在較低流速條件下與實驗結果十分接近,但較高流速時,產生了一定差異,主要原因是MARINTEK的管道模型軸向力是在靜水條件下施加的,但實驗過程中,隨著外部來流速度的增大,管道模型產生了一定的撓曲變形,直接導致軸向力大幅增加,而本文模型分析的軸向力取值相對偏小,導致了預報結果大于實驗觀測結果,但是整體的變化趨勢比較接近。所以,從整體分析結果判斷,該渦激振動預報模型的正確性和適用性是可以保證的。

圖2 MARINTEK實驗裝置示意圖Fig.2 Sketch of the experimental apparatus in MARINTEK

4 多跨管道渦激振動特性分析

本部分重點分析跨肩邊界條件、中間支撐條件和不同管跨的懸跨長度等對多跨管道渦激振動特性的影響,選擇實際工程中最常見的兩跨管道作為分析對象,基本的計算參數如表1所示。外部來流速度根據實際海流速度,選為0.4-3.5 m/s之間,并假定兩端跨肩處具有對稱的邊界條件(kt1=kt2=kt;ks1= ks2=ks)。

圖3 最大響應幅值結果對比結果Fig.3 Comparison of the results on maximum response amplitude

4.1 端部邊界條件影響

首先研究跨肩處拉伸彈簧系數足夠大(ks= 1010),扭轉彈簧系數變化時(kt=0、102、104、106、108和1010),討論其對海底兩跨管線渦激振動的最大響應振幅和頻率的影響,中間支撐假定為簡支(ktm= 1010,krm=0),令第一跨的懸跨長度為80 m,第二跨的懸跨長度為60 m。改變兩端扭轉彈簧彈性系數,模擬兩端邊界從簡支變化到固支的全過程。

圖4和圖5分別給出了扭轉彈簧彈性系數變化時,不同懸跨段管道結構的響應幅值和響應頻率。從圖4(a)中可以看出隨著扭轉彈簧彈性系數增加,最大響應幅值曲線的峰值幾乎沒有發生改變,僅是整體向右移動,扭轉彈簧彈性系數在0-106區間范圍時,扭轉彈簧對懸跨段L1的渦激振動影響很小;彈性系數介于106-1010時,彈性系數越大,懸跨段L1的響應幅值就越小,對響應幅值曲線的峰值影響不明顯。圖4(b)中有關L2管跨段最大響應幅值結果與L1段相似,其結論也相似。對比圖4(a)和圖4(b)能夠發現,L1管跨段結構響應幅值獲得極大值時,L2管跨段會出現極小值,反之亦然,這充分地證明了渦激振動能量存在在相鄰管跨之間的實時傳遞的現象。從圖5(a)中能夠得到扭轉彈簧系數在0-106時,其對懸跨段L1的渦激振動響應頻率幾乎沒有影響;扭轉彈簧系數在106-1010時,扭轉彈簧系數越大,懸跨段L1的渦激振動響應頻率就越大,高階模態也越難被激發。圖5(b)中的L2有相似結果。

表1 計算模型參數[10]Tab.1 Base set of system parameters[10]

圖4 扭轉彈簧彈性系數變化時,不同管跨的響應幅值Fig.4 Response amplitudes of different spans with torsional spring stiffness

圖5 扭轉彈簧彈性系數變化時,不同管跨的響應頻率Fig.5 Response frequencies of different spans with torsional spring stiffness

4.2 中間支撐條件影響

海底多跨管道中間一般是由土壤進行支撐的,扭轉彈簧和拉伸彈簧常被用來模擬中間支撐處的管土相互作用,在此,重點分析多跨管道中間支撐處扭轉彈簧和拉伸彈簧對兩跨管線渦激振動的影響。為了簡化起見,將兩端設定為簡支邊界。仍假設第一跨的懸跨長度為80 m,第二跨的懸跨長度為60 m。

4.2.1 中間扭轉彈簧

固定中間支撐處的拉伸彈簧系數(設為ksm=1010),變化中間扭轉彈簧系數(ktm=0、102、104、106、108和1010)。計算渦激振動響應的最大響應振幅和響應頻率,如圖6和圖7所示。

圖6 中間支撐處扭轉彈簧彈性系數變化時,不同管跨的響應幅值Fig.6 Response amplitudes of different spans with torsional spring stiffness at the middle support

圖7 中間支撐處扭轉彈簧彈性系數變化時,不同管跨的響應頻率Fig.7 Response frequencies of different spans with torsional spring stiffness at the middle support

分析圖6(a)最大響應幅值,扭轉彈簧系數處于0-106范圍內時,最大振幅響應曲線幾乎重合。扭轉彈簧系數介于106-108之間時,隨著中間扭轉彈簧系數的增加,L1的響應幅值減小,而且整個振幅曲線向右移動。這說明:扭轉彈簧系數在0-106時,中間扭轉彈簧系數對懸跨段L1渦激振動最大響應幅值幾乎沒有影響;扭轉彈簧系數在106-1010時,中間扭轉彈簧系數越大,L1渦激振動響應幅值就越小,高階模態激發就越困難。對于L2跨管道,圖6(b)的結果與圖6(a)大致相同。

圖7(a)給出了中間不同扭轉彈簧系數,不同來流速度下,L1管跨段的響應頻率。從圖中可以得到:扭轉彈簧系數在0-106時,懸跨段L1的響應頻率曲線重合;扭轉彈簧系數在106-1010時,隨著中間扭轉彈簧系數的增加,懸跨段L1的響應頻率略有增加。隨著來流速度的增加,較小的扭轉彈簧系數會提前激發高階振型。L2跨管道具有與L1跨管相近的特性,如圖7(b)所示。

4.2.2 中間拉伸彈簧

進一步分析中間拉伸彈簧對懸跨管線渦激振動響應的影響。設定兩端為簡支邊界,固定中間支撐處的扭轉彈簧系數(設為ktm=0),變化中間拉伸彈簧系數(ksm=0、102、104、106、108和1010)。計算兩跨管線渦激振動響應的最大振幅和頻率,如圖8和圖9所示。

圖8 中間支撐處線彈性彈簧彈性系數變化時,不同管跨的響應幅值Fig.8 Response amplitudes of different spans with stretching spring stiffness at the middle support

圖9 中間支撐處線彈性彈簧彈性系數變化時,不同管跨的響應頻率Fig.9 Response frequencies of different spans with stretching spring stiffness at the middle support

根據圖8(a)中懸跨段L1最大響應幅值可發現,隨著中間支撐處拉伸彈簧彈性系數增加,最大響應幅值曲線整體右移。同時明顯區分為三個振動狀態,中間支撐拉伸彈簧的彈簧系數介于0-104、108-1010之間時,響應幅值近似相同,而在104-108的過渡區間段,曲線產生了明顯的分離。圖8(b)中的懸跨段L2具有相似的結論。

圖9中兩組曲線的總體走勢基本相同,拉伸彈簧系數在0-104時,響應頻率曲線幾乎重合,拉伸彈簧系數在108-1010時,響應頻率曲線也十分接近,但是不同于過渡區間的情形。

4.3 懸跨長度的影響

懸跨長度是影響海底管道渦激振動特性的重要因素,尤其對于多跨管道。相鄰管跨之間、懸跨跨長對管道渦激振動影響機制尚不清晰。在此,懸跨段L1長度設定為80 m,變化懸跨段L2的長度(取為40、60、80和100m),探討研究懸跨長度對兩跨管道渦激振動的影響。繪制渦激振動響應最大響應幅值和頻率圖像,如圖10和圖11所示。

分析圖10(a)中的懸跨段L1最大響應幅值可知,隨著懸跨段L2長度的增加,響應幅值略有增加,曲線波谷和波峰都有向左移動的趨勢。這表明:懸跨段L2越長,懸跨段L1的渦激振動響應幅值就越大,而且越容易激發高階模態。圖10(b)中懸跨段L2最大響應幅值與圖10(a)中懸跨段L1最大響應幅值有明顯的不同。隨著懸跨段L2長度的增加,最大響應幅值曲線向左移動,峰值和谷值逐漸增加。這說明懸跨段L2越長,其渦激振動響應幅值就越大,而且高階模態容易被激發。

圖10 整體懸跨長度變化時,不同管跨的響應幅值Fig.10 Response amplitudes of different spans with the whole spanning length of pipeline

圖11 整體懸跨長度變化時,不同管跨的響應頻率Fig.11 Response frequencies of different spans with the whole spanning length of pipeline

圖11(a)懸跨段L1響應頻率可以發現相同外界流速的情況下,隨著懸跨段L2的增加,懸跨段L1的響應頻率逐漸降低;模態競爭區域向左移動。這表明:懸跨段L2長度越長,懸跨段L1的響應頻率就越大,越容易進入模態競爭區域,越容易激發高階模態。圖11(b)中懸跨段L2響應頻率和圖11(a)懸跨段L1響應頻率具有相同的現象和相似的結論。

5 結論

針對海底多跨懸空管道橫流向渦激振動問題,本文提出了適用于多跨管道渦激振動預報的尾流振子模型,通過將模型預報結果與DNV規范、VIVANA數值模擬和實驗結果的對比,驗證了該模型的適用性和精確性。運用該模型,分析了多跨管道的跨肩邊界條件、中間支撐條件以及懸跨長度對渦激振動的影響機制,通過對比分析最大響應幅值和響應頻率,可發現:

(1)兩端邊界處的扭轉彈簧彈性系數較低時,對多跨管道懸跨段的渦激振動幾乎沒有影響;但隨著彈性系數的增加,影響逐漸凸顯,懸跨段的最大響應幅值會略微變小,響應頻率會增加,高階模態被激發的難度變大。同時證實了渦激振動的能量會在相鄰管跨之間傳遞。

(2)中間支撐處的扭轉彈簧彈性系數處于較低區間范圍時,對多跨管道渦激振動的最大響應幅值和響應頻率幾乎都沒有影響;較高的扭轉彈簧彈性系數會讓中間支撐處的約束增強,導致高階模態不容易出現;中間支撐處的拉伸彈簧也可以得到類似的結論。

(3)隨著懸跨長度增加,多跨管道最大響應幅值會增加,響應頻率亦會增加,容易激發高階模態,也容易進入模態競爭區域。

[1]Bai Yong,Bai Qiang.Subsea pipelines and risers[M].Elsevier Science Ltd,2005.

[2]吳鈺驊,金偉良,龔順風,等.海底管道流固耦合振動數值模擬[J].浙江大學學報(工學版),2009,43(4):782-788. Wu Yuhua,Jin Weiliang,Gong Shunfeng,et al.Submarine pipeline vibration considering fluid-structure interaction[J]. Journal of Zhejiang University(Engineering Science),2009,43(4):782-788.(in Chinese)

[3]Srinil N,Zanganeh H.Modelling of coupled cross-flow/in-line vortex-induced vibrations using double Duffing and van der Pol oscillators[J].Ocean Engineering,2012,53:83-97.

[4]李小超,王永學,王國玉,等.固定海床上海底管跨順向和橫向耦合振動實驗研究[J].船舶力學,2013,17(5):547-558. Li Xiaochao,Wang Yongxue,Wang Guoyu,et al.Experimental investigation of in-line and cross-flow VIV of free spans near a plane boundary[J].Journal of Ship Mechanics,2013,17(5):547-558.(in Chinese)

[5]Nielsen F G,Soreide H T,Kvarme S O.VIV response of long free spanning pipelines[C]//The 21st International Conference on Offshore Mechanics and Arctic Engineering,OMAE.Oslo,Norway,2002.

[6]Soreide T H,Ilstad T,Paulsen G,Nielsen F G.Dynamic response of pipeline in long free spans or multi-spans[C]//The Sixth European Conference on Structural Dynamics,EURODYN.Paris,France,2005.

[7]DNV-RP-F105.Free spanning pipelines[M].HΦvik.Det Norske Veritas,2006.

[8]李明明.深水多跨管道VIV疲勞評估[D].哈爾濱:哈爾濱工程大學,2011. Li Mingming.VIV fatigue assessment of deepwater multi-span pipeline[D].Harbin:Harbin Engineering University,2011. (in Chinese)

[9]Soni P K.Hydrodynamic coefficients for vortex-induced vibrations of flexible beams[D].Norwegian:Norwegian University of Science and Technology(NTNU),2008.

[10]Koushan K.Vortex induced vibrations of free span pipelines[D].Norwegian:Norwegian University of Science and Technology(NTNU),2009.

[11]Facchinetti M L,Langre E,de Biolley F.Coupling of structure and wake oscillators in vortex-induced vibrations[J].Journal of Fluids and Structures,2004,19(2):123-140.

[12]Xü Wanhai,Wu Yingxiang,Zeng Xiaohui,Zhong Xingfu,Yü Jianxing.A new wake oscillator model for predicting vortex induced vibration of a circular cylinder[J].Journal of Hydrodynamics,Ser.B,2010,22(3):381-386.

[13]Bakhtiary AY,Ghaheri A,Valipour R.Analysis of offshore pipeline allowable free span length[J].International Journal of Civil Engineering,2007,5(1):84-91.

[14]Motaghian S E,Mofid M,Alanjari P.Exact solution to free vibration of beams partially supported by an elastic foundation [J].ScientiaIranica,2011,18(4):861-866.

Study on Vortex-Induced Vibrations(VIV)of multi-spans pipelines

Xü Wan-hai,XIE Wu-de,GAO Xi-feng,MA Ye-xüan
(State Key Laboratory of Hydraulic Engineering Simulation and Safety,Tianjin University,Tianjin 300072)

Multi-spans submarine pipelines are more normal than single-span ones at the bottom of deep sea.And their VIV characteristics are more complex.In this paper,the boundary condition in shoulders and pipe-soil interaction at the middle support were considered.A wake oscillator model for VIV prediction of multi-spans pipelines was introduced,and the model parameters were estimated.The influence of multispans length and boundary conditions in shoulders and middle support on the VIV of pipelines was analyzed.The results indicates that increasing the torsional spring stiffness and stretching spring stiffness may limit the occurrence of high mode,increasing the length of multi-spans submarine pipelines may cause high mode to easily happen.And it was also found that the energy could pass through the neighbor spans when the pipeline was suffered from VIV.

submarine pipelines;interaction of fluid and structure;Vortex-Induced Vibrations(VIV); multi-spans

TV312

A

10.3969/j.issn.1007-7294.2017.08.012

1007-7294(2017)08-1025-10

2016-12-28

國家自然科學基金(51209161、51379144、51479135、51679167);國家自然科學基金創新研究群體科學基金(51621092);國家重點基礎研究計劃-973計劃(2014CB046801)

徐萬海(1981-),男,副教授,碩士生導師,E-mail:xuwanhai@tju.edu.cn;謝武德(1989-),男,博士研究生,E-mail:tjuxwd@126.com。

猜你喜歡
模態振動
振動的思考
科學大眾(2023年17期)2023-10-26 07:39:14
噴水推進高速艇尾部振動響應分析
This “Singing Highway”plays music
振動攪拌 震動創新
中國公路(2017年18期)2018-01-23 03:00:38
中立型Emden-Fowler微分方程的振動性
車輛CAE分析中自由模態和約束模態的應用與對比
國內多模態教學研究回顧與展望
高速顫振模型設計中顫振主要模態的判斷
航空學報(2015年4期)2015-05-07 06:43:35
基于HHT和Prony算法的電力系統低頻振蕩模態識別
UF6振動激發態分子的振動-振動馳豫
計算物理(2014年2期)2014-03-11 17:01:44
主站蜘蛛池模板: 国产成年女人特黄特色大片免费| 国产欧美综合在线观看第七页| 在线欧美日韩国产| 国产美女丝袜高潮| 中文字幕日韩久久综合影院| 又黄又湿又爽的视频| 91亚洲精选| 欧美无遮挡国产欧美另类| av大片在线无码免费| 色一情一乱一伦一区二区三区小说| 国产在线观看人成激情视频| 久久青草精品一区二区三区| 国产精品欧美在线观看| 伊人AV天堂| 国产欧美日韩一区二区视频在线| 亚洲一区二区无码视频| 污污网站在线观看| 久久夜夜视频| 久热中文字幕在线| 2020国产在线视精品在| 午夜国产精品视频黄| 97超级碰碰碰碰精品| 青青青草国产| 日本色综合网| 色婷婷成人网| 内射人妻无码色AV天堂| 午夜爽爽视频| 国产伦精品一区二区三区视频优播 | 免费中文字幕在在线不卡| 99视频只有精品| 久久亚洲国产一区二区| 国产97公开成人免费视频| 欧美激情伊人| 国产精品露脸视频| 亚洲IV视频免费在线光看| 亚洲人成高清| 成人日韩欧美| 亚洲三级色| 毛片久久久| 国产成人免费视频精品一区二区| 成人毛片免费在线观看| 91人妻在线视频| 亚洲中文字幕精品| 欧美成人看片一区二区三区 | 香蕉99国内自产自拍视频| 国产精品19p| 波多野结衣视频网站| 精品丝袜美腿国产一区| 无码中文字幕精品推荐| 国产精品成人一区二区| 99这里只有精品免费视频| 一本色道久久88综合日韩精品| 手机精品视频在线观看免费| 日本三级欧美三级| 高清码无在线看| 亚洲a免费| 蜜桃臀无码内射一区二区三区 | 亚洲无码精彩视频在线观看| 波多野结衣AV无码久久一区| 国内精自视频品线一二区| www.99精品视频在线播放| 男女精品视频| 国产黄色爱视频| 久久国语对白| 国产精品香蕉在线| 日本精品αv中文字幕| 伊人久综合| 久久精品无码一区二区国产区| 国产免费福利网站| 亚洲综合久久一本伊一区| 91精品专区| 青青草91视频| 欧美在线视频不卡第一页| 国产精品成人一区二区不卡| 在线观看国产一区二区三区99| 久青草网站| 亚洲精品男人天堂| 精品久久久久久久久久久| 99资源在线| 精品国产Ⅴ无码大片在线观看81| 色哟哟色院91精品网站| 久久夜夜视频|