999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

分布式水文模型的構(gòu)建與應(yīng)用比較研究

2020-10-20 06:37:11菅浩然劉洪波童冰星
人民黃河 2020年5期

菅浩然 劉洪波 童冰星

摘 要:為了給分布式水文模型在陜西省中小流域的進(jìn)一步推廣與應(yīng)用提供參考,在板橋流域和馬渡王流域分別應(yīng)用TOPMODEL、TOPKAPI以及CASC2D等分布式水文模型對2000—2009年間的洪水進(jìn)行模擬,從兩流域的降雨徑流關(guān)系、地形指數(shù)特征以及土壤植被分布等方面進(jìn)行了分析。結(jié)果表明:在常發(fā)生超滲產(chǎn)流的板橋流域CASC2D模型的模擬效果較好,而在蓄滿產(chǎn)流發(fā)生概率較高的馬渡王流域采用蓄滿產(chǎn)流模式的TOPKAPI模型表現(xiàn)更優(yōu)。

關(guān)鍵詞:TOPMODEL模型;TOPKAPI模型;CASC2D模型;板橋流域;馬渡王流域

中圖分類號:P333 文獻(xiàn)標(biāo)志碼:A

doi:10.3969/j.issn.1000-1379.2020.05.005

Abstract: In this paper, TOPMODEL, TOPKAPI and CASC2D model were applied to simulate floods in Banqiao and Maduwang watersheds from 2000 to 2009. The results show that the CASC2D model has better effect in Banqiao watershed. While, the TOPKAPI model utilized a saturated runoff mechanism is more suitable for Maduwang watershed. In view of this phenomenon, this paper analyzed the relationship between rainfall and runoff, the spatial characteristics of topographic index and the distribution of soil vegetation in study watersheds. Through the analysis of the simulation results, the applicability of the hydrological model in these two basins was not only discussed, but also the hydrological runoff laws of these two basins were compared and analyzed. To a certain extent, the influence of comprehensive natural geographic conditions to the calculation results of distributed hydrological models was explained, which was beneficial to the further application of distributed hydrological models in medium and small river basins in Shaanxi Province.

Key words: TOPMODEL; TOPKAPI; CASC2D; Banqiao Watershed; Maduwang watershed

1 引 言

20世紀(jì)60年代,流域下墊面地理特征的空間分異性對降雨-徑流過程的影響逐漸得到水文工作者的重視[1-4]。由于概念性集總式模型難以對氣候時空變化和下墊面地理特性的空間分異性做出及時有效的反映,從而導(dǎo)致水文過程模擬精度降低,因此在實(shí)際的生產(chǎn)、研究過程中受到一定的限制。隨著計(jì)算機(jī)技術(shù)的蓬勃發(fā)展以及遙感信息技術(shù)的突飛猛進(jìn),越來越多的分布式水文模型相繼被提出[3-7],由于其能夠基于分布式的輸入數(shù)據(jù)以及邊界條件對流域內(nèi)部的水文過程進(jìn)行一定的描述,因此逐漸成為21世紀(jì)水文水資源領(lǐng)域的研究熱點(diǎn)之一[8-11]。目前已得到廣泛應(yīng)用的分布式水文模型包括TOPMODEL模型[12-14]、TOPKAPI模型[15-17]和CASC2D模型[18-20]等。

由Beven和Kirkby在1979年提出的TOPMODEL模型是基于地形的半分布式流域水文模型,該模型主要以地形指數(shù)的空間變化來模擬徑流產(chǎn)生的變動產(chǎn)流面積;2005年,張珂、李致家等提出了基于柵格數(shù)字高程和地形的GTOPMODEL模型[21-23];2002年,劉志雨等提出了TOPKAPI模型,基于運(yùn)動波原理對柵格單元內(nèi)地表水和壤中水的運(yùn)動進(jìn)行模擬計(jì)算;1991年,科羅拉多州立大學(xué)的Julien教授基于對二維地面徑流計(jì)算方法的研究提出了CASC2D模型[24-27]。

陜西省南部的板橋流域和馬渡王流域均屬于濕潤半濕潤地區(qū),由于兩流域地形起伏較大,土壤、植被等下墊面特征空間分異明顯,使用傳統(tǒng)的概念性集總式模型難以準(zhǔn)確地反映流域的降雨-徑流過程,因此筆者對兩流域分別構(gòu)建TOPMODEL模型、TOPKAPI模型和CASC2D模型,并應(yīng)用這些模型對兩流域的若干次洪水過程進(jìn)行模擬。通過對模擬結(jié)果的分析,不僅探討了各水文模型在兩流域的適用性,同時對兩流域的水文徑流規(guī)律進(jìn)行了比較分析,在一定程度上闡明了流域的綜合自然地理?xiàng)l件對分布式水文模型計(jì)算結(jié)果的影響,以期為分布式水文模型在陜西省中小流域的進(jìn)一步推廣應(yīng)用提供參考。

2 流域概況

板橋流域位于陜西省南部,流域面積約493 km2,流域內(nèi)地質(zhì)構(gòu)造以片麻巖、石灰?guī)r、碎屑巖為主。高山以森林為主,低山為荒山,耕地分布在河谷,地形為西北高東南低,流域的數(shù)字高程模型如圖1(a)所示。板橋流域多年平均降水量為729.0 mm,汛期(5—10月)降雨量占年降水量的80%~90%。流域中水系呈枝杈狀,主要的河流板橋河位于丹江左岸,屬丹江一級支流,發(fā)源于秦嶺主脊南側(cè),河長50.7 km,平均比降1.72%,流域入口為二龍山水庫,多年平均徑流量為8.752×107 m3。由于山體坡度大,因此暴雨發(fā)生時匯流快,很快形成洪峰,洪水暴漲暴落。

馬渡王流域位于陜西省南部,流域面積為2 577 km2,流域的數(shù)字高程模型如圖1(b)所示。該流域中的灞河發(fā)源于藍(lán)田縣灞源鄉(xiāng)麻家村以上的秦嶺北坡,灞河為不對稱水系,左岸支流少而長,右岸支流多而短,該河具有山溪性河流的特點(diǎn),河水陡漲陡落,多年平均降水量為630.9 mm,洪水均由夏季暴雨和初秋連陰雨造成,最早發(fā)生在4月上旬,最晚發(fā)生在10月下旬,7、8月洪水較大,洪峰漲落急劇。

3 模型構(gòu)建

結(jié)合研究流域的實(shí)際自然地理?xiàng)l件,分別對板橋流域和馬渡王流域構(gòu)建TOPMODEL模型、TOPKAPI模型以及CASC2D模型。

3.1 TOPMODEL模型

基于美國地質(zhì)調(diào)查局免費(fèi)提供的全球90 m×90 m的原始DEM數(shù)據(jù),將板橋流域和馬渡王流域分別劃分為若干個正交的網(wǎng)格單元。對于每一個網(wǎng)格單元,其水分運(yùn)動規(guī)律如圖2所示(P為降水量;E為蒸發(fā)量;S(i)為土壤飽和蓄水深;Srz為植被根系區(qū)的蓄水深;Srmax為植被根系區(qū)的最大蓄水深;Suz為土壤非飽和區(qū)蓄水深;qs為飽和坡面流量;qb為壤中流量;qv為下滲量)。

在進(jìn)行板橋流域和馬渡王流域的產(chǎn)流計(jì)算時,TOPMODEL模型采用流域變動產(chǎn)流理論[28-29],在整個計(jì)算過程中,源面積是不斷變化的,亦稱變動產(chǎn)流面積理論。研究流域中源面積的位置受流域地形和土壤水力特性兩個因素的影響,在一定意義上,變動產(chǎn)流面積可看作河道系統(tǒng)的延伸,如圖3所示。

TOPMODEL主要通過流域含水量(或缺水量)來確定源面積的大小和位置。而含水量的大小可由地形指數(shù)來計(jì)算。土壤相對含水量與地形指數(shù)呈正相關(guān)關(guān)系,土壤相對含水量和地形指數(shù)隨匯流面積的擴(kuò)大及坡度的減小而增加。通過計(jì)算得到的板橋和馬渡王流域的地形指數(shù)如圖4所示,以此為基礎(chǔ)描述兩流域中的徑流在重力排水作用下沿坡向的運(yùn)動。

3.2 TOPKAPI模型

TOPKAPI模型[30-33]將板橋流域和馬渡王流域的流域特性參數(shù)、降水和水文響應(yīng)的空間分布,在水平方向上用正交網(wǎng)格系統(tǒng),在垂直方向上用各網(wǎng)格所對應(yīng)的水平土柱進(jìn)行模擬。植被和土壤數(shù)據(jù)采用美國馬里蘭大學(xué)發(fā)布的全球土壤類型與土壤利用類型數(shù)據(jù)。其主要的計(jì)算模塊(如蒸散發(fā)、下滲、壤中流、地表徑流、地下徑流以及河道水流等)均采用相似的非線性水庫方程來進(jìn)行計(jì)算。非線性水庫方程的主要參數(shù)可以依據(jù)板橋流域和馬渡王流域的土壤空間分布估計(jì)得到,見表1。

3.3 CASC2D模型

在板橋流域和馬渡王流域,為了使CASC2D模型[34-35]計(jì)算更為快捷,通過重采樣等一系列操作,將90 m×90 m分辨率的DEM數(shù)據(jù)轉(zhuǎn)換為400 m分辨率的數(shù)據(jù),并基于美國馬里蘭大學(xué)的全球1 km精度的植被覆蓋數(shù)據(jù)描述兩流域內(nèi)植被的空間分布,如圖6所示,以0.16 km2為河道閾值提取流域中的水系。結(jié)合板橋流域和馬渡王流域內(nèi)實(shí)際狀況與匯流特征,設(shè)定模型的模擬計(jì)算步長為6 s,降雨輸入時段長為1 h。

結(jié)合CASC2D模型的使用經(jīng)驗(yàn)[36-37],參數(shù)率定過程遵循先調(diào)整水量再調(diào)整過程、先調(diào)整峰值再調(diào)整峰現(xiàn)時差的原則,然后結(jié)合實(shí)測徑流過程資料采用人工試錯法對參數(shù)進(jìn)行調(diào)整。CASC2D模型在使用過程中還需要根據(jù)板橋流域和馬渡王流域的DEM、土地利用類型、土壤類型以及河道特征來估算模型參數(shù),具體見表3和表4。

4 計(jì)算結(jié)果與分析

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

通過對兩流域數(shù)據(jù)資料較為完備的2000—2009年的洪水進(jìn)行分析模擬,分別對模擬結(jié)果的確定性系數(shù)、洪峰相對誤差和峰現(xiàn)時間誤差進(jìn)行了統(tǒng)計(jì)分析。具體結(jié)果見表5~表10。

4.1.1 板橋流域計(jì)算結(jié)果

基于表5~表7中的統(tǒng)計(jì)結(jié)果可知,在板橋流域中,TOPMODEL模型、TOPKAPI模型和CASC2D模型的平均確定性系數(shù)為-1.14、-0.42和0.58。結(jié)合洪峰相對誤差和峰現(xiàn)時間誤差可以看出,CASC2D模型應(yīng)用效果最好,洪峰相對誤差和峰現(xiàn)時間誤差最小,分別為19.7%和3 h,遠(yuǎn)低于TOPKAPI模型和TOPMODEL模型的模擬誤差。TOPMODEL模型和TOPKAPI模型在板橋流域的洪水模擬結(jié)果不是很好,確定性系數(shù)多為負(fù)值,其中TOPKAPI模型對半濕潤地區(qū)板橋流域的水文響應(yīng)不是很及時,這一點(diǎn)主要體現(xiàn)在模擬與實(shí)測的峰現(xiàn)時間大多不一致,這樣的模擬結(jié)果可能是由于TOPKAPI模型采用蓄滿產(chǎn)流機(jī)制,而板橋流域是半濕潤地區(qū),產(chǎn)流機(jī)制是蓄滿產(chǎn)流和超滲產(chǎn)流的綜合作用,模型產(chǎn)流機(jī)制與實(shí)際地區(qū)有出入,在一定程度上影響了模型模擬精度。

基于表8~表10中的統(tǒng)計(jì)結(jié)果可知,在馬渡王流域中,TOPMODEL模型、TOPKAPI模型和CASC2D模型的平均確定性系數(shù)為0.26、0.61和0.04。以相對誤差控制在20%以內(nèi)作為模擬合格的標(biāo)準(zhǔn),TOPKAPI模型的洪峰相對誤差僅為3.4%,明顯優(yōu)于其他模型,其確定性系數(shù)在0.7以上的有9場,11場洪水的洪峰相對誤差合格率為100%,峰現(xiàn)時間誤差基本控制在3 h以內(nèi),可以認(rèn)為TOPKAPI模型在馬渡王流域中水文響應(yīng)及時,模擬結(jié)果良好。

4.2 結(jié)果分析

為進(jìn)一步探究這3種模型的模擬結(jié)果差異產(chǎn)生的原因,筆者從板橋流域和馬渡王流域的降雨—徑流關(guān)系、地形指數(shù)特征以及土壤植被類型3個方面進(jìn)行了研究分析。

4.2.1 降雨—徑流關(guān)系

板橋流域和馬渡王流域的降雨與徑流的關(guān)系如圖7所示,圖中P為時段降雨量,Pa為時段前期降雨量,R為徑流深。降雨—徑流關(guān)系的趨勢線接近45°時,表明流域內(nèi)以蓄滿產(chǎn)流為主導(dǎo),即馬渡王流域更容易發(fā)生蓄滿產(chǎn)流。結(jié)合降雨歷史觀測數(shù)據(jù),雖然板橋的降雨量大于馬渡王的,但是實(shí)地勘察發(fā)現(xiàn),由于板橋流域內(nèi)多為巖石地質(zhì)構(gòu)造,多荒山,植被覆蓋程度較低,因此降雨量雖大,卻不易下滲,極易在表面發(fā)生超滲產(chǎn)流。

4.2.2 地形指數(shù)特征

利用TOPMODEL模型中計(jì)算地形指數(shù)的模塊計(jì)算板橋和馬渡王流域的地形指數(shù)及其對應(yīng)的面積比例,繪制流域地形指數(shù)—面積比例關(guān)系曲線,如圖8所示。

從圖8可以看出,占馬渡王流域面積比例最大的地形指數(shù)的值為8.7,小于板橋流域中的該值(約為9),曲線明顯偏右,這說明流域內(nèi)地形指數(shù)大的區(qū)域比例較高,坡度更為平緩,匯流面積更大。與板橋流域相比,在相同的時間段內(nèi),馬渡王流域在平緩的坡面上有更多的水流經(jīng)過。

4.2.3 土壤植被類型

由圖6可以明顯看出,板橋流域中存在大面積的莊稼地,相比于植物叢生、樹冠高大、遮陰面積較大的森林,莊稼地的地表較為裸露,水流速度較快,更容易發(fā)生超滲產(chǎn)流,加之板橋流域的土壤中石礫成分較多,下墊面涵養(yǎng)水源的能力較差,更不利于蓄滿產(chǎn)流的發(fā)生;馬渡王流域中針葉林、闊葉林等所占的面積比例較大,且實(shí)地勘察發(fā)現(xiàn),馬渡王流域中發(fā)育好的植被類型所占比例明顯低于板橋流域。

板橋流域地處半濕潤地區(qū),降雨量年際變化較大,降雨隨季節(jié)變化明顯,流域內(nèi)山體坡度大,匯流快,暴雨過后,河道流量迅速減小,洪水陡漲陡落。加之該流域內(nèi)巖石廣布,導(dǎo)致匯流時間較短,徑流成分以地面徑流為主,產(chǎn)流方式以超滲產(chǎn)流為主。超滲產(chǎn)流區(qū)突發(fā)性洪水居多,洪量較小,但洪峰相對較大,洪水過程線陡漲陡落。TOPKAPI與CASC2D模型擁有不同的產(chǎn)流機(jī)制,TOPKAPI模型的產(chǎn)流機(jī)制是蓄滿產(chǎn)流,與實(shí)際地區(qū)有出入,這在一定程度上影響模型模擬精度;而CASC2D模型具有超滲產(chǎn)流機(jī)制,對于時常發(fā)生陡漲陡落的尖瘦型大洪水的板橋流域而言,該模型的應(yīng)用效果良好,確定性系數(shù)最高,洪峰與洪量誤差較小。

馬渡王流域洪水類型多樣,屬于混合產(chǎn)流區(qū),既有緩漲緩落的蓄滿產(chǎn)流型洪水,又有陡漲陡落的超滲產(chǎn)流型洪水,既有大洪水,又有小洪水,流域內(nèi)地形指數(shù)大的區(qū)域所占比例較高,坡度更為平緩,而且土壤植被發(fā)育較好,林地所占比例高,相比之下蓄滿產(chǎn)流發(fā)生的幾率較高,這也可以解釋基于蓄滿產(chǎn)流原理的TOPKAPI模型在馬渡王流域應(yīng)用效果較好的原因。

TOPMODEL的表現(xiàn)不理想,這可能是因?yàn)門OPMODEL模型模擬流量大小與地下水位有關(guān),地下水位水深越淺,流量越大。由于地下水位不會陡漲陡落,因此TOPMODEL模型不太適用于馬渡王流域。

5 結(jié) 語

基于TOPMODEL模型、TOPKAPI模型以及CASC2D模型的原理,分別在板橋流域和馬渡王流域構(gòu)建了相應(yīng)的分布式水文模型,并對流域中的徑流過程進(jìn)行了模擬計(jì)算,最后對這3個模型在兩流域的應(yīng)用結(jié)果進(jìn)行了評價,結(jié)果表明:在常常發(fā)生超滲產(chǎn)流的板橋流域,CASC2D模型的應(yīng)用效果最好,TOPKAPI模型次之,主要考慮單一地形特征的TOPMODEL模型的應(yīng)用效果較差;在蓄滿產(chǎn)流發(fā)生概率較高的馬渡王流域,TOPKAPI模型的表現(xiàn)較好。針對這種現(xiàn)象,基于流域中的降雨—徑流關(guān)系、地形指數(shù)特征以及土壤植被類型進(jìn)行了研究分析,結(jié)果表明,馬渡王流域與板橋流域相比,降雨更為豐富,植被覆蓋條件更好,地形較為平緩,這些因素使得馬渡王流域的徑流系數(shù)更高,更容易發(fā)生蓄滿產(chǎn)流。如何利用TOPMODEL模型、TOPKAPI模型以及CASC2D模型對流域進(jìn)行洪水模擬與預(yù)報(bào)仍然是一個需要進(jìn)一步探索與研究的問題,而這也正是筆者所在研究團(tuán)隊(duì)下一步的研究方向。

參考文獻(xiàn):

[1] 劉昌明,李道峰,田英,等.基于DEM的分布式水文模型在大尺度流域應(yīng)用研究[J].地理科學(xué)進(jìn)展,2003,22(5):437-445,541-542.

[2] 王中根,劉昌明,左其亭,等.基于DEM的分布式水文模型構(gòu)建方法[J].地理科學(xué)進(jìn)展,2002,21(5):430-439.

[3] 宋建,莊瑋.新疆水資源分布的空間匹配性分析:基于洛倫茲曲線和基尼系數(shù)[J].水利規(guī)劃與設(shè)計(jì),2018(7):60-63,144.

[4] 石增福.城市節(jié)水與雨水資源化概述[J].水利規(guī)劃與設(shè)計(jì),2016(11):106,111.

[5] 楊大文,李翀,倪廣恒,等.分布式水文模型在黃河流域的應(yīng)用[J].地理學(xué)報(bào),2004,59(1):143-154.

[6] 王中根,劉昌明,吳險峰.基于DEM的分布式水文模型研究綜述[J].自然資源學(xué)報(bào),2003,18(2):168-173.

[7] 徐宗學(xué),程磊.分布式水文模型研究與應(yīng)用進(jìn)展[J].水利學(xué)報(bào),2010,41(9):1009-1017.

[8] 魏勝利,蘇偉杰.分布式地下水模擬模型MODFLOW介紹與應(yīng)用[J].水利規(guī)劃與設(shè)計(jì),2011(3):30-32.

[9] 王蕾,楊洋,趙彤彬.GIS技術(shù)在水利工程中的應(yīng)用展望[J].水利規(guī)劃與設(shè)計(jì),2018(2):174-176.

[10] ABBOTT M B, BATHURST J C, CUNGE J A, et al. An Introduction to the European Hydrological System-Systeme Hydrologique Europeen, “SHE”, 1: History and Philosophy of a Physically-Based, Distributed Modelling System[J].Journal of Hydrology,1986,87(1-2):45-59.

[11] BATHURST J C, OCONNELL P E. Future of Distributed Modelling: the Systeme Hydrologique Europeen[J].Hydrological Processes,1992,6(3):265-277.

[12] CHRISTIAENS K, FEYEN J. Use of Sensitivity and Uncertainty Measures in Distributed Hydrological Modeling with an Application to the MIKE SHE Model[J].Water Resources Research,2002,38(9):801-815.

[13] CHRISTIAENS K, FEYEN J. Constraining Soil Hydraulic Parameter and Output Uncertainty of the Distributed Hydrological MIKE SHE Model Using the GLUE Framework[J].Hydrological Processes,2002,16(2):373-391.

[14] 劉青娥,夏軍,王中根.TOPMODEL模型幾個問題的研究[J].水電能源科學(xué),2003,21(2):41-44.

[15] 凌峰,杜耘,肖飛,等.分布式TOPMODEL模型在清江流域降雨徑流模擬中的應(yīng)用[J].長江流域資源與環(huán)境,2010,19(1):48-53.

[16] WOLOCK D M, MCCABE G J. Comparison of Single and Multiple Flow Direction Algorithms for Computing Topographic Parameters in TOPMODEL[J].Water Resources Research,1995,31(5):1315-1324.

[17] 劉志雨,謝正輝.TOPKAPI模型的改進(jìn)及其在淮河流域洪水模擬中的應(yīng)用研究[J].水文,2003,23(6):1-7.

[18] CIARAPICA L, TODINI E. TOPKAPI: A Model for the Representation of the Rainfall-Runoff Process at Different Scales[J].Hydrological Processes,2002,16(2):207-229.

[19] LIU Z Y. ArcTOP: A Distributed Hydrologic Modeling System of Tight Coupling TOPKAPI with GIS[J].Hydrology,2005,25(3):325-339.

[20] LI Z J, QU C Y, HUANG P N, et al. Simulation of Runoff in Luanchuan Basin Using CASC2D and GSSHA Models[J].Journal of Hohai University,2017,45(1):1-6.

[21] MARSIK M, WAYLEN P. An Application of the Distributed Hydrologic Model CASC2D to a Tropical Montane Watershed[J].Journal of Hydrology,2006,330(3-4):481-495.

[22] 李致家,胡偉升,丁杰,等.基于物理基礎(chǔ)與基于柵格的分布式水文模型研究[J].水力發(fā)電學(xué)報(bào),2012,31(2):5-13,32.

[23] 張珂,李致家,包紅軍.GTOPMODEL模型與TOPMODEL模型比較[J].河海大學(xué)學(xué)報(bào)(自然科學(xué)版),2005,33(5),509-512.

[24] 姚成,韓從尚,李致家.基于柵格的新安江模型與GTOPMODEL模型[J].水力發(fā)電,2007,33(3):29-31.

[25] DENG H P, SUN S F. Extension of TOPMODEL Applications to the Heterogeneous Land Surface[J].Advances in Atmospheric Sciences,2010,27(1):164-176.

[26] 張漢辰,李致家,晁麗君,等.CASC2D模型和新安江模型的應(yīng)用比較[J].水力發(fā)電,2015,41(8):23-25,61.

[27] 李致家,屈晨陽,黃鵬年,等.CASC2D模型與GSSHA模型在欒川流域的徑流模擬[J].河海大學(xué)學(xué)報(bào)(自然科學(xué)版),2017,45(1):1-6.

[28] 晁麗君,李致家,李巧玲,等.基于CASC2D模型的縣北溝流域徑流模擬研究[J].水電能源科學(xué),2014,32(1):17-20.

[29] FORTIN J P, TURCOTTE R, MASSICOTTE S, et al. Distributed Watershed Model Compatible with Remote Sensing and GIS Data. Ⅱ: Application to Chaudière Watershed[J].Journal of Hydrologic Engineering,2001,6(2):91-99.

[30] BEVEN K, FREER J. A Dynamic TOPMODEL[J].Hydrological Processes,2001,15(10):1993-2011.

主站蜘蛛池模板: 麻豆国产精品视频| 久久a毛片| 国产小视频a在线观看| 2020国产在线视精品在| 久久国产高清视频| 亚洲热线99精品视频| 毛片大全免费观看| 妇女自拍偷自拍亚洲精品| 日韩大乳视频中文字幕| 在线高清亚洲精品二区| 麻豆精品在线| 国产综合在线观看视频| 四虎影视国产精品| 久久无码高潮喷水| 成年网址网站在线观看| 亚洲人成网站在线播放2019| Jizz国产色系免费| 国产亚洲欧美日韩在线一区二区三区| 国产在线观看91精品亚瑟| 日日碰狠狠添天天爽| 永久免费无码日韩视频| 九九热这里只有国产精品| 少妇精品在线| 色视频国产| 国产精品v欧美| 国内a级毛片| 国产精品福利尤物youwu| 国产高清色视频免费看的网址| 第九色区aⅴ天堂久久香| 久久国产V一级毛多内射| 1级黄色毛片| 亚洲大学生视频在线播放| 九色在线视频导航91| 成人午夜天| 亚洲区第一页| 久久久久久尹人网香蕉| 国产女人在线| 久久精品人妻中文系列| 亚洲天堂日韩av电影| 国产成人高清在线精品| 国产啪在线91| 国产精品区视频中文字幕| 国产专区综合另类日韩一区| 国产欧美在线观看精品一区污| 久久伊伊香蕉综合精品| 色噜噜综合网| 国产午夜福利在线小视频| 91日本在线观看亚洲精品| 成人年鲁鲁在线观看视频| 狠狠操夜夜爽| 999在线免费视频| 国产亚洲欧美在线专区| 亚洲码一区二区三区| 亚洲国产综合自在线另类| 午夜日b视频| 久久综合亚洲色一区二区三区 | 国产精品三区四区| 亚洲精品视频网| 精品久久777| 久久免费精品琪琪| 99热这里只有精品5| 综1合AV在线播放| 999精品免费视频| 蜜桃视频一区二区| 国产草草影院18成年视频| 天天综合天天综合| 亚洲区第一页| 四虎综合网| 9啪在线视频| 久久亚洲AⅤ无码精品午夜麻豆| 麻豆精品在线播放| 国产亚洲欧美在线人成aaaa| jizz国产在线| 在线国产91| 中文无码毛片又爽又刺激| 亚洲欧美h| 91视频精品| 谁有在线观看日韩亚洲最新视频| 成人年鲁鲁在线观看视频| 久青草免费在线视频| 熟妇无码人妻| 福利在线一区|