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

基于CEEMDAN-LSTM模型的大伙房水庫入庫流量中長期預報

2025-07-16 00:00:00王春羽張靜楊旭閆濱
人民珠江 2025年6期
關鍵詞:模型

Medium- and Long-term Forecast of Inflow of Dahuofang Reservoir Based on CEEMDANLSTMModel

WANGChunyu',ZHANGJing1*,YANGXu2,YANBin (1.CollgeofWaterResoures,ShenyangAgriculturalUniversityhenyangO866,China; 2.LiaoningDahuofangReir AdministrationLimitedLiabilityCompany,Fushun113oo7,China)

Abstract:Astocsisbcttinfallmpatufasanctivityteedultilectosf complexityitslwofhngesossigificantadomssnducertaintyreatlyiceasingtdifultyfedi-ndlote forecastandrestrictingitsaplicationectinproductionpractice.Therefore,howtobreakthroughthetechnicalbotleneckof medium-andlong-termflowforecasthasbecomeakeyproblemtobesolvedinthecurrnthydrologicalscienceresearchInviewof this,withDahuofang Reservoirastheresearchobject,the monthlyflow forecastof thereservoirwascariedoutbyusing gray correlationanalysisandprincipalcomponentanalysiscombinedwiththeCEEMDAN-LSTMmodel.Monthlyflow,rainfalland meteorologicaldatafrom1961to2008wereselectedtocalibratethemodelparameters.The monthlyflowdatafrom2009to2O20was usedtovalidatetemodel.Theidicatorsofeterinationcoefient,ootmansquareeandaverageelativeerorreapplied to evaluate theforecastresults.Aording totheresults,whenfourprincipalcomponentsareselectedbyusing principalcomponent analysis toreducethedimensionalityof theforecastfactorsetsaftertheadditionoftheprevious average temperatureandmaximum temperaturedata,theCEEMDAN-LSTMmodelcanefectivelyimprovetheforecastaccuracyandbecometheoptimalmodelforthe monthlyflowforecastofDahuofangReservoir.Technicalsupportisprovidedfortheformulationofthefuturemedium-andlong-te dispatching plan of Dahuofang Reservoir.

Keywords:iflow; medium-andlong-termforecast;principalcomponentanalysis;CEEMDAN-LSTMmodel; DahuofangReservoir

中長期流量預報在水資源規劃與管理中具有重要的理論價值與實踐意義,其預測結果可為流域水文情勢分析、水庫調度方案優化以及水資源合理配置提供科學依據。然而,由于流量形成過程具有顯著的隨機特性,且流量時間序列呈現復雜的非線性動力學特征,導致現有預測模型的準確性難以滿足實際應用需求,這也使得該領域成為當前水文科學研究中的重點攻關方向[2-6]。根據模型原理的差異,當前中長期流量預測方法可劃分為兩大類:一類是物理驅動模型,該方法通過耦合長期氣象預報產品與水文模型來實現預測;另一類數據驅動模型[7]。氣象預報產品提供的長預見期數據(如月尺度或季節尺度的降水量、氣溫等)通常具有較高的不確定性和較低的空間分辨率[8-10],難以直接滿足中小流域或局部區域的高精度水文預測需求。此外,長預見期氣象預報的系統性偏差和誤差隨預見期延長而增加,進一步限制了其在水文預測中的實用性,為提高水文預測結果的可靠性與應用價值,往往需要將水文機理模型與多源觀測數據進行融合分析。數據驅動模型是通過構建輸入輸出變量之間的統計關聯模型,從歷史水文資料中挖掘潛在的降雨-流量響應規律,進而實現流量過程的預測模擬[2]。人工智能技術與機器學習算法的快速發展催生了多種新型數據驅動建模方法,這些創新性技術為中長期預報開辟了新的方向[3。其中,基于深度學習的預報方法表現出強大的非線性映射能力和時序特征提取能力,如長短期記憶網絡(LongShort-TermMemory,LSTM)、門控循環單元(Gated RecurrentUnit,GRU)等模型已在多個流域得到成功應用。胡慶芳等[14]基于LSTM模型構建的安康站日徑流預測模型展現出優異的泛化能力和抗過擬合特性。Ding等15]提出基于VMD(Variational ModeDecomposition)分解的CNN-LSTM模型,驗證了其在不同時間尺度的預報穩定性,充分地挖掘了流量序列的隱含信息。李薇等提出的VMD-CNN-LSTM模型在預測精度和誤差控制方面比其他模型表現出明顯優勢。

然而,盡管已有多種水文預報模型和方法,但研究仍存在以下亟待解決的問題:首先,預報因子的篩選和降維方法尚待優化,難以充分挖掘多源數據的有效信息,導致模型輸入信息的冗余與不足;其次,單一模型在處理復雜水文過程時往往存在局限性,難以全面捕捉徑流序列的多尺度特征;此外模型的可解釋性和泛化能力仍需進一步提升,依賴于超參數的設定,以確保預報結果在實際應用中的可靠性。

針對上述問題,本研究以大伙房水庫1961年1月至2020年12月逐月入庫流量為研究對象,選取1961—2008年各月流量、逐月降雨量、12個氣象資料進行模型參數率定,利用2009—2020年各月入庫流量資料對模型進行檢驗,提出結合主成分分析法、自適應噪聲完全集合經驗模態分解、長短期記憶網絡的月組合預報模型。該模型在預報因子篩選與降維基礎上[17-18],利用CEEMDAN對原始流量序列進行分解降噪,以有效提取水文序列中的非線性和非平穩特征;結合LSTM模型對各分解分量分別建模預測并將分量重構得到最終預報結果。通過對比分析不同主成分個數的預報因子集在模型中的預報效果,探明主成分分析法和組合模型在大伙房水庫的適用性及精度,選取最優預報因子集作為大伙房水庫逐月最優預報模型,為大伙房水庫中長期流量預測提供參考。

1數據與模型方法

1.1 數據來源與處理

為滿足日常上下游站點降雨量、水位和人工水文站流量的日常傳輸工作,大伙房水庫管理局通過歷年建設,建立了一套水文自動測報系統,本文采用的大伙房水庫1961年1月至2020年12月各月入庫流量、降雨量、平均 2min 風速、平均水氣壓、平均氣溫、平均相對濕度、日照時數、最低氣壓、最低氣溫、最高氣壓、最高氣溫、最大風速、最小相對濕度數據資料均來源于此水文系統。此時段數據無缺失,質量可靠。氣象數據和水文數據在時間上完全匹配,為研究提供了堅實基礎。

1. 2 主成分分析法

主成分分析(Principal ComponentAnalysis,PCA)作為一種多元統計技術,其核心在于通過線性變換將一組存在相關性的觀測變量轉化為少數幾個相互獨立的新變量,這些新變量成為主成分。在數據降維技術,確保原始水文特征信息完整性的基礎上,實現對多維數據集維度的有效簡化[19]。設數據集有 n 個樣本 ?m 個變量,基本操作步驟如下。

步驟一對原始變量數據進行標準化處理:

1.3 CEEMDAN原理

自適應噪聲完全集合經驗模態分解CEEMDAN(Complete Ensemble Empirical Mode DecompositionwithAdaptiveNoise)是基于經驗模態分解(Empirical ModeDecomposition,EMD)算法改進的一種算法,其原理是將原始數據序列分解為一系列的振蕩數據序列,同時根據式(7)各振蕩數據序列添加自適應白噪聲,確保重構誤差(R)盡可能接近于零[20]

式中:""為原始變量 的平均值; σj"為原始變量的標準差。

步驟二 所有的系數 l 構成載荷矩陣 L, 設矩陣"并計算 L 的特征值及特征

向量:

求得 L 的特征值 λ1≥λ2≥…≥λm",以及相對應的單位正交向量u,u,…,um

步驟三 確定主成分貢獻率 τj"和累積貢獻 率nj

步驟四 得到主成分表達式:

S(t表示流量時間序列, Vi(t) 為第 i 次試驗中添加的具有標準正態分布的高斯白噪聲序列,第 i 次信號序列可表示為式(7):

Si(t)=S(t)+ε0Vi(t),i∈{1,2,…,M}

式中: Si(t) 為分解后新的時間序列; ε0 為噪聲系數;ε0Vi(t) 為自適應高斯白噪聲; M 為集成次數,一般為10\~20。

其算法公式具體如下:

步驟一對分解得到的第一個IMF分量IMF1i(t) 取均值得到第一個模態分量:

步驟二 將分解出的第一個分量 IMF1(t) 從原始時間序列 S(t) 中移除,得到第一個余量序列 R(t)

R1(t)=S(t)-IMF1(n)

步驟三 與步驟一相同,繼續進行分解,得到第二個IMF模態分量:

步驟四重復步驟一、二,計算其余各分解階段對應的的余量 RK 和 k+1 個模態分量,得到剩余的IMF模態分量,其中 K 為IMF模態總個數:

Rk(t)=Rk-1(t)-IMFk(t),k=(2,3,…,K)

步驟五分解過程的終止判定標準設定為:當余量序列中的極值點數量不超過2個時,即判定該序列已達到不可再分解狀態。此時,不可分解余量序列被定義為殘余分量,其表達式為:

故原始時間序列 S(t) 經過CEEMDAN分解后可用式(14)表示:

1.4LSTM模型原理

LSTM模型包含輸入門、遺忘門和輸出門[21-23]。結構見圖2,計算過程見式(15)。

圖2LSTM神經網絡結構

1.5 CEEMDAN-LSTM模型

通常的LSTM模型單層無法抓取較多數據信息,需構建多層隱藏層結構。通過構建CEEMDAN-LSTM模型,即搭建CEEMDAN分解與2層LSTM結構,可以使預測模型能更多地挖掘數據信息的映射關系。模型應用時,首先對預報因子集進行數據的分解,得到由CEEMDAN分解后的多個模態分量IMFk 和一個殘差項,然后將每個 IMFk 的統計特征組合成多維特征向量,將其劃分為訓練集和測試集,搭建LSTM模型,每層LSTM結束時鏈接一層dropout層,dropout層使每層LSTM的神經元按照設定的百分比隨機失活,不更新權重參數。引人dropout層的目的是防止雙層的LSTM訓練時發生過擬合。輸出層選擇用全連接層激活函數relu,使用均方誤差(Mean-squareError,MSE)作為模型訓練過程中的損失函數,Adam作為優化器,將各個模態分量和趨勢項預測結果疊加,最終輸出預測數據。CEEMDAN-LSTM模型結構見圖3。

Fig.2StructureofLSTMneuralnetwork圖3CEEMDAN-LSTM模型流程 Fig.3ProcessofCEEMDAN-LSTMmodel

1.6模型性能評價指標

為有效分析實驗結果,本文采用決定性系數(R2) )、均方根誤差(RootMean Square Error,ERMS) 、平均相對誤差(Mean Relative Error, EMR )3種指標對預報效果進行評估。

式中: xi 為預測值; 為預測平均值; yi 為實際值; 為實際平均值。 R2 指標越接近1,說明模型擬合效果好,預測結果誤差較小;反之,模型效果差,預測結果誤差大。 ERMS 越小,模型精度越高。

2 實例分析

大伙房水庫坐落于渾河流域中上游段,是一座綜合性大型水利工程設施,其主要功能包括防洪減災、城市供水、農業灌溉及水力發電等,也是跨流域調水工程的重要調節樞紐,承擔著遼寧省中部7座城市的供水任務[24]。因此該水庫的流量中長期預報非常重要。水庫流域概況見圖4。

圖4大伙房水庫流域

2.1關鍵預報因子篩選

為盡可能大范圍篩選出對入庫流量預報影響較大的因子,本文選取預報月的前14個月為影響預報時效。采用灰色關聯度法計算14個月內各月降雨量、平均 2min 風速、平均水氣壓、平均氣溫、平均水氣壓、平均相對濕度、日照時數、最低氣壓、最低氣溫、最高氣壓、最高氣溫、最大風速、最小相對濕度13個指標與預報流量之間的灰色關聯系數,選取灰色關聯系數大于0.95的指標為主要的預報因子,篩選結果見表1。

表1主要的預報因子Tab.1Principal forecastfactors

為進一步消除多個預報因子之間的數據冗余及相關性,采用主成分分析法將選定的平均氣溫、最高氣溫預報因子進行線性組合,通過線性變換將原始變量轉化為彼此正交的主成分,從而實現數據降維和特征提取。當主成分的累積方差貢獻率超過 90% 時,即可認為該成分能夠充分表征原始數據的主要特征。通過計算,前4個主成分的累積解釋方差達到99. 186% ,表明其幾乎完全保留原始變量的有效信息。各主成分對應的特征值及其方差貢獻率分布見圖5。

Fig.4DahuofangReservoirBasin圖5主成分分析結果Fig.5PCAresults

各影響因素對主成分之間的關系見表2,主成分得分系數矩陣表示了各影響因素與提取前4個主成分的線性相關性,基于該矩陣的數值分布特征,可計算出第一至第四主成分的數學表達式。

表2成分矩陣Tab.2Elementmatrix

代人式(6得到主成分表達式:

F1=0. 34W1+0. 39W2+0. 34W3+0. 19W4+0. 34Z1+ (20號 0.39Z2+0. 34Z3+0. 19Z4-0. 01Z5-0. 2Z6+0. 003Z11+ (2號 0.2Z12 (204號 (18

F2=-0.37W1+0.01W2+0.38W3+0.65W4-0.37Z1+ 0. 01Z+0. 38Z+0.65Z4+0. 74Z5+0.63Z6+0. 73Z11- (2 0.64Z12 (19)

F3=0.07W1-0.051W2-0.04W+0.08W4+0.07Z(24號 0. 08Z2-0. 07Z3+0. 09Z4+0. 19Z5+0. 14Z6+0. 17Z11- (20 0.19Z12 (204號 (20)

F4=0.35W1+0.19W2-0.22W3-0.29W4+0.47Z1+ 0.24Z2-0.28Z3-0.41Z4+0.09Z5+0.7Z6-0.61Z11- (2號 0.01Z12 (204號 (21)式中: Fi(i=1,2,3,4) 表示主成分。

計算各個樣本的主成分得分 F1?F2?F3?F4 ,基于上述分析,分別選定前2個、前3個和前4個主成分與前 t-1?t-2?t-11?t-12 時刻的入庫流量以及前 t- 1、t-2、t-3、t-4 時刻的降雨量組合作為模型的輸入,探討不同主成分組合對預報效果的影響,組成不同預報因子集(由降雨流量數據組成的預報因子集用XPQ 表示,由前2個主成分與 XPQ 組成的預報因子集用 Xs 表示,由前3個主成分組成與 XPQ 的預報因子集用 XT 表示,由前4個主成分組成與 XPQ 的預報因子集用 XF 表示),輸入模型進行對比分析。

2.2 CEEMDAN分解

將前述篩選的預報因子集(以 XPQ 為例)采用CEEMDAN分解法進行分解,結果見圖6。

圖6大伙房水庫預報因子集 XPQ CEEMDAN分解 Fig. 6 XPQ CEEMDANdecompositionofforecastfactorsets ofDahuofangReservoir

2.3 CEEMDAN-LSTM模型的應用

本文基于Python語言,在Jupyter環境中搭建CEEMDAN-LSTM模型,該環境既可在本地運行,也可部署于云服務器,實現遠程瀏覽器訪問。采用Anaconda進行配置管理,便于創建特定Python版本和庫函數的環境,避免版本沖突,支持快速部署與復制。神經網絡框架選用Tensorflow,這一開源機器學習框架庫豐富、更新迅速。開發過程中主要調用的庫包括 Tensorflow、Keras、Numpy、Matplotlib、Sklearn和Pandas等。

在水文預報模型的構建過程中,超參數的選取對于模型性能的優劣至關重要。為了確保模型在徑流預測中的準確性和穩定性,進行不同超參數組合的運算實驗,并綜合考慮了模型計算的參數需求以及訓練集和驗證集的性能指標。超參數優化選擇以驗證集均方誤差為目標函數,選用adam優化算法,初始學習率為0.0001,中間層dropout為0.2。LSTM記憶單元內的門控權重矩陣及偏置項的初始化用TensorFlow框架的seed隨機種子來確定。模型參數經試錯后運算結果最優值多在epoch為400左右時取到,將epoch設置為400。通過調整神經元數量及批量(batch_size)大小,尋找不同的預見期最優參數,見表3。

表3CEEMDAN-LSTM模型參數確定過程Tab.3ParameterdeterminationofCEEMDAN-LSTMmodel

在不同超參數的運算結果為依據的情況下,充分考慮模型計算的參數需求,以及訓練集和驗證集的 R2,ERMS ,不同預報因子集最終選擇參數:雙層LSTM神經網絡,batch_size :=8 ,epoch 1=400 。這一組合在多個性能指標上均表現出色,能夠滿足水文預報模型的精度要求,具體超參數組合見表4。

表4CEEMDAN-LSTM預測模型超參數設置

預測結果見表5、圖7,CEEMDAN-LSTM測試集 預報流量見圖8。

表5評價指標Tab.5Evaluationindicators
圖7CEMMDAN-LSTM預測Fig.7CEMMDAN-LSTMforecast

從預報結果可以看出:R2xp=0.837,R2x=0.842, , ERMS XP0=78.80 !ERMS X5=78.60 , ERMS XT=73 .35, ERMS XF=70.39 , R2 增大, ERMS 減小,分別下降了 0.20,5.25,2.96 。這一結果表明當預報因子集僅包含前期降雨和前期入庫流量時, R2 相對來說較低 ?ERMS 較高;在引入前期平均氣溫和最高氣溫的數據后采用主成分分析法進行降維處理、并依次選取前2、3、4個主成分作為模型的輸入變量時,分析表明 R2 值逐漸增大, ERMS 逐漸減小,這充分說明主成分分析法在提高預報精度方面發揮了積極作用。

從平均相對誤差上來看, XPQ?XS?XT?XF 的 EMR 分別為 19. 13% 18. 64% 16. 81% 15. 63%,Xpq 在4個模型中 EMR 相對較高。這表明,僅依靠前期降雨和前期入庫流量作為預報因子時,模型的預報準確性仍有提升空間。 XF 預報因子集精度最高,相較于前3個預報因子集 EMR 降低了3.50、3.01、1.18個百分點,這一變化趨勢明顯表明氣溫作為影響蒸發、降雨等水文過程的重要因素,納入預報模型顯著改善預報結果。

從實際值與測試值分布情況來看, Xs,XT,XF 預測結果大部分緊密圍繞 45° 線附近,明顯優于 XPQ 預報結果, XF 預報因子集相較于其他預報因子集更加緊密分布于 45° 線附近,整體貼近度最高,模擬效果最為理想。

圖8CEEMDAN-LSTM測試集預報流量 Fig.8ForecastflowsofCEEMDAN-LSTMtestset

3結論

以大伙房水庫1961年1月至2020年12月逐月入庫流量為研究對象,在篩選影響各月入庫流量預報主要影響因素的基礎上,分析不同主成分個數對各月入庫流量預報結果的影響,評價基于CEEMDAN-LSTM模型的月流量預測效果。

a)通過對CEEMDAN-LSTM模型進行不同預報因子集模擬預測效果綜合對比發現,分解算法能夠有效增強模型預測能力,總體上來說各個預報因子集預報結果均表現良好,證明此模型能夠提高中長期預報成果精度和穩定性,在處理復雜的非線性水文時間序列數據時具有較強的適應性和預測能力。

b)在預報因子集的篩選方面,采用主成分分析法對平均氣溫、最高氣溫因素進行降維分析,在確保原始水文數據信息無損的前提下,顯著降低了數據冗余度。通過主成分分析提取的4個關鍵特征分量,有效增強了模型的預測性能,其優化效果在中長期徑流預測領域具有重要的工程應用價值。

c)XF 預報因子集在各項指標上均表現出色,可作為大伙房水庫徑流逐月最優預報模型因子集。該模型不僅考慮了前期降雨和入庫流量等水文因素,還融入前期平均氣溫、最高氣溫氣象因素,為大伙房水庫精準調度提供了有力支持。

參考文獻:

[1]朱永英.水庫中長期徑流預報及興利調度方式研究[D].大連:大連理工大學,2008.

[2]張芳琴.多方法融合數據驅動模型在月徑流預報中的應用研究[D].楊凌:西北農林科技大學,2022.

[3]康艷,楊沁瑜,張芳琴,等.基于Box-Cox變換與Lasso回歸的支持向量回歸模型在月徑流預報中的應用[J].應用基礎與工程科學學報,2022,30(1):27-39.

[4]包苑村,解建倉,羅軍剛.基于VMD-CNN-LSTM模型的渭河流域月徑流預測[J].西安理工大學學報,2021,37(1):1-8.

[5]李福興.瑪納斯河徑流演變特性及其中長期徑流預報模型研究[D].石河子:石河子大學,2021.

[6]效文靜,周建中,楊建華,等.基于機器學習的長江干流中長期徑流預報[J].水電能源科學,2022,40(9):31-34,26.

[7]邢馨予.基于機器學習的中長期徑流預報方法研究[D].太原:太原理工大學,2022

[8]楊程.基于可解釋機器學習的黃河上游青海段中長期徑流預報研究[D].西安:西安理工大學,2024.

[9]楊靜欣.基于分解方法與機器學習組合的中長期徑流預報模型研究[D].鄭州:華北水利水電大學,2023.

[10]李保健.水電站群中長期徑流預報及發電優化調度的智能方法應用研究[D].大連:大連理工大學,2015.

[11]朱雙.流域中長期水文預報與水資源承載力評價方法研究[D].武漢:華中科技大學,2017.

[12」陳胥良.耦合中長期徑流預報的跨流域引水供水聯合調度研究[D].大連:大連理工大學,2021.

[13]李佳,牟時宇,朱艷軍,等.基于多維信息相似度的大渡河流域月徑流量預測研究[J].水電能源科學,2022,40(12):38-41.

[14]胡慶芳,曹士圯,楊輝斌,等.漢江流域安康站日徑流預測的LSTM模型初步研究[J].地理科學進展,2020,39(4):636-642.

[15]DINGYK,ZHUYL,FENGJ,etal.Interpretable spatio-temporalattention LSTMmodel forflood forecasting[J].Neurocomputing,2020,403:348-359.

[16]李薇,周建中,葉磊,等.基于主成分分析的三種中長期預報模型在柘溪水庫的應用[J].水力發電,2016,42(9):17-21.

[17]KRATZERTF,KLOTZD,BRENNERC,etal.Rainfall-runoffmodelling usingLongShort-TermMemory(LSTM)networks[J].Hydrologyand Earth SystemSciences,2018,22(11):6005-6022.

[18]農振學,王超,雷曉輝.基于主成分分析和BP神經網絡的贛江流域中長期徑流預報[J].水電能源科學,2018,36(1):16-19.

[19]解雪,陳軍鋒,鄭秀清,等.基于主成分分析和粒子群算法優化支持向量機的凍融土壤蒸發預報模型[J].節水灌溉,2020(1):61-65,72.

[20]胡斯曼.基于CEEMDAN-LSTM模型的中長期徑流預報研究及系統集成[D].武漢:華中科技大學,2020.

[21]SONG TY,DING W,WUJ,et al.Flash Flood ForecastingBasedon Long Short-Term Memory Networks[J].Water,2020,12(1).D0I:10.3390/w12010109.

[22]HOCHREITER,S,SCHMIDHUBER J.Long Short-TermMemory[J].Neural Computation,1997,9(8):1735-1780.

[23」汪海,沈延青,祁善勝,等.基于STL-CEEMDAN-LSTM模型的月徑流量預測[J].人民珠江,2025,46(4):39-46.

[24]王秀紅,王建民.大伙房水庫輸水洞改造工程開挖方法比較與分析[C]//遼寧省水利學會2016年學術年會論文集,2016:41-43.

(責任編輯:向飛)

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 欧美日韩精品一区二区视频| 操操操综合网| 天天爽免费视频| 日韩无码视频播放| 国产一级无码不卡视频| 99这里只有精品免费视频| 亚洲AV电影不卡在线观看| 四虎永久在线| 日韩欧美中文字幕在线精品| 丝袜无码一区二区三区| 国产精品网址在线观看你懂的| 亚洲一级毛片在线播放| 欧美亚洲日韩中文| 亚洲第一中文字幕| 一区二区三区四区日韩| 黄色网站不卡无码| 国产精品无码制服丝袜| 成人av手机在线观看| 欧美成人一级| 日韩无码白| 亚洲国产精品日韩欧美一区| 无码久看视频| 综合人妻久久一区二区精品| 国产精品久久久久无码网站| 亚洲日韩精品欧美中文字幕| 国产本道久久一区二区三区| AV无码一区二区三区四区| 国产精选小视频在线观看| 成人伊人色一区二区三区| 亚洲欧美激情小说另类| 亚洲欧美一级一级a| 精品人妻无码中字系列| 永久成人无码激情视频免费| 久久免费精品琪琪| 色婷婷在线播放| 97国产精品视频自在拍| 国产福利微拍精品一区二区| 亚洲性网站| 午夜老司机永久免费看片| 亚洲伊人天堂| 日韩性网站| 中文无码伦av中文字幕| 久久国产毛片| 成人在线观看不卡| 亚洲男人的天堂在线观看| 成人在线观看不卡| 欧美狠狠干| 激情无码视频在线看| 国产乱人伦精品一区二区| 色综合久久无码网| 亚洲精品免费网站| 久久综合色天堂av| 成人精品区| 高清国产va日韩亚洲免费午夜电影| 久久a级片| 四虎永久在线| 暴力调教一区二区三区| 亚洲一级色| 日本午夜精品一本在线观看| 国产精品流白浆在线观看| 亚洲av日韩av制服丝袜| 亚洲av成人无码网站在线观看| 亚洲国产成人麻豆精品| 91网站国产| 亚洲三级色| 无码免费视频| 国内黄色精品| 毛片最新网址| 无码免费试看| 亚洲一区二区成人| 欧美另类图片视频无弹跳第一页| 久草热视频在线| 中文无码精品A∨在线观看不卡 | 亚洲国产午夜精华无码福利| 亚洲精品自产拍在线观看APP| 午夜精品久久久久久久99热下载 | 精品福利一区二区免费视频| 亚洲va欧美va国产综合下载| 亚洲无码日韩一区| 91在线国内在线播放老师| 成人免费一级片| 国产精品视频猛进猛出|