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

    淺水方程源項處理的研究進展

    2021-08-11 04:35許仁義汪凌翔王遠見方浩川王俊賢戴惠東
    人民黃河 2021年7期

    許仁義 汪凌翔 王遠見 方浩川 王俊賢 戴惠東

    摘 要:淺水方程由納維-斯托克斯方程推導而來,是一種描述具有自由表面的淺水體在重力作用下流動的數(shù)學模型。實際淺水流動不可能是完全理想的情況,必須考慮底坡和摩阻等源項的存在,它們影響計算的穩(wěn)定性和精確性,如果處理不當還會造成格式的不和諧,因此源項處理是求解淺水流動方程的關鍵,也是國內外學者研究的一個重要方向。介紹目前求解淺水方程存在的4個主要困難,介紹淺水方程的各種離散方法以及優(yōu)缺點,闡述源項處理的重要性和源項的分類,回顧和總結國內外處理源項的方法以及目前存在的困難。根據目前的研究進展,源項處理依舊有廣闊的研究前景,如何平衡好源項處理方法的正效應和負效應是一個值得研究的課題。

    關鍵詞:淺水方程;源項處理;和諧性;有限體積法;不規(guī)則地形

    中圖分類號:TV131 文獻標志碼:A

    doi:10.3969/j.issn.1000-1379.2021.07.007

    引用格式:許仁義,汪凌翔,王遠見,等.淺水方程源項處理的研究進展[J].人民黃河,2021,43(7):35-40,83.

    Abstract: The shallow water equation, derived from the Navier-Stokes equation is a mathematical model to describe the flow of shallow water with a free surface under the action of gravity. Actual shallow water flow cant be a perfectly ideal situation, bottom slope and friction must always be considered such as the existence of the source term. It affects the calculations stability and precision, if handles not properly it will cause disharmony in the format. Thus, the source term processing is the key to solve the shallow water flow equation and also is an important direction in the research of scholars at home and abroad. This paper introduced four main difficulties in solving shallow water equation and various discrete methods of shallow water equation as well as their advantages and disadvantages, expounded the importance of source term treatment and the classification of source term and reviewed and summarized the methods of source term treatment at home and abroad as well as the existing difficulties. Finally, according to the current research progress, it put forward that the source term treatment still had a broad research prospect and how to balance the positive and negative effects of source term treatment would be a worthy topic to be studied.

    Key words: shallow water equation; source term processing; well-balance; finite volume method; irregular terrain

    水力學中的淺水是指水深尺度遠小于平面尺度且垂向流速小的水流[1],人們關心的潰壩問題、水環(huán)境污染問題、潮汐和涌浪等都可以用淺水流動描述,而淺水方程是描述各類淺水流動的數(shù)學形式。如今,隨著計算機性能大幅提升,淺水數(shù)值模擬以其無可替代的優(yōu)勢開始廣泛應用于河道、河口、水庫以及近海水域等環(huán)境。淺水方程是一種非線性雙曲型方程,如何更加精確且高效地求解是一個值得研究的課題,因此衍生出許多離散方式。由于現(xiàn)實環(huán)境的復雜性,為了得到更加真實的數(shù)值結果,必須考慮源項的存在,因此探求源項處理方法成為了許多學者研究的方向。

    1 控制方程

    控制方程為一維淺水方程守恒形式,滿足靜壓假定。

    式中:守恒變量向量U=hhu;通量向量F(U)=huhu2+12gh2;源項向量S(U)=Sb+Sf+Sw+Sc,如果忽略表面風應力和柯氏力,則S(U)=Sb+Sf=0gh(-Sbx-Sfx),其中Sbx=zx;Sfx為沿x方向的阻力;Sb、Sf、Sw、Sc分別為底坡源項、摩阻源項、表面風應力、柯氏力;z為河床曲面高程;u為x方向的平均流速;h為水深;g為重力加速度。

    2 求解淺水方程存在的困難

    淺水方程屬于非線性雙曲型偏微分方程,滿足質量守恒、動量守恒和能量守恒,這種物理量守恒的方程也可以被稱為雙曲守恒方程[2]。求解非線性偏微分方程組十分困難,很難獲得解析解,為了得到更精確的數(shù)值解,一般有兩種辦法:一種是提高計算機的算力,另一種就是不斷改良數(shù)值方法、提高計算精度和計算效率。

    目前數(shù)值求解淺水方程這類雙曲守恒方程存在以下困難:

    (1)解的間斷。擾動波傳播速度有限,可能產生間斷,在間斷處的導數(shù)無意義,無法求解微分方程。引入“弱解”的概念[3],即弱解在間斷點外的其他點上滿足微分方程,在間斷點上滿足一組跳躍條件(Rankine-Hugoniot條件)[4],它能夠將間斷點兩側的水力要素聯(lián)系起來。

    (2)弱解的非唯一性。淺水方程在推導的時候就已經做出前提假設,同時忽略了某些物理量帶來的影響,這就會導致模型與實際情況存在偏差[5]。

    (3)解的穩(wěn)定性。更高的模型精度以及求解時存在的誤差容易引起數(shù)值振蕩,因此需要優(yōu)化和改善高精度計算方法的數(shù)值振蕩問題[6]。

    (4)源項的處理。淺水方程如果不包含源項,則為齊次淺水流動方程,這類方程由于略去了地形變化、摩阻等計算項,形式上與空氣動力學中使用廣泛的歐拉方程類似[7],因此早期研究階段借鑒了空氣動力學中的求解方法[8]。比如基于近似Riemann解Godunov格式的Roe方法[9]就是一種最初用于空氣動力學的方法。但在水流數(shù)值模擬中,計算區(qū)域往往都是非平面的,對源項簡單處理會導致計算結果的不和諧,即靜水條件下,計算趨于穩(wěn)定后,無法滿足流速為零、水位為常數(shù)的結果[10]。下面將著重闡述源項的處理方法。本文提出了齊次淺水方程的解決方法,而非齊次淺水方程的問題有待解決,比如含源項的非齊次淺水方程的大時間步長格式現(xiàn)在仍未得到解決[11]。

    3 離散方法

    3.1 特征線法(MOC)

    特征線法是利用特征線以及特征相容關系進行離散得到數(shù)值解的方法,在計算機出現(xiàn)以前這還是CFD手工算法之一。它的物理意義明確,計算精度較高,但是不適合帶有源項的非齊次方程,目前很少直接采用該方法[12]。

    3.2 有限差分法(FDM)

    有限差分法就是在計算域上將原來控制方程中網格節(jié)點的微分項直接用差商近似代替,然后進行泰勒級數(shù)展開,可以在每個網格節(jié)點處得到離散方程,因此該方法簡單成熟,可以構造高精度的格式,但是對復雜網格的適應性不佳[13]。

    對方程離散時,常采用向前差商、向后差商和中心差商等格式。

    向前差商:

    3.3 有限元法(FEM)

    有限元法基本思想是將一個求解區(qū)域劃分成許多微小單元,根據極值原理,把所有微小單元極值之和作為總體的極值。由于微小單元可以是各種形狀,因此網格劃分較靈活。這種方法守恒性好,但是在處理復雜方程時比較困難[13]。在推導有限元格式時一般會采用加權余量法,它的基本思想就是選定一個試探函數(shù)φi代替原方程中的待求函數(shù),此時必然會產生誤差,可以稱為余量R,然后在計算域Ω內找到n個線性無關的權函數(shù)δWi(i=1,2…,n),滿足

    這就意味著余量R在加權平均的意義上為零。

    3.4 有限體積法(FVM)

    有限體積法基本思路是構造一系列控制體,對控制體中的微分方程積分,控制體內又滿足總質量、總動量和總能量守恒,因此方程守恒性好,可以處理復雜網格計算問題,結合了FDM和FEM的優(yōu)點[12],也是目前使用最廣泛的方法,本文所討論的源項處理大部分也是基于此方法。

    在網格單元Vi內對式(1)積分:

    通過Gauss-Green公式將式(6)轉化為沿單元邊界的線積分,可以寫成

    式中:ΔVi為網格單元的面積;Vi為網格單元的邊界;n為單元邊線外法線的單位向量;l為單元邊長;S為網格單元源項的積分值。

    3.5 間斷有限元法(DG)

    間斷有限元法很好地結合了FDM、FEM和FVM的長處,具有精度高、守恒性好和易于處理復雜網格的特點[14-15],但是計算量大、程序設計較復雜且難以處理間斷問題[2],因此目前一般都在光滑處用此方法。在式(5)的基礎上,將試探函數(shù)φi作為權函數(shù)δWi:

    3.6 粒子法(SPH)

    粒子法是最近發(fā)展起來的一種無網格方法,它的核心是一種插值方法,可以處理復雜的外形,但是精度并不易提高[16]。

    4 源項的處理

    淺水方程求解過程中,源項處理是數(shù)值模擬成敗的關鍵,處理不當可能導致在巨大的計算量下結果依舊失真和失穩(wěn),同時破壞求解格式和諧性,因此學者們提出了廣泛的解決方法。源項可以分為底坡引起的底坡源項、底摩擦引起的摩阻源項、地球自轉引起的柯氏力和表面風應力。李文俊等[17]建立的二維水沙模型充分考慮了以上4個源項的影響,在模擬實際的大潮中,模擬結果與實測值能良好吻合。下面重點闡述對方程影響最大的兩個源項。

    4.1 底坡源項的處理

    Zhao等[18]、譚維炎等[19]在20世紀90年代把斜底地形簡化為多級階梯、每個階梯視為平底,不需要考慮底坡的影響,可以用齊次淺水方程計算。這樣做會產生兩類誤差:一類是采用平均底高和平均水力要素代替實際情況導致,可以通過減小階梯的長度來減少誤差;另一類是計算通量時相鄰平底模型中法向動量通量的對流項出現(xiàn)變化所致,可以通過加修正值保持原值不變。此方法更適合于地形變化較緩的區(qū)域。

    Bermudez等[20]在同時期提出用迎風格式處理含有源項的雙曲守恒方程,并對Van Leer的Q格式進行了擴展。與以前的方法相比,此方法得到的數(shù)值解在穩(wěn)定性和精度上有明顯的提升。應用到淺水方程后,該方法在處理源項時并不能始終保證數(shù)值格式擁有良好的性能。后來Bermudez將同樣的思想用于解決更貼近實際的二維淺水問題。Vazquez[21]、Brufau等[22]在此基礎上又改進了對地形等源項的模擬方法,使其能夠適應復雜多變的水流環(huán)境,但在求解法向通量時只驗證了Roe格式有較好的性能。

    Rogers等[23-24]提出的數(shù)值平衡法用確保質量和動量守恒的代數(shù)方式拆分底坡源項。這種方法應用到淺水方程后消除了近似解算器中的數(shù)值振蕩,但是由于沒有考慮靜水壓力項的非線性分布,因此底坡源項和通量項的誤差不可避免,只要水位不為零就會出現(xiàn)虛假流動。

    Leveque[25]提出了一種水波傳播算法,在每個計算單元內部人為地引入一個不連續(xù)量來處理源項。該方法適用于準穩(wěn)態(tài)條件,但在預測含有沖擊的跨臨界流動時無法帶來穩(wěn)定的數(shù)值結果。Kurganov等[26]將中心迎風格式(CU)推廣到淺水方程組,提出了一種適用于復雜地形的自適應算法。但是這些方法在處理循環(huán)流時會出現(xiàn)巨大的數(shù)值耗散。

    Zhou等[27]先推導了適用于均質無黏方程的深度梯度法(DGM),為解決這種方法受地形影響很大的問題,又通過水位重構的方法得到了水面梯度法(SGM),這種分段線性重建的方法在沒有底坡源項的時候與DGM相同,最后用HLL-Riemann解算器進行了驗證。水面梯度法(SGM)更適合在結構網格下的水流計算,原因是在非結構網格上Bermudez提出的C-property守恒概念并不成立。類似的,Wang等[28]采用通量分裂修正技術分解源項,以滿足C-property守恒,在光滑區(qū)域保持其原始高階精度,并在強間斷附近保持基本上無振蕩的性質。

    潘存鴻等提出了水位方程法(WLF),該方法將淺水流動方程中水深變量替換為水位,處理后用經典方法求Riemann解,這在本質上與Zhou等提出的SGM方法是相同的,只是WLF法更易理解,更能說明理論上SGM法的誤差性質和誤差大小的由來。其中底坡源項在Riemann問題中的離散形式必須與壓力梯度項的離散形式保持一致,且水深也要采用Riemann解,不然解的和諧性將無法保證。潘存鴻用WLF法建立Godunov格式求解一維淺水方程的Riemann問題,并推廣到求解二維淺水方程的Riemann問題[29],該方法保持了計算的通用性、和諧性和高分辨率。

    Mohammadian等[30]提出用修正水流對底坡源項進行相容離散化,由于在所提出的方法中,不需要對源項進行額外的迎風求解或Riemann求解,而是直接對源項離散,因此不光能采用Roe方法,還可以采用已有的激波捕捉方法,如CU方法、HLL方法、HLLC方法[31]等。與現(xiàn)有的許多格式不同,該方法易于在非結構網格上實現(xiàn),并能靈活處理不規(guī)則邊界和局部網格重構問題。作者通過一系列算例證明提出的方法能夠準確地模擬復雜地形下的各種臨界流和跨臨界流動。

    Smolarkiewicz等[32]采用MPDATA法(multidimensional positive definite advection transport algorithm)處理帶有底坡源項的方程。MPDATA法原本是一種處理大氣平流的方法,經過很多學者的不懈努力,其已經有能力解決復雜的流體問題。Jenny等[33]提出一種改進的解算器Rankine-Hugoniot-Riemannsolver(簡稱RHR解算器),它的基本思想是將源項視為不連續(xù)項,得到單元內的雙曲條件,就可以應用Riemann不變量的思想來特征處理源項。RHR解算器所得結果比傳統(tǒng)的Riemann解算器得到的結果更精確,在處理更高維度的情況時也具有較好表現(xiàn)。

    王志力等[34]提出了特征分解底坡源項來平衡界面通量的方法。相比于簡單處理源項,該方法可以保證格式的和諧性。王昆等[10]把底坡源項表示成某個矩陣的散度形式,計算量小,概念也簡單,為了保證格式的和諧性,還對兩個相鄰單元水深使用均方根形式以抵消數(shù)值通量項。周浩瀾等[1]用靜水重構法重構后的水深替換原水深,也保證了計算格式的和諧性。魏紅艷等[35]也采用了類似的方法。于守兵[36]提出的靜水壓力項和底坡項的積分平衡法,可以準確計算兩者的值并保持嚴格的平衡,不產生虛假流動。

    求解這種雙曲問題常用基于有限體積法的數(shù)值格式,其目的是將區(qū)域預先離散成體積單元并在這些單元內集成信息和控制方程來提供問題的數(shù)值解[37]。用有限體積法處理淺水方程,把方程中水位梯度項分解成靜水壓力梯度項和底坡源項,靜水壓力隨著水深的增大是非線性變化的,因此在靜水中底坡源項和靜水壓力梯度項往往難以相互抵消,這就是產生虛假流動的原因。為了避免出現(xiàn)這種情況,學術界主要存在兩種方法:一種是通過對模型的修正消去底坡源項,比如譚維炎等[19]提出的平底模型;另一種是對底坡源項進行分解或者改造,使它能夠與靜水壓力梯度項相互抵消,比如Mohammadian等[30]、潘存鴻等[38]、王志力等[34]提出的方法。

    網格的選擇也對源項處理有著直接影響,比如在非結構三角形網格上就比四邊形網格更難滿足和諧性要求。王昆等[10]、宋利祥等[39]在三角形網格下建立了和諧性離散格式,并用算例驗證了其適用性。

    使用有限體積法時,兩個單元的界面處為一個局部間斷,底坡源項此時就變成階梯間斷,有限體積法在這里需要求解階梯Riemann問題(Step Riemann Problem,SRP)。盡管這不是一個新問題,也有很多研究成果,但仍然是有爭議的,其研究的方法大致可分為兩類:一類是使用質量守恒和能量守恒的方法,Alcrudo等[40]和Bukreev等[41]認為必須采用這種方法,因為在SRP中間斷處的河床底坡趨向于無窮大,所以此時依據動量守恒推導出的SWE沒有意義,然而眾所周知的是在階梯處發(fā)生流動的過程中能量其實并不守恒,由于湍流的存在會導致能量耗散,因此這兩人的理論基礎存在瑕疵。Gallouet等[42]、Seguin等[43]和Andrianov[44]認為既然間斷處的接觸間斷是駐波,那么在特征空間里Riemann不變量守恒,這必然使得質量和能量守恒,然而這也是從數(shù)學觀點去分析這個問題,在物理上或者說在實際工況中并不是這樣。還有一類是使用質量守恒和動量守恒的方法,比如Bernetti等[45]在其研究成果中把能量守恒僅僅作為一個約束條件來檢驗所得的結果是否正確。Rosatti等[46]通過研究得出,當上游河床高于下游河床時,目前所得出的研究成果都是錯誤的,這個問題至今未得到解決。

    4.2 摩阻源項的處理

    對于不規(guī)則地形的淺水模擬,在局部陡峭的底坡可能會出現(xiàn)水深過小而流速很大的情況,這時會引起摩阻源項的剛性問題[47]。很多基于無結構網格的數(shù)值方法都采用顯式,顯式在程序的編寫上相對簡單,但是局限性也很大[48]。若采用顯式數(shù)值方法離散摩阻源項,可能會出現(xiàn)明顯的物理錯誤,比如出現(xiàn)虛假流動、負水深等[1],或是影響時間步長、降低計算效率[39]。隱式在計算穩(wěn)定性上要優(yōu)于顯式,且不受Courant的限制,可以與自適應網格技術結合。Michel等[49]在幾個潰壩數(shù)值模擬中對比了顯式處理和隱式處理的差異,總體來說,隱式處理更接近真實情況,比如在模擬干濕邊界時顯式處理易出現(xiàn)偽振蕩,而隱式處理并未出現(xiàn)。因此,一般采用隱式或者半隱式處理摩阻源項。

    王黨偉等[50]、王鑫等[51]用三階龍格-庫塔法處理源項,并對摩阻項采用自適應步長法進行處理。Yoon等[52]通過算子分裂法完全隱式處理摩阻項,以防止干燥區(qū)附近小水深引起的數(shù)值不穩(wěn)定性。Martinez[53]采用半隱式離散方法,Hou等[54]提出了一種基于隱式概念的方法,并且能夠集成到常用的顯式方法中,避免了冗余迭代,并能夠使用更大的計算時間步長,從而提高了編程和計算的效率。Cea等[55]、Singh等[56]將摩阻源項用速度的乘積表示以獲得顯式格式,雖然可以有效避免由剛性問題引起的數(shù)值不穩(wěn)定,但是流動通常無法維持穩(wěn)態(tài),導致模擬結果失真。Xia等[57]提出一種在有限體積法下的新的隱式解法,可以同時保證數(shù)值穩(wěn)定性和模擬的準確性,該解法可以不使用迭代法直接計算摩擦項,而且理論上還可以用于其他數(shù)值方法,如有限差分法和有限元法。

    5 總結與展望

    淺水方程的源項處理問題在經過國內外學者的潛心研究下取得了豐碩的成果,未來繼續(xù)改進和開發(fā)更簡單有效且具有更好收斂性的格式將會是一個趨勢。應該注意到有些源項問題至今沒有得到解決,比如當上游河床高于下游河床時,目前所有的局部階梯Riemann問題的求解方法都是錯誤的,需要通過數(shù)值試驗和物理試驗的方法深入研究才能找到準確的答案,這也會是將來淺水方程源項問題研究的方向。對于新的數(shù)值格式,比如大時間步長格式,傳統(tǒng)的源項處理方法不再適用,新格式的源項處理方法也是一個值得探索的研究方向。

    淺水方程作為非線性的偏微分方程組,現(xiàn)在并沒有解析解,所有解決的方法都是數(shù)值近似解,目前所有研究都是讓數(shù)值解盡可能地逼近實際情況,而任何方法都不是完美的,每一種數(shù)值方法的改進都會帶來一定的負效應,需要權衡這種改進與負效應在實際工程中的應用。因此,源項的數(shù)值處理方法還有很大的發(fā)展空間。

    參考文獻:

    [1] 周浩瀾,陳洋波,任啟偉.不規(guī)則地形淺水模擬[J].水動力學研究與進展(A輯),2010,25(5):594-600.

    [2] 潘存鴻.淺水間斷流動數(shù)值模擬研究進展[J].水利水電科技進展,2010,30(5):77-84.

    [3] SINGER T, VESTBERG M. Local Boundedness of Weak Solutions to the Diffusive Wave Approximation of the Shallow Water equations[J].Journal of Differential Equations,2019,266(6):3014-3033.

    [4] ZIJLEMA M. The Role of the Rankine-Hugoniot Relations in Staggered Finite Difference Schemes for the Shallow Water Equations[J].Computers and Fluids,2019,192(10):104274.

    [5] BUCKMASTER T, VICOL V. Nonuniqueness of Weak Solutions to the Navier-Stokes Equation[J].Annals of Mathematics,2019,189(1):101-144.

    [6] DU H, LIU Y, LIU Y, et al. Well-Balanced Discontinuous Galerkin Method for Shallow Water Equations with Constant Subtraction Techniques on Unstructured Meshes[J].Journal of Scientific Computing,2019,81(3):2115-2131.

    [7] 汪德爟.計算水力學理論與應用[M].北京:科學出版社,2011:31-36.

    [8] TORO E F. Shock-Capturing Methods for Free-Surface Shallow Flows[M]. New Jersey: John Wiley,2001:324-339.

    [9] ROE P L. Approximate Riemann Solvers, Parameter Vectors, and Difference Schemes[J].Journal of Computational Physics,1997,135(2):250-258.

    [10] 王昆,金生,馬志強,等.基于和諧性離散格式求解帶源項的二維淺水方程[J].水動力學研究與進展(A輯),2009,24(5):535-542.

    [11] MORALES M, LACASTA A, MURILLO J, et al. A Large Time Step Explicit Scheme (CFL>1) on Unstructured Grids for 2D Conservation Laws: Application to the Homogeneous Shallow Water Equations[J].Applied Mathematical Modelling,2017,47:294-317.

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

    [13] 劉林,常福宣,肖長偉,等.潰壩洪水研究進展[J].長江科學院院報,2016,33(6):29-35.

    [14] 韓濤,逄勇,李一平,等.基于間斷有限元法求解三維NS水流方程[J].人民黃河,2008,30(1):67-69.

    [15] CALEFFI V, VALIANI A, LI G. A Comparison Between Bottom-Discontinuity Numerical Treatments in the DG Framework[J].Applied Mathematical Modelling,2016,40(17-18):7516-7531.

    [16] 劉謀斌,宗智,常建忠.光滑粒子動力學方法的發(fā)展與應用[J].力學進展,2011,41(2):217-234.

    [17] 李文俊,張慶河,李龍翔,等.基于無積分節(jié)點間斷有限元的二維水沙模型:(1)水動力[J].水道港口,2019,40(2):125-134.

    [18] ZHAO D H, SHEN H W, III G Q T, et al. Finite-Volume Two-Dimensional Unsteady-Flow Model for River Basins[J].Journal of Hydraulic Engineering,1994,120(7):863-883.

    [19] 譚維炎,胡四一.淺水流動計算中一階有限體積法Osher格式的實現(xiàn)[J].水科學進展,1994,5(4):262-270.

    [20] BERMUDEZ A, VAZQUEZ M E. Upwind Methods for Hyperbolic Conservation Laws with Source Terms[J].Computers & Fluids,1994,23(8):1049-1071.

    [21] VAZQUEZ C M E. Improved Treatment of Source Terms in Upwind Schemes for the Shallow Water Equations in Channels with Irregular Geometry[J].Journal of computational Physics,1999,148(2):497-526.

    [22] BRUFAU P, GARCIANAVARRO P, VAZQUEZCENDON M E. Zero Mass Error Using Unsteady Wetting-Drying Conditions in Shallow Flows over Dry Irregular Topography[J].International Journal for Numerical Methods in Fluids,2004,45(10):1047-1082.

    [23] ROGERS B D, BORTHWICK A G L, TAYLOR P H. Mathematical Balancing of Flux Gradient and Source Terms Prior to Using Roes Approximate Riemann Solver[J].Journal of Computational Physics,2003,192(2):422-451.

    [24] ROGERS B D, FUJIHARA M, BORTHWICK A G L. Adaptive Q-Tree Godunov-Type Scheme for Shallow Water Equations[J].International Journal for Numerical Methods in Fluids,2001,35(3):247-280.

    [25] LEVEQUE R J. Balancing Source Terms and Flux Gradients in High-Resolution Godunov Methods[J].Journal of Computational Physics,1998,146(1):346-365.

    [26] KURGANOV A, LEVY D. Central-Upwind Schemes for the Saint-Venant System[J].Mathematical Modelling and Numerical Analysis,2002,36(3):397-425.

    [27] ZHOU J G, CAUSON D M, MINGHAM C G, et al. 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.

    [28] WANG Z, ZHU J, ZHAO N. A New Fifth-Order Finite Difference Well-Balanced Multi-Resolution WENO Scheme for Solving Shallow Water Equations[J].Computers and Mathematics with Applications,2020,80(5):1387-1404.

    [29] 潘存鴻,林炳堯,毛獻忠.求解二維淺水流動方程的Godunov格式[J].水動力學研究與進展(A輯),2003,18(1):16-23.

    [30] MOHAMMADIAN A, ROUX D Y L. Simulation of Shallow Flows over Variable Topographies Using Unstructured Grids[J].International Journal for Numerical Methods in Fluids,2006,52(5):473-498.

    [31] 楊峰,倪玉芳,黃衛(wèi),等.出山店水庫潰壩洪水數(shù)值模擬研究[J].人民黃河,2020,42(1):27-31,36.

    [32] SMOLARKIEWICZ P K, MARGOLIN L G. MPDATA: a Finite-Difference Solver for Geophysical Flows[J].Journal of Computational Physics,1998,140(2):459-480.

    [33] JENNY P, MULLER B. Rankine-Hugoniot-Riemann Solver Considering Source Terms and Multidimensional Effects[J].Journal of Computational Physics,1998,146(2):575-610.

    [34] 王志力,耿艷芬,金生.具有復雜計算域和地形的二維淺水流動數(shù)值模擬[J].水利學報,2005,36(4):439-444.

    [35] 魏紅艷,梁艷潔,陳萌,等.基于Roe格式的不規(guī)則地形上淺水模擬[J].武漢大學學報(工學版),2019,52(1):7-12,82.

    [36] 于守兵.計算二維淺水方程中靜水壓力項與底坡項的積分平衡法[J].水利水電科技進展,2009,29(4):32-35.

    [37] NAVAS A, MURILLO J. Overcoming Numerical Shockwave Anomalies Using Energy Balanced Numerical Schemes. Application to the Shallow Water Equations with Discontinuous Topography[J].Journal of Computational Physics, 2017,340:575-616.

    [38] 潘存鴻,林炳堯,毛獻忠.一維淺水流動方程的Godunov格式求解[J].水科學進展,2003,14(4):330-336.

    [39] 宋利祥,周建中,王光謙,等.潰壩水流數(shù)值計算的非結構有限體積模型[J].水科學進展,2011,22(3):373-381.

    [40] ALCRUDO F, BENKHALDOUN F. Exact Solutions to the Riemann Problem of the Shallow Water Equations with a Bottom Step[J].Computers and Fluids,2001,30(6):643-671.

    [41] BUKREEV V I, GUSEV A V, OSTAPENKO V V. Breakdown of a Discontinuity of the Free Fluid Surface over a Bottom Step in a Channel[J].Fluid Dynamics,2003,38(6):889-899.

    [42] GALLOUET T, HRARD J M, SEGUIN N. Some Approximate Godunov Schemes to Compute Shallow-Water Equations with Topography[J].Computers and Fluids,2003,32(4):479-513.

    [43] SEGUIN N, CHINNAYYA A, LEROUX A Y. A Well-Balanced Numerical Scheme for the Approximation of the Shallow-Water Equations with Topography: the Resonance Phenomenon[J].International Journal of Finite Volumes,2004,1(1):1-33.

    [44] ANDRIANOV N. Performance of Numerical Methods on the Non-Unique Solution to the Riemann Problem for the Shallow Water Equations[J].International Journal for Numerical Methods in Fluids,2005,47(8-9):825-831.

    [45] BERNETTI R, TITAREV V A, TORO E F. Exact Solution of the Riemann Problem for the Shallow Water Equations with Discontinuous Bottom Geometry[J].Journal of Computational Physics,2008,227(6):3212-3243.

    [46] ROSATTI G, BEGNUDELLI L. The Riemann Problem for the One-Dimensional, Free-Surface Shallow Water Equations with a Bed Step: Theoretical Analysis and Numerical Simulations[J].Journal of Computational Physics,2009,229(3):760-787.

    [47] XIA X, LIANG Q, MING X, et al. An Efficient and Stable Hydrodynamic Model with Novel Source Term Discretization Schemes for Overland Flow and Flood Simulations[J].Water Resources Research,2017,53(5):3730-3759.

    [48] 唐岳灝.基于無結構化網格淺水方程的隱式解法[J].人民長江,2015,46(5):81-84,96.

    [49] MICHEL D V, BERTHON C, CLAIN S, et al. A Well-Balanced Scheme for the Shallow-Water Equations with Topography or Manning Friction[J].Journal of Computational Physics,2017,335:115-154.

    [50] 王黨偉,陳建國,吉祖穩(wěn).不規(guī)則地形上淺水模擬平衡性的實現(xiàn)[J].計算力學學報,2012,29(4):604-608,615.

    [51] 王鑫,曹志先,岳志遠.強不規(guī)則地形上淺水二維流動的數(shù)值計算研究[J].水動力學研究與進展(A輯),2009,24(1):56-62.

    [52] YOON T H, KANG S. Finite Volume Model for Two-Dimensional Shallow Water Flows on Unstructured Grids[J].Journal of Hydraulic Engineering,2004,130(7):678-688.

    [53] MARTINEZ V B. A Numerical Technique for Applying Time Splitting Methods in Shallow Water Equations[J].Computers & Fluids,2017,169:285-295.

    [54] HOU J, WANG T, LI P, et al. An Implicit Friction Source Term Treatment for Overland Flow Simulation Using Shallow Water Flow Model[J].Journal of Hydrology,2018,564:357-366.

    [55] CEA L, BLADE E. A Simple and Efficient Unstructured Finite Volume Scheme for Solving the Shallow Water Equations in Overland Flow Applications[J].Water Resources Research,2015,51(7):5464-5486.

    [56] SINGH J, ALTINAKAR M S, DING Y. Numerical Modeling of Rainfall-Generated Overland Flow Using Nonlinear Shallow-Water Equations[J].Journal of Hydrologic Engineering,2014,20(8):172-183.

    [57] XIA X, LIANG Q. A New Efficient Implicit Scheme for Discretising the Stiff Friction Terms in the Shallow Water Equations[J].Advances in Water Resources,2018,117:87-97.

    【責任編輯 張 帥】

    又紧又爽又黄一区二区| 成人午夜高清在线视频| 丁香欧美五月| 动漫黄色视频在线观看| 国产伦在线观看视频一区| 亚洲精品乱码久久久v下载方式 | 狂野欧美激情性xxxx| 精品熟女少妇八av免费久了| 成人性生交大片免费视频hd| 日韩亚洲欧美综合| 国产精品久久视频播放| 久久久精品欧美日韩精品| 舔av片在线| 69av精品久久久久久| 欧美性感艳星| 欧美性猛交黑人性爽| 亚洲国产欧美人成| 亚洲精华国产精华精| 操出白浆在线播放| 国产单亲对白刺激| 国产精品精品国产色婷婷| 一级毛片女人18水好多| 性欧美人与动物交配| 国产av在哪里看| 国产一级毛片七仙女欲春2| 两性午夜刺激爽爽歪歪视频在线观看| 在线看三级毛片| 高清日韩中文字幕在线| 亚洲精品日韩av片在线观看 | 全区人妻精品视频| 国内精品美女久久久久久| 精品一区二区三区视频在线观看免费| 一个人看视频在线观看www免费 | 国模一区二区三区四区视频| 国产午夜福利久久久久久| 99久久九九国产精品国产免费| 午夜两性在线视频| 亚洲一区高清亚洲精品| 91在线精品国自产拍蜜月 | 最近最新免费中文字幕在线| 综合色av麻豆| 精品久久久久久久毛片微露脸| 国产国拍精品亚洲av在线观看 | 香蕉丝袜av| 一区福利在线观看| 无遮挡黄片免费观看| 成人国产一区最新在线观看| 日韩人妻高清精品专区| 中文字幕av成人在线电影| 久久精品人妻少妇| 国产三级中文精品| 一个人看的www免费观看视频| 午夜日韩欧美国产| 真人做人爱边吃奶动态| 看黄色毛片网站| 欧美日本视频| 亚洲一区二区三区不卡视频| 婷婷亚洲欧美| 国产精品98久久久久久宅男小说| 一级毛片女人18水好多| 老司机福利观看| 婷婷丁香在线五月| 午夜福利在线观看吧| 天天添夜夜摸| 五月伊人婷婷丁香| 国产精品亚洲av一区麻豆| 国产成人欧美在线观看| 最近在线观看免费完整版| 亚洲精品在线观看二区| 国产精品99久久久久久久久| 少妇熟女aⅴ在线视频| 麻豆成人av在线观看| 神马国产精品三级电影在线观看| 亚洲自拍偷在线| 国产蜜桃级精品一区二区三区| 亚洲18禁久久av| 欧美日韩中文字幕国产精品一区二区三区| 观看美女的网站| 丝袜美腿在线中文| 欧美成狂野欧美在线观看| 一级黄片播放器| 九九热线精品视视频播放| 国产精品一区二区三区四区久久| 最近最新免费中文字幕在线| 国产色婷婷99| 无遮挡黄片免费观看| 狠狠狠狠99中文字幕| 欧美一级毛片孕妇| 激情在线观看视频在线高清| 日韩欧美三级三区| 热99re8久久精品国产| 色综合婷婷激情| 亚洲精品在线观看二区| 亚洲国产精品成人综合色| 熟女少妇亚洲综合色aaa.| 成年版毛片免费区| 国产成年人精品一区二区| 久9热在线精品视频| 免费观看精品视频网站| 亚洲精品亚洲一区二区| 69人妻影院| 欧美一级毛片孕妇| 成人特级黄色片久久久久久久| 久久久成人免费电影| 日本一二三区视频观看| netflix在线观看网站| 久久国产乱子伦精品免费另类| 在线观看美女被高潮喷水网站 | 嫁个100分男人电影在线观看| 夜夜躁狠狠躁天天躁| 国产免费一级a男人的天堂| 19禁男女啪啪无遮挡网站| 午夜激情欧美在线| 日韩精品中文字幕看吧| 99热只有精品国产| e午夜精品久久久久久久| 欧美日韩国产亚洲二区| 黄色女人牲交| 日韩亚洲欧美综合| 两性午夜刺激爽爽歪歪视频在线观看| 欧美区成人在线视频| 亚洲久久久久久中文字幕| 色综合站精品国产| 亚洲性夜色夜夜综合| 久久性视频一级片| 美女高潮的动态| 看片在线看免费视频| 欧美日韩黄片免| 亚洲激情在线av| av视频在线观看入口| 亚洲国产色片| 国产精品久久久久久久久免 | 国产精品久久久久久亚洲av鲁大| 久久精品国产99精品国产亚洲性色| 老司机午夜十八禁免费视频| 成人av在线播放网站| www.www免费av| 国产一区二区亚洲精品在线观看| 淫秽高清视频在线观看| 亚洲五月天丁香| 51国产日韩欧美| 亚洲第一电影网av| 波多野结衣巨乳人妻| 99国产精品一区二区三区| 精品一区二区三区视频在线观看免费| 一级a爱片免费观看的视频| 99精品欧美一区二区三区四区| 午夜亚洲福利在线播放| 真人一进一出gif抽搐免费| 亚洲avbb在线观看| 好看av亚洲va欧美ⅴa在| 波多野结衣高清无吗| 女生性感内裤真人,穿戴方法视频| 精品无人区乱码1区二区| 欧美日韩瑟瑟在线播放| 亚洲av成人av| 欧美另类亚洲清纯唯美| 日韩欧美 国产精品| 舔av片在线| 亚洲国产中文字幕在线视频| 精品电影一区二区在线| 一区二区三区免费毛片| 怎么达到女性高潮| 一a级毛片在线观看| 毛片女人毛片| 欧美激情久久久久久爽电影| 哪里可以看免费的av片| 内地一区二区视频在线| 不卡一级毛片| 亚洲av熟女| 两个人看的免费小视频| 精品久久久久久,| 欧美成人a在线观看| 精品国内亚洲2022精品成人| 国产精品av视频在线免费观看| 给我免费播放毛片高清在线观看| 搞女人的毛片| 久久精品影院6| 桃色一区二区三区在线观看| 欧美三级亚洲精品| 最近最新中文字幕大全免费视频| 国产激情欧美一区二区| 丰满的人妻完整版| 女人被狂操c到高潮| 在线天堂最新版资源| 琪琪午夜伦伦电影理论片6080| 亚洲精品456在线播放app | 久久久精品大字幕| 99国产极品粉嫩在线观看| 欧美成人a在线观看| 可以在线观看毛片的网站| 男女午夜视频在线观看| 一区二区三区免费毛片| 欧美中文综合在线视频| 午夜精品一区二区三区免费看| 夜夜躁狠狠躁天天躁| 在线观看舔阴道视频| 亚洲va日本ⅴa欧美va伊人久久| 老司机在亚洲福利影院| 国产高潮美女av| 男女做爰动态图高潮gif福利片| 精品久久久久久久毛片微露脸| 在线观看日韩欧美| 日本一本二区三区精品| 国产精品国产高清国产av| 人人妻人人澡欧美一区二区| 此物有八面人人有两片| 黄色女人牲交| 欧美成人a在线观看| 法律面前人人平等表现在哪些方面| 亚洲va日本ⅴa欧美va伊人久久| 亚洲欧美日韩无卡精品| 高清毛片免费观看视频网站| 国产在视频线在精品| 波野结衣二区三区在线 | 欧美日本亚洲视频在线播放| 搡老岳熟女国产| 色综合站精品国产| 国产精华一区二区三区| 欧美成人免费av一区二区三区| 日日摸夜夜添夜夜添小说| 国产一区二区三区视频了| 国产国拍精品亚洲av在线观看 | 国产综合懂色| 国产亚洲精品一区二区www| 99久久久亚洲精品蜜臀av| 国产私拍福利视频在线观看| 在线观看午夜福利视频| 国产在视频线在精品| 在线观看免费视频日本深夜| av欧美777| 国产精品久久久久久人妻精品电影| 亚洲精品456在线播放app | 国产一区二区亚洲精品在线观看| 国产亚洲精品av在线| 一个人看的www免费观看视频| 日本 av在线| 国模一区二区三区四区视频| 日韩欧美国产在线观看| 白带黄色成豆腐渣| 少妇熟女aⅴ在线视频| 法律面前人人平等表现在哪些方面| av专区在线播放| 精品国产超薄肉色丝袜足j| 2021天堂中文幕一二区在线观| 国产伦精品一区二区三区视频9 | 精品欧美国产一区二区三| 淫妇啪啪啪对白视频| 我的老师免费观看完整版| 欧美一区二区国产精品久久精品| 亚洲精品久久国产高清桃花| 国内揄拍国产精品人妻在线| 久久久精品欧美日韩精品| 亚洲精品456在线播放app | 丝袜美腿在线中文| 成年免费大片在线观看| 搡老熟女国产l中国老女人| 国产成人系列免费观看| 久久6这里有精品| 日韩欧美在线二视频| 日本黄色视频三级网站网址| 成人av在线播放网站| 啦啦啦免费观看视频1| 91久久精品国产一区二区成人 | 欧美日韩乱码在线| 两人在一起打扑克的视频| 亚洲欧美日韩东京热| 亚洲av免费在线观看| 可以在线观看的亚洲视频| 午夜免费男女啪啪视频观看 | 国产午夜精品久久久久久一区二区三区 | 久久久色成人| 国产精品久久久久久久久免 | 国产色婷婷99| 亚洲精品在线美女| 国产一区二区在线av高清观看| 一级黄片播放器| 国产伦一二天堂av在线观看| 久久99热这里只有精品18| 午夜老司机福利剧场| 99久久综合精品五月天人人| 99久久九九国产精品国产免费| 18禁国产床啪视频网站| 国产精品嫩草影院av在线观看 | 免费看美女性在线毛片视频| 午夜精品久久久久久毛片777| 少妇裸体淫交视频免费看高清| 两个人的视频大全免费| 又粗又爽又猛毛片免费看| 久久香蕉国产精品| 日韩亚洲欧美综合| 国内久久婷婷六月综合欲色啪| 亚洲国产色片| 欧美性感艳星| 国产成人aa在线观看| 99国产精品一区二区三区| av欧美777| 午夜福利在线观看吧| 男女那种视频在线观看| 2021天堂中文幕一二区在线观| 国产成人影院久久av| 国产亚洲精品综合一区在线观看| 亚洲精品一区av在线观看| 色在线成人网| 精品国产三级普通话版| 日本免费a在线| 可以在线观看的亚洲视频| 哪里可以看免费的av片| 好男人在线观看高清免费视频| 成人国产综合亚洲| 99国产精品一区二区蜜桃av| 日韩国内少妇激情av| 最近视频中文字幕2019在线8| 亚洲在线自拍视频| 欧美午夜高清在线| 中文在线观看免费www的网站| 国产一区在线观看成人免费| 少妇的丰满在线观看| 日韩精品青青久久久久久| 欧美成狂野欧美在线观看| 在线a可以看的网站| 日韩av在线大香蕉| 欧美精品啪啪一区二区三区| 亚洲欧美日韩卡通动漫| 激情在线观看视频在线高清| 国产精品久久久久久人妻精品电影| 特级一级黄色大片| 深爱激情五月婷婷| 亚洲国产日韩欧美精品在线观看 | 精品欧美国产一区二区三| 琪琪午夜伦伦电影理论片6080| 男女做爰动态图高潮gif福利片| 亚洲欧美一区二区三区黑人| av女优亚洲男人天堂| 日韩欧美在线乱码| 一区二区三区高清视频在线| 床上黄色一级片| 色av中文字幕| 久久精品亚洲精品国产色婷小说| 午夜免费成人在线视频| 午夜福利免费观看在线| 中出人妻视频一区二区| 欧美日韩黄片免| 一夜夜www| 中文字幕人妻熟人妻熟丝袜美 | 最近最新免费中文字幕在线| 怎么达到女性高潮| 久久这里只有精品中国| 特大巨黑吊av在线直播| 啦啦啦免费观看视频1| 国产精品亚洲av一区麻豆| 女人高潮潮喷娇喘18禁视频| 又爽又黄无遮挡网站| 亚洲自拍偷在线| 日本黄色片子视频| 美女免费视频网站| 亚洲国产中文字幕在线视频| 国产精品久久久人人做人人爽| 亚洲久久久久久中文字幕| 少妇丰满av| 一卡2卡三卡四卡精品乱码亚洲| 国产三级在线视频| 久久草成人影院| 色在线成人网| 色av中文字幕| 久久久成人免费电影| 国内久久婷婷六月综合欲色啪| 婷婷六月久久综合丁香| 国内久久婷婷六月综合欲色啪| 婷婷六月久久综合丁香| 国产又黄又爽又无遮挡在线| 亚洲国产欧洲综合997久久,| 国产欧美日韩精品一区二区| 一个人免费在线观看的高清视频| 欧美+日韩+精品| 18禁在线播放成人免费| www日本在线高清视频| av视频在线观看入口| 亚洲国产精品合色在线| 久久精品91无色码中文字幕| 久久精品91蜜桃| 给我免费播放毛片高清在线观看| 国产精品综合久久久久久久免费| 精品无人区乱码1区二区| 99国产极品粉嫩在线观看| aaaaa片日本免费| 女人高潮潮喷娇喘18禁视频| 精品午夜福利视频在线观看一区| 草草在线视频免费看| 在线播放国产精品三级| 午夜免费成人在线视频| 日韩欧美精品免费久久 | 熟女人妻精品中文字幕| 亚洲不卡免费看| 国产精品电影一区二区三区| 午夜福利在线在线| 国产三级黄色录像| 中文字幕精品亚洲无线码一区| 麻豆成人av在线观看| 国产精品女同一区二区软件 | 免费看美女性在线毛片视频| 欧美一级a爱片免费观看看| 国产精品一及| 国产色爽女视频免费观看| 高清日韩中文字幕在线| 成人18禁在线播放| 搡老岳熟女国产| 国产日本99.免费观看| 亚洲 欧美 日韩 在线 免费| 免费看美女性在线毛片视频| 亚洲第一欧美日韩一区二区三区| 亚洲av五月六月丁香网| 亚洲国产欧洲综合997久久,| 欧美色视频一区免费| 一区二区三区国产精品乱码| 精品国产三级普通话版| 欧美成人性av电影在线观看| 国产高清videossex| 久久精品人妻少妇| 高清毛片免费观看视频网站| 亚洲欧美日韩东京热| 成人特级av手机在线观看| 国产高清三级在线| 久久亚洲精品不卡| 欧美成人免费av一区二区三区| 婷婷精品国产亚洲av在线| 麻豆成人午夜福利视频| 欧美激情久久久久久爽电影| 男女床上黄色一级片免费看| 欧美色欧美亚洲另类二区| 精品一区二区三区视频在线观看免费| 最近最新中文字幕大全免费视频| 久久久久久久久中文| 国产主播在线观看一区二区| 精品欧美国产一区二区三| 在线十欧美十亚洲十日本专区| 亚洲无线在线观看| 亚洲在线观看片| 亚洲中文日韩欧美视频| 综合色av麻豆| 又黄又爽又免费观看的视频| 色精品久久人妻99蜜桃| 我的老师免费观看完整版| 少妇人妻精品综合一区二区 | 久久久久国内视频| 久久久久久久亚洲中文字幕 | 亚洲狠狠婷婷综合久久图片| 国产激情偷乱视频一区二区| 十八禁网站免费在线| 国产精品久久电影中文字幕| 日韩欧美 国产精品| 久久久久国内视频| 国产亚洲精品一区二区www| 亚洲狠狠婷婷综合久久图片| 中文字幕高清在线视频| 亚洲精品粉嫩美女一区| 久久久精品大字幕| 亚洲激情在线av| 丰满人妻熟妇乱又伦精品不卡| 国产精品av视频在线免费观看| 亚洲国产精品久久男人天堂| 久久这里只有精品中国| 熟女少妇亚洲综合色aaa.| 国产伦精品一区二区三区视频9 | 日日夜夜操网爽| 99精品欧美一区二区三区四区| 丰满人妻熟妇乱又伦精品不卡| 亚洲欧美日韩无卡精品| 老司机午夜十八禁免费视频| 女人高潮潮喷娇喘18禁视频| av片东京热男人的天堂| 97碰自拍视频| 老司机午夜十八禁免费视频| 午夜免费成人在线视频| 国产97色在线日韩免费| 精品人妻1区二区| x7x7x7水蜜桃| 麻豆久久精品国产亚洲av| 18美女黄网站色大片免费观看| 欧美乱色亚洲激情| 少妇裸体淫交视频免费看高清| 久久精品国产综合久久久| 久久精品国产清高在天天线| 成年免费大片在线观看| 三级国产精品欧美在线观看| 嫩草影院入口| 在线观看舔阴道视频| 色综合站精品国产| 国内精品久久久久精免费| 他把我摸到了高潮在线观看| 99久国产av精品| 国产三级在线视频| 久99久视频精品免费| 国产精品1区2区在线观看.| 观看美女的网站| 一区二区三区国产精品乱码| 国产精品久久久久久精品电影| 国产视频内射| 国产三级中文精品| 狂野欧美白嫩少妇大欣赏| 亚洲专区中文字幕在线| 久99久视频精品免费| 老熟妇仑乱视频hdxx| 亚洲欧美日韩高清专用| 老汉色∧v一级毛片| 真人一进一出gif抽搐免费| 91久久精品电影网| 91麻豆av在线| 国产成年人精品一区二区| 亚洲无线观看免费| 中国美女看黄片| 久久精品人妻少妇| 亚洲精品乱码久久久v下载方式 | 久久久久免费精品人妻一区二区| 日韩高清综合在线| 亚洲精品一区av在线观看| 国产精品永久免费网站| 国产一区二区激情短视频| 日韩精品中文字幕看吧| 一级作爱视频免费观看| 色噜噜av男人的天堂激情| 精品一区二区三区av网在线观看| 午夜亚洲福利在线播放| 欧美日韩黄片免| www国产在线视频色| 久久久久久久久久黄片| 国产伦一二天堂av在线观看| 国产精华一区二区三区| 日韩欧美精品v在线| 黑人欧美特级aaaaaa片| 亚洲国产欧美人成| 亚洲成av人片在线播放无| 偷拍熟女少妇极品色| 国产午夜精品论理片| 国产一区二区三区在线臀色熟女| 午夜视频国产福利| 中文字幕人妻丝袜一区二区| 中文在线观看免费www的网站| 久久99热这里只有精品18| 高清在线国产一区| 日韩欧美三级三区| 国产色婷婷99| 亚洲国产欧美人成| 免费av观看视频| 亚洲av一区综合| 免费看光身美女| 狠狠狠狠99中文字幕| 国产精品99久久99久久久不卡| bbb黄色大片| 美女免费视频网站| 日日夜夜操网爽| 一个人观看的视频www高清免费观看| 我的老师免费观看完整版| 三级毛片av免费| 欧美av亚洲av综合av国产av| 特大巨黑吊av在线直播| 亚洲国产精品sss在线观看| 少妇的丰满在线观看| 国产精品香港三级国产av潘金莲| 精品久久久久久久久久久久久| 亚洲国产欧洲综合997久久,| 中文字幕高清在线视频| 国产伦人伦偷精品视频| 欧美最新免费一区二区三区 | 免费观看人在逋| 婷婷精品国产亚洲av| 一本久久中文字幕| 99久久久亚洲精品蜜臀av| 国内精品久久久久久久电影| 久久人人精品亚洲av| 亚洲人成网站在线播| 色噜噜av男人的天堂激情| 最新在线观看一区二区三区| 亚洲精品亚洲一区二区| 日本撒尿小便嘘嘘汇集6| 91av网一区二区| 久久伊人香网站| 在线免费观看的www视频| 亚洲不卡免费看| 成人三级黄色视频| www国产在线视频色| 午夜a级毛片| 婷婷丁香在线五月| 午夜福利在线在线| 亚洲人与动物交配视频| 美女大奶头视频| 国产亚洲欧美98| 18禁黄网站禁片免费观看直播| 欧美成人免费av一区二区三区| avwww免费| 成人18禁在线播放| 高清日韩中文字幕在线| 久久久久久久亚洲中文字幕 | 国产精品自产拍在线观看55亚洲| 香蕉久久夜色| 欧美在线一区亚洲| 在线观看日韩欧美| 又紧又爽又黄一区二区| 亚洲最大成人手机在线| 日韩人妻高清精品专区| 女人被狂操c到高潮| 婷婷六月久久综合丁香| 国产成人a区在线观看| 琪琪午夜伦伦电影理论片6080| 日韩国内少妇激情av| 毛片女人毛片| 午夜影院日韩av| 国产av一区在线观看免费| 亚洲国产欧美人成| 欧美不卡视频在线免费观看| 日本一二三区视频观看| 最好的美女福利视频网| 一本精品99久久精品77| 中文字幕熟女人妻在线| 在线观看免费视频日本深夜|