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

    寒區(qū)水塘冰蓋生長和消融分析

    2021-05-08 01:32:14丁法龍茅澤育
    水利學(xué)報(bào) 2021年3期
    關(guān)鍵詞:冰溫冰蓋熱傳導(dǎo)

    丁法龍,茅澤育

    (清華大學(xué)水利水電工程系,北京 100084)

    1 研究背景

    靜冰生消過程是寒冷地區(qū)冬季普遍存在的一種自然現(xiàn)象,多發(fā)生于湖泊、水庫、水塘和流速較低的河流。冰蓋的形成,一方面削弱了水體與大氣之間的能量和物質(zhì)交換,從而影響水體生態(tài)環(huán)境[1];另一方面所產(chǎn)生靜冰壓力威脅水工建筑物的安全,如護(hù)坡工程的冰推破壞等。冰蓋厚度是冰工程研究中重要的特征參數(shù)之一[2-4],也是建立冰情災(zāi)害預(yù)報(bào)模式的關(guān)鍵物理參數(shù)之一[5-6]。另外,冰蓋強(qiáng)度與厚度有關(guān),因此冰蓋厚度也是估算冰對水工結(jié)構(gòu)物作用力、開展冰上交通以及冰期捕魚等活動(dòng)的重要參考指標(biāo)。

    冰蓋厚度計(jì)算的數(shù)學(xué)模型分析方法主要是依據(jù)冰層生消物理過程建立熱力及動(dòng)力學(xué)數(shù)學(xué)模型,通過公式簡化推導(dǎo)或者數(shù)值求解以確定冰厚。德國科學(xué)家Stefan基于冰層內(nèi)部一維穩(wěn)態(tài)熱傳導(dǎo)方程推導(dǎo)提出的冰凍度日法是最早的冰厚生長解析模型,此后又有學(xué)者根據(jù)冰凍度日法原理,通過實(shí)測數(shù)據(jù)擬合確定經(jīng)驗(yàn)參數(shù),研制出多種適用于當(dāng)?shù)貧庀蠛退臈l件的衍生模型[7-8]。冰凍度日法及其衍生形式,將冰厚的發(fā)展過程與累積負(fù)氣溫聯(lián)系到一起,具有應(yīng)用方便、計(jì)算簡單的優(yōu)點(diǎn),但同時(shí)也存在以下不足:(1)假定上表面冰溫等于氣溫,即一維單相Stefan 問題,這已被實(shí)測證明是不合理的;(2)僅考慮了冰內(nèi)熱傳導(dǎo)單一熱力過程對冰蓋厚度變化的影響;(3)當(dāng)用于冰蓋消融計(jì)算時(shí),會(huì)不可避免地引入更多的經(jīng)驗(yàn)參數(shù)或不確定項(xiàng)。

    根據(jù)冰凍度日法,冰厚只有在氣溫轉(zhuǎn)正后才開始減小,但實(shí)測資料顯示,當(dāng)氣溫尚低于0 ℃時(shí)冰蓋厚度就已經(jīng)開始減小,因此冰凍度日法是否適用于冰蓋消融計(jì)算尚值得商榷。對于這個(gè)問題,眾多學(xué)者進(jìn)行了研究,提出了不同的解決方法,其中比較有代表性的有:Biello[9]以累積正氣溫修正冰厚最大值的計(jì)算方法模擬冰蓋的消融過程,被稱為融冰度日法;文獻(xiàn)[10]通過修正Stefan公式中的經(jīng)驗(yàn)參數(shù),采用統(tǒng)一的公式對整個(gè)冰蓋的生長和消融過程進(jìn)行模擬;張學(xué)成等[11]在Stefan公式中增設(shè)了一個(gè)隨機(jī)項(xiàng)并利用現(xiàn)場實(shí)測資料進(jìn)行校核。這些方法或包含隨機(jī)項(xiàng),或包含經(jīng)驗(yàn)參數(shù),其預(yù)測精度顯然受制于實(shí)測資料,不利于廣泛應(yīng)用,它們沒有完善考慮影響冰蓋消融的各個(gè)熱力因素,僅側(cè)重于冰厚對氣溫的依賴性,實(shí)際上,氣溫是通過影響大氣與冰面之間的對流傳熱過程來影響冰蓋生消的。

    隨著計(jì)算機(jī)技術(shù)和數(shù)值方法的進(jìn)展,基于冰內(nèi)熱傳導(dǎo)方程建立冰厚熱力學(xué)數(shù)值模型,是一種便捷高效的研究手段,多位學(xué)者開展了這方面的研究工作[12-13],但這些數(shù)值模型在邊界條件設(shè)定時(shí),往往對大氣-冰面?zhèn)鳠?、輻射等熱力要素存在不同程度的近似或簡化,從而?dǎo)致其計(jì)算結(jié)果與實(shí)際情況存在一定的差異。

    靜冰的生長和消融是一個(gè)主要發(fā)生在垂直方向上多層熱交換、熱傳導(dǎo)的復(fù)雜過程,其熱力因素包括冰蓋上表面的太陽輻射、長波輻射、感熱和潛熱,冰蓋下表面的水體傳熱,以及冰蓋內(nèi)的熱傳導(dǎo),它們共同組成冰蓋生消的熱力條件。因此,冰蓋熱力學(xué)模型研究的是一個(gè)由大氣-冰-水組成的熱能耦連系統(tǒng),精確的冰蓋厚度計(jì)算應(yīng)該立足于分析影響冰蓋生消的各個(gè)熱力過程,忽略其中某些熱力要素甚至僅考慮冰厚對氣溫響應(yīng)關(guān)系的計(jì)算模型是不盡合理的。

    此外,截至目前,圍繞靜冰生消的研究多側(cè)重于湖泊和水庫,鮮有關(guān)于水塘的研究報(bào)道,相比于湖泊和水庫,水塘的面積和水深更小,岸邊和塘底對水體的熱交換作用更為顯著,在冬季冰蓋形成之后,水體與冰之間的傳熱作用也更為顯著,從而導(dǎo)致塘冰的生消表現(xiàn)出其特有的規(guī)律性。

    為明晰寒區(qū)靜冰生消過程,探究塘冰生消特性,本文通過分析影響靜冰生消的各熱力過程,建立生長期和消融期冰蓋厚度計(jì)算的熱力學(xué)模型,對黑龍江省大慶市青花湖8 號(hào)水塘進(jìn)行冰情原型觀測,采用塘冰觀測數(shù)據(jù)對冰厚計(jì)算模型進(jìn)行參數(shù)選定和模型驗(yàn)證,并基于塘冰溫度鏈觀測數(shù)據(jù),分析表面冰溫隨氣溫和風(fēng)速的響應(yīng)關(guān)系、冰溫垂向分布廓線、冰溫和水溫隨時(shí)間的變化等規(guī)律,以期為寒區(qū)靜冰生消分析及冰厚精確計(jì)算提供參考。

    2 冰蓋厚度計(jì)算模型

    2.1 冰蓋生長過程冰蓋的生長可簡化為一維熱傳導(dǎo)問題,描述為:

    式中:Ti為冰的溫度,℃;ρi為冰的密度,kg/m3;Ci為冰的比熱容,J/(kg·℃);t為時(shí)間;z為沿冰蓋厚度方向坐標(biāo);λi為冰的熱傳導(dǎo)系數(shù),W/(m·℃)。

    式中:Qi為冰蓋內(nèi)部熱傳導(dǎo)的熱通量,W/m2。

    邊界條件為:冰蓋下表面冰溫恒為冰的冰點(diǎn)Tm(即0℃);冰蓋上表面冰溫為Ts。因此,由冰內(nèi)熱傳導(dǎo)所決定的冰蓋生長過程可表示為:

    式中:hi為冰厚;L為冰的相變潛熱,kJ/kg。

    對式(3)進(jìn)行時(shí)間積分,可以得到:

    式(4)中假定冰蓋上表面冰溫Ts等于氣溫Ta,則得到著名的Stefan冰厚公式:

    式中:AFDD為冰凍度日,等于日均負(fù)氣溫乘以天數(shù)的累積總和,即1AFDD=-1 ℃·d。

    Stefan冰厚公式僅考慮了冰內(nèi)熱傳導(dǎo)過程的影響,將大氣-冰面對流傳熱和冰內(nèi)熱傳導(dǎo)過程統(tǒng)一考慮,根據(jù)傳熱學(xué)理論稱之為冰的傳熱過程,冰的傳熱過程兩個(gè)環(huán)節(jié)的熱通量表達(dá)式如下:

    式中,Q為冰的傳熱過程的熱通量,W/m2。

    將式(6)改寫成溫差的形式:

    將式(7)兩式相加,消去溫度Ts,整理后得:

    綜上,將冰的傳熱過程、水體-冰傳熱和輻射等熱力因素綜合考慮在內(nèi),則冰的生長過程可表達(dá)為:

    式中:Qwi為水體向冰蓋下表面?zhèn)鳠岬臒嵬?,W/m2;Qrad為冰面凈輻射通量密度,W/m2;Kia為大氣與冰面之間的對流傳熱系數(shù),W/(m2·℃)。

    將式(9)表達(dá)為時(shí)間增量的形式:

    式(10)即靜冰冰蓋厚度生長的熱力學(xué)模型。

    2.2 冰蓋消融過程研究表明[14],在冰蓋消融期間,氣溫在回升到0 ℃以前,輻射項(xiàng)占主導(dǎo)地位;當(dāng)氣溫上升到0 ℃以上后,冰面對流熱通量變得非常重要。基于該認(rèn)識(shí)和上文冰蓋生長期熱力過程的分析,建立冰蓋消融熱力學(xué)模型如下:

    將式(11)改寫為時(shí)間增量的形式:

    式中:Ta-mean為日平均氣溫,℃。

    式(12)即靜冰冰蓋消融的熱力學(xué)模型。

    2.3 模型參數(shù)取值上述分析提出了生長期和消融期靜冰冰蓋厚度變化的計(jì)算模型(式(10)和式(12)),包含了對各熱力過程的能量計(jì)算,模型中具體參數(shù)的取值如下。

    2.3.1 大氣-冰面對流傳熱系數(shù)Kia根據(jù)Ashton研究[15],大氣與冰面之間的對流傳熱系數(shù)Kia與風(fēng)力等級(jí)大小有關(guān),在無風(fēng)(0級(jí))與大風(fēng)(8級(jí))之間的常見風(fēng)力范圍內(nèi),Kia約為10~20 W/(m2·℃),其中無風(fēng)時(shí)為10 W/(m2·℃),風(fēng)級(jí)與Kia之間的關(guān)系列于表1,模型計(jì)算時(shí)的Kia取值結(jié)合現(xiàn)場觀測期間的風(fēng)力數(shù)據(jù),由表1查取。

    表1 Kia值與風(fēng)力級(jí)別的關(guān)系

    2.3.2 冰面凈輻射通量密度Qrad(1)短波輻射(太陽輻射)QS。晴空無云遮時(shí)到達(dá)地面的太陽輻射通量密度QS-clear與緯度和日期有關(guān),由氣象研究結(jié)果[16]獲得原型觀測所在地區(qū)的晴空太陽輻射通量密度的季節(jié)變化過程,通過擬合回歸得到下式:

    式中:QS-clear為晴空太陽輻射通量密度,W/m2;d為該日期在一年內(nèi)的日序數(shù),即1月1日的d為1,上一年的12月31的d為-1,本公式使用范圍為-91 ≤d≤152,即上年10月1日至本年6月1日。

    天空有云層時(shí),太陽輻射由于云層覆蓋而減少,云遮條件下的太陽輻射通量密度QS-cloud計(jì)算采用下式[16]:

    式中:QS-cloud為云遮條件下的太陽輻射通量密度,W/m2;C為云量。

    云遮條件下入射的太陽輻射又會(huì)在地表部分發(fā)生反射,因此,冰蓋實(shí)際接收到的輻射通量密度為:

    式中:QS為冰面實(shí)際接收到的凈短波輻射通量密度,W/m2;A為反照率,對于冰面取值0.5[17]。

    (2)長波輻射(地面輻射)QL。凈長波輻射通量密度QL等于來自天空的入射長波輻射通量密度QL-in與冰面發(fā)出的長波輻射通量密度QL-out之間的差值,計(jì)算公式如下:

    式中:QL為凈長波輻射通量密度,W/m2;QL-in為來自天空的入射長波輻射通量密度,W/m2;QL-out為冰面發(fā)出的長波輻射通量密度,W/m2;k為經(jīng)驗(yàn)系數(shù),取0.18[17];σ為Stefan-Boltzmann 常數(shù),5.67×10-8W/(m2·℃4);εi、εa分別為冰和大氣的黑度。

    2.3.3 水體-冰傳熱的熱通量Qwi水體向冰的傳熱發(fā)生在冰蓋下方極薄的水層內(nèi),其熱通量的計(jì)算表達(dá)式:

    式中:λw為水的導(dǎo)熱系數(shù),W/(m·℃);dT/dz為極薄水層內(nèi)的溫度梯度,該溫度梯度需要很高空間分辨率的精密測溫儀器才能測得,存在較大困難。

    針對湖水-湖冰傳熱過程的研究[18-19]表明,湖冰生長期及消融期前期的Qwi集中分布在5~7 W/m2值域內(nèi),考慮到水塘相較于湖泊,其水溫更容易受到塘底的傳熱影響,結(jié)合前人研究成果,本文分別選取5種不同的Qwi值(6、8、10、12和14 W/m2)進(jìn)行計(jì)算,結(jié)合冰厚觀測數(shù)據(jù),確定塘冰Qwi的合理取值。

    2.3.4 其他參數(shù) 冰蓋生消熱力學(xué)模型式(10)和式(12)中的其他參數(shù)取值為:冰密度ρi取916 kg/m3,冰相變潛熱L取335 kJ/kg,冰熱傳導(dǎo)系數(shù)λi取2.14 W/(m·℃)[18],冰厚hi采用人工鉆孔測量;冰蓋上表面冰溫Ts由溫度傳感器測得;氣溫Ta通過自行觀測獲得;日平均氣溫Ta-mean采用4 個(gè)定時(shí)平均法求得。

    3 塘冰原型觀測

    3.1 觀測場地概況靜冰原型觀測地點(diǎn)位于黑龍江省大慶市青花湖西側(cè),地處東經(jīng)124°45′,北緯45°56′,湖區(qū)總面積28 km2。青花湖每年10月下旬至11月初開始結(jié)冰,冰蓋到翌年4月中下旬融化,最大結(jié)冰厚度1.2~1.3 m。每年7—8月份給湖內(nèi)補(bǔ)水,冬季結(jié)冰期停止補(bǔ)水,基本沒有水位變動(dòng),因此青花湖的冬季環(huán)境接近靜水水域,提供了靜冰冰情觀測的理想條件。青花湖西端的部分水域被人為地以漿砌石擋墻劃分為若干水塘,本次觀測選取塘底高程起伏較小的青花湖8號(hào)水塘為研究目標(biāo)。青花湖8 號(hào)水塘呈東西長度120 m、南北寬度82 m的長方形布置,塘內(nèi)水深2.4~2.6 m。

    3.2 觀測裝置與方法根據(jù)對當(dāng)?shù)匾酝嗄隁鉁睾捅橘Y料的分析,選擇觀測時(shí)段為2017年10月1日—2018年5月1日,和2018年10月20日—2019年4月29日。觀測地點(diǎn)的氣象數(shù)據(jù)由A753 WS 無線自動(dòng)氣象站提供,包括氣溫、空氣濕度、風(fēng)速和風(fēng)向、氣壓等;冰溫及水溫測量采用羅漢姆公司生產(chǎn)的RH-8068熱電阻PT100鉑熱電阻溫度傳感器測定,測溫精度±0.05 ℃;溫度數(shù)據(jù)采集至CR1000X 數(shù)據(jù)采集儀,數(shù)據(jù)采集頻率為10 min/次。

    溫度測點(diǎn)位于水塘中心處,從冰面至冰下布設(shè)一條溫度觀測鏈,觀測鏈上安裝有18個(gè)PT100溫度傳感器探頭,以記錄冰溫和冰下水溫,冰下溫度探頭隨冰生長會(huì)凍結(jié)到冰內(nèi)以記錄冰溫。溫度探頭的垂向布設(shè)位置為:冰蓋上表面處(觀測冰蓋表面冰溫Ts);冰面以下10、20、30、40、50、60、70、80、90、100、120、140、160、180、200、220和240 cm。觀測鏈采用木質(zhì)框架固定,選擇木質(zhì)材料的原因是木頭的導(dǎo)熱系數(shù)低,對冰層的熱影響較小。觀測鏈導(dǎo)線連接在支架上,支架固定于預(yù)先埋設(shè)的混凝土管樁之間,溫度觀測鏈布置如圖1所示。冰厚采用冰鉆穿孔后用塔尺于每日12時(shí)進(jìn)行測量,在溫度測點(diǎn)附近選擇3個(gè)冰厚測點(diǎn),采用3個(gè)測點(diǎn)的平均值作為每個(gè)采集日期的冰蓋厚度。

    圖1 冰情原型觀測

    4 結(jié)果與分析

    4.1 氣象條件及冰厚歷程分析青花湖以往多年氣象數(shù)據(jù),青花湖每年10月份至次年4月份的日平均氣溫的多年平均值在-18.9~-12.3 ℃之間,冬季累積負(fù)氣溫在1458~2053 ℃之間,風(fēng)力級(jí)別多年平均集中在0~5級(jí)范圍內(nèi)。圖2展示了2017—2018年觀測期間的日平均氣溫、風(fēng)速和云量數(shù)據(jù)。觀測期間的日平均氣溫的平均值為-14.97 ℃,低于同期多年平均值-2.26%。平均風(fēng)速3.1 m/s,屬3級(jí)風(fēng)力(微風(fēng));云量主要集中在0(晴空)~3級(jí)之間。

    青花湖冰封期持續(xù)約5~6 個(gè)月,地面平均凍深1.63 m,最大凍深1.86 m,最大結(jié)冰厚度1.2~1.3 m之間。在觀測的2017—2018和2018—2019兩個(gè)年份:結(jié)冰起始日期分別為10月27日和10月22日;最大冰厚日期分別為3月16日(冰厚114 cm)和3月9日(冰厚106 cm);冰蓋消失日期分別為4月21日和4月15日。2017—2018年冰期時(shí)長為176 d,其中生長期139 d,平均生長速率0.82 cm/d;消融期37 d,平均消融速率3.08 cm/d。2018—2019年冰期時(shí)長為174 d,其中生長期137 d,平均生長速率0.77 cm/d;消融期37 d,平均消融速率2.86 cm/d。

    圖2 2017—2018年冬季氣象資料

    4.2 冰厚計(jì)算模型驗(yàn)證根據(jù)冰厚計(jì)算熱力學(xué)模型式(10)和式(12),代入日平均氣溫、風(fēng)速、云量等氣象數(shù)據(jù)和表面冰溫等觀測數(shù)據(jù),分別選取不同的水體-冰傳熱的熱通量,逐日計(jì)算2017—2018年冰蓋生長期和消融期的厚度變化,結(jié)果如圖3和圖4所示。需要說明的是:(1)水體-冰傳熱的熱通量Qwi是時(shí)間變量,本文不研究其具體的時(shí)間變化過程,本文及文獻(xiàn)[18-19]中Qwi的取值均為時(shí)間平均值,不同的是文獻(xiàn)[18-19]中將冰蓋生消全過程的Qwi進(jìn)行平均,本文則對生長期和消融期分別取平均值;(2)嚴(yán)格地說,冰的熱傳導(dǎo)系數(shù)λi也是時(shí)間變量,在封凍、穩(wěn)封及消融的不同階段,由于冰層組成類型、凍結(jié)密實(shí)度、冰層內(nèi)部溫度分布的變化,冰的熱傳導(dǎo)系數(shù)λi也隨之發(fā)生相應(yīng)變化。由于冰熱傳導(dǎo)系數(shù)的時(shí)空變化規(guī)律十分復(fù)雜,因此至今國內(nèi)外在進(jìn)行冰厚生消計(jì)算時(shí),多采用常數(shù)值代替,不同學(xué)者對這個(gè)取值亦有所差異,本文采用了Matti Lepparanta[17]在湖冰研究時(shí)的取值,即2.14 W/(m·℃)。

    圖3中,水體-冰傳熱熱通量Qwi分別取值6、8、10、12和14 W/m2時(shí),冰厚計(jì)算值與觀測值之間的均方根誤差RMSE分別為3.6、1.9、1.3、4.8 和5.9 cm,即當(dāng)Qwi取10 W/m2時(shí),生長期的冰厚計(jì)算值與觀測值最為接近,這比前人湖冰研究成果中的取值稍大,這是因?yàn)橄啾扔诤?,水塘的面積和水深更小,導(dǎo)致岸邊和塘底對水體的熱交換作用和水體向冰的傳熱作用更為顯著,從而使得其Qwi值更大。

    圖3 2017—2018年冰蓋生長期的冰厚變化過程

    圖4中Qwi分別取值6、8、10、12和14 W/m2時(shí),冰厚計(jì)算值與觀測值之間的均方根誤差RMSE分別為6.7、4.6、3.8、2.3 和5.2 cm,即冰蓋消融期Qwi的最佳取值是12 W/m2,比生長期的10 W/m2更大,這是因?yàn)樵诒w消融期間,氣溫回升及太陽輻射通量密度增強(qiáng),均使得塘內(nèi)水體溫度升高,導(dǎo)致水體與冰蓋交界面處的溫度梯度增大,因此水體向冰蓋的傳熱作用更強(qiáng)。上述結(jié)果表明,水體向冰的傳熱可顯著影響冰蓋的熱力生消過程,忽略該熱力因素或者選用不合適的熱通量Qwi值,都會(huì)造成冰厚計(jì)算的誤差。

    根據(jù)2017—2018年冰期的冰厚計(jì)算結(jié)果,分別選定10和12 W/m2作為生長期和消融期冰厚計(jì)算模型中Qwi的最佳值,并對2018—2019年的冰厚進(jìn)行計(jì)算和實(shí)測對比,結(jié)果如圖5所示。斯皮爾曼等級(jí)相關(guān)系數(shù)用于表征兩個(gè)定序變量之間的相關(guān)性,本文引入該相關(guān)系數(shù)(以P表示),來定量評估實(shí)測冰厚歷程曲線與計(jì)算冰厚歷程曲線之間的相關(guān)度。由圖5可見,冰厚計(jì)算值與觀測值吻合良好,冰蓋生長期和消融期的斯皮爾曼等級(jí)相關(guān)系數(shù)分別為0.96和0.92,表明本文提出的冰厚計(jì)算模型具有一定的合理性,結(jié)合場地氣象條件和水文條件選擇合理的模型參數(shù),可以較為精確地預(yù)測靜水冰冰蓋厚度的生消演變過程。

    圖4 2017—2018年冰蓋消融期的冰厚變化過程

    圖5 2018—2019年冰厚生消過程計(jì)算值與觀測值對比

    4.3 溫度鏈觀測數(shù)據(jù)分析溫度鏈對觀測期內(nèi)的表面冰溫、冰層內(nèi)冰溫和冰下水溫進(jìn)行了持續(xù)監(jiān)測。表面冰溫由大氣和冰面之間的熱交換所決定,并與冰內(nèi)的熱傳導(dǎo)過程相耦合,難以通過熱力過程分析確定,以往研究多是通過建立其與氣溫之間的經(jīng)驗(yàn)關(guān)系來估算表面冰溫的大小。由上文的分析可知,大氣與冰面之間的對流傳熱與風(fēng)力大小有關(guān),即相同負(fù)氣溫條件下,不同的風(fēng)力等級(jí)可能會(huì)形成不同的表面冰溫。根據(jù)塘冰觀測數(shù)據(jù),對表面冰溫Ts作兩因素方差分析,兩個(gè)預(yù)設(shè)影響因素分別為氣溫Ta和風(fēng)速Vw。方差分析結(jié)果表明,當(dāng)取顯著性水平為5%時(shí),表面冰溫Ts與氣溫Ta和風(fēng)速Vw均顯著相關(guān)。對表面冰溫Ts作多元回歸分析,得到Ts隨Ta和Vw的變化關(guān)系表達(dá)式:

    式中:Vw為風(fēng)速,m/s。

    式(20)回歸模型的統(tǒng)計(jì)量中,相關(guān)系數(shù)R2為0.937,殘差的方差為0.084,回歸結(jié)果的殘差離零點(diǎn)均較近,且殘差的置信區(qū)間均包含零點(diǎn),說明回歸模型能較好的符合原始數(shù)據(jù)。圖6為選取的典型氣象條件下的表面冰溫觀測值與式(20)計(jì)算值的對照,從圖6中可以看出,表面冰溫計(jì)算值與實(shí)測值集中均勻分布于1∶1正比例直線附近,計(jì)算值的平均相對誤差為4.72%,預(yù)測效果較為理想。

    冰層內(nèi)部的冰溫是冰的基本狀態(tài)參量,圖7為2018年3月1日至3月2日不同深度冰層的冰溫及氣溫變化。由圖7可見,冰溫變化與氣溫的波動(dòng)存在一定的響應(yīng)關(guān)系,但各層冰溫隨隨時(shí)間的波動(dòng)較氣溫更為平緩,表層(10 cm)冰溫變化與氣溫變化趨勢較為一致,隨著冰層深度的增加,冰溫波動(dòng)幅度逐漸減小,冰溫也逐步升高。氣溫對冰溫的影響深度至多可達(dá)30 cm,冰層深度30 cm以下的冰溫隨時(shí)間變化近似一組水平直線。冰溫與氣溫最小值出現(xiàn)于每日6時(shí)—7時(shí),冰溫比氣溫大約滯后1 h,最大值出現(xiàn)在下午12時(shí)—16時(shí),冰溫滯后約2 h。

    圖6 表面冰溫計(jì)算值與觀測值對比

    圖8為2018年12月15日不同時(shí)刻的冰溫垂向分布,由于冰層上部受氣溫變化影響顯著,因此冰溫垂向分布廓線的上半部分(30 cm以上)冰溫隨時(shí)間變動(dòng)較大;下半部分(30 cm以下)冰溫基本呈線性分布,這與圖7反映的規(guī)律基本相同;同一冰深位置的冰溫日平均值沿冰深呈良好的線性分布。

    圖7 不同深度的冰層溫度變化(2018.03.01—2018.03.02)

    圖8 2018年12月15日不同時(shí)刻的冰溫垂向分布

    冰下水溫是控制冰下熱通量的關(guān)鍵,它直接同冰底面融化速率有關(guān)[20],這一點(diǎn)也可由本文關(guān)于水體-冰傳熱的熱通量對冰厚計(jì)算影響的分析可見。圖9為2017—2018年冰期內(nèi)冰蓋下不同水深處水溫的時(shí)間變化過程。由于冰蓋厚度處于不斷的生消之中,這里的水深指當(dāng)日自冰蓋底部以下的垂向位置;水溫取當(dāng)日2、8、14和20時(shí)這4個(gè)不同時(shí)刻的平均值。

    由圖9可見,在冰蓋生長期,受負(fù)氣溫的持續(xù)影響,冰下水體持續(xù)降溫,越靠近冰蓋水溫越低,整個(gè)水體在垂向上呈逆溫分布,最高水溫出現(xiàn)在塘底,約為3.5~4.3℃。隨著水深的增加,水溫變化速率逐漸減小。在冰蓋消融期,3月19日之前水溫變化較為緩慢,3月20日之后,冰蓋逐層消融,水體溫度整體快速回升,最大升溫速率可達(dá)0.13 ℃/d。不同時(shí)期的水溫沿垂線方向均呈良好的線性分布,但每日不同時(shí)刻的垂向分布趨勢略有差異。

    圖9 2017—2018年冰期冰蓋下不同水深處的水溫變化

    對比冰蓋生長期和消融期的水溫變化:生長期水溫持續(xù)降低且變化速率較為緩慢,主要是受氣溫變化的影響;消融期水溫先是緩慢回升,隨后快速升高且變化速率高于同期的氣溫回升速率,因此,消融期水溫的升高應(yīng)是受到氣溫轉(zhuǎn)正、地溫升高以及太陽輻射增強(qiáng)等多個(gè)作用的綜合影響。

    5 結(jié)論

    (1)通過分析靜冰生消的各熱力要素,考慮大氣-冰面對流傳熱、冰面輻射(短波輻射、長波輻射)、水體-冰傳熱過程,建立了生長期和消融期靜冰冰蓋厚度計(jì)算熱力學(xué)模型,并給出了模型中相關(guān)參數(shù)的取值。(2)針對黑龍江省青花湖8號(hào)水塘,分別選取不同的水體-冰傳熱熱通量Qwi進(jìn)行了冰厚計(jì)算,實(shí)測資料顯示生長期和消融期的最佳Qwi值分別為10和12 W/m2,這比湖冰研究中的取值更大,是因?yàn)樗恋拿娣e和水深更小,導(dǎo)致岸邊和塘底對水體的熱交換作用和水體向冰的傳熱作用更為顯著,從而使得其Qwi值更大。(3)采用參數(shù)優(yōu)化后的冰厚計(jì)算模型對2018—2019年的冰厚歷程進(jìn)行了計(jì)算,計(jì)算值與觀測值吻合良好,生長期和消融期的斯皮爾曼等級(jí)相關(guān)系數(shù)分別為0.96 和0.92,表明本文提出的冰厚計(jì)算模型可以較為精確地預(yù)測靜冰厚度的生消過程。(4)對觀測到的表面冰溫?cái)?shù)據(jù)作方差分析,結(jié)果表明,表面冰溫Ts與氣溫Ta和風(fēng)速Vw均顯著相關(guān);通過多元回歸分析得到了Ts隨Ta和Vw變化的關(guān)系表達(dá)式。(5)冰層內(nèi)的冰溫變化與氣溫的波動(dòng)存在一定的響應(yīng)關(guān)系,且上層冰溫變化與氣溫變化趨勢較為一致,隨著冰深增加,冰溫波動(dòng)幅度減??;黑龍江省青花湖8號(hào)水塘的觀測資料表明,氣溫對冰溫的影響深度至多可達(dá)30 cm,冰溫與氣溫最小值出現(xiàn)于每日6 時(shí)—7時(shí),冰溫比氣溫大約滯后1 h,最大值出現(xiàn)在下午12時(shí)—16時(shí),冰溫滯后約2 h。(6)冰下水體溫度在垂向上呈逆溫分布;在冰蓋生長期,冰下水體持續(xù)降溫,隨著水深的增加,水溫變化速率逐漸減小;在冰蓋消融期,水體溫度先緩慢隨后快速回升,最大升溫速率可達(dá)0.13 ℃/d。

    致謝:感謝通威股份有限公司王小平在冰型觀測中給予的幫助。

    猜你喜歡
    冰溫冰蓋熱傳導(dǎo)
    軍事文摘(2024年6期)2024-04-30 03:13:59
    一類三維逆時(shí)熱傳導(dǎo)問題的數(shù)值求解
    格陵蘭島的冰蓋悄悄融化
    參花(下)(2022年1期)2022-01-15 00:45:01
    加強(qiáng)研究 提高冰溫技術(shù)在食品保鮮中的應(yīng)用
    中國食品(2020年19期)2020-10-27 09:30:20
    長距離輸水工程的冰期冰蓋數(shù)值模擬研究
    熱傳導(dǎo)方程解的部分Schauder估計(jì)
    一類非線性反向熱傳導(dǎo)問題的Fourier正則化方法
    冰溫真空干燥過程中維持冰溫的方法初探
    一類熱傳導(dǎo)分布參數(shù)系統(tǒng)的邊界控制
    冰溫結(jié)合納他霉素對綠蘆筍采后生理品質(zhì)的影響
    美女免费视频网站| 久久久久久久午夜电影| 免费人成视频x8x8入口观看| 亚洲欧美日韩高清专用| 国产精品永久免费网站| 国内精品一区二区在线观看| 午夜福利欧美成人| av视频在线观看入口| 99精品久久久久人妻精品| 一本精品99久久精品77| 欧美激情久久久久久爽电影| 国产一区在线观看成人免费| 免费观看人在逋| 日本一区二区免费在线视频| 亚洲在线自拍视频| 精品少妇一区二区三区视频日本电影| 少妇被粗大的猛进出69影院| 欧美成人性av电影在线观看| 婷婷丁香在线五月| 亚洲国产精品合色在线| 97人妻精品一区二区三区麻豆| 国产精品日韩av在线免费观看| 国产亚洲av嫩草精品影院| 国产激情久久老熟女| 熟妇人妻久久中文字幕3abv| 97人妻精品一区二区三区麻豆| 听说在线观看完整版免费高清| 久久久国产成人免费| 一本精品99久久精品77| 黑人巨大精品欧美一区二区mp4| 亚洲精品在线观看二区| 午夜福利在线在线| 好看av亚洲va欧美ⅴa在| 老熟妇仑乱视频hdxx| 宅男免费午夜| 午夜免费激情av| 亚洲av电影不卡..在线观看| 亚洲一区二区三区不卡视频| 亚洲国产精品sss在线观看| 国产精品爽爽va在线观看网站| tocl精华| 窝窝影院91人妻| 国产99白浆流出| 1024香蕉在线观看| 久久性视频一级片| 亚洲 欧美一区二区三区| 亚洲中文字幕日韩| 长腿黑丝高跟| 免费在线观看成人毛片| 欧美乱色亚洲激情| 99re在线观看精品视频| 18美女黄网站色大片免费观看| 亚洲avbb在线观看| 精品不卡国产一区二区三区| 久久这里只有精品中国| 又大又爽又粗| 欧美大码av| 久热爱精品视频在线9| 国产不卡一卡二| 精品免费久久久久久久清纯| 亚洲国产中文字幕在线视频| 小说图片视频综合网站| 免费电影在线观看免费观看| 亚洲一卡2卡3卡4卡5卡精品中文| 一级毛片女人18水好多| 黄色 视频免费看| 精品乱码久久久久久99久播| 操出白浆在线播放| 日日摸夜夜添夜夜添小说| 身体一侧抽搐| 久久久精品欧美日韩精品| 一进一出抽搐动态| 国产伦在线观看视频一区| 午夜福利在线在线| 大型av网站在线播放| 老司机午夜福利在线观看视频| 欧美日本视频| 午夜免费成人在线视频| 久久精品国产亚洲av高清一级| 免费看日本二区| 又黄又粗又硬又大视频| 一区福利在线观看| 午夜精品在线福利| 淫妇啪啪啪对白视频| 日本黄色视频三级网站网址| 床上黄色一级片| 日韩有码中文字幕| 欧美乱码精品一区二区三区| 日韩精品免费视频一区二区三区| bbb黄色大片| 久久久久久久久中文| 亚洲 国产 在线| 日韩三级视频一区二区三区| 夜夜夜夜夜久久久久| netflix在线观看网站| 大型黄色视频在线免费观看| 一夜夜www| 特大巨黑吊av在线直播| 亚洲最大成人中文| 国产91精品成人一区二区三区| 国产成人影院久久av| 国产人伦9x9x在线观看| 国产精品久久久久久久电影 | 亚洲电影在线观看av| 美女黄网站色视频| 国产成人一区二区三区免费视频网站| 久久久精品大字幕| 欧美日本亚洲视频在线播放| 久久天躁狠狠躁夜夜2o2o| www.熟女人妻精品国产| 午夜福利高清视频| 精品免费久久久久久久清纯| 午夜成年电影在线免费观看| 两个人的视频大全免费| 久久中文字幕一级| 正在播放国产对白刺激| 国产爱豆传媒在线观看 | 日本熟妇午夜| 久久婷婷人人爽人人干人人爱| 18美女黄网站色大片免费观看| 91大片在线观看| 国产午夜福利久久久久久| 国产午夜精品久久久久久| 亚洲人与动物交配视频| 在线十欧美十亚洲十日本专区| av欧美777| 色精品久久人妻99蜜桃| 日日摸夜夜添夜夜添小说| 国产av在哪里看| 亚洲七黄色美女视频| 久久伊人香网站| 亚洲熟女毛片儿| 久久精品影院6| 可以在线观看的亚洲视频| 成人18禁在线播放| 欧美国产日韩亚洲一区| 91成年电影在线观看| 亚洲av成人av| a在线观看视频网站| 熟女少妇亚洲综合色aaa.| 一级毛片高清免费大全| 国产99白浆流出| 久久婷婷人人爽人人干人人爱| 国产精品自产拍在线观看55亚洲| 成人av在线播放网站| 欧美大码av| 老熟妇乱子伦视频在线观看| 成人一区二区视频在线观看| 欧美一级a爱片免费观看看 | 国产三级中文精品| 亚洲成人免费电影在线观看| 精品电影一区二区在线| 成人特级黄色片久久久久久久| 日本五十路高清| 美女黄网站色视频| 亚洲七黄色美女视频| 黄色片一级片一级黄色片| √禁漫天堂资源中文www| 久久久久久久久中文| 全区人妻精品视频| 亚洲美女视频黄频| 亚洲精华国产精华精| 夜夜躁狠狠躁天天躁| 18禁观看日本| 亚洲五月天丁香| 日韩大尺度精品在线看网址| 欧美在线黄色| 久久香蕉激情| 亚洲午夜精品一区,二区,三区| 久久中文看片网| 一级作爱视频免费观看| 欧美性长视频在线观看| 99精品欧美一区二区三区四区| 日韩有码中文字幕| 亚洲精品一卡2卡三卡4卡5卡| 一二三四社区在线视频社区8| 国产麻豆成人av免费视频| 欧美乱色亚洲激情| 免费在线观看完整版高清| 国产麻豆成人av免费视频| 亚洲欧美日韩东京热| 欧美午夜高清在线| 成年免费大片在线观看| 熟女电影av网| 香蕉丝袜av| 天堂影院成人在线观看| 久久国产乱子伦精品免费另类| 一级片免费观看大全| 50天的宝宝边吃奶边哭怎么回事| 此物有八面人人有两片| 俺也久久电影网| 中文字幕av在线有码专区| 欧美在线黄色| 在线观看日韩欧美| 白带黄色成豆腐渣| 舔av片在线| 黄色丝袜av网址大全| 嫩草影视91久久| 亚洲av熟女| 九色国产91popny在线| 国产精品影院久久| 久久午夜亚洲精品久久| 在线观看舔阴道视频| 国产免费av片在线观看野外av| 久久久国产成人免费| 免费在线观看视频国产中文字幕亚洲| www.999成人在线观看| 狂野欧美白嫩少妇大欣赏| 亚洲最大成人中文| 午夜福利免费观看在线| 亚洲天堂国产精品一区在线| 久久国产乱子伦精品免费另类| 十八禁人妻一区二区| 色综合婷婷激情| 黄色女人牲交| 国产97色在线日韩免费| 日韩欧美 国产精品| 日本免费a在线| 久久热在线av| 99re在线观看精品视频| 99热6这里只有精品| 亚洲专区国产一区二区| 久久人妻福利社区极品人妻图片| 99精品在免费线老司机午夜| 国产av一区二区精品久久| 在线观看免费午夜福利视频| 国产真人三级小视频在线观看| 人成视频在线观看免费观看| 国产熟女xx| 一进一出好大好爽视频| 特大巨黑吊av在线直播| 99国产极品粉嫩在线观看| 久久精品91无色码中文字幕| 日本熟妇午夜| 一a级毛片在线观看| 欧美在线一区亚洲| 午夜日韩欧美国产| 国产私拍福利视频在线观看| 91麻豆av在线| 天天一区二区日本电影三级| 1024香蕉在线观看| 特大巨黑吊av在线直播| 变态另类成人亚洲欧美熟女| 老汉色av国产亚洲站长工具| 日韩欧美在线二视频| 精品电影一区二区在线| 毛片女人毛片| 成年版毛片免费区| 午夜亚洲福利在线播放| 免费在线观看日本一区| 久久精品国产清高在天天线| 老司机午夜福利在线观看视频| 三级国产精品欧美在线观看 | 99久久精品国产亚洲精品| 一边摸一边做爽爽视频免费| 69av精品久久久久久| 成人一区二区视频在线观看| 国产精品,欧美在线| 欧美3d第一页| 国产精品久久久久久亚洲av鲁大| 俺也久久电影网| 无限看片的www在线观看| 岛国视频午夜一区免费看| 成人国语在线视频| 久久欧美精品欧美久久欧美| 最近在线观看免费完整版| 亚洲av五月六月丁香网| 一个人免费在线观看电影 | 久久婷婷成人综合色麻豆| 大型av网站在线播放| 中文字幕高清在线视频| 他把我摸到了高潮在线观看| 精品久久久久久成人av| 又粗又爽又猛毛片免费看| 亚洲五月天丁香| 欧美乱码精品一区二区三区| videosex国产| 一级毛片女人18水好多| 在线免费观看的www视频| 欧美一级毛片孕妇| 久久久国产欧美日韩av| 国产麻豆成人av免费视频| 国产亚洲av高清不卡| 香蕉丝袜av| 成人欧美大片| 亚洲男人天堂网一区| 99热这里只有精品一区 | 国产又色又爽无遮挡免费看| 欧美3d第一页| 女人高潮潮喷娇喘18禁视频| 美女 人体艺术 gogo| 欧美成人午夜精品| 欧美黑人巨大hd| 国产亚洲精品av在线| 亚洲国产欧洲综合997久久,| 正在播放国产对白刺激| 久久这里只有精品中国| 一卡2卡三卡四卡精品乱码亚洲| 男人舔女人的私密视频| 午夜精品久久久久久毛片777| 一二三四在线观看免费中文在| 99国产综合亚洲精品| 夜夜躁狠狠躁天天躁| 久久精品国产亚洲av香蕉五月| 国产97色在线日韩免费| 日本五十路高清| 午夜精品久久久久久毛片777| 黑人巨大精品欧美一区二区mp4| 51午夜福利影视在线观看| 中文字幕高清在线视频| 琪琪午夜伦伦电影理论片6080| 精品第一国产精品| 宅男免费午夜| 成人18禁高潮啪啪吃奶动态图| 亚洲成人中文字幕在线播放| 成年免费大片在线观看| 男女午夜视频在线观看| 久久久国产成人免费| 国产亚洲精品久久久久久毛片| 岛国在线免费视频观看| 一本大道久久a久久精品| 国产一区二区在线av高清观看| 久久国产精品人妻蜜桃| 免费av毛片视频| 特大巨黑吊av在线直播| 无遮挡黄片免费观看| 午夜激情福利司机影院| 国产爱豆传媒在线观看 | 国产男靠女视频免费网站| 免费在线观看完整版高清| 亚洲男人天堂网一区| 日本一二三区视频观看| 中文字幕熟女人妻在线| e午夜精品久久久久久久| 欧美日韩瑟瑟在线播放| 91麻豆精品激情在线观看国产| 在线观看免费午夜福利视频| 欧美绝顶高潮抽搐喷水| 亚洲激情在线av| 午夜精品久久久久久毛片777| 成人永久免费在线观看视频| 婷婷六月久久综合丁香| 真人一进一出gif抽搐免费| 黑人操中国人逼视频| 99久久精品热视频| 99精品久久久久人妻精品| 成人永久免费在线观看视频| 1024香蕉在线观看| 99久久精品国产亚洲精品| 国产精品1区2区在线观看.| av在线播放免费不卡| 久久天躁狠狠躁夜夜2o2o| 曰老女人黄片| a级毛片在线看网站| 成年版毛片免费区| 国产精品久久久久久人妻精品电影| 欧美日本视频| 久久人妻福利社区极品人妻图片| 在线观看舔阴道视频| 国产一区二区三区在线臀色熟女| 色综合婷婷激情| 欧美性猛交黑人性爽| 精品久久久久久久人妻蜜臀av| 亚洲精品国产精品久久久不卡| 淫妇啪啪啪对白视频| 美女高潮喷水抽搐中文字幕| 免费观看精品视频网站| 国产一区二区在线av高清观看| 免费在线观看日本一区| av有码第一页| 成人精品一区二区免费| 国产精品精品国产色婷婷| 精品久久久久久久人妻蜜臀av| 国产一区二区激情短视频| 国产av麻豆久久久久久久| 国产高清视频在线观看网站| 男女那种视频在线观看| 十八禁网站免费在线| 最近最新中文字幕大全电影3| 91麻豆精品激情在线观看国产| 午夜免费成人在线视频| 两人在一起打扑克的视频| 男女午夜视频在线观看| 久久草成人影院| 巨乳人妻的诱惑在线观看| 老司机午夜福利在线观看视频| 国产在线观看jvid| 国产精品久久电影中文字幕| 亚洲激情在线av| 欧美大码av| 精品乱码久久久久久99久播| 欧美一级a爱片免费观看看 | 9191精品国产免费久久| 黄频高清免费视频| 激情在线观看视频在线高清| 最近最新中文字幕大全免费视频| 日韩国内少妇激情av| 日本免费a在线| 亚洲va日本ⅴa欧美va伊人久久| 欧美日韩亚洲综合一区二区三区_| 国产成+人综合+亚洲专区| 精品免费久久久久久久清纯| 日韩国内少妇激情av| 搡老妇女老女人老熟妇| 免费无遮挡裸体视频| 老鸭窝网址在线观看| 又紧又爽又黄一区二区| 精品欧美国产一区二区三| 成人av一区二区三区在线看| 很黄的视频免费| 村上凉子中文字幕在线| 国产精品久久电影中文字幕| 国产精品免费一区二区三区在线| 怎么达到女性高潮| 精品久久久久久久毛片微露脸| 日韩欧美在线乱码| 亚洲七黄色美女视频| 亚洲色图av天堂| 午夜福利在线观看吧| 麻豆成人av在线观看| 国产精品久久视频播放| 亚洲一区二区三区色噜噜| 欧美日韩福利视频一区二区| av中文乱码字幕在线| 夜夜夜夜夜久久久久| 一级a爱片免费观看的视频| 特大巨黑吊av在线直播| 香蕉av资源在线| 久久国产精品人妻蜜桃| 亚洲中文av在线| 亚洲午夜理论影院| 岛国在线免费视频观看| 欧美成人午夜精品| 久99久视频精品免费| 精品人妻1区二区| 亚洲第一欧美日韩一区二区三区| 日本 av在线| 久久人妻av系列| 久久久水蜜桃国产精品网| 国产成人精品久久二区二区免费| 日本a在线网址| 老司机在亚洲福利影院| 青草久久国产| 天天躁夜夜躁狠狠躁躁| 少妇的丰满在线观看| 一区福利在线观看| 夜夜躁狠狠躁天天躁| 午夜福利免费观看在线| 国内精品一区二区在线观看| 桃红色精品国产亚洲av| 97碰自拍视频| 成人亚洲精品av一区二区| 两个人免费观看高清视频| 亚洲天堂国产精品一区在线| 女生性感内裤真人,穿戴方法视频| 99国产精品99久久久久| 韩国av一区二区三区四区| 国产精品爽爽va在线观看网站| 欧美绝顶高潮抽搐喷水| 亚洲精品久久国产高清桃花| 久久久久久九九精品二区国产 | 日本 av在线| 精品第一国产精品| 国产精品爽爽va在线观看网站| 老司机深夜福利视频在线观看| av超薄肉色丝袜交足视频| 亚洲中文字幕一区二区三区有码在线看 | 久久精品aⅴ一区二区三区四区| 伦理电影免费视频| 久久久久久国产a免费观看| 国产精品亚洲美女久久久| 男女视频在线观看网站免费 | e午夜精品久久久久久久| 国产一区二区在线观看日韩 | 99riav亚洲国产免费| 12—13女人毛片做爰片一| 桃色一区二区三区在线观看| 久久婷婷人人爽人人干人人爱| 高潮久久久久久久久久久不卡| or卡值多少钱| 国语自产精品视频在线第100页| 午夜免费激情av| 又粗又爽又猛毛片免费看| 亚洲第一欧美日韩一区二区三区| 免费电影在线观看免费观看| 国产成人精品久久二区二区免费| 大型av网站在线播放| 亚洲精品在线美女| 97超级碰碰碰精品色视频在线观看| 男人舔女人下体高潮全视频| 精品国产乱子伦一区二区三区| 国产精品精品国产色婷婷| 777久久人妻少妇嫩草av网站| 午夜精品在线福利| 久久人人精品亚洲av| 成人18禁在线播放| 国产亚洲欧美在线一区二区| www日本在线高清视频| 日韩大码丰满熟妇| 日韩欧美国产一区二区入口| 日韩欧美在线二视频| 老汉色av国产亚洲站长工具| 99国产精品一区二区三区| 色精品久久人妻99蜜桃| 亚洲全国av大片| 成熟少妇高潮喷水视频| 久久精品影院6| 又黄又粗又硬又大视频| 久久久久国内视频| 午夜福利在线在线| www.熟女人妻精品国产| 欧美性猛交黑人性爽| 激情在线观看视频在线高清| 麻豆成人午夜福利视频| 国产一区在线观看成人免费| 啦啦啦观看免费观看视频高清| 香蕉丝袜av| 国产视频一区二区在线看| 热99re8久久精品国产| 美女大奶头视频| 国产精品一区二区免费欧美| 我要搜黄色片| 在线观看www视频免费| 国产亚洲欧美98| 99久久精品热视频| 国产成人精品久久二区二区免费| 久久天躁狠狠躁夜夜2o2o| 好男人电影高清在线观看| 蜜桃久久精品国产亚洲av| 黄色视频不卡| 色老头精品视频在线观看| 欧美人与性动交α欧美精品济南到| 亚洲国产精品sss在线观看| 久久人妻av系列| 露出奶头的视频| 变态另类成人亚洲欧美熟女| 18禁黄网站禁片午夜丰满| 婷婷精品国产亚洲av在线| 亚洲免费av在线视频| 两个人免费观看高清视频| 亚洲中文日韩欧美视频| 免费看十八禁软件| 免费看a级黄色片| 最近最新中文字幕大全电影3| 极品教师在线免费播放| 欧美三级亚洲精品| 中出人妻视频一区二区| 亚洲18禁久久av| 亚洲va日本ⅴa欧美va伊人久久| tocl精华| 99re在线观看精品视频| 日日摸夜夜添夜夜添小说| 夜夜爽天天搞| 久久欧美精品欧美久久欧美| 午夜久久久久精精品| 90打野战视频偷拍视频| 国产亚洲精品久久久久久毛片| 亚洲免费av在线视频| 亚洲精品色激情综合| 女人爽到高潮嗷嗷叫在线视频| 黄片大片在线免费观看| 97人妻精品一区二区三区麻豆| 啦啦啦观看免费观看视频高清| 亚洲aⅴ乱码一区二区在线播放 | 一卡2卡三卡四卡精品乱码亚洲| 19禁男女啪啪无遮挡网站| 免费无遮挡裸体视频| 午夜老司机福利片| 国产精品一及| 免费在线观看视频国产中文字幕亚洲| 99在线视频只有这里精品首页| 国产精品香港三级国产av潘金莲| 成人精品一区二区免费| a在线观看视频网站| 少妇人妻一区二区三区视频| 脱女人内裤的视频| 日本 欧美在线| 女人高潮潮喷娇喘18禁视频| 精品久久久久久久久久免费视频| 两个人视频免费观看高清| 亚洲专区中文字幕在线| 久久久久精品国产欧美久久久| 亚洲九九香蕉| 国产一级毛片七仙女欲春2| 老司机福利观看| 亚洲无线在线观看| 日韩中文字幕欧美一区二区| 91麻豆av在线| 欧美性猛交黑人性爽| 午夜两性在线视频| 亚洲欧美精品综合久久99| 国产欧美日韩精品亚洲av| 亚洲色图av天堂| 97超级碰碰碰精品色视频在线观看| 欧美黑人欧美精品刺激| 成年版毛片免费区| 国产视频内射| 国产欧美日韩精品亚洲av| 婷婷精品国产亚洲av| 国产97色在线日韩免费| 变态另类丝袜制服| 小说图片视频综合网站| 可以免费在线观看a视频的电影网站| 哪里可以看免费的av片| 变态另类成人亚洲欧美熟女| 欧美3d第一页| 亚洲中文av在线| 亚洲中文字幕日韩| 成熟少妇高潮喷水视频| 亚洲狠狠婷婷综合久久图片| 色哟哟哟哟哟哟| 国内久久婷婷六月综合欲色啪| 亚洲狠狠婷婷综合久久图片| 欧美成人午夜精品| 人妻久久中文字幕网|