• <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片| 五月玫瑰六月丁香| 亚洲熟妇熟女久久| 在线观看美女被高潮喷水网站 | 精品一区二区三区视频在线 | 久久精品国产综合久久久| xxx96com| 黑人欧美特级aaaaaa片| 啦啦啦观看免费观看视频高清| 日日夜夜操网爽| 色尼玛亚洲综合影院| 亚洲久久久久久中文字幕| 亚洲久久久久久中文字幕| 亚洲第一欧美日韩一区二区三区| 国产97色在线日韩免费| 日韩有码中文字幕| 久久九九热精品免费| 亚洲国产欧美人成| 亚洲av成人av| 90打野战视频偷拍视频| 国产高清有码在线观看视频| 成人av在线播放网站| 成人性生交大片免费视频hd| 亚洲人成伊人成综合网2020| 好男人在线观看高清免费视频| 天天躁日日操中文字幕| 成人无遮挡网站| 国产欧美日韩一区二区精品| 桃色一区二区三区在线观看| 蜜桃亚洲精品一区二区三区| 美女高潮的动态| 国产69精品久久久久777片| 国产高清有码在线观看视频| tocl精华| 老司机午夜福利在线观看视频| 亚洲成av人片免费观看| 国产伦人伦偷精品视频| 天堂√8在线中文| 19禁男女啪啪无遮挡网站| а√天堂www在线а√下载| 夜夜爽天天搞| 国产精品久久久久久人妻精品电影| 在线观看午夜福利视频| 免费在线观看成人毛片| 久久久久性生活片| 久久久久性生活片| 国产成人欧美在线观看| 蜜桃亚洲精品一区二区三区| 欧洲精品卡2卡3卡4卡5卡区| 精品国产超薄肉色丝袜足j| 欧美最新免费一区二区三区 | www国产在线视频色| www国产在线视频色| 成人一区二区视频在线观看| 9191精品国产免费久久| 老熟妇乱子伦视频在线观看| 久久国产精品影院| 国产高潮美女av| 中文字幕av成人在线电影| 2021天堂中文幕一二区在线观| 国产精品香港三级国产av潘金莲| 亚洲熟妇中文字幕五十中出| 国产精品久久久人人做人人爽| 听说在线观看完整版免费高清| 天美传媒精品一区二区| av视频在线观看入口| 黄色片一级片一级黄色片| 欧美丝袜亚洲另类 | 久久人人精品亚洲av| 国产精品综合久久久久久久免费| 波多野结衣高清无吗| 欧美av亚洲av综合av国产av| 国产蜜桃级精品一区二区三区| 国产精品一区二区三区四区久久| 哪里可以看免费的av片| 极品教师在线免费播放| 久久久精品欧美日韩精品| 欧美性感艳星| 免费观看的影片在线观看| 日本撒尿小便嘘嘘汇集6| 国产精品久久电影中文字幕| bbb黄色大片| 精品一区二区三区视频在线 | 婷婷精品国产亚洲av在线| 国产不卡一卡二| 亚洲熟妇中文字幕五十中出| 超碰av人人做人人爽久久 | 国内久久婷婷六月综合欲色啪| 窝窝影院91人妻| 欧美绝顶高潮抽搐喷水| 免费在线观看影片大全网站| 色在线成人网| 波多野结衣高清作品| 亚洲性夜色夜夜综合| 日本 欧美在线| 欧美色欧美亚洲另类二区| 亚洲av成人精品一区久久| 精品国产超薄肉色丝袜足j| 亚洲精品色激情综合| 男插女下体视频免费在线播放| 国产在视频线在精品| 日韩中文字幕欧美一区二区| 宅男免费午夜| 欧美乱妇无乱码| 亚洲人成网站在线播| 免费在线观看日本一区| 久久国产精品影院| 在线免费观看的www视频| 亚洲五月婷婷丁香| 国产aⅴ精品一区二区三区波| 国产黄片美女视频| 国产一级毛片七仙女欲春2| 午夜福利在线在线| 色吧在线观看| 亚洲狠狠婷婷综合久久图片| 精品久久久久久久末码| 欧美中文日本在线观看视频| xxx96com| 欧美不卡视频在线免费观看| 国产成年人精品一区二区| 亚洲成人久久爱视频| 久久精品国产综合久久久| 国产三级在线视频| 欧美zozozo另类| 亚洲国产中文字幕在线视频| 欧美三级亚洲精品| 男女之事视频高清在线观看| 一本精品99久久精品77| 99国产极品粉嫩在线观看| 女人十人毛片免费观看3o分钟| 老师上课跳d突然被开到最大视频 久久午夜综合久久蜜桃 | 亚洲电影在线观看av| 亚洲国产欧美人成| 法律面前人人平等表现在哪些方面| 国产亚洲精品久久久久久毛片| 长腿黑丝高跟| 久99久视频精品免费| 99久久综合精品五月天人人| 国产一区二区激情短视频| 欧美成人a在线观看| 男插女下体视频免费在线播放| 国产色婷婷99| 久久精品91蜜桃| 久久中文看片网| 天天躁日日操中文字幕| 蜜桃亚洲精品一区二区三区| 欧美极品一区二区三区四区| 日韩国内少妇激情av| 国产老妇女一区| 搡女人真爽免费视频火全软件 | 午夜日韩欧美国产| 一进一出抽搐gif免费好疼| 宅男免费午夜| 69av精品久久久久久| 国产69精品久久久久777片| 国产av麻豆久久久久久久| 久久欧美精品欧美久久欧美| 久久久久久国产a免费观看| 国产日本99.免费观看| 深夜精品福利| 99国产精品一区二区蜜桃av| 精品日产1卡2卡| 日本熟妇午夜| 国产真人三级小视频在线观看| 免费看美女性在线毛片视频| 在线看三级毛片| 国产欧美日韩精品一区二区| 国产午夜福利久久久久久| 国产高清激情床上av| 国产三级黄色录像| 最近在线观看免费完整版| 国产在线精品亚洲第一网站| 在线a可以看的网站| 婷婷六月久久综合丁香| 男女之事视频高清在线观看| 桃色一区二区三区在线观看| 亚洲欧美日韩东京热| 嫩草影视91久久| 一个人免费在线观看的高清视频| 国产麻豆成人av免费视频| 日日夜夜操网爽| 亚洲精品在线美女| 99精品欧美一区二区三区四区| 亚洲五月婷婷丁香| 亚洲人成网站高清观看| 亚洲 欧美 日韩 在线 免费| 国产真人三级小视频在线观看| 成人特级av手机在线观看| 亚洲成av人片免费观看| 亚洲国产中文字幕在线视频| 少妇高潮的动态图| 欧美在线黄色| 亚洲七黄色美女视频| 91麻豆av在线| 国产成+人综合+亚洲专区| 色在线成人网| 黄色视频,在线免费观看| 成人性生交大片免费视频hd| 可以在线观看毛片的网站| 99久久久亚洲精品蜜臀av| 97超视频在线观看视频| 一区二区三区免费毛片| 中文在线观看免费www的网站| 亚洲精品一卡2卡三卡4卡5卡| 久久亚洲精品不卡| 男女之事视频高清在线观看| 成年女人毛片免费观看观看9| 少妇熟女aⅴ在线视频| 日韩 欧美 亚洲 中文字幕| 成人午夜高清在线视频| 一卡2卡三卡四卡精品乱码亚洲| 免费无遮挡裸体视频| 国产主播在线观看一区二区| 女人十人毛片免费观看3o分钟| 亚洲不卡免费看| 国产精品亚洲一级av第二区| 日韩欧美三级三区| 男人和女人高潮做爰伦理| 99久久精品一区二区三区| h日本视频在线播放| 在线观看一区二区三区| 男女视频在线观看网站免费| 亚洲人成网站在线播| 欧美日韩黄片免| 久久久久久国产a免费观看| 亚洲第一欧美日韩一区二区三区| 亚洲专区国产一区二区| 亚洲精品粉嫩美女一区| 免费在线观看亚洲国产| 99国产精品一区二区三区| 99国产极品粉嫩在线观看| 欧美黄色片欧美黄色片| 久久性视频一级片| 精品一区二区三区av网在线观看| 夜夜看夜夜爽夜夜摸| 日韩 欧美 亚洲 中文字幕| 国产精品久久久久久人妻精品电影| 国产精品 欧美亚洲| 国产视频一区二区在线看| 久久精品影院6| 在线免费观看不下载黄p国产 | 久久久国产成人免费| 久久久久国产精品人妻aⅴ院| 女生性感内裤真人,穿戴方法视频| 国产v大片淫在线免费观看| 国产视频一区二区在线看| 男人舔女人下体高潮全视频| 亚洲 欧美 日韩 在线 免费| 国产精品av视频在线免费观看| 俺也久久电影网| 99热这里只有是精品50| 麻豆久久精品国产亚洲av| 小说图片视频综合网站| 亚洲av免费在线观看| 日韩高清综合在线| 色尼玛亚洲综合影院| 精品国产亚洲在线| 久久久国产精品麻豆| www国产在线视频色| 欧美最新免费一区二区三区 | 国产高清视频在线播放一区| 亚洲无线在线观看| 中出人妻视频一区二区| 女同久久另类99精品国产91| 长腿黑丝高跟| 国产97色在线日韩免费| 日韩精品中文字幕看吧| 亚洲激情在线av| 成人av在线播放网站| 亚洲av中文字字幕乱码综合| 日韩精品中文字幕看吧| 搡女人真爽免费视频火全软件 | 成年人黄色毛片网站| 国产v大片淫在线免费观看| 网址你懂的国产日韩在线| 日韩欧美三级三区| 51国产日韩欧美| 在线国产一区二区在线| 亚洲国产精品999在线| 成人午夜高清在线视频| 午夜视频国产福利| 99在线人妻在线中文字幕| 亚洲第一电影网av| 黄色女人牲交| www.熟女人妻精品国产| 久久久成人免费电影| 国产v大片淫在线免费观看| 99国产极品粉嫩在线观看| 别揉我奶头~嗯~啊~动态视频| 又黄又粗又硬又大视频| 国产日本99.免费观看| 国产精品久久视频播放| av黄色大香蕉| 欧美区成人在线视频| 精品福利观看| 看免费av毛片| 变态另类成人亚洲欧美熟女| 久久精品国产清高在天天线| 夜夜爽天天搞| 特级一级黄色大片| 99国产综合亚洲精品| 国产激情欧美一区二区| 亚洲精品色激情综合| 亚洲人成网站在线播| 国产三级在线视频| 成人鲁丝片一二三区免费| 神马国产精品三级电影在线观看| 很黄的视频免费| 999久久久精品免费观看国产| 嫩草影院入口| 国产黄片美女视频| 最近最新免费中文字幕在线| 网址你懂的国产日韩在线| 亚洲七黄色美女视频| 国产伦在线观看视频一区| 久久精品综合一区二区三区| 99久久精品一区二区三区| АⅤ资源中文在线天堂| 亚洲在线观看片| 在线免费观看的www视频| 在线播放国产精品三级| 18禁国产床啪视频网站| 最近最新中文字幕大全电影3| 亚洲片人在线观看| 欧洲精品卡2卡3卡4卡5卡区| 国产精品久久久人人做人人爽| 蜜桃久久精品国产亚洲av| 国产一区二区在线观看日韩 | 欧美激情久久久久久爽电影| 69av精品久久久久久| 亚洲最大成人中文| 色播亚洲综合网| 色老头精品视频在线观看| 一个人免费在线观看电影| 久久久久免费精品人妻一区二区| 少妇裸体淫交视频免费看高清| 亚洲国产日韩欧美精品在线观看 | 在线十欧美十亚洲十日本专区| 免费高清视频大片| 偷拍熟女少妇极品色| 久久精品影院6| 国产真实乱freesex| 一个人看视频在线观看www免费 | 日本免费a在线| 国内久久婷婷六月综合欲色啪| 国产三级中文精品| 一区二区三区免费毛片| 亚洲成av人片在线播放无| 国产91精品成人一区二区三区| 欧美日韩福利视频一区二区| 我的老师免费观看完整版| 1024手机看黄色片| 成人av一区二区三区在线看| 精品熟女少妇八av免费久了| 久久欧美精品欧美久久欧美| 欧美3d第一页| 少妇丰满av| 成熟少妇高潮喷水视频| 久久99热这里只有精品18| 高清在线国产一区| 欧美黑人欧美精品刺激| 国产精品野战在线观看| 天天添夜夜摸| 亚洲欧美日韩东京热| 日日干狠狠操夜夜爽| 国产探花极品一区二区| 亚洲av熟女| 白带黄色成豆腐渣| 性色avwww在线观看| 国产精品久久久久久久电影 | 亚洲精品亚洲一区二区| 国产69精品久久久久777片| 日韩成人在线观看一区二区三区| 校园春色视频在线观看| www.色视频.com| 制服人妻中文乱码| 国产av不卡久久| 免费一级毛片在线播放高清视频| 露出奶头的视频| 午夜日韩欧美国产| 亚洲av一区综合| 搞女人的毛片| 91久久精品国产一区二区成人 | 欧美性感艳星| 欧美激情久久久久久爽电影| 在线观看午夜福利视频| 麻豆久久精品国产亚洲av| 一级作爱视频免费观看| 女人高潮潮喷娇喘18禁视频| 人妻久久中文字幕网| 久久亚洲精品不卡| 国产在线精品亚洲第一网站| 成年女人永久免费观看视频| 韩国av一区二区三区四区| 18禁在线播放成人免费| 欧美激情久久久久久爽电影| 俺也久久电影网| 动漫黄色视频在线观看| 黄色片一级片一级黄色片| av中文乱码字幕在线| 欧美中文日本在线观看视频| 男女做爰动态图高潮gif福利片| 日本在线视频免费播放| 十八禁网站免费在线| 天堂av国产一区二区熟女人妻| or卡值多少钱| 亚洲va日本ⅴa欧美va伊人久久| 两个人视频免费观看高清| 亚洲av二区三区四区| 欧美成人一区二区免费高清观看| 人妻夜夜爽99麻豆av| 少妇人妻精品综合一区二区 | 五月玫瑰六月丁香| 久久久久精品国产欧美久久久| 黄片小视频在线播放| 脱女人内裤的视频| 中文字幕熟女人妻在线| 欧美绝顶高潮抽搐喷水| 丝袜美腿在线中文| 久久久久久久亚洲中文字幕 | 国产久久久一区二区三区| 老汉色∧v一级毛片| 国产精品精品国产色婷婷| 性欧美人与动物交配| 三级国产精品欧美在线观看| 欧美黑人巨大hd| 99国产精品一区二区三区| 18禁在线播放成人免费| 午夜视频国产福利| 国产av不卡久久| 国产综合懂色| 国产精品久久久久久精品电影| 中文字幕精品亚洲无线码一区| 欧美成人性av电影在线观看| 国产在视频线在精品| 波多野结衣高清无吗| 桃红色精品国产亚洲av| 99热精品在线国产| 91麻豆精品激情在线观看国产| 成熟少妇高潮喷水视频| 婷婷精品国产亚洲av在线| 嫩草影院入口| 岛国在线观看网站| 一个人免费在线观看的高清视频| 国内精品久久久久精免费| 国产欧美日韩一区二区三| 亚洲avbb在线观看| 国产真实伦视频高清在线观看 | 美女免费视频网站| 国产在视频线在精品| 床上黄色一级片| 国产伦一二天堂av在线观看| 欧美xxxx黑人xx丫x性爽| 12—13女人毛片做爰片一| 好男人在线观看高清免费视频| 亚洲精品在线观看二区| 国产成人啪精品午夜网站| 美女高潮喷水抽搐中文字幕| 啦啦啦观看免费观看视频高清| 日韩欧美三级三区| 丰满人妻熟妇乱又伦精品不卡| 亚洲成人精品中文字幕电影| 国产欧美日韩一区二区精品| h日本视频在线播放| 国产一区二区三区在线臀色熟女| 天堂网av新在线| 3wmmmm亚洲av在线观看| 非洲黑人性xxxx精品又粗又长| 老汉色∧v一级毛片| a在线观看视频网站| 国产精品99久久99久久久不卡| 欧美激情在线99| 欧美日韩瑟瑟在线播放| 亚洲av第一区精品v没综合| 欧美日韩国产亚洲二区| 天天添夜夜摸| 成年版毛片免费区| 成人性生交大片免费视频hd| 啦啦啦免费观看视频1| 最新美女视频免费是黄的| 九九热线精品视视频播放| 色老头精品视频在线观看| 男人和女人高潮做爰伦理| 2021天堂中文幕一二区在线观| 日本在线视频免费播放| 久久久久久久午夜电影| 岛国在线观看网站| 99在线视频只有这里精品首页| 亚洲自拍偷在线| 中文在线观看免费www的网站| 丰满乱子伦码专区| 国产精品99久久99久久久不卡| 99精品久久久久人妻精品| 搡老岳熟女国产| 亚洲人成网站高清观看| 观看美女的网站| 一区二区三区激情视频| 日本黄色片子视频| 精品乱码久久久久久99久播| 欧美日韩国产亚洲二区| 日本一二三区视频观看| 亚洲成人精品中文字幕电影| 国产色爽女视频免费观看| 露出奶头的视频| 国产精品99久久99久久久不卡| 在线观看免费午夜福利视频| 少妇丰满av| 窝窝影院91人妻| 女人十人毛片免费观看3o分钟| av在线蜜桃| 亚洲第一欧美日韩一区二区三区| 亚洲国产欧洲综合997久久,| 天天一区二区日本电影三级| 欧美色视频一区免费| 丁香欧美五月| 亚洲在线观看片| 国产黄色小视频在线观看| 国产亚洲精品久久久久久毛片| 亚洲欧美激情综合另类| 国产精品久久久人人做人人爽| 小蜜桃在线观看免费完整版高清| 亚洲精华国产精华精| 日韩精品青青久久久久久| 国内久久婷婷六月综合欲色啪| 久久久国产精品麻豆| 亚洲精品一卡2卡三卡4卡5卡| 无人区码免费观看不卡| 男女下面进入的视频免费午夜| 看免费av毛片| 欧美在线一区亚洲| 嫩草影视91久久| 国产av一区在线观看免费| 亚洲av电影在线进入| 国产精品98久久久久久宅男小说| 少妇人妻精品综合一区二区 | 国产一区二区亚洲精品在线观看| 国产乱人伦免费视频| 无人区码免费观看不卡| 精品国产美女av久久久久小说| 在线播放无遮挡| 黄片大片在线免费观看| 亚洲美女黄片视频| 亚洲欧美日韩高清在线视频| 狠狠狠狠99中文字幕| 老汉色∧v一级毛片| 欧美极品一区二区三区四区| 琪琪午夜伦伦电影理论片6080| av福利片在线观看| 两人在一起打扑克的视频| 黄片小视频在线播放| 无限看片的www在线观看| 精品欧美国产一区二区三| 国产不卡一卡二| 少妇的逼水好多| 性色avwww在线观看| 日韩成人在线观看一区二区三区| 中文字幕久久专区| 欧美黑人巨大hd| 精品人妻1区二区| 黑人欧美特级aaaaaa片| 国产高清有码在线观看视频| 日韩人妻高清精品专区| 天堂影院成人在线观看| 国产伦人伦偷精品视频| ponron亚洲| 日本三级黄在线观看| 亚洲狠狠婷婷综合久久图片| 熟女电影av网| 99久国产av精品| 少妇的逼好多水| 国产精品99久久99久久久不卡| 99久久精品国产亚洲精品| 五月玫瑰六月丁香| 成年版毛片免费区| xxxwww97欧美| 老汉色av国产亚洲站长工具| 国产精品一区二区三区四区免费观看 | 特大巨黑吊av在线直播| 床上黄色一级片| 一级黄色大片毛片| 桃色一区二区三区在线观看| 老汉色av国产亚洲站长工具| 国语自产精品视频在线第100页| 国产日本99.免费观看| 中文字幕人妻丝袜一区二区| 亚洲色图av天堂| 99热6这里只有精品| 国产精品一及| 国产精品乱码一区二三区的特点| 亚洲人与动物交配视频| 欧美午夜高清在线| 搞女人的毛片| 无人区码免费观看不卡| 深夜精品福利| 网址你懂的国产日韩在线| 色视频www国产| 19禁男女啪啪无遮挡网站| 欧美区成人在线视频| 特大巨黑吊av在线直播| 九色成人免费人妻av| 一本综合久久免费| 18美女黄网站色大片免费观看| 精品人妻偷拍中文字幕| 亚洲内射少妇av| 校园春色视频在线观看| 99热只有精品国产|