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

    潰壩洪水研究進(jìn)展

    2016-03-26 22:55:36常福宣肖長偉高延鴻程令社
    長江科學(xué)院院報 2016年6期
    關(guān)鍵詞:數(shù)值模擬

    劉 林,常福宣,肖長偉,高延鴻,程令社

    (1.長江科學(xué)院 水資源綜合利用研究所,武漢 430010;2.西藏自治區(qū)水利電力規(guī)劃勘測設(shè)計研究院,拉薩 850000 )

    潰壩洪水研究進(jìn)展

    劉林1,常福宣1,肖長偉2,高延鴻2,程令社2

    (1.長江科學(xué)院 水資源綜合利用研究所,武漢430010;2.西藏自治區(qū)水利電力規(guī)劃勘測設(shè)計研究院,拉薩850000 )

    摘要:潰壩洪水研究的目的是計算潰決壩址的流量、水位過程線,并向下游作洪水演進(jìn)得到沿程的流量、流速、水位、波前與洪峰的到達(dá)時間,評估下游洪水淹沒損失情況,以便于采取措施,降低洪水風(fēng)險。從潰壩水流理論研究、潰壩問題的試驗(yàn)研究、潰壩模擬及洪水在下游演進(jìn)這3個方面對潰壩洪水研究進(jìn)行了綜述,回顧和總結(jié)了國內(nèi)外潰壩洪水研究的發(fā)展歷程、已取得的成果和近些年的進(jìn)展,提出了將來要研究的重點(diǎn),并對研究前景進(jìn)行了展望。目前潰壩理論的數(shù)值求解發(fā)展迅速,由試驗(yàn)提出了潰壩機(jī)理,研究不斷模型化;但高強(qiáng)輸沙理論未建立,潰口沖刷過程未能準(zhǔn)確表達(dá),對梯級潰壩和冰湖潰決洪水研究較少。今后應(yīng)加強(qiáng)潰壩水流理論研究,開展大尺度、多庫潰壩模型試驗(yàn)研究,積極做好梯級潰壩和冰湖潰決洪水模擬,進(jìn)一步建立快捷可靠的區(qū)域潰壩洪水預(yù)報系統(tǒng)。

    關(guān)鍵詞:潰壩洪水;潰壩模型試驗(yàn);洪水演進(jìn);數(shù)值模擬;冰湖潰決洪水

    我國水庫眾多,據(jù)統(tǒng)計現(xiàn)有8.7萬余座。在大壩安全的前提下,水庫的合理運(yùn)行為區(qū)域社會經(jīng)濟(jì)發(fā)展和生態(tài)健康提供了保障,同時保證了區(qū)域的防洪安全。但因超標(biāo)準(zhǔn)洪水、大壩施工質(zhì)量欠佳、地基不良、水庫調(diào)度管理不合理、地震或戰(zhàn)爭等因素作用,一旦大壩發(fā)生潰決,巨大的潰壩水流瞬間高速宣泄而下,將對下游地區(qū)的生命財產(chǎn)和生態(tài)環(huán)境造成毀滅性災(zāi)害。1975年8月4日至7日,在河南省中部,因一場特大暴雨(75·8暴雨),導(dǎo)致板橋和石漫灘2座大型水庫垮壩[1],成為世界水庫垮壩慘劇。凡事預(yù)則立、不預(yù)則廢,人們在不斷積累對洪水的感性認(rèn)識的同時,對潰壩洪水進(jìn)行了大量的研究,主要包括潰壩水流理論研究、潰壩問題的試驗(yàn)研究、潰壩模擬及洪水在下游演進(jìn)等方面。

    1潰壩水流理論研究

    潰壩時庫區(qū)水體突然下泄,造成庫水位陡降與下游水位陡漲,并產(chǎn)生潰壩水流,這種水流的水位和流量在短時間內(nèi)劇烈變化,傳播過程中波形發(fā)生間斷,流態(tài)變化劇烈,通常用淺水方程來描述。對潰壩水流的理論研究主要是尋求淺水波方程的解析解,以精確描述潰壩水流運(yùn)動規(guī)律。

    1871年法國科學(xué)家圣維南提出描述具有自由表面的淺水體中漸變非恒定水流運(yùn)動規(guī)律的圣維南方程組,為潰壩水流理論分析奠定了基礎(chǔ)。

    在瞬間全潰的壩址峰頂流量計算方面,1892年Ritter采取圣維南方程組的特征形式,假定下游為干河床、無河床剪切摩擦作用,不計河床起伏,得出矩形斷面瞬時全潰Ritter解,該解適用于矩形河谷的自由出流情形[2]。1986年伍超[3]等用黎曼方程直接推導(dǎo)出U形斷面潰壩波的簡化解。1938年赫里斯佳洛維奇采用不連續(xù)波理論得出任意斷面的瞬時全潰解,但阿爾漢蓋里斯基指出若下游為干河床會導(dǎo)致流量無窮大的不合理結(jié)果[4];2001年張成林[5-6]發(fā)現(xiàn)該解產(chǎn)生不合理現(xiàn)象的原因是沒有嚴(yán)格依據(jù)物理學(xué)的規(guī)定構(gòu)建動力方程,使得方程中質(zhì)量與作用力的相互關(guān)系不匹配,為此重新選取研究對象構(gòu)建動力方程,推求出修正公式;但2010年趙明登等[7]認(rèn)為該修正公式的推導(dǎo)過程中并未考慮研究對象微小流束的水力要素的時變性,為此修正正負(fù)波相交公式導(dǎo)出適用范圍更廣的潰壩峰頂流量計算公式。1957年Stoker在前人研究的基礎(chǔ)上,不計河床剪切摩擦和河床起伏,但考慮下游有水情形,把壩址流態(tài)劃分為不連續(xù)波流、臨界流和連續(xù)波流3類,給出了矩形斷面與下游起始流速為0的解。1999年C.Wu等[8]將Stoker解的3個公式整合為1個公式,并給出了無量綱查用表格。謝任之[4]全面考慮壩址流態(tài)分區(qū),推導(dǎo)出適用于任意斷面,下游是干床或下游有水的統(tǒng)一公式。

    在瞬間部分潰壩的壩址峰頂流量計算方面,1949年Schoklitch由試驗(yàn)得到了經(jīng)驗(yàn)公式;1951年Frank不計橫向流速,在Ritter瞬時全潰解的基礎(chǔ)上,考慮水量平衡方程,進(jìn)行了公式推導(dǎo)[4]。1982年謝任之根據(jù)波堰流量相等原理,推求出適用于淹沒出流與自由出流的實(shí)用公式。2001年晏鄂川等[9]推求了平底無阻力河道大壩橫向局部潰決的峰頂流量計算公式。2012年方崇惠等[10]推導(dǎo)了計算大壩全潰,橫垂向局部潰壩的瞬時最大流量的同一公式,并給出了通式的上下限。在土石壩的逐漸潰決上,謝任之利用波匯流量相等原理,推導(dǎo)出計算公式,并且考慮庫區(qū)淤積和壩址口門沖深因素,對公式作了修正。

    在壩址流量過程線計算上,Dressler研究了平底矩形渠道考慮摩阻瞬時全潰的一階攝動解。1970年Su等把Dressler的攝動解應(yīng)用于不同形狀斷面明渠[6]。1980年林秉南等[11]采用特征線理論與黎曼解法推求了平底、無阻力,有限長棱柱體水庫下游小水深瞬間全潰問題的理論解。1984年謝任之[4]從不同途徑推導(dǎo)出類似于文獻(xiàn)[11]的解析解。1989年黃真理[12]在文獻(xiàn)[4]和[11]的基礎(chǔ)上,引入梯形斷面因子,以二次拋物線描述無限長梯形河渠瞬間全潰的水面曲線,得出了近似解,并給出了斷面形狀指數(shù)n隨梯形斷面因子變化的查用表格。1996年伍超等[13]采用形態(tài)參數(shù)分離法,克服決口粗糙帶來的計算困難,采用水力參量和一個無因次量關(guān)系的無因次過程表達(dá)水力參量與時間兩者的關(guān)系,給出了任意潰壩決口計算無因次過程線的方法與例子。2006年王立輝等[14]采用拋物型斷面的指數(shù)關(guān)系化簡非恒定明流控制方程,忽略底坡和摩阻影響,推導(dǎo)了無限長明渠瞬時全潰的潰壩涌波解,并分析了濕周斷面形狀系數(shù)對潰壩波各水流參數(shù)的影響。

    在計算逐漸潰壩過程線上,主要要認(rèn)清沖開口門大小的變化與發(fā)展規(guī)律,即要解決潰壩水流的輸沙計算問題,謝任之考慮口門變化,分峰前和峰后分別給出了斜底無阻力逐漸潰的計算公式。2001年鐘德鈺等[15]推導(dǎo)了挾沙水流內(nèi)懸浮顆粒的漂移擴(kuò)散方程,總結(jié)了影響顆粒濃度分布的因素為重力、虛擬質(zhì)量力、升力、Basset力和顆粒間相互作用。2004年王光謙等[16]分析了顆粒垂向脈動強(qiáng)度梯度和升力對沙粒擴(kuò)散的影響并總結(jié)了挾沙水流顆粒的垂向擴(kuò)散機(jī)理。

    潰壩波在下游河道演進(jìn)方面,Hunt將運(yùn)動波近似代替潰壩波,研究離壩有4倍水庫長度的下游河道陡坡河床上的潰壩問題。1985年潘漱芳[17]將Hunt運(yùn)動波解擴(kuò)展為水庫、大壩潰口和下游河渠不等寬時的近似解。2005年黃國富等[18]利用特征線理論研究了斜底,下游無水,梯形、矩形、三角形河道中一維瞬時全潰洪流演進(jìn)的解析模型,并得出正坡促流、逆坡阻流效應(yīng),這種效應(yīng)具有隨著時間增長的特性。2010年杜佐道等[19]探討了以河段槽蓄關(guān)系聯(lián)解圣維南方程組求解洪流演進(jìn)問題。2010年張琨等[20]結(jié)合李斯特萬公式和恒定流計算公式推求出一種計算潰壩洪水水面線的簡化方法,但指出該方法適用于下游出口斷面水位恒定的情況。

    對于淺水方程不考慮流速沿垂線的分布與實(shí)際水流有差別,2003年李榮輝[21]修正淺水方程組中動量方程的底部摩阻項與對流項,推求了適用于模擬潛壩平面流態(tài)的廣義淺水方程。考慮到潰壩洪水波的不連續(xù)性,2004年曾越[22]基于分子運(yùn)動能理論,研究用BGK波爾茲曼理論建立一維潰壩數(shù)值模型的理論基礎(chǔ),并求解了波爾茲曼方程。2005年鮑遠(yuǎn)林[23]研究了無碰撞Boltzmann方程在二維洪水演進(jìn)中的應(yīng)用問題。但是不管是對淺水方程的修正還是用新理論重新建立水流控制方程都還處于研究的初始階段。

    潰壩水流的間斷性和其控制方程的非線性使得推求方程解析解的難度較大。目前潰壩水流理論研究還是基于控制方程在簡化條件下的理論解。今后潰壩水流理論研究依然是潰壩洪水研究的一個重點(diǎn)和難點(diǎn)[24]。

    2 潰壩問題的試驗(yàn)研究

    潰壩試驗(yàn)?zāi)芴剿鳚嗡鞯奈锢頇C(jī)制,檢驗(yàn)數(shù)值模擬的計算成果,能在復(fù)雜地形中給出直觀的概念,但存在費(fèi)時長、投資大、不易制作的缺點(diǎn),在試驗(yàn)中主要存在庫區(qū)模擬問題、變態(tài)率問題、自動控制問題和糙率問題[4]。

    19世紀(jì)中葉法國就開始了潰壩的模型試驗(yàn)研究,20世紀(jì)50年代歐美嚴(yán)重的潰壩事件使?jié)⒖诹髁窟^程研究得到重視,20世紀(jì)50—60年代奧地利與美國分別做了大量室內(nèi)潰壩試驗(yàn),對比兩國的試驗(yàn)結(jié)果發(fā)現(xiàn)潰決時間尺度關(guān)系式接近。

    為配合三峽工程可行性研究和設(shè)計,20世紀(jì)60年代,中國水利水電科學(xué)研究院做了小比例(水平比例為1/30 000、垂直比例為1/300)的三峽變態(tài)模型潰壩試驗(yàn),按壩前庫水位和潰決形式的不同組合,確定了不同情況下的災(zāi)害路線和范圍;20世紀(jì)80年代,長江水利委員會進(jìn)行了大比例(水平比例為1/500、垂直比例為1/125)的三峽大壩變態(tài)模型潰壩試驗(yàn),按潰壩條件與入流的不同組合,進(jìn)行了200多次試驗(yàn)[24]。黃河水利委員會水利科學(xué)研究院于20世紀(jì)70年代就小浪底水庫進(jìn)行了水平比例為1/1 000和垂直比例為1/100的潰壩模型試驗(yàn),總結(jié)了潰壩流量過程經(jīng)驗(yàn)公式[4]。

    20世紀(jì)90年代,歐盟在前期開展的CADAM基礎(chǔ)上繼續(xù)開展IMPACT項目,考慮不確定性因素影響,就潰壩形成、洪水傳播、泥沙運(yùn)動3個方面進(jìn)行過5次大尺度的現(xiàn)場試驗(yàn)和22組小比例的室內(nèi)試驗(yàn),開展最不利洪水條件下的潰壩預(yù)測研究。

    在對潰壩實(shí)例的觀察和試驗(yàn)基礎(chǔ)之上,拉爾斯頓、瓦爾與漢森等研究者就土壩漫潰提出了陡坎沖刷理論[25];羅優(yōu)等[26]進(jìn)行了小尺度均質(zhì)土石壩的漫頂破壞室內(nèi)試驗(yàn),觀測到3種漫頂破壞模式:陡坎蝕退沖刷潰決模式、浸泡剝蝕破壞模式與剪切蝕退坍塌潰決模式,并指明有必要增加尺度和改變臨界條件作進(jìn)一步驗(yàn)證。2008年南京水利科學(xué)研究院與滁州市水利勘測設(shè)計院在滁州南譙區(qū)的大洼水庫就不同壩體土材、不同壩體坡度、壓實(shí)度和含水量進(jìn)行過多組次現(xiàn)場潰壩試驗(yàn),試驗(yàn)結(jié)果驗(yàn)證了“陡坎”理論,在試驗(yàn)結(jié)果的基礎(chǔ)上確定了在壩體開挖泄流槽引導(dǎo)泄流、降低壩前水位的應(yīng)急搶險方案[27]。

    2011年Y.Xue等[28]在陡坡矩形水槽內(nèi)研究梯級水庫瞬時潰壩情形下的水流特點(diǎn),得出水庫初始水深比大壩間距對下游洪峰水深的影響大。2014年C.Y.Chen等[29]通過水槽試驗(yàn)研究不同回水影響下,上游潰壩洪水給下游壩施加的壓力荷載,應(yīng)用線性回歸方法總結(jié)了預(yù)測最大壓力荷載的經(jīng)驗(yàn)公式,并引進(jìn)正弦衰減函數(shù)描述下游壩壓力荷載的時變過程。冰湖潰決也是一個關(guān)注點(diǎn),2014年黃金輝等[30]以光謝錯冰磧湖潰決為原型,利用1∶1 000的水槽系統(tǒng)按不同涌浪條件設(shè)計18組實(shí)驗(yàn),研究由冰崩引起的涌浪沖刷誘發(fā)的潰壩的發(fā)展過程和潰決流量特征。李云等[31]從研究內(nèi)容、模型相似理論、模型比例尺、試驗(yàn)測試手段4個方面對潰壩模型試驗(yàn)研究發(fā)展動向作了詳細(xì)的論述。

    以上大都是基于單庫不同因素組合的潰壩試驗(yàn),深入理解潰壩的機(jī)理還有待進(jìn)一步試驗(yàn)研究,同時,在流域梯級開發(fā)的背景下,需要開展多庫潰壩試驗(yàn)研究。

    3 潰壩模擬及洪水在下游演進(jìn)

    3.1 模擬方法

    潰壩模擬及洪水在下游的演進(jìn)以往是根據(jù)經(jīng)驗(yàn)公式求取,隨著計算機(jī)和存儲技術(shù)的發(fā)展,近些年主要應(yīng)用數(shù)值模擬來求解??蒲泄ぷ髡唛L期研究潰壩洪水,開發(fā)了一些代表性模型,比如DAMBRK模型、BREACH模型、FLDWAV模型、BEED模型、DELFT23D模型、DHI模型系列等。2003年朱勇輝等[25]對10個主要的土壩潰壩模型從水動力學(xué)、泥沙輸移、潰口形狀、所需參數(shù)4個方面以表格形式作了總結(jié)。2007年王立輝等[24]對DAMBRK,BREACH,SMPDBK模型進(jìn)行了介紹。

    據(jù)統(tǒng)計顯示,世界范圍內(nèi)土壩占絕大多數(shù)[1],而土壩的失事率很高,因此對土壩的研究比較深入。土壩的潰壩模型分為基于參數(shù)的模型和基于物理過程的模型。前者利用潰口歷時、潰口寬度等參數(shù)的簡單時變過程或者建立參數(shù)和潰口發(fā)展速度、潰壩峰頂流量間的回歸方程來模擬潰壩過程,此類模型簡單方便,但較粗略。后者通過學(xué)科交叉構(gòu)建時變過程來模擬潰壩過程與洪水過程線,此類模型涉及實(shí)際潰壩機(jī)理,精度較高[25]。2類模型大多采用了牽引應(yīng)力分析和泥沙輸移理論,但T.L.Wahl[32]認(rèn)為采用牽引應(yīng)力分析法應(yīng)該考慮區(qū)分不同類型土質(zhì)對分析結(jié)果的影響,而現(xiàn)有的模型基本上沒能區(qū)分水流的輸沙能力與水流對大壩的沖刷速度2個概念。

    洪水波的數(shù)值模擬主要采用的數(shù)值方法為有限差分法、有限元法、有限體積法。有限差分法先對定義域劃分網(wǎng)格,再在網(wǎng)格點(diǎn)上,利用數(shù)值微分公式將方程用差商代替微商做進(jìn)一步求解,該方法靈活、通用性強(qiáng)、易實(shí)現(xiàn),但其卻常常不能適應(yīng)復(fù)雜的幾何求解域;有限元法基于變分原理和微分方程的弱解形式,其網(wǎng)格劃分能處理幾何形狀復(fù)雜的求解區(qū)域,但須假定網(wǎng)格點(diǎn)間值的變化規(guī)律且計算量比較大;有限體積法結(jié)合了有限差分法靈活和有限元法適應(yīng)復(fù)雜求解區(qū)域的優(yōu)點(diǎn)。王立輝等[24]對以上3種方法在潰壩問題中的應(yīng)用作了分析與比較。此外,馮民權(quán)等[33]采用隨機(jī)游動法求解河渠非恒定流問題。S.C.Chang[34]提出計算守恒方程的時空守恒元和解元方法,此新方法把時空統(tǒng)一起來,局部和全局都遵守守恒律,不需求解黎曼問題且易推廣到多維情形,但是該方法比較難處理源項與求解變量關(guān)系不明確的情況。張增產(chǎn)等[35-37]對該方法作了改進(jìn)與優(yōu)化,張永祥等[38]應(yīng)用改進(jìn)的時空守恒元和解元方法模擬瞬間全潰以及局部潰壩洪水的演進(jìn)與繞流過程。

    3.2 單庫潰壩模擬

    潰壩洪水?dāng)?shù)值模擬的任務(wù)是需要計算出壩址的流量、水位過程線,并向下游作洪水演進(jìn)得到沿程的水位、流量、流速、波前與洪峰到達(dá)時間、洪水持續(xù)時間,其實(shí)質(zhì)是求解動床不恒定流擬線性雙曲偏微分方程組的間斷問題[4]。

    2003年朱勇輝等[25]介紹并分析了土壩潰決機(jī)理學(xué)說——“陡坎”沖刷理論,指明“陡坎”的移動速度dX/dt的預(yù)測是一個難點(diǎn)并且支持該理論的模型試驗(yàn)都只局限于潰口發(fā)展到壩底之前的階段,實(shí)際中潰口還會繼續(xù)擴(kuò)寬。2008年隆文非等[39]有機(jī)結(jié)合物理成因模型、參數(shù)模型與瞬時全潰模型來研究潰壩洪水。2012年沈洋等[40]對金牛山水庫,用MIKE11和MIKE21軟件耦合的方法,同時考慮不同水庫水位與入庫設(shè)計洪水組合情況,計算潰壩洪水及其演進(jìn)過程。

    鑒于水力計算中,齊次形式的淺水方程和歐拉方程在數(shù)學(xué)上有相同的特性,引入計算氣動力學(xué)的研究成果到淺水流動的數(shù)值模擬中[41-42]。這些研究成果能較好處理水流間斷,其中關(guān)鍵是構(gòu)建守恒逆風(fēng)格式,包括3類:①從代數(shù)特征分析出發(fā)的通量向量分裂法(FVS);②應(yīng)用黎曼解的格式,如通量差分裂法(FDS);③對反擴(kuò)散項或者黏性項給以限制,如通量校正輸運(yùn)法(FCT)[43]。胡四一等[44]分別對以上3種高性能格式應(yīng)用到一維非恒定明渠流情況作了介紹和推導(dǎo),并指出對于實(shí)際明渠水流,這種移植須作出一定的修改和處理。從空氣動力學(xué)中引進(jìn)的單純高精度無振蕩捕捉格式雖然可以解決復(fù)雜流態(tài)的過渡和間斷問題,但不能應(yīng)用于地形變化較大的情況,否則會造成計算結(jié)果的不和諧[45]。但在有限體積法里面結(jié)合有限差分法構(gòu)建新的離散格式,對地形變化的處理能力能得到較大改善。如2007年向波等[46]應(yīng)用有限體積法離散二維淺水波方程,并在非結(jié)構(gòu)網(wǎng)格中結(jié)合有限差分法求解界面通量,重力項也用差分法處理,使得界面通量達(dá)到二階精度,并且適用于地形變化較大的情況。另外,1994年Zhao等[47]采用特殊單元水力模型、2001年Zhou等[48]采用水面梯度法、2003年潘存鴻等[49-50]采用水位方程法處理河床的地形起伏問題。

    對笛卡爾坐標(biāo)系進(jìn)行坐標(biāo)變換,可以較好地適應(yīng)河道的復(fù)雜邊界,使模擬更精準(zhǔn)。王船海等采用正交曲線擬合坐標(biāo)解決了自然河道邊界形狀復(fù)雜和長與寬的尺度不協(xié)調(diào)給有限差分帶來的模擬困難。2008年于曰晻等[51]也研究了平面二維流場數(shù)值模擬中正交曲線網(wǎng)格的生成方法。正交曲線變換的計算網(wǎng)格要嚴(yán)格保證正交,而非正交曲線變換卻不受此限制,網(wǎng)格生成比較靈活。2014年M.Guerra等[52]采用非正交曲線坐標(biāo)系模擬非規(guī)則邊界和劇烈變化地形下的水流。另外,隨著計算條件的改善,應(yīng)用間斷有限元方法能處理定解域的不規(guī)則問題[53]。韓濤等[54-55]探討了間斷有限元在淺水方程中的應(yīng)用,并指出間斷有限元與有限體積法的耦合值得進(jìn)一步研究。

    潰壩波的漲落引起水陸邊界變化,產(chǎn)生淺水動邊界問題,數(shù)值模擬中常用的處理方法有凍結(jié)法[56]、切削法[57]、窄縫法[58]與干濕法[59]。宋志堯等[60]對以上4種方法作了評價,并提出了線邊界法。潘存鴻等[61]認(rèn)為如何在動邊界條件下準(zhǔn)確模擬間斷流動,且保證水量平衡是動邊界計算中的最大困難,并將模擬間斷流的干底黎曼解應(yīng)用到求解動邊界問題。2008年張大偉等[62]利用Roe格式的近似黎曼解計算界面通量、用Godunov格式作空間離散、運(yùn)用特征分解與迎風(fēng)處理底坡源項、用半隱式離散方法處理摩阻項,設(shè)定一個較小水深值hmin處理干濕邊界,基于非結(jié)構(gòu)網(wǎng)格構(gòu)建潰壩水流的干濕變化過程數(shù)值模擬。2010年夏軍強(qiáng)等[63]認(rèn)為在干濕界面處理上,對于簡單地形如試驗(yàn)水槽,最小水深值可取的相對較小,而對復(fù)雜地形最小水深宜取相對較大點(diǎn)。

    鑒于有些有限體積模型當(dāng)中,時間離散常采用較簡單的顯格式使計算效率低,而采用隱格式算法復(fù)雜、不易編程,Jameson[64]提出雙時間步方法,雙時間步指虛擬時間步和物理時間步。Helenbrook等[65]曾就一階精度算法對雙時間步法在淺水方程中應(yīng)用的預(yù)處理加速問題進(jìn)行了理論分析。2014年喻海軍等[66]應(yīng)用雙時間步方法模擬了4個典型算例,并與傳統(tǒng)顯式算法比較,結(jié)果顯示,雙時間步方法的時間步長能取得顯式格式10倍或以上。

    基于網(wǎng)格的數(shù)值方法已經(jīng)廣泛應(yīng)用到潰壩洪水模擬當(dāng)中,但是出于對網(wǎng)格的依賴性,如何保證動邊界網(wǎng)格質(zhì)量一直是個難點(diǎn)。無網(wǎng)格粒子法能較好適應(yīng)自由面大變形情況,如張馳等[67]應(yīng)用常用的2種粒子方法,SPH(Smoothed Particle Hydrodynamics)方法和MPS(Moving Particle Semi-Implicit)方法在實(shí)驗(yàn)室進(jìn)行潰壩模擬,并對這2種方法作了比較。姜偉等研究了基于不完全黎曼解的SPH方法。陳善群等[68]用最小二乘無網(wǎng)格有限差分法求解非守恒形式的二維淺水方程,通過圓形潰壩模型檢驗(yàn)該算法能處理對稱間斷流問題。但是無網(wǎng)格粒子法在潰壩洪水?dāng)?shù)值模擬中的應(yīng)用還并不多見,有待進(jìn)一步研究。

    3.3 梯級壩潰決模擬

    梯級潰壩洪水由于上游潰壩洪水和下游再生洪水波的疊加而出現(xiàn)多峰性。2014年黃衛(wèi)等[69]運(yùn)用一維淺水動力學(xué)模型模擬不同工況下梯級大壩潰決洪水。對比分析單壩潰決洪水和梯級壩潰決洪水模擬結(jié)果,發(fā)現(xiàn)后者存在漸進(jìn)增強(qiáng)機(jī)制。具體表現(xiàn)為:峰值流量和水位增大,洪水到達(dá)時間提前。梯級潰壩模擬需在單庫潰壩模擬的基礎(chǔ)上,考慮不同大壩類型、梯級水庫的連鎖響應(yīng)和更加復(fù)雜的洪水組成等問題。如劉慶紅等[70]以某河流上7個梯級水庫按不同大壩類型進(jìn)行了庫群洪水潰壩模擬。

    3.4 冰湖潰決模擬

    冰川活動或者退縮產(chǎn)生的冰雪融水在冰川的前部或者側(cè)部匯集形成冰湖。冰湖可分為多種類型,其中以冰川終磧湖規(guī)模最大、分布數(shù)量最多、災(zāi)害風(fēng)險較高,當(dāng)湖體的沖擊力大于堆積體的抗?jié)⒘r易發(fā)生潰決。2003年Cenderelli等[71]采用一維恒定流模型計算冰湖潰決洪水洪峰流量,2006年Miyamoto等[72]考慮水流的非恒定性,應(yīng)用擴(kuò)散波方程研究冰湖潰決洪水特征,目前冰湖潰決洪水一般按潰壩洪水計算。2007年岳志遠(yuǎn)等[73]基于淺水方程構(gòu)建二維水動力學(xué)模型計算冰湖潰決洪水過程,并指出冰湖潰決洪水容易演變?yōu)槟嗍?,需發(fā)展水沙耦合模型作進(jìn)一步研究。2014年吳秀山[74]對西藏的年楚河流域潛在危險冰湖——黃湖,分6種不同潰決模式,模擬二維潰壩水流及其演進(jìn)。結(jié)果表明:從上游近潰口端到下游河道,不同潰壩模式對水力要素如洪水最大流速、洪峰流量和洪水最大水深的影響逐漸減小,其中近潰口端較敏感;瞬時潰壩工況下的洪峰流量比逐漸潰壩工況大,小潰口工況下的洪峰流量比大潰口工況大,即大潰口不一定對應(yīng)著高洪峰流量。

    4結(jié)論與展望

    潰壩水流理論表達(dá)潰壩水流運(yùn)動的物理規(guī)律,模型試驗(yàn)還原和再現(xiàn)潰壩水流形態(tài),數(shù)值模擬對潰壩水流進(jìn)行模擬和預(yù)測。潰壩洪水研究為潰壩洪水災(zāi)害應(yīng)對提供依據(jù),目前研究的主要困難是高強(qiáng)輸沙理論未建立,口門沖刷過程未能準(zhǔn)確表達(dá)。另外,目前潰壩模型眾多,選用哪個模型更合適,仍是需要研究的問題。由于直接推求潰壩水流控制方程解析解存在困難,目前潰壩水流理論研究還是基于控制方程在簡化條件下的理論解,利用數(shù)學(xué)的新理論新方法,逼近潰壩洪水實(shí)際情況,繼續(xù)加強(qiáng)潰壩水流理論研究依然是一個長期過程。

    現(xiàn)行的潰壩試驗(yàn)多是基于中小尺度下進(jìn)行的,對于潰口發(fā)展的真實(shí)過程不能清晰地反映,需繼續(xù)進(jìn)行單庫潰壩模型試驗(yàn)如大尺度試驗(yàn),由于壩的破壞涉及多相非均質(zhì)問題,需結(jié)合流體力學(xué)、土力學(xué)和計算機(jī)等學(xué)科,加強(qiáng)潰壩機(jī)理研究。此外,現(xiàn)行的潰壩試驗(yàn)大多是基于單庫不同因素組合的潰壩試驗(yàn),開展多庫潰壩模型試驗(yàn)研究也是一個方向。

    現(xiàn)行潰壩洪水?dāng)?shù)值模擬主要存在水流間斷處理、底坡源項和摩阻項處理、干濕動邊界、不規(guī)則邊界擬合、起伏地形處理等技術(shù)問題,針對這些問題不同研究者提出了不同的解決方法,但各有優(yōu)缺點(diǎn),有必要研究新方法更好地解決這類問題,以實(shí)現(xiàn)潰壩水流的精細(xì)模擬。同時,在流域梯級開發(fā)和氣候變化的背景下,積極做好梯級潰壩洪水和冰湖潰決洪水的模擬研究顯得尤為重要。由于潰壩洪水災(zāi)害的應(yīng)急性,建立快捷可靠的區(qū)域潰壩洪水預(yù)報系統(tǒng)將是未來應(yīng)該著重考慮的。

    參考文獻(xiàn):

    [1]劉寧.國內(nèi)外大壩失事分析研究[M].武漢:湖北科學(xué)技術(shù)出版社,2001.

    [2]李煒.水力計算手冊[M].北京:中國水利水電出版社,2006.

    [3]伍超,譚振宏.U形河槽潰壩波的簡化計算[J].重慶交通學(xué)院學(xué)報,1986,(3):166-171.

    [4]謝任之.潰壩水力學(xué)[M].濟(jì)南:山東科學(xué)技術(shù)出版社,1993.

    [5]張成林.對常用的計算潰壩最大流量和波速公式的分析及修正[J].水利水電技術(shù),2001,32(5):7-10.

    [6]張成林,朱春英.對常用河渠立波公式問題的分析及修正[J].水利學(xué)報,1997,(4):79-83.

    [7]趙明登,李靚亮,周湘靈.潰壩最大流量計算公式的問題與修正[J].中國農(nóng)村水利水電,2010,(6):66-68.

    [8]WU Chao, HUANG Guo-fu, ZHENG Yong-hong. Theoretical Solution of Dam-break Shock Wave[J]. Journal of Hydraulic Engineering, 1999, 125(11): 1210-1215.

    [9]晏鄂川,鄭萬模,唐輝明,等.滑坡堵江壩潰決洪水及其演進(jìn)的理論分析[J].水文地質(zhì)工程地質(zhì),2001,(6):15-17.

    [10]方崇惠,方堃.瞬時潰壩最大流量計算新通式推導(dǎo)及驗(yàn)證[J].水科學(xué)進(jìn)展, 2012, 23(5):721-727.

    [11]林秉南,龔振瀛,王連祥.突泄壩址過程線簡化分析[J].清華大學(xué)學(xué)報,1980,20(1):17-31.

    [12]黃真理. 梯形河渠瞬時全潰問題的近似解[J].水動力學(xué)研究與進(jìn)展,1989,4(2):42-52.

    [13]伍超,趙文謙,鄭永紅,等. 任意形狀潰口突泄壩無因次過程線計算方法[J].水利學(xué)報,1996,(3):76-84.

    [14]王立輝,胡四一.拋物型斷面明渠潰壩波的簡化分析[J].水電能源科學(xué),2006,24(1):58-61.

    [15]鐘德鈺,王光謙,王士強(qiáng).挾沙水流中懸浮顆粒的漂移—擴(kuò)散方程[J].水動力學(xué)研究與進(jìn)展,2001,16(1):51-55.

    [16]王光謙,傅旭東.挾沙水流顆粒垂向擴(kuò)散機(jī)理[J].科學(xué)通報,2004,49(4):321-324.

    [17]潘漱芳.潰壩問題的解法[J].東北水利水電,1985,(6):44-48.

    [18]黃國富,張翼飛,伍超,等.一維潰壩波在斜底坡河道中演進(jìn)的解析解[J].水動力學(xué)研究與進(jìn)展,2005,20(5):597-603.

    [19]杜佐道,向德明,馬軍,等.用槽蓄關(guān)系聯(lián)解圣維南方程組模擬洪水演進(jìn)的方法探討[J].水利規(guī)劃與設(shè)計,2010,(6):14-15.

    [20]張琨,陳樹山,華家鵬,等.潰壩洪水水面線的簡化計算方法[J].水電能源科學(xué),2010,28(2):41-43.

    [21]李榮輝.廣義淺水方程及工程應(yīng)用初步研究[D].南京:河海大學(xué),2003.

    [22]曾越.GBK波爾茲曼理論在一維潰壩水流中的應(yīng)用[D].武漢:武漢大學(xué),2004.

    [23]鮑遠(yuǎn)林.有限體積流矢量分裂法在洪水演進(jìn)中的應(yīng)用[D].武漢:華中科技大學(xué),2005.

    [24]王立輝,胡四一.潰壩問題研究綜述[J].水利水電科技進(jìn)展,2007,27(1):80-85.

    [25]朱勇輝,寥鴻志,吳中如.國外土壩潰壩模擬綜述[J].長江科學(xué)院院報,2003,(4):26-29.

    [26]羅優(yōu),陳立,郝婕妤,等.均質(zhì)土石壩不同因素與漫頂破壞模式的內(nèi)在聯(lián)系[J].武漢大學(xué)學(xué)報,2014,4(5):610-614.

    [27]沈登樂.土壩潰壩試驗(yàn)研究與應(yīng)用[J].江淮水利科技,2009,(5):35-37.

    [28]XUE Yang, XU Wei-lin, LUO Shu-jing,etal. Experimental Study of Dam-break Flow in Cascade Reservoirs with Steep Bottom Slope[J]. Journal of Hydrodynamics, Serial B, 2011, 23(4): 491-497.

    [29]CHEN C Y, XU W L, DENG J,etal. Experimental Investigation of Pressure Load Exerted on a Downstream Dam by Dam-break Flow[J]. Journal of Hydraulic Engineering, 2014, 140(2): 199-207.

    [30]黃金輝,劉建康,程尊蘭,等. 涌浪規(guī)模對冰磧湖潰決的影響實(shí)驗(yàn)[J].山地學(xué)報, 2014, 32(2):241-248.

    [31]李云,李君.潰壩模型試驗(yàn)研究綜述[J].水科學(xué)進(jìn)展,2009, 20(2):304-310.

    [32]WAHL T L. Prediction of Embankment Dam Breach Parameters: A Literature Review and Needs Assessment. DSO-98-004, Dam Safety Research Report[R]. US: Bureau of Reclamation, 1998.

    [33]馮民權(quán),趙明登,鄭邦民.河渠非恒定流及其物質(zhì)輸運(yùn)的數(shù)值模擬[M].北京:科學(xué)出版社,2012.

    [34]CHANG S C. The Method of Space-time Conservation Element and Solution Element:A New Approach for Solving the Navier-stokes and Euler Equations[J]. Computational Physics, 1995, 119(2): 295-324.

    [35]張增產(chǎn),沈孟育.改進(jìn)的時空守恒元和解元方法[J].清華大學(xué)學(xué)報,1997,37(8):65-68.

    [36]張增產(chǎn),沈孟育.一種嚴(yán)格保證時空守恒率的數(shù)值方法[J]. 計算物理, 1997,14(6):835-841.

    [37]張增產(chǎn),沈孟育.求解二維Euler方程的時空守恒格式[J].力學(xué)學(xué)報,1999,(2):1-7.

    [38]張永祥,陳景秋.用守恒元和解元法數(shù)值模擬二維潰壩洪水波[J].水利學(xué)報,2005,36(10):1224-1229.

    [39]隆文非,張新華,黃金池,等.水庫潰壩洪水預(yù)測方法研究及應(yīng)用[J].四川大學(xué)學(xué)報,2008,40(1):21-26.

    [40]沈洋,王佳妮. 基于MIKE軟件的潰壩洪水?dāng)?shù)值模擬[J].水電能源科學(xué), 2012, 30(6):56-58.

    [41]汪德爟.計算水力學(xué)理論與應(yīng)用[M].北京:科學(xué)出版社,2011.

    [42]譚維炎.計算淺水動力學(xué)-有限體積法的應(yīng)用[M].北京:清華大學(xué)出版社,1998.

    [43]譚維炎.淺水動力學(xué)的回顧和當(dāng)代前沿問題[J].水科學(xué)進(jìn)展,1999,10(3),296-303.

    [44]胡四一,譚維炎.一維不恒定明流計算的三種高性能差分格式[J].水科學(xué)進(jìn)展, 1991,2(1):11-21.

    [45]潘存鴻.三角形網(wǎng)格下求解二維淺水方程的和諧Godunov格式[J].水科學(xué)進(jìn)展, 2007, 18(2):204-209.

    [46]向波,藍(lán)霄峰,紀(jì)昌明,等. 基于二階差分法和非結(jié)構(gòu)網(wǎng)格的有限體積法的潰壩模擬[J].水動力學(xué)研究與進(jìn)展(A輯),2007,22(6):737-743.

    [47]ZHAO D H, SHEN W H, LAI J S,etal. Finite-volume Two-dimensional Unsteady-flow Model for River Basins[J]. Journal of Hydraulic Engineering, 1994, 120(7): 863-882.

    [48]ZHOU J G, CAUSON D M, MINGHAM C G,etal. The Surface Gradient Method for the Treatment of Source Terms in the Shallow Water Equations[J]. Journal of Computational Physics, 2001, 168(1): 1-25.

    [49]潘存鴻,林炳堯,毛獻(xiàn)忠.一維淺水流動方程的Godunov格式求解[J].水科學(xué)進(jìn)展,2003,14(4):430-436.

    [50]潘存鴻,林炳堯,毛獻(xiàn)忠.求解二維淺水流動方程的Godunov格式[J].水動力學(xué)研究與進(jìn)展,2003,18(1):1-8.

    [51]于曰晻,沈永明,吳修廣,等.正交數(shù)值網(wǎng)格的生成及平面二維流場的數(shù)值模擬[J].計算力學(xué)學(xué)報,2008,25(1):90-93.

    [52]GUERRA M, CIENFUEGOS R, ESCAURIAZA C,etal. Modeling Rapid Flood Propagation Over Natural Terrains Using a Well-balanced Scheme[J]. Journal of Hydraulic Engineering, 2014. DOI: 10.1061/(ASCE)HY.1943-7900.0000881.

    [53]蔚喜軍,周鐵.流體力學(xué)方程的間斷有限元方法[J].計算物理,2005,22(2):109-116.

    [54]韓濤,逄勇,翟金波,等.淺水方程的間斷有限元解法[J].中國農(nóng)村水利水電,2007,(9):23-25.

    [55]韓濤,逄勇,安婷,等.基于間斷有限元求解淺水方程[J].西安交通大學(xué)學(xué)報,2007,41(3):377-379.

    [56]程文輝,王船海.用正交曲線網(wǎng)格及“凍結(jié)”法計算河道流速場[J].水利學(xué)報,1988,(6):18-24.

    [57]彭凱,方鐸,曹叔尤.在二維流動計算中應(yīng)用“河床切削”技術(shù)處理動邊界問題[J].水動力學(xué)研究與進(jìn)展(A輯),1992,7(2):200-205.

    [58]何少苓,王連祥.窄縫法在二維邊界變動水域計算中的應(yīng)用[J].水利學(xué)報,1986,(12):11-19.

    [59]史峰巖,孔亞珍,丁平興.潮灘海域邊界自適應(yīng)網(wǎng)格潮流數(shù)值模型[J].海洋工程,1998,16(3):68-75.

    [60]宋志堯,薛鴻超,嚴(yán)以新.線邊界法在潮流模擬中的應(yīng)用[J].海洋工程,2000,18(4):49-54.

    [61]潘存鴻,于普兵,魯海燕.淺水動邊界的干底 Riemann 解模擬[J].水動力學(xué)研究與進(jìn)展(A輯),2009,24(3):305-312.

    [62]張大偉,李丹勛,王興奎.基于非結(jié)構(gòu)網(wǎng)格的潰壩水流干濕變化過程數(shù)值模擬[J].水力發(fā)電學(xué)報,2008,27(5):98-102.

    [63]夏軍強(qiáng),王光謙,LIN Bin-lian,等.復(fù)雜邊界及實(shí)際地形上潰壩洪水流動過程模擬 [J].水科學(xué)進(jìn)展,2010,21(3):289-298.

    [64]JAMESON A. Time Dependent Calculations Using Multigrid, with Applications to Unsteady Flows Past Airfoils and Wings[C]∥Proceedings of AIAA 10th Computational Fluid Dynamics Conference. AIAA, Honolulu, June 24-27, 1991: 1-13.

    [65]HELENBROOK B T, COWLES G W. Preconditioning for Dual-time-stepping Simulations of the Shallow Water Equations Including Coriolis and Bed Friction Effects[J]. Journal of Computational Physics, 2008, 227(9): 4425-4440.

    [66]喻海軍,黃國如,武傳號. 雙時間步法在二維淺水方程求解中的應(yīng)用[J].水科學(xué)進(jìn)展,2014,25(4):542-549.

    [67]張馳,張雨新,萬德成.SPH方法和MPS方法模擬潰壩問題的比較分析[J].水動力學(xué)研究與進(jìn)展,2011,26(6):736-746.

    [68]陳善群,廖斌,李海峰. 用最小二乘無網(wǎng)格有限差分方法求解二維淺水方程[J].長江科學(xué)院院報,2012,29(6):36-39.

    [69]黃衛(wèi),曹志先.梯級大壩潰決洪水漸進(jìn)增強(qiáng)機(jī)制數(shù)值模擬[J].武漢大學(xué)學(xué)報,2014,47(2):160-164.

    [70]劉慶紅,付湘,何海源,等.梯級水庫群洪水潰壩模擬[J].中國農(nóng)村水利水電,2011,(3):173-177.

    [71]CENDERELLI D A, WOHL E E. Flow Hydraulics and Geomorphic Effects of Glacial-lake Outburst Floods in the Mount Everest Region, Nepal[J]. Earth Surface Processes and Landforms, 2003, 28(4): 385-407.

    [72]MIYAMOTO H, ITOH K, KOMATSU G,etal. Numerical Simulations of Large-scale Cataclysmic Flood Water: A Simple Depth-averaged Model and an Illustrative Application[J]. Geomorphology, 2006, 76(1/2): 179-192.

    [73]岳志遠(yuǎn),曹志先,車濤,等.冰湖潰決洪水的二維水動力學(xué)數(shù)值模擬[J].冰川凍土,2007,29(5):756-763.

    [74]吳秀山.不同潰決模式下冰湖潰壩洪水演進(jìn)模擬[D].杭州:浙江大學(xué),2014.

    (編輯:劉運(yùn)飛)

    LIU Lin1,CHANG Fu-xuan1,XIAO Chang-wei2,GAO Yan-hong2,CHENG Ling-she2

    (1.Water Resources Department,Yangtze River Scientific Research Institute,Wuhan430010,China;

    2.Hydropower Planning Survey and Design Institute of Tibet Autonomous Region,Lhasa850000,China)

    Research Progress on Dam-break Flood

    Abstract:The purpose of research on dam-break flood is to calculate the dam-break flow, water level hydrograph and the downstream flood flow, flow velocity, water level, arrival time of the wave front and peak along the flow path, to evaluate the loss caused by downstream flood and take measures to reduce flood risk. A review was made on the study of dam-break flood from aspects of theoretical analysis, physical experiment and numerical simulation on dam failure and flood routing in the downstream. The developments, achievements and progresses of researches on dam-break flood in recent years in China and abroad are summarized, and the prospects and focus were prospected. At present, the numerical solution of the dam break theory is developing rapidly, the dam-break mechanism has been proposed according to dam-break experiments, the research has become increasingly modelled. But it’s hard to accurately express the dyke erosion process while high strength sediment transport theory has not been established yet, and also the research about cascade dam-break floods and glacier lake outburst floods is rare. In the future, we should strengthen the theoretical research, conduct large-scale model test for dam breach of cascade reservoirs, improve simulation of cascade dam-break floods and glacier lake outburst floods, and further build a fast and reliable forecasting system for regional dam-break flood.

    Key words:dam-break flood; dam-break experiment; flood routing; numerical simulation; glacier lake outburst flood

    收稿日期:2015-03-27;修回日期:2015-04-23

    基金項目:水利部公益性行業(yè)科研專項( 201301037)

    作者簡介:劉林( 1991 - ) ,男,江西萍鄉(xiāng)人,碩士研究生,主要從事水文水資源研究,(電話) 18627196150( 電子信箱) linzi3720@163.com。

    doi:10.11988/ckyyb.20150225

    中圖分類號:TV122.4

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

    文章編號:1001-5485(2016)06-0029-07

    2016,33(06):29-35

    猜你喜歡
    數(shù)值模擬
    基于AMI的雙色注射成型模擬分析
    錐齒輪精密冷擺輾成形在“材料成型數(shù)值模擬”課程教學(xué)中的應(yīng)用
    基于氣象信息及風(fēng)場信息的風(fēng)機(jī)輪轂處風(fēng)速預(yù)測
    鉆孔灌注樁樁底沉渣對樁體承載特性影響的模擬分析
    西南地區(qū)氣象資料測試、預(yù)處理和加工研究報告
    科技資訊(2016年18期)2016-11-15 08:01:18
    張家灣煤礦巷道無支護(hù)條件下位移的數(shù)值模擬
    科技視界(2016年18期)2016-11-03 23:14:27
    張家灣煤礦開切眼錨桿支護(hù)參數(shù)確定的數(shù)值模擬
    科技視界(2016年18期)2016-11-03 22:57:21
    跨音速飛行中機(jī)翼水汽凝結(jié)的數(shù)值模擬研究
    科技視界(2016年18期)2016-11-03 20:38:17
    姚橋煤礦采空區(qū)CO2防滅火的數(shù)值模擬分析
    雙螺桿膨脹機(jī)的流場數(shù)值模擬研究
    科技視界(2016年22期)2016-10-18 14:53:19
    国产成人精品无人区| videosex国产| 18禁动态无遮挡网站| 欧美 亚洲 国产 日韩一| 黄色怎么调成土黄色| 久久ye,这里只有精品| 丝袜脚勾引网站| 天堂俺去俺来也www色官网| 99热全是精品| 观看美女的网站| 日韩欧美精品免费久久| 777久久人妻少妇嫩草av网站| 国产一级毛片在线| 成年动漫av网址| 国产又爽黄色视频| 日韩欧美一区视频在线观看| 日韩电影二区| 久久av网站| 深夜精品福利| av女优亚洲男人天堂| 精品一区二区三区四区五区乱码 | 国产黄色视频一区二区在线观看| 777久久人妻少妇嫩草av网站| 老司机在亚洲福利影院| 久久狼人影院| 黄频高清免费视频| 极品人妻少妇av视频| 亚洲色图综合在线观看| 久久久久久久大尺度免费视频| 欧美黄色片欧美黄色片| 久久性视频一级片| 国产精品女同一区二区软件| 久久久久久免费高清国产稀缺| 久久天堂一区二区三区四区| 美女大奶头黄色视频| av网站在线播放免费| 中文字幕精品免费在线观看视频| 日本欧美视频一区| 中文字幕人妻丝袜制服| 久久久久国产一级毛片高清牌| 国产高清国产精品国产三级| 亚洲综合色网址| 亚洲精品国产av蜜桃| 亚洲熟女毛片儿| 麻豆乱淫一区二区| 制服人妻中文乱码| 日本黄色日本黄色录像| 亚洲视频免费观看视频| 日本欧美视频一区| 欧美老熟妇乱子伦牲交| 亚洲成av片中文字幕在线观看| 亚洲七黄色美女视频| 街头女战士在线观看网站| xxx大片免费视频| 如何舔出高潮| 美女视频免费永久观看网站| 青春草视频在线免费观看| 午夜免费男女啪啪视频观看| 国产精品av久久久久免费| 免费女性裸体啪啪无遮挡网站| 叶爱在线成人免费视频播放| 亚洲成人av在线免费| 久久久久久久久久久久大奶| 少妇被粗大猛烈的视频| 亚洲美女视频黄频| 亚洲国产中文字幕在线视频| 在线观看国产h片| 激情五月婷婷亚洲| 黄片播放在线免费| 久久97久久精品| 日韩 亚洲 欧美在线| 天天影视国产精品| 亚洲av日韩在线播放| 精品午夜福利在线看| 国产成人a∨麻豆精品| 日韩中文字幕欧美一区二区 | 国产一区二区三区av在线| 日韩一本色道免费dvd| 国产精品 国内视频| 最新的欧美精品一区二区| 人人妻人人澡人人看| 国产 一区精品| 国产精品秋霞免费鲁丝片| 各种免费的搞黄视频| 高清av免费在线| 日日摸夜夜添夜夜爱| 91精品伊人久久大香线蕉| 亚洲婷婷狠狠爱综合网| 久久国产亚洲av麻豆专区| 一区二区三区四区激情视频| 一本大道久久a久久精品| 热99国产精品久久久久久7| 久久久久久免费高清国产稀缺| 丁香六月天网| h视频一区二区三区| 在线观看免费视频网站a站| 亚洲成人国产一区在线观看 | bbb黄色大片| 91国产中文字幕| 69精品国产乱码久久久| 日本wwww免费看| 成年女人毛片免费观看观看9 | 久久久国产欧美日韩av| 亚洲精华国产精华液的使用体验| 欧美在线一区亚洲| 9色porny在线观看| 九九爱精品视频在线观看| 天天添夜夜摸| 777米奇影视久久| 老司机深夜福利视频在线观看 | 黄色怎么调成土黄色| 成人漫画全彩无遮挡| 免费观看性生交大片5| 久久精品久久精品一区二区三区| 我要看黄色一级片免费的| 99热全是精品| av在线app专区| 亚洲国产日韩一区二区| 日本色播在线视频| 一区二区三区乱码不卡18| 精品一区二区三卡| 欧美日韩成人在线一区二区| 嫩草影院入口| 国产免费又黄又爽又色| 这个男人来自地球电影免费观看 | 久久精品人人爽人人爽视色| 久久久国产一区二区| 久久人人爽人人片av| 中文字幕人妻丝袜一区二区 | 久久狼人影院| 亚洲天堂av无毛| 久久韩国三级中文字幕| 青春草视频在线免费观看| 丝袜人妻中文字幕| 91精品伊人久久大香线蕉| 国产成人系列免费观看| 婷婷色麻豆天堂久久| 日韩成人av中文字幕在线观看| 久久久久久免费高清国产稀缺| 久久狼人影院| 欧美精品亚洲一区二区| 久久国产亚洲av麻豆专区| 2021少妇久久久久久久久久久| 国产精品熟女久久久久浪| 国产亚洲午夜精品一区二区久久| 久久 成人 亚洲| 丝袜美足系列| 天天操日日干夜夜撸| 十八禁人妻一区二区| 视频区图区小说| 日本猛色少妇xxxxx猛交久久| 欧美日韩视频高清一区二区三区二| 2021少妇久久久久久久久久久| 伊人亚洲综合成人网| 午夜免费观看性视频| 久久精品亚洲av国产电影网| 最近最新中文字幕免费大全7| 欧美另类一区| 欧美xxⅹ黑人| av免费观看日本| 一区二区av电影网| 2018国产大陆天天弄谢| av在线播放精品| 一边摸一边抽搐一进一出视频| 久久人人97超碰香蕉20202| 欧美日韩综合久久久久久| xxx大片免费视频| 韩国高清视频一区二区三区| 欧美激情高清一区二区三区 | 亚洲美女搞黄在线观看| 日韩人妻精品一区2区三区| 我的亚洲天堂| 在线观看一区二区三区激情| 欧美精品一区二区免费开放| 2021少妇久久久久久久久久久| 宅男免费午夜| 女人精品久久久久毛片| 一边摸一边做爽爽视频免费| 欧美激情高清一区二区三区 | 午夜福利免费观看在线| 一边亲一边摸免费视频| 18禁观看日本| 国产精品 国内视频| 啦啦啦在线免费观看视频4| 在线观看一区二区三区激情| 欧美日韩一级在线毛片| av女优亚洲男人天堂| 精品人妻熟女毛片av久久网站| 桃花免费在线播放| 亚洲三区欧美一区| 欧美日韩一区二区视频在线观看视频在线| 综合色丁香网| 在线观看免费视频网站a站| 观看av在线不卡| 热99国产精品久久久久久7| 亚洲男人天堂网一区| 菩萨蛮人人尽说江南好唐韦庄| 日韩 亚洲 欧美在线| av有码第一页| 国产精品三级大全| 女人高潮潮喷娇喘18禁视频| 亚洲精品,欧美精品| 欧美日韩国产mv在线观看视频| 久热这里只有精品99| 久久精品久久久久久噜噜老黄| 在线 av 中文字幕| 成人免费观看视频高清| 成人漫画全彩无遮挡| 精品人妻一区二区三区麻豆| 成人手机av| 欧美人与性动交α欧美软件| 国产成人精品福利久久| 久久久久精品国产欧美久久久 | 王馨瑶露胸无遮挡在线观看| 亚洲久久久国产精品| 成年女人毛片免费观看观看9 | 日本91视频免费播放| 在线观看免费日韩欧美大片| 中文乱码字字幕精品一区二区三区| 不卡视频在线观看欧美| 黑人猛操日本美女一级片| 999久久久国产精品视频| 交换朋友夫妻互换小说| 丰满少妇做爰视频| 亚洲成av片中文字幕在线观看| 中文字幕色久视频| 熟妇人妻不卡中文字幕| 黄片无遮挡物在线观看| 国产精品一区二区精品视频观看| 精品国产乱码久久久久久男人| 久久人人爽av亚洲精品天堂| 天天影视国产精品| 少妇人妻 视频| 少妇 在线观看| 国产av精品麻豆| 两个人看的免费小视频| 汤姆久久久久久久影院中文字幕| 欧美日韩成人在线一区二区| 欧美国产精品va在线观看不卡| 久久鲁丝午夜福利片| 中文字幕色久视频| 老司机在亚洲福利影院| 青春草视频在线免费观看| 大陆偷拍与自拍| 精品一区二区三区四区五区乱码 | 精品国产乱码久久久久久小说| 多毛熟女@视频| 精品久久久精品久久久| videos熟女内射| 色网站视频免费| 中国三级夫妇交换| 久久久久精品性色| 午夜免费男女啪啪视频观看| 欧美xxⅹ黑人| 久久久久视频综合| 老司机亚洲免费影院| 人人妻人人爽人人添夜夜欢视频| av片东京热男人的天堂| a 毛片基地| 午夜福利,免费看| 在线观看免费日韩欧美大片| 亚洲美女视频黄频| 一本久久精品| 午夜影院在线不卡| 欧美在线一区亚洲| 欧美另类一区| 国产无遮挡羞羞视频在线观看| 深夜精品福利| 日韩伦理黄色片| 又黄又粗又硬又大视频| 亚洲欧美色中文字幕在线| 久久99一区二区三区| 欧美少妇被猛烈插入视频| 韩国av在线不卡| 激情视频va一区二区三区| 亚洲精品第二区| 久久狼人影院| 麻豆av在线久日| 久久久精品94久久精品| 777久久人妻少妇嫩草av网站| 久久毛片免费看一区二区三区| 十八禁人妻一区二区| 日本av手机在线免费观看| 国产深夜福利视频在线观看| 伦理电影大哥的女人| 亚洲国产av影院在线观看| 人妻人人澡人人爽人人| 国产日韩欧美在线精品| 人成视频在线观看免费观看| 一级片'在线观看视频| 一二三四在线观看免费中文在| 欧美变态另类bdsm刘玥| 91国产中文字幕| 天天躁日日躁夜夜躁夜夜| 考比视频在线观看| 国产一区二区三区av在线| 成年人午夜在线观看视频| 亚洲成人av在线免费| 日日爽夜夜爽网站| 啦啦啦 在线观看视频| 精品一区在线观看国产| 亚洲人成网站在线观看播放| 大码成人一级视频| 国产成人免费观看mmmm| 亚洲一级一片aⅴ在线观看| 中文字幕制服av| 国产精品偷伦视频观看了| 丝袜脚勾引网站| 欧美xxⅹ黑人| 老司机亚洲免费影院| 美女脱内裤让男人舔精品视频| 久久 成人 亚洲| 国产精品二区激情视频| 国产精品一区二区在线不卡| kizo精华| 国产日韩一区二区三区精品不卡| 国产午夜精品一二区理论片| 亚洲欧美成人综合另类久久久| 校园人妻丝袜中文字幕| 女人精品久久久久毛片| 免费在线观看视频国产中文字幕亚洲 | 韩国精品一区二区三区| 国产免费福利视频在线观看| 熟女少妇亚洲综合色aaa.| 少妇 在线观看| 国产在线免费精品| 又粗又硬又长又爽又黄的视频| 热re99久久精品国产66热6| 精品国产露脸久久av麻豆| av网站在线播放免费| 一级黄片播放器| 久久久久久久久久久免费av| 最近最新中文字幕免费大全7| 少妇 在线观看| 伦理电影大哥的女人| 亚洲精品国产一区二区精华液| 亚洲熟女毛片儿| 国产有黄有色有爽视频| 国产精品女同一区二区软件| 一个人免费看片子| 只有这里有精品99| 欧美 亚洲 国产 日韩一| 亚洲精品美女久久av网站| 亚洲欧洲日产国产| 97在线人人人人妻| 一区二区三区精品91| 老司机深夜福利视频在线观看 | 晚上一个人看的免费电影| 精品福利永久在线观看| 嫩草影院入口| 亚洲成av片中文字幕在线观看| 美女扒开内裤让男人捅视频| 国产欧美日韩一区二区三区在线| 亚洲精品日韩在线中文字幕| 高清欧美精品videossex| 亚洲,欧美,日韩| 一区二区三区精品91| 好男人视频免费观看在线| 99re6热这里在线精品视频| 色婷婷av一区二区三区视频| 一二三四在线观看免费中文在| 精品人妻熟女毛片av久久网站| 日本wwww免费看| 美女视频免费永久观看网站| 这个男人来自地球电影免费观看 | 日韩一本色道免费dvd| 午夜91福利影院| 一本一本久久a久久精品综合妖精| av又黄又爽大尺度在线免费看| 免费看av在线观看网站| 男女之事视频高清在线观看 | 视频区图区小说| 97精品久久久久久久久久精品| 国产成人精品在线电影| 欧美中文综合在线视频| 亚洲av在线观看美女高潮| 亚洲精品自拍成人| 欧美 亚洲 国产 日韩一| 久久久久精品人妻al黑| 免费观看性生交大片5| 国产一区二区三区综合在线观看| 一二三四中文在线观看免费高清| 天天添夜夜摸| 一本久久精品| 欧美在线一区亚洲| 欧美成人午夜精品| 哪个播放器可以免费观看大片| 美女视频免费永久观看网站| 亚洲欧洲国产日韩| av在线app专区| 少妇的丰满在线观看| 国产精品久久久av美女十八| 久久久国产一区二区| 亚洲精品久久成人aⅴ小说| 夫妻午夜视频| a级片在线免费高清观看视频| 搡老乐熟女国产| 国产成人一区二区在线| 高清不卡的av网站| 夫妻性生交免费视频一级片| 亚洲欧美一区二区三区久久| 久久鲁丝午夜福利片| 国产亚洲午夜精品一区二区久久| www.自偷自拍.com| 亚洲精品久久久久久婷婷小说| 中文乱码字字幕精品一区二区三区| 色94色欧美一区二区| 黑人猛操日本美女一级片| 男女午夜视频在线观看| 最新的欧美精品一区二区| 我要看黄色一级片免费的| 岛国毛片在线播放| 激情五月婷婷亚洲| 国产精品一二三区在线看| 岛国毛片在线播放| 国产高清国产精品国产三级| 久久亚洲国产成人精品v| 精品一区二区三区av网在线观看 | 国产成人精品久久久久久| 久久精品久久精品一区二区三区| 尾随美女入室| 毛片一级片免费看久久久久| 天堂中文最新版在线下载| a级毛片黄视频| 麻豆av在线久日| 久久久精品区二区三区| 婷婷色综合大香蕉| 黄网站色视频无遮挡免费观看| av一本久久久久| 亚洲精品国产区一区二| 精品国产乱码久久久久久小说| 亚洲av在线观看美女高潮| 一本—道久久a久久精品蜜桃钙片| 国产又爽黄色视频| 国产亚洲最大av| 亚洲,欧美,日韩| 国产1区2区3区精品| 久久国产精品大桥未久av| av不卡在线播放| a级毛片在线看网站| 亚洲欧美激情在线| 少妇被粗大猛烈的视频| 在线观看人妻少妇| 欧美少妇被猛烈插入视频| 高清视频免费观看一区二区| 亚洲婷婷狠狠爱综合网| 亚洲一码二码三码区别大吗| 亚洲精品成人av观看孕妇| 黄片无遮挡物在线观看| 乱人伦中国视频| 日韩中文字幕视频在线看片| 亚洲欧美成人精品一区二区| 亚洲精品av麻豆狂野| 99精国产麻豆久久婷婷| 日韩 亚洲 欧美在线| 如日韩欧美国产精品一区二区三区| 自线自在国产av| 国产成人欧美在线观看 | 午夜福利,免费看| 久久 成人 亚洲| 肉色欧美久久久久久久蜜桃| 亚洲五月色婷婷综合| 午夜av观看不卡| 欧美精品av麻豆av| 国产成人一区二区在线| 国产免费一区二区三区四区乱码| 国产福利在线免费观看视频| 国产乱来视频区| 日韩伦理黄色片| 少妇精品久久久久久久| 欧美精品av麻豆av| 看十八女毛片水多多多| 中文乱码字字幕精品一区二区三区| 汤姆久久久久久久影院中文字幕| 成年av动漫网址| 成人国语在线视频| av视频免费观看在线观看| 我的亚洲天堂| 在线天堂最新版资源| 日韩 亚洲 欧美在线| 国产熟女午夜一区二区三区| 亚洲成人免费av在线播放| 男女边吃奶边做爰视频| 日本黄色日本黄色录像| 超碰成人久久| 亚洲人成电影观看| 人人妻,人人澡人人爽秒播 | 亚洲精品av麻豆狂野| 免费女性裸体啪啪无遮挡网站| 一级a爱视频在线免费观看| 欧美人与性动交α欧美精品济南到| 国产亚洲精品第一综合不卡| 国产亚洲最大av| 国产精品.久久久| 99久国产av精品国产电影| 最新的欧美精品一区二区| 久久人人97超碰香蕉20202| 色综合欧美亚洲国产小说| 欧美激情 高清一区二区三区| 亚洲成人手机| 母亲3免费完整高清在线观看| 老司机靠b影院| 亚洲精品乱久久久久久| 久久99热这里只频精品6学生| 性色av一级| 国产精品秋霞免费鲁丝片| 精品视频人人做人人爽| 久久精品久久精品一区二区三区| 国产成人精品久久久久久| 久久久精品区二区三区| 国产成人91sexporn| 午夜免费观看性视频| 亚洲,一卡二卡三卡| 国产免费福利视频在线观看| 亚洲精品国产区一区二| 少妇的丰满在线观看| 久久av网站| 人妻人人澡人人爽人人| 只有这里有精品99| videos熟女内射| 国产激情久久老熟女| 秋霞伦理黄片| 精品一区在线观看国产| 国产99久久九九免费精品| 国产精品久久久人人做人人爽| 一级片'在线观看视频| 国产精品久久久av美女十八| 中国三级夫妇交换| av在线老鸭窝| 国产麻豆69| 日日啪夜夜爽| e午夜精品久久久久久久| 新久久久久国产一级毛片| 精品一区二区三区四区五区乱码 | 国语对白做爰xxxⅹ性视频网站| 亚洲成人国产一区在线观看 | 伦理电影免费视频| 久久久精品免费免费高清| 亚洲国产中文字幕在线视频| 国产精品麻豆人妻色哟哟久久| 国产成人免费观看mmmm| 国产亚洲av片在线观看秒播厂| 色吧在线观看| 亚洲欧美色中文字幕在线| 777米奇影视久久| 成年动漫av网址| 欧美国产精品一级二级三级| 国产成人精品福利久久| 国产黄色视频一区二区在线观看| 亚洲美女视频黄频| 国产不卡av网站在线观看| 99久久99久久久精品蜜桃| 日韩大码丰满熟妇| 精品人妻一区二区三区麻豆| 国产一区二区 视频在线| 中文字幕另类日韩欧美亚洲嫩草| 精品国产超薄肉色丝袜足j| 最新在线观看一区二区三区 | 精品国产一区二区三区久久久樱花| 国产精品国产三级国产专区5o| 国产97色在线日韩免费| 男女床上黄色一级片免费看| 五月天丁香电影| 免费在线观看黄色视频的| 成人午夜精彩视频在线观看| 午夜91福利影院| 国产精品成人在线| 日韩视频在线欧美| 啦啦啦视频在线资源免费观看| 电影成人av| 女人被躁到高潮嗷嗷叫费观| 多毛熟女@视频| 久久国产精品大桥未久av| 激情视频va一区二区三区| 久久久久国产一级毛片高清牌| 啦啦啦 在线观看视频| 欧美精品亚洲一区二区| 菩萨蛮人人尽说江南好唐韦庄| 久久亚洲国产成人精品v| 不卡av一区二区三区| 日韩一本色道免费dvd| 丰满饥渴人妻一区二区三| 母亲3免费完整高清在线观看| 国产精品国产三级专区第一集| 色播在线永久视频| 亚洲成色77777| 婷婷色综合www| 久久久亚洲精品成人影院| 亚洲欧美激情在线| 免费观看a级毛片全部| 亚洲精品国产av蜜桃| 桃花免费在线播放| 久久狼人影院| 制服人妻中文乱码| 九草在线视频观看| 我要看黄色一级片免费的| 国产精品久久久久久人妻精品电影 | 久久久久久久久免费视频了| 国产福利在线免费观看视频| 黑人巨大精品欧美一区二区蜜桃| 丰满迷人的少妇在线观看| 亚洲一码二码三码区别大吗| 乱人伦中国视频| 欧美日韩视频高清一区二区三区二| 人成视频在线观看免费观看| 97在线人人人人妻| 最近中文字幕高清免费大全6| 一级,二级,三级黄色视频| av在线老鸭窝| 久久国产亚洲av麻豆专区| 亚洲精品,欧美精品| 精品视频人人做人人爽| 一级毛片我不卡| 亚洲色图 男人天堂 中文字幕| 两个人免费观看高清视频|