• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    流冰對(duì)引水隧洞撞擊破壞力學(xué)特性數(shù)值分析與驗(yàn)證

    2018-08-10 07:19:04李雅嫻靳春玲
    關(guān)鍵詞:撞擊力隧洞流速

    貢 力,李雅嫻,靳春玲

    ?

    流冰對(duì)引水隧洞撞擊破壞力學(xué)特性數(shù)值分析與驗(yàn)證

    貢 力,李雅嫻,靳春玲

    (蘭州交通大學(xué)土木工程學(xué)院,蘭州 730070)

    寒旱地區(qū)冬季寒冷,晝夜溫差大,長(zhǎng)距離輸水工程在解凍期易形成流冰,不同尺寸的流冰對(duì)隧洞形成不同的擠壓力或撞擊力,甚至導(dǎo)致襯砌破損,工程失效。該文通過(guò)深入研究流冰與引水隧洞碰撞時(shí)的相互作用問(wèn)題,并利用有限元接觸-碰撞算法的對(duì)稱罰函數(shù)算法,進(jìn)行了流冰撞擊引水隧洞襯砌的接觸-碰撞算法的理論分析。應(yīng)用ANSYS/LS-DYNA軟件,建立隧洞模型和流冰模型,選用LS-DYNA SOLVER求解器進(jìn)行求解,分析流冰對(duì)隧洞的破壞程度;并通過(guò)物理模型試驗(yàn)采用幾何比尺為C為28進(jìn)行驗(yàn)證,揭示流冰對(duì)引水隧洞的撞擊破壞機(jī)理。結(jié)果表明:流冰在不同流速、不同平面尺寸、不同厚度等工況下,流冰與隧洞襯砌碰撞時(shí),隧洞襯砌表面會(huì)產(chǎn)生不同程度的變形和破壞;隨著流速增大,撞擊應(yīng)力值也相應(yīng)增大,兩者呈線性關(guān)系;當(dāng)流冰平面尺寸變化時(shí),其撞擊應(yīng)力隨著流冰平面尺寸的增大而增大,兩者呈非線性關(guān)系;當(dāng)流冰厚度增加時(shí),流冰厚度小于0.5 m時(shí),撞擊力隨著流冰厚度的增大而增大;其厚度超過(guò)0.5 m時(shí),撞擊應(yīng)力值變化不大,其流冰平面尺寸和最大應(yīng)力呈現(xiàn)近似線性關(guān)系。同時(shí),通過(guò)軟件模擬和試驗(yàn)觀測(cè)得出的計(jì)算結(jié)果基本一致,流冰與隧洞襯砌碰撞時(shí),隧洞襯砌表面會(huì)產(chǎn)生不同程度變形,變形對(duì)隧洞穩(wěn)定性影響不大。但是,流冰沖刷會(huì)導(dǎo)致隧洞襯砌表面破碎,長(zhǎng)期會(huì)影響結(jié)構(gòu)的強(qiáng)度與穩(wěn)定性。其結(jié)果可為寒旱地區(qū)冬季輸水工程安全提供理論支撐和技術(shù)保障。

    力學(xué)特性;數(shù)值模擬;模型;冰工程;引水隧洞;撞擊力

    0 引 言

    中國(guó)西北地區(qū)位于高寒、干寒地帶,冬季氣溫低,冰期長(zhǎng),為了緩解西部地區(qū)嚴(yán)重缺水的現(xiàn)狀,中國(guó)先后在西部地區(qū)投資建設(shè)了引大入秦工程、引洮工程等大批引調(diào)水工程,用于農(nóng)田灌溉、人畜飲水,目前籌建的南水北調(diào)西線項(xiàng)目更是投資空前[1]。然而如此規(guī)模龐大的投資建設(shè),建成的項(xiàng)目卻因遭受寒旱氣候的影響,產(chǎn)生了諸如引水隧洞凍融破壞、引水隧洞地質(zhì)缺陷滲漏等一系列問(wèn)題,對(duì)結(jié)構(gòu)安全性與供水效率產(chǎn)生嚴(yán)重影響。尤其是高寒地區(qū)冬季寒冷,晝夜溫差大,在氣溫降低到0 ℃以下后,水流內(nèi)首先會(huì)出現(xiàn)水內(nèi)冰,水內(nèi)冰經(jīng)過(guò)發(fā)展逐漸形成冰花,繼而堆積形成冰蓋,至第二年解凍期,由于高寒地區(qū)溫差大,冰蓋融化過(guò)程中,易形成武開(kāi)河或半文半武開(kāi)河,產(chǎn)生的流冰,對(duì)引水隧洞造成不同的擠壓力或撞擊力,導(dǎo)致襯砌破損,直至工程失效,達(dá)不到規(guī)劃建設(shè)的預(yù)期目標(biāo),進(jìn)而造成受水區(qū)的供水分配壓力,影響春季引水灌溉和春耕生產(chǎn)。因此,流冰對(duì)引水隧洞的破壞機(jī)理研究迫在眉睫。

    國(guó)內(nèi)外學(xué)者對(duì)冰期的冰情演變進(jìn)行了研究,國(guó)內(nèi)1921年在黃河下游對(duì)冰凌進(jìn)行了首次科學(xué)觀測(cè)[2];Gilberto等[3]進(jìn)行了彎道冰塞堆積的水槽試驗(yàn),對(duì)冰塞堆積現(xiàn)象的一些特性開(kāi)展了研究;茅澤育等[4-6]建立了基于河流動(dòng)力學(xué)和熱力學(xué)原理的冰水力學(xué)數(shù)學(xué)模型的基本流動(dòng)路線圖和相應(yīng)程序;肖建民等[7]開(kāi)展了冰蓋的形成與消融原理方面的研究;王軍等[8]通過(guò)冰塞堆積試驗(yàn),使得冰塞的形成機(jī)理和規(guī)律方面的研究向前推進(jìn)了一步;徐國(guó)賓等[9]在天津大學(xué)低溫冰工程實(shí)驗(yàn)室進(jìn)行了冰力學(xué)模型試驗(yàn),探索了冰的膨脹力、冰蓋穩(wěn)定性和流冰對(duì)橋墩等水工建筑物的撞擊力等冰力學(xué)問(wèn)題;貢力等[10]對(duì)西部地區(qū)引水工程病害特點(diǎn)進(jìn)行了研究和分析。綜上所述,國(guó)內(nèi)外的研究者對(duì)河流冰情演變規(guī)律研究較多,經(jīng)歷了從早期的原型觀測(cè)、試驗(yàn)研究到現(xiàn)在廣泛應(yīng)用的數(shù)值模擬幾個(gè)階段,但國(guó)內(nèi)外專門針對(duì)流冰對(duì)寒旱地區(qū)長(zhǎng)距離輸水工程中引水隧洞的影響研究較少。

    為了研究寒旱地區(qū)長(zhǎng)距離輸水工程中解凍期流冰對(duì)引水隧洞的影響,本文采用理論分析、數(shù)值模擬和試驗(yàn)研究的方法,開(kāi)展流冰對(duì)引水隧洞的襯砌撞擊影響力的研究。運(yùn)用ANSYS LS-DYNA建立流冰與引水隧洞之間發(fā)生碰撞的有限元模型,搭建流冰撞擊引水隧洞試驗(yàn)單邊坡試驗(yàn)裝置,模擬隧洞中流冰與引水隧洞碰撞的全演變過(guò)程,發(fā)現(xiàn)流冰在不同工況下對(duì)引水隧洞的撞擊力影響規(guī)律,為寒旱地區(qū)冬季輸水工程安全提供理論支撐和技術(shù)保障。

    1 流冰對(duì)引水隧洞撞擊的數(shù)值模擬理論基礎(chǔ)

    一般描述各種非線性物質(zhì)運(yùn)動(dòng)和變形以及碰撞等問(wèn)題的控制方程有Lagrange法、Euler法兩類。本文研究的接觸-碰撞問(wèn)題一般釆用Lagrange描述法。Lagrange法質(zhì)點(diǎn)隨著物體運(yùn)動(dòng)的過(guò)程中,其質(zhì)量保持不變。流冰撞擊隧洞邊壁的過(guò)程中,質(zhì)點(diǎn)在整個(gè)運(yùn)動(dòng)系統(tǒng)中必須保持能量守恒方程、質(zhì)量守恒方程、動(dòng)量守恒方程[11-14]。

    1.1 質(zhì)量守恒方程

    式中0(, 0)為模型=0時(shí)刻的初始模型中的介質(zhì)密度(kg/m3);(,)為模型在時(shí)刻的構(gòu)形介質(zhì)密度(kg/m3);(,)為模型的Jacobin行列式。

    式中e為排列張量,X為物質(zhì)坐標(biāo),x為空間坐標(biāo)。

    1.2 動(dòng)量守恒方程

    式中v(,)為質(zhì)點(diǎn)在時(shí)刻坐標(biāo)為的瞬時(shí)速度(m/s)。

    1.3 能量守恒方程

    式中v為質(zhì)點(diǎn)在方向的瞬時(shí)速度,m/s;v為質(zhì)點(diǎn)在方向的瞬時(shí)速度,m/s。

    1.4 邊界條件

    隧洞模型現(xiàn)時(shí)構(gòu)型的面力邊界條件為:=n,其中為面元上的應(yīng)力矢量,為現(xiàn)時(shí)構(gòu)型有向面元,為現(xiàn)時(shí)構(gòu)型的現(xiàn)時(shí)構(gòu)型的應(yīng)力張量。計(jì)算要求流冰與隧洞之間的整個(gè)運(yùn)動(dòng)過(guò)程當(dāng)中,在求解范圍內(nèi)全部滿足動(dòng)量守恒方程式,但在實(shí)際工程問(wèn)題中幾乎不可能直接求解出結(jié)果。而數(shù)值計(jì)算方法可以從微分方程的弱形式角度考慮,只需在內(nèi)積條件下得到基本動(dòng)量方程,進(jìn)而對(duì)模型的虛功率方程式進(jìn)行推導(dǎo),再經(jīng)過(guò)有限元離散化后,得到模型的節(jié)點(diǎn)位移方程。通過(guò)加權(quán)余量法,選取虛速度作為加權(quán)系數(shù),則動(dòng)量方程的弱形式見(jiàn)式(7)。

    式中v()為虛速度,其值根據(jù)下式(8)進(jìn)行求解。

    式中0,0分別為0時(shí)刻的向量和阻尼矩陣,A為模型的邊界面,運(yùn)用分部積分原理,接觸面力的平衡式如式(9)所示。

    那么,模型中任意質(zhì)點(diǎn)的速度、加速度、虛速度及變形率可以寫成方程組

    把方程組(10)中的各方程寫成矩陣形式,然后代入虛功率公式(9)中。整理后可得

    求解式(14),便可得出當(dāng)下時(shí)刻的質(zhì)點(diǎn)位移,進(jìn)一步解得該時(shí)刻的結(jié)構(gòu)應(yīng)變與應(yīng)力。

    2 流冰對(duì)隧洞撞擊模型的建立

    2.1 工程實(shí)例

    盤道嶺3號(hào)隧洞全長(zhǎng)15.723 km,縱坡1/1 000,設(shè)計(jì)流量為29 m3/s,最大流量為34 m3/s。該隧道為無(wú)壓引水隧洞,洞身為混合式襯砌型式,一次襯砌為錨桿、噴射混凝土、鋼筋網(wǎng)片和鋼拱架,二次襯砌為現(xiàn)澆混凝土和鋼筋混凝土,圓拱直墻、底板為反拱的斷面,凈寬4.2 m,凈高4.4 m,頂為半圓,半徑為2.1 m,反拱底板半徑9.75 m,側(cè)墻與反拱交接處加設(shè)貼角,貼角高0.404 m,水平寬0.337 m。隧洞前期支護(hù)隧洞斷面圖及模型圖如圖1所示。

    注:C20為抗凍強(qiáng)度,F(xiàn)150為抗?jié)B強(qiáng)度。

    2.2 模型建立

    ANSYS LS-DYNA具有ALE和Euler算法、熱分析和流體-固體耦合分析功能、靜力分析功能和隱式分析功能等,可以快速地求解平面或空間內(nèi)高速碰撞、爆炸等動(dòng)態(tài)非線性問(wèn)題。采用ANSYS/LS-DYNA模擬流冰對(duì)隧洞的撞擊。因此,流冰位于水面之上,所以在ANSYS/LS-DYNA軟件模擬模型的荷載時(shí),應(yīng)忽略豎向荷載,如重力、浮力等,只考慮風(fēng)、水流拖曳力等水平方向的荷載[15-16]。水平荷載主要通過(guò)流冰的初速度表達(dá),體積不等、厚度不一的流冰分別以不同的速度作用于隧洞內(nèi)壁時(shí),各種不同工況下的撞擊力也是不同的,通過(guò)建模來(lái)模擬計(jì)算流冰對(duì)隧洞的撞擊力。

    本文選取單面自動(dòng)接觸(automatic single surface contact, ASSC)接觸類型。假設(shè)初速度,定義接觸、邊界條件后進(jìn)行輸出控制,利用動(dòng)力學(xué)分析命令流文件(ANSYS軟件中稱為K文件)進(jìn)行輸入計(jì)算。在計(jì)算中,流冰和隧洞模型采用3D Solid164實(shí)體單元模擬;流冰和隧洞材料為線彈性材料;單元?jiǎng)澐植捎糜成渚W(wǎng)格劃分。隧洞襯砌面為主體,流冰接觸面為從界面。參數(shù)見(jiàn)表1,其中冰的參數(shù)根據(jù)于天來(lái)等[17-18]對(duì)冰的彈性模量與冰的關(guān)系確定。

    表1 隧洞和流冰物理性能參數(shù)

    3 流冰對(duì)隧洞撞擊模型的應(yīng)用

    由于流冰與隧洞在碰撞過(guò)程中,牽扯諸如動(dòng)態(tài)、流固耦合等不穩(wěn)定因素較多,是一個(gè)很復(fù)雜的冰結(jié)構(gòu)之間相互碰撞的問(wèn)題,影響因素主要有特征冰體的速度、平面尺寸、厚度、幾何形狀、抗拉強(qiáng)度、抗壓強(qiáng)度、彈性模量、結(jié)構(gòu)的柔性、摩擦力、碰撞角度等[19-21]。本文主要研究流冰在不同流速、不同流冰平面尺寸、不同流冰厚度等工況下流冰對(duì)引水隧洞邊壁的撞擊力。方便計(jì)算,計(jì)算中碰撞角度采用正碰,幾何形狀采用長(zhǎng)方形,不考慮二次碰撞等問(wèn)題,其隧洞流冰組合模型網(wǎng)格劃分如圖2所示。

    圖2 隧洞流冰組合模型網(wǎng)格劃分圖

    3.1 流速對(duì)撞擊力的影響

    數(shù)值模擬中選取水深為2 m,實(shí)際碰撞過(guò)程中流冰的尺寸不規(guī)則,參照徐國(guó)賓等[9]的試驗(yàn),模擬流冰尺寸設(shè)置為2 m×0.5 m×2 m,根據(jù)引大入秦工程運(yùn)營(yíng)實(shí)際水流的控制流速,流速分別取0.5、0.8、1.0、1.3、1.5、1.8、2.0、2.3、2.5、2.8、3.0、3.5 m/s。在流冰允許破壞范圍內(nèi),流冰在=0.5 m/s時(shí)撞擊力最大應(yīng)力云圖如圖3所示。

    注:時(shí)間t=0.004 984 8 s,接觸面取1.0 m2。v為流速,下同。

    圖3為流冰對(duì)隧洞的撞擊力最大應(yīng)力云圖,由于撞擊過(guò)程中撞擊力為動(dòng)態(tài)變化,當(dāng)=0.004 984 8 s時(shí),可以得到撞擊力瞬間最大值為0.53×103kN。同理,可以得到不同流速相應(yīng)的撞擊力最大應(yīng)力云圖。通過(guò)云圖可知,最大撞擊力分別0.94、1.18、1.47、1.77、2.12、2.36、2.71、2.93、3.27、3.42、4.09×103kN??梢缘玫阶畲笈鲎擦﹄S著流冰流速的變化關(guān)系曲線,如圖4所示。

    圖4 最大撞擊力-流速關(guān)系圖

    由圖4可以看出,在其他工況不變的情況下,只改變流冰速度,流速對(duì)撞擊應(yīng)力極值有明顯的影響,其撞擊力隨著流冰流速的增大而增大,其變化關(guān)系呈現(xiàn)出近似線性的關(guān)系。由于該模擬過(guò)程中接觸面積為1.0 m2,可以發(fā)現(xiàn),流冰流速與最大撞擊力的關(guān)系近似如公式(12)所示。

    =1.2(0.5m/s≤≤3.5 m/s) (12)

    式中為流冰的流速,m/s;為流冰與隧洞邊壁碰撞時(shí)所產(chǎn)生的最大撞擊力值,kN。

    由動(dòng)能定理和動(dòng)量定理可知,當(dāng)同一質(zhì)量物體的速度越大,它所具有的動(dòng)能就越大,產(chǎn)生的撞擊力就越大。但是,冰與普通物體(鋼、鐵、石頭等)不同,冰在一定條件下由于其韌、脆轉(zhuǎn)變典型特性的存在,超過(guò)一定大小的撞擊力作用下,冰會(huì)發(fā)生破碎等現(xiàn)象,當(dāng)冰破碎后,流冰尺寸變小,但仍然形成二次碰撞,由于涉及因素復(fù)雜,且影響小,本文不作計(jì)算。

    3.2 流冰平面尺寸對(duì)撞擊力的影響

    為了研究流冰體積對(duì)隧洞襯砌作用力的影響,在不改變接觸面的前提下,僅改變流冰平面尺寸。水深選取2 m,流冰流速選取=5 m/s,冰厚選取0.5 m,在其他參數(shù)固定不變的條件下,只改變流冰平面尺寸,分別取0.5 m×0.5 m×2 m、1.0 m×0.5 m×2 m、1.5 m×0.5 m× 2 m、2.0 m×0.5 m×2 m、2.5 m×0.5 m×2 m。平面尺寸分別為0.5 m×0.5 m的1、2、3、4、5倍。其實(shí)質(zhì)是碰撞接觸面積不發(fā)生變化,為1.0 m2,而流冰的體積和質(zhì)量發(fā)生變化。在流冰允許破壞范圍內(nèi),圖5為流冰尺寸為1.0 m×0.5 m×2 m對(duì)隧洞的撞擊力最大應(yīng)力云圖。

    Note: t= 0.002 466 2 s, v=5 m×s-1

    從圖5可以看出,由于撞擊過(guò)程中撞擊力為動(dòng)態(tài)變化,當(dāng)=0.002 466 2 s時(shí),可以得到撞擊力瞬間最大值為4.34×103kN。同樣的方法得到0.5 m×0.5 m×2 m、1.5 m×0.5 m×2 m、2.0 m×0.5 m×2 m、2.5 m×0.5 m×2 m的撞擊力最大應(yīng)力云圖。對(duì)應(yīng)的撞擊力分別1.10×103、4.99×103、5.69×103、6.23×103kN。由接觸面為1.0 m2,可以得到不同流速下的最大撞擊應(yīng)力分別為1.10、4.99、5.69、6.23 MPa。即可以得到在不同面積尺寸下,最大撞擊力與流冰平面尺寸如圖6所示。

    由圖6可以看出,在不改變流冰其他參數(shù)的情況下,只改變流冰平面尺寸時(shí),撞擊應(yīng)力極值的變化明顯,說(shuō)明質(zhì)量和體積的變化對(duì)撞擊力有較大影響。由于在不同工況下接觸面積均取1.0 m2,所以撞擊應(yīng)力關(guān)系曲線與撞擊力關(guān)系曲線一致。當(dāng)流冰平面尺寸較小時(shí),其撞擊應(yīng)力隨著流冰平面尺寸的增大而增大較為明顯,但是當(dāng)質(zhì)量增大到一定值時(shí),撞擊力增加緩慢,其變化關(guān)系呈非線性的關(guān)系(流冰在發(fā)生韌、脆轉(zhuǎn)變前)。

    3.3 流冰厚度對(duì)撞擊力的影響

    為了研究流冰厚度對(duì)隧洞襯砌作用力的影響,水深選取2 m,選取引大入秦工程運(yùn)營(yíng)中的較大流速,取流冰流速=3 m/s,流冰平面尺寸選取2.0 m×2.0 m,在不改變其他參數(shù)的情況下,只改變流冰厚度,厚度分別取0.3、0.4、0.5、0.6、0.7、0.8、0.9、1.0、1.1、1.2、1.3、1.5 m共12種工況。在流冰允許的破壞范圍內(nèi),得到流冰厚度為0.3 m對(duì)隧洞的撞擊力最大應(yīng)力云圖,如圖7所示。

    Note: t=0.007 451 s. v=3 m×s-1

    由圖7可以看出,由于撞擊過(guò)程中撞擊力是動(dòng)態(tài)變化的,當(dāng)=0.007 45 s時(shí),可以得到撞擊力瞬間最大值為1.58×103kN。同理,可以得到厚度為0.3、0.4、0.5、0.6、0.7、0.8、0.9、1.0、1.1、1.2、1.3、1.5m的撞擊力最大應(yīng)力云圖和局部應(yīng)力云圖。通過(guò)云圖可知,撞擊力分別1.58、2.26、3.42、4.02、5.10、5.49、6.34、7.08、7.37、8.40、9.15、10.47×103kN。當(dāng)接觸面取0.6、0.8、1.0、1.2、1.4、1.6、1.8、2.0、2.2、2.4、2.6、3.0 m2時(shí),得到不同流速下的最大撞擊應(yīng)力分別為2.63、2.82、3.42、3.35、3.64、3.43、3.52、3.54、3.35、3.50、3.52、3.49 MPa。隨之得到在不同厚度下,最大應(yīng)力值與流冰厚度的關(guān)系、最大撞擊力與流冰厚度關(guān)系如圖8所示。

    圖8 流冰厚度與最大應(yīng)力及最大撞擊力關(guān)系圖

    由圖8可知,當(dāng)流冰厚度小于0.5 m時(shí),撞擊應(yīng)力隨著流冰厚度的增大而增大;其厚度超過(guò)0.5 m時(shí),撞擊應(yīng)力值變化不大,最大應(yīng)力呈現(xiàn)線性關(guān)系。撞擊力與流冰厚度之間,當(dāng)流冰厚度較小時(shí),其撞擊力隨著冰厚度的增大而增大,其變化關(guān)系呈現(xiàn)出線性的關(guān)系。

    綜上結(jié)果看出:在不同流冰流速、不同流冰平面尺寸、不同流冰厚度等工況下,流冰與隧洞襯砌碰撞時(shí),隧洞襯砌表面會(huì)產(chǎn)生不同程度的變形。這種變形對(duì)隧洞穩(wěn)定性影響不大,但是,隧洞襯砌會(huì)產(chǎn)生表面破碎,經(jīng)過(guò)長(zhǎng)時(shí)間水流的沖刷,將會(huì)造成隧洞襯砌表面脫落等現(xiàn)象,破壞結(jié)構(gòu)的強(qiáng)度和穩(wěn)定性。

    4 流冰對(duì)隧洞撞擊的物理模型試驗(yàn)

    4.1 試驗(yàn)設(shè)計(jì)

    本文利用室內(nèi)模型試驗(yàn)研究解凍期流冰對(duì)引水隧洞的碰撞作用。對(duì)不同流冰平面尺寸、不同流冰厚度、流冰速度等工況下流冰對(duì)引水隧洞的撞擊力展開(kāi)系統(tǒng)的室內(nèi)模型試驗(yàn)[22-26]。本文中試驗(yàn)裝置與原型的幾何比尺C取為28,試驗(yàn)中液體采用普通水, 材料密度比尺取為1.0,在常重力場(chǎng)條件下進(jìn)行試驗(yàn), 加速度比尺C取為1.0[27-30]。模型試驗(yàn)裝置圖如圖9所示,其中a,b,c,d,e點(diǎn)為應(yīng)變片一側(cè)的位置,另一側(cè)對(duì)稱布置,應(yīng)變片位置高度為17 cm。

    1.循環(huán)水管 2.固定支柱 3.水泵 4.水箱 5.測(cè)試儀器 6.應(yīng)變片(a, b, c, d, e) 7.水槽 8.測(cè)試線 9.旋轉(zhuǎn)螺旋

    4.2 試驗(yàn)內(nèi)容

    根據(jù)試驗(yàn)設(shè)計(jì),其具體試驗(yàn)步驟如下:

    1)研究流速對(duì)撞擊力的影響試驗(yàn)。在水箱中充滿水;打開(kāi)總電源,利用泵使得水箱中的水沖入水槽;將準(zhǔn)備好的若干規(guī)格為7 cm×3.5 cm×7 cm冰塊放入水槽中;利用旋轉(zhuǎn)螺旋來(lái)調(diào)節(jié)流速;用流速計(jì)測(cè)得流速為0.5 m/s;利用水位探針確定水深,根據(jù)幾何比尺,水深為17.2 cm;在隧洞模型兩邊邊壁水面線位置分別貼10個(gè)有機(jī)玻璃應(yīng)變片;使瞬態(tài)應(yīng)變測(cè)試儀處于關(guān)閉狀態(tài);利用測(cè)試線將應(yīng)變片與瞬態(tài)應(yīng)變測(cè)試儀鏈接,每個(gè)應(yīng)變片與瞬態(tài)應(yīng)變測(cè)試儀上面的接頭對(duì)應(yīng)相連;打開(kāi)瞬態(tài)應(yīng)變測(cè)試儀,調(diào)成手動(dòng)裝置;隨著冰塊碰撞隧洞邊壁,依次在瞬態(tài)應(yīng)變測(cè)試儀上度數(shù),讀出應(yīng)變值;記錄數(shù)據(jù),處理數(shù)據(jù);接著改變流速大小使得流速分別為0.5、0.8、1.0、1.2、1.5、1.8、2.0、2.3、2.5、2.8、3.0、3.5 m/s,重復(fù)上述步驟,觀察流速對(duì)撞擊應(yīng)力的影響。

    2)流冰平面尺寸大小對(duì)撞擊力的影響試驗(yàn)。改變工況,水深為17.2 cm;流速為3.0 m/s,改變冰塊尺寸分別取1.8 cm×1.8 cm×7 cm、3.5 cm×1.8 cm×7 cm、5.5 cm× 1.8 cm×7 cm、7.0 cm×1.8 cm×7 cm、9.0 cm×1.8 cm×7 cm;重復(fù)上述步驟,測(cè)試流冰平面尺寸對(duì)撞擊應(yīng)力的影響。

    3)流冰厚度對(duì)撞擊力的影響試驗(yàn)。流速和冰塊厚度一定,冰塊平面尺寸為變量進(jìn)行試驗(yàn),具體為:選定水深為17.2 cm,流速為3.0 m/s,流速和平面尺寸一定,改變冰塊尺寸分別取1.1、1.4、1.8、2.14、2.5、2.9、3.2、3.6、3.9、4.3、4.6、5.4 cm;重復(fù)上述步驟,測(cè)試流冰平面尺寸對(duì)撞擊應(yīng)力的影響。

    4.3 流速對(duì)撞擊力的影響

    按照試驗(yàn)步驟1)開(kāi)展試驗(yàn),得到不同流速下,隧洞邊壁所受的撞擊力及應(yīng)力值,將上述試驗(yàn)數(shù)據(jù)處理后,按照模型相似比進(jìn)行計(jì)算,其軟件模擬計(jì)算結(jié)果、試驗(yàn)觀測(cè)計(jì)算結(jié)果對(duì)比圖如圖10 a所示??梢园l(fā)現(xiàn),在流速不同的條件下,數(shù)值模擬值與試驗(yàn)值有大致相同的趨勢(shì),吻合性較好,總體來(lái)說(shuō),隨著流速的增大,流冰對(duì)隧洞的撞擊力力增大,且基本程線性分布關(guān)系。模擬值與試驗(yàn)值的最大誤差為1.83%,最小誤差為0.04%,平均相對(duì)誤差為0.94%,在允許范圍內(nèi)[28]。由此可知,本文采用的計(jì)算模型準(zhǔn)確可靠,模擬結(jié)果可信。分析可知,誤差主要是試驗(yàn)值波動(dòng)導(dǎo)致,可能與試驗(yàn)中流冰的二次碰撞有關(guān)。

    圖10 試驗(yàn)結(jié)果與模型計(jì)算結(jié)果對(duì)比圖

    4.4 流冰平面尺寸大小對(duì)撞擊力的影響

    按照試驗(yàn)步驟2)開(kāi)展試驗(yàn),得到不同平面尺寸下,隧洞邊壁所受的撞擊力及應(yīng)力值,將上述試驗(yàn)數(shù)據(jù)處理后,按照模型相似比進(jìn)行計(jì)算,可以得到試驗(yàn)結(jié)果與模型計(jì)算結(jié)果對(duì)比圖10b,觀察圖可以發(fā)現(xiàn),當(dāng)流冰平面尺寸較小時(shí),其撞擊應(yīng)力隨著流冰平面尺寸的增大而明顯增大。但是當(dāng)平面尺寸大于1.2 m2以后,撞擊力的增量趨緩,說(shuō)明當(dāng)面積增大時(shí),冰水粘滯力和拖曳力增加,使得撞擊力變緩。通過(guò)計(jì)算可知,試驗(yàn)值與模型值最大誤差為2.24%,最小誤差為0.05%,平均相對(duì)誤差為1.15%,在允許范圍內(nèi)[29],吻合性較好。

    4.5 流冰厚度對(duì)撞擊力的影響

    按照試驗(yàn)步驟3)開(kāi)展試驗(yàn),根據(jù)試驗(yàn)的測(cè)算,流冰的厚度不同,隧洞邊壁所受的撞擊力及應(yīng)力值也不同,將上述試驗(yàn)數(shù)據(jù)處理后,按照模型相似比進(jìn)行計(jì)算,可以得到試驗(yàn)結(jié)果與模型計(jì)算結(jié)果對(duì)比圖10c,觀察圖可以發(fā)現(xiàn),在流冰厚度不同的條件下,模擬值與試驗(yàn)值有大致相同的趨勢(shì),吻合性較好,總體來(lái)說(shuō),隨著流冰厚度的增大,流冰對(duì)隧洞的撞擊力增大,且基本程線性分布關(guān)系。模擬值和試驗(yàn)值結(jié)果吻合性較好,最大誤差為1.98%,最小誤差為0.06%,平均相對(duì)誤差為1.02%,在允許范圍內(nèi)[30]。由此可知,本文采用的計(jì)算模型準(zhǔn)確可靠,模擬結(jié)果可信。

    5 結(jié) 論

    本研究應(yīng)用有限元接觸-碰撞算法建立流冰撞擊引水隧洞襯砌作用時(shí)的數(shù)值計(jì)算仿真模型,并相應(yīng)地開(kāi)展室內(nèi)模型試驗(yàn),數(shù)值模擬結(jié)果和試驗(yàn)結(jié)果良好吻合,并根據(jù)流冰的流速、平面尺寸、厚度等參數(shù)變化對(duì)引水隧洞破壞力學(xué)特性得到以下結(jié)論:

    1)隨著流冰流速的增大,流冰對(duì)隧洞襯砌的碰撞力增大,且流冰流速對(duì)隧洞最大撞擊力之間存在線性關(guān)系(=1.2),試驗(yàn)值與計(jì)算誤差在5%以內(nèi)。因此,在滿足供水要求的情況下,解凍期運(yùn)營(yíng)要控制流速,減小流冰對(duì)隧洞襯砌的撞擊力。

    2)流冰在平面尺寸增大時(shí),其實(shí)質(zhì)是體積和質(zhì)量增大,當(dāng)流冰平面尺寸小于1.0 m2時(shí),其撞擊應(yīng)力隨著流冰平面尺寸的增大而明顯增大,但是當(dāng)流冰平面尺寸超過(guò)1.0 m2時(shí),撞擊力增加緩慢,二者之間呈非線性的關(guān)系。說(shuō)明流冰質(zhì)量增加后,水流的挾冰能力下降,撞擊力變緩。

    3)流冰厚度小于0.5 m時(shí),撞擊力隨著流冰厚度的增大而增大;其厚度超過(guò)0.5m時(shí),撞擊應(yīng)力值變化不大,說(shuō)明冰水粘滯力和拖曳力增加,使得撞擊力變緩。

    實(shí)際工程中冰水兩相之間的運(yùn)動(dòng)是非常復(fù)雜的,其間存在著粘滯力、拖曳力等多因素耦合影響的問(wèn)題;流冰與隧洞襯砌之間的碰撞是不規(guī)則的,在模擬流冰與隧洞的碰撞時(shí),碰撞角度選取的形式等問(wèn)題,需要做進(jìn)一步的研究。

    [1] 楊開(kāi)林. 長(zhǎng)距離輸水水力控制的研究進(jìn)展與前沿科學(xué)問(wèn)題[J]. 水利學(xué)報(bào),2016,47(3):424-435.

    Yang Kailin. Review and frontier scientific issues of hydraulic control for long distance water diversion[J]. Journal of Hydraulic Engineering, 2016, 47(3): 424-435. (in Chinese with English abstract)

    [2] 張成,王開(kāi). 冰期輸水研究進(jìn)展[J]. 南水北調(diào)與水利科技,2006,4(6):59-63.

    Zhang Cheng, Wang Kai. The relational research and experience of water transfer during freezing period[J]. South-to-North Water Transfers and Water Science& Technology, 2006, 4(6): 59-63. (in Chinese with English abstract)

    [3] Gilberto E U, Robert E. Bend ice jams: Laboratory observations[J]. Canadian Journal of Civil Engineering, 1992, 19: 855-864.

    [4] 茅澤育,吳劍疆,張磊,等. 天然河道冰塞演變發(fā)展的數(shù)值模擬[J]. 水科學(xué)進(jìn)展,2003,14(6):700-705.

    Mao Zeyu, Wu Jianjiang, Zhang Lei, et al. Numerical simulation of river ice jam[J]. 水科學(xué)進(jìn)展, 2003, 14(6): 700-705. (in Chinese with English abstract)

    [5] 茅澤育,趙升偉,相鵬,等. 冰蓋下水流流動(dòng)的摻混特性[J]. 水利學(xué)報(bào),2005,36(3):291-297.

    Mao Zeyu, Zhao Shengwei, Xiang Peng, et al. Turbulentmixing characteristics of ice-covered river flow[J]. Journal of Hydraulic Engineering, 2005, 36(3): 291-297. (in Chinese with English abstract)

    [6] 茅澤育,許昕,王愛(ài)民,等. 基于適體坐標(biāo)變換的二維河冰模型[J]. 水科學(xué)進(jìn)展,2008,19(2):214-223.

    Mao Zeyu, Xu Xin, Wang Aimin, et al. 2D numerical model for river-ice processes based upon body-fitted coordinate[J]. Advances in Water Science, 2008, 19(2): 214-223. (in Chinese with English abstract)

    [7] 肖建民,金龍海,謝永剛,等. 寒區(qū)水庫(kù)冰蓋形成與消融機(jī)理分析[J]. 水利學(xué)報(bào),2004,35(6):80-85.

    Xiao Jianmin, Jin Longhai, Xie Yonggang, et al. Study on mechanism of formation and melting of reservoir ice cover in cold area[J]. Journal of Hydraulic Engineering, 2004, 35(6): 80-85. (in Chinese with English abstract)

    [8] 王軍,陳胖胖,江濤,等. 冰蓋下冰塞堆積的數(shù)值模擬[J]. 水利學(xué)報(bào),2009,40(3):348-354+363.

    Wang Jun, Chen Pangpang, Jang Tao, et al. Parameter model of water-conducting device specification for indirect subsurface drip irrigation[J]. Journal of Hydraulic Engineering, 2009, 40(3): 348-354+363. (in Chinese with English abstract)

    [9] 徐國(guó)賓,李大冉,黃焱,等. 南水北調(diào)中線輸水工程若干冰力學(xué)問(wèn)題試驗(yàn)研究[J]. 水科學(xué)進(jìn)展,2010,21(6):808-815.

    Xu Guobin, Li Daran, Huang Yan, et al. Laboratory study of problems in ice mechanics encountered in the Middle Route of South-to-North Water Transfer Project[J]. Advances in Water Science, 2010, 21(6): 808-815. (in Chinese with English abstract)

    [10] 貢力,靳春玲. 西部干寒地區(qū)引水明渠病害特點(diǎn)及治理措施[J]. 建設(shè)監(jiān)理,2014(12):66-68, 73.

    Gong Li, Jin Chunling. The characteristics and treatment measures of open diversion canal diseases in dry and cold regions of western China [J]. Journal of Construction Supervision, 2014(12): 66-68, 73. (in Chinese with English abstract)

    [11] 李釗. 石頭河水庫(kù)灌區(qū)東干渠水流數(shù)值模擬研究[D]. 西安:西安理工大學(xué),2009.

    Li Zhao. Research on the Stream Numerical Simulation of Shitou River Reservoir Irrigation District East Main Channel[D]. Xi’an: Xi’an University of Technology, 2009. (in Chinese with English abstract)

    [12] Zufelt J E, Ettema R. Fully coupled model of ice-jamdynamics[J]. Journal of Cold Regions Engineering, 2000, 14(1): 24-41.

    [13] Shen Hungtao, Sun Junshan, Liu Lianwu. SPH simulation of river ice dynamics[J]. Journal of Computational Physics, 2000, 165(2): 752-770.

    [14] 陳明千. 西藏高寒地區(qū)引水渠道冰花生消規(guī)律研究[D]. 成都:四川大學(xué),2006.

    Chen Mingqian. Study on the Process of Ice Formation and Melting in Diversion Channel of Tibet[D]. Chengdu: Sichuan University, 2006. (in Chinese with English abstract)

    [15] 吳素杰,宗全利,鄭鐵剛,等. 高寒區(qū)多口融冰井引水渠道水溫變化三維模擬及井群優(yōu)化布置[J]. 農(nóng)業(yè)工程學(xué)報(bào),2017,33(14):130-137.

    Wu Sujie, Zong Quanli, Zheng Tiegang, et al. 3D simulation on water temperature change of diversion channel and optimal arrangement of multi-wells at high altitude and cold regions[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2017, 33(14): 130-137. (in Chinese with English abstract)

    [16] Betchelor G K. Mass transfer from small particles suspended in turbulent fluid[J]. Journal of Fluid Mechanics, 1980, 98(3): 609-623.

    [17] 于天來(lái),袁正國(guó),黃美蘭. 河冰力學(xué)性能試驗(yàn)研究[J]. 遼寧工程技術(shù)大學(xué)學(xué)報(bào):自然科學(xué)版,2009,28(6):937-940.

    Yu Tianlai, Yuan Zhengguo, Huang Meilan. Experimental study on mechanical behavior of river ice[J]. Journal of Liaoning Technical University: Natural Science, 2009, 28(6): 937-940. (in Chinese with English abstract)

    [18] 張鳳德,李廣一. 基于LS-DYNA的流冰撞擊壩體仿真計(jì)算[J]. 水利建設(shè)與管理,2013,33(2):19-21.

    Zhang Fengde, Li Guangyi. Simulation of flow ice impact on dam body based on LS-DYNA[J]. Water Conservancy Construction and Management, 2013, 33(2): 19-21. (in Chinese with English abstract)

    [19] 李抗彬,沈冰,李智錄,等. 基于非恒定水流模擬的灌區(qū)明渠水力響應(yīng)特征分析[J]. 農(nóng)業(yè)工程學(xué)報(bào),2015,31(10):107-114.

    Li Kangbin, Shen Bing, Li Zhilu, et al. Open channel hydraulic response characteristics in irrigation area based on unsteady flow simulation analysis[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2015, 31(10): 107-114. (in Chinese with English abstract)

    [20] 王亞玲,劉應(yīng)中,繆國(guó)平. 圓柱繞流的三維數(shù)值模擬[J]. 上海交通大學(xué)學(xué)報(bào),2001,35(10):1464-1469.

    Wang Yaling, Liu Yingzhong, Miao Guoping. Three-dimensional numerical simulation of viscous flow around circular cylinder[J]. Journal of Shanghai Jiaotong University, 2001, 35(10): 1464-1469. (in Chinese with English abstract)

    [21] 張寬地,呂宏興,陳俊英. 馬蹄形過(guò)水?dāng)嗝媾R界水深的直接計(jì)算法[J]. 農(nóng)業(yè)工程學(xué)報(bào),2009,25(4):15-18.

    Zhang Kuandi, Lü Hongxing, Chen Junying. Direct calculation of critical depth of horseshoe section tunnel[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2009, 25(4): 15-18. (in Chinese with English abstract)

    [22] 李煒. 水力計(jì)算手冊(cè)[M]. 北京:中國(guó)水利水電出版社,2006.

    [23] 張成. 南水北調(diào)中線工程非恒定輸水響應(yīng)及運(yùn)行控制研究[D]. 北京:清華大學(xué),2008.

    Zhang Cheng. Research on Response of Unsteady Water Transport and Operation Control in the Middle Route of the South-to-North Water Diversion Project[D]. Beijing: Tsinghua University, 2008. (in Chinese with English abstract)

    [24] 韓延成,高學(xué)平. 長(zhǎng)距離自流型渠道輸水控制的二步法研究[J]. 水科學(xué)進(jìn)展,2006,17(3):414-418.

    Han Yancheng, Gao Xueping. Two-step optimal operation and control method for long distance gravity-flow delivery in canals[J]. Advance in Water Science, 2006, 17(3): 414-418. (in Chinese with English abstract)

    [25] 高霈生,靳國(guó)厚,呂斌秀. 南水北調(diào)中線工程輸水冰情的初步分析[J]. 水利學(xué)報(bào),2003,34(11):96-102.

    Gao Peisheng, Jin Guohou, Lü Binxiu. Preliminary study on ice regime in the middle route of south to north water transfer project[J]. Journal of Hydraulic Engineering, 2003, 34(11): 96-102. (in Chinese with English abstract)

    [26] 李家春. 現(xiàn)代流體力學(xué)發(fā)展的回顧與展望[J]. 力學(xué)進(jìn)展,1995,25(4):442-450.

    Li Jiachun. Retrospects and prospects of fluid mechanics[J]. Advances in Mechanic, 1995, 25(4): 442-450. (in Chinese with English abstract)

    [27] 劉孟凱,王長(zhǎng)德,馮曉波. 長(zhǎng)距離控制渠系結(jié)冰期的水力響應(yīng)分析[J]. 農(nóng)業(yè)工程學(xué)報(bào),2011,27(2):20-27.

    Liu Mengkai, Wang Changde, Feng Xiaobo. Analysis on the hydraulic response of long distance canal control system during ice period[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2011, 27(2): 20-27. (in Chinese with English abstract)

    [28] 王曉玲,周正印,蔣志勇,等. 考慮氣溫變化影響的引水渠道水內(nèi)冰演變數(shù)值模擬[J]. 天津大學(xué)學(xué)報(bào):自然科學(xué)與工程技術(shù)版,2010,43(6):515-522.

    Wang Xiaoling, Zhou Zhengyin, Jiang Zhiyong, et al. Numerical simulation of frazil ice evolution in diversion channel considering effect of temperature variation[J]. Journal of Tianjin University, 2010, 43(6): 515-522. (in Chinese with English abstract)

    [29] 成都科學(xué)技術(shù)大學(xué)水力學(xué)教研室. 水力學(xué)(下冊(cè))[M]. 北京:人民教育出版社,1979.

    [30] 王曉玲,張自強(qiáng),李濤,等. 引水流量對(duì)引水渠道中水內(nèi)冰演變影響的數(shù)值模擬[J]. 水利學(xué)報(bào),2009,40(11):1307-1312.

    Wang Xiaoling, Zhang Ziqiang, Li Tao, et al. Numerical simulation of diversion water flux effect on frazil ice evolution in diversion channel[J]. Journal of Hydraulic Engineering, 2009, 40(11): 1307-1312. (in Chinese with English abstract)

    Numerical simulation and verification on impact damage mechanical property of drift ice on diversion tunnel

    Gong Li, Li Yaxian, Jin Chunling

    (School of Civil Engineering, Lanzhou Jiaotong University, Lanzhou 730070, China)

    In high latitude region of western China, the environment is in harsh condition that it is cold and dry in winter and has long ice period, which causes many water resource problems. In order to relieve the serious water shortage condition in cold and dry region, a large number of long distance water diversion projects were established to improve the water resource condition, such as increasing farm irrigation, human and animal drinking. While the ice damage occurs frequently under severe ice conditions in cold and dry region, especially in ice period in winter and thawing period in spring, it is easy to form drift ice with different velocities, different plan sizes and different thicknesses, which produces different extrusion forces or impact forces to damage tunnel lining, causing project failure. The failure project could not realize the original planning and construction goal, giving rise to the water allocation pressure. The water allocation would cause water shortage which influences diversion irrigation and farming production in spring. Based on the intense researches on the collision simulation problem of the interaction between drift ice and diversion tunnel, this paper used the symmetric penalty function in the finite element contact-impact algorithm to conduct the theoretical study on collision simulation problem between drift ice and water diversion tunnel. ANSYS/LS-DYNA was adopted as the platform to establish tunnel model and drift ice model. LS-DYNA SOLVER was used as the solver to solve and analyze the damage degrees of drift ice on tunnel. The physical model tests were conducted to verify and reveal the impact damage mechanism of drift ice on diversion tunnel. The physical model was constructed by the geometric scale of 28, which is the ratio of the experiment facility to the prototype in the test. The results show that tunnel lining surface will form varying degrees of deformation and failure when the tunnel lining is impacted by the drift ice with different velocities, different plane sizes and different thicknesses. It is also discovered that the impact stress increases with the flow velocity and their relationship presents linear variation. The impact stress also increases with the drift ice’s plane size and their relationship presents nonlinear variation. The impact stress increases with the drift ice thickness when the drift ice thickness is less than 0.5 m. While the drift ice thickness is greater than 0.5 m, the maximum stress value shows little change. The relationship between drift ice’s plane size and maximum stress shows approximately linear variation. Meanwhile, the software simulation and test observation results are almost the same. The impact of drift ice on the tunnel lining would cause the deformation of lining, but the deformation has little influence on the tunnel stability. The drift ice’s long time erosion would cause the tunnel lining surface to fall off, and further break the strength and stability of the tunnel structures. The study supplies theoretical support and technical guarantee for water diversion project security in cold and dry region of western China.

    mechanical performance; numerical simulation; models; ice engineering; diversion tunnel; impact force

    貢 力,李雅嫻,靳春玲. 流冰對(duì)引水隧洞撞擊破壞力學(xué)特性數(shù)值分析與驗(yàn)證[J]. 農(nóng)業(yè)工程學(xué)報(bào),2018,34(13):144-151. doi:10.11975/j.issn.1002-6819.2018.13.017 http://www.tcsae.org

    Gong Li, Li Yaxian, Jin Chunling. Numerical simulation and verification on impact damage mechanical property of drift ice on diversion tunnel[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2018, 34(13): 144-151. (in Chinese with English abstract) doi:10.11975/j.issn.1002-6819.2018.13.017 http://www.tcsae.org

    2018-03-12

    2018-06-01

    國(guó)家自然科學(xué)基金項(xiàng)目(51669010,51541902);甘肅省自然基金(17JR5RA105,17JR5RA101)

    貢力,教授,博士,主要從事輸水工程安全的研究。Email:gongli@mail.lzjtu.cn

    10.11975/j.issn.1002-6819.2018.13.017

    TV672

    A

    1002-6819(2018)-13-0144-08

    猜你喜歡
    撞擊力隧洞流速
    “流體壓強(qiáng)與流速的關(guān)系”知識(shí)鞏固
    『流體壓強(qiáng)與流速的關(guān)系』知識(shí)鞏固
    隧洞止水帶安裝質(zhì)量控制探討
    滇中引水工程大理段首條隧洞順利貫通
    山雨欲來(lái)風(fēng)滿樓之流體壓強(qiáng)與流速
    愛(ài)虛張聲勢(shì)的水
    接觸面對(duì)駁船撞擊橋墩動(dòng)力響應(yīng)的影響
    自密實(shí)混凝土在水工隧洞襯砌中的應(yīng)用
    受撞橋梁結(jié)構(gòu)撞擊力仿真分析研究
    樁基布置對(duì)高樁碼頭撞擊力分配的影響
    水道港口(2014年1期)2014-04-27 14:14:40
    午夜福利视频精品| 亚洲av日韩精品久久久久久密 | 男女之事视频高清在线观看 | 欧美中文综合在线视频| 美女国产高潮福利片在线看| 涩涩av久久男人的天堂| 国产熟女欧美一区二区| 日本wwww免费看| 久久青草综合色| 激情视频va一区二区三区| 成人午夜精彩视频在线观看| 国产极品天堂在线| 制服人妻中文乱码| 97精品久久久久久久久久精品| 精品一品国产午夜福利视频| 国产成人精品无人区| 一级毛片 在线播放| 亚洲成人手机| 99久久人妻综合| 岛国毛片在线播放| 久久这里只有精品19| 丁香六月欧美| 亚洲av成人不卡在线观看播放网 | 18禁国产床啪视频网站| 青草久久国产| 久久久久久久久免费视频了| 黄网站色视频无遮挡免费观看| 在线观看免费视频网站a站| 天堂中文最新版在线下载| 精品久久久久久电影网| 一级毛片我不卡| 国产在线一区二区三区精| 亚洲国产日韩一区二区| 性少妇av在线| 老司机亚洲免费影院| 少妇精品久久久久久久| 制服丝袜香蕉在线| 夫妻午夜视频| 水蜜桃什么品种好| 亚洲欧美成人精品一区二区| av免费观看日本| 91老司机精品| 激情视频va一区二区三区| 成年动漫av网址| 国产一区二区三区综合在线观看| a级片在线免费高清观看视频| 老司机深夜福利视频在线观看 | 国产淫语在线视频| 亚洲久久久国产精品| 不卡视频在线观看欧美| 国产高清不卡午夜福利| 免费观看av网站的网址| 国产黄色免费在线视频| 亚洲美女黄色视频免费看| 日本午夜av视频| 亚洲精品乱久久久久久| 自拍欧美九色日韩亚洲蝌蚪91| 成人18禁高潮啪啪吃奶动态图| 久久 成人 亚洲| 亚洲欧美精品综合一区二区三区| 热99久久久久精品小说推荐| 久久性视频一级片| 久久ye,这里只有精品| 国产av国产精品国产| 精品一区二区免费观看| 欧美日韩综合久久久久久| 黄频高清免费视频| 亚洲熟女精品中文字幕| 日韩一卡2卡3卡4卡2021年| 欧美精品av麻豆av| 日本av手机在线免费观看| 久久97久久精品| av不卡在线播放| 日韩 亚洲 欧美在线| 久久久久久免费高清国产稀缺| 97在线人人人人妻| 亚洲成人免费av在线播放| 亚洲视频免费观看视频| 日韩 欧美 亚洲 中文字幕| 人妻一区二区av| 国产精品二区激情视频| 日韩熟女老妇一区二区性免费视频| 国产熟女午夜一区二区三区| 91精品国产国语对白视频| 亚洲精品国产av蜜桃| 美女国产高潮福利片在线看| 亚洲一区中文字幕在线| 婷婷色av中文字幕| 精品视频人人做人人爽| 国产探花极品一区二区| 欧美在线一区亚洲| 观看av在线不卡| 制服人妻中文乱码| 国产 精品1| 精品亚洲成国产av| 欧美激情极品国产一区二区三区| 欧美日韩成人在线一区二区| 少妇人妻久久综合中文| 亚洲国产精品一区三区| 午夜激情av网站| 国产在线一区二区三区精| 午夜福利乱码中文字幕| 久久精品国产亚洲av涩爱| 满18在线观看网站| 国产极品粉嫩免费观看在线| 国产精品三级大全| 黄色毛片三级朝国网站| 女性被躁到高潮视频| 欧美成人午夜精品| 久久久国产欧美日韩av| 人人妻人人爽人人添夜夜欢视频| 黑丝袜美女国产一区| 久久99一区二区三区| 午夜老司机福利片| 久久国产亚洲av麻豆专区| 美女脱内裤让男人舔精品视频| 亚洲七黄色美女视频| 极品人妻少妇av视频| 在线观看免费高清a一片| 少妇人妻精品综合一区二区| 涩涩av久久男人的天堂| 亚洲成av片中文字幕在线观看| 成人毛片60女人毛片免费| 国产国语露脸激情在线看| 男女高潮啪啪啪动态图| 啦啦啦视频在线资源免费观看| 老司机深夜福利视频在线观看 | 婷婷色综合大香蕉| 国产高清国产精品国产三级| 免费日韩欧美在线观看| 亚洲精品av麻豆狂野| a级毛片在线看网站| 色吧在线观看| 精品第一国产精品| 久久久久精品性色| 美女午夜性视频免费| videosex国产| 国产精品国产三级专区第一集| 欧美少妇被猛烈插入视频| 欧美日韩精品网址| 午夜福利乱码中文字幕| 国产成人精品无人区| 女人精品久久久久毛片| 热re99久久国产66热| 久久久久久久大尺度免费视频| 岛国毛片在线播放| 丝袜人妻中文字幕| 啦啦啦在线观看免费高清www| 99久久99久久久精品蜜桃| 飞空精品影院首页| 久久精品aⅴ一区二区三区四区| 欧美日韩综合久久久久久| 亚洲色图 男人天堂 中文字幕| www.精华液| 国产在视频线精品| av在线老鸭窝| 免费观看av网站的网址| 亚洲av男天堂| 色播在线永久视频| 美女视频免费永久观看网站| 蜜桃在线观看..| 你懂的网址亚洲精品在线观看| 69精品国产乱码久久久| 国产日韩一区二区三区精品不卡| 国产亚洲av高清不卡| 久久婷婷青草| 国产精品久久久久久人妻精品电影 | 十八禁高潮呻吟视频| 亚洲精品久久成人aⅴ小说| 一级片免费观看大全| 国产精品亚洲av一区麻豆 | 国产在线视频一区二区| 纵有疾风起免费观看全集完整版| 亚洲av男天堂| 亚洲欧美清纯卡通| 婷婷色av中文字幕| 尾随美女入室| 高清黄色对白视频在线免费看| 亚洲欧美成人综合另类久久久| 一级毛片电影观看| 女性被躁到高潮视频| 少妇精品久久久久久久| 在线天堂最新版资源| 亚洲欧美成人综合另类久久久| 欧美中文综合在线视频| 在现免费观看毛片| 一级爰片在线观看| 性高湖久久久久久久久免费观看| 国产精品成人在线| 日韩制服骚丝袜av| 久久 成人 亚洲| 国产深夜福利视频在线观看| 99国产精品免费福利视频| 国产黄频视频在线观看| 色综合欧美亚洲国产小说| 丁香六月欧美| 国产又爽黄色视频| 国产有黄有色有爽视频| 女人精品久久久久毛片| 亚洲av福利一区| 人人澡人人妻人| 性色av一级| 久久久欧美国产精品| 日韩人妻精品一区2区三区| av免费观看日本| 国产亚洲欧美精品永久| 精品卡一卡二卡四卡免费| 美女国产高潮福利片在线看| 丝瓜视频免费看黄片| 精品人妻熟女毛片av久久网站| 成人手机av| 精品少妇久久久久久888优播| 咕卡用的链子| 伊人久久大香线蕉亚洲五| 欧美激情极品国产一区二区三区| 免费高清在线观看视频在线观看| 在线观看三级黄色| 99久久99久久久精品蜜桃| av在线app专区| 肉色欧美久久久久久久蜜桃| 中国三级夫妇交换| 老司机在亚洲福利影院| 人妻一区二区av| 久久性视频一级片| 亚洲欧洲精品一区二区精品久久久 | 国产精品女同一区二区软件| 亚洲精品中文字幕在线视频| 亚洲久久久国产精品| 久久久精品国产亚洲av高清涩受| 国产一区二区 视频在线| 日韩 亚洲 欧美在线| 精品一品国产午夜福利视频| 免费观看av网站的网址| 久久久精品国产亚洲av高清涩受| 欧美日韩国产mv在线观看视频| 巨乳人妻的诱惑在线观看| 国产熟女欧美一区二区| 精品一区二区三区四区五区乱码 | 亚洲国产av新网站| 亚洲国产欧美网| 婷婷色av中文字幕| 日本黄色日本黄色录像| 精品国产乱码久久久久久小说| 国产成人av激情在线播放| 九九爱精品视频在线观看| 亚洲成人av在线免费| 操出白浆在线播放| 美女福利国产在线| videos熟女内射| 永久免费av网站大全| 色综合欧美亚洲国产小说| 在线免费观看不下载黄p国产| 人妻一区二区av| 在线观看免费日韩欧美大片| 久久ye,这里只有精品| 亚洲一区中文字幕在线| 亚洲精品一区蜜桃| 尾随美女入室| 日韩人妻精品一区2区三区| 国产亚洲午夜精品一区二区久久| 一二三四在线观看免费中文在| 中文字幕人妻丝袜一区二区 | 亚洲五月色婷婷综合| 亚洲精品美女久久久久99蜜臀 | 精品第一国产精品| 国产亚洲av高清不卡| 一区二区日韩欧美中文字幕| 成人亚洲欧美一区二区av| 日本一区二区免费在线视频| 99久久综合免费| 精品国产一区二区三区四区第35| 婷婷色麻豆天堂久久| 中文字幕最新亚洲高清| 丝袜脚勾引网站| 亚洲男人天堂网一区| 欧美日韩一级在线毛片| av国产精品久久久久影院| 久久久久久人妻| 久久久国产精品麻豆| 亚洲欧美成人综合另类久久久| 夜夜骑夜夜射夜夜干| 日韩电影二区| 国产黄频视频在线观看| 校园人妻丝袜中文字幕| 免费看av在线观看网站| 亚洲,欧美,日韩| 亚洲美女黄色视频免费看| 久久精品国产综合久久久| 男人添女人高潮全过程视频| 香蕉国产在线看| 欧美国产精品va在线观看不卡| 久久狼人影院| 1024视频免费在线观看| 国产片内射在线| 亚洲av成人不卡在线观看播放网 | 天堂8中文在线网| 免费日韩欧美在线观看| 欧美精品高潮呻吟av久久| 十八禁人妻一区二区| 99久久综合免费| 免费黄频网站在线观看国产| 男女之事视频高清在线观看 | 日本午夜av视频| 久久精品亚洲熟妇少妇任你| 91精品三级在线观看| 午夜激情久久久久久久| 激情视频va一区二区三区| 精品酒店卫生间| 亚洲欧洲日产国产| 性色av一级| 国产一区二区三区综合在线观看| 极品人妻少妇av视频| 亚洲精品在线美女| 最近中文字幕2019免费版| √禁漫天堂资源中文www| 悠悠久久av| 国产成人精品久久久久久| 亚洲国产精品一区三区| 国产成人系列免费观看| 国产黄色免费在线视频| 久久久久久久久免费视频了| 午夜福利视频在线观看免费| 一级片免费观看大全| 伦理电影免费视频| 十八禁高潮呻吟视频| 久久久久久久久久久久大奶| 一级毛片电影观看| 亚洲国产精品国产精品| 国产精品成人在线| 91精品伊人久久大香线蕉| 波多野结衣av一区二区av| 久久久久久久久久久久大奶| 亚洲精品av麻豆狂野| 男人添女人高潮全过程视频| 日韩精品有码人妻一区| 婷婷色综合www| 国产在线免费精品| 在线观看国产h片| 免费少妇av软件| 精品一区二区三区四区五区乱码 | 国产精品亚洲av一区麻豆 | 制服丝袜香蕉在线| 2018国产大陆天天弄谢| 国产精品免费视频内射| 国产一区二区三区综合在线观看| 2018国产大陆天天弄谢| 999精品在线视频| 香蕉丝袜av| 亚洲国产精品一区三区| av在线app专区| 国产精品一国产av| 女人爽到高潮嗷嗷叫在线视频| 另类亚洲欧美激情| 日本一区二区免费在线视频| 老司机深夜福利视频在线观看 | 在现免费观看毛片| 国产又爽黄色视频| 两性夫妻黄色片| 亚洲婷婷狠狠爱综合网| 国产色婷婷99| 午夜老司机福利片| 曰老女人黄片| 国产一区二区三区综合在线观看| 色婷婷av一区二区三区视频| 少妇人妻精品综合一区二区| 女人久久www免费人成看片| 久久久国产一区二区| 亚洲精品国产区一区二| 你懂的网址亚洲精品在线观看| 精品第一国产精品| 男女午夜视频在线观看| 最新的欧美精品一区二区| 欧美日韩国产mv在线观看视频| 久久精品aⅴ一区二区三区四区| 最黄视频免费看| 亚洲欧美清纯卡通| h视频一区二区三区| 亚洲欧美成人综合另类久久久| 久久久久精品性色| 亚洲精品美女久久久久99蜜臀 | 国产精品欧美亚洲77777| www.自偷自拍.com| 日韩av在线免费看完整版不卡| 欧美日韩视频高清一区二区三区二| 国产日韩欧美亚洲二区| 国产精品国产三级专区第一集| 欧美人与性动交α欧美软件| 黄色一级大片看看| 国产人伦9x9x在线观看| av在线老鸭窝| 日韩伦理黄色片| 国产精品国产三级国产专区5o| 黄频高清免费视频| 老司机亚洲免费影院| 国产女主播在线喷水免费视频网站| a 毛片基地| 视频在线观看一区二区三区| 久久久久久免费高清国产稀缺| 亚洲熟女精品中文字幕| 99久久人妻综合| e午夜精品久久久久久久| 天天躁狠狠躁夜夜躁狠狠躁| 男女边吃奶边做爰视频| 一本一本久久a久久精品综合妖精| 欧美日韩视频高清一区二区三区二| 十八禁网站网址无遮挡| 亚洲国产欧美网| 一级毛片黄色毛片免费观看视频| 老司机影院成人| 国产女主播在线喷水免费视频网站| 欧美精品亚洲一区二区| 国产精品成人在线| 久久人人97超碰香蕉20202| 亚洲美女视频黄频| 在线观看一区二区三区激情| 欧美少妇被猛烈插入视频| 尾随美女入室| 91精品国产国语对白视频| 大陆偷拍与自拍| 99久国产av精品国产电影| 永久免费av网站大全| 少妇人妻精品综合一区二区| 一区二区三区精品91| 欧美 日韩 精品 国产| 亚洲国产av新网站| 亚洲国产看品久久| 久久热在线av| 性高湖久久久久久久久免费观看| 亚洲自偷自拍图片 自拍| 亚洲av欧美aⅴ国产| 侵犯人妻中文字幕一二三四区| 一区二区三区精品91| 晚上一个人看的免费电影| 精品一区二区三区四区五区乱码 | 宅男免费午夜| 丁香六月欧美| 亚洲第一区二区三区不卡| 黄色视频在线播放观看不卡| 亚洲美女黄色视频免费看| 自线自在国产av| 爱豆传媒免费全集在线观看| 免费看不卡的av| 校园人妻丝袜中文字幕| 国产精品av久久久久免费| 国产片特级美女逼逼视频| 国产极品天堂在线| 天天躁夜夜躁狠狠躁躁| 亚洲国产欧美在线一区| 一区二区日韩欧美中文字幕| xxxhd国产人妻xxx| 久久国产精品大桥未久av| 老司机影院成人| 亚洲精品久久久久久婷婷小说| 欧美 日韩 精品 国产| 日本vs欧美在线观看视频| 欧美在线黄色| 亚洲情色 制服丝袜| 美女高潮到喷水免费观看| 一本久久精品| 欧美精品人与动牲交sv欧美| 午夜福利网站1000一区二区三区| 国产日韩欧美视频二区| 亚洲av成人精品一二三区| 搡老岳熟女国产| 久久久精品94久久精品| 精品酒店卫生间| 2021少妇久久久久久久久久久| 亚洲av中文av极速乱| 欧美日韩亚洲高清精品| 成年女人毛片免费观看观看9 | 亚洲成色77777| 国产精品嫩草影院av在线观看| 亚洲精品日韩在线中文字幕| 欧美日韩亚洲高清精品| 日韩中文字幕视频在线看片| 国产又色又爽无遮挡免| 中文字幕av电影在线播放| av.在线天堂| 最近2019中文字幕mv第一页| 日本午夜av视频| 午夜福利视频在线观看免费| 国产精品一区二区在线不卡| 亚洲av电影在线进入| 久久久久精品性色| 别揉我奶头~嗯~啊~动态视频 | 熟妇人妻不卡中文字幕| 日韩av免费高清视频| 深夜精品福利| 午夜福利影视在线免费观看| 日韩大码丰满熟妇| 伊人久久国产一区二区| 男女国产视频网站| 性高湖久久久久久久久免费观看| 亚洲国产精品成人久久小说| 日本午夜av视频| 国产精品久久久av美女十八| 制服人妻中文乱码| 精品人妻一区二区三区麻豆| 日韩精品有码人妻一区| 狠狠精品人妻久久久久久综合| av免费观看日本| 国产欧美日韩一区二区三区在线| 国产熟女午夜一区二区三区| 国产精品久久久久久久久免| 99香蕉大伊视频| 丝瓜视频免费看黄片| 亚洲国产av新网站| 极品少妇高潮喷水抽搐| 亚洲一级一片aⅴ在线观看| tube8黄色片| xxx大片免费视频| avwww免费| 丰满乱子伦码专区| 欧美老熟妇乱子伦牲交| 91老司机精品| 国产成人精品福利久久| 色播在线永久视频| 国产有黄有色有爽视频| 欧美国产精品va在线观看不卡| 亚洲国产精品一区二区三区在线| 一级黄片播放器| 精品国产一区二区久久| 考比视频在线观看| 99久久99久久久精品蜜桃| 亚洲激情五月婷婷啪啪| 亚洲精品一区蜜桃| 久久午夜综合久久蜜桃| 国产av码专区亚洲av| 亚洲熟女毛片儿| 亚洲免费av在线视频| 尾随美女入室| 国产精品久久久久成人av| 大香蕉久久成人网| 亚洲国产毛片av蜜桃av| 日本av手机在线免费观看| 亚洲欧美一区二区三区黑人| 如何舔出高潮| 午夜久久久在线观看| 国产精品.久久久| 久久久久精品性色| 午夜老司机福利片| 亚洲视频免费观看视频| 午夜福利影视在线免费观看| 国产精品无大码| 午夜影院在线不卡| 亚洲色图综合在线观看| 男女之事视频高清在线观看 | 少妇精品久久久久久久| 日韩成人av中文字幕在线观看| 亚洲第一区二区三区不卡| 老汉色∧v一级毛片| 久久久久国产一级毛片高清牌| 十八禁网站网址无遮挡| 下体分泌物呈黄色| 精品国产一区二区久久| 宅男免费午夜| 老司机深夜福利视频在线观看 | 亚洲av中文av极速乱| 男女午夜视频在线观看| 天堂8中文在线网| 亚洲av日韩在线播放| 国产欧美日韩一区二区三区在线| 欧美国产精品一级二级三级| 国产av精品麻豆| 妹子高潮喷水视频| 大码成人一级视频| 国产成人精品福利久久| 97精品久久久久久久久久精品| av片东京热男人的天堂| 在线精品无人区一区二区三| 国产片特级美女逼逼视频| 人人妻人人澡人人看| 国产不卡av网站在线观看| 婷婷色综合www| 波多野结衣一区麻豆| 亚洲精品久久成人aⅴ小说| 99久久99久久久精品蜜桃| 成人手机av| 男女无遮挡免费网站观看| 精品少妇内射三级| 大香蕉久久成人网| 欧美97在线视频| 男男h啪啪无遮挡| 亚洲人成电影观看| 国产xxxxx性猛交| 中文字幕高清在线视频| 在线观看免费日韩欧美大片| 日韩视频在线欧美| 久久国产亚洲av麻豆专区| 丝袜人妻中文字幕| 亚洲第一av免费看| 热re99久久国产66热| 一区二区三区精品91| 日韩视频在线欧美| 亚洲国产欧美网| 如日韩欧美国产精品一区二区三区| 久久精品熟女亚洲av麻豆精品| 亚洲国产欧美网| 国产精品一国产av| 成人午夜精彩视频在线观看| 午夜福利影视在线免费观看| 婷婷色综合大香蕉| 久久久欧美国产精品| 日本欧美国产在线视频| 国产国语露脸激情在线看| 欧美激情 高清一区二区三区| 五月天丁香电影| 精品久久久久久电影网| 十八禁网站网址无遮挡| 在线免费观看不下载黄p国产| 女人精品久久久久毛片| 国产成人午夜福利电影在线观看| 精品少妇久久久久久888优播| 精品国产超薄肉色丝袜足j| 肉色欧美久久久久久久蜜桃|