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

近水面水平圓柱在波浪作用下的水動力系數

2020-07-01 05:06:58毛鴻飛赫巖莉袁劍平潘新祥賈寶柱
廣東海洋大學學報 2020年4期
關鍵詞:水平

毛鴻飛,赫巖莉,袁劍平,潘新祥,賈寶柱

近水面水平圓柱在波浪作用下的水動力系數

毛鴻飛1,赫巖莉1,袁劍平1,潘新祥2,賈寶柱2

(1. 廣東海洋大學海洋工程學院,廣東 湛江 524088;2. 廣東海洋大學海運學院,廣東 湛江 524088)

【目的】研究大幅非線性波浪作用下近水面水平圓柱上波浪力和水動力系數的特征。【方法】基于黏性流理論和有限體積方法開發數值波浪水槽模型,分析入射波波幅和波頻率對波浪力特征的影響,提出可以考慮圓柱完全暴露于空氣情況的改進Morison方程,通過大量的數據擬合得到水動力系數,并討論改進的Morison方程對波浪力預測的適用性。【結果】圓柱上波浪力隨波幅和頻率增大逐漸增大,而由于浮力占重要比例,正垂向力受波頻率的影響并不顯著;所得近水面水平圓柱慣性力系數和速度力系數最大值比淹沒圓柱在線性波浪作用下的水動力系數大。【結論】所提出改進的Morison方程對波浪力的預測有較好的適用性。

改進的Morison方程;數值模擬;水平圓柱;波浪力;水動力系數

海洋工程中,水平圓柱結構的應用非常廣泛,如海上輸油管線等和海上采油平臺導管架等。為保證近水面水平圓柱結構的安全,避免坍塌和疲勞破壞,需要對其在波浪作用下的受力特征進行研究,在設計階段對波浪力進行準確預測。早期對近水面水平圓柱上的波浪力和水動力系數的研究大多采用物理試驗方法和Morison方程預測方法[1]。Dixon等[2]和Easson等[3]先后對小幅波浪作用下部分淹沒水平圓柱上的波浪力和水動力系數開展了物理試驗研究,提出了在考慮慣性力、速度力和浮力后改進的Morison方程,在垂向方向上擬合水動力系數。其結果表明,慣性力系數小于2.0,且速度力系數分散性強。Chaplin[4]對近水面淹沒水平圓柱上的波浪力開展了物理試驗研究,發現慣性力系數隨著波幅增大的變化趨勢為先減小后增大,并分析其原因是有旋流動產生的升力以及邊界層分離的影響。Contento等[5]也對淹沒水平圓柱上的波浪力開展了物理試驗研究,考察了低Kc數(Keulegan–Carpenter number)下圓柱淹沒深度對波浪力的影響,結果表明淹沒深度對波浪力的高倍頻分量影響比較顯著。

隨著計算流體動力學(Computational Fluid Dynamics,CFD)的發展,對近水面水平圓柱上波浪力相繼開展了數值模擬研究。Guerber等[6-7]基于勢流理論和邊界元方法建立了二維數值水槽模型,對非線性波浪作用下近水面淹沒水平圓柱上的波浪力進行計算,獲得了波浪力及其高倍頻分量。Hu等[8]基于黏性流理論建立數值水槽模型,對線性波浪作用下半淹沒水平圓柱上的波浪力進行了數值計算,其數值結果與前人實驗數據吻合較好。Westphalen等[9]應用基于黏性流理論采用四種數值方法建立水槽模型,對線性波浪作用下半淹沒水平圓柱上的波浪力進行計算,數值結果與前人實驗數據總體上吻合較好。Bai等[10]基于勢流理論,采用邊界元方法建立了完全非線性三維數值水槽模型,計算了近水面淹沒水平圓柱上波浪力。Chen等[11]應用基于黏性流理論所建立的二維數值水槽模型,對部分淹沒水平圓柱和方柱上的波浪力進行了數值計算,并分析了水平和垂向波浪力的特征。研究發現入射波長對應慣性力有影響。Teng等[12]基于黏性流理論,采用有限體積法建立了數值水槽模型,對近水面淹沒水平圓柱上的波浪力開展了計算研究,解釋了黏性流與勢流理論下慣性力差距隨波幅增大而增大的原因。毛鴻飛等[13-14]應用二維數值波浪水槽模型,分別對波浪作用下非完全淹沒水平圓柱所受沖擊力和水動力系數開展了計算,研究發現大波幅下沖擊力在波浪力中占主導,但研究中未討論所提出簡化計算方程的適用性。

上述對波浪作用下近水面水平圓柱上波浪力的研究尚未全面系統地考察大幅波浪作用下圓柱所受波浪力的特征,改進Morison方程對波浪力的預測無法考慮大幅波浪下圓柱完全暴露于空氣中的情況,對水動力系數的擬合僅在垂向單一方向。借鑒Gao等[15-17]采用黏性流CFD模型對窄縫共振問題研究和毛鴻飛等[18-19]對水面附近平板受波浪沖擊力研究的成功經驗,本研究基于黏性流理論,采用有限體積法建立數值波浪水槽模型,并應用該模型,針對上述前人研究的局限性,對大幅非線性波浪作用下近水面水平圓柱上的波浪力特征和水動力系數進行數值計算研究。首先,通過與物理實驗數據的對比驗證數值模型的準確性。進而,通過不同入射波波幅和波頻率下波浪力的對比,分析波浪參數對波浪力特征的影響。進一步,提出可以考慮大幅波浪作用下,圓柱完全暴露于空氣中情況的改進的Morison方程,并在水平和垂向方向上對水動力系數聯立求解,獲得其參考值。最后,討論改進的Morison方程對波浪力預測的適用性。

1 模型原理和改進的Morsion方程

1.1 數學模型

數值模型控制方程以Navier-Stokes方程,連續性方程和動量守恒方程可表示為:

式中,u為流體質點速度在方向上的分量,為時間,為流體的密度,為流體壓強,為流體動力學黏性系數,f為體積力。

兩相流自由面采用VOF方法進行捕捉。將定義為體積分數,氣相中= 0,液相中= 1,自由面內= 0 ~ 1。

進而,單元密度和黏性系數表示為:

數值波浪水槽模型的造波和消波功能通過采用速度邊界法和松弛方法相結合的形式來實現,邊界和松弛區定義如圖1所示。圖中,“入口邊界”按初始流體速度所需波浪理論給定,壓強為?*/?= 0,其中*=?為動態壓強,為某點到自由面的距離;“出口邊界”和“底邊界”均給定不可滑移邊界條件,即初始速度和壓強為u= 0,?*/?= 0;“頂邊界”給定可自由進出邊界條件,即初始速度和壓強為?u/?= 0,*= 0。水槽中結構物為固定不動條件下,物面邊界按照不可滑移邊界條件給定。松弛區I作用是波浪生成和吸收從結構物反射波浪,松弛區II作用是消除消波邊界反射波浪[21]。松弛區內,根據解析形式對數值求解過程進行修正,從而實現穩定長時間的造波和消波功能,該區域內流體速度、壓強和體積分數的解析修正形式為:

圖1 邊界條件

松弛函數的分布如圖1所示。

在獲得壓強和速度場分布后,作用在結構物上的波浪力可通過下式進行計算,

1.2 改進的Morison方程

傳統的Morison方程主要針對完全淹沒的柱體結構受力簡化計算,而對部分淹沒柱體上波浪力計算,需將方程中速度投影面積和圓柱淹沒體積轉變為以時間變量的函數形式,Dixon[2]和Easson等[3]相繼改進了Morison方程,但僅針對柱體始終處于部分淹沒狀態,無法考慮大波幅情況下柱體完全暴露于空氣中的情況。為了解決這樣的局限性問題,針對大幅入射波作用下水平圓柱在的波浪力預測,本研究提出了改進的Morison方程,即

式中,、FFF分別為總波浪力,慣性力,速度力和浮力,下標和表示水平分量和垂向分量;CC分別為慣性力和速度力系數;() 為圓柱在任意時刻淹沒體積,0為圓柱初始時刻淹沒體積;()為任意時刻圓柱淹沒部分的速度投影面積。在波浪波幅和波長相對于圓柱尺度較大條件下,簡化波浪與圓柱作用的幾何關系,如圖2所示。

圖2 圓柱和波面幾何關系

將結構附近自由面簡化為一條水平線,有如下幾何關系:自由面高度可表示為

圓柱底部到自由面垂向距可表示為:

圓柱和自由面形成的圓心角度可表示為:

從而,將結構淹沒體積和速度投影面積以時間相關函數表達為

2 模型驗證計算

2.1 計算設置

為了驗證數值模型準確性,以經典物理問題“淹沒水平圓柱上的波浪力”進行數值驗證。參考Chaplin[4]的物理實驗,驗證計算的計算域和參數設置如圖3所示。圖中,為圓柱半徑,為圓柱淹沒深度,即圓柱軸心距靜水面距離,為水深,為入射波波長。圓柱半徑為5.1 cm, 水深為0.80 m,圓柱的淹沒深度為0.255 m。圓柱軸心所在水平方向位置,松弛區I和松弛區II的長度均根據入射波波長的一定比例進行設置:圓柱位置距離造波和消波區分別為1.65和4.0,造波和消波區長度分別設置為2.5和3.0。以5階Stokes理論[22]為造波理論生成入射波,周期考慮為= 1.20 s,選擇多個入射波波幅。

圖3 驗證計算域示意圖

2.2 波浪力驗證對比

圖4為垂向波浪力隨/(波幅的無量綱化形式)變化情況,本研究數值結果與前人實驗結果的對比。為了便于對比分析,波浪力無量綱化采用形式為F/(π22e-kz)。從對比結果可見,本研究數值結果總體上與前人實驗結果吻合較好,這驗證了對于波浪作用下結構物受力的計算,本研究數值模具備正確性,并具有較好的計算精度。

圖4 垂向波浪力隨波幅變化情況

3 波浪力數值計算和結果分析

3.1 算例設置

波浪作用于近水面水平圓柱計算域示意圖見圖5。圖中,為圓柱軸心至靜水面的距離,圓柱半徑為5.1 cm,水深為1.785 m。圓柱水平位置距離造波和消波區分別為1.65和4,造波和消波區長度分別設置為2.5和3,圓柱垂向位置確定為其軸心至于靜水面及其以上,圓柱垂向位置為= 0.00 ~ 0.28 m。選取多組入射波幅和波頻率,波幅為= 1.28 ~ 30.60 cm;頻率的無量綱化用表示,其中為波數,選取值為0.05、0.10、0.15、0.20、0.30、0.40。

3.2 波浪力特征

圖6為不同波浪頻率下,波浪力隨波幅變化的對比。圖中,用正負極值來描述波浪力大小,并為了便于對比分析,波浪力極值轉化為無量綱形式/(2)表示。從圖中對比可見,由于波浪場中流體加速度和速度與波幅成正比例關系,無量綱波浪力總體上呈現出隨著波幅增大而增大的現象。在部分工況下,正垂向力隨著波幅增大而增大的速率變快(如/= 4.50 ~ 6.00,/= 0.00,= 0.05,其中/為波幅的無量綱形式)。這種現象可解釋為,波幅較大情況下的流體運動速度較大,圓柱受到了一定的沖擊作用,即在波浪接觸到圓柱的短時間內波浪力突然增大,出現一個峰值,甚至在一定條件下,沖擊力在總波浪力中占主導作用,即沖擊力峰值大小為波浪力極值[13]。正、負水平力呈現出隨頻率增大而變大的現象;而正垂向力則呈現出不同頻率下波浪力接近,甚至一些情況下出現頻率較小時波浪力較大的現象。這種現象可以解釋為,波浪場中流速與頻率成正比關系,加速度與頻率平方成正比關系;波浪頻率較小情況下波長較長,自由面平緩,且當圓柱垂向位置較低時,浮力的影響比較顯著。

圖5 波浪作用于近水面水平圓柱計算域示意圖

圖6 振蕩方向上的一倍頻流體作用力

Fig. 6 The first harmonic components of the hydrodynamic forces in the oscillation direction

4 水動力系數特征及方程適用性

4.1 水動力系數

擬合水動力系數,需要排除沖擊力對其影響。首先,利用數據點時間間隔特征,將數值結果中時間間隔小于0.002時間步長的沖擊力數據點刪除;其次,采用數據點差分補充方法填補波浪力完整數據,形成連續數據化,即時間歷程不間斷;再次,將淹沒體積和速度投影面積時間函數帶入改進的Morison方程,即將式(15)、(16)、(17)代入式(18)等號左側項,將波浪力數值結果也代入式(18)等號右側項;最后,采用最小二乘法對方程組求解。

在深水條件下,對于波浪與近水面水平圓柱作用問題,Kc數可定義為

式中,為水質點流速,為圓柱直徑。

圖7為不同條件下,水動力系數隨Kc數的分布情況,圖7-a為慣性力系數,圖7-b為速度力系數。從圖中系數分布特征可見,慣性力系數總體上隨著Kc數增大逐漸減小,而速度力系數隨著Kc數增大逐漸增大;隨著增大,慣性力系數變化范圍不大,而速度力系數逐漸減小。在本研究計算范圍內,獲得的慣性力系數和速度力系數最大值分別為C= 2.20和C= 2.54,大于淹沒圓柱在線性波浪作用下的水動力系數參考值C= 1.2,C= 2.0。

圖7 水動力系數分布

4.2 改進的Morison方程適用性

將得到的水動力系數帶入改進的Morison方程,并以/= 1.00,= 0.10,/= 3.00為例,將方程預測的波浪力和數值結果進行對比,如圖8所示。可見,改進的Morison方程預測的波浪力與數值結果總體上吻合較好,表明本研究改進的Morison方程對近水面水平圓柱所受波浪力的預測有較好的適用性。

圖8 波浪力時間歷程

5 結論

本研究基于黏性流理論,采用有限體積方法建立數值波浪水槽模型,并應用該模型對大幅非線性波浪作用下近水面水平圓柱上的波浪力進行數值計算,提出了可以考慮圓柱完全暴露于空氣情況的改進Morison方程,通過數據擬合獲得了水動力系數,并討論了改進的Morison方程對波浪力預測的適用性。研究結論總結如下:

(1)波浪力總體上呈現隨波幅增大而逐漸增大的趨勢;水平力和負垂向力體現出隨波浪頻率增大而逐漸增大的特點,而受浮力的影響,正垂向力受頻率的影響并不顯著。

(2)近水面水平圓柱在大幅波浪作用下的水動力系數最大值比淹沒圓柱在線性波浪作用下的系數參考值大,本研究范圍內獲得了水動力系數參考值,慣性力系數最大值為2.20,速度力系數最大值為2.54。

(3)以分段函數形式表達近水面水平圓柱在大幅波浪作用下的淹沒體積和速度投影面積,進而提出的改進的Morison方程,對圓柱上波浪力的預測具有較好的適用性。

[1] MORISON J R. The force exerted by surface waves on piles[J]. Petroleum Transaction, American Institute of Mining, Metallurgical and Petroleum Engineers, 1950, 189: 149-154.

[2] GRAHAM DIXON A, SALTER S, GREATED C. Wave forces on partially submerged cylinders[J]. Journal of the Waterway, Port, Coastal and Ocean Division, 1979, 105(4): 421-438.

[3] EASSON W J, GREATED C A, DURANNI T S. Force spectra from partially submerged circular cylinders in random seas[J]. Journal of Waterway Port Coastal and Ocean Engineering-case, 1985, 111(5): 856-879.

[4] CHAPLIN J R. Nonlinear forces on a horizontal cylinder beneath waves[J]. Journal of Fluid Mechanics, 1984, 147: 449-464.

[5] CONTENTO G, CODIGLIA R. Non-linear free surface induced pressure on a submerged horizontal circular cylinder at low Keulegan-Carpenters numbers[J]. Applied Ocean Research, 2001, 23(3): 175-185.

[6] GUERBER E, BENOIT M, GRILLI S T, et al. Modeling of fully nonlinear wave interactions with moving submerged structures[C]. Proceedings of the International Offshore and Polar Engineering Conference, 2010, 529-536.

[7] GUERBER E, BENOIT M, GRILLI S T, et al. A fully nonlinear implicit model for wave interactions with submerged structures in forced or free motion[J]. Engineering Analysis with Boundary Elements, 2012, 36(7): 1151-1163.

[8] HU Z Z, CAUSON D M, MINGHAM C G, et al. Numerical simulation of floating bodies in extreme free surface waves[J]. Natural Hazards and Earth System Sciences, 2011, 11(2): 519-527.

[9] WESTPHALEN J, GREAVES D M, RABY A, et al. Investigation of wave-structure interaction using state of the art CFD techniques[J]. Open Journal of Fluid Dynamics, 2014, 4(1): 18-43.

[10] BAI W, HANNAN M A, ANG K K. Numerical simulation of fully nonlinear wave interaction with submerged structures: Fixed or subjected to constrained motion[J]. Journal of Fluids and Structures, 2014, 49: 534-553.

[11] CHEN B, LU L, GREATED C A, et al. Investigation of wave forces on partially submerged horizontal cylinders by numerical simulation[J]. Ocean Engineering, 2015, 107: 23-31.

[12] TENG B, MAO H F, LU L. Viscous effects on wave forces on a submerged horizontal circular cylinder[J]. China Ocean Engineering, 2018, 32(3): 245-255.

[13] 毛鴻飛, 陳洪洲. 非完全淹沒水平圓柱上波浪力特征的數值模擬[J]. 水科學進展, 2019, 30(5): 749-759.

[14] 毛鴻飛, 滕斌. 波浪作用下近水面水平圓柱的水動力系數研究[C]. 重慶: 第十九屆中國海洋(岸)工程學術探討會論文集, 2019, 402-408.

[15] GAO J L, HE Z W, ZANG J, et al. Topographic effects on wave resonance in the narrow gap between fixed box and vertical wall[J]. Ocean Engineering, 2019, 180: 97-107.

[16] GAO J L, ZANG J, CHEN L F, et al. On hydrodynamic characteristics of gap resonance between two fixed bodies in close proximity[J]. Ocean Engineering, 2019, 173: 28-44.

[17] GAO J L, HE Z W, ZANG J, et al. Numerical investigations of wave loads on fixed box in front of vertical wall with a narrow gap under wave actions[J]. Ocean Engineering, 2020, 206: 107323.

[18] 毛鴻飛, 李芳成, 吳光林, 等. 基于黏性流理論對平板受波浪沖擊的兩相流數值研究[J]. 廣東海洋大學學報, 2019, 39(4): 73-80.

[19] 毛鴻飛. 基于粘性流理論的水平圓柱上波浪力數值研究[D]. 大連: 大連理工大學, 2018.

[20] WELLER H G, TABOR G, JASAK H, et al. A tensorial approach to computational continuum mechanics using object oriented techniques[J]. Computers in Physics, 1998, 12(6): 620.

[21] JACOBSEN N G, FUHRMAN D R, FREDSOE J. A wave generation toolbox for the open‐source CFD library: OpenFoam?[J]. International Journal for Numerical Methods in Fluids, 2012, 70(9): 1073-1088.

[22] FENTON J D. A fifth-order Stokes theory for steady waves[J]. Journal of Waterway, Port, Coastal and Ocean Engineering, 1985, 111: 216-234.

Hydrodynamic Coefficients of a Horizontal Circular Cylinder near Free Surface under Wave Action

MAO Hong-fei1, HE Yan-li1, YUAN Jian-ping1, PAN Xin-xiang2, JIA Bao-zhu2

(1.,,524088,; 2.,,524088,)

【Objective】This study investigates wave forces and hydrodynamic coefficients on a horizontal circular cylinder near the free surface under large-amplitude wave actions.【Method】A numerical wave tank model is developed based on the viscous fluid theory and finite volume method. By analyzing the characteristics of the wave forces on the circular cylinder, the influences of wave amplitude and wave frequency on the wave forces are investigated. A modified Morison’s equation for predicting the wave forces is derived, consider the condition of the circular cylinder is completely exposed to the air phase under the action of large-amplitude wave. The inertia and drag coefficients are obtained by solving the modified Morison’s equation. The applicability of the modified Morison’s equation for predicting the wave forces on the circular cylinder is discussed.【Result】The wave forces on the circular cylinder increase gradually with the increase of wave amplitude and wave frequency. The effects of wave frequency on the positive vertical waves are not significant because the buoyancy forces take up an important proportion. The maximums of the hydrodynamic coefficients obtained are larger than those for a submerged circular cylinder under linear wave action.【Conclusion】The modified Morison’s equation works well for the prediction of the wave forces.

modified Morison's equation; hydrodynamic coefficient; numerical simulation; horizontal circular cylinder; wave force

TV139.2+6

A

1673-9159(2020)04-0109-07

10.3969/j.issn.1673-9159.2020.04.015

2020-04-16

廣東省自然科學基金(51979045);國家國防科工局穩定支持課題(JCKYS2019604SXJQR-02);廣東省教育廳高校青年創新人才項目(2019KQNCX045);廣東海洋大學科研啟動費項目(120602-R19024)

毛鴻飛(1985-),男,博士,講師,研究方向為波浪與結構物的相互作用。E-mail:maohongfei-gdou@qq.com

袁劍平(1979-),男,碩士,高級工程師,研究方向為船舶智能控制系統開發。E-mail: yjp_103@163.com

毛鴻飛,赫巖莉,袁劍平,等. 近水面水平圓柱在波浪作用下的水動力系數[J].廣東海洋大學學報,2020,40(4):109-115.

(責任編輯:劉朏)

猜你喜歡
水平
張水平作品
作家葛水平
火花(2019年12期)2019-12-26 01:00:28
深化精神文明創建 提升人大工作水平
人大建設(2019年6期)2019-10-08 08:55:48
加強上下聯動 提升人大履職水平
人大建設(2019年12期)2019-05-21 02:55:32
水平有限
雜文月刊(2018年21期)2019-01-05 05:55:28
加強自身建設 提升人大履職水平
人大建設(2017年6期)2017-09-26 11:50:44
老虎獻臀
中俄經貿合作再上新水平的戰略思考
建機制 抓落實 上水平
中國火炬(2010年12期)2010-07-25 13:26:22
做到三到位 提升新水平
中國火炬(2010年8期)2010-07-25 11:34:30
主站蜘蛛池模板: 国产激爽大片高清在线观看| 伊人久久青草青青综合| h视频在线观看网站| 国产xxxxx免费视频| 欧美 亚洲 日韩 国产| 国产精品永久免费嫩草研究院| 天天做天天爱夜夜爽毛片毛片| 亚洲色图另类| 超碰aⅴ人人做人人爽欧美 | 国禁国产you女视频网站| 亚洲床戏一区| 亚洲欧美成人在线视频| 在线观看视频一区二区| 超级碰免费视频91| 无码专区国产精品第一页| 91精品啪在线观看国产60岁 | 亚洲第一视频网| 青草国产在线视频| 国产97区一区二区三区无码| 日本一区二区不卡视频| 日韩免费成人| 国产一级妓女av网站| 天天躁夜夜躁狠狠躁图片| 午夜福利视频一区| 欧美爱爱网| 99热这里只有免费国产精品 | 日韩欧美国产三级| 香蕉蕉亚亚洲aav综合| 午夜国产精品视频| 国产精品无码翘臀在线看纯欲| 国内精品视频在线| 日韩第一页在线| 久久婷婷国产综合尤物精品| 在线永久免费观看的毛片| 日韩欧美国产区| 97国产在线观看| 婷婷色婷婷| 久久久久无码精品国产免费| 伊人久久大香线蕉aⅴ色| 免费精品一区二区h| 97人人做人人爽香蕉精品| 精品国产aⅴ一区二区三区| 第一区免费在线观看| 久久精品波多野结衣| 97色伦色在线综合视频| 韩国自拍偷自拍亚洲精品| 少妇精品网站| 黄色三级网站免费| 91香蕉视频下载网站| 国产美女人喷水在线观看| 欧美成人精品一级在线观看| 狠狠亚洲五月天| 亚洲三级网站| 亚洲第一黄色网址| 五月婷婷导航| h网站在线播放| 国产麻豆另类AV| 欧美成人在线免费| 久久网欧美| 欧美日一级片| 亚洲,国产,日韩,综合一区 | 中文字幕波多野不卡一区| 国产一区二区精品福利| 久久精品娱乐亚洲领先| 女人18毛片一级毛片在线 | 在线观看国产精品日本不卡网| 久久免费视频6| 99精品视频在线观看免费播放| 欧美亚洲国产日韩电影在线| 久久永久精品免费视频| 精品少妇三级亚洲| 国产亚洲欧美在线人成aaaa | AⅤ色综合久久天堂AV色综合| 久久综合五月婷婷| 欧美性色综合网| 99性视频| 国产成年无码AⅤ片在线| 尤物在线观看乱码| 激情五月婷婷综合网| 一本色道久久88综合日韩精品| 亚洲精品亚洲人成在线| 亚洲一欧洲中文字幕在线|