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

    基于LBM-DEM 耦合方法的突水潰砂運移規(guī)律研究

    2021-04-17 09:12:00郭世儒劉德俊許軍策
    煤炭科學技術(shù) 2021年2期
    關(guān)鍵詞:砂粒突水砂層

    浦 海,郭世儒,劉德俊,許軍策,王 健

    (1.中國礦業(yè)大學 深部巖土力學與地下工程國家重點實驗室,江蘇 徐州 221116;2. 新疆工程學院 采礦工程與地質(zhì)學院,新疆 烏魯木齊 830091;3.中國礦業(yè)大學 力學與土木工程學院,江蘇 徐州 221116)

    0 引 言

    隨著我國煤炭開發(fā)基地的戰(zhàn)略性西移,西部礦區(qū)煤層開采時因其獨特的地質(zhì)條件和高強度開采方式涌現(xiàn)出許多工程問題,其中在采掘工作面頻發(fā)的突水潰砂災害已成為影響煤礦安全生產(chǎn)的關(guān)鍵問題之一。 高強度開采過程中的擾動促使上覆含水層底部產(chǎn)生裂隙,在裂隙發(fā)育形成導水通道后進一步貫通厚松散砂層使水砂混合物涌入工作面是突水潰砂災害形成的主要方式之一。 在形成導通裂隙后,水砂運移的方式將決定礦井的危害程度,迫切需要研究水砂運移規(guī)律,尤其是在形成裂隙通道時,水砂在其內(nèi)部的流動規(guī)律及對通道的擴張影響[1-3]。

    近年來,眾多學者針對突水潰砂問題在理論分析、相似試驗模擬、基礎(chǔ)試驗驗證、數(shù)值計算模擬等方面做出大量的研究成果,為深入了解突水潰砂發(fā)生機理與災害防控提供了有意義的參考[4-12]。 基于工程分析角度,張桂民等[13]分析探討了煤礦鉆孔誘發(fā)的突水潰砂機理,提出了防治和治理突水潰砂的措施。 許延春等[14]研究了工作面突水潰砂機制,并提出相應的預防措施。 張勃陽等[15]采用流體力學方法,分析煤礦潰砂事故,提出了鉆孔潰砂災害的通用計算方法。 基于試驗研究角度,張凱等[16]根據(jù)實際工程情況進行了涌砂判據(jù)的試驗研究,通過滲流試驗研究了砂巖的水力特性,得出涌砂發(fā)生的臨界條件。 李宏杰等[17]進行了突水防砂的實例研究,根據(jù)試驗結(jié)果提出了分層開采、降低上下分層開采厚度等優(yōu)化設(shè)計和技術(shù)措施。 袁奇[18]分析了陜北礦區(qū)潰砂主要來源,并設(shè)計整套試驗裝置分析突水潰砂現(xiàn)象,攝錄了不同砂粒粒徑、砂厚和水頭高度情況下潰砂剖面特征,直觀展示了突水潰砂過程。 楊鑫等[19]采用自主研發(fā)的水砂兩相高速滲流試驗設(shè)備,對西部榆橫礦區(qū)風積砂含水層的起動、潰砂規(guī)律及水砂兩相高速滲流進行室內(nèi)試驗研究。 提出西部典型礦區(qū)風積砂含水層突水潰砂災害的臨界判據(jù),并分析潰砂過程中的能量傳遞機制。 楊偉峰等[20]設(shè)計了水砂混合涌流裝置,模擬了覆巖裂隙通道水砂混合突水啟動過程及遷移過程,分析了裂縫通道中水砂混合流的突水特征。 高炳倫[21]針對工程中裂隙潰砂問題研制了裂隙潰砂試驗系統(tǒng),實現(xiàn)了單一平直裂隙條件下裂隙潰砂過程的試驗模擬。 基于模擬研究角度,王迎超等[22]采用DEM 數(shù)值模擬方法建立了固-水兩相流動的三維模型,直觀地展現(xiàn)了突水潰砂時水砂運移過程。 梁艷坤等[23-24]采用分形理論分析了垮落帶破碎巖體的塊度分布特征,并利用離散元方法(DEM) 根據(jù)破碎巖體的塊度分布生成不規(guī)則碎石,建立了垮落帶破碎巖體潰砂數(shù)值模型。

    目前針對突水潰砂問題的研究并不透徹,主要原因是地下空間為封閉不可觀測系統(tǒng),地層水砂分布、承壓含水層壓力、巖層裂隙網(wǎng)絡(luò)分布等因素難以觀測。 加之水砂性質(zhì)的多樣,進行相似模擬試驗時選取相似材料過程繁瑣,通過試驗手段對突水潰砂危害性評估有很大的難度。 數(shù)值模擬因其可重復性和直觀性等特點,在研究突水潰砂災害時有明顯的優(yōu)勢。

    突水潰砂數(shù)值模擬中,重點在于針對砂粒及流體選取合適的計算模型,能夠?qū)烧哌M行耦合以得出符合實際的運動過程。 傳統(tǒng)的研究方法中將砂看成另一種形式的流相,水砂混合流動看成兩相流問題,或者將砂簡化為質(zhì)點,忽略顆粒之間的碰撞及翻轉(zhuǎn),在分析突水潰砂下水砂在裂隙網(wǎng)絡(luò)中的運移規(guī)律時達不到真實效果。 為使模擬更符合實際問題,筆者采用格子Boltzmann 方法(Lattice Boltzmann Method,LBM) 模 擬 液 相 流 動, 用 離 散 單 元 法(Discrete Element Method,DEM)構(gòu)建球形砂粒,并采用浸入移動邊界法(Immersed Moving Boundary Method, IMB)方法實現(xiàn)流體和砂粒間的耦合,由此模擬砂粒在流場中的流動情況,有效地分析突水潰砂下水砂兩相流動運移規(guī)律。

    1 突水潰砂分析模型

    西部礦區(qū)煤層開采時突水潰砂災害的演變過程如圖1 所示。 西部礦區(qū)煤層開采時劇烈的擾動使上層覆巖形成大量裂隙網(wǎng)絡(luò),當擾動加劇時將促使裂隙網(wǎng)絡(luò)貫通巷道和上層含砂水層,水砂混合流將通過裂隙大量涌入生產(chǎn)工作面,形成突水潰砂災害,帶來巨大的經(jīng)濟損失和人員傷亡。 由于西部礦區(qū)地表存在風積砂,若擾動加劇上層含砂水層形成與地表的貫通裂隙,將有源源不斷的砂供給,進一步加劇突水潰砂的災害性。 為明確突水潰砂的成災原因及危害程度,重點研究圖1 中紅圈位置水砂混合流在突發(fā)裂隙開口情況下的運移情況。 將實際突水潰砂時含砂水層中水砂進入破裂巖體這一情況簡化為恒定壓力下的單裂隙開口通道水砂運移數(shù)值計算模型如圖2 所示。

    圖1 西部礦區(qū)突水潰砂示意Fig.1 Schematic of water and sand inrush in western mining area

    圖2 突水潰砂簡化計算模型Fig.2 Simplified analysis model of water and send inrush

    矩形通道中寬度為d的喉道口表征突水潰砂發(fā)生時突然增大的裂隙開口,邊界壓力P表征地下含砂水層中的水力壓力情況。 水分布在開口上方的矩形區(qū)域,密度為ρw、動力黏度為υ。 砂粒用半徑為r、密度為ρs、厚度為h、數(shù)量為N的球形顆粒群表征。 由于砂粒分布區(qū)域固定,且分析時以自然堆積方式排布,砂層厚度和砂粒數(shù)量呈正比,僅取1 項作為變量即可。 以P、d、h作為變量,從數(shù)值模擬角度分析水砂混合流在單裂隙下的運移規(guī)律。

    2 LBM-DEM 耦合方法及模擬方案

    LBM-DEM 耦合方法在模擬中對顆粒進行直接建模,顆粒擁有物理邊界和確切的體積,流體假定為連續(xù)相并構(gòu)建小于流場中最小渦長尺寸的網(wǎng)格進行直接數(shù)值模擬。 隨著越來越多的學者不斷完善基于LBM-DEM 方法下的液固耦合理論,此方法能夠很好的分析液固耦合時顆粒流動和流場分布情況,并開始廣泛應用于生物、醫(yī)療等多領(lǐng)域中[25]。

    2.1 格子Boltzmann 方法(LBM)

    格子Boltzmann 方法是應用非連續(xù)介質(zhì)思想求解流體力學問題的方法。 其基本思想是把流體看成由許多只有質(zhì)量沒有體積的微小粒子組成,在時間、空間上離散,按照符合特定基本物理規(guī)律的演化規(guī)則碰撞,并沿網(wǎng)格線在節(jié)點之間遷移,利用分布函數(shù)的統(tǒng)計平均得到系統(tǒng)的宏觀物理參數(shù)[26]。

    借助于離散的單松弛時間格子Boltzmann -BGK 方程描述流體相互作用,即

    其中;x為粒子空間位置矢量;ci為粒子速度矢量;t為時間;δt為時間步長;fi x,t( ) 為粒子分布函數(shù);Fi x,t( ) 為離散的外力項;Ωi為離散的碰撞算子,可表達式為

    其中:τc為2 次碰撞間隔時間;u為流體宏觀速度矢量;fi,equ,ρ( ) 為離散后的平衡態(tài)速度分布函數(shù),表達式為

    其中:ωi為權(quán)重系數(shù);cs為格子聲速,表達式為

    式中:R為氣體常數(shù);T為熱力學溫度;D為空間維數(shù)。

    考慮到實際應用,選取單松弛時間的D3Q19(3是空間維數(shù),19 是離散速度數(shù))三維模型,如圖3 所示。 并根據(jù)三維模型各點坐標表征ci,i=0,1,2,…,18,指可能遷移方向數(shù)。

    圖3 D3Q19 模型Fig.3 D3Q19 model

    若知道粒子分布函數(shù),可計算模型的宏觀密度、速度,其表達式為

    2.2 離散單元法(DEM)

    離散單元法被廣泛應用于巖石力學和顆粒流等離散系統(tǒng)模擬中,是基于非連續(xù)介質(zhì)力學的數(shù)值方法[27]。 計算各個時間步上顆粒間所受的接觸力及顆粒運動,分析時將運動分為平動和轉(zhuǎn)動,平動在全局坐標系下計算,由牛頓方程給出,即

    轉(zhuǎn)動在局部坐標系下計算,需要計算出繞顆粒質(zhì)心的合力矩,顆粒轉(zhuǎn)動的控制方程為

    式中:Ip為主慣性矩;ω為角速度;rk為顆粒上第k個接觸處的矢徑;Md為阻尼力矩。

    2.3 浸入移動邊界法(IMB)

    浸入移動邊界法是LBM 和DEM 的耦合方法,其具備碰撞算子的局部性和遷移操作的簡單性,并在此基礎(chǔ)上克服了傳統(tǒng)耦合方法中的動量不連續(xù)性問題,具有足夠代表性的網(wǎng)格非一致邊界[28-30]。

    所有離散的格子被分為純流體格子、邊界格子和純固體格子3 類。 引入格子固含率εs定義每個格子內(nèi)流固相互作用的影響,并通過格子內(nèi)顆粒占據(jù)面積與格子占據(jù)總面積之比定量表征格子內(nèi)流場和顆粒分布情況。 流體格子的固含率εs=0,固體格子固含率εs=1,邊界格子的固含率εs為0~1,計算時用砂粒所占面積除以格子總面積,如圖4 所示。

    圖4 格子分類示意Fig.4 Schematic of grid classification

    為表征顆粒對流體的影響,在忽略外力項時將Boltzmann - BGK 方程添加源項變?yōu)?/p>

    β(εs,τ) 是與格子固含率有關(guān)的附加碰撞項的加權(quán)函數(shù)。 表示為格子固含率εs和松弛因子τ的函數(shù)為

    為考慮流體對顆粒的影響,通過加入所有被當前顆粒覆蓋的格子上對顆粒的力和力矩表示周圍流體作用于該顆粒的力和力矩,該力即顆粒受到的曳力,曳力F和力矩T的表達式分別為

    其中,n為該顆粒覆蓋邊界格子的數(shù)目;m為離散速度數(shù);xj為第j個邊界格子的位置;xc為顆粒中心位置。 根據(jù)受力和力矩,顆粒對位置進行更新,其顆粒平動速度和轉(zhuǎn)動角速度的更新式為

    式中:ρL為流體密度;ρp為顆粒密度;mL為流體質(zhì)量;g為重力加速度;Δt為計算時間步長。 最終通過顆粒和流體作用方程的不斷迭代計算,實現(xiàn)水砂混合流的液固耦合計算。

    2.4 計算流程

    基于LBM-DEM 實現(xiàn)了水砂混合流的液固耦合,大量學者已通過對比單顆粒沉降及兩顆粒DKT(追擊、接觸、翻轉(zhuǎn))試驗過程驗證LBM-DEM 耦合方法的合理性,可用于進行水砂混合流的液固耦合計算[31]。 其計算流程為在初始t0時刻輸入流場和顆粒信息,判斷各格點上的格子固含率。 通過LBM計算下一時刻t1時流場信息,并根據(jù)DEM 循環(huán)計算顆粒上的碰撞,將其信息更新為t1時刻。 依次對時間步循環(huán),直至判定精度達到要求結(jié)束計算,實現(xiàn)水砂混合流的液固耦合計算過程如圖5 所示。

    2.5 模型構(gòu)建及模擬方案

    考慮到實際突水潰砂問題,構(gòu)建單裂隙通道三維數(shù)值計算模型如圖6 所示。

    圖5 LBM-DEM 耦合計算過程Fig.5 Calculation process of LBM-DEM

    圖6 單裂隙通道三維數(shù)值計算模型Fig.6 Three-dimensional numerical model of single fracture passage

    以裂隙開口檔板最左側(cè)為O 點、建立YOZ 坐標系,三維空間尺寸為0.015 m×0.400 m×0.600 m(X×Y×Z),水流區(qū)域和初始砂粒區(qū)域均位于裂隙開口上方,上下兩面作用恒定壓力驅(qū)使水流流動,砂粒計算前先進行自由堆積,再通過坐標導入計算模型,使其均勻分布在開口上方。 模型參數(shù)設(shè)置:ρw=1 000 kg/m3;υ=0.001 Pa·s;ρs=2 500 kg/m3;r=0.002 m;E=4.9 GPa。

    為探究突水潰砂成災的主導因素,改變不同邊界壓力P、裂隙開口寬度d和砂層厚度h分析水砂混合流在模型中的流動情況,統(tǒng)計裂隙開口處截面流量及潰砂速率的變化,方案見表1。

    表1 模擬方案Table 1 Simulation schemes

    3 計算結(jié)果分析

    3.1 不同邊界壓力下水砂運移規(guī)律

    模擬不同邊界壓力下水砂混合流的運移情況如圖7 所示,隨著邊界壓力的增加,潰砂量和裂隙開口處最大流速明顯增加。 為定量研究潰砂增加量及流場變化,對裂隙開口處截面的流量Sw和單位時間潰砂速率vs與時間t關(guān)系進行分析(圖8、圖9)。

    隨著邊界壓力的增加,最終的穩(wěn)定流量逐漸增加(圖8)。 邊界壓力較小時流量先增加后保持恒定,而高壓力時在達到早期穩(wěn)定區(qū)持續(xù)一段時間后流量會再次升高,以P-5 最為明顯。 對比圖9 中潰砂速率變化可發(fā)現(xiàn)同樣的規(guī)律,壓力較小時潰砂速率先增至峰值后保持恒定,而隨著壓力的增加,潰砂速率從初始到峰值的增幅變大,但在峰值平臺區(qū)持續(xù)的時間變短,之后出現(xiàn)明顯降低,且隨著壓力增加降幅逐漸增大。 考慮到流量二次增加的時間點和潰砂速率從峰值開始降低的時間基本一致。 以P-5為例,分析了5 個等距時間步上截面流量、潰砂速率及潰砂形貌(圖10),并通過紅色虛線標記不同壓力時顆粒中的主要受力區(qū)。

    圖7 不同邊界壓力下水砂運移Fig.7 Migration of sediment at different boundary pressures

    圖8 不同邊界壓力下截面流量變化Fig.8 Flow changes of cross section under different boundary pressures

    圖9 不同邊界壓力下潰砂速率變化Fig.9 Changes of sand breaking rate underdifferent boundary pressures

    從圖10 所知,在①—③步中截面流量和潰砂速率基本不發(fā)生變化,此時在水流驅(qū)動作用下砂粒將在裂隙出口處形成臨時的拱狀受力結(jié)構(gòu),顆粒間相互作用使砂粒更加密集,對流場有一定的阻礙作用。隨著砂粒的不斷流出,出口處的拱狀結(jié)構(gòu)逐漸變化,其中心處不斷凹陷。 至第④步時拱狀結(jié)構(gòu)變?yōu)轭愃啤癕”狀結(jié)構(gòu),仍對水流及流砂有阻礙作用。 直至第⑤步時,中心砂粒完全流出喉道口,拱狀結(jié)構(gòu)完全消失后不再對流場起到阻礙作用,此時流量會快速上升,而由于砂粒已大量流出其潰砂速率會快速下降。

    圖10 P-5 截面流量和潰砂速率對比Fig.10 Comparison of cross section flow rate and sediment burst rate of Orifice P-5

    為表征不同邊界壓力對水砂流動的影響程度,統(tǒng)計單位時間潰砂速率最大值及初始階段受砂粒阻礙作用的平均流量隨壓力的變化規(guī)律如圖11 所示。其中隨著邊界壓力的增加,最大潰砂速率和平均流量均單調(diào)遞增。 計算最大潰砂速率及初始平均孔口流量隨壓力的區(qū)間增幅,并通過區(qū)間增幅除以區(qū)間范圍得出單位壓力增幅以表征不同壓力數(shù)值下的影響情況,計算結(jié)果見表2。

    圖11 截面流量和潰砂速率隨邊界壓力變化Fig.11 Variations of cross section flow rate and sediment burst rate with boundary pressures

    表2 截面流量和潰砂速率隨邊界壓力的增幅Table 2 Increase of cross section flow rate and sediment burst rate with boundary pressure

    從表3 發(fā)現(xiàn),低壓力時隨壓力的增加最大潰砂速率和平均流量的增幅均大于高壓力時。 說明僅考慮邊界壓力情況下,隨著邊界壓力的增加,同等時間下突水潰砂災害程度會逐漸加劇,但低壓力下由于壓力增加引發(fā)的災害加劇程度較大,此時應注意由壓力增加引起的劇烈波動。 由于模型中設(shè)置恒定壓差表征地下承壓含水層作用,未考慮實際流動中壓力降低的作用,結(jié)論僅能在一定程度上適用于工程實際問題,后續(xù)將針對復雜模型繼續(xù)研究。

    針對干砂單裂隙開口通道下流出孔口時單位質(zhì)量流量計算,1961 年BEVERLOO 等[32]最早給出單位質(zhì)量流量W的計算式為

    式中:C為與顆粒間摩擦性質(zhì)相關(guān)的無因次量;ρB為顆粒堆積密度,通過顆粒密度ρs與空隙度ε計算,表示為ρB=ρs1-ε( ) ;D0為孔口直徑;dp為顆粒直徑;k為無因次系數(shù)。

    但式(1) 未考慮流體與顆粒并存情況,BULSARA 等[33]引入壓力梯度概念得出含流體時顆粒流出孔口時單位質(zhì)量流量計算式為

    其中:P1和P2分別為孔板下方及上方壓力。 郭帥等[34]通過試驗分析驗證了其合理性,模擬時設(shè)置不同恒定壓差以驅(qū)動水流,參照式(2)給出最大潰砂量隨壓差的計算式為

    其中:ΔP為模型恒定壓差;W0為干砂在孔口的單位質(zhì)量流量,能有效表征壓差為零時干砂受自重下落的情況。 圖12 給出式(3)擬合曲線與數(shù)值模擬數(shù)據(jù)點對比情況,曲線擬合情況較好。

    圖12 單位質(zhì)量流量數(shù)據(jù)點與擬合曲線對比Fig.12 Comparison of data points per unit mass flow and fitting curve

    3.2 不同開口寬度下水砂運移規(guī)律

    模擬不同開口寬度下的水砂混合流運移情況表明:隨開口寬度的增加,潰砂量及開口處最大流速明顯增加,且流速高速區(qū)逐漸貫穿整個流場。 為明確潰砂增加量及流場變化情況,圖13、圖14 中結(jié)合裂隙開口處相同截面的流量Sw和單位時間潰砂速率vs進行分析。

    圖13 不同開口寬度下截面流量變化Fig.13 Variation of flow rate of orifice section with different opening widths

    隨著開口寬度的增加,最終的穩(wěn)定流量逐漸增加,如圖13 所示。 開口較小時出口截面流量在達到初始流量峰值后保持恒定一段時間后再逐漸增加,開口較大時,達到初始流量峰值后將出現(xiàn)明顯的下降階段,之后又由于顆粒大量流出開口使流量逐漸增至最終穩(wěn)定流量,以K-4、K-5 最為明顯。

    隨著開口寬度的增加,潰砂速率從初始到峰值的增幅變大,且在峰值平臺區(qū)持續(xù)的時間變短,如圖14 所示。 開口寬度較小時潰砂速率增加到峰值將保持恒定一段時間,后由于砂粒逐漸流出裂隙出現(xiàn)緩慢降低。 開口寬度較大時,潰砂速率達到峰值保持一定時間后由于砂粒流出開口急速降低。 在開口進一步擴大后,甚至在峰值持續(xù)階段出現(xiàn)一段小幅降低。

    圖14 不同開口寬度下潰砂速率變化Fig.14 Change of sediment discharge rate under different opening widths

    圖15 K-4 截面流量和潰砂速率對比Fig.15 Comparison of cross section flow rate and sediment burst rate of K-4 orifice

    因此,以K-4 為例,結(jié)合流量、潰砂速率、小幅降低時刻的水砂分布情況進行分析,并用紅色虛線標示出大量砂粒受擠壓聚集在開口處的現(xiàn)象,如圖15 所示。 分析可知砂粒受擠壓在開口處形成密實堆積結(jié)構(gòu)影響水砂運移,造成潰砂速率和截面流量小幅下降,但由于結(jié)構(gòu)中的砂粒具有流動性,形成的密實結(jié)構(gòu)會隨著砂粒流出開口逐漸失去其阻礙效果,直至砂粒中心平面降至孔口下方時該結(jié)構(gòu)消失。

    圖16 中選取潰砂速率峰值時不同孔口下水砂分布形態(tài)進行分析,如圖16 所示,以紅色虛線畫出不同開口上方的密實結(jié)構(gòu)區(qū)域,隨著開口寬度的增加,在孔口處形成的密實結(jié)構(gòu)體積不斷增加,對流場造成的影響程度將不斷增大。 開口較小時,由于密實結(jié)構(gòu)范圍小,其受擠壓砂粒相對總體砂粒占比較小,影響效果不明顯。 開口較大時密實區(qū)范圍較大,影響砂粒數(shù)量占比較多,將對整體水砂運移產(chǎn)生巨大影響。

    圖16 不同開口寬度下水砂運移Fig.16 Migration diagram of sediment with different opening widths

    為表征不同開口寬度對水砂流動的影響程度,圖17 中給出單位時間潰砂速率最大值及初始階段受砂粒阻礙作用的平均流量隨開口寬度的變化規(guī)律。

    圖17 截面流量和潰砂速率隨開口寬度的變化Fig.17 Variation of cross section flow rate and sand breaking rate of orifice vary with opening width

    隨著邊界壓力的增加,最大潰砂速率和初始階段平均流量均單調(diào)遞增(圖17)。 為量化開口寬度的影響,表3 給出最大潰砂速率和初始平均流量隨開口寬度的區(qū)間增幅。

    通過表3 發(fā)現(xiàn),開口寬度較小時隨寬度增加對最大潰砂速率和初始平均流量的影響大于開口寬度較大時。 說明固定的砂含量下,隨著裂隙開口寬度的增加,突水潰砂災害程度的增幅逐漸減小,在寬度大于一定值后,由于砂含量的限制災害逐漸趨于恒定。 突增現(xiàn)象最容易出現(xiàn)在細小裂隙開口擴張的過程中。

    表3 截面流量和潰砂速率隨開口寬度的增幅Table 3 Variation of cross section flow rate and sand breaking rate increase with width of opening

    3.3 不同砂層厚度下水砂運移規(guī)律

    密實填充下砂層厚度與砂粒數(shù)目呈正比,不同砂層厚度可表征不同砂粒數(shù)目對水砂運移的影響。模擬不同砂層厚度下水砂運移結(jié)果表明:隨著砂層厚度的增加,砂粒數(shù)目增加,最終的潰砂量明顯增加,且開口兩側(cè)堆積的未流出砂粒也逐漸增加,但流速未產(chǎn)生明顯變化。 為明確潰砂增加量及流場變化情況,圖18、圖19 結(jié)合開口處的截面流量Sw和單位時間潰砂速率vs進行分析。

    圖18 不同砂層厚度下截面流量變化Fig.18 Variation of cross section flow of orifice with sand thickness

    圖19 不同砂層厚度下潰砂速率變化Fig.19 Changes of sediment discharge rate under different sand thickness

    從圖18 所知,最終穩(wěn)定流量不隨砂層厚度的增加發(fā)生變化。 而隨著砂層厚度的增加,截面流量在達到初始峰值后均出現(xiàn)一定的下降,用對應曲線顏色的虛線記錄最終時間,①表示初始流量峰值點,②③④⑤分別對應N-2、N-3、N-4、N-5 時最終下落結(jié)束的時間點,發(fā)現(xiàn)下降持續(xù)的時間隨著砂層厚度的增加逐漸增長。 說明隨著砂層厚度的增加,下落過程中對水流的阻礙作用越明顯。

    如圖19 所示,隨著砂層厚度的增加,潰砂速率峰值點位置不斷延后。 在砂層厚度較大時,均呈現(xiàn)先急速增至初始高速率后持續(xù)降低一段時間,再升至峰值速率后快速下降的趨勢。 同樣用①表示初始速率峰值點,②③④⑤分別對應N-2、N-3、N-4、N-5 時速率持續(xù)下降結(jié)束的時間點,發(fā)現(xiàn)隨著砂層厚度的增加,從初始高速率降至最小速率的持續(xù)時間不斷增長,說明砂層厚度的增加對砂層整體潰突阻礙作用逐漸增加。 為分析砂層厚度較大時對水砂運移的阻礙形式,圖20 以N-5 為例,結(jié)合流量、潰砂速率及不同時刻水砂分布情況進行分析。

    圖20 N-5 截面流量和潰砂速率對比Fig.20 Comparison of cross section flow rate and sediment burst rate of Orifice N-5

    從圖20 可知,對水流及整體砂粒流速起阻礙作用的原因主要為砂粒在裂隙開口處聚集形成密實區(qū)域,對比各階段可以發(fā)現(xiàn)①②③④時刻大量砂粒受力后聚集在裂隙開口,使得水砂流動減緩,在⑤時刻砂粒形成的密實結(jié)構(gòu)崩塌,對流體的阻礙消失,截面流量逐漸開始上升。 而此時由于密實結(jié)構(gòu)崩塌后大量砂粒傾瀉將造成此時的潰砂速率達到峰值。 為明確不同砂層厚度對水砂運移的影響,圖21 分析不同砂層高度下單位時間潰砂速率最大值及初始階段受砂粒阻礙作用階段平均流量的變化規(guī)律。

    圖21 截面流量和潰砂速率隨砂層厚度的變化Fig.21 Variation of cross section flow rate of orifice and rate of sediment burst change with thickness of sand layer

    隨著砂層厚度的增加,最大潰砂速率和平均流量均單調(diào)遞減(圖21),在砂層厚度較小時出現(xiàn)均快速降低的現(xiàn)象。 當砂層厚度較大時,隨砂層厚度增加對最大潰砂速率和平均流量的影響較小。 為定量研究砂層厚度對水砂運移的影響程度,表4 計算給出不同砂層厚度下最大潰砂速率和平均流量的區(qū)間降幅。

    表4 截面流量和潰砂速率隨砂層高度的降幅Table 4 Decrease of cross section flow rate and sediment burst rate with thickness of sand layer

    由表5 知,在砂層厚度較小時隨砂層厚度增加最大潰砂速率和平均流量的降幅遠大于砂層厚度較大時,在N-2 時的降幅基本為N-5 時的10 倍。 說明在邊界壓力、砂層范圍及裂隙開口恒定的情況下,隨著砂層厚度的增加,包含的砂粒數(shù)量增多,同等時間下由于砂粒阻礙作用增強其突水潰砂災害程度會逐漸減小。 當砂層厚度較小時,水驅(qū)動砂快速流出裂隙孔口,災害最為嚴重。 而砂層越厚砂粒越容易在裂隙口擁堵阻礙砂粒及流體運動,此時在砂粒密實結(jié)構(gòu)失穩(wěn)造成砂粒傾瀉時發(fā)生最嚴重的災害。 分析僅考慮同等時間下的災害影響,未研究不同砂層厚度砂粒全部流出情況,此情況明顯為砂層厚度越大,可流出砂粒增加造成災害愈加嚴重。

    4 結(jié)論及展望

    1)隨邊界壓力增加,截面流量及最大潰砂速率均單調(diào)增加,同等時間下突水潰砂災害程度會逐漸加劇。 由于高壓力下砂粒更加容易在喉道口堆積成臨時性密實結(jié)構(gòu)阻礙水砂運移,升高單位壓力時引起的最大潰砂量及平均截面流量的增幅逐漸降低。在Bulsara 式基礎(chǔ)上考慮流體與顆粒并存情況下給出的計算式能較好地擬合最大單位質(zhì)量流量曲線。

    2)隨開口寬度增加,截面流量及最大潰砂速率均單調(diào)增加,同等時間下突水潰砂災害程度會逐漸加劇。 砂粒在開口處形成的密實結(jié)構(gòu)體積逐漸增加單增速不斷降低,截面流量及最大潰砂速率增幅逐漸降低,小開口時由寬度增加引發(fā)的災害加劇大于大開口情況。 發(fā)現(xiàn)砂粒在開口處形成的密實結(jié)構(gòu)體積增加速率對災害程度起關(guān)鍵影響,細小裂隙開口擴張的過程中更易出現(xiàn)水砂流速突增,誘發(fā)突水潰砂災害。

    3)隨砂層厚度增加,截面流量及最大潰砂速率均單調(diào)減少,且降幅在小厚度時達最大后逐級減小至趨于平穩(wěn),小厚度區(qū)間降幅可達大厚度的10 倍,在開口處形成的砂粒密實結(jié)構(gòu)持續(xù)時間增長,對水砂運移的阻礙效果逐漸增加,同等時間下突水潰砂災害程度會逐漸減弱。 不同厚度下災害特點不同,厚度較小時,災害發(fā)生在水驅(qū)動砂快速流出裂隙孔口階段,厚度較大時,災害發(fā)生在砂粒密實結(jié)構(gòu)失穩(wěn)造成砂粒傾瀉階段。

    實際突水潰砂問題不僅需要考慮裂隙開口處的水砂流動問題,還需要分析在裂隙網(wǎng)絡(luò)中的流動及對裂隙網(wǎng)絡(luò)的侵蝕等相關(guān)問題。 后續(xù)將針對平面裂隙網(wǎng)絡(luò)中的流動問題、空間破裂巖體通道中的流動問題及裂隙通道中的滲流侵蝕問題進行研究,為突水潰砂災害防控提供理論指導。

    猜你喜歡
    砂粒突水砂層
    礦井突水水源的判別方法
    第四系膠結(jié)砂層水文地質(zhì)鉆進技術(shù)研究初探
    礦井滯后突水機理探討
    下降管蓄熱器中沙漠砂流動性數(shù)值分析
    主動出擊
    新少年(2020年10期)2020-10-30 02:04:05
    用于粒子分離器的砂粒反彈特性實驗研究
    上海航天(2018年3期)2018-06-25 02:57:48
    李雅莊礦滯后突水主控因素分析
    西安主城區(qū)砂層工程地質(zhì)特征分析
    大相嶺隧道高壓突水機理與預測分析
    汽車制動檢驗臺滾筒粘砂層脫落的原因及維護方法
    久久精品亚洲熟妇少妇任你| 一个人免费在线观看的高清视频| 久久九九热精品免费| 成人国语在线视频| 成人精品一区二区免费| 伊人久久大香线蕉亚洲五| av视频免费观看在线观看| 香蕉丝袜av| 久久香蕉激情| 国产在线精品亚洲第一网站| 黄色成人免费大全| 在线永久观看黄色视频| 最近最新免费中文字幕在线| www.精华液| 又黄又粗又硬又大视频| 亚洲精品在线观看二区| 大香蕉久久网| 亚洲人成电影免费在线| 黄片播放在线免费| 亚洲精品国产区一区二| 国产一区二区三区综合在线观看| 欧美人与性动交α欧美精品济南到| 国产亚洲午夜精品一区二区久久| 一边摸一边抽搐一进一出视频| 亚洲精品在线观看二区| 狠狠婷婷综合久久久久久88av| 欧美日韩福利视频一区二区| 美女福利国产在线| 久久精品国产a三级三级三级| 欧美黑人欧美精品刺激| 亚洲色图 男人天堂 中文字幕| 自拍欧美九色日韩亚洲蝌蚪91| 成年女人毛片免费观看观看9 | 国产精品美女特级片免费视频播放器 | 亚洲国产av影院在线观看| 黄色怎么调成土黄色| 老熟妇乱子伦视频在线观看| 美女午夜性视频免费| 国产精品一区二区在线观看99| 国产欧美日韩一区二区精品| xxxhd国产人妻xxx| 亚洲国产成人一精品久久久| 欧美日韩中文字幕国产精品一区二区三区 | 操出白浆在线播放| 丝袜美足系列| 国产男靠女视频免费网站| 99精品欧美一区二区三区四区| 电影成人av| 自拍欧美九色日韩亚洲蝌蚪91| 国产精品久久久久成人av| 亚洲精品在线观看二区| 人人妻人人澡人人爽人人夜夜| 麻豆国产av国片精品| av天堂久久9| 国产精品久久久久久精品电影小说| 午夜福利视频在线观看免费| 国产成人精品久久二区二区91| 99热网站在线观看| 在线观看免费高清a一片| 少妇粗大呻吟视频| 黄色视频不卡| 精品一区二区三区四区五区乱码| 国产精品美女特级片免费视频播放器 | 下体分泌物呈黄色| 午夜福利乱码中文字幕| 性色av乱码一区二区三区2| 两人在一起打扑克的视频| 中文字幕av电影在线播放| 超碰97精品在线观看| 精品乱码久久久久久99久播| 中文字幕人妻熟女乱码| 免费在线观看日本一区| 老司机靠b影院| 高清视频免费观看一区二区| 1024视频免费在线观看| 精品午夜福利视频在线观看一区 | 中文字幕人妻丝袜一区二区| 成人手机av| 两个人看的免费小视频| 久久久久久久大尺度免费视频| 久久久久精品人妻al黑| 欧美日韩福利视频一区二区| 国产成人精品久久二区二区免费| 午夜福利视频在线观看免费| 母亲3免费完整高清在线观看| 欧美另类亚洲清纯唯美| 亚洲三区欧美一区| 国产一区二区三区综合在线观看| 性少妇av在线| 成年人午夜在线观看视频| 国产日韩欧美在线精品| 免费观看av网站的网址| 国产片内射在线| 亚洲专区国产一区二区| 天堂动漫精品| 无人区码免费观看不卡 | 亚洲av日韩在线播放| 成人国语在线视频| 久久国产精品影院| 老汉色av国产亚洲站长工具| 国产色视频综合| 黄色怎么调成土黄色| 在线亚洲精品国产二区图片欧美| 99在线人妻在线中文字幕 | 久久午夜综合久久蜜桃| 国产亚洲一区二区精品| 视频区欧美日本亚洲| 一区二区av电影网| 五月天丁香电影| 激情视频va一区二区三区| 波多野结衣av一区二区av| 美女高潮喷水抽搐中文字幕| 国产欧美亚洲国产| 亚洲七黄色美女视频| 欧美亚洲 丝袜 人妻 在线| 热99re8久久精品国产| 法律面前人人平等表现在哪些方面| 青青草视频在线视频观看| 又黄又粗又硬又大视频| 每晚都被弄得嗷嗷叫到高潮| 久久久国产欧美日韩av| 亚洲av美国av| 国产av精品麻豆| 69精品国产乱码久久久| 国产成人精品久久二区二区免费| 免费黄频网站在线观看国产| 十八禁人妻一区二区| 9色porny在线观看| 国产亚洲精品第一综合不卡| 一级,二级,三级黄色视频| av天堂久久9| 黄色 视频免费看| 国产精品九九99| 国产精品美女特级片免费视频播放器 | 欧美久久黑人一区二区| 在线观看免费高清a一片| 一区二区三区激情视频| 男女床上黄色一级片免费看| 视频区图区小说| 国产在视频线精品| 久久精品国产亚洲av香蕉五月 | 中文字幕最新亚洲高清| 国产欧美日韩一区二区三区在线| 亚洲av第一区精品v没综合| 美女扒开内裤让男人捅视频| 国产日韩欧美视频二区| 成人免费观看视频高清| 一本—道久久a久久精品蜜桃钙片| 欧美精品啪啪一区二区三区| 国产欧美日韩一区二区三| 国产精品成人在线| 在线十欧美十亚洲十日本专区| 久久久欧美国产精品| 最黄视频免费看| 午夜91福利影院| 国产又色又爽无遮挡免费看| 男人操女人黄网站| 在线十欧美十亚洲十日本专区| a级片在线免费高清观看视频| 别揉我奶头~嗯~啊~动态视频| 国产亚洲精品一区二区www | 中国美女看黄片| 脱女人内裤的视频| av在线播放免费不卡| 国产黄频视频在线观看| 无限看片的www在线观看| 久久精品aⅴ一区二区三区四区| 超色免费av| 啦啦啦视频在线资源免费观看| 无遮挡黄片免费观看| 中亚洲国语对白在线视频| 色视频在线一区二区三区| 亚洲少妇的诱惑av| 99国产精品一区二区蜜桃av | 亚洲色图 男人天堂 中文字幕| 久久久久视频综合| 国精品久久久久久国模美| av不卡在线播放| 老熟女久久久| 国产精品二区激情视频| 啪啪无遮挡十八禁网站| 国产区一区二久久| 窝窝影院91人妻| 一本—道久久a久久精品蜜桃钙片| 久久久久视频综合| 不卡av一区二区三区| kizo精华| 成人手机av| 日本av手机在线免费观看| 亚洲三区欧美一区| 亚洲人成77777在线视频| 老司机午夜十八禁免费视频| 欧美精品啪啪一区二区三区| 免费在线观看完整版高清| 国产午夜精品久久久久久| 久久久国产成人免费| 日本av手机在线免费观看| 黄片小视频在线播放| 女同久久另类99精品国产91| 国产日韩欧美视频二区| 久久人人97超碰香蕉20202| 在线观看一区二区三区激情| 欧美黑人精品巨大| 大香蕉久久网| 丝袜美腿诱惑在线| 精品国产一区二区久久| netflix在线观看网站| 亚洲少妇的诱惑av| 久久99一区二区三区| 成人影院久久| 亚洲成人免费电影在线观看| 黄片小视频在线播放| 久久久久久亚洲精品国产蜜桃av| 不卡av一区二区三区| 国产xxxxx性猛交| 制服诱惑二区| 国产亚洲精品一区二区www | 美女国产高潮福利片在线看| 大香蕉久久网| 性少妇av在线| 侵犯人妻中文字幕一二三四区| 欧美黄色片欧美黄色片| 黄色片一级片一级黄色片| a级毛片黄视频| 18禁国产床啪视频网站| 国产成人一区二区三区免费视频网站| 精品国产国语对白av| 又紧又爽又黄一区二区| 国产深夜福利视频在线观看| 国产一卡二卡三卡精品| videosex国产| 男女之事视频高清在线观看| 人人妻人人澡人人看| 国产精品欧美亚洲77777| 色播在线永久视频| 午夜激情av网站| 国产精品影院久久| 欧美变态另类bdsm刘玥| 母亲3免费完整高清在线观看| 老司机午夜十八禁免费视频| 丰满人妻熟妇乱又伦精品不卡| 免费在线观看视频国产中文字幕亚洲| 日韩欧美一区视频在线观看| 亚洲人成伊人成综合网2020| 精品少妇一区二区三区视频日本电影| 国产av精品麻豆| 自拍欧美九色日韩亚洲蝌蚪91| 久热爱精品视频在线9| 亚洲,欧美精品.| 别揉我奶头~嗯~啊~动态视频| 午夜激情av网站| 99在线人妻在线中文字幕 | 首页视频小说图片口味搜索| 亚洲中文字幕日韩| 国产精品久久久久久人妻精品电影 | aaaaa片日本免费| 国产精品熟女久久久久浪| 精品一品国产午夜福利视频| 99精品久久久久人妻精品| 一二三四在线观看免费中文在| 中文字幕人妻丝袜一区二区| 欧美黄色淫秽网站| 日韩一区二区三区影片| 国产三级黄色录像| av天堂久久9| 伊人久久大香线蕉亚洲五| 男女免费视频国产| 最新美女视频免费是黄的| 亚洲成人免费av在线播放| 精品福利永久在线观看| 一级毛片女人18水好多| 国产一区有黄有色的免费视频| av福利片在线| 久久精品亚洲精品国产色婷小说| 国产一卡二卡三卡精品| 桃花免费在线播放| 色婷婷久久久亚洲欧美| 久久中文看片网| 午夜精品久久久久久毛片777| 91九色精品人成在线观看| 91字幕亚洲| 国产精品成人在线| 香蕉久久夜色| 天堂动漫精品| 高清av免费在线| 亚洲伊人色综图| 欧美精品一区二区大全| 免费一级毛片在线播放高清视频 | 欧美激情 高清一区二区三区| 久久久国产成人免费| 亚洲五月色婷婷综合| 岛国毛片在线播放| www.熟女人妻精品国产| 亚洲成人免费av在线播放| 中文字幕人妻熟女乱码| av免费在线观看网站| 国产伦理片在线播放av一区| 免费不卡黄色视频| 亚洲国产欧美日韩在线播放| 中文字幕人妻丝袜制服| 我的亚洲天堂| 久久久久久久大尺度免费视频| 免费高清在线观看日韩| 9191精品国产免费久久| 欧美午夜高清在线| 国产一区二区在线观看av| 热re99久久精品国产66热6| 黄色视频,在线免费观看| 搡老乐熟女国产| 女人被躁到高潮嗷嗷叫费观| 国产精品成人在线| 久久精品亚洲熟妇少妇任你| av有码第一页| 欧美成人免费av一区二区三区 | e午夜精品久久久久久久| 手机成人av网站| 欧美日韩亚洲综合一区二区三区_| 成人影院久久| 在线观看人妻少妇| 久9热在线精品视频| 一级黄色大片毛片| 少妇裸体淫交视频免费看高清 | 91成人精品电影| 久久亚洲真实| 丰满少妇做爰视频| 久久国产精品大桥未久av| 法律面前人人平等表现在哪些方面| 欧美日韩亚洲高清精品| 视频在线观看一区二区三区| 久热爱精品视频在线9| 露出奶头的视频| 国产亚洲午夜精品一区二区久久| 丁香六月天网| 狂野欧美激情性xxxx| 成人特级黄色片久久久久久久 | 国产成人av教育| 久久精品亚洲精品国产色婷小说| 亚洲人成77777在线视频| 黑人猛操日本美女一级片| 午夜精品国产一区二区电影| a在线观看视频网站| 欧美日韩av久久| 久久久国产精品麻豆| 在线播放国产精品三级| 日韩制服丝袜自拍偷拍| 亚洲av国产av综合av卡| 搡老熟女国产l中国老女人| 另类亚洲欧美激情| 久久精品国产99精品国产亚洲性色 | 欧美乱码精品一区二区三区| kizo精华| 高清av免费在线| 妹子高潮喷水视频| 欧美精品人与动牲交sv欧美| 久热爱精品视频在线9| 久久人妻福利社区极品人妻图片| 桃花免费在线播放| 午夜福利影视在线免费观看| 乱人伦中国视频| 午夜福利在线免费观看网站| 99热国产这里只有精品6| 成年版毛片免费区| 亚洲国产精品一区二区三区在线| 香蕉久久夜色| 国内毛片毛片毛片毛片毛片| 日日摸夜夜添夜夜添小说| 一二三四社区在线视频社区8| 国产精品久久电影中文字幕 | 美女视频免费永久观看网站| a级毛片在线看网站| 咕卡用的链子| 亚洲精品粉嫩美女一区| 亚洲av成人一区二区三| 捣出白浆h1v1| 搡老乐熟女国产| 18禁国产床啪视频网站| 亚洲精品乱久久久久久| 蜜桃国产av成人99| 纵有疾风起免费观看全集完整版| 亚洲av片天天在线观看| 日韩视频一区二区在线观看| 男男h啪啪无遮挡| 亚洲成人国产一区在线观看| 亚洲情色 制服丝袜| h视频一区二区三区| 精品久久久久久电影网| 亚洲国产欧美日韩在线播放| 一区二区av电影网| 久久久久久久大尺度免费视频| 亚洲精品自拍成人| 国产1区2区3区精品| 大香蕉久久网| 丰满迷人的少妇在线观看| 搡老岳熟女国产| 国产麻豆69| 一边摸一边做爽爽视频免费| 欧美亚洲 丝袜 人妻 在线| 国产老妇伦熟女老妇高清| 国产男女超爽视频在线观看| 久久精品亚洲熟妇少妇任你| 俄罗斯特黄特色一大片| 日日摸夜夜添夜夜添小说| 久久天堂一区二区三区四区| 成人18禁高潮啪啪吃奶动态图| 久久国产精品男人的天堂亚洲| 黄色片一级片一级黄色片| 亚洲中文字幕日韩| 波多野结衣一区麻豆| 一级a爱视频在线免费观看| 中文字幕人妻丝袜制服| 高潮久久久久久久久久久不卡| 啦啦啦 在线观看视频| 欧美黄色淫秽网站| 国产精品亚洲一级av第二区| 大型av网站在线播放| 一区二区三区乱码不卡18| 久久人人97超碰香蕉20202| 黄色视频不卡| 99热网站在线观看| 无限看片的www在线观看| 啦啦啦免费观看视频1| 国产欧美日韩一区二区精品| 亚洲中文字幕日韩| 国产黄色免费在线视频| 超碰成人久久| 啪啪无遮挡十八禁网站| 欧美 亚洲 国产 日韩一| 国产三级黄色录像| 亚洲av片天天在线观看| 久久精品亚洲av国产电影网| 97人妻天天添夜夜摸| 精品国产乱码久久久久久男人| 真人做人爱边吃奶动态| 精品视频人人做人人爽| 国产人伦9x9x在线观看| 两个人免费观看高清视频| 成年女人毛片免费观看观看9 | 亚洲人成电影免费在线| 十八禁高潮呻吟视频| 国产欧美日韩一区二区精品| 热re99久久精品国产66热6| 免费女性裸体啪啪无遮挡网站| 在线 av 中文字幕| 狠狠狠狠99中文字幕| 亚洲午夜理论影院| 欧美精品高潮呻吟av久久| 亚洲国产欧美一区二区综合| 91麻豆av在线| 99国产精品一区二区三区| 一区二区日韩欧美中文字幕| tube8黄色片| 亚洲欧美精品综合一区二区三区| 热99久久久久精品小说推荐| av片东京热男人的天堂| 亚洲全国av大片| a在线观看视频网站| 曰老女人黄片| 免费女性裸体啪啪无遮挡网站| 欧美精品亚洲一区二区| cao死你这个sao货| 色婷婷久久久亚洲欧美| 久久性视频一级片| 女同久久另类99精品国产91| 久久久久国内视频| 悠悠久久av| 久久久国产成人免费| 久久国产精品影院| 亚洲av国产av综合av卡| 日韩免费av在线播放| 久久午夜亚洲精品久久| 日韩一区二区三区影片| 一级,二级,三级黄色视频| 不卡av一区二区三区| 亚洲av日韩在线播放| 超碰成人久久| 老司机影院毛片| 两个人看的免费小视频| 91成年电影在线观看| 日本av手机在线免费观看| 精品乱码久久久久久99久播| 中文字幕av电影在线播放| av有码第一页| 国产精品电影一区二区三区 | 他把我摸到了高潮在线观看 | 涩涩av久久男人的天堂| 精品久久蜜臀av无| 欧美精品高潮呻吟av久久| 日韩一卡2卡3卡4卡2021年| 精品视频人人做人人爽| 国产精品久久电影中文字幕 | 国产野战对白在线观看| 高清黄色对白视频在线免费看| 久久影院123| av天堂在线播放| 91国产中文字幕| 999久久久精品免费观看国产| 久久精品成人免费网站| 无限看片的www在线观看| 久久久久久久精品吃奶| 亚洲一区中文字幕在线| 亚洲七黄色美女视频| 国产在视频线精品| 桃红色精品国产亚洲av| 热re99久久精品国产66热6| 两性夫妻黄色片| 亚洲少妇的诱惑av| 精品午夜福利视频在线观看一区 | 久久精品国产亚洲av高清一级| 麻豆成人av在线观看| 久久久久国内视频| 午夜两性在线视频| 国产色视频综合| 涩涩av久久男人的天堂| 精品久久蜜臀av无| 国产一区有黄有色的免费视频| 久久精品国产99精品国产亚洲性色 | 欧美国产精品va在线观看不卡| 久久国产精品大桥未久av| 欧美 亚洲 国产 日韩一| 少妇裸体淫交视频免费看高清 | 亚洲精品在线美女| 两人在一起打扑克的视频| 免费黄频网站在线观看国产| 后天国语完整版免费观看| 在线观看人妻少妇| 首页视频小说图片口味搜索| 99re在线观看精品视频| 久久国产精品人妻蜜桃| 国产在线免费精品| 在线观看舔阴道视频| 交换朋友夫妻互换小说| 在线观看66精品国产| 国产色视频综合| 欧美av亚洲av综合av国产av| 亚洲avbb在线观看| 久久毛片免费看一区二区三区| 热re99久久国产66热| 丝袜美腿诱惑在线| 一边摸一边抽搐一进一出视频| 国产欧美日韩一区二区精品| 18禁观看日本| 精品国产一区二区久久| 一区二区日韩欧美中文字幕| 我的亚洲天堂| 一个人免费在线观看的高清视频| 在线观看66精品国产| 国产精品成人在线| 在线观看免费视频日本深夜| 极品人妻少妇av视频| 久久久国产欧美日韩av| 美女扒开内裤让男人捅视频| avwww免费| av天堂久久9| 久久毛片免费看一区二区三区| 18禁黄网站禁片午夜丰满| 三上悠亚av全集在线观看| 精品久久久久久电影网| 日韩欧美国产一区二区入口| 91国产中文字幕| 老熟妇乱子伦视频在线观看| 99久久国产精品久久久| 波多野结衣av一区二区av| 视频区图区小说| 成人特级黄色片久久久久久久 | 国内毛片毛片毛片毛片毛片| 久久人妻av系列| 国产男靠女视频免费网站| 国产熟女午夜一区二区三区| 亚洲精品一卡2卡三卡4卡5卡| 久久精品国产亚洲av高清一级| 精品福利观看| 最近最新中文字幕大全免费视频| 国产精品99久久99久久久不卡| 成在线人永久免费视频| 日日爽夜夜爽网站| 18禁国产床啪视频网站| 色综合欧美亚洲国产小说| 极品少妇高潮喷水抽搐| 麻豆av在线久日| 午夜日韩欧美国产| 黄色视频不卡| 50天的宝宝边吃奶边哭怎么回事| 视频区欧美日本亚洲| 两性午夜刺激爽爽歪歪视频在线观看 | 老司机午夜十八禁免费视频| 国产欧美日韩一区二区三| 亚洲精品av麻豆狂野| 一级,二级,三级黄色视频| 老熟妇仑乱视频hdxx| 久久久久精品国产欧美久久久| 人人妻人人爽人人添夜夜欢视频| 99国产精品一区二区三区| av超薄肉色丝袜交足视频| 涩涩av久久男人的天堂| 色尼玛亚洲综合影院| 精品国产乱码久久久久久男人| 一级a爱视频在线免费观看| 天堂中文最新版在线下载| 日韩免费高清中文字幕av| videosex国产| www日本在线高清视频| 肉色欧美久久久久久久蜜桃| 国精品久久久久久国模美| 午夜福利在线免费观看网站| 国产成人免费无遮挡视频| 亚洲人成电影观看| 亚洲欧美精品综合一区二区三区| 肉色欧美久久久久久久蜜桃| 亚洲专区字幕在线| 纵有疾风起免费观看全集完整版| 岛国在线观看网站| 男女之事视频高清在线观看| 精品少妇一区二区三区视频日本电影| 91av网站免费观看| 精品少妇内射三级|