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

丁法龍，茅澤育

（清華大學(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。

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