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

墜物撞擊海底管道損傷的BP神經網絡預測模型

2023-11-09 01:37:44張春會趙文豪田英輝牛習現佘虹宇齊曉亮
河北科技大學學報 2023年5期
關鍵詞:模型

張春會,趙文豪,田英輝,牛習現,王 樂,黃 鑫,佘虹宇,齊曉亮

(1.河北科技大學建筑工程學院,河北石家莊 050018;2.河北省巖土與結構體系防災減災技術創新中心,河北石家莊 050018; 3.墨爾本大學基礎設施工程系,維多利亞墨爾本 3010;4.河北青年管理干部學院信息系,河北石家莊 050031;5.天津大學水利工程仿真與安全國家重點實驗室,天津 300354;6.中海油能源發展股份有限公司采油服務分公司,天津 300451)

海底管道是海洋油氣輸送的生命線,在服役期內除遭受波流和溫度荷載作用外,還經常遭受拖網漁具、拋錨、集裝箱與其他墜落物等偶然荷載的撞擊,威脅海底管道安全[1-5]。已有統計數據結果表明,47%的海底管道故障是由這些偶然荷載撞擊引起的[6]。在墜物等偶然荷載撞擊下,海底管道的損傷預測是海底管道安全防護技術研究的基礎。目前,國內外對于墜物撞擊下海底管道損傷預測開展了很多研究[7-21]??傮w來看,目前國內外主要采用DNV RP-F107規范[7]公式預測海底管道損傷,但該公式假定海床為剛性,這與實際海底海床條件并不符。TIAN等[16]采用數值模擬方法分析了海床土力學特性對海底管道損傷的影響,給出了類似于DNV規范的海底管道損傷計算經驗公式。這些計算公式在計算管道損傷時都需要先確定墜物撞擊下海底管道的吸收能,然而海底管道吸收能受海床土力學特性、管道特性及墜物特性等諸多因素的影響,目前尚無有效的方法對其進行估算,針對具體工程只能通過數值模擬確定相應的管道吸收能,這給海底管道損傷預測帶來了極大不便。

BP神經網絡是一種典型的前饋型網絡,適合在大量數據基礎上建立復雜的非線性映射[22],為直接通過墜物、管道和海床土特性預測受撞擊海底管道損傷提供了新思路。本文在ABAQUS下建立墜物撞擊海底管道的數值模型,開展大量數值模擬,分析墜物撞擊能量、管道直徑與壁厚、鋼材等級、內壓和海床土力學特性對管道損傷的影響,將模擬結果作為訓練樣本建立管道損傷BP神經網絡預測模型,為海底管道損傷預測提供一種新的方法。

1 墜物撞擊海底管道的數值模型

本文在有限元軟件ABAQUS中建立墜物撞擊下海底管道與海床土相互作用的動力數值模型。

1.1 數值模型建立

海底管道為X60型鋼,管道直徑D為1 m,質量密度為7 850 kg/m3,楊氏模量為2×1011Pa,泊松比為0.3,管道為Mises材料,屈服強度為450 MPa,采用S4R殼體單元進行模擬。海床為飽和黏土,土體有效質量密度為600 kg/m3,墜物撞擊時間短,土體為不排水狀態,Tresca材料,內摩擦角與剪脹角均為0,泊松比為0.49,海床土不排水抗剪強度Su=10 kPa,彈性模量為500Su,土體采用C3D8R實體單元模擬。已有研究結果表明[16],在撞擊能量相同的條件下球形墜物造成的管道損傷最大。本文將墜物簡化為球形進行計算。球形墜物直徑為1 m,質量為1.9 t,將球形墜物視為剛體,采用C3D8R實體單元模擬,墜物撞擊點為管道中央。根據球形墜物撞擊海底管道的對稱性,沿球形墜物中心和撞擊點,取半模型進行數值分析,如圖1所示,其中xoz面為對稱面,y軸方向為沿管道軸向,z軸方向為海床高度方向。

圖1 墜物撞擊管道有限元模型圖Fig.1 Finite element model of subsea pipe striked by falling object

海床長寬高分別為16D,16D和10D,海底管道平鋪于海床上,長度為16D。墜物受到海水阻力,下落一定深度后速度不變[7],數值模型中半球形墜物從管道上方勻速豎直下落。海底管道與飽和黏土海床之間摩擦系數為0.3,海床土底部邊界固定,前后左右4個表面邊界只允許垂直移動。在管-土接觸長寬高方向1D范圍內網格加密。經過試算,加密區域網格尺寸為D/40,建立的球形墜物動力撞擊海底管道數值模型如圖1所示。

1.2 損傷定義

管道損傷示意圖見圖2。墜落物撞擊海底管道,Q1點為撞擊前管道的上部中心點,即球形落物撞擊管道的撞擊點,P點為撞擊前管道下底點,Q1與P之間的距離為管道直徑D,管道上部中心點受到豎向撞擊,Q1點產生豎向變形,向下移動至Q2點,Q1和Q2之間的距離即為撞擊引起的豎向凹陷深度δ,本文中管道損傷定義為

圖2 管道損傷示意圖Fig.2 Diagram of pipe damage

(1)

1.3 模型驗證

利用本文數值模型,模擬PINHEIRO等[23]開展的剛性床上支撐管道加載實驗,同時將本文數值模擬結果與DNV RP-F107[7]的規范公式計算結果進行對比,驗證本文數值模型的正確性。

1.3.1 剛性床上管道承受橫向荷載實驗

文獻[23]利用半球形壓頭對支撐在剛性底座上的管道(D=73 mm,t=3.05 mm)進行靜態橫向加載實驗,在距離中心點17 mm和25 mm處分別布設4個環向和軸向應變片,編號為SG1,SG2,SG3和SG4,如圖3所示,測試壓頭加載過程中4個測點的應變值。利用本文數值模型模擬實驗過程,在數值模擬中,壓頭和海床均為剛體,靠近管道-壓頭相互作用區域使用細小網格(0.05D),其余區域使用0.2D的單元尺寸,數值模型如圖4所示。

圖3 實驗模型圖Fig.3 Experimental model diagram

圖4 數值模型圖Fig.4 Numerical model diagram

數值結果和實驗結果對比見圖5,其中橫軸為應變,縱軸為壓頭位移。SG1和SG3沿軸向和環向測得的應變負值表示拉伸,正值表示收縮。由圖5可見,在加載初始階段,壓頭加載應變增長較緩;壓頭位移超過約1.7 mm后,應變增長速率變快。這是由于壓頭位移超過約1.7 mm后,管道初始圓形橫截面發生顯著變形,管道變得更加柔韌,呈現出更快的應變增長速率。由圖5可知,數值模擬結果與實驗結果基本吻合,表明本文數值模型可靠。

圖5 本文數值模擬結果與實驗結果對比Fig.5 Comparison of numerical simulation results with experimental results in this paper

1.3.2 數值模型與DNV規范對比驗證

DNV RP-F107規范[7]給出了剛性海床上海底管道受楔形墜物撞擊的凹痕損傷計算公式為

(2)

式中:E為管道吸收能量;σy為管道屈服應力;δ為管道凹陷深度;t為管道壁厚。

將本文數值模型中的海床改為剛性海床,墜物形狀修改為楔形,將楔形墜物視為剛體,楔形墜物采用C3D8R實體單元進行模擬,撞擊點為管道中央;海底管道為X60型鋼,管道直徑為1 m,質量密度為7 850 kg/m3,屈服強度為450 MPa,楊氏模量為2×1011Pa,泊松比為0.3,管道為Mises材料,采用S4R殼體單元進行模擬,取半模型進行數值分析。管道與海床之間光滑無摩擦,楔形墜物位于管道上方0.1 m處,海床剛性面的邊界條件被完全限制,管道的遠端被完全固定,管道的近端約束法向移動,楔形墜物只沿豎向運動,如圖6所示。算例中楔形墜物撞擊能量從3.5 kJ增加到197 kJ,數值計算獲得的管道損傷與DNV RP-F107規范計算結果對比如圖7所示。由圖7可知,數值模擬結果與DNV規范計算結果相差不超過5%,表明本文數值模型計算可靠。

圖6 數值模擬驗證模型Fig.6 Numerical simulation to validate the model

圖7 本文數值結果與DNV公式結果對比Fig.7 Comparison between the numerical results in this paper and the results of the DNV formula

1.4 計算方案

墜物撞擊能量、管道直徑、壁厚、鋼材等級及內壓、飽和黏土海床強度及模量是海底管道損傷的主要影響因素。為建立墜物撞擊下海底管道損傷預測的神經網絡模型,使用本文數值模型對飽和黏土海床上海底管道墜物撞擊開展數值模擬分析。改變墜物速度,控制墜物能量在42~198 kJ范圍內,墜物能量覆蓋了98%的能量帶[7]。依據API Spec 5L[24]管線規范,常用的海底管道實際管徑范圍通常為0.5~1 m,工程中管道壁厚通常為17~25 mm,鋼材等級范圍X52~X80。本文分別選取海底管道直徑為0.70,0.75,0.80,0.85和0.90 m,海底管道壁厚分別為19,20,21,22,和23 mm,管道鋼材等級選取X52,X56,X60,X65和X70。管道輸送油氣資源時,內部通常處于高壓狀態,常見壓力一般在10 MPa左右,最高不超過70 MPa[25-26]。本文分別取0(輸油管道關閉)和3,6,9,12 MPa共5種管道內壓值,常見海洋飽和黏土的不排水抗剪強度Su在10~50 kPa之間[27-29],飽和黏土模量為(200~1 000)Su。已有研究[16]表明,飽和黏土彈性模量對海底管道損傷影響微小,本文僅考慮不排水抗剪強度的影響,不排水抗剪強度取Su=10,20,30,40,50 kPa,飽和黏土的彈性模量取500Su。計算方案中的影響因素涵蓋了絕大多數實際工程中的飽和黏土海床海底管道工作條件。本文制定了105組算例,分別研究管道直徑、壁厚、鋼材等級、管道內壓以及海床不排水抗剪強度Su對海底管道損傷的影響,計算方案如表1所示。

表1 墜物撞擊海底管道計算方案Tab.1 Falling object impact calculation program for submarine pipelines

2 計算結果及分析

2.1 管道特性參數的影響

按照計算方案1—方案25開展計算,獲得不同管徑與墜物撞擊能量下海底管道損傷結果如圖8所示。從圖8可見,管道直徑為0.7 m時,墜物能量42.5 kJ時引起的管道損傷為5.01%,墜物能量197.5 kJ時引起的管道損傷為17.67%。相同管道直徑下,隨著墜物能量的增加,管道損傷增大。從圖8還可知,墜物能量為170 kJ時,管道直徑0.70,0.75,0.80,0.85,0.90 m引起的海底管道損傷依次為16.01%,15.36%,15.00%,14.70%和14.36%。這表明在相同墜物能量下,隨著管道直徑的增加,管道損傷減小。這是由于管道直徑增加,墜物和管道相接觸作用面的曲率減小,沖擊面積增大,對管道的撞擊作用減小,相應的管道損傷也減小。

圖8 海底管道直徑與損傷關系Fig.8 Relationship between submarine pipeline diameter and damage

根據計算方案21—方案45,得到不同壁厚與墜物撞擊能量下海底管道損傷結果如圖9所示。由圖9可見,管道壁厚為19 mm時,墜物能量42.5 kJ引起的管道損傷為6.66%,墜物能量為197.5 kJ時對管道引起的損傷為18.04%,相同管道壁厚下,墜物能量增加,管道損傷增大。由圖9還可以看出,墜物能量為170.0 kJ、管道壁厚分別為19,20,21,22,23 mm時,海底管道的損傷依次為16.52%,16.12%,15.24%,14.52%和14.07%,表明在相同墜物能量下,管道壁厚增加,管道損傷減小。這是由于管道壁厚增加會提升管道的剛度,管道抵抗沖擊荷載作用下塑性變形的能力也相應增加,從而使管道損傷相應減小。

圖9 海底管道壁厚與損傷關系Fig.9 Relationship between wall thickness and damage in submarine pipelines

依照計算方案21—方案25、方案46—方案65,獲得了管道不同鋼材等級與撞擊能量作用下海底管道損傷結果,如圖10所示。文中管道鋼材等級X52,X56,X60,X65和X70分別代表管道屈服強度為370,400,450,480,和520 MPa。由圖10可知,管道鋼材等級為X52時,墜物能量為42.5 kJ所引起的管道損傷為6.05%;墜物能量為197.5 kJ時,對管道引起的損傷為17.18%。相同管道鋼材等級下,墜物能量增加,管道損傷越大。由圖10還可知,墜物能量為170 kJ時,管道鋼材等級分別為X52,X56,X60,X65和X70所引起的管道損傷依次為15.73%,15.12%,14.37%,13.98%和13.43%,這表明相同墜物能量下,管道鋼材等級越高,管道損傷越小。這是由于管道鋼材等級提高會提升管道強度,與增加管道壁厚情形類似,管道在沖擊荷載作用下抵抗塑性變形的能力增加,從而使管道損傷減小。

圖10 海底管道鋼材等級與損傷關系Fig.10 Subsea pipeline steel grade and damage relationship

參照計算方案21—方案25和方案66—方案85,獲得了不同管道內壓值與墜物撞擊能量下海底管道損傷情況,如圖11所示。從圖11可見,管道內壓為0 MPa時,墜物能量42.5 kJ所引起管道損傷為4.51%,墜物能量197.5 kJ時對管道引起的損傷為15.46%,由此可見,在相同管道內壓情況下,墜物能量增大,管道損傷呈負指數函數增長且減小速率趨緩。從圖11還可知,墜物能量為170.0 kJ時,管道內壓為0,3,6,9和12 MPa所引起的管道損傷依次為14.37%,10.77%,8.66%,6.23%和4.20%,表明在相同墜物能量下,管道內壓增加,管道損傷減小。這是由于管道內壓增大,有助于提高管道抵抗外界撞擊的能力,減小管道損傷。

圖11 海底管道內壓與損傷關系Fig.11 Relationship between internal pressure and damage in submarine pipelines

2.2 海床土不排水抗剪強度的影響

土體的不排水抗剪強度Su是反映土體特性的重要力學指標,也決定了土體強度及極限承載能力。在沖擊過程中,不排水抗剪強度的大小會在海底管道受墜物撞擊時影響接觸土體引起破壞及破壞的程度,并且在墜物撞擊海底管道這一沖擊過程中,土體的不排水抗剪強度Su增大,海床土吸收能減小,從而影響海底管道的損傷程度。依照計算方案21—方案25和方案86—方案105,得到如圖12所示的管道損傷與不排水抗剪強度關系圖。由圖12可以看出,不排水抗剪強度Su=10 kPa時,墜物能量42.5 kJ所引起的管道損傷為4.51%,墜物能量197.5 kJ時對管道引起的損傷為15.46%,表明在相同不排水抗剪強度下,墜物能量增加,管道損傷增加。圖13分析了墜物能量為42.5 kJ時,不排水抗剪強度與管道吸收能量及損傷關系圖。不排水抗剪強度Su=10,20,30,40,50 kPa,所產生的管道損傷依次為4.51%,4.96%,5.22%,5.30%和5.37%,表明相同墜物能量時,土體不排水抗剪強度增加,管道損傷增大。這是由于海床土強度增加,土體變形和塑性區減小,土體吸收能量減小,管道吸收的能量增大,管道損傷增大。

圖12 海床土不排水抗剪強度與損傷關系Fig.12 Relationship between undrained shear strength and damage in seabed soils

圖13 海床土不排水抗剪強度與能量及損傷關系Fig.13 Undrained shear strength of seabed soils in relation to energy and damagea

3 海底管道損傷預測的BP神經網絡模型

本文建立的墜物撞擊下海底管道損傷預測的BP神經網絡模型是一種將誤差逆向傳播的多層前饋神經網絡,包括輸入層、中間層(也稱為隱藏層)和輸出層,同層單元之間互不連接。模型的輸入層包括墜物撞擊能量、管道直徑、壁厚、內壓和屈服強度、飽和黏土海床不排水抗剪強度6個參數,管道損傷為輸出層。輸入層6個參數經過中間隱藏層節點作用,與輸出層建立非線性映射關系,輸出管道損傷預測值。模型結構如圖14所示。

圖14 BP神經網絡結構圖Fig.14 BP neural network structure diagram

本文模型的BP神經網絡學習過程由正向傳播和誤差反向傳遞組成。在正向傳播過程中,輸入層不進行任何處理,隱藏層對輸入數據進行線性與非線性運算,并且通過 Sigmoid 函數將不同大小的輸出值控制在(0,1)的范圍內,最后將結果傳送到輸出層,Sigmoid函數為[30]

(3)

式中:x為輸入層的6個特征參數;f(x)為處理后的特征參數。

當預測誤差較大,即誤差達不到神經網絡模型設置的精度時,將正向傳播的預測值進行反向傳播。在誤差反向傳遞階段,按原來的通路逆向工作,采用梯度下降法對BP神經網絡的連接權值和誤差進行調整優化,使預測值和真實值更為接近;同時,通過調整學習因子與步長來調節權重矩陣和誤差,使誤差達到最小。其公式如下[31]:

(4)

式中:ΔWij為輸入層與隱藏層之間權重矩陣的誤差;n為時間步;η1為主要學習因子;η2為自適應學習因子;φi為輸入層輸出節點的計算誤差;Oj為隱藏層節點輸出結果;μ為動量因子;φ為步長;Hp(n)為時間步n的管道損傷值;Hp(n-1)為n-1時間步管道損傷值。

本文采用數值模擬的105組算例結果作為訓練樣本,以管道直徑、壁厚、鋼材等級、內壓、墜物撞擊能量和飽和黏土不排水抗剪強度為輸入參數,代入神經網絡模型,通過模型計算,輸出海底管道損傷,模型的允許誤差為10-5。經過模型訓練迭代52次,預測結果滿足精度要求,由此建立墜物撞擊下飽和黏土海床海底管道損傷預測的神經網絡模型。

4 BP神經網絡模型與已有研究結果的對比驗證

文獻[16]建立了墜物撞擊海底管道數值模型,計算了平鋪于黏土海床上海底管道受球形墜物撞擊引起的損傷,數值模擬中管道直徑0.4 m,管道壁厚12.7 mm,管道屈服強度360 MPa,楊氏模量為2×1011Pa,管道無內壓,海床土的彈性模量為5 MPa,墜物能量在0.1~40.7 kJ范圍,其10組數值計算結果如表2所示。利用本文神經網絡預測模型獲得相應的海底管道損傷預測值,也如表2所示。

表2 與文獻[16]海底管道撞擊模型驗證案例對比Tab.2 With documents [16] submarine pipeline impact model validation case

由表2可知,本文神經網絡模型預測結果與文獻[16]的數值結果最小偏差為0.31%,最大偏差為9.41%。可見,本文模型預測結果與文獻[16]的數值結果基本一致。文獻[16]數值計算中,部分墜物能量在DNV規范范圍之外,這可能是一些結果偏差較大的原因。總體上,本文神經網絡模型預測結果可靠。

文獻[32]建立的墜物撞擊海底管道數值模型,計算了墜物能量范圍在10~400 kJ、管道直徑0.8~1 m、壁厚17~25 mm、鋼材等級X52~X70、管道內壓0~8 MPa以及海床土不排水抗剪強度10~200 kPa條件下球形墜物撞擊引起的海底管道損傷。本文從文獻[32]的數值結果中選取10組作為本文神經網絡模型的驗證算例,10組算例結果如表3所示。該驗證算例的主要參數如下:海床土不排水抗剪強度為10 kPa、土的彈性模量為5 MPa,管道屈服強度為450 MPa,管道無內壓。

表3 與文獻[32]海底管道撞擊模型驗證案例對比Tab.3 With documents [32] submarine pipeline impact model validation case

由表3可知,利用本文神經網絡模型預測海底管道損傷,預測結果與文獻[32]的數值結果偏差在10%以內,數值結果基本相符,表明本文構建的BP神經網絡計算結果可靠。

5 結 語

本文在ABAQUS下使用顯式動力非線性分析方法建立了墜物撞擊海底管道損傷的有限元模型,利用該模型開展數值分析,結合實際海底管道工作條件參數變化范圍,研究了管道直徑、壁厚、鋼材等級、內壓、土體不排水抗剪強度和墜物撞擊能量6個因素對海底管道損傷的影響規律;將6個因素作為輸入層參數,以海底管道損傷作為輸出參數,以數值模擬結果作為訓練樣本,通過學習、訓練,構建形成了海底管道損傷預測的BP神經網絡模型,并對模型預測結果的可靠性進行了驗證,得到如下結論。

1)墜物撞擊能量越大,管道損傷越大,管道損傷增長速率隨墜物撞擊能量的增大而趨緩,墜物能量由42.5 kJ增加到197.5 kJ,管道損傷由5.01%增加至17.67%。

2)管道直徑從0.7 m增加到0.9 m,管道損傷由16.01%降低到14.36%;管道壁厚由19 mm增加至23 mm,管道損傷從16.52%降低到14.07%。管道直徑、壁厚越大,管道損傷越小;鋼材等級越高,管道內壓越大,管道損傷越小。

3)飽和黏土海床不排水抗剪強度增加,海底管道損傷增加,這主要是由于不排水抗剪強度越大,土體吸收的撞擊能量減小,管道吸收撞擊變形的能量增加。

4)建立的海底管道損傷BP神經網絡預測模型能夠較好地預測飽和黏土海床海底管道受墜物撞擊的損傷情況。另外,本文數值算例涵蓋了常見飽和黏土海床海底管道的工作條件,建立的神經網絡模型預測海底管道損傷僅需墜物撞擊能量、管道直徑、壁厚、鋼材等級、內壓和海床土不排水抗剪強度6個參數,模型簡單、便捷,適用性好,為海底管道損傷預測提供了新思路。

下一步將進一步探究管道埋深、雙層管道結構等對海底管道損傷的影響及預測方法。

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 国模沟沟一区二区三区| 天天躁夜夜躁狠狠躁图片| 久久综合伊人 六十路| 国产中文在线亚洲精品官网| 欧美日韩国产系列在线观看| 亚洲无码精彩视频在线观看| 欧美成人在线免费| аⅴ资源中文在线天堂| 亚洲黄色激情网站| 黄色不卡视频| 2020亚洲精品无码| 伊人天堂网| 欧美日韩第二页| 激情国产精品一区| 午夜少妇精品视频小电影| 亚洲最大看欧美片网站地址| 午夜丁香婷婷| 99国产在线视频| 国产精品片在线观看手机版| 国产成人AV综合久久| 国产自产视频一区二区三区| 国产视频自拍一区| 五月婷婷欧美| 98超碰在线观看| 91无码网站| 成人在线观看不卡| 精品福利一区二区免费视频| 中文字幕欧美日韩高清| 无码人中文字幕| 日本91在线| 91在线激情在线观看| 成人精品视频一区二区在线| 欧美成人国产| 亚洲国产精品不卡在线| 亚洲一区色| 国产自在线拍| 亚洲侵犯无码网址在线观看| 国产人人乐人人爱| 国产在线小视频| 婷婷亚洲最大| 免费A∨中文乱码专区| 亚洲精品人成网线在线| 亚洲天堂精品视频| 中文字幕波多野不卡一区| 成人综合久久综合| 免费A级毛片无码免费视频| 91亚洲国产视频| 国产欧美中文字幕| 亚国产欧美在线人成| 尤物特级无码毛片免费| 重口调教一区二区视频| 国产精品亚洲片在线va| 日韩小视频网站hq| 久草视频中文| 国产成人精品视频一区视频二区| 日韩精品成人在线| 国产白浆在线观看| 精品久久高清| 亚洲综合激情另类专区| 亚洲AV无码久久天堂| 日本欧美中文字幕精品亚洲| 亚洲美女一区| 久久综合九色综合97婷婷| 国产成人凹凸视频在线| 亚洲一区二区三区在线视频| 九色国产在线| 精品国产Ⅴ无码大片在线观看81 | 国产乱人伦精品一区二区| 在线亚洲小视频| 成年A级毛片| julia中文字幕久久亚洲| 99久久精彩视频| 国产精品亚洲一区二区三区z| 精品国产香蕉在线播出| 九一九色国产| 高清精品美女在线播放| 亚洲码一区二区三区| 久久人人97超碰人人澡爱香蕉 | 色综合五月婷婷| 色呦呦手机在线精品| 欧美一区二区精品久久久| 亚洲精品成人片在线播放|