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

    基于水循環(huán)的分布式水資源調(diào)配模型開發(fā)與應(yīng)用

    2022-05-19 08:31:04楊明智許繼軍桑連海
    水利學(xué)報(bào) 2022年4期
    關(guān)鍵詞:水循環(huán)水文水源

    楊明智,許繼軍,桑連海,劉 強(qiáng)

    (1.流域水資源與生態(tài)環(huán)境科學(xué)湖北省重點(diǎn)實(shí)驗(yàn)室 長江水利委員會(huì)長江科學(xué)院,湖北 武漢 430010;2.長江水利委員會(huì) 長江經(jīng)濟(jì)帶保護(hù)與發(fā)展戰(zhàn)略研究中心,湖北 武漢 430010)

    1 研究背景

    長期以來,國內(nèi)外對水資源配置的實(shí)踐和研究一般是通過還原實(shí)測徑流、評估天然水資源量的方式進(jìn)行水量分配,而隨著社會(huì)取用水量的加大,需要還原的水資源量占比逐漸加大,還原法所產(chǎn)生的誤差已經(jīng)難以被忽略[1-2]。并且,傳統(tǒng)的水資源配置一般采用集總式配置方法。在空間上,忽略計(jì)算分區(qū)內(nèi)自然地理?xiàng)l件和經(jīng)濟(jì)社會(huì)參數(shù)的空間差異性,均化了分區(qū)的下墊面和人類活動(dòng)強(qiáng)度空間特性,不利于基層單元的水資源管理[3];在時(shí)間上,一般采用逐月配置的模式,不能實(shí)時(shí)響應(yīng)水資源的動(dòng)態(tài)變化對水資源配置的影響[4]。

    隨著分布式水文模型的發(fā)展日漸成熟,一些學(xué)者提出了以水循環(huán)為基礎(chǔ)的水資源動(dòng)態(tài)配置模式。趙勇等[5]將平原區(qū)分布式水循環(huán)模型(PDWCM)與廣義水資源合理配置模型(GWRAM)耦合,采用分解和聚合的方式進(jìn)行信息交互,實(shí)現(xiàn)了“自然-社會(huì)”水循環(huán)的全過程模擬。張洪剛等[6]利用NAM模型模擬流域水文過程,并將徑流轉(zhuǎn)化到MIKEBASIN 模型所劃分的各分區(qū)單元,實(shí)現(xiàn)了集總式水文模型與水資源配置模型的耦合,并對漢江流域進(jìn)行水資源配置研究。陳強(qiáng)等[7]通過改進(jìn)SWAT 模型,建立水文響應(yīng)單元與配置分區(qū)的空間拓?fù)潢P(guān)系,對水資源配置結(jié)果進(jìn)行時(shí)空展布,實(shí)現(xiàn)了SWAT 模型與水資源配置模型的松散耦合,計(jì)算了天津13 個(gè)區(qū)(縣)的農(nóng)業(yè)灌溉水量。張守平[8]通過時(shí)空展布方式,將以水資源三級區(qū)套地市為單位的配置分區(qū)數(shù)據(jù)傳輸?shù)揭宰恿饔蚋叱處澐值挠?jì)算單元上,實(shí)現(xiàn)了集總式水資源配置模型與分布式水循環(huán)模型WEP-L 的耦合,并進(jìn)行了渭河流域水量水質(zhì)聯(lián)合配置研究。呂偉[9]將水資源配置信息進(jìn)行尺度轉(zhuǎn)換,實(shí)現(xiàn)與MIKE SHE 灌溉模塊的信息交互,構(gòu)建了灌區(qū)農(nóng)業(yè)水資源優(yōu)化調(diào)配耦合模型,并在人民勝利渠灌區(qū)進(jìn)行了應(yīng)用。曾思棟等[10]針對傳統(tǒng)水資源配置模型在水循環(huán)聯(lián)系方面的不足,通過插件式架構(gòu)實(shí)現(xiàn)分布式水文模型DTVGM 與水資源配置的信息交互,開發(fā)了分布式水資源配置模型DTVGM-WEAR,并在宜昌市主城區(qū)及東風(fēng)灌區(qū)進(jìn)行了應(yīng)用。桑學(xué)鋒等[11-12]通過水文模塊與水資源調(diào)配模塊緊密耦合,構(gòu)建了概念性分布式水資源綜合模擬與調(diào)配模型WAS,實(shí)現(xiàn)了自然-社會(huì)水循環(huán)過程的實(shí)時(shí)互饋模擬。

    水資源配置模型與分布式水文模型的耦合研究至今已有十幾年時(shí)間,模型的發(fā)展也比較迅速。但是,其耦合方式多以松散耦合為主,即將水資源配置信息和水文模擬信息進(jìn)行尺度轉(zhuǎn)換,采用公用參數(shù)或者輸出文件單向傳遞的方式實(shí)現(xiàn)集總式水資源配置模型與分布式水文模型的耦合[3]。這種耦合方式雖然充分考慮了區(qū)域的空間異質(zhì)性和水循環(huán)過程,但主要是進(jìn)行信息的單向傳遞,忽略了自然水循環(huán)與社會(huì)水循環(huán)之間的動(dòng)態(tài)互饋?zhàn)饔?。以降?產(chǎn)流-入滲-蒸發(fā)為代表的自然水循環(huán)與以取水-用水-耗水-排水為代表的社會(huì)側(cè)枝水循環(huán)兩大過程相互聯(lián)系、密不可分,前一時(shí)段的水資源開發(fā)利用活動(dòng)直接影響著下一個(gè)時(shí)段的水資源數(shù)量及其空間分布[11]?;谒缮Ⅰ詈系乃Y源配置方法割裂了兩大水循環(huán)過程的聯(lián)系,不能很好地反映人類活動(dòng)用水過程與天然水文過程之間的相互影響,容易導(dǎo)致水資源配置結(jié)果出現(xiàn)偏差。在水文模型的結(jié)構(gòu)設(shè)計(jì)中嵌入水資源配置模塊,也即是二者的雙向耦合,是實(shí)現(xiàn)自然-人工水資源復(fù)雜系統(tǒng)精準(zhǔn)模擬的有效手段。盡管如此,當(dāng)前與分布式水文模型耦合的研究仍然較少,而且大多屬于半松散耦合的方式,只有少部分實(shí)現(xiàn)了雙向耦合,如WAS 模型。

    因此,本文在分布式水文模型SWAT 的基礎(chǔ)上,通過改進(jìn)其內(nèi)置模塊,并將之與嵌入的水資源配置模塊相銜接,實(shí)現(xiàn)二者的雙向耦合,構(gòu)建基于SWAT 的分布式水資源調(diào)配模型(Water resources allocation and regulation model based on the SWAT, SWAT-WARM)。模型采用具有多元屬性的嵌套式單元?jiǎng)澐址椒?,便于水文模擬單元與水資源配置模擬單元之間的實(shí)時(shí)信息交互,實(shí)現(xiàn)了水文模型與水資源配置模型的雙向耦合,使SWAT-WARM 模型具有“自然-社會(huì)”二元水循環(huán)與水資源調(diào)配一體化模擬功能,可以在模擬過程中時(shí)刻保持自然水循環(huán)與社會(huì)側(cè)枝水循環(huán)之間的動(dòng)態(tài)耦合關(guān)系,詳細(xì)描述自然-社會(huì)水資源復(fù)雜互饋系統(tǒng),既可以從宏觀上描述流域整體水循環(huán)狀況和水資源配置結(jié)果,也可從微觀上探究日尺度模擬單元的取用水結(jié)構(gòu)。

    2 模型開發(fā)

    2.1 模型結(jié)構(gòu)SWAT-WARM 模型是采用FORTRAN 語言在Microsoft Visual Studio2012 平臺(tái)下開發(fā)完成的。模型由分布式水文模型SWAT[13]和水資源配置模型兩大子模型聚合而成,其中水資源配置模型是為了便于與SWAT 模型結(jié)合而開發(fā)。通過改進(jìn)SWAT 模型中的相關(guān)模塊,將水資源配置子模型嵌入SWAT 模型,實(shí)現(xiàn)雙方模擬數(shù)據(jù)之間的實(shí)時(shí)交互。

    水資源配置子模型主要由數(shù)據(jù)空間展布模塊、多源供水模塊、配置模塊、耗水模塊和排水模塊構(gòu)成,以SWAT 模型劃分的水文模擬單元作為基本配置單元,依托SWAT 模型驅(qū)動(dòng)運(yùn)行。在模擬開始時(shí),數(shù)據(jù)空間展布模塊將輸入的需水?dāng)?shù)據(jù)(或歷史用水?dāng)?shù)據(jù))展布到計(jì)算單元上;多源供水模塊根據(jù)水源類型、取水規(guī)則等信息和由水文模型提供的水資源信息進(jìn)行逐日的取水計(jì)算;配置模塊根據(jù)每個(gè)計(jì)算單元取水量、行業(yè)用水需求和用水規(guī)則信息進(jìn)行逐日的水量分配,通過耗水和排水模塊輸出逐日取水、排水等信息,傳遞給水文模型;

    水文子模型(SWAT)刻畫水文循環(huán)過程,實(shí)時(shí)模擬水循環(huán)變化對社會(huì)用水的影響,以及模擬社會(huì)用水對下一階段的水資源、供用水變化的影響,為水資源配置模型提供河道徑流量、湖庫蓄水量、地下含水層儲(chǔ)水量等實(shí)時(shí)水資源邊界條件信息。

    兩大子模型互為一體、緊密聯(lián)系,又分擔(dān)各自的功能,通過輸入、輸出數(shù)據(jù)實(shí)時(shí)互饋,構(gòu)成了一個(gè)完整的“自然-社會(huì)”水資源模擬系統(tǒng),實(shí)現(xiàn)二元水循環(huán)和水資源配置一體化模擬。模型結(jié)構(gòu)框架見圖1 所示。

    圖1 SWAT-WARM 模型結(jié)構(gòu)框架

    2.2 耦合方法水資源配置子模型與水文子模型的信息雙向交互方式如下:在運(yùn)行開始時(shí),(1)水資源配置模型的多源供水模塊將每個(gè)計(jì)算單元的需水信息和取水優(yōu)先序信息傳遞給水文模型,調(diào)用各水源模塊進(jìn)行取水量計(jì)算;(2)取水?dāng)?shù)據(jù)反向傳遞給水資源配置模型中的配置模塊,對每一個(gè)計(jì)算單元內(nèi)的各類用水部門進(jìn)行水量分配,之后將數(shù)據(jù)依次傳遞給耗水模塊和排水模塊,計(jì)算耗水量和排水量;(3)將排水?dāng)?shù)據(jù)返回水文模型的點(diǎn)面源模塊,進(jìn)行點(diǎn)源和面源排水計(jì)算;(4)重復(fù)(1)—(3)步驟,進(jìn)行逐日計(jì)算,得到整個(gè)模擬時(shí)段的區(qū)域水循環(huán)動(dòng)態(tài)變化和水資源實(shí)時(shí)配置模擬結(jié)果。通過模塊銜接,SWAT-WARM 模型可以在模擬過程中實(shí)現(xiàn)水文模型與水資源配置模型之間逐日信息的雙向傳遞。模塊之間的數(shù)據(jù)傳遞關(guān)系見圖2 所示。

    圖2 模塊的數(shù)據(jù)傳遞關(guān)系

    為實(shí)現(xiàn)以上所描述的信息雙向交互,SWAT 模型與水資源配置模型的耦合設(shè)計(jì)包括計(jì)算單元?jiǎng)澐帧⒛K修改與銜接兩大部分內(nèi)容。

    2.2.1 計(jì)算單元?jiǎng)澐?水文模型以自然分水嶺線的延伸方向劃分天然子流域,但是灌區(qū)邊界和行政區(qū)邊界一般不與自然地貌重合,導(dǎo)致所劃分的計(jì)算單元在空間位置上與行政區(qū)、灌區(qū)難以很好地對應(yīng),影響模擬結(jié)果的準(zhǔn)確性[14]。本研究采用基于“子流域-行政區(qū)-灌區(qū)-土地利用-土壤”多元屬性的嵌套式坡面離散方法進(jìn)行計(jì)算單元?jiǎng)澐郑▓D3(a))。首先,提取河網(wǎng)水系,劃分天然子流域;然后,疊加土地利用信息、土壤類型信息和坡度信息,劃分水文響應(yīng)單元(HRU);最后,將HRU 依次疊加行政區(qū)邊界和灌區(qū)邊界,將其進(jìn)一步切割成多個(gè)水文模擬單元,最終完成計(jì)算單元?jiǎng)澐?。每一個(gè)計(jì)算單元都具有子流域?qū)傩浴⑿姓^(qū)屬性、灌區(qū)屬性、土地利用類型屬性和土壤類型屬性。

    通過計(jì)算單元?jiǎng)澐郑顾哪M單元與基本配置單元一一對應(yīng),從而建立水資源配置模型計(jì)算單元與水文模型計(jì)算單元之間的空間拓?fù)潢P(guān)系。每一個(gè)計(jì)算單元既是水文模型的水文模擬單元,又是水資源配置模型的基本配置單元。計(jì)算單元內(nèi)的用水部門類型劃分視水文模擬單元的土地利用類型而定:如果一個(gè)水文模擬單元的土地利用類型是城鎮(zhèn)用地,則對應(yīng)的用水部門包括城鎮(zhèn)居民、工業(yè)和服務(wù)業(yè)等;如果是農(nóng)村用地,則對應(yīng)的用水部門包括農(nóng)村居民和畜禽等;如果是農(nóng)業(yè)用地,則對應(yīng)的用水部門包括農(nóng)田、林地和草地等(圖3(b))。

    圖3 計(jì)算單元?jiǎng)澐?/p>

    這樣的劃分方式,既能體現(xiàn)分布式水文模型的單元?jiǎng)澐痔攸c(diǎn),也滿足了區(qū)域水資源管理的需求,為兩個(gè)子模型模擬信息的雙向傳遞奠定基礎(chǔ),進(jìn)而實(shí)現(xiàn)不同尺度的水資源配置模擬。在系統(tǒng)運(yùn)行時(shí),水資源配置模型可以根據(jù)計(jì)算單元的類型識別其各用水部門的取水源及其取水優(yōu)先序,進(jìn)而傳遞和調(diào)用水文模型中的各類水源模塊,計(jì)算各個(gè)單元的行業(yè)用水量和水源取水量,從而進(jìn)一步統(tǒng)計(jì)各行政分區(qū)和整個(gè)研究區(qū)的水資源配置結(jié)果。

    2.2.2 模塊修改與銜接 為實(shí)現(xiàn)水資源配置模型與SWAT 模型的耦合,需要對SWAT 模型的部分模塊進(jìn)行代碼修改,并與水資源配置模型中的相關(guān)模塊進(jìn)行銜接。

    (1)耗用水模塊修改。SWAT 模型中的耗用水模塊包括rchuse 模塊(河道取水)、res 模塊(水庫取水)和watuse 模塊(淺水層取水、深水層取水和坑塘取水),但耗用水的設(shè)計(jì)僅能反映其逐月變化特征。因此,對三個(gè)耗用水模塊進(jìn)行代碼修改,將月取水參數(shù)修改為日取水參數(shù),實(shí)現(xiàn)逐日用水模擬。此外,在watuse 模塊中添加外調(diào)水取水和再生水處理與使用程序,實(shí)現(xiàn)外調(diào)水和再生水在生活、工業(yè)和服務(wù)業(yè)等用水部門中的使用。

    (2)灌溉水源模塊修改。SWAT 模型中的灌溉水源模塊包括irr_rch 模塊(河道灌溉)、irr_res 模塊(水庫灌溉)和irrsub 模塊(淺層地下水灌溉、深層地下水灌溉和外調(diào)水灌溉)。在irrsub 模塊中添加坑塘取水程序,實(shí)現(xiàn)坑塘水在農(nóng)業(yè)灌溉供水方面的使用。SWAT 模型采用將水源取水量直接添加到HRU 上的方式進(jìn)行灌溉模擬,渠系蒸發(fā)、渠系滲漏和渠系退水都作為系統(tǒng)損失,這會(huì)影響水文模擬的準(zhǔn)確性。因此,對三個(gè)灌溉水源模塊進(jìn)行修改,相應(yīng)的計(jì)算公式如下:

    式中:ETcan為渠系蒸發(fā)損失量,mm;Irrcan為進(jìn)入渠道的灌溉水量,mm;Lscan為渠系滲漏損失量,mm;Surpcan為渠系退水量,mm;φ為渠系水有效利用系數(shù);α為渠系蒸發(fā)系數(shù);β為渠系滲漏系數(shù)。

    (3)水源模塊的銜接。多源供水模塊首先識別計(jì)算單元的土地利用類型。若是居工地,則程序進(jìn)入城鄉(xiāng)供水計(jì)算流程;若是農(nóng)田用地,則程序進(jìn)入灌溉供水計(jì)算流程;若是其他土地利用類型,則該計(jì)算單元不進(jìn)行供用水模擬。然后,根據(jù)設(shè)定的取水優(yōu)先序規(guī)則識別計(jì)算單元的水源識別碼和取水優(yōu)先序代碼,依次調(diào)用各水源模塊(圖4)。河道、水庫、淺層地下水、深層地下水、外調(diào)水、坑塘水和再生水的水源識別碼分別是1、2、3、4、5、6、7。

    圖4 多源供水模塊中的水源模塊調(diào)用流程

    (4)地下水模塊修改。對gwmod 模塊相關(guān)代碼進(jìn)行修改,描述渠系滲漏水量、城鎮(zhèn)供水管網(wǎng)漏失水量補(bǔ)給潛水這一過程。修改公式如下:

    式中:rht為第t天的潛水補(bǔ)給量,mm;rht-1為第t-1 天的潛水補(bǔ)給量,mm;prc為補(bǔ)給地下水的土壤滲漏水量,mm;Lscan為渠系漏損量,mm;Lspip為城鎮(zhèn)供水管網(wǎng)漏失量,mm;delay為地下水補(bǔ)給延遲系數(shù),mm;WSP為水源取水量,m3;pip為城鎮(zhèn)供水管網(wǎng)漏損率;Area為計(jì)算單元的面積,m2。

    (5)地表產(chǎn)流模塊修改。對surface 模塊進(jìn)行修改,在地表產(chǎn)流計(jì)算公式中添加渠道退水參數(shù)。這樣,渠道退水隨地表徑流進(jìn)入河道,參與后續(xù)的河道匯流過程模擬。

    (6)點(diǎn)源、面源模塊銜接。為實(shí)現(xiàn)水資源配置模型的排水模塊與SWAT 模型的點(diǎn)源模塊(recday 模塊和recmon 模塊)的銜接,將點(diǎn)源模塊中的排水參數(shù)置換為排水模塊計(jì)算出的排水參數(shù),采用錯(cuò)日傳遞法[14]將排水?dāng)?shù)據(jù)傳遞給點(diǎn)源模塊,將產(chǎn)生的污水安排在下一日模擬時(shí)排放,進(jìn)而實(shí)現(xiàn)點(diǎn)源的逐日銜接模擬。

    式中:PWy,i和WPy,i分別為第y年第i日的點(diǎn)源排放量和生活工業(yè)污水排放量,m3;dlt為一年的天數(shù),dlt= 365 或366。

    城鎮(zhèn)生活污水和工業(yè)廢水產(chǎn)生量由排水模塊模擬計(jì)算獲得,計(jì)算公式如下:

    式中:Wrd為城鎮(zhèn)生活污水/工業(yè)廢水產(chǎn)生量,m3,是在耗水模塊中采用行業(yè)實(shí)際用水量乘以耗水系數(shù)計(jì)算獲得;vd為廢污水處理率;re為再生水利用率。

    農(nóng)村居民生活污水和畜禽養(yǎng)殖污水排放量計(jì)入其所在計(jì)算單元的地表產(chǎn)流計(jì)算中,相關(guān)參數(shù)在subbasin 模塊中修改。

    2.3 配置規(guī)則考慮到水資源系統(tǒng)的復(fù)雜性,采用目標(biāo)函數(shù)和約束條件構(gòu)建的優(yōu)化配置模型容易出現(xiàn)求解困難或優(yōu)化結(jié)果偏離實(shí)際的情況[15]。SWAT-WARM 模型采用基于規(guī)則的水資源配置方法,通過對水資源系統(tǒng)各過程的模擬與控制,在設(shè)定的規(guī)則下進(jìn)行水量分配,比優(yōu)化模型具有更強(qiáng)的適應(yīng)性和仿真性,可以針對水資源復(fù)雜系統(tǒng)有效地進(jìn)行多水源向多用戶的水量分配[16]。所設(shè)定的配置規(guī)則主要包括水量分配規(guī)則、供水優(yōu)先序規(guī)則、用水優(yōu)先序規(guī)則、水量平衡規(guī)則、水源利用規(guī)則、生態(tài)保護(hù)規(guī)則等。在配置求解方法上,模型以水源可利用水量、工程供水能力、取水總量限制等為約束條件,按照設(shè)定的水源優(yōu)先順序依次計(jì)算各水源供水量,以計(jì)算單元的編號順序,按照設(shè)定的用水優(yōu)先級進(jìn)行逐日循環(huán)計(jì)算,盡可能地滿足受水單元。

    2.3.1 水量分配規(guī)則 水量分配計(jì)算以最大程度地滿足計(jì)算單元的用水需求為原則,計(jì)算公式如下

    式中:TWf為計(jì)算單元的缺水量;I為計(jì)算單元內(nèi)的用水部門個(gè)數(shù);WDi為第i個(gè)用水部門的需水量,m3,由輸入的總需水量數(shù)據(jù)展布得來;WUi為第i個(gè)用水部門的實(shí)際用水量,m3。

    2.3.2 用水優(yōu)先序規(guī)則 計(jì)算單元內(nèi)的城鄉(xiāng)生活、工業(yè)、服務(wù)業(yè)、農(nóng)業(yè)灌溉等用水部門對水量的分配存在優(yōu)先順序要求。SWAT-WARM 模型允許每一個(gè)計(jì)算單元根據(jù)自身的管理需求,對單元內(nèi)的用水部門設(shè)置或調(diào)整用水次序方案。程序根據(jù)計(jì)算單元當(dāng)日的取水量依次分配給各用水部門,直至分配完畢為止。用水優(yōu)先計(jì)算流程見圖5(a)所示。行業(yè)用水計(jì)算表達(dá)式如下:

    式中:WDi和WUi分別為計(jì)算單元內(nèi)用水優(yōu)先序?yàn)閕的用水部門當(dāng)日的需水量和用水量,m3;I為用水部門個(gè)數(shù);WSP為計(jì)算單元當(dāng)日的取水量,m3;pip為供水管網(wǎng)漏損率(城鎮(zhèn)生活、工業(yè)等)或灌溉水利用率(農(nóng)業(yè)灌溉)。

    2.3.3 供水優(yōu)先序規(guī)則 計(jì)算單元對設(shè)定的外調(diào)水、水庫水、河道水、坑塘水、地下水、再生水等多種水源存在供水優(yōu)先順序要求。SWAT-WARM 模型允許每一個(gè)計(jì)算單元根據(jù)自身的水源偏好,設(shè)定或調(diào)整水源供水次序方案。程序從指定的水源進(jìn)行順序取水計(jì)算,直至滿足該計(jì)算單元的當(dāng)日需水量或者到最后一個(gè)供水源供水完畢為止。供水優(yōu)先計(jì)算的流程如圖5(b)所示。水源取水計(jì)算表達(dá)式如下:

    圖5 用水、供水優(yōu)先計(jì)算流程圖

    式中:WSPi為計(jì)算單元中供水優(yōu)先序?yàn)閕的水源當(dāng)日的取水量,m3;I為供水源個(gè)數(shù);WD為計(jì)算單元當(dāng)日的需水量,m3;WF為取水工程的供水能力,m3;Wsc為水源當(dāng)日的可利用水量,m3;WMX為區(qū)域取水限制總量,m3。

    2.3.4 水量平衡規(guī)則 包括河流水量平衡、湖庫水量平衡、土壤水量平衡、地下水量平衡等內(nèi)容。這些在SWAT 模型內(nèi)置的相關(guān)模塊中都有體現(xiàn)。

    2.3.5 取水總量控制規(guī)則 模型中以設(shè)置的區(qū)域地表水利用總量和地下水開采總量為取水總量限制條件,通過分階段配置的方式實(shí)現(xiàn)地表、地下水資源的聯(lián)合調(diào)配。取水總量可以在年內(nèi)滾動(dòng)提取,不能跨年度累積。

    式中:I、J分別為年總天數(shù)和計(jì)算單元總個(gè)數(shù);Wsij為第i天、計(jì)算單元j的地表水利用量,m3;WSM為區(qū)域年度地表水允許利用量,m3;Wgij為第i天、計(jì)算單元j的地下水開采量,m3;WGM為區(qū)域年度地下水允許開采量,m3。

    2.3.6 生態(tài)保護(hù)規(guī)則 為了維持水體自凈能力、保護(hù)水生生態(tài)系統(tǒng)等要求,需要滿足河道內(nèi)的生態(tài)流量需求。對此,模型允許設(shè)置生態(tài)流量控制閾值,當(dāng)河道流量低于該閾值時(shí),將不能在河道取水,以保護(hù)河道生態(tài)流量。

    式中:Qeco為河道生態(tài)流量控制閾值,m3/s;I為從河道取水的用水戶個(gè)數(shù);WSPi為用水戶i的河道取水量,m3;Qrch為河道流量,m3/s。

    通過以上的設(shè)計(jì),SWAT-WARM 模型具有了“自然-社會(huì)”二元水循環(huán)與水資源實(shí)時(shí)配置一體化模擬功能。

    3 案例應(yīng)用

    3.1 研究區(qū)概況唐河流域地處漢江流域中部,是唐白河流域的一條支流流域,北起伏牛山地區(qū),西與白河為鄰,東靠桐柏山脈,與淮河隔山相望,南抵翟灣兩河口,流域總面積約為8473 km2(圖6)。地勢由東北向西南傾斜,唐河水系主要由干流唐河及其分支毗河、泌陽河、桐河、澗河等,涉及社旗縣、沁陽縣、桐柏縣、方城縣、宛城區(qū)、臥龍區(qū)、新野縣、襄州區(qū)和棗陽市等地。唐河流域境內(nèi)河道縱橫,渠系密布,南水北調(diào)中線渠在流域北部的方城縣境內(nèi)穿過,供水范圍涉及方城縣、社旗縣和唐河縣三縣。唐河流域氣候溫和、雨量充沛、地形平坦、土壤肥沃,盛產(chǎn)小麥、水稻、花生等糧油作物。但是,流域降水變率大,具有基流小、洪峰高的特點(diǎn),是漢江徑流的低值區(qū)[17],年內(nèi)分配極為不均,地表水資源年際變化大、豐枯非常懸殊,常發(fā)生洪、旱災(zāi)害。為防洪抗旱、保障農(nóng)業(yè)灌溉用水,流域內(nèi)建有宋家場水庫、虎山水庫等水庫。

    圖6 研究區(qū)地理位置

    3.2 模型構(gòu)建

    3.2.1 計(jì)算單元?jiǎng)澐?從SRTM 數(shù)據(jù)下載平臺(tái)(http://srtm.csi.cgiar.org/)下載研究區(qū)DEM 數(shù)據(jù);從資源環(huán)境科學(xué)與數(shù)據(jù)中心下載土地利用數(shù)據(jù)和土壤類型數(shù)據(jù),并通過查詢中國土壤數(shù)據(jù)庫獲取土壤相關(guān)參數(shù)?;谝陨线b感數(shù)據(jù)劃分天然子流域和HRU,再依次疊加行政區(qū)劃信息及灌區(qū)分布信息,將HRU 進(jìn)一步細(xì)分,共獲得19 個(gè)子流域和760 個(gè)計(jì)算單元(圖7)。這樣,每個(gè)計(jì)算單元都具有子流域、行政區(qū)、灌區(qū)、土壤類型和土地利用類型屬性,從而構(gòu)建水文模擬單元與基本配置單元之間的一一對應(yīng)關(guān)系。

    圖7 計(jì)算單元?jiǎng)澐纸Y(jié)果

    3.2.2 數(shù)據(jù)輸入 本研究的模型構(gòu)建需要輸入的數(shù)據(jù)主要有:

    (1)氣象觀測數(shù)據(jù):包括唐河流域內(nèi)及鄰近的4 個(gè)氣象站的1981—2016年氣象數(shù)據(jù),包括日降水量、日最高及最低氣溫、相對濕度、日照時(shí)數(shù)風(fēng)速等氣象要素;

    (2)經(jīng)濟(jì)社會(huì)用水?dāng)?shù)據(jù):通過查詢2006—2016年南陽市水資源公報(bào)、襄陽市水資源公報(bào)等資料,獲得歷年各行政區(qū)行業(yè)用水?dāng)?shù)據(jù)和地表水、地下水資源利用量數(shù)據(jù);

    (3)水利工程信息:主要包括水庫庫容、灌溉渠道過水能力、電機(jī)井日提水能力等數(shù)據(jù);

    (4)耗用水效率信息:主要包括渠系水有效利用系數(shù)、行業(yè)耗水率、污水處理率等內(nèi)容。考慮到唐河流域的再生水利用量很少,不到用水總量的0.1%,故本次不考慮流域再生水利用情況,不設(shè)置再生水利用率參數(shù)。

    3.2.3 配置規(guī)則設(shè)置 根據(jù)唐河流域供用水實(shí)際情況及水庫、渠道等水利工程空間分布情況,明晰輸配水去向,采用導(dǎo)入信息表的方式設(shè)置取用水優(yōu)先序規(guī)則,主要包括供水源類型、水源供水優(yōu)先次序方案、用水部門類型和行業(yè)用水優(yōu)先次序方案等規(guī)則信息。由此構(gòu)建用水部門與取水源之間的空間拓?fù)潢P(guān)系。

    3.3 模型校驗(yàn)本文參照文獻(xiàn)[12]的多指標(biāo)校驗(yàn)方法,按照先水文校驗(yàn)、后水量校驗(yàn)的次序,分別對河道徑流、特征頻率年的徑流總量、國民經(jīng)濟(jì)行業(yè)用水量和水源取水量進(jìn)行校驗(yàn)。首先調(diào)整徑流過程,然后對校驗(yàn)?zāi)陱搅骺偭?,確定水文參數(shù);最后,對行業(yè)用水量和水源取水量進(jìn)行校驗(yàn),確定社會(huì)水循環(huán)參數(shù)。在評價(jià)指標(biāo)選擇方面,采用相關(guān)系數(shù)R2和納什效率系數(shù)Ens兩個(gè)性能指標(biāo)評價(jià)河道斷面徑流模擬效果,采用偏差百分比指標(biāo)PBIAS 評價(jià)年徑流總量、行業(yè)用水量、地表和地下水資源利用量的模擬效果。一般認(rèn)為,當(dāng)PBIAS 值在±25%范圍內(nèi)、R2>0.60、Ens>0.50 時(shí),模擬精度是令人滿意的[20-21]。

    3.3.1 徑流校驗(yàn) 通過對徑流曲線數(shù)、蒸發(fā)補(bǔ)償系數(shù)、河床水力傳導(dǎo)度等主要水文參數(shù)的調(diào)整,郭灘水文站的模擬徑流過程與實(shí)測徑流過程對比結(jié)果見表1和圖8 所示。在率定期(1983—2002年),郭灘站月徑流模擬值與實(shí)測值的相關(guān)系數(shù)為0.741,納什效率系數(shù)為0.734。在驗(yàn)證期(2003—2016年),郭灘站月徑流模擬值與實(shí)測值的相關(guān)系數(shù)為0.724,納什系數(shù)為0.685。可以看出,模型對水文過程的模擬精度達(dá)到了要求值。

    圖8 郭灘站的實(shí)測與模擬月徑流過程對比

    表1 郭灘水文站的月徑流模擬結(jié)果

    3.3.2 徑流總量校驗(yàn) 通過1983—2016年長系列模擬,獲得郭灘水文站多年平均、25%、50%、75%和90%水平年的年徑流量結(jié)果,與實(shí)測年徑流量的偏差百分比結(jié)果見表2。由表可知,郭灘站多年平均、25%、50%、75% 和90%頻率年的偏差百分比分別為7.82%、23.96%、2.17%、-19.34 和-7.58%。可以看出,模型對徑流總量的模擬精度達(dá)到了要求值。

    表2 不同水平年郭灘水文站年徑流總量的偏差百分比

    3.3.3 供用水校驗(yàn) 由于水資源公報(bào)中的供用水?dāng)?shù)據(jù)是以行政區(qū)為單位進(jìn)行統(tǒng)計(jì)的,本次基于《南陽市水資源公報(bào)》(2010—2016),選擇行政范圍基本全部位于唐河流域境內(nèi)的唐河縣和社旗縣兩個(gè)縣的歷年實(shí)際用水?dāng)?shù)據(jù)和水資源開發(fā)利用量數(shù)據(jù)進(jìn)行取用水量校驗(yàn)。通過對渠系水有效利用系數(shù)、渠系輸水蒸發(fā)率、生活工業(yè)耗水率等主要社會(huì)水循環(huán)參數(shù)的調(diào)整,供用水量模擬結(jié)果的偏差百分比見表3所示,模擬結(jié)果對比見圖9 所示。

    表3 2010—2016年社旗縣和唐河縣供用水模擬結(jié)果的偏差百分比

    圖9 唐河縣和社旗縣2010—2016年供用水實(shí)際值與模擬值對比

    社旗縣的供用水量偏差百分比幾乎全部為0%,模擬效果很好。唐河縣用水總量偏差百分比為5.6%,地表水、地下水供水量偏差百分比分別為3.1%、6.6%,均在10%以內(nèi)。這是因?yàn)樘坪涌h經(jīng)濟(jì)社會(huì)用水以地下水為主,長系列模擬造成了模型中設(shè)置為供水源的子流域內(nèi)的地下水量采補(bǔ)失衡,使唐河縣地下水利用量模擬誤差相對偏高,也致使行業(yè)用水量的模擬值在整體上略小于實(shí)際值??傮w上,模型對取用水的模擬結(jié)果在可接受的范圍內(nèi),可以較好地反映實(shí)際情況。

    3.4 結(jié)果分析

    3.4.1 流域年度水資源配置結(jié)果分析 SWAT-WARM 模型可以統(tǒng)計(jì)、分析并輸出流域年度水資源配置結(jié)果。以75%特征頻率年為例,唐河流域的行業(yè)用水和水源供水結(jié)果見表4所示。唐河流域用水總量為9.19 億m3,其中生活用水(含生態(tài)用水)、工業(yè)用水和農(nóng)業(yè)用水分別為1.18 億m3、1.18 億m3和6.83 億m3。其中,南陽市用水量最大(5.48 億m3),其次為駐馬店市(1.96 億m3),襄陽市用水量最?。?.76 億m3)。在各供水水源中,地下水供水量最高(4.91 億m3),占總供水量的53.37%,其次為水庫水(3.64 億m3),占總供水量的39.65%,外流域供水量最?。s0.02 億m3),僅占總供水量的0.2%。唐河流域內(nèi)南陽市的地下水利用量最高,達(dá)3.73 億m3,占地下水用水總量的76.04%;其次是駐馬店市(1.03 億m3),占地下水用水總量的20.91%;襄陽市地下水利用量最小(0.15億m3),僅占用水總量的3.05%。

    表4 唐河流域75%頻率年的水資源配置結(jié)果

    3.4.2 流域年度水循環(huán)轉(zhuǎn)化分析 以75%頻率年為例,通過匯總唐河流域的水循環(huán)模擬結(jié)果,獲得流域水循環(huán)轉(zhuǎn)化關(guān)系見圖10 所示。全流域降水總量65.04 億m3,地表產(chǎn)流量9.86 億m3,總徑流量11.13 億m3,年末土壤水蓄變量-10.31 億m3、地表水蓄變量1.78 億m3、地下水蓄變量為-4.34 億m3。全流域耗水總量為66.62 億m3,其中土壤蒸發(fā)18.56 億m3,植被散發(fā)44.20 億m3,截留蒸發(fā)2.10 億m3,積雪升華0.37 億m3,水面蒸發(fā)0.47 億m3,生活工業(yè)用水消耗0.92 億m3。全流域的經(jīng)濟(jì)社會(huì)總用水量為9.19 億m3,地表水利用量4.29 億m3(含外調(diào)水0.02 億m3),地下水利用量4.91 億m3,經(jīng)濟(jì)社會(huì)耗水量和排水量分別為5.74 億m3和3.45 億m3,流域出境總水量8.17 億m3。

    圖10 唐河流域75%頻率年的自然-社會(huì)水循環(huán)轉(zhuǎn)化關(guān)系圖(單位:億m3)

    3.4.3 行政區(qū)月供需平衡分析 SWAT-WARM 模型可以統(tǒng)計(jì)、分析并輸出每個(gè)行政單元的逐月水量平衡分析結(jié)果。以社旗縣為例, 75%頻率年的水量供需平衡分析結(jié)果見表5所示。由表可知,75%頻率年社旗縣經(jīng)濟(jì)社會(huì)需水量為12 305 萬m3,供水量為11 178 萬m3,全年缺水1127 萬m3。缺水月份分別為3月、10—12月,其中,11月份缺水量最大(754 萬m3),其次是10月份(326 萬m3),3月份和12月份為輕微缺水(分別是36 萬m3和11 萬m3)。

    表5 2016年社旗縣在75%頻率年的水量供需平衡結(jié)果

    3.4.4 計(jì)算單元逐日取用水分析 SWAT-WARM 模型可以模擬并輸出每一個(gè)計(jì)算單元的逐日水量分配結(jié)果。本次選取了第365 號計(jì)算單元(農(nóng)田)和第757 號計(jì)算單元(城市)作為分析對象。兩個(gè)計(jì)算單元在75%頻率年都利用了河道、水庫和地下含水層三種水源。其中,第365 號計(jì)算單元的水源供水次序?yàn)椋核畮欤竞拥溃镜叵潞畬樱坏?57 號計(jì)算單元的水源供水次序?yàn)椋汉拥溃舅畮欤镜叵潞畬?。兩個(gè)計(jì)算單元在75%頻率年的行業(yè)用水和水源供水逐日模擬結(jié)果見圖11 所示。

    圖11 75%頻率年唐河流域第365 號和757 號計(jì)算單元的逐日水量分配結(jié)果圖

    由圖可知,SWAT-WARM 模型根據(jù)指定的配置規(guī)則對兩個(gè)計(jì)算單元進(jìn)行了多源互濟(jì)供水,滿足了計(jì)算單元內(nèi)各用水部門每日的用水需求。第365 號計(jì)算單元的水庫水利用量最大(324.33 萬m3),占總用水量的67.35%,河道水利用水量最?。?6.47 萬m3),占比僅為3.42%。在6—8月份,地下水成為了最重要的補(bǔ)充灌溉水源。這是因?yàn)樗畮煸谇皫讉€(gè)月的供水,蓄水量不足,而這段時(shí)間正是唐河流域的降雨集中期,較多的降雨有利于潛水補(bǔ)給。因此,水庫供水量減小,地下水利用量增加。第757號計(jì)算單元以河道取水為主,全年取水量占總量的74.53%(2241.43 萬m3)。在1月12—23日和6月25日—7月9日兩段時(shí)間因河道水量不足,進(jìn)行了水庫補(bǔ)充供水,全年供水量為139.62 萬m3,僅占總量的4.64%。在進(jìn)入10月份后,河道徑流量降低,難以繼續(xù)支撐該計(jì)算單元的用水需求。因此,在10月16日開始使用地下水進(jìn)行補(bǔ)充。

    4 結(jié)論與展望

    4.1 結(jié)論

    (1)本文通過修改SWAT 模型源代碼,設(shè)計(jì)并嵌合水資源配置模塊,構(gòu)建了基于水循環(huán)的分布式水資源調(diào)配模型SWAT-WARM,提出了包含計(jì)算單元?jiǎng)澐?、模塊修改與銜接以及配置規(guī)則的雙向耦合方法,實(shí)現(xiàn)了“自然-社會(huì)”水循環(huán)動(dòng)態(tài)互饋與水資源實(shí)時(shí)配置一體化模擬,具有了刻畫水資源系統(tǒng)多水源、多工程、多用戶的功能,能夠?qū)⑺磁c用水戶之間的拓?fù)潢P(guān)系、水量在水資源系統(tǒng)中的傳遞關(guān)系得到客觀、清晰的描述,不僅能夠反映自然水循環(huán)過程對社會(huì)取用水的影響,還能夠反映經(jīng)濟(jì)社會(huì)活動(dòng)對自然水循環(huán)過程的實(shí)時(shí)干預(yù),充分體現(xiàn)經(jīng)濟(jì)社會(huì)用水過程與天然水文過程之間的動(dòng)態(tài)反饋?zhàn)饔谩?/p>

    (2)模型采用嵌套式單元?jiǎng)澐址椒?,既能體現(xiàn)分布式水文模型的單元?jiǎng)澐痔攸c(diǎn),也滿足了區(qū)域水資源管理的需求,拓展了水資源配置的時(shí)空尺度:宏觀上,能夠提供區(qū)域的年尺度水循環(huán)轉(zhuǎn)化及水資源配置結(jié)果;中觀上,能夠提供行政分區(qū)的月尺度供需平衡分析;微觀上,能提供計(jì)算單元的日尺度行業(yè)用水過程和水源取水過程。模型采用基于規(guī)則的模擬方法,仿真性強(qiáng),可以通過設(shè)計(jì)或調(diào)整區(qū)域用水需求、水源取水限制條件及供、用水優(yōu)先順序等多種規(guī)則,制定出代表不同利益的配置方案,為區(qū)域水資源規(guī)劃與精細(xì)化管理提供良好的決策。

    (3)以唐河流域?yàn)檠芯繉ο螅捎冒ê拥罃嗝鎻搅鬟^程、特征頻率徑流總量、國民經(jīng)濟(jì)用水量和水資源開發(fā)利用量等多指標(biāo)校驗(yàn)方法對模型的模擬效果進(jìn)行校驗(yàn)。模型在郭灘站率定期(1983—2002年)的月徑流模擬的相關(guān)系數(shù)為0.741,納什效率系數(shù)為0.734,在驗(yàn)證期(2003—2016年)月徑流模擬的相關(guān)系數(shù)為0.724,納什效率系數(shù)為0.685;流域多年平均、25%、50%、75%和95%水平年的年徑流量模擬值與實(shí)際值的偏差百分比分別為7.82%、23.96%、2.17%、-19.34 和-7.58%;社旗縣的供用水量模擬值與實(shí)際值的偏差百分比均接近0%,唐河縣供用水量的偏差百分比均小于10%。以上結(jié)果均達(dá)到了模擬精度要求。在此基礎(chǔ)上,分別分析了年尺度全流域水資源配置結(jié)果、月尺度行政分區(qū)供需平衡結(jié)果和日尺度計(jì)算單元取用水過程模擬結(jié)果,并給出了全流域自然-社會(huì)水循環(huán)轉(zhuǎn)化關(guān)系圖。綜上可見,模型能夠反映流域水循環(huán)過程和水資源開發(fā)利用情況,可以支撐水資源精細(xì)化管理需求。

    4.2 展望

    (1)本研究所構(gòu)建的SWAT-WARM 模型是一種基于規(guī)則的水資源配置模擬模型,由于難以有效控制數(shù)量眾多的約束條件及參數(shù),一般采用經(jīng)驗(yàn)法組建多種情景,容易忽略系統(tǒng)最優(yōu)解,致使模擬方案的目標(biāo)針對性不強(qiáng)。如何將優(yōu)化算法與SWAT-WARM 模型進(jìn)行結(jié)合,綜合優(yōu)化模型和模擬模型的優(yōu)點(diǎn),還需要做進(jìn)一步的探討研究。

    (2)在模型中,計(jì)算單元需水量是按照計(jì)算單元所在行政區(qū)的面積占比進(jìn)行時(shí)空展布獲得的,忽略了不同城鎮(zhèn)之間人口分布、產(chǎn)業(yè)結(jié)構(gòu)的空間差異,也未考慮各用水部門對水資源的需求強(qiáng)度的時(shí)間差異,當(dāng)模擬面積較大區(qū)域時(shí),可能會(huì)影響模擬精度。此外,模型沒有設(shè)計(jì)水源對不同供水對象的分水比例,當(dāng)在有明確分水比例的地區(qū)進(jìn)行應(yīng)用時(shí),可能會(huì)造成一定的誤差。這些設(shè)計(jì)方法都有待進(jìn)一步的改善。

    (3)SWAT-WARM 模型內(nèi)嵌合了分布式水文模型(SWAT)的水質(zhì)相關(guān)模塊,為精細(xì)描述水質(zhì)的動(dòng)態(tài)變化過程提供了良好的基礎(chǔ)。但本次研究僅限于對水量的分配,沒有考慮水質(zhì)對水資源配置的影響,后續(xù)工作可以考慮在模擬規(guī)則中增加分質(zhì)供水條件,實(shí)現(xiàn)區(qū)域水量水質(zhì)聯(lián)合模擬與配置,進(jìn)一步完善模型的綜合性能與實(shí)用性能。

    猜你喜歡
    水循環(huán)水文水源
    用袋裝水模擬水循環(huán)
    2022年《中國水文年報(bào)》發(fā)布
    保護(hù)水源
    品牌研究(2022年18期)2022-06-29 05:33:24
    硫酸法鈦白粉水洗及水循環(huán)膜處理工藝
    化工管理(2021年7期)2021-05-13 00:46:24
    水文
    水文水資源管理
    南水源頭清如許
    水文
    尋找水源
    水循環(huán)高效礦井乏風(fēng)熱泵系統(tǒng)分析與應(yīng)用
    同煤科技(2015年4期)2015-08-21 12:51:02
    av免费在线看不卡| 亚洲怡红院男人天堂| 亚洲欧洲日产国产| 国产片特级美女逼逼视频| 99久久人妻综合| 欧美高清成人免费视频www| 亚洲av日韩在线播放| 成人国产麻豆网| 亚洲av电影在线观看一区二区三区 | 欧美高清性xxxxhd video| 亚洲欧美精品自产自拍| 久久久久久久亚洲中文字幕| 亚洲国产成人一精品久久久| 人妻系列 视频| 女人被狂操c到高潮| 欧美成人精品欧美一级黄| 国产伦精品一区二区三区四那| 亚洲中文字幕一区二区三区有码在线看| 日本色播在线视频| 国产高清国产精品国产三级 | 插逼视频在线观看| 在线观看一区二区三区| 天天躁日日操中文字幕| 国产乱人偷精品视频| 欧美激情久久久久久爽电影| 久久精品夜夜夜夜夜久久蜜豆| 搡老妇女老女人老熟妇| 欧美另类亚洲清纯唯美| av天堂中文字幕网| 色吧在线观看| 能在线免费看毛片的网站| 国产黄色视频一区二区在线观看 | 免费av不卡在线播放| 波野结衣二区三区在线| av福利片在线观看| av.在线天堂| 高清午夜精品一区二区三区| 少妇人妻一区二区三区视频| av.在线天堂| 精品久久久噜噜| a级一级毛片免费在线观看| 色5月婷婷丁香| 女人久久www免费人成看片 | 久久草成人影院| av专区在线播放| 亚洲av日韩在线播放| 精品午夜福利在线看| 亚洲欧美清纯卡通| 欧美潮喷喷水| 亚洲欧美精品专区久久| 99国产精品一区二区蜜桃av| 亚洲国产精品久久男人天堂| 国产毛片a区久久久久| 亚洲国产精品sss在线观看| 欧美性猛交╳xxx乱大交人| 秋霞在线观看毛片| 精品少妇黑人巨大在线播放 | 精品久久久久久久久久久久久| 大香蕉久久网| 草草在线视频免费看| 久久久久精品久久久久真实原创| 精品久久久久久成人av| 中文字幕久久专区| 国产国拍精品亚洲av在线观看| av在线观看视频网站免费| 麻豆久久精品国产亚洲av| 亚洲激情五月婷婷啪啪| av在线观看视频网站免费| 国产男人的电影天堂91| 亚洲国产色片| 成人美女网站在线观看视频| 亚洲电影在线观看av| 亚洲综合色惰| 毛片女人毛片| 少妇的逼水好多| ponron亚洲| 久久久久久久久久黄片| 精品一区二区三区视频在线| 亚洲,欧美,日韩| av播播在线观看一区| 校园人妻丝袜中文字幕| 亚洲国产成人一精品久久久| 免费观看在线日韩| 亚洲精品成人久久久久久| 久久久久久国产a免费观看| 欧美成人免费av一区二区三区| 精品99又大又爽又粗少妇毛片| 十八禁国产超污无遮挡网站| 日韩国内少妇激情av| 欧美潮喷喷水| 伊人久久精品亚洲午夜| av天堂中文字幕网| 国产免费福利视频在线观看| 成人三级黄色视频| 最近视频中文字幕2019在线8| 国产成人aa在线观看| 国产色婷婷99| 日本黄大片高清| 免费看a级黄色片| 一区二区三区高清视频在线| 青春草视频在线免费观看| 国产高清三级在线| 精品久久久久久成人av| 国产av一区在线观看免费| av专区在线播放| 国产精品久久久久久久电影| 2021天堂中文幕一二区在线观| 有码 亚洲区| 欧美成人一区二区免费高清观看| 啦啦啦啦在线视频资源| 一级爰片在线观看| 成人高潮视频无遮挡免费网站| 日韩av在线免费看完整版不卡| 长腿黑丝高跟| 国模一区二区三区四区视频| 日本黄色片子视频| 国产色婷婷99| 伦理电影大哥的女人| ponron亚洲| 丰满人妻一区二区三区视频av| 国产亚洲最大av| 好男人视频免费观看在线| 精品人妻熟女av久视频| 国产精品一二三区在线看| 大话2 男鬼变身卡| 永久免费av网站大全| ponron亚洲| 国产不卡一卡二| 久久精品影院6| 色播亚洲综合网| 国产乱来视频区| 特级一级黄色大片| 婷婷六月久久综合丁香| 免费观看性生交大片5| 男人和女人高潮做爰伦理| 大又大粗又爽又黄少妇毛片口| 纵有疾风起免费观看全集完整版 | 国产精品,欧美在线| 午夜爱爱视频在线播放| 少妇熟女aⅴ在线视频| 成人国产麻豆网| 大香蕉97超碰在线| 色5月婷婷丁香| 久久婷婷人人爽人人干人人爱| 18禁在线无遮挡免费观看视频| 亚洲av二区三区四区| 97超碰精品成人国产| 欧美区成人在线视频| 久久99蜜桃精品久久| 97超视频在线观看视频| 最近中文字幕高清免费大全6| 欧美xxxx性猛交bbbb| 国产精品电影一区二区三区| 2021天堂中文幕一二区在线观| 91久久精品国产一区二区成人| 国产在视频线在精品| 欧美日韩一区二区视频在线观看视频在线 | 99久国产av精品| 熟女电影av网| 激情 狠狠 欧美| 欧美成人免费av一区二区三区| 国产熟女欧美一区二区| 国产精品女同一区二区软件| 观看免费一级毛片| 国产精品福利在线免费观看| 男插女下体视频免费在线播放| 可以在线观看毛片的网站| 国产真实伦视频高清在线观看| 亚洲国产欧美在线一区| 精品免费久久久久久久清纯| 搡老妇女老女人老熟妇| 熟女电影av网| 婷婷色麻豆天堂久久 | 99久久精品国产国产毛片| 99热这里只有是精品在线观看| 免费人成在线观看视频色| 精品少妇黑人巨大在线播放 | 成年免费大片在线观看| 国产免费又黄又爽又色| 国产极品精品免费视频能看的| 国产一区二区在线观看日韩| 欧美zozozo另类| 亚洲av电影不卡..在线观看| 91精品国产九色| 最近中文字幕高清免费大全6| 一本一本综合久久| 啦啦啦啦在线视频资源| av在线亚洲专区| 神马国产精品三级电影在线观看| 亚洲18禁久久av| 亚洲va在线va天堂va国产| 99久国产av精品| 国产毛片a区久久久久| 99热网站在线观看| 三级毛片av免费| 亚洲国产精品合色在线| 国产欧美另类精品又又久久亚洲欧美| 91精品国产九色| 人妻夜夜爽99麻豆av| 嫩草影院新地址| 欧美一级a爱片免费观看看| 日韩欧美精品v在线| 免费观看的影片在线观看| 亚洲精华国产精华液的使用体验| 日日撸夜夜添| 午夜精品国产一区二区电影 | 欧美成人免费av一区二区三区| 大话2 男鬼变身卡| 午夜福利在线观看吧| 久久久欧美国产精品| 亚洲国产精品专区欧美| 三级国产精品片| 国产精品一区二区在线观看99 | 久久精品久久精品一区二区三区| 国产午夜精品一二区理论片| av又黄又爽大尺度在线免费看 | 日韩人妻高清精品专区| 欧美成人一区二区免费高清观看| 免费看光身美女| 日韩三级伦理在线观看| 精品久久久久久成人av| 91aial.com中文字幕在线观看| 日韩高清综合在线| 99久久九九国产精品国产免费| 久久精品影院6| 久久久久久久国产电影| 日本与韩国留学比较| 六月丁香七月| 黄色配什么色好看| 日韩强制内射视频| 国产淫语在线视频| 亚洲伊人久久精品综合 | 亚洲av日韩在线播放| 岛国在线免费视频观看| 国产成人一区二区在线| 国产成人福利小说| 亚洲在线自拍视频| 国产成人免费观看mmmm| 国产精品蜜桃在线观看| 神马国产精品三级电影在线观看| 亚洲av免费高清在线观看| 国产伦理片在线播放av一区| 哪个播放器可以免费观看大片| 亚洲精品日韩在线中文字幕| 少妇熟女aⅴ在线视频| 免费看a级黄色片| 免费观看人在逋| 超碰97精品在线观看| 亚洲,欧美,日韩| 欧美bdsm另类| 最近2019中文字幕mv第一页| 亚洲第一区二区三区不卡| 国产黄色小视频在线观看| 天堂网av新在线| 精华霜和精华液先用哪个| 最近2019中文字幕mv第一页| 午夜福利网站1000一区二区三区| 国产精华一区二区三区| 99久久九九国产精品国产免费| 淫秽高清视频在线观看| 日韩一区二区视频免费看| 久久这里只有精品中国| 久热久热在线精品观看| 日产精品乱码卡一卡2卡三| 亚洲av免费高清在线观看| 少妇高潮的动态图| 免费观看人在逋| 亚洲av中文字字幕乱码综合| 天堂√8在线中文| 亚洲人成网站在线播| 亚洲精品乱码久久久久久按摩| 亚洲欧美清纯卡通| 99视频精品全部免费 在线| 免费看美女性在线毛片视频| 免费av观看视频| 最近最新中文字幕免费大全7| 婷婷色av中文字幕| 国产精品国产三级专区第一集| 男人的好看免费观看在线视频| 大香蕉久久网| 在线播放国产精品三级| 校园人妻丝袜中文字幕| 免费黄色在线免费观看| 免费黄网站久久成人精品| 九九在线视频观看精品| 久久人人爽人人爽人人片va| 午夜激情欧美在线| 国产精品99久久久久久久久| 在线观看一区二区三区| 男人的好看免费观看在线视频| 97超碰精品成人国产| 草草在线视频免费看| 少妇人妻精品综合一区二区| 免费一级毛片在线播放高清视频| 国产成人aa在线观看| 亚洲av电影不卡..在线观看| 国产午夜精品一二区理论片| 男女下面进入的视频免费午夜| 久久国内精品自在自线图片| 日韩在线高清观看一区二区三区| av国产免费在线观看| 国模一区二区三区四区视频| 可以在线观看毛片的网站| 亚洲精品成人久久久久久| 日本欧美国产在线视频| 免费一级毛片在线播放高清视频| 美女大奶头视频| 欧美xxxx黑人xx丫x性爽| 国产日韩欧美在线精品| 国产熟女欧美一区二区| 黑人高潮一二区| 国产真实乱freesex| 99热这里只有是精品50| 日韩av在线免费看完整版不卡| 乱人视频在线观看| 久久久国产成人免费| 黄色欧美视频在线观看| 国产免费福利视频在线观看| 日韩强制内射视频| 内地一区二区视频在线| 天堂中文最新版在线下载 | 欧美高清性xxxxhd video| 午夜福利高清视频| 69人妻影院| 18禁在线无遮挡免费观看视频| 色综合站精品国产| 99久久精品国产国产毛片| 99热这里只有是精品50| 国产精品国产三级专区第一集| 国产精品久久视频播放| 大又大粗又爽又黄少妇毛片口| 国产三级在线视频| 成人无遮挡网站| 又黄又爽又刺激的免费视频.| 欧美日本亚洲视频在线播放| 三级国产精品片| 国产三级中文精品| 国产老妇伦熟女老妇高清| 熟女电影av网| 日本欧美国产在线视频| 两个人的视频大全免费| 亚洲成人久久爱视频| 一区二区三区高清视频在线| 最近2019中文字幕mv第一页| 久久久午夜欧美精品| 午夜激情福利司机影院| 国产三级中文精品| 国产亚洲一区二区精品| 成人三级黄色视频| 成人亚洲精品av一区二区| 亚洲乱码一区二区免费版| 免费黄色在线免费观看| 免费看a级黄色片| 久久精品国产亚洲网站| 国产一区二区三区av在线| 51国产日韩欧美| 国产精品综合久久久久久久免费| 久久精品夜色国产| 日本与韩国留学比较| 久久99热6这里只有精品| av在线观看视频网站免费| 久久精品国产亚洲网站| 精品人妻视频免费看| 一级av片app| 精品人妻视频免费看| 又爽又黄无遮挡网站| 亚洲欧美日韩无卡精品| 热99re8久久精品国产| 精品不卡国产一区二区三区| 一级二级三级毛片免费看| 毛片女人毛片| 久久久色成人| 亚洲最大成人中文| 噜噜噜噜噜久久久久久91| 高清毛片免费看| 午夜日本视频在线| 国产真实乱freesex| 伦理电影大哥的女人| 中文字幕熟女人妻在线| 国产色爽女视频免费观看| 精华霜和精华液先用哪个| 亚洲人与动物交配视频| 夫妻性生交免费视频一级片| 九九热线精品视视频播放| 国产中年淑女户外野战色| 国产av码专区亚洲av| 久久精品久久精品一区二区三区| 在线播放国产精品三级| 国产真实伦视频高清在线观看| 美女cb高潮喷水在线观看| 久久久久国产网址| 国产精品日韩av在线免费观看| 亚洲国产精品久久男人天堂| 国产三级中文精品| 亚洲高清免费不卡视频| 亚洲av中文字字幕乱码综合| 日韩欧美三级三区| 特大巨黑吊av在线直播| 蜜臀久久99精品久久宅男| 亚洲va在线va天堂va国产| 白带黄色成豆腐渣| 亚洲av电影在线观看一区二区三区 | 国产伦理片在线播放av一区| 爱豆传媒免费全集在线观看| 成人亚洲精品av一区二区| 国产三级中文精品| 国产精华一区二区三区| 国产精品.久久久| 日韩强制内射视频| 伦精品一区二区三区| 岛国在线免费视频观看| 99视频精品全部免费 在线| 国产亚洲精品久久久com| 午夜福利网站1000一区二区三区| 久久久久性生活片| 精品久久久久久久久亚洲| 久久精品国产亚洲av天美| 我的女老师完整版在线观看| 国产精品麻豆人妻色哟哟久久 | 水蜜桃什么品种好| 亚洲成人精品中文字幕电影| 欧美另类亚洲清纯唯美| 亚洲一级一片aⅴ在线观看| 免费观看精品视频网站| 丰满乱子伦码专区| 天天一区二区日本电影三级| 最后的刺客免费高清国语| 午夜福利在线观看免费完整高清在| 久久久色成人| 国产av在哪里看| 久久久国产成人精品二区| 久久久久精品久久久久真实原创| 男女那种视频在线观看| 久久久久久大精品| 九草在线视频观看| 99久久精品国产国产毛片| 在线观看av片永久免费下载| 日韩一区二区三区影片| 国内揄拍国产精品人妻在线| 国产成人a区在线观看| 老女人水多毛片| 免费不卡的大黄色大毛片视频在线观看 | 看黄色毛片网站| 视频中文字幕在线观看| 伦理电影大哥的女人| 最近手机中文字幕大全| 亚洲激情五月婷婷啪啪| 中文字幕亚洲精品专区| 99久久中文字幕三级久久日本| 最近最新中文字幕大全电影3| 一边亲一边摸免费视频| 国产免费男女视频| 日韩欧美精品免费久久| 国产精品久久久久久精品电影| 99久久成人亚洲精品观看| 亚洲国产欧洲综合997久久,| 波多野结衣巨乳人妻| 男女国产视频网站| 国产精品国产高清国产av| 亚洲精品国产av成人精品| 22中文网久久字幕| 非洲黑人性xxxx精品又粗又长| 亚洲精品日韩av片在线观看| 国产精品久久电影中文字幕| 久久久久久久久久黄片| 欧美日本视频| 老女人水多毛片| 九色成人免费人妻av| 亚洲性久久影院| 69av精品久久久久久| 小说图片视频综合网站| 国产麻豆成人av免费视频| 亚洲成人精品中文字幕电影| 日韩欧美国产在线观看| 久久精品久久久久久噜噜老黄 | 边亲边吃奶的免费视频| 午夜日本视频在线| 男插女下体视频免费在线播放| 伦理电影大哥的女人| 久久精品综合一区二区三区| 级片在线观看| 免费播放大片免费观看视频在线观看 | 美女内射精品一级片tv| 精品少妇黑人巨大在线播放 | 伊人久久精品亚洲午夜| 国产在视频线精品| 又粗又爽又猛毛片免费看| 如何舔出高潮| 亚洲国产欧美在线一区| 欧美激情国产日韩精品一区| 亚洲综合精品二区| 亚洲精品亚洲一区二区| 欧美成人午夜免费资源| 男女边吃奶边做爰视频| 99热这里只有精品一区| 人人妻人人澡欧美一区二区| 精华霜和精华液先用哪个| 99在线视频只有这里精品首页| 一个人看的www免费观看视频| av在线播放精品| 日韩成人伦理影院| or卡值多少钱| 在线播放无遮挡| av在线天堂中文字幕| 久久亚洲国产成人精品v| 久久久精品欧美日韩精品| 日韩在线高清观看一区二区三区| or卡值多少钱| 尾随美女入室| or卡值多少钱| 老司机影院毛片| 亚洲人与动物交配视频| 搡女人真爽免费视频火全软件| 国产色婷婷99| 岛国在线免费视频观看| 日韩欧美 国产精品| 特大巨黑吊av在线直播| 一边亲一边摸免费视频| 少妇高潮的动态图| 久99久视频精品免费| 精华霜和精华液先用哪个| 精品久久久噜噜| 嫩草影院新地址| 寂寞人妻少妇视频99o| 久久久国产成人精品二区| 亚洲精品乱久久久久久| 国产免费男女视频| 老司机影院成人| 国产精品无大码| 99热这里只有是精品50| 日本爱情动作片www.在线观看| 国产成人精品久久久久久| 嘟嘟电影网在线观看| 又粗又硬又长又爽又黄的视频| 伊人久久精品亚洲午夜| 99在线视频只有这里精品首页| 精华霜和精华液先用哪个| 成人综合一区亚洲| 少妇裸体淫交视频免费看高清| 日本-黄色视频高清免费观看| 久久人人爽人人片av| 久久久久久伊人网av| 日韩成人伦理影院| 国产亚洲av嫩草精品影院| 亚洲欧美日韩无卡精品| 亚洲三级黄色毛片| 国产免费视频播放在线视频 | 亚洲成人中文字幕在线播放| 99热这里只有是精品50| 久久久精品欧美日韩精品| 人妻制服诱惑在线中文字幕| 一个人观看的视频www高清免费观看| 亚洲精华国产精华液的使用体验| 日本黄色视频三级网站网址| 一个人观看的视频www高清免费观看| 少妇熟女aⅴ在线视频| 国产高清有码在线观看视频| 看免费成人av毛片| 精品久久久久久久末码| 精品国产一区二区三区久久久樱花 | 亚洲成人精品中文字幕电影| 国产久久久一区二区三区| 精品久久久久久久久av| 国产精品一区二区在线观看99 | 久久精品国产鲁丝片午夜精品| 久久精品综合一区二区三区| 色网站视频免费| 成人鲁丝片一二三区免费| 亚洲国产欧美在线一区| 99久久九九国产精品国产免费| 狂野欧美白嫩少妇大欣赏| 日本wwww免费看| 日日摸夜夜添夜夜添av毛片| 亚洲精品国产av成人精品| 国产精品无大码| 啦啦啦韩国在线观看视频| 欧美成人a在线观看| 久久99精品国语久久久| 日本一本二区三区精品| 日韩成人av中文字幕在线观看| 边亲边吃奶的免费视频| 亚洲国产精品国产精品| 亚洲av一区综合| 日本黄大片高清| 九九爱精品视频在线观看| 听说在线观看完整版免费高清| 亚洲欧美精品自产自拍| 免费黄网站久久成人精品| 中国国产av一级| 亚洲国产高清在线一区二区三| 超碰97精品在线观看| 精品久久久久久久久av| 欧美日韩综合久久久久久| 国产单亲对白刺激| 国国产精品蜜臀av免费| 全区人妻精品视频| 免费看日本二区| 欧美一区二区国产精品久久精品| 国产高清不卡午夜福利| 一级av片app| 色综合站精品国产| 久久精品91蜜桃| 如何舔出高潮| 国产精品爽爽va在线观看网站| 日本-黄色视频高清免费观看| 中文字幕免费在线视频6| 欧美性猛交╳xxx乱大交人| 我的老师免费观看完整版| 亚洲久久久久久中文字幕| 深爱激情五月婷婷| 久久久久久九九精品二区国产| 乱码一卡2卡4卡精品| 久久99热这里只频精品6学生 | 久久久色成人| 精品一区二区免费观看|