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

    基于流固耦合的平板閘門動(dòng)水關(guān)閉過程數(shù)值模擬

    2021-09-18 11:48:40戴冰清茅澤育
    人民長(zhǎng)江 2021年8期
    關(guān)鍵詞:開度閘門水流

    戴冰清 茅澤育

    摘要: 閘門動(dòng)水關(guān)閉是一個(gè)復(fù)雜的物理過程,在此過程中的水動(dòng)力特性及閘門本身結(jié)構(gòu)特性易受到多種因素的綜合影響。針對(duì)平板閘門動(dòng)態(tài)關(guān)閉過程,采用動(dòng)網(wǎng)格技術(shù)、RNG kε紊流模型及VOF自由水面處理技術(shù)三者相結(jié)合的方法建立了三維非定常水氣兩相流模型,在此基礎(chǔ)上,進(jìn)一步結(jié)合流固耦合方法針對(duì)不同因素綜合影響下閘門的形變、等效應(yīng)力的變化規(guī)律進(jìn)行數(shù)值模擬研究,綜合探究來流流量、閘門底緣型式、閘門關(guān)閉速度等因素對(duì)平板閘門閉門過程的水動(dòng)力特性及閘門結(jié)構(gòu)特性的影響。研究結(jié)果表明:閘門最大形變出現(xiàn)在底緣處,最大等效應(yīng)力出現(xiàn)在門頂處;合理減小閘門關(guān)閉速度及上游來流流量、采用底緣前傾角較大的閘門等措施,有利于改善閘孔出流條件,減小閘門門體所受負(fù)壓的區(qū)域,并減小閘門形變及等效應(yīng)力。研究成果可為閘門的優(yōu)化設(shè)計(jì)及確定合理運(yùn)行工況提供科學(xué)依據(jù)。

    關(guān) 鍵 詞: 平板閘門; 水動(dòng)力特性; 結(jié)構(gòu)特性; 數(shù)值模擬; 流固耦合; 動(dòng)網(wǎng)格

    中圖法分類號(hào): ?TV135.4

    文獻(xiàn)標(biāo)志碼: ?A

    DOI: 10.16232/j.cnki.1001-4179.2021.08.033

    0 引 言

    閘門是水利樞紐安全運(yùn)行的關(guān)鍵部分之一,它能否安全運(yùn)行直接關(guān)系著泄水建筑物運(yùn)行的技術(shù)可行性和安全可靠性[1],也是一直受到關(guān)注的問題。閘門動(dòng)水關(guān)閉是一個(gè)復(fù)雜的物理過程,同時(shí)受到上游水頭、水流流速、下閘速度、閘門體型和底緣型式等多因素的綜合影響。國內(nèi)外對(duì)閘門失事原因的調(diào)查結(jié)果表明,橫向流、折沖水流及水躍沖擊等水流會(huì)對(duì)閘門的結(jié)構(gòu)造成一定的破壞。因此需要綜合探究閘門閉門過程的水動(dòng)力特性及閘門結(jié)構(gòu)特性的變化規(guī)律,為工程上閘門的優(yōu)化設(shè)計(jì)及確定合理運(yùn)行工況提供科學(xué)依據(jù)。

    近年來,學(xué)者們分別采用物理模型試驗(yàn)、原型觀測(cè)、數(shù)值計(jì)算等方法研究了閘門關(guān)閉過程中的水動(dòng)力學(xué)特性[2]以及影響水流水力特性的因素。如Willi等[3]采用半經(jīng)驗(yàn)法研究了溢洪道形狀及閘門開度對(duì)溢洪道中水流水力特征的影響;Marcou等[4]采用Lattice Boltzmann模型,研究了閘門關(guān)閉過程中的明渠水流特性;Masoud[5]認(rèn)為水閘的出水性能與水流的弗勞德數(shù),和上游水深與閘門開度的比值相關(guān);Marzieyh[6]等采用多普勒測(cè)速計(jì)測(cè)得了閘門附近水流速度的分布情況。南水北調(diào)工程建設(shè)期間,方神光等[7]模擬計(jì)算了同時(shí)改變渠首流量、閘門開度和分水口流量時(shí)渠道中的非恒定過渡過程;吳保生等[8]根據(jù)渠道水位、流量等信息來優(yōu)化渠道自動(dòng)化控制系統(tǒng)。還有一些學(xué)者采用數(shù)值模擬方法進(jìn)行了初步探索,劉昉等[9]將閘門關(guān)閉過程簡(jiǎn)化為閘門的單向運(yùn)動(dòng),再采用動(dòng)網(wǎng)格方法進(jìn)行數(shù)值模擬;文林森等[10]基于VOF方法,對(duì)某水電站事故閘門閉門過程進(jìn)行了初步數(shù)值模擬研究;郭桂禎等[11]就流體對(duì)結(jié)構(gòu)自振特性的影響規(guī)律進(jìn)行了初步研究。

    可以看出,前人對(duì)閘門動(dòng)水關(guān)閉過程的研究大多僅限于水流的水動(dòng)力學(xué)特性本身,而很少或沒有考慮水流對(duì)閘門門體本身結(jié)構(gòu)特性的影響。實(shí)際工程中,閘門動(dòng)水關(guān)閉過程中水氣兩相流交替變化復(fù)雜,水在流經(jīng)閘門底緣過程中由于固體邊界的影響,必然形成繞流彎曲流場(chǎng),水流分離現(xiàn)象不可避免造成閘門底部壓力梯度變化大,致使整個(gè)流域具有復(fù)雜的動(dòng)邊界問題。

    鑒于此,本文運(yùn)用FLUENT軟件,結(jié)合動(dòng)網(wǎng)格技術(shù),綜合考慮了不同入口流量、不同底緣型式、不同閘門關(guān)閉速度等因素的綜合影響,采用有限元數(shù)值模擬計(jì)算方法分析平板閘門閉門過程中的水動(dòng)力學(xué)特性變化。同時(shí)結(jié)合流固耦合方法,對(duì)閘門門體的形變、等效應(yīng)力變化規(guī)律進(jìn)行研究分析。

    1 數(shù)值模型及理論基礎(chǔ)

    有限元方法具有計(jì)算精度高、能適應(yīng)較為復(fù)雜的幾何形狀條件等優(yōu)點(diǎn),被視為一種具有廣闊應(yīng)用前景的工程分析手段。本文采用平板閘門有限元法展開分析,整個(gè)過程可以分為前處理(包括建立模型及劃分網(wǎng)格)、求解、后處理(包括采集處理并分析結(jié)果)3個(gè)部分。

    1.1 水流控制方程

    采用雷諾時(shí)均模型[12],不可壓縮水流的連續(xù)性方程和動(dòng)量方程分別如下:

    ρ t +div ρ u? =0 (1)

    ρu? t +div ρu u? =- p x + τxx x + τyx y + τzx z +Fx (2)

    ρv? t +div ρv u? =- p y + τxy x + τyy y + τzy z +Fy (3)

    ρw? t +div ρw u? =- p z + τxz x + τyz y + τzz z +Fz (4)

    式中:p是壓強(qiáng),τxx等是因分子黏性作用及湍流黏性影響而產(chǎn)生的作用在微元體表面上的黏性應(yīng)力τ的分量;Fx、Fy、Fz是微元體上的體積力,本文模型中體積力只有重力,z軸豎直向上,F(xiàn)x=0、Fy=0、Fz=-ρg。

    1.2 紊流模型

    采用RNG k-ε紊流模型計(jì)算結(jié)果具有更高的可信度[13]。其中k方程及ε方程分別為

    t? ρk + ?xi? ρkui = ?xj? αkueff k xj? +Gk-ρε? (5)

    t? ρε + ?xi? ρεui = ?xj? αεueff ε xj? +C*1ε ε k Gk-C2ερ ε2 k?? (6)

    μt=ρCμ k2 ε

    Gk=μt? ui xj + uj xi?? ui xj

    C*1ε=C1ε- η 1-η/η0? 1+βη3?? (7)

    η= 2Eij·Eij 1/2 k ε? (8)

    Eij= 1 2?? ui xj + uj xi?? (9)

    式中:Gk為壓力生成項(xiàng);ueff為有效黏性系數(shù);Eij為時(shí)均應(yīng)變率。

    Cμ=0.084 5,C1ε=1.42,C2ε=1.68

    ,αk=αε=1.393,η0=4.377,β=0.012。其他物理量含義參見文獻(xiàn)[14]。

    1.3 VOF模型

    水氣界面的自由表面跟蹤采用VOF方法[14]。假設(shè)在一個(gè)單元中,第q相流體體積分?jǐn)?shù)為αq,則有3種情況:αq=0表示單元網(wǎng)格中無第q相流體;αq=1表示單元網(wǎng)格中全為第q相流體;0<αq<1時(shí),表示通過求解水的容積分?jǐn)?shù)來追蹤水氣交界面。

    1.4 動(dòng)網(wǎng)格技術(shù)

    閘門動(dòng)水關(guān)閉過程的復(fù)雜性體現(xiàn)在水氣兩相流交替變化。當(dāng)水流流經(jīng)閘門底緣時(shí),由于固體邊界的影響,將形成繞流彎曲流場(chǎng)。閘門底緣處產(chǎn)生水流分離現(xiàn)象,從而導(dǎo)致閘門底部壓力梯度變大,整個(gè)流域具有復(fù)雜的動(dòng)邊界問題。對(duì)于這類問題,需要?jiǎng)泳W(wǎng)格模型實(shí)現(xiàn)網(wǎng)格劃分[15]。

    在進(jìn)行初始模型建立及網(wǎng)格劃分時(shí),同時(shí)采用邊界型函數(shù)或者UDF(用戶自定義函數(shù)),用以定義指定的邊界運(yùn)動(dòng)方式。將運(yùn)動(dòng)函數(shù)定義在所劃分的局部網(wǎng)格面或網(wǎng)格區(qū)域上,如果流場(chǎng)中包含運(yùn)動(dòng)與不運(yùn)動(dòng)兩種區(qū)域,則需要將它們組合在初始網(wǎng)格中以進(jìn)行識(shí)別。用于動(dòng)網(wǎng)格的更新有彈簧光順模型、動(dòng)態(tài)層模型及局部重劃模型。

    1.5 流固耦合理論

    依據(jù)數(shù)據(jù)傳遞方式的不同,流固耦合可以分為單向流固耦合和雙向流固耦合。由于平板閘門整體剛度較大、位移小,相對(duì)于流場(chǎng)中的變化可忽略不計(jì)[16],因此采用單向流固耦合方法進(jìn)行不同工況下平板閘門變形及等效應(yīng)力分析。單向流固耦合適用于相互作用中結(jié)構(gòu)位移場(chǎng)變化較小、對(duì)流場(chǎng)影響可以忽略的耦合計(jì)算,并且僅將流場(chǎng)結(jié)果數(shù)據(jù)傳遞給固體結(jié)構(gòu)仿真計(jì)算。先進(jìn)行流場(chǎng)計(jì)算,然后通過耦合面把流場(chǎng)結(jié)果數(shù)據(jù)傳遞給結(jié)構(gòu)場(chǎng)作為荷載或邊界條件,最后求解結(jié)構(gòu)位移場(chǎng)[17]。

    2 模型驗(yàn)證

    本文采用加拿大Mica電站進(jìn)水口閘門水力學(xué)試驗(yàn)資料對(duì)本文構(gòu)建的數(shù)值模型進(jìn)行驗(yàn)證[18]。Mica電站進(jìn)水口采用上游底緣傾角為30°的平板閘門,閘門門寬6.7 m,門高7.188 m,最大門厚1.143 m,啟閉機(jī)設(shè)計(jì)容量為248.1 t,閘門事故緊急閉門速度為6.1 m/min,事故閘門上有長(zhǎng)為14.5 m的斜喇叭進(jìn)口,縱向高程差為67 m。電站進(jìn)水口及引水管道系統(tǒng)的水工模型和事故閘門的幾何比尺λl=18,流量比尺λq=1 374.6,壓力比尺λp=18,時(shí)間比尺λt=4.24。在閘門門體上設(shè)置了20個(gè)測(cè)點(diǎn)。在試驗(yàn)工況H=71.5 m,Q=335 m3/s且閘門關(guān)閉速度為6.1 m/min時(shí),試驗(yàn)得到閘門底緣所受壓強(qiáng)隨閘門相對(duì)開度變化的曲線。

    圖1是考慮流固耦合作用后數(shù)值模擬計(jì)算結(jié)果與模型結(jié)果的對(duì)比。圖1(a)表示閘門底緣測(cè)點(diǎn)10所受壓強(qiáng)的試驗(yàn)值和數(shù)值模擬值,兩者變化規(guī)律基本一致。在水流阻滯作用下,閘門底緣所受壓強(qiáng)快速上升;隨著相對(duì)開度減小到0.8,水流逐步脫離閘門底緣,底緣所受壓強(qiáng)持續(xù)減小并在相對(duì)開度為0.2~0.3時(shí)出現(xiàn)最小值;隨著開度繼續(xù)減小,底緣所受壓強(qiáng)逐漸升高,關(guān)閉至底部時(shí)接近上游靜水壓強(qiáng)。圖1(b)為根據(jù)20個(gè)測(cè)點(diǎn)實(shí)測(cè)得到的壓強(qiáng)值繪制的閘門壓強(qiáng)分布圖,與圖1(c)中模擬閘門壓強(qiáng)等值線圖對(duì)比發(fā)現(xiàn),閘門典型開度下的試驗(yàn)與模擬所得的壓強(qiáng)分布規(guī)律基本一致,靠近底緣處開始呈逆壓梯度分布,下游面板處壓強(qiáng)降低為零值附近。驗(yàn)證結(jié)果表明本文采用的數(shù)值模型可以滿足平板閘門閉門過程的水動(dòng)力特性三維模擬。

    3 數(shù)值試驗(yàn)的模型建立及邊界條件

    3.1 計(jì)算區(qū)域網(wǎng)格劃分

    將模型簡(jiǎn)化為進(jìn)口段、閘室段和出口段組成的三維模型。計(jì)算尺寸:寬2 m,高2 m,閘前區(qū)域長(zhǎng)3 m,閘后區(qū)域長(zhǎng)15 m。閘門門寬2.5 m,門高2.2 m,最大厚度0.2 m,采用鋼制平面閘門型式(見圖2)。計(jì)算區(qū)域以四面體為單元進(jìn)行初始網(wǎng)格劃分,根據(jù)接近度和曲率,選擇最高關(guān)聯(lián)度,設(shè)置局部網(wǎng)格加密。以底緣前傾角為60°型式的閘門為例,初始狀態(tài)計(jì)算區(qū)域有限元計(jì)算網(wǎng)格總計(jì)約43萬個(gè),平均網(wǎng)格質(zhì)量為0.755,如圖3所示。

    在初始網(wǎng)格劃分基礎(chǔ)上,首先應(yīng)用ICEM CFD軟件進(jìn)行初始模型建立及網(wǎng)格劃分,同時(shí)采用UDF(用戶自定義函數(shù))定義指定的邊界運(yùn)動(dòng)方式為閘門以不同速度勻速下落,并且將運(yùn)動(dòng)函數(shù)定義在所劃分的局部網(wǎng)格面或網(wǎng)格區(qū)域上。

    本文采用動(dòng)網(wǎng)格中的動(dòng)態(tài)分層法及局部重劃法,以適應(yīng)閘門幾何和網(wǎng)格結(jié)構(gòu)的復(fù)雜性[19]。其中,動(dòng)態(tài)分層法的中心思想是:在邊界發(fā)生運(yùn)動(dòng)時(shí),如果緊鄰邊界的網(wǎng)格層高度增加到一定程度,就將其劃分為2個(gè)網(wǎng)格層;如果網(wǎng)格層高度降低到一定程度,就將緊鄰邊界的2個(gè)網(wǎng)格層合并為1層。而局部重劃法,則是用于將計(jì)算過程中畸變率或尺寸變化過大的網(wǎng)格集中在一起進(jìn)行局部網(wǎng)格的重新劃分。若重劃后的網(wǎng)格可以滿足畸變率及尺寸要求,則用新網(wǎng)格代替原來的網(wǎng)格,否則放棄重新劃分的結(jié)果[20]。將這兩種方法進(jìn)行結(jié)合,可以更大程度地保證網(wǎng)格質(zhì)量及計(jì)算精度,使閘門區(qū)網(wǎng)格更新后不發(fā)生明顯的畸變現(xiàn)象。

    3.2 數(shù)值計(jì)算方法

    設(shè)置計(jì)算域左端進(jìn)口條件為壓力進(jìn)口,水深為1.5 m;出口為壓力出口,參考大氣壓位置設(shè)置在自由水面處,并設(shè)置參考操作工質(zhì)密度為1.225 kg/m3。流場(chǎng)與閘門界面設(shè)置單向流固耦合邊界條件,壁面默認(rèn)為無滑移邊界條件。采用控制變量法進(jìn)行工況計(jì)算方案設(shè)計(jì),計(jì)算工況如表1所列。

    計(jì)算方法:選擇基于壓力求解器的瞬態(tài)計(jì)算方法;選擇用于非穩(wěn)態(tài)可壓縮或不可壓縮流體流場(chǎng)中求解壓力速度耦合關(guān)系的PISO算法;選擇PRESTO!格式空間離散化壓力方程;采用二階迎風(fēng)格式動(dòng)量方程。其余參數(shù)保持默認(rèn)值,并設(shè)置時(shí)間步長(zhǎng)為0.005 s,計(jì)算1 000步,瞬態(tài)場(chǎng)每一時(shí)間步的最大迭代次數(shù)為10。先將靜態(tài)閘門全開時(shí)的狀態(tài)按上述設(shè)置進(jìn)行計(jì)算,直到整個(gè)自由水面以下的計(jì)算流域充滿水(見圖4),再將最后時(shí)刻的case和data值導(dǎo)入動(dòng)網(wǎng)格計(jì)算case中作為初始化數(shù)據(jù),開始計(jì)算。

    3.3 流固耦合計(jì)算

    閘門材料為Q345C鋼,其彈性模量為2.1×105 MPa,泊松比為0.3,密度為7.85×103 kg/m3。由材料力學(xué)知識(shí)可知,由于金屬材料Q345C鋼為延性材料,通常以屈服形式失效,符合米塞斯屈服準(zhǔn)則,一般采用第四強(qiáng)度理論(又稱形狀改變比能理論,不論什么應(yīng)力狀態(tài),只要形狀改變比能達(dá)到極限值,便引起屈服)進(jìn)行校核。根據(jù)第四強(qiáng)度理論[20],等效應(yīng)力計(jì)算公式為

    σ = 1? 2??? σx-σy 2+ σy-σz 2+ σz-σx 2+6 τ2xy+τ2yz+τ2zx???? (10)

    實(shí)際工程中閘門通過啟閉機(jī)啟閉,僅在空間z軸上有位移,因此將閘門上端設(shè)置為固定端。將平底底緣及前傾角底緣型式閘門以四面體為單元?jiǎng)澐志W(wǎng)格。在設(shè)置的耦合面上加載來自流體的作用力,求解閘門總變形及等效應(yīng)力。

    4 ?平板閘門關(guān)閉過程的數(shù)值模擬及水動(dòng)力特性分析

    將閘門上端設(shè)置為固定端,底緣為自由端,閘門可視為懸臂梁模型,分析研究流固耦合作用下閘門關(guān)閉過程中的水動(dòng)力特性。

    4.1 結(jié)果計(jì)算

    對(duì)于閘底底坎為平頂堰,當(dāng)相對(duì)開度e/H≤0.65時(shí),為閘孔出流;e/H>0.65時(shí),為堰流。對(duì)于底坎為平頂堰的閘孔出流,有以下計(jì)算公式[21]:

    Qv=μbe 2gH0? (11)

    式中:e為閘門開度;H為閘門底坎頂?shù)阶杂杀砻嫔疃?μ為閘孔自由出流流量系數(shù),μ=0.60-0.76 e H ;b為閘孔寬度;H0為閘孔全水頭,H0=H+α0v0/(2g)。質(zhì)量流量與體積流量的換算,Qm=ρQv。

    理論計(jì)算值與數(shù)值模擬(工況6)得到的閘孔出流質(zhì)量流量值隨閘門開度的變化曲線如圖5所示。顯然,計(jì)算值與數(shù)值模擬值誤差較小。

    4.2 不同工況下的壓強(qiáng)分布規(guī)律

    本節(jié)控制閘門關(guān)閉速度為0.3 m/s,探究其他因素對(duì)閘門關(guān)閉過程的影響規(guī)律。圖6是平底底緣與前傾角60°底緣的閘門在相對(duì)開度分別為0.65和0.30的情況下計(jì)算區(qū)域內(nèi)的壓強(qiáng)分布。顯然,閘門底緣處,壓強(qiáng)呈大梯度逆壓梯度分布。開度較大時(shí)閘后的回流現(xiàn)象并不明顯,形成的漩渦對(duì)閘門底緣的壓強(qiáng)分布影響較小,底緣沒有形成明顯的負(fù)壓區(qū);當(dāng)開度較小時(shí),水流與閘底分離產(chǎn)生脫流現(xiàn)象,閘門后水流漩渦區(qū)明顯增大。這一區(qū)域處于水氣兩相流混合交替狀態(tài),閘門底緣產(chǎn)生射流擾動(dòng)作用,需要補(bǔ)入氣體,吸氣作用強(qiáng)烈,造成圖示閘門底緣處的負(fù)壓分布,從而使閘門底緣壓強(qiáng)呈明顯的逆梯度分布特征。

    由圖6還可以看到:當(dāng)Q=30 m3/s,v=0.3 m/s時(shí),動(dòng)水壓強(qiáng)峰值出現(xiàn)在閘門上端,負(fù)壓出現(xiàn)在閘門底緣。平底底緣型式閘門的關(guān)閉過程中,較小開度時(shí)水流分離點(diǎn)一般出現(xiàn)在上游,在閘門下部底緣產(chǎn)生空腔時(shí)無法及時(shí)從下游補(bǔ)充空氣,閘門底緣產(chǎn)生極不穩(wěn)定的負(fù)壓,工程中此時(shí)的動(dòng)水垂直力會(huì)呈現(xiàn)下吸力特性。而前傾角底緣型式的閘門設(shè)計(jì)更為合理,流線型輪廓的水流條件較好。閘門底緣的水流分離現(xiàn)象不明顯,幾乎沒有形成負(fù)壓區(qū)域,工程上可以減少底緣的空化現(xiàn)象。前傾角底緣閘門所受到的壓強(qiáng)峰值較小,且底緣傾角為60°及45°時(shí),在閘門處基本不產(chǎn)生負(fù)壓。因此,可以通過改變閘門底緣型式以改善閘孔出流條件,減小過流對(duì)閘門產(chǎn)生的影響,并且宜選用較大傾角。

    由圖7(a)和圖8(a)可見,最大壓強(qiáng)隨開度減小而增大。來流流量分別為18,30,36 m3/s時(shí),最大負(fù)壓分別為29.80,61.60,83.16 kPa。

    由圖7(b)和圖8(b)可見,閘門關(guān)閉前期,所受到的負(fù)壓處于較為穩(wěn)定的小范圍不規(guī)律變動(dòng)狀態(tài)。原因主要有:① 在來流流量不變的情況下,隨著閘門的關(guān)閉,過流斷面面積減小,流速增大,弗勞德數(shù)增大,水流本身處于強(qiáng)烈的紊動(dòng)狀態(tài),各個(gè)物理量及補(bǔ)氣條件受影響,則閘門底緣所受負(fù)壓產(chǎn)生不規(guī)律變動(dòng);② 由于計(jì)算機(jī)內(nèi)存及計(jì)算時(shí)間限制,網(wǎng)格的數(shù)量及質(zhì)量對(duì)計(jì)算精度產(chǎn)生了一些影響;③ 在模型建立時(shí),對(duì)邊界條件及初始條件的簡(jiǎn)化在一定程度上也影響著計(jì)算精度。相對(duì)開度減小到0.4時(shí),閘門所受最大負(fù)壓(圖中的最小壓強(qiáng))急劇增大。

    來流流量越大,閘門受到的負(fù)壓也越大,負(fù)壓峰值分別為-8.7,-16.9,-23.4 kPa。原因?yàn)椋簛砹髁髁吭酱?,斷面流速越大,在閘孔處斷面收縮時(shí)水流分離形成的渦區(qū)較大,水流分離現(xiàn)象越劇烈,而閘后補(bǔ)氣越困難,則造成底緣處承受負(fù)壓越大。

    4.3 不同閘門關(guān)閉速度對(duì)閘門的影響

    表2是在3個(gè)典型開度下,以平底底緣閘為研究對(duì)象,采用不同關(guān)閉速度時(shí)閘門所受的最大及最小壓強(qiáng)值。同一開度下,關(guān)閉速度為0.4 m/s時(shí),閘門所受壓強(qiáng)最大,0.2 m/s時(shí)最小。開度越大,關(guān)閉速度對(duì)壓強(qiáng)大小影響越大。各個(gè)關(guān)閉速度下的閘門所受負(fù)壓峰值相差不大。因此減小閘門的關(guān)閉速度,可以減小閘門承受的最大壓強(qiáng),但對(duì)閘門底緣承受的負(fù)壓影響不大。

    5 平板閘門結(jié)構(gòu)特性數(shù)值模擬及分析

    5.1 平板閘門變形情況分析

    由圖9~10可知:當(dāng)Q=30 m3/s,v=0.3 m/s,e/H=0.2時(shí),從頂端至底部閘門總變形逐漸增大,在底緣處結(jié)果偏離初始狀態(tài)距離最遠(yuǎn)。閘門的最大形變量隨開度減小而增大。開度相同時(shí),來流流量越大,閘門最大形變量越大。例如平底底緣閘門在開度為0.2時(shí),3種工況下閘門的最大形變值分別達(dá)到了3.25×10-4,7.96×10-4,10.96×10-4m。

    同一工況下,前傾角底緣型式閘門的最大變形區(qū)域及最大變形量均小于平底底緣型式閘門,說明前傾角底緣型式閘門具有更好的結(jié)構(gòu)穩(wěn)定性,且傾角為60°時(shí),閘門產(chǎn)生的變形量相對(duì)較小。

    5.2 平板閘門等效應(yīng)力分析

    由圖11~12可見:當(dāng)Q=30 m3/s,v=0.3 m/s,e/H=0.2時(shí),隨著開度的減小,閘門等效應(yīng)力值呈增大趨勢(shì)。主要原因有:① 動(dòng)水壓力隨開度的減小逐漸增大;② 動(dòng)水壓力對(duì)閘門的作用面積逐漸增大。最大等效應(yīng)力值出現(xiàn)在閘門頂部,并且往閘門底部逐漸減小。在相同開度時(shí),來流流量越大,閘門最大等效應(yīng)力值越大。例如,在開度0.2時(shí),3種工況下,平底底緣閘門的最大等效應(yīng)力值分別達(dá)到了5.85,14.93,20.68 MPa。

    60°前傾角底緣型式閘門的最大等效應(yīng)力值始終小于平底底緣型式閘門,更不容易被破壞。這與前文中前傾角底緣型式閘門的水流條件較好、不易形成負(fù)壓區(qū)域等優(yōu)點(diǎn)相結(jié)合,在工程中建議使用較大前傾角底緣型式閘門而非平底式。

    6 結(jié) 論

    閘門動(dòng)水關(guān)閉過程中涉及到復(fù)雜的水氣兩相流,整個(gè)計(jì)算域具有復(fù)雜的動(dòng)邊界問題。為了研究平板閘門閉門過程的水動(dòng)力特性及閘門結(jié)構(gòu)特性,本文采用動(dòng)網(wǎng)格技術(shù),應(yīng)用RNG kε紊流模型、VOF方法及流固耦合方法,建立了平板閘門動(dòng)水關(guān)閉過程的三維非定常水氣兩相流模型并應(yīng)用試驗(yàn)資料進(jìn)行了驗(yàn)證,進(jìn)一步進(jìn)行了數(shù)值模擬,主要研究結(jié)果如下。

    (1) 閘后出現(xiàn)的漩渦區(qū)可以減小水霧對(duì)閘門的侵蝕作用。在閘門底緣處,呈大壓強(qiáng)梯度變化并呈逆壓梯度分布現(xiàn)象。

    (2) 閘門最大形變出現(xiàn)在底緣處,最大等效應(yīng)力出現(xiàn)在門頂處。因此工程上可以針對(duì)閘門的不同部位采取措施,進(jìn)一步提高閘門的性能。

    (3) 合理減小閘門的關(guān)閉速度可以減小閘門承受的最大壓強(qiáng),但對(duì)閘門底緣承受的負(fù)壓無明顯影響。入口流量對(duì)流域及閘門也有一定影響,入口流量越大,閘門承受更大的負(fù)壓值,產(chǎn)生較大形變及等效應(yīng)力值。因此也需要合理控制上游來流流量使閘門安全運(yùn)行。

    (4) 前傾角底緣型式閘門具有水流條件較好、不易形成負(fù)壓區(qū)域、閘門形變及等效應(yīng)力較小等優(yōu)點(diǎn),在工程中建議使用較大前傾角(如60°)底緣型式閘門而非平底式。可改善閘孔出流條件,減小過流對(duì)閘門產(chǎn)生的影響。

    參考文獻(xiàn):

    [1] ?黃金林.平面閘門底緣型式及選擇[J].中國水利,2004(14):44-45.

    [2] 趙夢(mèng)麗.泄洪洞事故閘門動(dòng)水閉門水力及爬振特性研究[D].天津:天津大學(xué),2017.

    [3] WILLI H H,ROGER B.Plane gate on standard spillway[J].Journal of Hydraulic Engineering,1988,114(11):1390-1397.

    [4] MARCOU O,CHOPARD B,S EL YACOUBI,et al.Lattice boltzmann model for the simulation of flows in open channels with application to flows in a submerged sluice gate[J].Journal of Irrigation and Drainage Engineering,2010,136(12):809-822.

    [5] MASOUD G.Flow through side sluice gate[J].Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE),2003,129(6):458-463.

    [6] MARZIEYH E,MANOUCHEHR H,SAYED SAEID E.Flow characteristics of a sharp-crested side sluice gate[J].Journal of Irrigation and Drainage Engineering,2014,141(7):6014007.

    [7] 方神光,李玉榮,吳保生.大型輸水渠道閘前常水位的研究[J].水科學(xué)進(jìn)展,2008(1):68-71.

    [8] 吳保生,尚毅梓,崔興華,等.渠道自動(dòng)化控制系統(tǒng)及其運(yùn)行設(shè)計(jì)[J].水科學(xué)進(jìn)展,2008(5):746-755.

    [9] 劉昉,趙夢(mèng)麗,冷東升,等.不同底緣形式的平板閘門水力特性數(shù)值模擬[J].水利水電科技進(jìn)展,2017,37(5):46-50,77.

    [10] ?文林森,王才歡,楊偉,等.水工附環(huán)閘門閉門過程水力特性數(shù)值模擬研究[J].長(zhǎng)江科學(xué)院院報(bào),2017,34(10):68-73.

    [11] 郭桂禎.平板閘門垂向流激振動(dòng)特性與數(shù)值計(jì)算研究[D].天津:天津大學(xué),2011.

    [12] 王福軍.計(jì)算流體動(dòng)力學(xué)分析[M].北京:清華大學(xué)出版社,2004.

    [13] 胡揭玄.橫風(fēng)下公路橋梁與車輛系統(tǒng)的氣動(dòng)特性的數(shù)值模擬及試驗(yàn)研究[D].長(zhǎng)沙:長(zhǎng)沙理工大學(xué),2012.

    [14] 王福軍.計(jì)算流體動(dòng)力學(xué)分析—CFD軟件原理與應(yīng)用[M].北京:清華大學(xué)出版社,2004.

    [15] 郭澤宇.基于自適應(yīng)動(dòng)網(wǎng)格技術(shù)的類魚游動(dòng)水動(dòng)力及流場(chǎng)特性的數(shù)值研究[C]∥中國力學(xué)學(xué)會(huì)慶祝中國力學(xué)學(xué)會(huì)成立60周年大會(huì)論文集,2017.

    [16] 李明.弧形閘門動(dòng)力特性及流激振動(dòng)數(shù)值模擬[D].長(zhǎng)沙:長(zhǎng)沙理工大學(xué),2013.

    [17] 聞德蓀,李兆年,黃正華.工程流體力學(xué)(水力學(xué))[M].北京:高等教育出版社,2004.

    [18] 中國水利水電科學(xué)研究院.Mica電站進(jìn)水口閘門水力學(xué)試驗(yàn)研究[R].北京:中國水利水電科學(xué)研究院,2011.

    [19] 張斌,楊濤,豐志偉,等.網(wǎng)格質(zhì)量反饋的彈性體動(dòng)網(wǎng)格改進(jìn)[J].國防科技大學(xué)學(xué)報(bào),2018,40(1):10-16.

    [20] 甘慶明,付芳琴,張琴,等.基于動(dòng)網(wǎng)格的往復(fù)式壓縮機(jī)進(jìn)氣閥數(shù)值模擬[J].石油機(jī)械,2019,47(9):105-110,117.

    [21] 吳持恭.水力學(xué)[M].北京:高等教育出版社,2005.

    (編輯:胡旭東)

    引用本文:

    戴冰清,茅澤育.基于流固耦合的平板閘門動(dòng)水關(guān)閉過程數(shù)值模擬

    [J].人民長(zhǎng)江,2021,52(8):214-221.

    Numerical simulation of plane gate closing under running water based on

    fluid-structure interaction

    DAI Bingqing,MAO Zeyu

    ( Department of Hydraulic Engineering,Tsinghua University,Beijing 100084,China )

    Abstract:

    The gate closing under running water is a complex physical process.In this process,the hydrodynamic characteristics and the structural characteristics of the gate itself are affected by many factors.Aiming at the dynamic closing process of the plane gate,a 3D unsteady water-gas two-phase flow model was established by combining dynamic mesh technology,turbulence k-ε model and VOF free surface treatment technology.On this basis,combined with the method of fluid-structure interaction,the variation law of gate deformation and equivalent stress under the influence of different factors was numerically simulated,and the influence of inflow flow,gate bottom edge type,gate closing speed and other factors on the hydrodynamic characteristics and gate structure characteristics of plane gate closing process were comprehensively explored.The results showed that the maximum deformation of the gate occurred at the bottom edge,and the maximum equivalent stress occurred at the top of the gate.Reasonable reduction of gate closing speed and upstream inflow flow,and the gates with large bottom rake angle were conducive to improving the outflow conditions downstream,reducing the negative pressure area of gate body,and reducing the gate deformation and equivalent stress.The above research results can provide a scientific basis for the optimal design of the gate and the determination of reasonable operating conditions.

    Key words:

    plane gate;hydrodynamic characteristics;structural characteristics;numerical simulation;fluid-structure interaction;dynamic mesh

    猜你喜歡
    開度閘門水流
    我終于會(huì)過地鐵閘門了
    哪股水流噴得更遠(yuǎn)
    能俘獲光的水流
    掘進(jìn)機(jī)用截止閥開度對(duì)管路流動(dòng)性能的影響
    增大某車型車門開度的設(shè)計(jì)方法
    北京汽車(2021年2期)2021-05-07 03:56:26
    燃燒器二次風(fēng)擋板開度對(duì)爐內(nèi)燃燒特性的影響
    我只知身在水中,不覺水流
    文苑(2020年6期)2020-06-22 08:41:56
    把住醫(yī)?;鹚亻l門
    找準(zhǔn)入口,打開思路的閘門
    我國專家發(fā)現(xiàn)藥物成癮新機(jī)制,找到“閘門”開關(guān)
    欧美三级亚洲精品| 久久久久久久久久久免费av| 免费av观看视频| 日韩视频在线欧美| 亚洲三级黄色毛片| 午夜免费男女啪啪视频观看| 欧美日本视频| 女人久久www免费人成看片| 美女脱内裤让男人舔精品视频| av福利片在线观看| 在线看a的网站| 精品国产乱码久久久久久小说| 欧美性猛交╳xxx乱大交人| 丝袜脚勾引网站| 99热这里只有是精品50| .国产精品久久| 街头女战士在线观看网站| 国产精品人妻久久久影院| 丰满少妇做爰视频| 最近最新中文字幕大全电影3| 精品酒店卫生间| 少妇人妻一区二区三区视频| av国产精品久久久久影院| 啦啦啦啦在线视频资源| 亚洲,欧美,日韩| freevideosex欧美| 在线观看美女被高潮喷水网站| 欧美xxxx性猛交bbbb| 久久精品熟女亚洲av麻豆精品| 免费观看在线日韩| 各种免费的搞黄视频| 丰满少妇做爰视频| 日韩视频在线欧美| 亚洲精品乱码久久久久久按摩| 2018国产大陆天天弄谢| 欧美xxxx性猛交bbbb| 久久这里有精品视频免费| 国内精品美女久久久久久| 国产精品人妻久久久影院| 男插女下体视频免费在线播放| 亚洲国产精品国产精品| 全区人妻精品视频| 久久久国产一区二区| 禁无遮挡网站| 国产欧美日韩精品一区二区| 亚洲精品aⅴ在线观看| 中文精品一卡2卡3卡4更新| 亚洲精品国产av成人精品| 尤物成人国产欧美一区二区三区| 国产精品一区二区在线观看99| 又大又黄又爽视频免费| 一区二区三区免费毛片| 久久久久久久亚洲中文字幕| 成人毛片a级毛片在线播放| 五月玫瑰六月丁香| 国产精品久久久久久精品电影小说 | 国产老妇女一区| 亚洲欧美成人综合另类久久久| 色播亚洲综合网| 18禁动态无遮挡网站| 九草在线视频观看| 晚上一个人看的免费电影| 亚洲在久久综合| 国产精品无大码| 久久人人爽av亚洲精品天堂 | 欧美高清成人免费视频www| 久久精品人妻少妇| 国产综合精华液| 干丝袜人妻中文字幕| 91狼人影院| 黄片无遮挡物在线观看| 午夜福利视频精品| 别揉我奶头 嗯啊视频| 国产黄片美女视频| 亚洲国产精品专区欧美| 毛片女人毛片| 18+在线观看网站| 亚洲精品久久久久久婷婷小说| 91久久精品国产一区二区三区| 秋霞在线观看毛片| 欧美+日韩+精品| 少妇人妻一区二区三区视频| 一区二区三区乱码不卡18| 乱码一卡2卡4卡精品| 一级毛片黄色毛片免费观看视频| 啦啦啦在线观看免费高清www| 亚洲四区av| 黄色欧美视频在线观看| 国产亚洲5aaaaa淫片| 亚洲av一区综合| 国产一级毛片在线| 国产有黄有色有爽视频| 激情五月婷婷亚洲| 国产v大片淫在线免费观看| 少妇高潮的动态图| 久久99热这里只有精品18| 亚洲国产最新在线播放| 成人综合一区亚洲| 国产精品福利在线免费观看| 一个人观看的视频www高清免费观看| 亚洲不卡免费看| 国产男女超爽视频在线观看| 国产成人精品婷婷| 青春草国产在线视频| 欧美人与善性xxx| 99热国产这里只有精品6| 王馨瑶露胸无遮挡在线观看| tube8黄色片| 大话2 男鬼变身卡| 我要看日韩黄色一级片| 亚洲不卡免费看| 肉色欧美久久久久久久蜜桃 | av在线亚洲专区| 深夜a级毛片| 久久6这里有精品| 国产在线一区二区三区精| 欧美区成人在线视频| 爱豆传媒免费全集在线观看| 国产精品一二三区在线看| 伊人久久精品亚洲午夜| 七月丁香在线播放| 建设人人有责人人尽责人人享有的 | 99热这里只有是精品50| 国产精品一区二区性色av| 精品久久久久久久久亚洲| 少妇丰满av| 大话2 男鬼变身卡| 亚洲精品久久久久久婷婷小说| 在线观看av片永久免费下载| 视频区图区小说| 欧美性猛交╳xxx乱大交人| 2021少妇久久久久久久久久久| 久久久久久久国产电影| 亚洲欧美成人精品一区二区| 欧美区成人在线视频| 国产精品偷伦视频观看了| 大片电影免费在线观看免费| 91精品伊人久久大香线蕉| 久久久久久久大尺度免费视频| 国产成人午夜福利电影在线观看| 亚洲精品日韩av片在线观看| 99久国产av精品国产电影| 99久久精品一区二区三区| 亚洲欧美精品专区久久| 一区二区三区免费毛片| 尾随美女入室| 在线观看一区二区三区| 狂野欧美激情性bbbbbb| 天天一区二区日本电影三级| 国产男女内射视频| 中文字幕免费在线视频6| 99精国产麻豆久久婷婷| 亚洲av一区综合| 久久精品国产鲁丝片午夜精品| 纵有疾风起免费观看全集完整版| 91在线精品国自产拍蜜月| 日日啪夜夜爽| 中文欧美无线码| 18禁在线无遮挡免费观看视频| 天堂网av新在线| 午夜日本视频在线| 久久97久久精品| 国产高清三级在线| 亚洲精品自拍成人| 禁无遮挡网站| www.av在线官网国产| 一本—道久久a久久精品蜜桃钙片 精品乱码久久久久久99久播 | 三级男女做爰猛烈吃奶摸视频| 日韩大片免费观看网站| 亚洲成人一二三区av| av福利片在线观看| 综合色丁香网| 联通29元200g的流量卡| 精品久久久久久电影网| 麻豆精品久久久久久蜜桃| 国产精品.久久久| 亚洲自偷自拍三级| 国产爽快片一区二区三区| 三级国产精品欧美在线观看| 91久久精品国产一区二区成人| 亚洲av成人精品一区久久| 自拍偷自拍亚洲精品老妇| 好男人视频免费观看在线| 国产伦精品一区二区三区四那| 亚洲精品成人久久久久久| 亚洲最大成人手机在线| 国产毛片a区久久久久| 亚州av有码| 天堂俺去俺来也www色官网| 久久女婷五月综合色啪小说 | 老师上课跳d突然被开到最大视频| 国产永久视频网站| 欧美成人一区二区免费高清观看| a级毛片免费高清观看在线播放| 免费看av在线观看网站| 国产免费福利视频在线观看| 午夜免费鲁丝| av专区在线播放| 日日啪夜夜爽| 亚洲精品视频女| 免费观看无遮挡的男女| 欧美一级a爱片免费观看看| 国产精品无大码| 日本一二三区视频观看| 99久久九九国产精品国产免费| 自拍欧美九色日韩亚洲蝌蚪91 | 久久久久精品性色| 成人亚洲欧美一区二区av| 97在线视频观看| 国产久久久一区二区三区| 22中文网久久字幕| 丰满少妇做爰视频| 国产一区二区三区综合在线观看 | 久久久久网色| 国产精品一区二区三区四区免费观看| a级毛片免费高清观看在线播放| 别揉我奶头 嗯啊视频| 亚洲av中文av极速乱| 91精品伊人久久大香线蕉| 久热这里只有精品99| 国内揄拍国产精品人妻在线| 国产亚洲一区二区精品| 欧美精品一区二区大全| 性色avwww在线观看| 欧美高清成人免费视频www| 国产成人免费观看mmmm| 欧美极品一区二区三区四区| 国产精品成人在线| 少妇 在线观看| 色网站视频免费| 嫩草影院精品99| 国产高清三级在线| 3wmmmm亚洲av在线观看| 亚洲精品第二区| 国产熟女欧美一区二区| 国产黄片视频在线免费观看| 另类亚洲欧美激情| 亚洲国产欧美在线一区| 国产乱人偷精品视频| 蜜臀久久99精品久久宅男| 18禁裸乳无遮挡免费网站照片| 免费av观看视频| 18禁裸乳无遮挡动漫免费视频 | 搡女人真爽免费视频火全软件| 国产中年淑女户外野战色| 国产有黄有色有爽视频| 性插视频无遮挡在线免费观看| 在线播放无遮挡| 成人国产av品久久久| 国产黄色免费在线视频| 简卡轻食公司| 亚洲电影在线观看av| 青青草视频在线视频观看| 国产日韩欧美在线精品| 精品久久久久久电影网| 欧美3d第一页| 亚洲精品视频女| 午夜福利视频1000在线观看| 国产精品一区二区在线观看99| 欧美性感艳星| 免费av观看视频| 亚洲国产欧美在线一区| 观看免费一级毛片| 欧美bdsm另类| 97精品久久久久久久久久精品| 国产老妇女一区| 日本爱情动作片www.在线观看| 成人午夜精彩视频在线观看| 黑人高潮一二区| 成人亚洲精品av一区二区| 亚洲av中文av极速乱| 日本欧美国产在线视频| 日本免费在线观看一区| 能在线免费看毛片的网站| 人人妻人人爽人人添夜夜欢视频 | 天天躁夜夜躁狠狠久久av| 插阴视频在线观看视频| 亚洲av二区三区四区| 欧美日韩视频高清一区二区三区二| 久久精品熟女亚洲av麻豆精品| 少妇人妻久久综合中文| 久久精品国产亚洲av天美| 国产精品无大码| 乱系列少妇在线播放| 最近中文字幕2019免费版| 久久综合国产亚洲精品| 日日摸夜夜添夜夜添av毛片| 久久久久久久午夜电影| 国产探花在线观看一区二区| 亚洲av福利一区| 99久久人妻综合| av免费观看日本| 在线免费观看不下载黄p国产| 夜夜看夜夜爽夜夜摸| 国产免费视频播放在线视频| 国产毛片a区久久久久| 欧美性感艳星| 91狼人影院| av一本久久久久| 欧美一级a爱片免费观看看| 少妇人妻久久综合中文| 看十八女毛片水多多多| 男女边吃奶边做爰视频| 日韩欧美一区视频在线观看 | 在线观看一区二区三区| 午夜精品国产一区二区电影 | 国产精品久久久久久精品电影小说 | 国产高潮美女av| 亚洲高清免费不卡视频| 亚洲经典国产精华液单| 亚洲婷婷狠狠爱综合网| 日韩,欧美,国产一区二区三区| 中国国产av一级| 人妻 亚洲 视频| 国产成人精品婷婷| 日韩电影二区| 国产欧美另类精品又又久久亚洲欧美| 男女那种视频在线观看| av免费在线看不卡| 丝袜喷水一区| av又黄又爽大尺度在线免费看| 成人亚洲精品av一区二区| 黄片无遮挡物在线观看| 国产精品久久久久久精品古装| 午夜福利在线在线| 亚洲精品成人av观看孕妇| 久久久久久九九精品二区国产| 在现免费观看毛片| 亚洲精品aⅴ在线观看| 国内精品宾馆在线| 我的老师免费观看完整版| 春色校园在线视频观看| 亚洲电影在线观看av| 97热精品久久久久久| 黄色配什么色好看| 国产成人a∨麻豆精品| 日韩国内少妇激情av| 免费大片18禁| 欧美激情国产日韩精品一区| 色视频在线一区二区三区| av在线亚洲专区| 我要看日韩黄色一级片| 亚洲欧美成人精品一区二区| 偷拍熟女少妇极品色| 久久久欧美国产精品| 国产av不卡久久| 免费少妇av软件| 免费电影在线观看免费观看| 在线观看人妻少妇| 亚洲国产欧美在线一区| 欧美97在线视频| 亚州av有码| av在线老鸭窝| 大片电影免费在线观看免费| a级一级毛片免费在线观看| av国产精品久久久久影院| 日韩在线高清观看一区二区三区| 欧美性猛交╳xxx乱大交人| 国产精品人妻久久久久久| 91久久精品国产一区二区三区| 80岁老熟妇乱子伦牲交| 国产一区二区亚洲精品在线观看| 免费看不卡的av| 久久ye,这里只有精品| 亚洲人成网站在线观看播放| 成人漫画全彩无遮挡| kizo精华| 日本av手机在线免费观看| 搞女人的毛片| 亚洲,一卡二卡三卡| 亚洲最大成人中文| 91久久精品国产一区二区成人| 国产精品一区二区三区四区免费观看| 嫩草影院入口| 免费观看av网站的网址| 老司机影院成人| 寂寞人妻少妇视频99o| 99久久精品国产国产毛片| 欧美xxⅹ黑人| 国产精品久久久久久精品电影小说 | 亚洲精品亚洲一区二区| 久久精品国产亚洲网站| av播播在线观看一区| 国产亚洲一区二区精品| 欧美bdsm另类| 欧美高清性xxxxhd video| 26uuu在线亚洲综合色| 欧美极品一区二区三区四区| 美女国产视频在线观看| 日本一本二区三区精品| 亚洲国产欧美人成| 色5月婷婷丁香| 国产精品精品国产色婷婷| 久久久久久久久久成人| 91精品伊人久久大香线蕉| av福利片在线观看| 热99国产精品久久久久久7| 精品熟女少妇av免费看| 亚洲国产色片| 只有这里有精品99| 婷婷色av中文字幕| 日本wwww免费看| 精品亚洲乱码少妇综合久久| 18禁裸乳无遮挡动漫免费视频 | 久久人人爽av亚洲精品天堂 | 视频区图区小说| 伊人久久国产一区二区| 国产亚洲91精品色在线| 97超碰精品成人国产| 成人漫画全彩无遮挡| 免费看光身美女| 久久97久久精品| 人人妻人人澡人人爽人人夜夜| 丰满少妇做爰视频| 国产高潮美女av| 久久精品夜色国产| 免费高清在线观看视频在线观看| 极品教师在线视频| 又大又黄又爽视频免费| 国产高清三级在线| 丰满乱子伦码专区| 伦精品一区二区三区| 亚洲色图综合在线观看| 九九在线视频观看精品| 国产精品一区www在线观看| 日韩av在线免费看完整版不卡| 少妇熟女欧美另类| 精品国产乱码久久久久久小说| 尤物成人国产欧美一区二区三区| 国产亚洲91精品色在线| 国产精品不卡视频一区二区| 久久精品久久精品一区二区三区| 99视频精品全部免费 在线| 五月天丁香电影| 男女边摸边吃奶| 又爽又黄无遮挡网站| 看非洲黑人一级黄片| 别揉我奶头 嗯啊视频| 国产精品人妻久久久影院| 成人午夜精彩视频在线观看| 久久99热这里只频精品6学生| 国产毛片a区久久久久| 日韩欧美精品免费久久| 我的老师免费观看完整版| 一级黄片播放器| 黄片无遮挡物在线观看| 精品亚洲乱码少妇综合久久| 精品久久久久久久久av| 亚洲在久久综合| 性插视频无遮挡在线免费观看| 欧美三级亚洲精品| 麻豆乱淫一区二区| 亚洲欧美成人综合另类久久久| 男女那种视频在线观看| 搡老乐熟女国产| 麻豆成人av视频| 午夜福利视频精品| 三级国产精品欧美在线观看| 色视频www国产| 日韩一本色道免费dvd| 亚洲精品日本国产第一区| 久久久久久九九精品二区国产| 99热这里只有是精品在线观看| 丰满少妇做爰视频| 久久久久久久久久成人| 搡老乐熟女国产| 免费av观看视频| 日韩强制内射视频| 人妻 亚洲 视频| 亚洲内射少妇av| 国产日韩欧美在线精品| 高清日韩中文字幕在线| 亚洲图色成人| 精品国产乱码久久久久久小说| 成人二区视频| 国产精品国产三级专区第一集| 亚洲人与动物交配视频| 亚洲av欧美aⅴ国产| 国内精品宾馆在线| 欧美日韩一区二区视频在线观看视频在线 | 免费观看性生交大片5| 欧美激情国产日韩精品一区| 久久久久久久午夜电影| 国产探花在线观看一区二区| 亚洲av福利一区| 国产精品麻豆人妻色哟哟久久| 九九久久精品国产亚洲av麻豆| 亚洲精品乱码久久久v下载方式| 少妇人妻久久综合中文| 80岁老熟妇乱子伦牲交| 卡戴珊不雅视频在线播放| 国内精品美女久久久久久| 韩国高清视频一区二区三区| 久久久久久久久大av| 精品久久久久久久人妻蜜臀av| 亚洲,一卡二卡三卡| 国产成人精品久久久久久| 性色av一级| 少妇丰满av| 男女无遮挡免费网站观看| 综合色丁香网| 天堂中文最新版在线下载 | 交换朋友夫妻互换小说| 18禁动态无遮挡网站| 大码成人一级视频| 久久精品国产亚洲av涩爱| 美女cb高潮喷水在线观看| 免费高清在线观看视频在线观看| 尤物成人国产欧美一区二区三区| 综合色av麻豆| 亚洲经典国产精华液单| 少妇 在线观看| 亚洲av成人精品一二三区| 美女被艹到高潮喷水动态| 国产高清有码在线观看视频| 我的老师免费观看完整版| 嘟嘟电影网在线观看| 各种免费的搞黄视频| 国产成人一区二区在线| 日韩免费高清中文字幕av| av在线观看视频网站免费| 精品人妻一区二区三区麻豆| 91久久精品国产一区二区三区| 成人亚洲精品一区在线观看 | 蜜桃久久精品国产亚洲av| 日本猛色少妇xxxxx猛交久久| 国产探花极品一区二区| 五月天丁香电影| 亚洲成人中文字幕在线播放| 精品一区二区三卡| 最近中文字幕高清免费大全6| 国产免费福利视频在线观看| 国产精品久久久久久久久免| 中文精品一卡2卡3卡4更新| 十八禁网站网址无遮挡 | 国模一区二区三区四区视频| 汤姆久久久久久久影院中文字幕| 老师上课跳d突然被开到最大视频| 我要看日韩黄色一级片| 在线天堂最新版资源| 天堂俺去俺来也www色官网| 中国三级夫妇交换| 在线精品无人区一区二区三 | 黄色欧美视频在线观看| av在线播放精品| 纵有疾风起免费观看全集完整版| 欧美极品一区二区三区四区| 国产一区二区亚洲精品在线观看| 免费av不卡在线播放| 免费在线观看成人毛片| 最近最新中文字幕免费大全7| 亚洲欧美成人精品一区二区| 男插女下体视频免费在线播放| 黄片无遮挡物在线观看| 夜夜爽夜夜爽视频| 日本爱情动作片www.在线观看| 成人毛片60女人毛片免费| 如何舔出高潮| 一边亲一边摸免费视频| av免费观看日本| 国产成人精品一,二区| av在线播放精品| 天天躁夜夜躁狠狠久久av| 久久久久久久久久久丰满| 欧美性猛交╳xxx乱大交人| 成年女人看的毛片在线观看| 香蕉精品网在线| 亚洲国产欧美人成| 久久久久久久久久人人人人人人| 韩国av在线不卡| 97超视频在线观看视频| 久久久久国产网址| 久久久久国产精品人妻一区二区| 精品一区二区免费观看| 人人妻人人爽人人添夜夜欢视频 | 又粗又硬又长又爽又黄的视频| 国产一区二区三区综合在线观看 | 久久久欧美国产精品| 一区二区三区乱码不卡18| 精品人妻视频免费看| 黄色欧美视频在线观看| 久久韩国三级中文字幕| 亚洲精品视频女| 真实男女啪啪啪动态图| 国产av国产精品国产| 国产乱人视频| 国产乱人偷精品视频| 少妇裸体淫交视频免费看高清| 纵有疾风起免费观看全集完整版| av线在线观看网站| 青青草视频在线视频观看| 边亲边吃奶的免费视频| 高清在线视频一区二区三区| 视频中文字幕在线观看| 22中文网久久字幕| 国产综合懂色| 国精品久久久久久国模美| 亚洲激情五月婷婷啪啪| videos熟女内射| 国产黄片视频在线免费观看| 男女边摸边吃奶| 免费黄网站久久成人精品| 国产成年人精品一区二区| 国产午夜精品一二区理论片| 青青草视频在线视频观看| 黄色日韩在线| 男人和女人高潮做爰伦理| 午夜日本视频在线| 在线亚洲精品国产二区图片欧美 | 男插女下体视频免费在线播放| 在线观看av片永久免费下载| 嘟嘟电影网在线观看| 亚洲一级一片aⅴ在线观看| 午夜免费男女啪啪视频观看| 免费观看无遮挡的男女|