亚洲免费av电影一区二区三区,日韩爱爱视频,51精品视频一区二区三区,91视频爱爱,日韩欧美在线播放视频,中文字幕少妇AV,亚洲电影中文字幕,久久久久亚洲av成人网址,久久综合视频网站,国产在线不卡免费播放

        ?

        參數(shù)率定方案對水溫數(shù)值模擬結(jié)果的影響

        2014-06-09 14:24:23陳青生魏小旺
        水利水運(yùn)工程學(xué)報 2014年4期
        關(guān)鍵詞:實(shí)測值庫區(qū)水溫

        劉 肖,陳青生,魏小旺

        (1.河海大學(xué) 水利水電學(xué)院,江蘇 南京 210098;2.浙江省錢塘江管理局勘測設(shè)計(jì)院,浙江 杭州 310016)

        參數(shù)率定方案對水溫數(shù)值模擬結(jié)果的影響

        劉 肖1,陳青生1,魏小旺2

        (1.河海大學(xué) 水利水電學(xué)院,江蘇 南京 210098;2.浙江省錢塘江管理局勘測設(shè)計(jì)院,浙江 杭州 310016)

        建立了立面二維水動力-水溫模型,運(yùn)用建庫前天然河道及鄰近區(qū)域同類型水庫兩組不同實(shí)測資料對模型進(jìn)行參數(shù)率定。根據(jù)這兩組計(jì)算參數(shù)分別模擬了驗(yàn)證時段內(nèi)D水庫水溫分布,以及A水庫庫區(qū)水溫變化過程。數(shù)值模擬結(jié)果表明,在其他條件保持一致的情況下,天然河道及鄰近區(qū)域同類型水庫資料均可用于數(shù)值模擬參數(shù)的率定,將兩組驗(yàn)證參數(shù)用于實(shí)際計(jì)算時,所得結(jié)果表現(xiàn)出相同的規(guī)律,但具體數(shù)值存在一定差異。分析認(rèn)為,鄰近區(qū)域同類型水庫資料驗(yàn)證更符合實(shí)際情況,故在此類水庫水溫的數(shù)值模擬中,應(yīng)盡量選擇鄰近區(qū)域同類型水庫資料作為驗(yàn)證資料,若在條件不足的情況下,相關(guān)的天然河道資料也可以作為參考。

        參數(shù)率定;水庫水溫;數(shù)值模擬

        大型水庫建成后,庫區(qū)內(nèi)水深增加,水體流速減小,水體熱量輸運(yùn)過程隨之發(fā)生變化,使得水溫分布結(jié)構(gòu)發(fā)生變化,庫區(qū)水溫分布及下泄水體的水溫不同于天然河道,對下游的水環(huán)境會產(chǎn)生一定影響[1-2]。因此,國內(nèi)外學(xué)者對模擬及預(yù)測水庫水溫進(jìn)行了大量的研究。張士杰等[3]總結(jié)了我國水庫水溫研究現(xiàn)狀及存在問題,并提出了今后的發(fā)展趨勢;胡平等[4]介紹了目前工程中預(yù)測水庫水溫分布的主要方法,并采用水庫水溫數(shù)值分析軟件,對二灘水庫和錦屏一級水庫水溫進(jìn)行數(shù)值預(yù)測。數(shù)值模擬方法是預(yù)測水溫的一種常用方法,而對所建數(shù)學(xué)模型的參數(shù)率定是數(shù)值模擬過程中很重要的一個步驟,參數(shù)的選擇是否合適,直接關(guān)系到計(jì)算結(jié)果的準(zhǔn)確性和可靠性[5-9]。然而,水溫預(yù)測研究工作一般是在水庫尚未建成前進(jìn)行,往往難以獲得水庫實(shí)測資料,一般可以采用鄰近區(qū)域同類型水庫或者天然河道的相關(guān)實(shí)測資料對模型參數(shù)進(jìn)行率定。目前對兩種模型參數(shù)的率定方式的對比分析研究較少。本文以西北地區(qū)某A水庫為例,采用立面二維水動力-水溫數(shù)學(xué)模型[10-11],利用天然河道水溫(方案1)及鄰近區(qū)域D水庫水溫(方案2)兩種不同的驗(yàn)證方案對模型分別進(jìn)行了驗(yàn)證,經(jīng)過參數(shù)率定,得到兩組不同的模型計(jì)算參數(shù),在保證其他條件一致的情況下利用這兩組計(jì)算參數(shù)模擬計(jì)算A水庫水溫,并對數(shù)值計(jì)算結(jié)果進(jìn)行了對比分析。

        1 水庫水溫數(shù)學(xué)模型的建立

        水庫水溫分層主要發(fā)生在垂向,水面風(fēng)速、風(fēng)向作為影響水體熱交換重要因素,其分布在河道縱向上需要加以考慮,結(jié)合水庫規(guī)模(壩前水深、水面面積、水面寬度和回水長度等)及調(diào)節(jié)性能,對于本文所研究的狹長型水庫確定采用沿橫向積分的立面二維的水動力水溫數(shù)學(xué)模型。立面二維水動力-水溫數(shù)學(xué)模型主要控制方程如下,方程定義河道走向?yàn)閤軸,垂向垂直于河道為z軸,橫向垂直于河道為y軸。

        連續(xù)方程:

        動量方程:

        式中:Tw為水溫;U,W為流速在水平方向和垂直方向上的分量;B為水面在y方向上的寬度;P為壓力;g為重力加速度;子x,子z為x方向和z方向上的切應(yīng)力;ρ為密度;η為水位;α為底坡;Ux為支流流速的x分量; q為單位寬度上的流量;qΦ為每單位體積橫向入流/出流中的溫度或組分質(zhì)量流量;Φ為溫度或組分濃度, Dx為縱向溫度或組分分散系數(shù),Dz為垂向溫度或組分分散系數(shù);SΦ為橫向平均的源/匯項(xiàng)(本文計(jì)算中無橫向入流/出流及源/匯項(xiàng),故取值為0)。

        上游入流邊界給定逐月來流流量、水溫,出流邊界給定水庫逐月下泄流量;水氣交界面給定月平均氣溫、風(fēng)速風(fēng)向、云量等。

        初始條件分兩次設(shè)定,首次計(jì)算初始條件給定計(jì)算起始時刻水位,初始溫度分布為計(jì)算起始時刻平均來流水溫,經(jīng)一個完整周期的計(jì)算得到穩(wěn)定結(jié)果后,采用首次計(jì)算的結(jié)束時刻結(jié)果作為正式計(jì)算的初始條件。

        本文選用有限差分法對控制方程進(jìn)行離散求解,得出式(1)~(6)中的U,W,TW,P,ρ和η共6個變量。采用交錯網(wǎng)格系統(tǒng)解決壓力梯度和連續(xù)方程[6]的離散困難。

        2 模型驗(yàn)證

        模型驗(yàn)證的目的主要在于校核相關(guān)氣象條件與計(jì)算參數(shù)的選擇是否恰當(dāng),能否滿足計(jì)算需要。根據(jù)現(xiàn)有實(shí)測資料,本文模型驗(yàn)證分為方案1水溫驗(yàn)證和方案2水溫驗(yàn)證2個步驟進(jìn)行,其中方案1還包括水面線的驗(yàn)證。

        2.1 方案1水溫驗(yàn)證

        根據(jù)A水庫壩址、庫區(qū)河段大斷面及實(shí)測水面線資料進(jìn)行水面線驗(yàn)證。選擇A水庫的天然河道作為水面線驗(yàn)證的計(jì)算驗(yàn)證區(qū)域,該區(qū)域河道長11.9 km,落差123.8 m。

        圖1 A庫區(qū)水面線驗(yàn)證Fig.1 Verification of the flow profile in Reservoir A

        根據(jù)選定的計(jì)算區(qū)域?qū)W(wǎng)格進(jìn)行劃分,從30 m精度DEM地形圖上確定該河段的走勢,以便更好地模擬實(shí)際河道水流狀態(tài),主流方向的網(wǎng)格節(jié)點(diǎn)根據(jù)實(shí)測河道深泓點(diǎn)位置來布置,尺度為15~228 m,垂向網(wǎng)格尺度均為1 m。根據(jù)實(shí)測水位與河道斷面資料由謝才公式反推來流流量,求得Q=360 m3/s,以此作為入流流量,下游邊界設(shè)為水位邊界,通過調(diào)整河道糙率使得計(jì)算值與實(shí)測值基本吻合,驗(yàn)證結(jié)果如圖1所示。

        由圖1可知,計(jì)算值與實(shí)測值二者基本吻合,表明各項(xiàng)水動力參數(shù)及條件選擇恰當(dāng),所建立數(shù)學(xué)模型可以準(zhǔn)確模擬該河段的水動力計(jì)算。

        以水面線驗(yàn)證為基礎(chǔ)進(jìn)行水溫驗(yàn)證。通過該河道上游S1及下游S2兩個水文站某年同一日期所測水溫值的分析比較,可以得出逐月的單位千米溫升值。由于A水庫庫區(qū)縱向長度僅為12 km,河段溫升不顯著,故根據(jù)所給庫區(qū)中的30個實(shí)測地形斷面資料,概化出一個平均斷面,以此為基礎(chǔ),按照此段河道的河道底坡,向上下游分別對計(jì)算區(qū)域進(jìn)行了一定長度的延伸,得到25 km的概化河道,為排除上下游邊界對河道水溫的影響,取該段河道中間的20 km為有效計(jì)算區(qū)域進(jìn)行分析。計(jì)算區(qū)域沿主流方向網(wǎng)格尺度為200 m,垂向尺度為0.4 m。入流條件采用S2水文站該年逐日實(shí)測流量資料,上游邊界給定S2水文站實(shí)測天然水溫,并給定氣溫、風(fēng)速風(fēng)向、云量等氣象資料,計(jì)算結(jié)果如圖2所示??梢?計(jì)算最大誤差為4%左右,計(jì)算值與實(shí)測值基本吻合,表明各項(xiàng)參數(shù)及條件選擇恰當(dāng),所建立數(shù)學(xué)模型可以用于水溫預(yù)測分析。

        2.2 方案2水溫驗(yàn)證

        根據(jù)現(xiàn)有實(shí)測地形及水文資料,對下游D水庫壩前垂向水溫進(jìn)行模擬驗(yàn)證。

        D水庫位于A水庫下游16 km處,壩高110 m,壩頂高程1 654 m,正常蓄水位1 646 m,正常蓄水位下總庫容1.25億m3。根據(jù)國際科學(xué)數(shù)據(jù)服務(wù)平臺提供的30 m精度DEM地形圖,對D水庫建模,沿河道方向網(wǎng)格尺度為100 m,網(wǎng)格數(shù)為122,垂向網(wǎng)格尺度為1 m,網(wǎng)格數(shù)為103。計(jì)算時采用的邊界條件見圖3。

        通過模擬計(jì)算,對相關(guān)參數(shù)進(jìn)行率定,得到D水庫6月3日(典型年)壩前垂向水溫分布,計(jì)算水溫與實(shí)測水溫的比較見圖4,可見,計(jì)算結(jié)果與實(shí)測值吻合良好。

        圖2 水溫驗(yàn)證Fig.2 Verification of the water temperature

        圖3 D水庫入流出流過程及氣溫和水溫分布Fig.3 Inflow and outflow,air temperature and inflow temperature of Reservoir D

        圖4 D水庫6月3日壩前水溫驗(yàn)證Fig.4 Water temperature verification of Reservoir D in front of the dam

        3 數(shù)值模擬結(jié)果及分析

        圖5 兩種方案下D水庫6月3日壩前水溫對比Fig.5 Water temperature comparison in front of the dam under two schemes

        3.1 運(yùn)用兩組方案分別模擬D水庫水溫

        為比較兩組參數(shù)對水庫水溫數(shù)值模擬的影響,首先采用方案1(天然河道)率定好的參數(shù)對D水庫水溫進(jìn)行了數(shù)值模擬,計(jì)算中其他初始條件、邊界條件、計(jì)算時段與方案2保持一致。圖5給出了同一時間節(jié)點(diǎn)上的數(shù)值模擬結(jié)果。從圖5可以看出,由方案2驗(yàn)證資料得到的計(jì)算模型參數(shù)對D水庫的水溫模擬與實(shí)測值擬合程度較好,而由方案1得到的D水庫的水溫模擬與實(shí)測值有一定的偏差,其偏差在水庫上層尤為明顯,與實(shí)測資料相比,方案2表層水溫相差約0.5℃,而方案1表層水溫與實(shí)測值相差約6.5℃。另一方面,兩種方案計(jì)算結(jié)果在水庫中下層均與實(shí)測結(jié)果相近,且?guī)靺^(qū)水溫垂向分布呈現(xiàn)相同的規(guī)律,存在明顯的水庫分層結(jié)構(gòu),表溫層、溫躍層、低溫層高度基本一致。

        可見,由于水庫建成后,與天然河道相比,水動力條件明顯改變,水深增加亦使得底部邊界熱動力條件發(fā)生變化,故而采用方案1(天然河道)率定的參數(shù)用于水庫水溫計(jì)算時,應(yīng)用條件存在明顯差異,導(dǎo)致計(jì)算結(jié)果偏差較大,且偏差主要體現(xiàn)在水庫表層。故采用鄰近區(qū)域水庫實(shí)測資料作為模型驗(yàn)證資料更為合理,其結(jié)果準(zhǔn)確性與可靠性更高。

        3.2 模擬A水庫庫區(qū)及下泄水溫

        A水庫正常蓄水位為1 820 m,正常蓄水位下總庫容約2.19億m3。以A水庫實(shí)測地形資料及樞紐特性為基礎(chǔ),構(gòu)造計(jì)算區(qū)域的網(wǎng)格,主流方向的網(wǎng)格節(jié)點(diǎn)根據(jù)實(shí)測河道深泓點(diǎn)的位置來布置,尺度為15~228 m,垂向網(wǎng)格尺度均為1 m。另外從30 m精度DEM地形圖上確定該河段的走勢情況,以便更好地模擬水庫水溫分布情況。

        根據(jù)S2水文站實(shí)測水溫利用沿程增溫率法推算A水庫庫尾的逐月水溫值作為A水庫的入庫水溫,模擬計(jì)算A水庫平水年(P=50%)正常蓄水位1 820 m下的水庫水溫分布及下泄水溫過程,A水庫上層水溫及電站引水口溫度計(jì)算結(jié)果如圖6所示。由圖6可見,兩組不同參數(shù)下的上、下層水體及電站引水口溫度變化規(guī)律一致,12月,1月,2月處于冬季,氣溫低于0℃,水面會出現(xiàn)結(jié)冰現(xiàn)象,河道及水庫上層溫度接近0℃。在其他月份,由于水庫內(nèi)水流流動較天然狀態(tài)下變緩,而相應(yīng)的氣溫值高于水溫值,故水庫上層水體的溫升比天然河道條件下大,水庫表層溫度高于天然河道的水溫。由此可知,兩組不同參數(shù)對A水庫庫區(qū)及下泄水溫的計(jì)算結(jié)果遵循同樣的變化規(guī)律,但在結(jié)果上存在一定差別。

        圖6 上層水體及電站引水口水溫對比Fig.6 Comparison between the water temperature of upper water body and that of intake of hydropower station

        4 結(jié) 語

        本文建立了立面二維水動力-水溫模型,運(yùn)用建庫前天然河道及鄰近區(qū)域同類型水庫兩組不同實(shí)測資料對模型進(jìn)行參數(shù)率定。根據(jù)這兩組計(jì)算參數(shù)分別模擬了驗(yàn)證時段內(nèi)D水庫庫區(qū)水溫分布,以及A水庫庫區(qū)水溫變化完整過程,得出以下結(jié)論:

        (1)由方案2(鄰近水庫)驗(yàn)證資料得到的計(jì)算模型參數(shù)對D水庫的水溫模擬與實(shí)測值擬合程度較好;而由方案1(天然河道)得到的D水庫的水溫模擬與實(shí)測值有一定的偏差,其偏差在水庫上層尤為明顯。

        (2)兩種方案計(jì)算結(jié)果在水庫中下層均與實(shí)測結(jié)果相近,且?guī)靺^(qū)水溫垂向分布呈現(xiàn)相同的規(guī)律,存在明顯的水庫分層結(jié)構(gòu),表溫層、溫躍層、低溫層高度基本一致。對A水庫的完整水溫變化過程數(shù)值模擬結(jié)果表明,兩種方案對A水庫庫區(qū)及下泄水溫的計(jì)算結(jié)果遵循同樣的變化規(guī)律,但在數(shù)值上稍有不同。

        (3)在其他條件保持一致的情況下,天然河道及鄰近區(qū)域同類型水庫資料均可用于數(shù)值模擬參數(shù)的率定,鄰近區(qū)域同類型水庫資料驗(yàn)證更符合實(shí)際情況,故在此類水庫水溫的數(shù)值模擬中,應(yīng)盡量選擇鄰近區(qū)域同類型水庫資料作為驗(yàn)證資料,若在條件不足的情況下,相關(guān)的天然河道資料也可以作為參考,仍能反映出水溫變化的大致規(guī)律。

        [1]蔡為武.水庫及下游河道的水溫分析[J].水利水電科技進(jìn)展,2001,21(5):20-23.(CAI Wei-wu.Analysis of water temperature of reservoir and downstream river[J].Advances in Science and Technology of Water Resources,2001,21(5):20-23.(in Chinese))

        [2]薛聯(lián)芳,顧洪賓,崔磊,等.紅水河干流梯級開發(fā)對水溫累積影響的調(diào)查研究[J].水力發(fā)電,2010,36(11):5-8. (XUE Lian-fang,GU Hong-bin,CUI Lei,et al.Study on cumulative effects to water temperature caused by the hydropower cascade development in Hongshui River[J].Water Power,2010,36(11):5-8.(in Chinese))

        [3]張士杰,劉昌明,王紅瑞,等.水庫水溫研究現(xiàn)狀及發(fā)展趨勢[J].北京師范大學(xué)學(xué)報:自然科學(xué)版,2011,47(3):316-320.(ZHANG Shi-jie,LIU Chang-ming,WANG Hong-rui,et al.A review on water temperature of reservoirs[J].Journal of Beijing Normal University(Natural Science Edition),2011,47(3):316-320.(in Chinese))

        [4]胡平,劉毅,唐忠敏,等.水庫水溫數(shù)值預(yù)測方法[J].水利學(xué)報,2010,41(9):1045-1053.(HU Ping,LIU Yi,TANG Zhong-min,et al.Numerical prediction of water temperature in reservoir[J].Journal of Hydraulic Engineering,2010,41(9): 1045-1053.(in Chinese))

        [5]董壯,羅龍洪,鄭福壽.岔管流動的數(shù)值模擬[J].河海大學(xué)學(xué)報:自然科學(xué)版,2007,35(1):15-17.(DONG Zhuang, LUO Long-hong,ZHENG Fu-shou.Numerical simulation of flow in bifurcated pipes[J].Journal of Hohai University(Natural Science),2007,35(1):15-17.(in Chinese))

        [6]潘俊,高維春,田海龍,等.潛流式人工濕地的數(shù)值模擬與參數(shù)率定[J].水利學(xué)報,2010,41(11):1280-1285.(PAN Jun,GAO Wei-chun,TIAN Hai-long,et al.Numerical simulation and parameter calibration for subsurface flow wetland[J]. Journal of Hydraulic Engineering,2010,41(11):1280-1285.(in Chinese))

        [7]董壯.三維水流數(shù)值模擬研究進(jìn)展[J].水利水運(yùn)工程學(xué)報,2002(3):66-73.(DONG Zhuang.Advances in 3D flow numerical simulation[J].Hydro-Science and Engineering,2002(3):66-73.(in Chinese))

        [8]SHANKAR N J,CHAN E S,ZHANG Q Y.Three dimensional simulation for an open channel flow with a constriction[J]. Journal of Hydraulic Reserach,2001,39(2):187-201.

        [9]FISCHER A T,STOESSER T,BATES P,et al.3D numerical modelling of open channel flow with submerged vegetation[J]. Journal of Hydraulic Reserach,2001,39(3):303-310.

        [10]郝紅升,李克鋒,梁瑞峰,等.支流影響下的水庫水溫預(yù)測模型[J].水利水電科技進(jìn)展,2006,26(5):7-9.(HAO Hong-sheng,LI Ke-feng,LIANG Rui-feng,et al.Water temperature prediction model for the reservoirs under the influence of branches[J].Advances in Science and Technology of Water Resources,2006,26(5):7-9.(in Chinese))

        [11]王冠,韓龍喜,常文婷.基于立面二維水動力-水溫耦合模型的水庫水溫分布[J].水資源保護(hù),2009,25(2):59-63. (WANG Guan,HAN Long-xi,CHANG Wen-ting.Modeling water temperature distribution in reservoirs with 2D laterally averaged flow-temperature coupled model[J].Water Resources Protection,2009,25(2):59-63.(in Chinese))

        Effects of parameter calibration schemes on numerical simulation results of water temperature

        LIU Xiao1,CHEN Qing-sheng1,WEI Xiao-wang2
        (1.College of Water Conservancy and Hydropower,Hohai University,Nanjing 210098,China;2.Surveying and Design Institute of Qiantang River Administration of Zhejiang Province,Hangzhou 310016,China)

        A vertical two dimensional hydrodynamic-temperature model is established for calculating variation in water temperature in reservoirs.The parameter calibration is carried out by comparing two sets of the observation data from the natural channel and the same types of reservoirs located nearby.According to the two sets of calibrated parameters,the water temperature in Reservoir D during the verification period and the complete process of the water temperature distribution in Reservoir A have been simulated respectively.The numerically simulated results indicate that the data from natural channels and the same types of the reservoirs located nearby can be used in the model parameter calibration when the other conditions remain the same.The numerically simulated results show the same tendency but different values applying those two sets of calibrated parameters.It is confirmed that the parameters obtained from the nearby reservoirs are more reasonable and close to realities.Therefore,it is recommended that in the reservoir water temperature simulation,the data from the same types of reservoirs located nearby are preferred,and the data from the natural channels can be used as a reference for analysis and studies of variation in the reservoir water temperature when the data are insufficient.

        parameter calibration;reservoir water temperature;numerical simulation

        TV697.2+1

        A

        1009-640X(2014)04-0082-05

        2013-12-09

        劉 肖(1988-),女,湖南岳陽人,碩士研究生,主要從事工程水力學(xué)研究。E-mail:liuxiao_1224@163.com

        猜你喜歡
        實(shí)測值庫區(qū)水溫
        丹江口庫區(qū)加強(qiáng)文物保護(hù)(之一)
        ±800kV直流輸電工程合成電場夏季實(shí)測值與預(yù)測值比對分析
        常用高溫軸承鋼的高溫硬度實(shí)測值與計(jì)算值的對比分析
        哈爾濱軸承(2020年1期)2020-11-03 09:16:22
        市售純牛奶和巴氏殺菌乳營養(yǎng)成分分析
        中國奶牛(2019年10期)2019-10-28 06:23:36
        基于PLC的水溫控制系統(tǒng)設(shè)計(jì)
        電子制作(2019年7期)2019-04-25 13:18:10
        一種基于實(shí)測值理論計(jì)算的導(dǎo)航臺電磁干擾分析方法
        電子制作(2018年23期)2018-12-26 01:01:22
        基于DS18B20水溫控制系統(tǒng)設(shè)計(jì)
        電子制作(2018年17期)2018-09-28 01:56:38
        突變理論在庫區(qū)移民后期扶持成效評估中的應(yīng)用
        庫區(qū)防護(hù)工程社會穩(wěn)定風(fēng)險識別
        廣東河源萬綠湖庫區(qū)的“雙音話”
        亚洲精品成人网久久久久久| 高清日韩av在线免费观看| 婷婷色综合视频在线观看| 美女又色又爽视频免费| 国产欧美一区二区成人影院| 99精品国产成人一区二区在线| 人妻少妇艳情视频中文字幕| 欧美牲交videossexeso欧美| 久久精品久久久久观看99水蜜桃| 欧美激情国产一区在线不卡| 亚洲精品成人一区二区三区 | 国产成人精品中文字幕| 大香焦av一区二区三区| 午夜福利试看120秒体验区| 九九在线视频| 国产免费99久久精品| 国产精品无码翘臀在线观看| 久久精品人人做人人综合| 丁香婷婷色| 国产91九色免费视频| 国产激情无码视频在线播放性色| 色噜噜狠狠色综合成人网| 手机AV片在线| 色播视频在线观看麻豆| 少妇厨房愉情理伦bd在线观看| 国产精品久久国产精品99gif| 精品日本一区二区视频| 亚洲精品国产av成人精品| 欧洲女人性开放免费网站| 四虎影视久久久免费| 国产久久久自拍视频在线观看| 欧美国产亚洲日韩在线二区| 国产mv在线天堂mv免费观看| 老肥熟女老女人野外免费区| 麻豆69视频在线观看| 国产免国产免费| 色噜噜狠狠色综合中文字幕| 国产亚洲精品一区二区在线观看 | 亚洲精品国产熟女久久久| 夜夜爽日日澡人人添| 女人与牲口性恔配视频免费|