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

    珠江磨刀門河口潮波振幅梯度與上下游動力邊界的關(guān)系異變研究*

    2022-05-30 01:15:40歐素英蔡華陽楊清書
    海洋與湖沼 2022年3期
    關(guān)鍵詞:分潮馬口河口

    李 博 楊 昊 歐素英 蔡華陽①劉 鋒 楊清書

    (1. 中山大學(xué)海洋工程與技術(shù)學(xué)院 河口海岸研究所 廣東廣州 510275; 2. 河口水利技術(shù)國家地方聯(lián)合工程實驗室 廣東廣州510275; 3. 廣東省海岸與島礁工程技術(shù)研究中心 廣東廣州 510275; 4. 南方海洋科學(xué)與工程廣東省實驗室(珠海) 廣東珠海519000)

    河口是陸海徑潮動力耦合的獨特區(qū)域, 既有上游徑流匯入, 又有外海潮波上溯, 還受人類活動疊加效應(yīng)的影響, 導(dǎo)致其動力結(jié)構(gòu)復(fù)雜多變。河口動力的典型特征包括徑潮相互作用、斜壓效應(yīng)、波流耦合等,其中, 徑潮相互作用對洪水下泄、鹽淡水混合、最大渾濁帶形成、營養(yǎng)鹽和污染物輸移等物理、化學(xué)和生物過程均有直接影響。潮波振幅梯度, 即潮波振幅沿程的平均變化率, 綜合反映了徑流、潮汐和底床摩擦等多因子影響在河口空間上的差異, 是探討河口復(fù)雜徑潮動力結(jié)構(gòu)的有效切入點。因此, 揭示潮波振幅梯度與上下游動力邊界(即上游徑流和外海潮汐)的耦合過程及機制是河口動力學(xué)研究的基礎(chǔ)科學(xué)前沿問題, 可為河口海岸區(qū)域自然災(zāi)害(如咸潮上溯、洪澇、赤潮等)的防治和水資源的高效開發(fā)利用等提供重要科學(xué)依據(jù)(Schuttelaarset al, 2013; Caoet al, 2020; 楊昊等, 2020)。

    國內(nèi)外針對河口徑潮相互作用的研究由來已久。研究表明, 當潮波受到徑流和地形的調(diào)制作用時呈現(xiàn)出明顯的非線性變化(Hoitinket al, 2017; Jiet al,2019), 河口潮波衰減主要通過流量和底床的摩擦效應(yīng)實現(xiàn)(Godin, 1985; Savenije, 1992; Horrevoetset al,2004)。當流量增大時, 全日分潮和半日分潮的振幅以不同幅度衰減, 而河道底床下切, 水深增大可導(dǎo)致底床摩擦效應(yīng)減弱, 引起沿程潮波振幅不同程度的增加(Caiet al, 2012)。對河口地形進行概化并考慮下泄徑流影響, 采用一維潮波傳播解析模型可通過潮波振幅梯度和余水位梯度等參數(shù)揭示潮波傳播的時空變化特征及其動力學(xué)機制(Caiet al, 2016)。當河口動力以徑流作用為主時, 傳統(tǒng)調(diào)和分析(Pawlowiczet al,2002)的水位預(yù)報誤差明顯增大。針對流量影響下河口潮波傳播的非定常和非線性問題, 對水位序列采用連續(xù)小波變換(Guoet al, 2015; Moftakhariet al,2016)、經(jīng)驗?zāi)B(tài)分解(Panet al, 2018a)、經(jīng)驗正交函數(shù)(Panet al, 2019)、非平穩(wěn)潮汐調(diào)和分析(Matteet al,2013, 2014; Panet al, 2018b)等處理方法探討流量對潮波傳播的影響過程及機制均取得較好效果。

    近年來, 流量對潮波傳播的衰減及其閾值效應(yīng)研究得到廣泛關(guān)注。研究表明, 在長江河口感潮河段,流量與潮波振幅梯度之間存在閾值效應(yīng), 即流量小于臨界閾值時, 隨著流量增大潮波衰減效應(yīng)增強, 但當潮波振幅梯度達到極小值后, 隨著流量增大衰減效應(yīng)反而減弱(Caiet al, 2019; 張先毅等, 2019)。在珠江磨刀門河口也觀察到類似現(xiàn)象, 即在甘竹-馬口河段潮波振幅梯度隨流量變化有先減小后上升的特點,在人類活動干預(yù)較弱的自然演變階段, 潮波振幅梯度在流量為5 000~10 000 m3/s 時達到極小值, 而在人類活動強烈影響后的恢復(fù)調(diào)整階段, 該流量閾值上升至15 000~20 000 m3/s, 且其對應(yīng)的潮波振幅梯度極小值增大(Yanget al, 2020; 張先毅等, 2020), 反映因采砂導(dǎo)致的河槽容積增大對潮波振幅梯度閾值的增大效應(yīng)。綜上所述, 對于河口徑潮相互作用及潮波振幅梯度與流量動力邊界關(guān)系的研究已有較多成果,但是對于定量辨識流量驅(qū)動下潮汐調(diào)和分析結(jié)果的階段性演變以及潮波振幅梯度與下游動力邊界關(guān)系的研究等相關(guān)科學(xué)問題還尚待深入。

    位于我國粵港澳大灣區(qū)國家戰(zhàn)略核心經(jīng)濟圈的珠江三角洲是全球受到人類活動影響程度最大的區(qū)域之一。高強度的人類活動(如水庫建設(shè)、口門圍墾、無序采砂、航道疏浚等)導(dǎo)致來水來沙和地形地貌發(fā)生明顯變化, 進而顯著改變珠江河網(wǎng)的徑潮動力格局。水庫建設(shè)使河流挾沙量銳減, 三角洲河道下切加劇; 地形變化導(dǎo)致西江、北江的水沙分配格局發(fā)生變化, 河網(wǎng)中上部地區(qū)余水位梯度變緩(Liuet al, 2019);口門圍墾、口門整治疊加海平面上升的影響則導(dǎo)致磨刀門河口的形態(tài)向窄深化趨勢演變(Tanet al, 2019)。以上變化導(dǎo)致珠江主要出??? 即磨刀門河口的潮流界和潮區(qū)界顯著向上游移動, 咸潮上溯災(zāi)害風(fēng)險增大, 對粵港澳大灣區(qū)城市的居民生活、工農(nóng)業(yè)生產(chǎn)用水的水質(zhì)安全形成嚴峻挑戰(zhàn)。因此, 本文針對強人類活動驅(qū)動下的徑潮動力異變問題, 采用數(shù)據(jù)驅(qū)動模型探究流量驅(qū)動下珠江磨刀門河口潮波傳播特征的異變過程, 并分析潮波調(diào)和常數(shù)的階段性演變, 揭示潮波振幅梯度對上下游動力邊界的響應(yīng)過程及機制, 研究成果可為受到高強度人類活動影響的河口三角洲的治理和保護等提供科學(xué)依據(jù)。

    1 研究區(qū)域概況

    珠江河網(wǎng)(圖1)流域面積達4.5×105km2, 具有三江(東江、西江、北江)匯流, 八口(自東向西依次為虎門、蕉門、洪奇門、橫門、磨刀門、雞啼門、虎跳門、崖門)入海的特點, 河網(wǎng)縱橫交錯導(dǎo)致徑潮動力相互作用過程復(fù)雜。珠江每年約有2 823×108m3淡水流量及72.4×106t 懸沙量流入南海(Liuet al, 2018), 平均流量約為8 952 m3/s, 口門附近潮差為1.0~1.7 m。磨刀門河口為西江入??? 是八大口門中輸水輸沙量最大的河口, 該河口屬典型的河控型河口。據(jù)三角洲頂端馬口水文站1959~2016 年月均流量資料, 其多年平均流量為7 078 m3/s, 月均流量最大值為29 000 m3/s(6 月), 最小值為1 210 m3/s (1 月), 洪枯季差異明顯(楊昊等, 2020)。潮汐以不規(guī)則半日潮為主, 日不等現(xiàn)象顯著, 燈籠山站平均漲落潮歷時分別為5.37 h 和7.25 h, 口門處三灶站潮型系數(shù)為1.30。

    圖1 研究區(qū)域及所用潮位站或水文站位置Fig.1 The study area and the deployment of the tidal gauging stations and the hydrological station

    隨著經(jīng)濟社會高速發(fā)展, 強人類活動已成為河口三角洲系統(tǒng)演變的第三驅(qū)動力(陳吉余等, 2008)。口門圍墾、聯(lián)圍筑閘、河道采砂及上游流域水庫建設(shè)等高強度人類活動極大影響磨刀門河口的河道形態(tài)和徑潮動力, 其“動力-沉積-地貌”體系發(fā)生顯著變化(吳超羽, 2018)。20 世紀80 年代前, 磨刀門地區(qū)主要人類活動為大規(guī)模灘涂圍墾, 至 1 9 9 7 年已有16.5 km2的淺灘被圍填成陸地(Tanet al, 2019)。圍墾、疏浚、口門整治等一系列工程措施使河口口門區(qū)水域面積減小, 入海水道水深增大, 至20 世紀90 年代初入??陂T向海延伸約13 km (何用等, 2018)。其次, 20世紀50~70 年代中期, 大規(guī)模的聯(lián)圍筑閘導(dǎo)致三角洲腹部水位壅高, 中上游水面比降變緩, 流速減小, 水流挾沙能力降低, 河床產(chǎn)生淤積(劉鋒等, 2011)。20世紀80 年代中期開始, 由于用沙需求激增, 在珠江河網(wǎng)區(qū)出現(xiàn)大規(guī)模河道采砂活動。磨刀門采砂活動主要集中在珠海大橋以北河段, 1991~2000 年共挖沙420×104m3(韓志遠等, 2010)。截至2006 年, 西江流域大壩總庫容達5.38×1010m3, 占其入海淡水流量的24.7% (Liuet al, 2018)。大壩建成后流域內(nèi)大部分沉積物被截留在上游地區(qū), 西江干支流輸沙量均呈下降趨勢, 下游泥沙供給量銳減(Wanget al, 2011)。在上述強人類活動的綜合影響下, 珠江磨刀門河口沿程河床大幅下切, 河槽容積增大, 導(dǎo)致近年來潮汐動力顯著增強, 咸潮上溯加劇(韓志遠等, 2010; 張先毅等, 2020)。

    2 數(shù)據(jù)與方法

    2.1 數(shù)據(jù)及其處理方法

    本文所用數(shù)據(jù)取自《廣東省水文年鑒》第八卷和廣東省水文局, 潮位原始數(shù)據(jù)高程為凍結(jié)基面, 已統(tǒng)一轉(zhuǎn)換至珠江基面。數(shù)據(jù)包括磨刀門河口沿程4 個潮位站(三灶、燈籠山、竹銀、甘竹站)逐日高、低潮位數(shù)據(jù)和相應(yīng)時段馬口水文站日均流量及逐日高、低潮位數(shù)據(jù), 并對潮位和流量數(shù)據(jù)進行分段三次Hermite插值處理后獲得的逐時數(shù)據(jù)。站點具體信息見表1。

    表1 磨刀門河口各站點數(shù)據(jù)信息Tab.1 Information of the adopted tidal gauging stations and hydrological station in the Modaomen estuary

    2.2 R_TIDE 數(shù)據(jù)驅(qū)動模型

    本文采用徑潮耦合的數(shù)據(jù)驅(qū)動模型(river-tidal harmonic analysis, R_TIDE)對逐時潮位數(shù)據(jù)進行調(diào)和分析(歐素英等, 2017), 模型可輸出重構(gòu)的逐時水位、特定分潮的振幅和遲角等。在徑流動力占主導(dǎo)的河口,潮波上溯受到徑流季節(jié)性周期變化的調(diào)制影響, 傳統(tǒng)的潮汐調(diào)和分析模型對河口區(qū)潮波特征的分析和預(yù)報誤差較大, 不能滿足徑潮信號分離的要求。為有效辨識潮波在流量影響下的傳播特征, 基于Matte 等(2013, 2014)提出的非平穩(wěn)調(diào)和分析(nonstationary tidal harmonic analysis, NS_TIDE)思路, 假設(shè)河口任意位置潮波振幅和遲角主要受上游流量和地形變化的非線性調(diào)制影響, 在此基礎(chǔ)上進一步假定分析時段內(nèi)河口地形邊界不變, 將信號分成由流量Q引起的水位變化和徑流調(diào)制影響下的潮水位變化, 通過信號分離定量分析感潮河段徑潮動力相互作用。

    無徑流調(diào)制時, 對于月球和太陽引起的周期性潮汐現(xiàn)象可看作是許多假想天體引起的簡單的潮汐簡諧波動的總和(Pawlowiczet al, 2002), 即潮位曲線可近似看作是由許多振幅、周期、相位不同的分潮組成。潮位可表達成以下多個余弦函數(shù)的疊加形式:

    其中,z(t)為站點的實測潮位;t為時間;z0為平均海面高度;N為分潮個數(shù);fi為分潮振幅的訂正因子, 為時間的函數(shù), 常取一年的中值;Hi為分潮平均振幅;σi為分潮角頻率;Vi為分潮初相角;ui為天文相角的交角訂正角,gi為由于海底摩擦、海水慣性引起的遲角。其中Hi和gi合稱為分潮調(diào)和常數(shù)。

    Kukulka 等(2003)基于一維圣維南方程組推導(dǎo)得出河口三角洲內(nèi)任意位置x的潮波振幅變化:

    其中,j=1,2,3,…,m,m為逐時序列的數(shù)據(jù)個數(shù);r為衰減系數(shù);ε0(tj)為河口口門位置x=0 處的潮波振幅;a'為常數(shù)。衰減系數(shù)r為潮波傳播速度、水深、流量和河道摩擦系數(shù)的函數(shù)(Jayet al, 1997)。若摩擦系數(shù)采用曼寧公式, 且摩擦系數(shù)與潮波傳播速度均為水深的函數(shù), 則衰減系數(shù)r可簡寫為水深h(tj)和流量Q(tj)的函數(shù):

    或:

    其中,h(tj)和Q(tj)表示隨時間變化的水深和流量;p0、p1、p2、β、γ均為待定參數(shù)。式(3)未考慮外海非潮動力邊界比如風(fēng)暴等因素的影響, 這是與Matte 等(2013,2014)提出的非定常調(diào)和分析模型的不同之處。假定研究時間段河道底床高程不隨時間變化, 則式(4)中水深項可歸到常數(shù)項, 衰減系數(shù)r可簡寫為

    即河口潮汐非線性變化僅受流量驅(qū)動影響。結(jié)合式(1)、(2)、(5)可得流量影響下任意位置x的潮水位變化:

    式(6)的待定系數(shù)bi、ci、γ可通過粒子群優(yōu)化算法(Kennedyet al, 1995)或求解帶約束的非線性多變量方程組來確定該參數(shù)的最優(yōu)值。由于流量對不同分潮的影響程度存在差異(Guoet al, 2015), 半日分潮的衰減效應(yīng)大于全日分潮, 將流量對每一個分潮族的調(diào)制影響均采用不同的待定系數(shù), 則式(6)可表示為:

    其中,

    式(7)即為流量驅(qū)動下的R_TIDE 模型,S(tj)為底床高程或水深變化、海平面、流量等引起的平均水面變化, 稱為余水位項;F(tj)代表不同流量條件下徑潮耦合引起的水位變化, 稱為徑潮耦合項。d0、d1、v0、v1、r0、r1為模型的待定參數(shù), 通過求解帶約束的非線性多變量方程組來確定。同時, 由于磨刀門河口徑流動力占優(yōu)勢導(dǎo)致徑潮動力非線性因素作用較強, 洪水流量時導(dǎo)致三角洲上段為徑流控制, 無潮汐波動。因此, 在R_TIDE 模型中引入臨界流量QC, 當連續(xù)2 d 潮差小于某個值(如0.001 m)時所對應(yīng)的最小流量即為臨界流量, 當流量大于QC時, 表示該站點及其以上河段不存在潮汐信號。

    假設(shè)河口三角洲任意位置采用時間間隔為Δt的觀測潮位y(tj)(j=1,2,3,…,m)序列, 根據(jù)實測資料序列長度mΔt、分潮頻率差和 Rayleigh 準則判據(jù) Δσ=max[(mΔt)-1,σi-σj], 選擇待提取的分潮, 將參數(shù)tj、σj、Q(tj)代入式(9), 要求模型計算結(jié)果y(tj)和實測水位z(tj)的誤差平方和達到最小, 即:

    選擇對應(yīng)資料長度的全日分潮、半日分潮等分潮族, 對實測潮位進行回歸擬合, 回歸模型效果以均方根誤差ERMS和相關(guān)指數(shù)R2表示, 即:

    2.3 徑潮動力的特征指標

    3 結(jié)果

    3.1 徑潮動力格局的階段性演變

    前期研究表明, 基于燈籠山-馬口河段月均余水位梯度(或潮波振幅梯度絕對值)與流量的雙累積曲線變化趨勢, 可知曲線斜率在1990 年與2000 年發(fā)生明顯改變(張先毅等, 2020)。根據(jù)上述結(jié)果, 可將磨刀門河口徑潮動力格局的演變分為三個階段: 1990 年前為自然演變和人類活動共同作用的階段, 但此階段人類活動尚未對河口徑潮動力格局產(chǎn)生明顯影響, 河口處在動態(tài)平衡期, 以自然演變?yōu)橹? 稱為“平衡期”(簡稱P1); 2000 年后, 河道中下游采砂活動基本已停止, 灘涂圍墾面積銳減, 河口在人類活動影響下已發(fā)生徑潮動力格局轉(zhuǎn)換, 為強人類活動干預(yù)后河口的自適應(yīng)調(diào)整階段, 稱為“調(diào)整期”(簡稱P2); 1991~2000年為人類活動影響最強烈的時期, 河口從自然演變?yōu)橹鞯钠胶鈶B(tài)向強人類活動干預(yù)后的平衡態(tài)轉(zhuǎn)變,為兩個平衡態(tài)的過渡階段, 稱為“過渡期”。

    3.2 R_TIDE 數(shù)據(jù)驅(qū)動模型率定效果

    圍繞珠江磨刀門河口徑潮動力異變問題, 采用流量驅(qū)動的R_TIDE 數(shù)據(jù)驅(qū)動模型對徑潮信號進行分離。表2 為磨刀門河口兩個階段沿程站點的臨界流量QC、模型主要率定參數(shù)及相關(guān)的模型評價指標(ERMS和R2)結(jié)果。C1~C9為不同分潮族所使用的參數(shù), 其中C1代表長周期氣象分潮, 其周期約為15 d, 反映一個大小潮周期變化;C2、C3分別代表全日和半日分潮族,為主導(dǎo)該區(qū)域潮汐動力的分潮族;C4~C9分別代表短周期波動, 由分潮淺水變形效應(yīng)(如M4、MS4、M6、M8)產(chǎn)生, 反映地形、底床摩擦等邊界條件的影響。由表2 可知, 越往下游ERMS值越大, 表明模型對潮流優(yōu)勢區(qū)域水位的重構(gòu)效果略低于上游區(qū)域, 但是其R2均大于0.79, 模型結(jié)果合理。對于馬口站, 由于其為流量驅(qū)動邊界, 除受洪峰影響外還受邊界效應(yīng)影響, 因此其ERMS相對較大。調(diào)整期階段的ERMS比平衡期普遍增大,R2減小, 表明強人類活動對地形的改造導(dǎo)致水面線為適應(yīng)地形變化, 水位對流量的響應(yīng)敏感度有所降低, 由此一定程度上影響了模型的準確性。

    圖2 為R_TIDE 模型在平衡期和調(diào)整期重構(gòu)得到的不同站點日均水位與實測值的對比結(jié)果, 其中黑色虛線表示實測值與重構(gòu)值完全一致。當流量大于臨界流量QC時潮波基本消失, 率定過程中流量大于QC時間段所對應(yīng)的潮波信號為虛假信號, 結(jié)果分析部分已剔除這部分信息。由圖2 可知, 馬口站下游站點模型擬合的效果較好(ERMS介于0.14~0.22 m,R2介于0.79~0.99)。但由于夏季洪峰出現(xiàn)次數(shù)多且持續(xù)時間長, 對徑潮信號分離產(chǎn)生較大影響, 模型誤差較大,因此當流量接近臨界流量QC時, 馬口站重構(gòu)水位比實測值偏大。

    圖2 R_TIDE 數(shù)據(jù)驅(qū)動模型率定結(jié)果Fig.2 Calibration of the adopted R_TIDE data-driven model

    由表2 可見, 調(diào)整期的QC和γ均大于平衡期, 主要原因為調(diào)整期階段潮汐動力增強, 潮區(qū)界向上游移動。流量指數(shù)γ反映地形邊界條件變化, 其值增大表明河床下切, 過水斷面增大。其中三灶站變化最大,反映入海水道水深劇增, 而上游各站點水深均有不同程度增大。竹銀、燈籠山站主要受采砂影響, 而馬口、甘竹站點除受采砂影響外, 還受上游流域水庫調(diào)蓄等影響, 導(dǎo)致馬口水文站的分流、分沙比發(fā)生變化, 進而顯著改變底床形態(tài)。竹銀及其下游站點為潮流優(yōu)勢河段,QC均一致, 而甘竹、馬口則受強烈徑流影響, 流量越大, 潮區(qū)界越往下游移動, 因此,QC逐漸減小。

    表2 R_TIDE 數(shù)據(jù)驅(qū)動模型的部分參數(shù)及效果評價Tab.2 Calibrated parameters derived from R_TIDE data-driven model and its model performance

    3.3 M2 分潮振幅及其梯度的階段性變化特征

    由于磨刀門水位及其梯度發(fā)生階段性異變, 導(dǎo)致潮波振幅及其梯度亦發(fā)生顯著變化。磨刀門河口潮汐動力主要由M2分潮主導(dǎo), 因此選用其作為代表性分潮進行分析。M2分潮振幅及其梯度的季節(jié)及年均變化如表3 所示。潮波向上游傳播過程中振幅沿程衰減, 洪季振幅小于枯季(三灶站除外), 表明徑流對潮波傳播存在明顯的衰減作用, 尤其是在磨刀門河口上游。平衡期燈籠山站振幅的洪枯季差異僅為0.03 m,甘竹站則為0.07 m, 調(diào)整期潮波振幅有明顯增大(三灶站除外), 且上游馬口站振幅的洪枯季差異增大(約為0.02 m), 其下游的甘竹、竹銀站則減小, 燈籠山、三灶站振幅的洪枯季差異基本不變。

    表3 磨刀門河口各潮位站M2 分潮振幅及其梯度的季節(jié)變化及年均變化Tab.3 Seasonal and annual alterations in amplitude and its gradient of M2 tidal constituent observed at each tidal gauging station along the Modaomen estuary

    為探究潮波傳播過程的時空不均勻性, 對比分析強人類活動前后各站點日均M2分潮振幅及相鄰站點間日均M2分潮潮波振幅梯度的階段性演變, 繪制時間過程線如圖3 所示, 圖中黑色實線表示階段性演變?yōu)?。圖3a 為強人類活動前后各站點日均M2分潮振幅變化值的年內(nèi)分布, 其中日均M2分潮振幅變化值定義為調(diào)整期的M2分潮振幅的日平均值減去平衡期所對應(yīng)的值, 記為Δη(以下日均M2分潮潮波振幅梯度變化值的定義類似, 記為Δδ)。由圖3a 可見, 站點距離口門越遠其值越大, 其中三灶站振幅年均減小0.01 m, 而馬口站則變化0.11 m。季節(jié)變化也有類似規(guī)律(馬口站夏季和甘竹站冬季除外), 表明振幅變化量的季節(jié)性差異有向上游累積的效應(yīng)。日均潮波振幅梯度的時空變化如圖3b 所示。潮波振幅梯度最小值出現(xiàn)在夏季, 最大值出現(xiàn)在冬季。三灶-甘竹河段的潮波振幅梯度在調(diào)整期階段的洪枯季差異減小(各河段洪枯季差異的階段性變化分別為 7.14×10-7、4.86×10-6、5.11×10-6m-1), 且越往上游減小量越大, 表明強人類活動對潮波振幅梯度的增大效應(yīng)在甘竹以下河段都有累積效應(yīng), 但在甘竹-馬口河段, 洪枯季差異雖減小但幅度小于竹銀-甘竹河段(洪枯季差異的階段性變化為4.24×10-6m-1), 表明強人類活動干預(yù)后甘竹-馬口河段潮波振幅梯度對流量響應(yīng)減弱。從潮波振幅梯度變化量來看, 三灶-燈籠山河段在秋季變化最大,燈籠山-甘竹河段則為夏季變化最大, 甘竹-馬口河段則是夏季。上述結(jié)果表明三灶-燈籠山河段水動力主要受海平面變化控制(秋季三灶余水位最大); 而甘竹-馬口河段受徑流控制為主, 由于該區(qū)域在調(diào)整期明顯受水庫調(diào)蓄影響(即洪季蓄水、枯季放水), 雖然潮波振幅梯度夏季階段性變化最大, 冬季最小, 但是其洪枯季差異減小, 是流量季節(jié)性差異減小的直接反映。

    圖3 強人類活動前后磨刀門河口日均M2 分潮潮波振幅及其梯度的階段性差異Fig.3 Stepwise alteration in tidal amplitude and gradient under intensive human activities in the Modaomen estuary

    3.4 潮波振幅梯度與上下游動力邊界的關(guān)系

    圖4 為三灶-馬口河段M2分潮日均潮波振幅梯度與馬口站日均流量的關(guān)系曲線(圖4a)及其季節(jié)變化(圖4b~4e)。平衡期階段三灶-馬口河段潮波振幅梯度極小值(記為δM)為-1.41×10-5m-1, 對應(yīng)流量閾值為10 500 m3/s, 而在調(diào)整期階段則分別變?yōu)?1.38×10-5m-1和19 500 m3/s (其中流量閾值定義為當潮波振幅梯度達到極小值所對應(yīng)的流量值, 記為QM)。日均尺度下, 閾值效應(yīng)在調(diào)整期的冬季基本消失。從季節(jié)變化來看, 平衡期的秋季潮波振幅梯度極小值最小, 調(diào)整期則為夏季最小, 冬季最大, 大流量下潮波衰減效應(yīng)更明顯。調(diào)整期階段潮波振幅梯度極小值(約為-1.38×10-5m-1)和對應(yīng)流量閾值(約為19 500 m3/s)均大于平衡期。

    圖4 三灶-馬口河段M2 分潮日均潮波振幅梯度(δ)與馬口站日均流量(Q)的閾值關(guān)系(a)及其季節(jié)變化(b~e)Fig.4 Relationship between the daily averaged M2 tidal amplitude gradient over SZ-MK reach (δ) and the river discharge (Q) observed at the MK hydrological station to illustrate the threshold effect (a) and its seasonal alteration (b~e)

    潮波振幅梯度是徑潮動力相互作用的直接結(jié)果,因此類似流量閾值效應(yīng), 外海潮汐動力變化(即周期約為15 d 的大小潮)亦會導(dǎo)致潮波振幅梯度產(chǎn)生振幅閾值效應(yīng)(記為ηM)。圖5 為磨刀門河口三灶-馬口河段日均M2分潮潮波振幅梯度與三灶站日均M2分潮振幅的關(guān)系曲線。由圖5 可見, 日均潮波振幅梯度極小值在調(diào)整期略有增大, 約為2.49×10-7m-1, 對應(yīng)的振幅閾值基本不變(約增大0.14 cm)。圖5b~5e 為日均M2分潮潮波振幅梯度與三灶站日均M2分潮振幅閾值關(guān)系的季節(jié)變化。與圖4 類似, 除調(diào)整期的冬季振幅閾值效應(yīng)消失外(圖5e), 其余季節(jié)均存在較明顯的振幅閾值效應(yīng), 其中夏季振幅梯度極小值最小, 其次為春、秋兩季, 最大為冬季。而平衡期階段秋季的振幅閾值最小, 冬季最大。調(diào)整期冬季閾值效應(yīng)消失的主要原因是下泄流量較小, 即便經(jīng)過流域水庫調(diào)蓄,調(diào)整期口門振幅依舊未能達到使衰減效應(yīng)達到最大的臨界值(約0.43 m)。達到振幅閾值前, 由于大潮流速振幅較大導(dǎo)致其底床摩擦較大, 衰減效應(yīng)增強, 潮波振幅梯度減小。而達到振幅閾值后, 流速振幅增大導(dǎo)致的潮能耗散不足以平衡地形輻聚增強(大潮時水深增大)導(dǎo)致的潮能增加, 因此, 潮波振幅梯度略微增大。圖6 為相鄰站點間沿程潮波振幅梯度隨口門三灶站振幅的變化曲線。結(jié)果表明, 振幅閾值效應(yīng)僅存在于徑流動力占優(yōu)勢的河段, 主要存在于甘竹-馬口河段(僅調(diào)整期的冬季消失), 另外, 在竹銀-甘竹河段(春、夏季和平衡期的秋季)及燈籠山-竹銀河段(調(diào)整期的夏季)也有類似現(xiàn)象, 而在三灶-燈籠山河段, 該現(xiàn)象則不存在。對于甘竹-馬口河段, 平衡期階段, 其秋季潮波振幅梯度極小值最小, 為-3.55×10-5m-1,冬季值最大, 為-3.52×10-5m-1。雖然冬季余水位減小、底床摩擦增大, 但由于流量減小, 對潮能的耗散作用亦減弱, 潮波振幅梯度反而增大, 表明在年內(nèi)變化中甘竹-馬口河段主要受徑流影響。而在調(diào)整期階段, 冬季振幅閾值效應(yīng)消失, 這是由于冬季衰減效應(yīng)減小, 尚未達到閾值(約為-3.23×10-5m-1); 其余季節(jié)均存在閾值效應(yīng), 表明要在更大的外海潮汐動力驅(qū)動下才能使該河段發(fā)生閾值效應(yīng), 這與強人類活動驅(qū)動下各站點(三灶站除外)振幅均有不同程度的抬升有關(guān), 即兩個站點間振幅差值需要達到臨界值, 潮波振幅梯度極小值才會出現(xiàn)。

    圖5 三灶-馬口河段M2 分潮日均潮波振幅梯度(δ)與三灶站日均M2 分潮振幅(ηSZ)的閾值關(guān)系(a)及其季節(jié)變化(b~e)Fig.5 Relationship between the daily averaged M2 tidal amplitude gradient (δ) in SZ-MK sector and the M2 tidal amplitude (ηSZ)observed at the SZ tidal gauging station to illustrate the threshold effect (a) and its seasonal alteration (b~e)

    圖6 相鄰站點間河段M2 分潮日均潮波振幅梯度(δ)與三灶站日均振幅(ηSZ)關(guān)系曲線的季節(jié)變化Fig.6 Seasonal alteration in the relationship between the daily averaged M2 tidal amplitude gradient of the sub-sector along neighboring tidal gauging stations and tidal amplitude observed in SZ tidal gauging station

    對于下游河段(竹銀以下), 該區(qū)域主要受口門圍墾、河道整治工程等地形邊界影響, 對徑流非線性作用的響應(yīng)不明顯, 潮波振幅梯度隨口門振幅的變化基本呈線性關(guān)系, 而燈籠山-竹銀河段出現(xiàn)閾值效應(yīng)(平衡期的夏季)主要是由于夏季衰減效應(yīng)已達到調(diào)整期階段的閾值(約為-3.92×10-5m-1), 而平衡期振幅梯度極值更小, 潮波衰減尚未達到該臨界值。竹銀-甘竹河段閾值效應(yīng)出現(xiàn)的原因與甘竹-馬口河段類似。

    表 4 為三灶-馬口河段及相鄰站點間河段的潮波振幅梯度閾值及其對應(yīng)的振幅閾值。沿程潮波振幅梯度極小值變化規(guī)律與該河段M2分潮潮波振幅梯度的空間變化規(guī)律相似, 竹銀-甘竹河段衰減效應(yīng)減弱,到甘竹-馬口河段衰減效應(yīng)又增大, 表明甘竹站以上河段受強烈徑流調(diào)制, 潮波衰減效應(yīng)增強。在上游河段振幅沿程減小, 若要達到潮波振幅梯度閾值已不需要更大的口門振幅, 因此, 對應(yīng)振幅閾值沿程略有減小。

    表4 三灶-馬口河段M2 分潮潮波振幅梯度極小值(δM)及其對應(yīng)的振幅閾值(ηM)的季節(jié)變化Tab.4 Seasonal alterations in the minimum of daily averaged M2 tidal amplitude gradient (δM) in SZ-MK sector and the critical tidal amplitude (ηM) observed at SZ tidal gauging station

    4 討論

    上述結(jié)果表明, 磨刀門河口流量-余水位及其梯度關(guān)系在調(diào)整期階段發(fā)生顯著變化, 同流量下沿程各站點余水位及其梯度明顯下降(張蔚等, 2008; Yanget al, 2020)。然而三角洲頂端的馬口水文站年均流量僅減小816 m3/s (約為11%), 下游口門處三灶站M2分潮振幅年均值僅變化0.01 m, 表明上下游動力邊界的階段性變化不是導(dǎo)致潮波振幅梯度與上下游關(guān)系異變的主要原因, 地形異變及其導(dǎo)致的底床摩擦效應(yīng)變化才是導(dǎo)致徑潮動力格局轉(zhuǎn)換的主導(dǎo)因素。

    余水位梯度作為表征徑潮動力的重要參數(shù), 其形成演變是探究地形異變對潮波傳播影響機制的重要依據(jù)和有效切入點, 其解析表達式可由一維動量守恒方程得到(Savenijeet al, 2005)

    式中,x為距離口門處的距離;ρ為水體密度;U為斷面平均流速;Z為自由水面高程;D為水深;g為重力加速度;K為曼寧摩擦系數(shù)的倒數(shù)。忽略密度效應(yīng)和對流加速度項并假定流速具有周期性, 式(15)在一個潮周期內(nèi)積分可得式(16)(Caiet al, 2014)。

    式(16)為余水位梯度的解析表達式, 其中Z為一個潮周期內(nèi)的余水位。表明余水位梯度主要與有效摩擦項相平衡(Buschmanet al, 2009; Sassiet al, 2013), 即余水位梯度越大, 潮波所受有效摩擦越大。圖7 為磨刀門河口全河段及各河段M2分潮日均潮波振幅梯度和日均余水位梯度的關(guān)系曲線。三灶-馬口河段的曲線均為下凹形, 表明余水位梯度增大過程中, 潮波振幅梯度減小速度減緩并趨近于其極小值, 達到極小值后變大, 表明潮波振幅梯度的變化已不受底床摩擦主導(dǎo)。竹銀以下河段曲線為上凸形, 表明余水位梯度增大時有效摩擦增大, 在該區(qū)域中, 靠近口門的三灶-燈籠山河段(圖7b)在調(diào)整期變化速率更快, 該變化主要由口門圍墾導(dǎo)致的水深變淺主導(dǎo), 而燈籠山-竹銀河段(圖7c)則受采砂和圍墾的共同影響。竹銀-甘竹河段(圖7d)曲線由平衡期的下凹轉(zhuǎn)變?yōu)檎{(diào)整期的上凸形, 根據(jù)該河段余水位梯度和潮波振幅梯度曲線,相同的余水位梯度條件下, 潮波振幅梯度增大, 表明有效摩擦對其影響減弱, 這與航道疏浚后水深增大有關(guān)。甘竹-馬口河段(圖7e)曲線為下凹形, 這表明影響該河段潮波傳播的多種因素耦合作用增強, 調(diào)整期余水位梯度顯著減小而水深增大, 當兩種效應(yīng)的耦合作用達到平衡時即達到潮波振幅梯度閾值, 因此調(diào)整期階段潮波振幅梯度閾值大于平衡期。

    圖7 相鄰站點間河段M2 分潮日均潮波振幅梯度(δ)與日均余水位梯度(S)的變化Fig.7 Stepwise alteration in the relationship between the daily averaged M2 tidal amplitude gradient (δ) in SZ-MK sector and the corresponding residual water level gradient (S) observed in each sub-sector

    余水位梯度的階段性變化是水面線適應(yīng)地形變化的直接反映。在人類活動干預(yù)前, 河床緩慢抬升,底坡降變化不大。1962~1977 年為緩慢淤積狀態(tài), 磨刀門河口竹洲頭-大排沙河段、大排沙-燈籠山河段平均水深分別減小0.49, 0.28 m; 1977~1999 年, 水道受到明顯沖刷, 河道深泓下切深度為0.59~1.70 m 不等,平均水深增加1.13 m (劉鋒等, 2011); 磨刀門河口沿程平均寬深比從1960 年的6.25 減小到1999 年的4.73(河道趨向窄深化趨勢發(fā)展), 同期河槽容積增大34.15×106m3(Liuet al, 2019)。表5 統(tǒng)計了三灶-竹銀河段1964、1977、1999、2016 年各水深區(qū)間的容積和面積。按水深分成0~2、2~5、5~10 和>10 m 四個區(qū)間, 分別計算各區(qū)間的河槽容積和面積。1964~1977 年, 各水深段容積和面積均顯著減小, 反映以圍墾為主導(dǎo)的淤積作用。1977~1999 年, 口門圍墾和采砂共存, 前者導(dǎo)致口門大部分水域被圍填成陸地, 0~2 m 水深區(qū)間容積和面積均銳減, 然而水深5 m 以下容積和面積顯著增大, 表明河床大幅下切,兩岸侵蝕效應(yīng)增強, 河道被沖刷, 容積和面積增加量分別為每年3.78×106m3和0.90×106m2。1999~2016年, 強人類活動逐漸減弱, 但潮波傳播過程已發(fā)生異變, 加上上游水庫建設(shè)導(dǎo)致泥沙減少, 導(dǎo)致河道水沙分配格局發(fā)生變化, 潮汐動力顯著增強, 沖刷效應(yīng)進一步加劇, 5 m 以下水深容積持續(xù)增大, 增加量為每年1.06×106m3, 面積僅增大2.93×106m2。當沖淤態(tài)勢由淤積轉(zhuǎn)變?yōu)闆_刷時, 河道容積增大, 底床摩擦減小, 潮波振幅梯度增大。

    表5 三灶-竹銀河段1964、1977、1999、2016 年各水深區(qū)間的河槽容積和面積對比Tab.5 Comparison in water volume and surface area of the channel in different water depths in SZ-ZY sector observed in 1964, 1977,1999, 2016

    5 結(jié)論

    本文基于珠江磨刀門河口沿程4 個潮位站的高、低潮數(shù)據(jù)及相應(yīng)時段馬口水文站高、低潮數(shù)據(jù)和流量數(shù)據(jù), 采用基于流量驅(qū)動的R_TIDE 數(shù)據(jù)驅(qū)動模型對河口潮波傳播的異變過程進行分析, 根據(jù)余水位及其梯度、M2分潮振幅及其梯度的階段性變化重點探討潮波振幅梯度與上下游動力邊界的關(guān)系異變過程及機制, 主要結(jié)論如下:

    (1) 與平衡期相比, 在強人類活動影響后的調(diào)整期階段, 各站點日均M2分潮振幅及其梯度增大, 其中振幅在秋季變化最大(燈籠山-馬口分別為0.06、0.10、0.10、0.13 m), 冬季最小(馬口除外, 燈籠山、竹銀、甘竹分別為0.05、0.08、0.07 m), 年內(nèi)差異不大; 潮波振幅梯度變化則較為復(fù)雜, 三灶-甘竹河段沿程減小, 甘竹-馬口河段階段性變化增大。

    (2) 日均M2分潮潮波振幅梯度與下游動力邊界之間存在閾值效應(yīng), 主要出現(xiàn)在甘竹-馬口河段。強人類活動驅(qū)動下潮波振幅梯度閾值增大(三灶-馬口河段增量約為2.50×10-7m-1), 振幅閾值基本不變(僅增大約0.14 cm)。由于底床有效摩擦減小導(dǎo)致沿程潮能耗散減弱, 因此當口門振幅相同時, 調(diào)整期的潮波振幅梯度大于平衡期, 隨著振幅增大, 底床摩擦效應(yīng)越明顯, 導(dǎo)致調(diào)整期潮波衰減速度減緩, 當兩個階段的口門振幅均達到振幅閾值時(此處假定振幅閾值相同),則調(diào)整期對應(yīng)的潮波振幅梯度大于平衡期。

    (3) 對于潮波振幅梯度與大小潮的關(guān)系異變分析需綜合考慮地形、摩擦、流量等多因素影響。由于上下游動力邊界的階段性變化不足以導(dǎo)致上述關(guān)系的異變, 因此地形及其引起的底床摩擦變化則是徑潮動力異變的根本原因。從1964~2016 年, 10 m 以下水深的容積增大15.26×106m3, 0~2 m 水深容積則減小366.94×106m3, 反映河道窄深化的演變趨勢。由于地形的異變, 在口門振幅增大過程中, 達到振幅閾值前,由于底床摩擦增大, 小潮時的潮波振幅梯度大于大潮, 而達到振幅閾值后, 此時以地形輻聚為主導(dǎo)影響因子導(dǎo)致潮波振幅梯度增大。

    猜你喜歡
    分潮馬口河口
    馬口魚
    垂釣(2023年11期)2024-01-21 16:07:04
    大亞灣雙峰水位的形成條件及準調(diào)和分量應(yīng)用的分析
    山東鄰海長周期分潮對深度基準面的影響分析
    非遺視野下湖北馬口窯的保護、繼承與開發(fā)
    流行色(2019年5期)2019-12-13 18:38:50
    來自馬口窯的對話
    ——馬口窯文獻與當代陶藝創(chuàng)作研究展
    馬口煤礦分層開采工作面上覆采空區(qū)自燃治理技術(shù)研究
    中國煤炭(2016年9期)2016-06-15 20:29:53
    他們?yōu)槭裁催x擇河口
    河口,我們的家
    特殊的河口水
    河口
    伦理电影大哥的女人| 欧美少妇被猛烈插入视频| 天堂中文最新版在线下载| 日韩,欧美,国产一区二区三区| 91成人精品电影| 国产亚洲精品久久久com| 国产精品久久久久久av不卡| 久久精品久久精品一区二区三区| 国产精品免费大片| 亚洲av男天堂| 日本午夜av视频| 伊人亚洲综合成人网| 成年人免费黄色播放视频| 日韩在线高清观看一区二区三区| 国产淫语在线视频| 日韩制服丝袜自拍偷拍| 老熟女久久久| av在线老鸭窝| 激情视频va一区二区三区| 国语对白做爰xxxⅹ性视频网站| 亚洲 欧美一区二区三区| 亚洲欧洲日产国产| 亚洲综合色网址| 免费播放大片免费观看视频在线观看| 亚洲色图 男人天堂 中文字幕 | 久久久久国产网址| 考比视频在线观看| 人妻 亚洲 视频| 日韩精品免费视频一区二区三区 | 亚洲av日韩在线播放| 黑人高潮一二区| 汤姆久久久久久久影院中文字幕| 亚洲av电影在线进入| 亚洲国产日韩一区二区| 女性生殖器流出的白浆| 插逼视频在线观看| 永久网站在线| 天天躁夜夜躁狠狠久久av| 亚洲美女黄色视频免费看| 成年美女黄网站色视频大全免费| 91午夜精品亚洲一区二区三区| 日韩人妻精品一区2区三区| 永久免费av网站大全| 免费高清在线观看日韩| 中国美白少妇内射xxxbb| 中国美白少妇内射xxxbb| 在线天堂中文资源库| 亚洲成av片中文字幕在线观看 | 一本大道久久a久久精品| 中文乱码字字幕精品一区二区三区| 日韩欧美精品免费久久| 一级毛片我不卡| 国产成人免费无遮挡视频| 国产免费视频播放在线视频| 校园人妻丝袜中文字幕| 免费大片18禁| 亚洲情色 制服丝袜| 99热网站在线观看| 91国产中文字幕| 久久久久视频综合| 男女国产视频网站| 国产精品女同一区二区软件| 精品少妇内射三级| 性色avwww在线观看| 日韩欧美精品免费久久| 成人无遮挡网站| 欧美 亚洲 国产 日韩一| av国产久精品久网站免费入址| 男男h啪啪无遮挡| 性色avwww在线观看| 国产 精品1| 午夜精品国产一区二区电影| 热re99久久精品国产66热6| 内地一区二区视频在线| 成人二区视频| 晚上一个人看的免费电影| av一本久久久久| 久久精品久久久久久噜噜老黄| 精品一区二区免费观看| 国产熟女午夜一区二区三区| 亚洲欧洲国产日韩| 一级毛片 在线播放| 少妇人妻久久综合中文| 自拍欧美九色日韩亚洲蝌蚪91| 国产精品欧美亚洲77777| 色哟哟·www| 伦理电影免费视频| 黑人欧美特级aaaaaa片| 天堂8中文在线网| 久久99热这里只频精品6学生| 两个人看的免费小视频| 天堂8中文在线网| 国产免费一区二区三区四区乱码| 高清不卡的av网站| 午夜福利视频精品| 日本欧美国产在线视频| 欧美日韩av久久| 999精品在线视频| 亚洲伊人色综图| 午夜视频国产福利| 亚洲国产精品999| 国产高清不卡午夜福利| 国产成人精品久久久久久| 亚洲激情五月婷婷啪啪| 在线观看免费日韩欧美大片| 中文字幕人妻熟女乱码| 亚洲国产av影院在线观看| 午夜激情久久久久久久| 丰满迷人的少妇在线观看| 中文字幕av电影在线播放| 99九九在线精品视频| 热99国产精品久久久久久7| 在线 av 中文字幕| 国产精品一区www在线观看| 免费黄色在线免费观看| 秋霞在线观看毛片| 久久精品国产亚洲av天美| 成人综合一区亚洲| 国产精品国产av在线观看| 亚洲国产日韩一区二区| 国产免费福利视频在线观看| 日本av免费视频播放| 人人妻人人澡人人看| 亚洲精品色激情综合| 国产亚洲最大av| 久久影院123| av国产久精品久网站免费入址| 五月天丁香电影| 久久精品久久久久久久性| 黄色 视频免费看| videos熟女内射| 免费看不卡的av| 国产男女超爽视频在线观看| 国产精品 国内视频| av免费观看日本| 最近最新中文字幕免费大全7| 夫妻午夜视频| 一二三四中文在线观看免费高清| 亚洲四区av| 一二三四在线观看免费中文在 | 亚洲,欧美,日韩| 亚洲成色77777| 涩涩av久久男人的天堂| 久久久久久久久久久免费av| 我的女老师完整版在线观看| 最近2019中文字幕mv第一页| 午夜福利网站1000一区二区三区| 国产一区二区三区综合在线观看 | 看免费av毛片| 一级毛片我不卡| 婷婷成人精品国产| 新久久久久国产一级毛片| 国产成人午夜福利电影在线观看| 人人妻人人澡人人爽人人夜夜| 男女啪啪激烈高潮av片| 国产精品秋霞免费鲁丝片| 日韩电影二区| 成人综合一区亚洲| 欧美日本中文国产一区发布| 国产有黄有色有爽视频| 亚洲精品久久久久久婷婷小说| 欧美bdsm另类| 欧美另类一区| 午夜老司机福利剧场| 51国产日韩欧美| 国产精品久久久av美女十八| 巨乳人妻的诱惑在线观看| 男女午夜视频在线观看 | 寂寞人妻少妇视频99o| 成人手机av| 日韩精品有码人妻一区| 搡老乐熟女国产| 嫩草影院入口| 久久人人爽人人片av| 亚洲精品一区蜜桃| 久久久精品免费免费高清| 免费看光身美女| 在线看a的网站| 99香蕉大伊视频| 欧美少妇被猛烈插入视频| 在线观看免费高清a一片| www日本在线高清视频| 九九爱精品视频在线观看| 国产av一区二区精品久久| 日韩熟女老妇一区二区性免费视频| 亚洲国产精品一区二区三区在线| 五月伊人婷婷丁香| 97在线视频观看| 毛片一级片免费看久久久久| 日韩 亚洲 欧美在线| 日韩视频在线欧美| 国产精品.久久久| av播播在线观看一区| 日韩av免费高清视频| 色视频在线一区二区三区| 黄色视频在线播放观看不卡| 丁香六月天网| 国产伦理片在线播放av一区| 国产成人精品无人区| 爱豆传媒免费全集在线观看| 日日爽夜夜爽网站| 99热6这里只有精品| 国产亚洲精品第一综合不卡 | 一级a做视频免费观看| 成人无遮挡网站| 亚洲国产成人一精品久久久| 80岁老熟妇乱子伦牲交| 黑丝袜美女国产一区| 两性夫妻黄色片 | 丝袜美足系列| 韩国av在线不卡| 丝袜喷水一区| 亚洲精品美女久久久久99蜜臀 | 日本黄大片高清| 国产精品一国产av| 丰满乱子伦码专区| 亚洲四区av| 午夜福利视频在线观看免费| 国产成人一区二区在线| 永久免费av网站大全| 日韩制服骚丝袜av| 热re99久久国产66热| 久久久久久久国产电影| 亚洲精品一二三| 欧美精品一区二区大全| 日韩制服骚丝袜av| 久久精品久久精品一区二区三区| 亚洲精品第二区| 美女中出高潮动态图| 在线天堂最新版资源| 国产精品一国产av| 成人影院久久| 最近的中文字幕免费完整| 欧美97在线视频| 肉色欧美久久久久久久蜜桃| 国产一区二区在线观看av| 综合色丁香网| 少妇 在线观看| 日韩一区二区视频免费看| 日本-黄色视频高清免费观看| 女人精品久久久久毛片| 91精品伊人久久大香线蕉| 久热久热在线精品观看| 人人妻人人爽人人添夜夜欢视频| 日韩av在线免费看完整版不卡| 777米奇影视久久| 人妻人人澡人人爽人人| 久久精品国产亚洲av天美| 欧美+日韩+精品| 久久久亚洲精品成人影院| 国产男人的电影天堂91| 国产免费一区二区三区四区乱码| 热re99久久精品国产66热6| 两个人免费观看高清视频| 免费观看在线日韩| 国产麻豆69| 丝袜美足系列| 久久久久久人妻| 久久久久久久久久久免费av| 成人黄色视频免费在线看| 狂野欧美激情性xxxx在线观看| 美国免费a级毛片| 熟妇人妻不卡中文字幕| 人妻一区二区av| 亚洲精品视频女| 有码 亚洲区| 国产精品久久久久成人av| 乱人伦中国视频| 久久久久精品人妻al黑| 亚洲精品视频女| 久久人人爽av亚洲精品天堂| 一边摸一边做爽爽视频免费| 26uuu在线亚洲综合色| 国产精品一区二区在线不卡| 久久精品夜色国产| 国产色婷婷99| 欧美xxⅹ黑人| 国产欧美另类精品又又久久亚洲欧美| 黑人欧美特级aaaaaa片| 亚洲国产精品一区三区| 国产精品熟女久久久久浪| 国产日韩欧美在线精品| 波野结衣二区三区在线| 日韩中字成人| 色视频在线一区二区三区| 日韩一区二区三区影片| 日韩电影二区| 中文天堂在线官网| 精品久久蜜臀av无| av不卡在线播放| 亚洲国产毛片av蜜桃av| 黑丝袜美女国产一区| 色婷婷av一区二区三区视频| 国产av国产精品国产| 亚洲国产欧美日韩在线播放| 欧美激情国产日韩精品一区| av在线观看视频网站免费| 成人影院久久| 国产精品99久久99久久久不卡 | 看非洲黑人一级黄片| 一区二区三区四区激情视频| 成人二区视频| 看十八女毛片水多多多| 精品一区二区三区四区五区乱码 | 中国三级夫妇交换| 美国免费a级毛片| 国产高清三级在线| 插逼视频在线观看| 黑人高潮一二区| 亚洲四区av| 午夜视频国产福利| 亚洲精品日韩在线中文字幕| 青春草亚洲视频在线观看| 免费看av在线观看网站| 欧美激情极品国产一区二区三区 | 人妻系列 视频| 国产欧美亚洲国产| 最新中文字幕久久久久| 欧美日韩国产mv在线观看视频| 新久久久久国产一级毛片| 国产精品久久久久久久电影| 久久久久视频综合| 五月天丁香电影| 成人午夜精彩视频在线观看| av网站免费在线观看视频| 少妇的逼好多水| 国产 精品1| 国产精品人妻久久久久久| 欧美日韩视频精品一区| 看免费av毛片| av在线app专区| 男的添女的下面高潮视频| 亚洲综合色网址| 建设人人有责人人尽责人人享有的| 人妻一区二区av| 99热国产这里只有精品6| 又黄又粗又硬又大视频| 精品亚洲成a人片在线观看| 成人毛片60女人毛片免费| 搡老乐熟女国产| 国产成人欧美| 波野结衣二区三区在线| 久久人人爽人人爽人人片va| 九色亚洲精品在线播放| 国产不卡av网站在线观看| 如日韩欧美国产精品一区二区三区| 久久精品久久精品一区二区三区| 久久精品国产亚洲av天美| 成年美女黄网站色视频大全免费| 2018国产大陆天天弄谢| 一二三四在线观看免费中文在 | 欧美精品人与动牲交sv欧美| 免费av中文字幕在线| 久久久精品区二区三区| 国产成人精品福利久久| 成人黄色视频免费在线看| 午夜老司机福利剧场| 日韩欧美一区视频在线观看| 男女下面插进去视频免费观看 | 午夜福利视频在线观看免费| 制服丝袜香蕉在线| 丰满迷人的少妇在线观看| 99久久人妻综合| 亚洲欧美精品自产自拍| 欧美精品av麻豆av| 激情视频va一区二区三区| 性高湖久久久久久久久免费观看| 亚洲国产av新网站| 国产成人aa在线观看| 蜜桃在线观看..| 国产熟女午夜一区二区三区| 国产在视频线精品| 欧美变态另类bdsm刘玥| 日本色播在线视频| 咕卡用的链子| 成人午夜精彩视频在线观看| 在线观看免费日韩欧美大片| 国精品久久久久久国模美| 男人舔女人的私密视频| 在线天堂最新版资源| 2018国产大陆天天弄谢| 五月天丁香电影| 精品一品国产午夜福利视频| 蜜臀久久99精品久久宅男| 国产免费视频播放在线视频| tube8黄色片| 永久免费av网站大全| 日本av手机在线免费观看| 精品亚洲成国产av| 久久久亚洲精品成人影院| av国产精品久久久久影院| 美女中出高潮动态图| 1024视频免费在线观看| 国产免费视频播放在线视频| 国产成人精品一,二区| 国产精品一国产av| 一本久久精品| kizo精华| 又大又黄又爽视频免费| 亚洲精品日本国产第一区| 国产一区二区激情短视频 | 欧美日本中文国产一区发布| 国产精品国产三级国产av玫瑰| 欧美成人精品欧美一级黄| 国产av码专区亚洲av| 中国国产av一级| 久久久久久人妻| 边亲边吃奶的免费视频| 国产成人av激情在线播放| 99久久人妻综合| 日本-黄色视频高清免费观看| 国产成人av激情在线播放| 亚洲伊人色综图| av网站免费在线观看视频| 亚洲精品久久午夜乱码| 国产色婷婷99| 久久99蜜桃精品久久| 国语对白做爰xxxⅹ性视频网站| 久久久欧美国产精品| 午夜av观看不卡| 国产精品熟女久久久久浪| 不卡视频在线观看欧美| 中文精品一卡2卡3卡4更新| 免费看av在线观看网站| 中文字幕另类日韩欧美亚洲嫩草| www.熟女人妻精品国产 | 两性夫妻黄色片 | 涩涩av久久男人的天堂| 2022亚洲国产成人精品| 中文字幕精品免费在线观看视频 | 在线观看国产h片| 丁香六月天网| 免费观看a级毛片全部| 下体分泌物呈黄色| 久久久国产精品麻豆| 久久影院123| 亚洲欧美色中文字幕在线| 免费观看a级毛片全部| 男人爽女人下面视频在线观看| 99re6热这里在线精品视频| 最近最新中文字幕免费大全7| 在线观看www视频免费| 亚洲欧美精品自产自拍| 日日摸夜夜添夜夜爱| 80岁老熟妇乱子伦牲交| 国产一区二区在线观看av| 国产乱来视频区| 午夜福利在线观看免费完整高清在| 啦啦啦中文免费视频观看日本| 青春草视频在线免费观看| 又大又黄又爽视频免费| 99热这里只有是精品在线观看| 高清黄色对白视频在线免费看| 人人妻人人添人人爽欧美一区卜| 成年动漫av网址| 精品国产一区二区久久| 亚洲五月色婷婷综合| 在线天堂最新版资源| 黄色怎么调成土黄色| 成年动漫av网址| 日韩大片免费观看网站| 久久人人97超碰香蕉20202| 日韩av不卡免费在线播放| 久久韩国三级中文字幕| 国产不卡av网站在线观看| 777米奇影视久久| 男人舔女人的私密视频| av在线观看视频网站免费| 色婷婷久久久亚洲欧美| 人人妻人人澡人人爽人人夜夜| 免费大片18禁| 精品人妻在线不人妻| 精品少妇内射三级| 亚洲五月色婷婷综合| 欧美精品人与动牲交sv欧美| 男男h啪啪无遮挡| 欧美激情 高清一区二区三区| 日韩大片免费观看网站| 晚上一个人看的免费电影| 免费高清在线观看视频在线观看| 丝袜脚勾引网站| 国产熟女欧美一区二区| 亚洲欧美日韩卡通动漫| 考比视频在线观看| 一本色道久久久久久精品综合| 欧美丝袜亚洲另类| 亚洲精品视频女| 精品少妇内射三级| 亚洲,欧美精品.| xxx大片免费视频| 乱人伦中国视频| 两性夫妻黄色片 | 国产精品久久久av美女十八| 成人黄色视频免费在线看| 国产黄色视频一区二区在线观看| 免费av不卡在线播放| 9色porny在线观看| 亚洲人成网站在线观看播放| xxx大片免费视频| 97在线人人人人妻| 亚洲熟女精品中文字幕| 免费黄色在线免费观看| 久久人妻熟女aⅴ| 9191精品国产免费久久| 亚洲性久久影院| 狂野欧美激情性xxxx在线观看| 午夜福利乱码中文字幕| 免费不卡的大黄色大毛片视频在线观看| 自线自在国产av| 一区在线观看完整版| 久久精品久久久久久久性| 亚洲激情五月婷婷啪啪| www日本在线高清视频| 国产免费现黄频在线看| 国产综合精华液| 在线精品无人区一区二区三| 夫妻午夜视频| 美女中出高潮动态图| 爱豆传媒免费全集在线观看| 亚洲伊人色综图| 亚洲欧洲国产日韩| 国产爽快片一区二区三区| 大香蕉久久成人网| 人妻少妇偷人精品九色| 肉色欧美久久久久久久蜜桃| av国产久精品久网站免费入址| 成年人免费黄色播放视频| 久久久久久久精品精品| 欧美+日韩+精品| 久久影院123| 少妇被粗大的猛进出69影院 | 人妻人人澡人人爽人人| 国产免费现黄频在线看| 精品少妇内射三级| 99视频精品全部免费 在线| 搡老乐熟女国产| 久久99蜜桃精品久久| 最近的中文字幕免费完整| 精品福利永久在线观看| 一本大道久久a久久精品| 国产精品.久久久| 一二三四中文在线观看免费高清| 亚洲在久久综合| 色视频在线一区二区三区| 亚洲成av片中文字幕在线观看 | av播播在线观看一区| 婷婷色av中文字幕| 日本与韩国留学比较| 热re99久久精品国产66热6| 亚洲天堂av无毛| 香蕉精品网在线| 欧美bdsm另类| 免费黄网站久久成人精品| 久久精品久久久久久噜噜老黄| 国产精品蜜桃在线观看| 欧美bdsm另类| 日韩制服丝袜自拍偷拍| 亚洲av成人精品一二三区| 自线自在国产av| 亚洲伊人久久精品综合| 男女午夜视频在线观看 | 自线自在国产av| 极品人妻少妇av视频| 中文字幕制服av| 91aial.com中文字幕在线观看| 欧美丝袜亚洲另类| 91aial.com中文字幕在线观看| 亚洲av综合色区一区| 观看av在线不卡| 啦啦啦啦在线视频资源| 国产精品偷伦视频观看了| 热re99久久精品国产66热6| 午夜av观看不卡| 久久久久久久国产电影| 菩萨蛮人人尽说江南好唐韦庄| av线在线观看网站| 一本—道久久a久久精品蜜桃钙片| 日本91视频免费播放| 啦啦啦在线观看免费高清www| 国产又爽黄色视频| 大片电影免费在线观看免费| 日韩一区二区视频免费看| 久久影院123| 亚洲精品美女久久久久99蜜臀 | 一级毛片我不卡| 久久人妻熟女aⅴ| 在现免费观看毛片| 日韩av不卡免费在线播放| 麻豆精品久久久久久蜜桃| 制服诱惑二区| 男人舔女人的私密视频| av免费在线看不卡| 男女无遮挡免费网站观看| 亚洲欧洲日产国产| 日韩成人伦理影院| 久久精品久久精品一区二区三区| 日韩熟女老妇一区二区性免费视频| 18禁在线无遮挡免费观看视频| 日本-黄色视频高清免费观看| 亚洲欧美中文字幕日韩二区| 国产精品国产三级专区第一集| 999精品在线视频| 日日撸夜夜添| 青春草视频在线免费观看| 亚洲欧洲日产国产| 人人妻人人爽人人添夜夜欢视频| 一区二区三区四区激情视频| 欧美最新免费一区二区三区| 夜夜骑夜夜射夜夜干| 天天操日日干夜夜撸| 三级国产精品片| 免费在线观看黄色视频的| 涩涩av久久男人的天堂| 国产亚洲欧美精品永久| 亚洲欧美成人精品一区二区|