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

基于SWAT模型的淮河上游土地利用變化情景的水文響應

2016-12-27 12:09:09董曉華方燕琴李英海喻丹蘇歡李璐
南水北調與水利科技 2016年5期

董曉華 方燕琴 李英海 喻丹 蘇歡 李璐 姚著喜

摘要:淮河流域不合理的土地利用是導致水土流失的重要因素之一,探究土地利用變化下的水文響應規律是開展水土保持工作的重要前提。以淮河王家壩水文站以上流域為研究區域構建SWAT模型,并驗證模型的適用性;對研究區耕地進行適宜性評價,結合2000年土地利用現狀構造最佳土地利用情景;利用驗證好的SWAT模型分析最佳土地利用情景、生態保護情景、土地利用現狀情景的水文響應。結果表明,相比2000年土地利用,生態保護情景水文響應變幅較小,而最佳土地利用情景下4129%的耕地轉化為林地,年、月均地表徑流和年、月均泥沙量減小,年、月均地下徑流增加,年、月均蒸發量變化不大。說明最佳土地利用狀態有助于提高研究區域涵養水源、保持水土的能力,且在汛期效果較為明顯

關鍵詞:淮河流域;土地利用情景;SWAT模型;適宜性評價;水文響應

中圖分類號:P339;P343.9 文獻標志碼:A 文章編號:16721683(2016)05003209

全球氣候變化及人類活動造成了環境變化的影響[1]。在長時期尺度上,氣候變化是造成環境變化的主要因素[2],但在短期尺度上,人類活動造成的土地利用方式的變化對水文的影響更加明顯[34]。土地利用/覆被變化(LUCC)通過改變地表蒸發、地表徑流、土壤水分的入滲能力等因素,影響了流域的水文循環及水量平衡,進而對生態環境產生影響[1]。目前,水文模型方法被廣泛應用于LUCC水文響應研究中,這一研究多以歷史土地利用或現狀土地利用為基準,建立多種土地利用情景,如全耕地、全林地的極端土地利用方式[4]以及根據地區工農業發展設置不同比例的土地利用[5]、根據水土保持政策優化土地利用空間分布[6]等情景,進而分析不同土地利用情景對水文過程的影響。然而這種假設的情景并不能充分代表實際情況下土地利用的演變結果。

SWAT模型作為一種發展和應用都相對較為成熟的水文模型[710],以水量平衡為驅動力,能夠定量描述不同土地利用/覆被以及氣候條件下水文循環各組分的過程,具有較強的物理機制。本文以淮河流域王家壩以上部分流域為研究對象,結合研究區土地利用變化驅動因子進行耕地適宜性評價,對研究區土地利用現狀中耕地不適宜的區域進行退耕還林,以此得到符合研究區土地利用發展方向的最優土地利用狀態。并以SWAT模型作為研究工具,比較最佳土地利用、現狀土地利用以及生態保護情景下土地利用的水文響應,為土地管理及規劃提供了有效依據。

1 研究區概況與數據來源

1.1 研究區概況

淮河流域地處我國東部,介于長江和黃河兩流域之間,位于東經111°55′-121°25′,北緯30°55′-36°36′,面積約為27萬km2。淮河干流發源于桐柏山太白頂北麓,依次流經河南、湖北、安徽、江蘇。淮河流域地處我國南北氣候過渡帶,淮河以北屬暖溫帶區,淮河以南屬北亞熱帶區,氣候溫和,年平均氣溫為11 ℃~16 ℃。氣溫變化由北向南,由沿海向內陸遞增。極端最高氣溫達445 ℃,極端最低氣溫達-241 ℃。淮河流域的氣候特點是冬春干旱少雨,夏秋悶熱多雨,冷暖和早澇轉變急劇。流域多年平均降水量約為920 mm,其分布狀況大致是由南向北遞減,山區多于平原,沿海大于內陸。流域西部、西南部及東北部為山區、丘陵區,其余為廣闊的平原。山丘區面積約占總面積的三分之一,平原面積約占總面積的三分之二。隨著人類社會經濟的高速發展與生態環境問題的日益突出,淮河流域的生態脆弱問題逐漸受到重視[11],尤其是淮河流域水土流失的形勢更加不容樂觀。究其原因,淮河流域土地利用方式單一、林地的掠奪式開發以及城市用地的擴張等不合理的土地利用情況是導致水土流失的重要因素之一。

研究區位于淮河流域上游,總面積約30 630 km2,研究區內分布著138個雨量站,4個流量站(圖1),本文不考慮蓄、滯洪區設施的應用對水文的影響。對研究區有影響的水庫有6個,分別是南灣水庫、石山口水庫、五岳水庫、潑河水庫、板橋水庫、宿鴨湖水庫。由于宿鴨湖水庫有桂莊和夏屯兩個出流點,因此流域中水庫共有7個出流點。

1.2 數據來源

本文所需要的主要數據包括地形、土壤、氣象和水文、土地利用數據等,具體數據類型及來源見表 1。其中,坡度數據是利用ArcGIS中的slope工具從DEM中提取得到;降雨數據是根據研究區內138個雨量站的日實測數據算術平均處理得到;由于無研究區氣象站點的相關資料數據,[HJ1.95mm]因此通過建立氣象數據庫來模擬研究區的氣候特征;徑流數據采用研究區內4個水文站(息縣站、潢川站、班臺站、王家壩站)的日徑流數據。

2 淮河上游流域SWAT模型構建

2.1 模型基本原理

SWAT(Soil and Water Assessment Tool)[1214]模型是由美國農業部(USDA)開發的一個基于物理機制的、連續性的流域尺度水文模型。該模型是基于GIS的分布式水文模型,可模擬流域中多種不同的水文物理化學過程,并可模擬和預測下墊面因素及氣候因素變化下流域水文循環的影響[15]。

式中:SWt表示土壤最終含水量;SW0表示土壤初始含水量;Rday表示第i天的降水量; Qsurf表示第i天的地表徑流;Ea表示第i天的蒸發量;Wseep表示第i天透過土壤剖面底部的滲透量和測流量;Qgw表示第i天的地下水量;單位均為mm。t表示時間步長(d)。

模型在模擬過程中,基于DEM數據,確定流域總出口位置為王家壩水文站。考慮到水文站點及水庫,共將流域劃分為11個子流域(圖 1),并設置土地利用類型、土壤類型、坡度的閾值,將子流域劃分多個水文響應單元(HRU)。[HJ2.05mm]由于本文研究的是土地利用變化對水文過程的影響,設土地利用類型的閾值會導致面積的重新分配,最終無法確定土地利用變化的水文響應是重新分配產生還是真實變化的結果,因此土地利用類型、土壤類型、坡度的閾值分別設置為0%、10%、20%,將11個子流域劃分為221個水文響應單元。

2.3 模型的率定和驗證

采用SWAT自帶的參數率定模塊,使用SCEUA[17]算法對以上選取的10個敏感參數進行率定。選取2000年為模型預熱期,2001-2005年為率定期,2006-2010年為驗證期對研究區進行月徑流模擬,并以流域出口王家壩水文站流量數據為觀測流量,采用NashSutcliffe(Ens)和R2決定系數對模型模擬精度進行評價,計算公式如下:

3 淮河上游流域土地適宜性評價方法

土地適宜性評價是對于某種土地用途和利用方式的適宜性程度進行評估[18],以獲得適用于研究目的的最佳土地利用類型[19]。由于流域水土流失主要原因是耕地的不合理利用,因此本次適宜性評價的對象為耕地。

3.1 評價因子及評價單元的確定

評價因子的合理選擇對評價起主導限制作用,選擇比較穩定且可度量的因素作為評價因子,是評價工作的關鍵[20]。考慮研究區的實際情況以及資料的完備程度,本次研究從地形條件、土壤條件、氣象條件和水源條件四大方面分別選取高程、坡度、地貌、土壤有機質含量、土壤質地、土壤pH值、多年平均降雨量、≥10°積溫、距水源的距離等因子為研究區耕地適宜性評價因子。

評價單元是土地適宜性評價的基本單位,常采用疊置法、網格法、地塊法、多邊形法等方法進行劃分。本次研究采用網格法劃分,每一個柵格網格即為一個評價單元。確定合適的網格大小是網格法劃分評價單元的關鍵[21]。因所選研究區的評價因子中網格單元最小的是90×90 m高程因子圖,因此以該網格大小為標準,使用ArcGIS的Resample工具將所有的因子重采樣成90×90 m的柵格大小。

3.2 評價因子權重

權重表示評價因子對評價對象影響程度的貢獻率。結合德爾菲法(專家打分法)和層次分析法確定各評價因子的權重。首先將所選取的因子分成幾個層次(見圖 2),判斷出各層次兩兩指標的重要程度,并建立判斷矩陣,其中階數大于2的矩陣需進行一致性檢驗,一致性比率CR需小于01,否則重新構建判斷矩陣進行計算。經計算得出耕地權重值,見表 3。

3.3 評價因子的分級量化

采用等級賦分法[22]將研究區各因子依據不同的適宜范圍進行分級劃分,然后對每級賦予相應的分值(0~100分),并分成5個等級,分值越大,則其適宜程度越高。根據國家農用地分等規程[23]及研究區域的特征,并參考文獻[2426]及專家意見,評價因子分級及其作用分值見表 4。

3.4 土地適宜性評價模型

根據所選因子及計算的權重值,建立研究區土地適宜性評價模型——加權指數和模型如下:

采用ArcGIS中的Raster Calculator工具,根據式計算研究區各評價單元的適宜性總分值,并自動生成評價等級圖。運用ArcGIS重分類中的“自然間斷分類方法[22]”劃分研究區耕地土地適宜性級別,依據FAO《土地評價綱要》中的評價原則,將耕地適宜性總分值劃分為高度適宜(S1)、中等適宜(S2)、勉強適宜(S3)、不適宜(S4)四個級別。其中高度適宜是指土地對該種土地利用方式沒有限制或限制性很小,經濟效益好,能持續利用;中等適宜類指土地對該種土地利用方式的持續利用有中等程度的限制;勉強適宜類指土地對該種土地利用方式的持續利用有較大的限制,經濟效益較差,利用不當會導致水土流失;不適宜類是指在當前的條件下,這類土地對所特定用途不能利用或不能持續利用,否則會加快水土流失的速度。經計算,土地適宜性劃分級別見表 5。

4 計算結果與分析

4.1 SWAT模型模擬結果分析

SWAT模型月徑流模擬結果見圖 3,由圖 3可以看出模擬的流量過程線與實測較為吻合。經計算,率定Ens和R2分別為089、089,驗證期Ens和R2分別089、095。Ens和R2均在08以上,說明SWAT模型能比較準確的模擬該流域的徑流,具有較好的適用性,可用于研究區土地利用變化水文響應的相關研究。

4.2 淮河上游流域耕地適宜性評價

通過對研究區耕地適宜性評價得到耕地適宜性評價等級圖(見圖 4),評價結果統計見表 6。

4.3 最佳土地利用狀態

將淮河上游流域耕地適宜性評價結果柵格圖(圖4)與2000年土地利用現狀圖(圖 5)疊加,即可得到基于土地利用現狀的耕地適宜性等級圖(圖 6)。具體做法是:利用ArcGIS的選擇功能分別從現狀圖的屬性表中評價選取耕地、林地數據并導出,然后用raster calculator分別與耕地適宜性評價結果進行柵格疊加,即可得到基于土地利用現狀耕地適宜性等級圖。

由表 7可知,2000年現狀的耕地以中等適宜和勉強適宜為主,高度適宜區所占總面積為2518%,勉強適宜和不適宜區占3814%。

在土地利用現狀基礎上,耕地適宜性高的地區則保留耕地,適宜性低的地區則退耕還林[27]。對耕地不適宜的地區如果不進行適當整改,可能會加快水土流失的速度,從而逐步破壞到適宜類區域。根據圖 6結果,將評價單元中 S1(高度適宜類)、S2(中等適宜類)的區域予以保留, S3(勉強適宜類)、 S4(不適宜類)的區域進行退耕還林,這些不適宜區域主要是淮河以北地區,如平輿、新蔡、項城、臨泉、新縣東部縣城的部分區域,將這些區域進行退耕還林,即可得到研究區最佳土地利用圖,如圖 7。

4.4 水文響應分析

為分析不同土地利用狀態對水文過程的影響,設置三種情景:一是基準情景,即2000年土地利用圖;二是按照國家相關政策設定的生態保護情景,即坡度15°~25°退耕還草,坡度>25°的退耕還林;三是最佳土地利用狀態。各情景下土地利用變化情況見表 8。

利用已經率定好的模型,分別輸入三種情景的土地利用,保持其他變量和參數不變(固定參數最優值,如表 2所示,除了HRU級別以外參數,其他參數均不變),模擬各情景下地表、地下徑流深、蒸發量及泥沙量的變化情況。部分參數值會因為HRU內要素的變化而變化,如從CN2和SolAwc,土地利用變化導致HRU的重新分配,則每個子流域內的各個HRU的參數值會發生相應的變化。土地利用變化實質導致的是HRU的重新分配以及HRU級別的部分參數(如產流參數CN2)的變化,從而導致最終的產流變化[28],這是模型原理所決定的。

由此可知,退耕還林后,林地對徑流、蒸發量、泥沙量的影響在枯水年表現的不明顯,而在豐水年或平水年,影響較大。當土地朝著最佳利用方向發展時,耕地減少,林地增加,冠層截留及林地枯落物截留使得地表徑流減少;林地的土壤孔隙度較大,增大了下滲率[29],耕地的土壤較為密實,下滲率較低,因此林地地下水的補給較耕地多;林地涵養水源作用較強,且林地增加,植被覆蓋度增大,所以最佳土地利用狀態下泥沙量減少,蒸發量增加。由于研究區耕地基本在平原區,坡度大于15°的面積較小,因此在生態保護情景下,面積變化較小,耕地、林地、草地面積變化率分別為-068%、022%、2568%,所以對水文響應影響不大。

根據上述結果可知,最佳土地利用狀態下,徑流、泥沙量比生態保護情景及現狀少,地下徑流補給增加,但蒸發量變化不明顯,說明在該種狀態下,土地保持水土能力較強,能明顯改善該流域水土流失的現狀。

5 結論

本文通過在淮河流域王家壩水文站以上流域建立SWAT模型,結合研究區土地利用變化驅動因子進行耕地適宜性評價,對研究區現狀中耕地不適宜的區域進行退耕還林,得到土地最佳利用狀態,并以此作為一種土地利用情景,對比生態保護情景及2000年土地利用現狀情景,模擬不同情景下水文響應的變化。研究結果表明,研究區中大部分淮河以北地區的耕地不利于耕作,將不適宜的耕地退耕還林后,流域涵養水源作用、保持水土能力增強,遠超過生態保護情景下的土地利用狀態,能大大改善流域內的生態問題,為土地管理及規劃提供了依據。

本文目前僅對現狀進行比較分析,而對未來自然發展土地利用進行預測,并與優化的土地利用進行水文響應比較,從而為流域土地合理的管理和規劃提供更有針對性的依據,將是下一步研究的方向。

參考文獻(References):

[1] 歐春平,夏軍,王中根.土地利用/覆被變化對SWAT模型水循環模擬結果的影響研究——以海河流域為例[J].水力發電學報, 2009,(04):124129.(OU Chunping,XIA Jun,WANG Zhonggen.LUCC influence on SWAT hydrological simulation—a case study of Haihe river basin[J].Journal Of Hydroelectric Engineering,2009,(04):124129.(in Chinese))

[2] Sanjay K Jain,Jaivir Tyagi,Vishal Singh.Simulation of runoff and sediment yield for a Himalayan watershed using SWAT model[J].Journal of Water Resource and Protection,2010,2(03):267.

[3] Bfw Croke,Ws Merritt,Aj Jakeman.A dynamic model for predicting hydrologic response to land cover changes in gauged and ungauged catchments[J].Journal of Hydrology,2004,291(1):115131.

[4] 張青.基于水文效應分析的土地利用優化模式研究[D].西安:陜西師范大學,2013.(ZHANG Qing.The research on effects of land use optimization pattern based on the analysis of the hydrological[D].Xi′an:Shaanxi Normal University,2013.(in Chinese))

[5] A.Montenegro,R.Ragab.Hydrological response of a Brazilian semiarid catchment to different land use and climate change scenarios:a modelling study[J].Hydrol Process,2010,24(19):27052723.doi:10.1002/hyp.7825.

[6] 羅巧,王克林,王勤學.基于SWAT模型的湘江流域土地利用變化情景的徑流模擬研究[J].中國生態農業學報, 2011, (06):14311436.( LUO Qiao,WANG Kelin,WANG Qinxue.Using SWAT to simulate runoff under different land use scenarios in Xiangjiang River Basin[J].Chinese Journal of EcoAgriculture,2011,(06):14311436.(in Chinese))

[7] Zhixiang Lu,Songbing Zou,Zuodong Qin,et al.Hydrologic Responses to Land Use Change in the Loess Plateau:Case Study in the Upper Fenhe River Watershed[J].Advances in Meteorology,2015,2015.

[8] Hadi Memarian,Siva K.Balasundram,Karim C.Abbaspour,et al.SWATbased hydrological modelling of tropical landuse scenarios[J].Hydrological Sciences Journal,2014,59(10):18081829.doi:10.1080/02626667.2014.892598.

[9] 劉玉明,張靜,武鵬飛,等.北京市媯水河流域人類活動的水文響應[J].生態學報,2012,(23):75497558.(LIU Yuming, ZHANG Jing,WU Pengfei,et al.The hydrological response to human activities in Guishui River Basin, Beijng[J].Acta Ecologica Sinica, 2012,(23):75497558.(in Chinese))[ZK)]

[10] [ZK(#]王學,張祖陸,寧吉才.基于SWAT模型的白馬河流域土地利用變化的徑流響應[J].生態學雜志,2013,(01):186194.(WANG Xue, ZHANG Zulu, NING Jicai.Runoff response to land use change in Baimahe River of China based on SWAT model[J].Chinese Journal of Ecology,2013,(01):186194.(in Chinese))

[11] 葉正偉.基于生態脆弱性的淮河流域水土保持策略研究[J].水土保持通報,2007,(03):141145,156.(YE Zhengwei.Countermeasures for soil and water conservation in the Huaihe River Basin based on ecological vulnerability[J].Bulletin of Soil and Water Conservation,2007,(03):141145,156.(in Chinese))

[12] Jeffrey G Arnold,Peter M Allen,Gilbert Bernhardt.A comprehensive surfacegroundwater flow model[J].Journal of hydrology, 1993,142(1):4769.

[13] Jimmy R Williams,Ad Nicks,Jeffrey G Arnold.Simulator for water resources in rural basins[J].Journal of Hydraulic Engineering, 1985,111(6):970986.

[14] Susan L Neitsch,Jeffrey G Arnold,Jim R Kiniry,et al.Soil and water assessment tool theoretical documentation version 2009[R].Texas Water Resources Institute,2011.

[15] Jeffrey G Arnold,Raghavan Srinivasan,Ranjan S Muttiah,et al.Large area hydrologic modeling and assessment part I:Model development1.Wiley Online Library,1998.

[16] Sl Neitsch,Jg Arnold,Jr E A1 Kiniry,et al.Soil and water assessment tool users manual version 2000[J].GSWRL report,2002, 202(0206).

[17] K Eckhardt,Jg Arnold.Automatic calibration of a distributed catchment model[J].Journal of hydrology, 2001, 251(1):103109.

[18] YanSui Liu,JieYong Wang,LiYing Guo.GISBased Assessment of Land Suitability for Optimal Allocation in the Qinling Mountains,China1[J].Pedosphere,2006,16(5):579586.doi:http://dx.doi.org/10.1016/S10020160(06)60091X.

[19] Fazel Amiri,Abm Shariff.Application of geographic information systems in landuse suitability evaluation for beekeeping:A case study of Vahregan watershed (Iran)[J].African Journal of Agricultural Research,2012,7(1):8997.

[20] 任小寧,郭新成.基于GIS的新疆克拉瑪依市土地適宜性評價[J].安徽農業科學, 2007 (08):23952396,2427.(REN Xiaoning,GUO Xincheng.Land suitability assessment of Kelamayi City based on GIS[J].Journal of Anhui Agri.Sci, 2007 (08):23952396,2427.(in Chinese))

[21] 王令超.農用土地定級方法初探[J].國土資源科技管理, 2001 (01):49.(WANG Lingchao.The research on agricultural land grading method[J].Science and technology of land and resources management,2001 (01):49.(in Chinese))

[22] 馬春霞.基于GIS的土地適宜性評價及開發戰略研究[D].烏魯木齊:新疆大學,2010.(Ma Chunxia.The Research on Suitabiliy for land and the Strategy of Exploitation Based on GIS[D].Urumchi :Xinjiang University,2010.(in Chinese))

[23] Td/T 農用地分等規程[S].2003.(Td/T The Regulations of Farmland Grading[S].2003.(in Chinese))

[24] David G.Rossiter.A theoretical framework for land evaluation[J].Geoderma,1996,72(34):165190.doi:http://dx.doi.org/10.1016/00167061(96)000316.

[25] Bakhtiar Feizizadeh,Thomas Blaschke.Land suitability analysis for Tabriz County,Iran:a multicriteria evaluation approach using GIS[J].Journal of Environmental Planning and Management,2013,56(1):123.

[26] Farag F.Abushnaf,Kevin J.Spence,Ian D.Rotherham.Developing a Land Evaluation Model for the Benghazi Region in Northeast Libya using a Geographic Information System and Multicriteria Analysis[J].APCBEE Procedia,2013,5(0):6975.doi:http://dx.doi.org/10.1016/j.apcbee.2013.05.013.

[27] 陳軍鋒,李秀彬.土地覆被變化的水文響應模擬研究[J].應用生態學報,2004(05):833836.(CHEN Junfeng,LI Xiubing.Simulation of hydrological response to landcover changes[J].Chinese Journal of Appled Ecology,2004(05):833836.(in Chinese))

[28] 胡彩虹,劉麗娜,吳澤寧.未來氣候情景下覆被變化對徑流量的影響[J].人民黃河,2012(01):1921,24.(HU Caihong,LIU Lina,WU Zening.Effect of land cover change to runoff under future climate scenarios[J].Yellow River,2012(01):1921,24.(in Chinese))

[29] 盛前麗.香溪河流域土地利用變化徑流效應研究[D].北京: 北京林業大學,2008.(SHENG Qianli.Research on response of runoff to landuse change in Xiangxi River Watershed[D].Beijing:Beijing Forestry University.2008.

由于研究區內有6個水庫,因此需修改模型的水庫特性。SWAT模型中計算水庫出流量有四種方法:實測日出流量、實測月出流量、無控制水庫的年均泄流量、控制水庫的目標泄流量。其中南灣水庫、石山口水庫、五岳水庫、潑河水庫數據較全,采用實測日出流量方法計算,而板橋水庫、宿鴨湖水庫缺乏實測日出流數據,因此采用控制水庫的目標泄流量計算。通過對水庫數據的準確輸入,使SWAT模型在該流域的模擬更有準確性。

2.2 參數敏感性分析

模型中關于徑流模擬的參數有26個,但是只有少數幾個參數對模擬結果起決定性作用。為了使模型更好的反應研究區特征,提高模型運行效率,需對模型進行敏感性分析。本文選取研究區2000年-2010年的數據,采用SWAT 2005模型自帶的LHOAT[16]參數敏感性分析方法對參數進行敏感性分析,選取敏感度排在前10位的參數,結果見表2。

主站蜘蛛池模板: 国产精品久久自在自2021| 57pao国产成视频免费播放| 暴力调教一区二区三区| 亚洲精品无码成人片在线观看| 成人免费视频一区| 玖玖精品在线| 再看日本中文字幕在线观看| 午夜久久影院| 精品人妻无码中字系列| 热伊人99re久久精品最新地| 久久黄色一级视频| 国产产在线精品亚洲aavv| 婷婷激情亚洲| 国产电话自拍伊人| 国产黄色爱视频| 成人av手机在线观看| 久久国产精品无码hdav| 成人国产小视频| 在线观看免费人成视频色快速| 国产丝袜无码一区二区视频| 有专无码视频| 亚洲无码在线午夜电影| 无码福利视频| 亚洲中文字幕在线观看| 日韩欧美国产另类| 国产最新无码专区在线| 亚洲AV无码一区二区三区牲色| 中国美女**毛片录像在线| 中文字幕久久波多野结衣| 国产精品女人呻吟在线观看| 欧美日韩中文字幕在线| 性网站在线观看| 中文国产成人久久精品小说| 亚洲国产看片基地久久1024| 亚洲高清中文字幕| 她的性爱视频| 91麻豆国产精品91久久久| 亚洲国产综合精品一区| 69精品在线观看| 欧美黄网在线| 麻豆国产在线不卡一区二区| 国产日韩欧美一区二区三区在线| 福利一区在线| 亚洲综合精品第一页| 色欲色欲久久综合网| 国产精品综合久久久| 亚洲资源站av无码网址| 亚洲男人的天堂网| 久久香蕉欧美精品| 美女国内精品自产拍在线播放| 久久精品电影| 青青草一区二区免费精品| 毛片a级毛片免费观看免下载| 久久综合色天堂av| 婷婷五月在线| 麻豆精品国产自产在线| 国产成人亚洲精品色欲AV| 亚洲男人天堂2020| 婷婷丁香在线观看| 国产导航在线| 亚洲欧美日韩成人在线| 国产精品亚洲一区二区三区z | 91福利免费| 日韩午夜福利在线观看| 亚洲第一黄色网址| 波多野结衣第一页| 一级毛片免费观看不卡视频| 色天天综合| 亚洲成av人无码综合在线观看| 免费激情网址| 激情六月丁香婷婷四房播| 97成人在线视频| 精品少妇人妻av无码久久| 亚洲va在线∨a天堂va欧美va| 99精品视频播放| 99福利视频导航| 黄色片中文字幕| 欧美精品成人一区二区视频一| 国产黄网站在线观看| 最新国语自产精品视频在| 有专无码视频| 在线中文字幕网|