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

海底管道軸向移動數值模擬研究

2016-11-11 05:45:00陳志華黃金超趙思玥何永禹劉紅波
哈爾濱工程大學學報 2016年9期

陳志華,黃金超,趙思玥,何永禹,劉紅波

(1.天津大學 水利工程仿真與安全國家重點實驗室,天津 300072;2.天津大學 建筑工程學院,天津 300072)

?

海底管道軸向移動數值模擬研究

陳志華1,2,黃金超2,趙思玥2,何永禹2,劉紅波1,2

(1.天津大學 水利工程仿真與安全國家重點實驗室,天津 300072;2.天津大學 建筑工程學院,天津 300072)

非線性的溫度荷載是造成管道產生軸向移動的主要原因,針對海底管道在多次開啟/關閉的循環作用下軸向移動現象,采用有限元分析的方法,通過ABAQUS軟件建立循環溫度作用下平坦海床上短直管道有限元模型,模型中考慮了非線性的溫度加載梯度和高溫對管道材料的影響,揭示短管道在多次開啟/關閉循環荷載作用下有效軸向力分布以及管道各點軸向移動現象的形成原因及發展規律。最后通過參數化分析,得到了管土摩擦系數與管道軸向移動速率之間的關系,并對比兩種管道內壓的加載方式對管道軸向移動速率的影響,為工程中減緩管道軸向移動速率的方法提供借鑒和幫助。

短直海底管道;數值模擬;軸向移動;非線性循環溫度;參數化分析

對于海底管道軸向移動研究,Torne等[1]應用庫倫管土摩擦模型建立有限元模型,得出管道內有效軸力的分布和管道軸向移動形態;Carr等[2]提出海底管道軸向移動的3個條件:立管連接端的張力、管道長度方向上海床的傾斜度、管道開啟/關閉過程中的溫度效應,進而提出管道軸向移動速率的解析方程;David等[3]研究了管土作用對管道軸向移動和側向屈曲的影響;劉羽霄等[4]研討了管道軸向移動和側向屈曲的相互作用;Rong等[5]利用ABAQUS有限元軟件,考慮了非線性管土摩擦與材料非線性建立埋地管道軸向移動模型,提出短管道更容易發生軸向移動現象的結論。

對于深海管道來說,在高溫和高壓環境下軸向膨脹現象和整體屈曲[6]尤為突出,并且由于經濟技術原因,深海管道多是直接鋪設于海床之上,其邊界條件對于軸向移動的阻抗作用相比于埋地管道也更小,多種不利因素高溫高壓、直接鋪設于平坦海床之上都大大加大軸向移動現象發生,因此本文研究背景為處于高溫高壓直接鋪設于海床上環境條件下管道軸向移動現象的研究。

其次海底管道的開啟/關閉的循環作用對于軸向移動起控制作用,但長直管道其沿管道長度方向溫度梯度較小,且在使用壽命期內管道開啟/關閉的次數較少,開啟/關閉循環引起的管道軸向位移較小,可以控制在設計范圍之內。因此,在深海高溫高壓平坦海床環境條件下,本文著重對短直管在多次開啟/關閉循環軸向有效應力分布及各點軸向移動進行了分析詳盡分析,最后對管土摩擦以及管道內壓進行參數化分析,得到軸向移動現象的形成原因及發展規律。

1 管道膨脹產生原因

引起海底管道膨脹的主要原因為以下3種:溫度、內外壓和材料泊松比,主要由以下幾部分組成[7]:

1)端帽效應即由管道內部流體介質作用在管道端部面上從而產生的應變εp;

2)泊松比效應引起的管道壓縮應變εv是由管道內壓產生的環向應力引起的,導致管道在軸向方向上發生收縮;

3)管道安裝和運營過程中溫度的變化會引起管道溫度應力的產生,在自由條件下,溫度變化引起的溫度應變為εt;

4)阻抗應變是由于管道和河床之間的摩擦阻力阻礙著管道軸向膨脹產生的應變εf。

管道的總的凈應變εnet如下式:

(1)

根據上式可以看出沿著管道長度方向有一個特定點凈應變為零,此點稱為錨固點,對于長直管道會出現兩個錨固點,位于兩個錨固點之間的管道處于完全約束狀態,在這些管道的軸向應變為零的區域稱為錨固段。對于短直管道,管道長度沒有達到完全約束的位置,通常在管道靠近中點的位置有效軸力達到最大,同時形成虛擬錨固點。與長直管道不同的是,虛擬錨固點的位移為零,但是軸向應變不為零,管道的最大有效軸力遠小于全約束軸向力。

2 短直管道軸向移動解析計算

假設管道為均勻對稱模型,且海床為剛性海床,考慮小變形和線彈性假設,溫度引起的管道軸向荷載[8]fθ為:

(2)

(3)

式中:As為管道凈截面面積,E為彈性模量,αT為熱膨脹系數,T為管道內外溫差,L為管道長度。

溫度荷載作用下管道的軸向移動速率計算為

(4)

(5)

式中:f管土摩擦阻力,μa為管道與海床的摩擦系數,Wsub為管道的淹沒重量。

由式(5)可以看出管道軸向移動現象與管道海床之間的摩擦系數、管道重量等有關。本節對管道軸向移動速度的解析結果進行計算,以便驗證后續的數值模擬結果,計算得到管道的軸向移動速度為0.059 m/循環,采取的主要計算參數見表1。

表1 管道參數

3 短直管道軸向移動解析計算

3.1有限元模型的建立

在海底管道有限元分析中,由于建立的管道模型一般較長,一個方向的尺寸(長度)遠大于另外兩個尺寸,并且以縱向應力為主,所以選擇PIP31梁單元模擬海底管道,具體參數見表1。管道材料采用Ramberg-Osgood本構關系來模擬管道的應力應變關系,其形式為:

(6)

式中:n為材料的硬化系數,σ0.7為斜率為,0.7E的直線與材料應力應變曲線交點所對應的應力。

邊界條件的模擬過程中,將海床模擬成剛性地基,由于海底管道發生軸向膨脹產生位移,管道與土體之間將會產生摩阻力,所克服的摩擦阻力與其軸向位移呈非線性關系,本文采用二折線模型模擬管土軸向相對作用[9]如圖1所示。

圖1 管土相互作用兩折線模型Fig.1 Double line model of pipeline-soil interaction

對于二折線模型中Ubk為管道滑移距離,進行管土側向分析時取0.1倍管徑,軸向分析時取0.01倍管徑或0.005 m。H為摩擦阻力,V為管道淹沒重量,兩者有如下關系:

(7)

在建模過程中,通過添加非線性彈簧單元來模擬管道與海床的摩擦接觸,非線性彈簧單元含有兩個節點,一個節點定義在管道外壁上,另一個節點定義在剛性海床上。通過定義彈簧變形量與彈簧受力之間的非線性關系,實現管道與海床間的非線性摩擦力模擬。

3.2溫度加載曲線

管道在加熱過程中,內外溫差會引起內部流體溫度損失,導致管道兩端溫度不一致,所以管道升溫是一個非線性漸變的過程。與之相反,冷卻過程發生在管道關閉之后,管道在失去溫度荷載之后緩慢恢復到環境溫度,是一個均勻線性降溫的過程。在一個開啟/關閉循環過程中,管道軸向移動現象通常發生在管道開啟(非線性升溫)過程中,并且在管道關閉(降溫)過程中變形不可逆轉,所以非線性的溫度分布是造成管道軸向變形的主要因素。海底管道在多次開啟/關閉循環過程中軸向變形產生積累,所以在每次循環時,管道所處的狀態是不同的。圖2為管道一次開啟過程中的溫度變化曲線[10],每次開啟過程分15步逐步升溫。

圖2 管道加熱溫度曲線Fig.2 Pipeline heating temperature curve

3.3有效軸向力分布

圖3為第一次開啟/關閉循環過程中管道有效軸向力分布圖。STEP 1~15代表15步升溫加載過程,STEP-16代表降溫過程。

管道在加熱過程中膨脹產生軸向力,在第一次循環加載過程中,管道的虛擬錨固點伴隨加熱過程從管道加熱端向中點移動,直到第8個加載步,管道虛擬錨固點位于管道中點,軸向力被完全激發。之后隨著溫度荷載的增加,管道會繼續受熱膨脹,但是管道有效軸向力沿管道分布形態幾乎沒有變化。在第16荷載步,管道降溫到環境溫度,管道收縮,導致摩擦力反向,管道的有效軸向力變為拉力,虛擬錨固點仍位于管道中部。

第二次和第三次開啟/關閉循環過程中管道的有效軸向力分布如圖4、圖5所示,由于第一次降溫過程中軸向摩擦產生管道殘余軸向拉力,第二次循環及之后的循環的有效軸向力沿管長的分布與第一次循環是不同的。管道進行第一步加載后,產生兩個虛擬錨固點,分別位于26 m和1 004 m處。隨著管道持續升溫,虛擬錨固點的位置分別向管道中點(1 000 m)和非加熱端(2 000 m)移動。直至第9個加熱步,管道膨脹被完全激發,虛擬錨固點位于管道中點,有效軸向力達到最大值(551 kN)并且分布趨于穩定。

圖3 第一次循環加載有效軸向力分布圖Fig.3 Effective axial force distribution graph under the first cyclic loading

圖4 第二次循環加載有效軸向力分布圖Fig.4 Effective axial force distribution graph under the second cyclic loading

圖5 第三次循環加載有效軸向力分布圖Fig.5 Effective axial force distribution graph under the third cyclic loading

3.4管道軸向位移分布

圖6~8所示為管道在前三次開啟/關閉過程中管道膨脹的累計位移分布圖,隨著溫度的逐漸上升,管道會發生非均勻的膨脹作用,并傾向于由管道中部向兩端膨脹。并且,在管道完全激發之后,管道非加載端的位移隨著溫度荷載的升高而急速增加。

圖6 第一次循環加載管道位移分布圖Fig.6 Pipeline displacement distribution under the first cyclic loading

圖7 第二次循環加載管道位移分布圖Fig.7 Pipeline displacement distribution under the second cyclic loading

圖8 第三次循環加載管道位移分布圖Fig.8 Pipeline displacement distribution under the third cyclic loading

由前三次開啟/關閉循環位移對比可看出,由于管道溫度的不均勻升高,管道的整體趨勢是由加熱端向非加熱端移動。

3.5管道軸向移動現象

3.5.1 管道中點軸向移動

由前三次開啟/關閉循環過程中管道中點的軸向位移數據,繪制出管道中點軸向移動趨勢圖,見圖9。由圖9可見,隨著溫度荷載的升高,管道中點位移值逐步增加,并且是由管道加載端向非加載端移動。但是在管道被完全激發之后,由于管道非加熱端溫度的大幅度升高,管道中點位移會短暫下降,之后趨于穩定。管道中點軸向移動在第二次管道循環加熱過程中開始發生,并且隨著循環加載的進行,中點的軸向位移也逐漸增加。管道的軸向移動速度約為0.06 m/循環,與式(7)數值計算結果相似。

圖9 管道中點軸向移動趨勢圖Fig.9 Mid-point axial movement trend of the pipeline

3.5.2 管道端點軸向移動

圖10為管道加熱端在前三次升溫降溫循環過程中的軸向移動趨勢圖。如圖10所示,在首次循環荷載作用下,管道加熱端的位移隨著溫度荷載的升高而快速增加,但在管道被完全激發之后,管道軸向位移值增加較緩,直至管道降溫冷卻收縮,管道位移迅速下降,但會存在部分殘余變形。在隨后的循環過程中,會繼承之前管道的殘余變形,重復之前的變形過程,進而發生管道加熱端的軸向移動現象。

圖10 管道加熱端軸向移動趨勢圖Fig.10 Axial movement of heated end of pipeline

圖11 管道非加熱端軸向移動趨勢圖Fig.11 Axial movement of non-heated end of pipeline

圖11為管道非加熱端在前三次升溫降溫循環過程中的軸向移動趨勢圖。

如圖11所示,在首次循環荷載作用下,管道被完全激發之前,管道非加熱端的位移隨著溫度荷載的升高而增加,但位移變化很小。直至管道被完全激發,非加載端的軸向位移迅速增加,且由加載端向非加載端移動。在管道冷卻過程中,端部位移會緩慢下降,但會有部分殘余變形存在。在隨后的循環過程中,會繼承之前管道的殘余變形,重新發生位移循環過程,進而發生管道非加熱端的軸向移動現象。

圖12為管道軸向移動趨勢圖。如圖所示,在管道每次升溫降溫循環過程之后,管道都位于一個新的位置。這意味著,管道的軸向移動現象是多次循環荷載作用下軸向位移累積的結果。

圖12 管道軸向移動趨勢圖Fig.12 Axial movement trend of the pipeline

綜上所述,在非線性循環溫度荷載作用下,管道有效軸向力將產生循環變化,從而發生管道軸向移動現象,且其方向為由管道加熱端向非加熱端移動。

4 管土摩擦對管道軸向移動的影響

4.1有效軸向力分布

采用前文所述二折線管土相互作用模型,圖13~16為不同的管土摩擦系數下,管道的有效軸向力變化趨勢,升溫過程中管道中存在軸向壓力,而在降溫冷卻過程中管道存在軸向拉力。如圖所示,隨著管道與海床摩擦力的增加,管道的有效軸向力逐漸增大。

圖13 有效軸力分布圖(μ=0.3)Fig.13 Distribution of effective axial force(μ=0.3)

當管土摩擦系數較小時,管道能夠被完全激發,且在完全激發之后,管道的虛擬錨固點都位于管道中部,如圖13、14。

圖14 有效軸力分布圖(μ=0.9)Fig.14 Distribution of effective axial force(μ=0.9)

圖15 有效軸力分布圖(μ=1.2)Fig.15 Distribution of effective axial force(μ=1.2)

圖16 有效軸力分布圖(μ=1.8)Fig.16 Distribution of effective axial force(μ=1.8)

當管土摩擦系數增加到某一值時,此時由于管道升溫產生的軸向力不足以抵消管土摩擦力,管道會在中間段產生虛擬錨固區域,并隨著摩擦系數的增加,錨固區域擴大。在管道的虛擬錨固區域內管道沒有被完全激發,這意味著在虛擬錨固區內,管土相互作用是發生在管土作用模型的彈性區域內的。因此,管道的有效軸向力不會從最大的壓縮力變為最大的拉伸力。

4.2管道軸向移動速度

根據不同的管土摩擦系數,建立模型計算管道前三次荷載循環作用下的軸向移動現象。得出不同摩擦系數下管道中點軸向位移值,如圖17所示。當管土摩擦系數較低時,隨著管土摩擦系數的提高,管道中點處的位移也增加。但隨著摩擦系數的繼續提高,管道的中點位移反而下降。

(a) 摩擦系數較小

(b)摩擦系數較大圖17 不同摩擦系數下管道中點軸向位移Fig. 17 Displacement of pipeline mid-point with different friction coefficient

圖18顯示了不同管土摩擦系數對管道軸向移動速率的影響。當管土摩擦系數較小時,管道再升溫過程中會被完全激活,且隨著管土摩擦系數的增大,管道的軸向移動速度也隨之增加。但隨著管土摩擦系數繼續增大,管道的虛擬錨固區域逐漸增加,同時管道軸向移動速度開始下降。

圖18 管道軸向移動速度隨管土摩擦系數變化曲線Fig.18 Relationship between the pipeline-soil friction coefficient and velocity of axial movement of pipeline

綜上所述,管土之間軸向摩擦力的提高會使管道在中間段產生虛擬錨固區域,降低管道軸向移動速率。

5 管道內壓對軸向移動的影響

考慮恒定內壓與循環內壓兩種方式下管道軸向移動。恒定內壓加載方式為海底管道所受的壓力荷載在整個溫度循環過程種是恒定的,并在管道降溫過程中下降為零;循環內壓加載方式為在管道每一次溫度循環過程中內壓荷載也循環,并在降溫荷載下,管道內壓降低為零。

圖19 恒定內壓下管道有效軸力分布圖Fig.19 Distribution of effective axial force under constant inner pressure

圖20 管道中點軸向位移Fig.20 Distribution of axial displacement of pipeline mid-point

圖19為恒定內壓條件下管道有效軸力分布,管道加熱共分為15個加載步,第16步為降溫加載步。由圖中可以看出,在第一個加熱步為管道施加內壓荷載,在管道全長范圍內有效軸向力增加,但在管道兩端逐漸釋放。此后加載步中有效軸向力的分布形態與無內壓管道相同,但管道所達到的最大有效軸向力相較無內壓管道有所增加。在降溫加載步中,管道溫度降低內壓下降,管道收縮,摩擦力反向,管道有效軸力轉變為拉力;循環內壓條件下管道有效軸力分布與圖3相似,管道加熱共分為15個加載步,每個荷載步中管道內壓逐漸上升,第16步為降溫加載步,此加載步中管道內壓為零。在管道的每個升溫加載步中,由于管道內壓逐漸上升,管道內的有效軸向力也均勻上升。到第8個加載步時,管道被完全激發,虛擬錨固點位于管道中點。在此之后,雖然管道的虛擬錨固點不再變化,但管道有效軸向力還存在小幅升高。在降溫加載步中,管道溫度和內壓逐漸下降,管道收縮,摩擦力反向,管道有效軸力轉變為拉力。

圖20所示為設計內壓為10 MPa管道和無內壓管道在循環溫度荷載的作用下管道中點的軸向位移值。由圖可知,在管道內壓為10 MPa時,管道軸向移動速率約為0.081 m/循環,相較于無壓力管道0.060 m/循環,增長了26%。由此可見,管道內壓會加快管道的軸向移動速率。

由圖20可知,內壓循環管道的軸向位移值明顯低于恒定內壓管道。恒定內壓管道和循環內壓管道的中點軸向移動速率對比如表2所示。由表可知,與恒定內壓管道相比,循環內壓管道對應的中點軸向移動速率降低了14.8%。由此可見,管道內壓循環有益于降低管道的軸向移動速率。

表2恒定內壓管道與循環內壓管道對比

Table 2Comparison between pipeline under constant inner pressure and cyclic inner pressure

恒定內壓管道循環內壓管道軸向位移速率0.081m/循環0.069m/循環

綜上所述,管道內壓會提高管道在循環溫度荷載下的軸向移動速率,但內壓加載方式,管道軸向移動速率也不相同,相比于恒定內壓的加載方式,循環內壓方式有益于降低管道的軸向移動速率。

6 結論

本文建立了循環溫度荷載下平坦海床上直管道的非線性有限元模型,模型中考慮了非線性的溫度加載梯度和高溫對管道材料的影響。分析了短直管道軸向運動機理,考慮了管道內壓加載方式對管道軸向移動現象的影響。通過研究得到如下結論:

1)短直管道在非線性循環溫度荷載作用下更容易發生軸向移動現象,這是因為管道長度不足以使其達到完全約束狀態,管道虛擬錨固點位于管道中點處。

2)管道的非線性溫度循環使其產生由管道加熱端向非加熱端的軸向移動現象,因此管道端部的連接件可靠性尤為重要。

3)管土摩擦對管道軸向移動速率有很大影響,當管土摩擦系數大于某一范圍后,管道軸向移動速率隨著摩擦系數的增加而下降。

4)管道內壓會提高管道在非線性循環溫度荷載下的軸向移動速率,但內壓加載方式,管道軸向移動速率也不相同,相比于恒定內壓的加載方式,循環內壓方式有益于降低管道的軸向移動速率。

[1]TORNES K, JURY J, OSE B. Axial creeping of high temperature flowlines caused by soil ratcheting[C]//Proceedings of the Conference on Offshore Mechanics and Arctic Engineering. ASME, (s.l.),2000.

[2]CARR M, SINCLAIR F, BRUTON D. Pipeline walking-understanding the field layout challenges, and analytical solutions developed for the SAFEBUCK JIP[C]//Offshore Technology Conference, 2006. Houston, Texas: OTC, 2006.

[3]BRUTON D A S, Bolton M, CARR M, et al. Pipe-soil interaction with flowlines during lateral buckling and pipeline walking-The SAFEBUCK JIP[C]//Offshore Technology Conference, 2008. Houston, Texas: OTC, 2008.

[4]ZHOU Jing, LIU Yuxiao, LI Xin. Pipe walking-lateral buckling interaction[C]//Earth and Space 2010: Engineering, Science, Construction, and Operations in Challenging Environments. Honolulu, Hawaii, United States: American Society of Civil Engineers, 2010: 3318-3327.

[5]RONG Haicheng, INGLIS R, BELL G, et al. Evaluation and mitigation of axial walking with a focus on deep water flowlines[C]//Offshore Technology Conference 2009. Houston, Texas: OTC, 2009.

[6]王哲, 馬克儉, 陳志華, 等. 深海管道整體屈曲研究綜述[J]. 天津大學學報: 自然科學與工程技術版, 2014, 47(S): 17-23.

WANG Zhe, MA Kejian, CHEN Zhihua, et al. Overview of global buckling of deep-sea pipeline[J]. Journal Tianjin university: science and technology, 2014, 47(S): 17-23.

[7]《海洋石油工程設計指南》編委會. 海洋石油工程安裝設計[M]. 北京: 石油工業出版社, 2007.

[8]TAYLOR N, TRAN V. Experimental and theoretical studies in subsea pipeline buckling[J]. Marine structures, 1996, 9(2): 211-257.

[9]American Petroleum Institute (API). API 2A, Recommended practice for planning, designing and constructing fixed offshore platforms-working stress design[S]. 21st ed. Washington: American Petroleum Institute, 2000.

[10]CARR M, SINCLAIR F, BRUTON D. Pipeline walking-understanding the field layout challenges, and analytical solutions developed for the SAFEBUCK JIP[C]//Offshore Technology Conference 2006.Houston, Texas, 2006.

本文引用格式:

陳志華,黃金超,趙思玥,等. 海底管道軸向移動數值模擬研究[J]. 哈爾濱工程大學學報, 2016, 37(9): 1197-1203.

CHEN Zhihua, HUANG Jinchao, ZHAO Shiyue,et al. Numerical simulation of axial travel of submarine pipelines[J]. Journal of Harbin Engineering University, 2016, 37(9): 1197-1203.

Numerical simulation of axial travel of submarine pipelines

CHEN Zhihua1,2, HUANG Jinchao2, ZHAO Shiyue2, HE Yongyu2, LIU Hongbo1,2

(1.State Key Laboratory of Hydraulic Engineering Simulation and Safety, Tianjin University; 2.Department of Civil Engineering, Tianjin University, Tianjin 300072,China)

Non-linear temperature load is a main factor causing the axial travel of pipelines. To investigate the axial travel of pipelines on the seabed under frequent on/off cyclic action, in this paper, we used ABAQUS software to build a finite element model of a short-straight pipeline located on the flat seabed under a cyclic temperature environment. This model took into consideration the effects of the non-linear temperature application gradient and high temperature on the pipe material, and revealed the distribution of the effective axial force on the short-straight pipe under an on/off cyclic load as well as the formation mechanism and development pattern of the axial travel of the pipe. Lastly, by parametric analysis, we obtained the relationship between the pipe-seabed friction coefficient and the axial moving speed of the pipeline. Also, we compared the effects of two different loading means of internal pressure in the pipeline on the axial travel speed. Our study conclusions can provide assistance in developing methods to reduce the axial travel speed of pipelines.

short straight submarine pipeline; numerical simulation; axial travel;nonlinear cyclic temperature; parametric analysis

2015-07-18.

時間:2016-09-05.

國家重點基礎研究發展計劃(2014CB046801).

陳志華(1966-), 男, 教授,博士生導師;

劉紅波,E-mail:hb_liu2008@163.com.

10.11990/jheu.201507047

TE973

A

1006-7043(2016)09-1197-07

劉紅波(1983-), 男, 副教授,博士.

網絡出版地址:http://www.cnki.net/kcms/detail/23.1390.u.20160905.0910.004.html

主站蜘蛛池模板: 国产精品黑色丝袜的老师| 日本妇乱子伦视频| 91视频区| 激情网址在线观看| 国产丝袜丝视频在线观看| 日韩第一页在线| 久久亚洲综合伊人| 欧美性猛交一区二区三区| 免费国产小视频在线观看| 国产成人久久综合一区| 亚洲成年人片| 91麻豆国产视频| 毛片基地美国正在播放亚洲| 亚洲成aⅴ人片在线影院八| 99热这里都是国产精品| 欧美精品v日韩精品v国产精品| 中文字幕久久精品波多野结| 欧美日韩国产精品va| 91无码视频在线观看| 精品国产免费人成在线观看| 欧美黄网在线| 2022国产无码在线| 天堂岛国av无码免费无禁网站 | 黄色福利在线| 免费不卡视频| a级毛片毛片免费观看久潮| 久久综合丝袜日本网| 一级毛片网| 综合人妻久久一区二区精品| 一级做a爰片久久毛片毛片| 麻豆精品在线视频| 一边摸一边做爽的视频17国产| 日韩a级片视频| 久久精品人人做人人综合试看| 成人午夜网址| 九九九久久国产精品| 亚洲中文字幕日产无码2021| 亚洲伊人久久精品影院| 91在线免费公开视频| 一区二区自拍| 青青草国产精品久久久久| 欧美在线国产| 国产玖玖玖精品视频| 国内精品视频| 美女无遮挡被啪啪到高潮免费| 久久精品日日躁夜夜躁欧美| 国产综合精品一区二区| 国产成人高清精品免费软件| 欧美精品不卡| 亚洲国产中文欧美在线人成大黄瓜 | 久久国产精品无码hdav| 91福利国产成人精品导航| 国产尤物视频在线| www中文字幕在线观看| 无码精品国产VA在线观看DVD| 免费无码AV片在线观看国产| 激情综合婷婷丁香五月尤物| 欧美丝袜高跟鞋一区二区| 欧美日韩精品一区二区视频| 亚洲成a人片| 久久伊人操| 久久久久九九精品影院| 久久精品国产精品一区二区| 亚洲精品男人天堂| 国产欧美日韩综合一区在线播放| 无码'专区第一页| 欧美日韩中文字幕在线| 亚洲无码视频一区二区三区 | 国产黑人在线| 91精品日韩人妻无码久久| 97在线碰| 亚洲天堂视频网站| 日本黄色不卡视频| 成人亚洲天堂| 日韩欧美成人高清在线观看| 日本一本在线视频| 在线观看热码亚洲av每日更新| 国产成人免费| 国产一级妓女av网站| 精品久久蜜桃| 国产原创演绎剧情有字幕的| 国产麻豆va精品视频|