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

山丘區小流域水文模型適用性研究

2018-03-21 03:03:28常凊睿
中國農村水利水電 2018年2期
關鍵詞:模型

王 璐,葉 磊,吳 劍,常凊睿,張 弛

(1.大連理工大學水利工程學院,遼寧 大連 116000;2.中國水利水電科學研究院,北京 100038)

我國地處歐亞大陸東南部,屬東亞季風氣候區,山丘區暴雨頻發,地質地貌條件復雜,加之日益嚴重的災害性天氣和人類活動的強烈影響,山洪災害發生頻繁,是世界上山洪災害最嚴重的國家之一。山洪突發性強,陡漲陡落,一次山洪過程歷時短且成災范圍分散,使得防洪工程措施無法根本上阻擋山洪對危險區的侵襲,因此,開展山洪災害的預警預報工作,使處于山洪災害危險區的居民能夠提前采取預防措施,對于提高山丘區應急管理水平以及降低人民生命財產損失具有重要的現實意義[1]。

基于水文模型的山洪預警預報研究是山洪災害防治工作中重要的非工程措施體系之一,是在在山洪機理分析及水文計算的基礎上,提前判斷山洪發生的時間、地點、規模、危害范圍以及可能造成的損失[2],進而為災害防治工作提供決策支持。葉金印等[3]針對濕潤地區中小河流,建立基于新安江模型和API-Nash模型的2套山洪預報方案,達到良好的應用效果。鄒大偉等[4]利用多源降雨數據,構建了以TOPMODEL模型為核心的暴雨山洪預報系統。薛霞[5]在山東省開展了基于HEC分布式水文模型的臨界雨量計算和山洪預警指標確定的研究,提高了水文模型在山洪預警預報工作中的實用性。Braud[6]采用分布式水文模型分析了法國東南部流域2002年的山洪過程,取得了較好的模擬結果。目前,國內外山洪預警預報研究大多是在某個特定流域選取某種水文模型進行模擬,缺少多個流域、多種模型的對比分析,且水文模型的選取主要依據研究人員的經驗,存在較大的隨機性。

基于此,本文在全國范圍內選擇多個具有代表性的典型山丘區小流域,采用多種產匯流機制不同、復雜程度各異的常見水文模型模擬山丘區小流域暴雨洪水過程,并對比分析模型的模擬結果,總結不同水文模型的適用條件,以期在實際山洪災害防治工作中,能夠根據小流域特點選擇適宜的水文模型,提高山洪預警預報水平,減少山洪災害損失。

1 水文模型概述與參數率定

1.1 水文模型概述

目前用于洪水預報的水文模型很多,有傳統的集總式水文模型,也有較為復雜的半分布式、分布式水文模型。我國地域遼闊,不同山丘區小流域的地理特征、氣候條件、下墊面情況差異顯著,且不同小流域的水文氣象資料掌握情況也不盡相同。因此,為了深入探討水文模型在山丘區小流域的適用性,本文選用5種常用的產匯流機制不同、復雜程度各異的水文模型,包括API水文模型、新安江水文模型和大伙房水文模型3種集總式水文模型、TOPMODEL半分布式水文模型及HEC-HMS分布式水文模型。

1.1.1 API水文模型

API水文模型根據經驗性的降雨徑流關系和單位線模擬洪水的產匯流過程,其核心是在成因分析和數理統計的基礎上,利用多場洪水的流域平均降水量和相應產生的徑流總量以及影響洪水的主要因素(最常用的是前期影響雨量)建立的P~Pa~R相關圖[7],模型原理簡單、參數調試容易,在濕潤地區應用廣泛。

山丘區小流域大多位于偏遠地區,缺乏長期可靠的觀測資料,很難得到規律性強的P~Pa~R相關圖,一定程度上限制了P~Pa~R相關圖的應用。據此,本文對傳統較為復雜的P~Pa~R關系曲線進行了簡化處理。山丘區暴雨洪水突發性強,起漲較快,故忽略只下滲不產流階段,并以是否發生全流域蓄滿產流為轉折點,將產流過程分為2段折線(部分產流和全流域產流)以近似代替光滑的P~Pa~R關系曲線,從而拓展API水文模型在觀測資料較少的山丘區的應用。

1.1.2 新安江水文模型

新安江水文模型采用流域蓄水曲線考慮下墊面不均勻對產流面積變化的影響,將每個單元面積預報的流量過程演算至流域出口并疊加,即為整個流域的預報流量過程,其產匯流計算的重要特點是分水源和分階段[8]。分水源是指以流域自由水蓄水曲線控制徑流的組成和分布,進行水源劃分,即地表徑流、壤中流和地下徑流,3種水源的匯流速度不同,地表最快、地下最慢;分階段是指將匯流過程分為坡面匯流階段和河網匯流階段,坡地計算時,地表徑流、壤中流和地下徑流采用線性水庫方式,最后采用馬斯京根法進行河道匯流。新安江水文模型是典型的蓄滿產流模型,在我國濕潤、半濕潤地區得到了廣泛的應用。

1.1.3 大伙房水文模型

大伙房水文模型產流部分采用8參數超滲產流計算模型,引入雙層入滲曲線方程進行下滲損失計算,以拋物線方程的形式描述土壤上層蓄水量和雙層下滲率的分布情況;匯流部分則采用8參數變強度和變匯流速度的經驗單位線進行計算,各個計算時段信息不同,得到的單位線形狀和匯流速度均不同,錯時段疊加得到流域總出流過程[9];蒸散發量按照2層模式進行計算,表層按照流域平均蒸散發能力蒸發,待表層需水量蒸發完以后,下層開始蒸發。該模型主要適應于北方干旱和半干旱地區,近些年在北方一些流域已經得到廣泛的應用,并取得了良好的效果[10,11]。

1.1.4 TOPMODEL水文模型

TOPMODEL水文模型以地形為基礎,利用每個單元網格的地形指數反映水文過程。模型的核心是變動產流面積的概念[12],當包氣帶中的含水量達到飽和含水量,即滿足完全重力排水含水量時,土壤中的水都變成自由水完全在重力的作用下流動。由于垂直排水及流域內的側向水分運動,一部分流域面積地下水位抬升至地表面產生飽和坡面流。在整個降水過程中,飽和源面積是不斷變化的[13],其位置受流域地形和土壤水力特性2個因素的影響。模型結構簡單,對資料要求不高,在國內濕潤、半濕潤半干旱地區應用效果較好[14]。

1.1.5 HEC-HMS水文模型

HEC-HMS水文模型主要由HEC-GeoHMS、HEC-DSSVue和HEC-HMS 3部分組成[15],3者分工不同,采用松散耦合模式整合。HEC-GeoHMS模塊通過處理研究流域DEM、土地利用、土壤類型等分布式數據,能夠更為真實地反映出山丘區小流域下墊面的特點[16];HEC-DSS是HEC模型系統中的數據管理系統,可以有效存儲和更新降雨徑流數據;HEC-HMS水文模型對實際產流過程進行簡化,不考慮河道與地面和地下含水層的相互作用,把壤中流和地面徑流匯合為直接徑流。模型結構簡單,可綜合考慮氣候和下墊面條件,同時,模型可以以根據不同的流域、資料情況以及計算要求選擇不同的計算方法,具有廣泛的適應性[17]。本文選擇初損后損法計算產流,Snyder單位線計算匯流,運動波進行河道演算,該組合方法原理簡單,參數較少,在缺資料或無資料地區易于推廣。

1.2 目標函數選取及參數率定

傳統的水文模型參數率定主要采用單一目標函數,然而水文模型的實際應用表明,基于單目標的水文模型參數優化率定僅僅考慮了水文過程某一方面的特征,不能充分反映水文系統的不同動力學行為特征。為了使優化的參數能更好地反映流域水文特征,本文基于洪水預報精度評定標準選擇4個目標函數,從洪量、洪峰、洪水過程線3方面對模擬結果進行評價。指標的具體計算方法如下。

(1)徑流深平均相對誤差:

(1)

式中:Rsim、Robs分別為場次徑流深的模擬值與實測值,mm;Rrel為場次相對誤差,%,所有場次取其平均值即為流域徑流深平均相對誤差。

(2)洪峰流量平均相對誤差:

(2)

式中:Qsim、Qobs分別為場次洪峰流量的模擬值與實測值,m3/s;Qrel為場次洪峰流量相對誤差,%,所有場次取其平均值即為流域洪峰流量平均相對誤差。

(3)峰現時間差:

ΔT=Tsim-Tobs

(3)

式中:Tsim、Tobs分別為峰現時間的模擬值與實測值;ΔT為峰現時間差,h,2 h內視為合格。

(4)平均確定性系數:

(4)

將所有場次洪水相接得到總流量過程線,此流量過程線的確定性系數即為流域平均確定性系數。

在上述目標函數選取的基礎上,本文采用多目標優化算法eNSGA-II進行水文模型的參數率定,以獲得更為準確的模型參數。eNSGA-II是目前應用較為廣泛的多目標遺傳算法之一,它降低了非劣排序遺傳算法的復雜性,并保證優良種群個體在進化過程中不會被丟棄[18,19],從而提高優化結果的精度,具有運行速度快、解集的收斂性好等優點。

2 模型應用與結果對比

2.1 研究流域選取

我國地域遼闊,地貌類型復雜,季風氣候特征多樣,不同區域的產匯流機制不同導致山洪災害地域差異明顯。因此,應該在全國范圍內盡可能選擇具有代表性的山丘區小流域。本文從河南、河北、吉林、甘肅、福建、江西、浙江7省共選取14個典型山丘區小流域開展水文模型對比研究,流域空間分布見圖1。

圖1 研究流域空間分布Fig.1 Spatial distribution of study watersheds

選取的研究流域整體呈南北分布,橫跨我國地形的第2、3階梯,涵蓋了不同的徑流帶、溫度帶,充分代表了我國自然地理特征的空間分布和差異性,流域具體信息見表1。

表1 研究流域基本情況Tab.1 Watersheds for simulation and comparion

2.2 模型結果對比與討論

在選取的14個山丘區小流域采用eNSGA-II優化算法率定5種水文模型的模型參數,將優選的參數進行水文模擬,結果見表2。

(1)API水文模型。API水文模型效果主要受P~Pa~R相關圖影響,在福建、江西及浙江等濕潤小流域,包氣帶缺水量少,蓄滿產流型洪水居多,P~Pa~R實測點分布規律較為明顯,易于擬合曲線,故模型在濕潤流域表現較好。然而,由表2可知,API水文模型在福建永春流域的模擬效果要劣于其余模型,這主要由于永春流域降雨資料不連續,洪水過程的前期實測資料不完整導致Pa的計算存在明顯偏差,見圖2(a),P~Pa~R曲線擬合程度低,嚴重影響API水文模型的精度。類似的,在福建省小安下流域和吉林省周戶屯流域,缺失洪水過程的前期降雨數據且實測資料年限較短,用于擬合曲線的場次數據點不足,見圖2(b)~(c),P~Pa~R曲線不具可靠性,導致小安下流域模擬出現確定性系數為負的情況,周戶屯流域模擬結果也要明顯劣于其余模型。對于半濕潤半干旱小流域,超滲與蓄滿產流模式并存,尤其在短歷時強降雨條件下,極易發生超滲產流型洪水,其P~Pa~R實測點相關性較差,且大多集中在未蓄滿產流階段,影響模擬精度,例如在河南下河村、河北王岸及甘肅華亭流域,API水文模型效果明顯劣于其余模型,確定性系數為負。

由此可見,API水文模型作為典型的蓄滿產流模型,對實測降雨徑流資料要求較高,并不是所有的濕潤區小流域都適用,若降雨徑流資料不連續或洪水場次較少,導致P~Pa~R曲線擬合較差,API水文模型則不能達到較好的模擬效果。

表2 5種水文模型在14個流域的模擬結果Tab.2 Simulation results for fourteen basins of five models

(2)新安江水文模型。本文選取的山丘區小流域大多年平均降雨量在800 mm以上,地處濕潤地區,基于蓄滿產流機制的新安江水文模型適用性較高,同時模型采用的張力水蓄水容量曲線和表層自由水蓄水容量曲線呈拋物線分布,充分考慮了流域下墊面的空間分布特征,即使前期土壤較為干旱,也不會發生只降雨不產流的情況,因此在大多數流域均能達到較好的模擬效果。

圖2 P~Pa~R關系曲線Fig.2 P~Pa~R curve

然而在洪水過程前期實測資料短缺的福建省小安下流域和吉林省周戶屯流域,新安江水文模型的產流評價指標(徑流深平均相對誤差)為眾模型最差,平均相對誤差比其余模型高出近20%到30%,分析認為由于洪水前期降雨數據的缺失導致前期土壤含水量的計算嚴重偏低,進而影響了水文模擬精度。圖3和4分別為小安下及周戶屯流域的典型洪水過程,流域平均張力水蓄水容量Wm取值已在100 mm左右(正常范圍120~200 mm[20]),但仍有大部分降雨損失,產流較小,造成模擬結果的嚴重偏低。在前期降雨資料同樣短缺的永春流域,盡管新安江模擬結果與其余模型精度相近,但其Wm取值為95 mm,模型參數的物理意義也不切合實際。

(3)大伙房水文模型。基于超滲產流機制的大伙房水文模型在東北地區的干旱半干旱流域應用廣泛,而本文選取的研究流域大多植被覆蓋率高且水系發育程度較好,故模型產流模擬效果一般,但大伙房水文模型對于洪峰流量和峰現時間的模擬效果較為突出。對于濕潤小流域,其洪峰流量相對誤差基本在20%左右,對于半濕潤半干旱小流域,其洪峰模擬精度處于眾模型的中等水平,尤其在河北省王岸及西臺峪流域,徑流深模擬最差的情況下,洪峰流量的模擬效果處于中等。分析認為大伙房水文模型匯流計算采用變參數、變匯流速度的經驗單位線,具有起漲迅速、退水快的特征,適用于山丘區小流域突發性洪水的峰值模擬。由此可見,大伙房水文模型在山洪預報中仍有較大的應用潛力,但后續需要改進其產流模塊以提高洪量模擬精度。李致家等[21]研究發現,即使在超滲產流區,也存在蓄滿產流因素,超滲產流模型能準確模擬洪峰,但難以準確模擬洪量,這與本文大伙房模型的結果一致。

圖3 小安下流域典型洪水過程模擬Fig.3 Simulation of flood in Xiaoanxia basin

圖4 周戶屯流域典型洪水過程模擬Fig.4 Simulation of flood in Zhouhutun basin

(4)TOPMODEL水文模型。TOPMODEL水文模型在濕潤小流域應用較為成功,而在半濕潤半干旱小流域適用性一般。這是因為TOPMODEL模擬的流量大小與地下水位有關,地下水位越淺,模擬流量越大,而地下水位往往不會陡漲陡落[21]。在濕潤小流域,基流通常維持在相對較高的水平,流域平均地下水深較淺,大部分場次洪水起漲較緩,陡漲陡落型洪水較少,故模型整體產匯流模擬效果好,相對誤差基本保持在20%以內,確定性系數也在0.7以上;而在經常出現無降雨期的半濕潤半干旱流域,流域基流較小,洪水前期土壤常常處于干旱狀態,易發生陡漲陡落型洪水,模型應用存在局限。

圖5為前期無雨的4場典型陡漲陡落洪水,其起漲流量均在1 m3/s左右,模型認為初始地下水深由初始流量決定,因而認為這4場的初始地下水位很深。由結果可看出,對此類型洪水,洪峰流量模擬值遠低于實測值,且實測流量越高,模擬效果越差。黃鵬年等[22]通過分析TOPMODEL模型在半濕潤半干旱流域的模擬結果也得出流域平均地下水深太深,模型難以模擬洪水峰值的類似結論。

(5)HEC-HMS水文模型。HEC-HMS水文模型在14個小流域均表現突出,山洪模擬適用性強。在半濕潤半干旱小流域,模型產流模塊中的參數初損考慮植被截流、地面填洼等損失,后損則代表流域的最大下滲能力,計算方法屬超滲產流機制,相比于長歷時多峰型洪水,更適用于局地暴雨條件下的短歷時山洪模擬[23],因而在此類型洪水居多的欒川、下河村、王岸等半濕潤半干旱小流域的水文模擬結果要優于其余4種水文模型。王莉莉等[24]對比新安江、TOPMODEL及HEC-HMS在半濕潤半干旱研究流域的模擬結果,也表明HEC-HMS水文模型更具優勢。在濕潤小流域,產流參數初損可理解為流域初始缺水量,后損為土壤飽和下滲率,產流計算損失整體比半濕潤半干旱區偏大,在一定程度上體現蓄滿產流模式[25],故達到與蓄滿產流模型相近的模擬效果。

圖5 陡漲陡落型洪水過程模擬Fig.5 Simulation of rapid flood

表3 欒川流域“20050817”洪水模擬評價指標Tab.3 Simulation index of “20050817” flood in Luanchuan basin

圖6 欒川流域“20050817”洪水過程模擬Fig.6 Simulation of “20050817” flood in Luanchuan basin

同時,在甘肅省華亭流域,見圖7,共分布3個雨量站,其中趙莊和王寨并列位于流域中游,華亭位于出口處,站點分布極為不均。對此,由表2可知,HEC-HMS分布式水文模型對徑流深及洪峰的模擬效果優于其余集總式降雨輸入模型。

圖7 華亭流域Fig.7 Huating basin

3 結 語

本文將5種常用水文模型應用于濕潤及半濕潤半干旱地區的14個典型山丘區小流域,開展模型適用研究。研究發現,在選取的5種水文模型當中,大伙房水文模型對山丘區暴雨洪水的適用性最低,盡管模型對洪峰模擬精度較高,但其產流模塊區域特征鮮明,在植被覆蓋良好、水系發育程度較高的山丘區小流域產流模擬效果較差,影響模型的整體應用。

在濕潤小流域,API、新安江、TOPMODEL及分布式水文模型HEC-HMS模擬效果較好,但API和新安江水文模型對降雨徑流資料要求高,TOPMODEL對超滲產流比重大的陡漲陡落洪水應用較為局限。因此,若小流域實測資料充 足,在濕潤小流域可優先選擇建模簡單的API、新安江集總式水文模型,但山丘區小流域大多缺少實測資料,而HEC-HMS分布式水文模型可根據下墊面數據確定大部分參數,對水文資料依賴性較低,更適用于缺資料或無資料小流域。

在半濕潤半干旱小流域,短歷時強降雨條件下超滲產流型洪水居多,HEC-HMS水文模型表現突出,在眾模型中適用性最強。但在部分小流域,表現最好的HEC-HMS水文模型也難以滿足洪水預警預報要求。這主要由于半濕潤半干旱小流域內超滲與蓄滿產流機制并存,暴雨洪水過程復雜多變,導致現有水文模型應用存在較大的局限。摸清半濕潤半干旱區流域混合產流時空轉換機制,完善降雨-徑流預報研究基礎理論,是開發改進山丘區半濕潤半干旱小流域水文模型亟需解決的問題。

[1] 姬 晶, 劉 攀, 江炎生, 等. 一種無資料地區山洪災害臨界雨量計算方法研究[J]. 中國農村水利水電, 2017,(1):76-79.

[2] 王文川, 和 吉, 邱 林. 我國山洪災害防治技術研究綜述[J]. 中國水利,2011,(13).

[3] 葉金印, 吳勇拓, 李致家,等. 濕潤地區中小河流山洪預報方法研究與應用[J]. 河海大學學報(自然科學版), 2012,40(6):615-621.

[4] 鄒大偉, 張興強, 張 艷. 基于TOPMODEL模型的大汶河流域暴雨誘發山洪的預報研究[C]∥ 第33屆中國氣象學會年會S 水文氣象災害預報預警. 2016.

[5] 薛 霞. 山洪防治區小流域分類與預警指標確定研究[D]. 濟南:山東大學, 2016.

[6] Braud I. The use of distributed hydrological models for the Gard 2002 flash flood event: analysis of associated hydrological processes[J]. Journal of Hydrology, 2010,394(1-2):162-181.

[7] 袁作新. 流域水文模型[M]. 北京: 水利電力出版社, 1990.

[8] 趙人俊. 流域水文模擬——新安江模型和陜北模型[M]. 北京:水利電力出版社, 1984.

[9] 大連理工大學, 國家防汛抗旱總指揮部辦公室. 水庫防洪預報方法與應用[M]. 北京: 中國水利水電出版社, 1996.

[10] 王 猛, 彭 勇, 梁國華. 大伙房模型在石佛寺流域洪水預報方案中的應用[J]. 南水北調與水利科技, 2012,(2):93-97.

[11] 丁 陽. 大伙房水庫洪水預報模型的研究[J]. 東北水利水電, 1999,(4):7-9.

[12] Vijaz P Singh. Computer models of watershed hydrology[M]. Water Resources Publications, 1999.

[13] 劉青娥. 流域水文模擬的TOPMODEL模型及其應用研究[D]. 武漢:武漢大學, 2003.

[14] 彭 偉. 基于三種水文模型的流域徑流模擬和土壤含水量模擬應用研究[D]. 成都:四川農業大學, 2009.

[15] Hydrologic Engineering Center. Hydrologic modeling system HEC-HMS users manual[Z]. Hydrologic Engineering Center, 2000.

[16] US Army Corps of Engineers Hydrologic Engineering Center. HEC-GeoHMS Users Manual: Version 1.0[Z]. US Army Corps of Engineers Hydrologic Engineering Center, 2000.

[17] 董小濤, 李致家. HEC模型在洪水預報中的應用[J]. 東北水利水電, 2004,(11):43-44.

[18] Srinivas N, Kalyanmoy Deb. Multi-objective optimization using non-dominated sorting in genetic algorithms[J]. Evolutionary Computation, 1994,2(3):221-248.

[19] Kalynanmoy Deb, Amrit Pratap, Sameer Agarwal, et al. A fast and elitist multi-objective genetic algorithm: NSGA-II[J]. IEEE Transactions on Evolutionary Computation, 2002,6(2):182-197.

[20] 徐宗學. 水文模型[M]. 北京: 科學出版社, 2009.

[21] 李致家, 黃鵬年, 姚 成,等. 靈活架構水文模型在不同產流區的應用[J]. 水科學進展, 2013,(11).

[22] 黃鵬年, 李致家, 姚 成,等. 半干旱半濕潤流域水文模型應用與比較[J]. 水力發電學報, 2013,32(4):4-9.

[23] 梁 睿. HEC-HMS水文模型在北張店流域的應用研究[D]. 太原:太原理工大學, 2012.

[24] 王莉莉, 李致家, 包紅軍. 三個水文模型在沂河洪水模擬中的分析[J]. 水電能源科學, 2008,26(4):47-49.

[25] 董小濤, 李致家, 李利琴. 不同水文模型在半干旱地區的應用比較研究[J]. 河海大學學報(自然科學版), 2006,34(2):132-135.

[26] 馬玉軍. 碧流河水庫流域洪水預報方案及設計洪水過程修勻研究[D]. 大連: 大連理工大學, 2013.

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 久久亚洲国产视频| 爽爽影院十八禁在线观看| 精品人妻无码中字系列| 大学生久久香蕉国产线观看| 91在线播放免费不卡无毒| 666精品国产精品亚洲| 欧美日本在线观看| 国产成人AV男人的天堂| 国产特级毛片| 91美女在线| 东京热高清无码精品| 色网在线视频| 一级爆乳无码av| AV不卡无码免费一区二区三区| 在线看AV天堂| 制服丝袜国产精品| 99久久亚洲精品影院| 成人a免费α片在线视频网站| 99999久久久久久亚洲| 波多野结衣一区二区三区AV| 免费在线看黄网址| 尤物亚洲最大AV无码网站| 女同国产精品一区二区| 波多野结衣久久高清免费| 国产网站黄| 91日本在线观看亚洲精品| 中文字幕乱码中文乱码51精品| 国产精品无码AV片在线观看播放| 在线播放国产99re| 538精品在线观看| 国产区在线观看视频| 再看日本中文字幕在线观看| 久久国产拍爱| 国产美女丝袜高潮| 黄色网页在线播放| 国产成人h在线观看网站站| 午夜国产不卡在线观看视频| 中文字幕中文字字幕码一二区| 巨熟乳波霸若妻中文观看免费| 亚洲人成日本在线观看| 欧美日韩福利| 欧美一区二区三区国产精品| 久久亚洲国产最新网站| 欧美一区精品| 亚洲美女操| 91精品小视频| 伊人久久精品无码麻豆精品| 欧美成人日韩| 国产人成在线观看| 欧美天堂在线| 国产麻豆精品在线观看| 99热国产在线精品99| 五月天久久综合| 国产美女叼嘿视频免费看| 国产日本一区二区三区| 亚洲欧美不卡中文字幕| 国产精品专区第1页| 亚洲免费福利视频| 亚洲国产成人久久77| 欧美一级黄片一区2区| 亚洲国产黄色| 免费高清毛片| 日本黄色不卡视频| 国产99久久亚洲综合精品西瓜tv| 精品久久久久无码| 精品偷拍一区二区| 日本91视频| 欧美精品aⅴ在线视频| 午夜福利亚洲精品| 欧美伦理一区| 国产精品无码久久久久AV| 日本免费一级视频| 亚洲AV无码久久天堂| 亚洲综合18p| 99中文字幕亚洲一区二区| 亚洲男人的天堂在线观看| 91麻豆精品国产高清在线| 成人午夜天| 国产麻豆aⅴ精品无码| 亚洲第一在线播放| 免费在线播放毛片| 国产精品3p视频|