方 晴，馬劍波，沈文龍，熊 絲，陸寶宏

(1.河海大學(xué)水文水資源學(xué)院，江蘇 南京 210098；2.江蘇省秦淮河水利工程管理處，江蘇 南京 210022；3.江西省水利規(guī)劃設(shè)計(jì)研究院，江西 南昌 330029)

流域內(nèi)水庫的建設(shè)會(huì)形成流動(dòng)相對(duì)緩慢的大體積水域，使天然河道的水溫分布規(guī)律發(fā)生變化，加上與太陽輻射、入流、出流、庫底河床等的熱量傳輸，形成了特有的水庫溫度場(chǎng)。調(diào)節(jié)性能強(qiáng)、庫水交換周期長的大型水庫，沿水深方向多出現(xiàn)分層現(xiàn)象。

國外對(duì)水庫水溫的研究工作開展較早，美國在20世紀(jì)30年代便開始了水溫監(jiān)測(cè)，并在60、70年代提出了數(shù)值模擬方法，使水庫水溫研究步入了新時(shí)代。我國在20世紀(jì)50年代起步，對(duì)部分水庫進(jìn)行了監(jiān)測(cè)和研究，70年代提出一些經(jīng)驗(yàn)公式，80年代后引入并改進(jìn)國外數(shù)學(xué)模型，開展了大量關(guān)于水庫水溫的模擬研究[1]。但研究區(qū)域多集中在亞熱帶和溫帶地區(qū)，海南等熱帶地區(qū)的水庫水溫分布規(guī)律研究甚少。海南島屬典型的海洋性熱帶季風(fēng)氣候，太陽輻射強(qiáng)、日照時(shí)間長，氣溫較高且年變幅小，異于我國大多數(shù)地區(qū)的氣候特征。因氣候條件是影響水溫分布的重要因素，故該地區(qū)的水溫變化規(guī)律與其他地區(qū)存在較大差異。

本文以海南省松濤水庫為例，采用垂向一維模型模擬了水庫的水溫分布，分析了該水庫的水溫結(jié)構(gòu)和壩前垂向水溫分布規(guī)律，對(duì)當(dāng)?shù)赜嘘P(guān)部門制定保護(hù)措施、進(jìn)行流域規(guī)劃的環(huán)境影響評(píng)價(jià)具有一定的科學(xué)指導(dǎo)意義。

1 研究區(qū)域概況

松濤水庫位于海南省南渡江上游段，儋州市和白沙縣境內(nèi)，屬南渡江流域(見圖1)。控制流域面積1 496 km2，屬多年調(diào)節(jié)大(1)型水庫，總庫容33.45億m3，多年平均入庫水量14.9億m3，年庫水替換次數(shù)為0.45，庫面寬闊，平均水面寬2 460 m，最寬處達(dá)7 500 m，縱向尺度相對(duì)較小，庫區(qū)流速緩慢，整體形成湖泊型水庫。

圖1 松濤水庫位置示意

2 研究方法

垂向一維模型忽略溫度、密度等特性在橫向和縱向上的變化，有輸入條件簡(jiǎn)單、計(jì)算穩(wěn)定的特點(diǎn)，對(duì)于調(diào)節(jié)性能較強(qiáng)、水溫在橫縱向上變化甚微的湖泊型水庫，模擬垂向水溫效果較好[2]。松濤水庫為湖泊型水庫，水面寬闊，垂向變化明顯，符合一維的假定，故采用垂向一維水溫模型——DYRESM模擬水溫分布，探究熱帶地區(qū)湖泊型水庫水溫演變規(guī)律。

2.1 模型原理

DYRESM運(yùn)用拉格朗日分層理論，將水庫簡(jiǎn)化為一個(gè)具有多水平層且各層屬性一致的系統(tǒng)。模型可模擬、預(yù)測(cè)水庫不同時(shí)間尺度的熱量特征值，對(duì)環(huán)境因子和水庫屬性特征的長期變化進(jìn)行敏感性檢測(cè)，已成功應(yīng)用于多個(gè)湖泊及水庫的溫度、鹽度、密度預(yù)測(cè)[3- 6]。

模型以一維假定為前提，引入湖泊指數(shù)LN[7]和地球自轉(zhuǎn)的比率參數(shù)R[8]來校核一維假定的有效性。模型運(yùn)行主要遵守水量平衡、熱量平衡兩大定律，考慮因水體垂向移流、水氣界面熱交換、分子和紊動(dòng)擴(kuò)散、入流、出流引起的熱傳遞及水體吸收的太陽輻射。

2.2 模型構(gòu)建

模型需輸入地形資料、動(dòng)力變量、經(jīng)驗(yàn)參數(shù)和狀態(tài)變量。其中，地形資料包括水庫水位-水面面積關(guān)系、入流河口處詳細(xì)情況、取水口、大壩及溢洪道高程；動(dòng)力變量包括氣象條件、入流流量及水溫、出流流量等；經(jīng)驗(yàn)參數(shù)多為學(xué)者野外觀測(cè)、實(shí)驗(yàn)室試驗(yàn)或率定后確定[9]；狀態(tài)變量如溫度、鹽度、水位，用于模型初始化和驗(yàn)證。

2.2.1氣象資料

本次模擬時(shí)間為2016年3月2日～2017年3月1日，時(shí)間步長取一日，故氣象數(shù)據(jù)取日均值或日累積量。氣象數(shù)據(jù)選自距水庫最近的儋州氣象臺(tái)站同期觀測(cè)資料，輻射數(shù)據(jù)選自同緯度的?？谂_(tái)站同期觀測(cè)資料。

2.2.2地形、初始剖面資料

地形文件主要輸入水庫的地理位置和形態(tài)特征。

初始剖面資料需輸入一組自庫底至水面的高度及相對(duì)應(yīng)的水溫關(guān)系數(shù)組。采用松濤水庫2016年3月2日實(shí)測(cè)水溫(見圖2)作為初始輸入資料，由圖2可知，水庫于3月初已出現(xiàn)水溫分層，溫躍層厚8.81 m，水溫遞減率為0.13 ℃/m。

圖2 模型輸入初始水溫分布

2.2.3入、出庫流量資料

入庫流量選取2016年3月2日至2017年3月1日的松濤水庫實(shí)測(cè)入庫資料(見圖3a)，其中8月17日～19日和10月13日～19日期間入庫流量突增。入流水溫為南溪河入庫斷面處福才水文站多年平均水溫。出流資料為松濤水庫實(shí)測(cè)出流過程，見圖3b。

2.2.4參數(shù)、配置文件資料

根據(jù)已有的研究成果及模型經(jīng)驗(yàn)設(shè)置各參數(shù)，見表1和表2。

圖3 模擬時(shí)段內(nèi)日平均入、出庫流量過程

表1 模型參數(shù)推薦值

表3 模型率定采用的參數(shù)

表4 水溫垂向分布特征值統(tǒng)計(jì)結(jié)果

表2 配置文件參數(shù)初始值

2.3 模型率定

DYRESM是一個(gè)基于過程的參數(shù)化模型，涉及參數(shù)較多，且參數(shù)具有一定通用性，采用模型默認(rèn)參數(shù)即可得到較滿意的模擬結(jié)果[9]。研究表明，該模型模擬中小型水庫水溫結(jié)構(gòu)時(shí)精度較高，對(duì)于較大、深水庫，一些參數(shù)需重新率定，如垂向混合系數(shù)C、允許最大層厚、最小層厚等[10-12]?？紤]到國內(nèi)應(yīng)用DYRESM模型模擬水溫結(jié)構(gòu)的研究較少，尤其是熱帶地區(qū)，故模型默認(rèn)參數(shù)不一定完全適用，本次基于實(shí)測(cè)水溫?cái)?shù)據(jù)對(duì)參數(shù)重新率定。

率定時(shí)以模擬值與實(shí)測(cè)值的絕對(duì)誤差、均方根誤差最小為原則，所有參數(shù)在其適用范圍內(nèi)反復(fù)調(diào)試，最后采用參數(shù)值見表3。

采用該組參數(shù)對(duì)2016年4月26日和6月13日水溫進(jìn)行模擬，結(jié)果見圖4和表4。

圖4 松濤水庫水溫結(jié)構(gòu)分布

由圖4可知，模型基本能反映松濤水庫的水溫變化趨勢(shì)和結(jié)構(gòu)特點(diǎn)，兩組模擬值與實(shí)測(cè)水溫?cái)?shù)據(jù)總體吻合，擬合度較好，確定性系數(shù)分別為1.00和0.97，均方根誤差分別為0.46 ℃和1.25 ℃，最大絕對(duì)誤差出現(xiàn)于溫躍層，分別為1.15 ℃和3.25 ℃。

由表4可知，實(shí)測(cè)水溫與模擬值相比，滯溫層的擬合度最高，溫躍層擬合度最低。滯溫層的溫差不超過0.05 ℃，且位置基本吻合；表溫層水溫模擬值偏高，厚度偏大，但根據(jù)已有研究成果[3]，該模擬差異在可接受范圍內(nèi)；模擬的溫躍層中心位置上移，溫度梯度偏大，4月26日溫躍層中心上移2.25 m，溫度梯度偏大0.33 ℃/m，6月13日溫躍層中心上移0.81 m，溫度梯度偏大0.36 ℃/m。

綜上，率定后的DYRESM模型模擬效果較好，總體能揭示水庫水溫分布特征及變化趨勢(shì)，但模擬結(jié)果與實(shí)際值仍存在一定差異，且主要存在于溫躍層。經(jīng)分析，產(chǎn)生差異的主要原因?yàn)棰倌M出的溫躍層偏陡是DYRESM模型研究水庫水溫普遍存在的一個(gè)現(xiàn)象[11-13]，有些學(xué)者認(rèn)為可能是模型結(jié)構(gòu)存在一定問題。②DYRESM模型假定消光系數(shù)為恒定值，但實(shí)際卻隨時(shí)間和水深變化，其取常數(shù)時(shí)對(duì)模擬水溫有一定影響[10]。③溫躍層深度受內(nèi)部波運(yùn)動(dòng)及晝夜溫差的影響，內(nèi)部波運(yùn)動(dòng)會(huì)造成溫躍層的暫時(shí)性偏移。④研究表明[11]，DYRESM模擬中小型水庫水溫結(jié)構(gòu)效果更佳，因?yàn)槠渌綔囟炔町愋孕?，易滿足一維假定，風(fēng)應(yīng)力、熱源輸入的空間變化對(duì)大型水庫均會(huì)造成水溫在水平面上的差異。⑤本文利用的氣象資料源于儋州氣象站，氣象監(jiān)測(cè)數(shù)據(jù)尤其是風(fēng)速能否代表湖面實(shí)際氣象條件存在一定爭(zhēng)議[12]。

2.4 模型驗(yàn)證

松濤水庫僅有三日實(shí)測(cè)水溫?cái)?shù)據(jù)，已用于模型輸入和率定，且該地區(qū)同類型水庫水溫資料稀缺，故無法移用其他水庫驗(yàn)證模型，本文利用松濤水庫模擬水位過程與實(shí)際水位過程的擬合程度驗(yàn)證模型，見圖5。由圖5可知，模擬結(jié)果的擬合度較高，均方根誤差僅為0.18 m，說明模型能較準(zhǔn)確地模擬松濤水庫水量、水位變化過程，用于模擬該庫水溫變化是合理的。

圖5 松濤水庫模擬期內(nèi)實(shí)測(cè)、模擬水位過程

3 水庫水溫結(jié)構(gòu)分析

3.1 經(jīng)驗(yàn)法判別水庫水溫結(jié)構(gòu)

采用以下3種經(jīng)驗(yàn)法初步判別松濤水庫的水溫結(jié)構(gòu)。

(2)密度佛汝德數(shù)法。密度佛汝德數(shù)定義為慣性力與浮力的比值，即

(1)

式中，L為水庫縱向長度；H為平均水深；V為總庫容。松濤水庫回水長度約53 km，總庫容為33.45億m3，平均水深51.3 m，多年平均入庫流量51.4 m3/s，計(jì)算得Fr=0.005 1，可知松濤水庫水溫結(jié)構(gòu)為穩(wěn)定分層型。

(3)寬深比法。該方法的表達(dá)式為

R=B/H

(2)

式中，B為水庫水面平均寬度；H為水庫平均水深。松濤水庫平均寬度2 460 m，平均深度51.3 m，計(jì)算得R=47.95，可知松濤水庫的水溫結(jié)構(gòu)為過渡型。

綜上，不同經(jīng)驗(yàn)方法得出的結(jié)論不同。相比較而言，α-β指數(shù)法和密度佛汝德數(shù)法綜合考慮了水文條件和水庫特征，而寬深比法僅考慮了水庫的水深和水面寬，故前者更優(yōu)，且SL278—2002《水利水電工程水文計(jì)算規(guī)范》[14]推薦使用α-β指數(shù)法。但經(jīng)驗(yàn)法判定水溫結(jié)構(gòu)時(shí)僅考慮了水庫特征、入庫徑流量，導(dǎo)致判定結(jié)果可能有所偏差，而數(shù)值模型法理論嚴(yán)密，能較準(zhǔn)確地模擬出水庫水溫的時(shí)空變化。

3.2 水溫分布年內(nèi)變化規(guī)律

模擬期內(nèi)松濤水庫的水溫垂向分布變化過程見圖6、7，雙溫躍層結(jié)構(gòu)如圖8所示。

圖6 松濤水庫模擬期內(nèi)水溫垂向變化過程

圖7 水庫各月水溫垂向分布

圖8 水庫雙溫躍層結(jié)構(gòu)分布

由圖6～8可知，松濤水庫水溫呈顯著的穩(wěn)定分層結(jié)構(gòu)，垂向水溫變化有明顯的季節(jié)性：6、7月平均水溫最高，水面最高溫度達(dá)33.83 ℃；2月最低，庫表與庫底水體摻混均勻，平均水溫19.20 ℃；3月～翌年1月為水庫的分層期，2月為混合期；全年內(nèi)40 m以下至庫底水體水溫不隨時(shí)間改變，常年保持在19 ℃左右。

水庫水溫在年內(nèi)的具體變化規(guī)律為：

(1)升溫期。2月份進(jìn)入升溫初期，水溫為全年最低，上下層水體由于摻混水溫呈均勻分布，垂向溫差僅0.40 ℃；3月份隨著氣溫與太陽輻射的逐漸上升，庫表開始從外界吸收熱量，表層水體水溫上升，出現(xiàn)厚度較小(約3 m)、深度較淺、溫度梯度僅為0.55 ℃/m的溫躍層；4、5月進(jìn)入升溫末期，水體攝入的熱能大幅上升，且水庫已吸收儲(chǔ)存了大量熱能，分層強(qiáng)度持續(xù)加劇，表層溫度升高，同溫層厚度增大，溫躍層下移，溫度梯度達(dá)1.7 ℃/m；6、7月氣溫與太陽輻射達(dá)到頂峰，此時(shí)溫度結(jié)構(gòu)為全年最穩(wěn)定狀態(tài)，表層平均水溫達(dá)33 ℃，溫度梯度為1.6～1.88 ℃/m。

(2)降溫期。8、9月進(jìn)入降溫初期，水庫從外界吸收的熱量開始減少，但氣溫、太陽輻射仍處于較高水平，水溫穩(wěn)定分層，僅表層水溫較6、7月略低；8月18日出現(xiàn)雙溫躍層，19日中部均溫層厚度變大，該雙溫躍層穩(wěn)定存在至9月20日，此后上層溫躍層逐漸減弱，直至10月初完全消失，水庫恢復(fù)單溫躍層結(jié)構(gòu)。10月14日水庫再次出現(xiàn)雙溫躍層，15日中層均溫層加厚，至19日達(dá)到最厚，此時(shí)水庫不斷向外界釋放熱量，表層水溫逐漸下降，加之中部同溫層因入流減少、湍流擴(kuò)散減弱而逐漸變薄，上部溫躍層結(jié)構(gòu)受到擾動(dòng)，11月5日上部溫躍層消失，水庫恢復(fù)單溫躍層結(jié)構(gòu)。隨后氣溫與太陽輻射大幅下降，表層水體冷卻加快、密度增加，上層冷水與下層暖水不斷摻混，表溫層厚度增大，逐漸破壞30 m深的溫躍層結(jié)構(gòu)，溫度梯度減小，溫躍層下移，至1月底該溫躍層完全消失，水庫恢復(fù)上下層水體同溫狀態(tài)。由于該溫躍層深度大(水面以下30 m)，而水庫水溫與外界氣溫溫差逐漸減小，水體向外界散熱速率減慢，故該溫躍層存在時(shí)間長。

圖9為模擬期內(nèi)水庫表層水溫、入流水溫和氣溫的變化過程，可知全年內(nèi)三者關(guān)系為氣溫<入流水溫<表層水溫，不符合水溫滯后于氣溫的規(guī)律?？梢姛釒У貐^(qū)的水庫水體長時(shí)間處于強(qiáng)太陽輻射和高氣溫的環(huán)境下，吸熱多，水溫常年維持在一個(gè)較高水平，同時(shí)也反映了除氣溫外，其他條件如太陽輻射、風(fēng)速、水汽壓、入流水溫等也可能會(huì)對(duì)水溫產(chǎn)生影響。

圖9 水庫表層水溫、入流水溫與氣溫變化過程

3.3 雙溫躍層結(jié)構(gòu)原因分析

松濤水庫分別在8月18日和10月14日出現(xiàn)雙溫躍層，原因分析如下：兩次雙溫躍層的出現(xiàn)均與前一天入庫流量(見圖3a)猛增有關(guān)。模擬期水庫的年平均入庫流量為68.62 m3/s ，8月17日入庫流量1 950 m3/s，入庫水溫小于庫表水溫，入流水體沿水庫向下潛入與其密度相同的水層，該部分水體紊動(dòng)擴(kuò)散能力大，摻混作用強(qiáng)，形成一定厚度的中層同溫層，與頂層同溫層、下層同溫層形成雙溫躍層結(jié)構(gòu)。8月18日入流流量再次加大，達(dá)2 860 m3/s，19日中層同溫層厚度進(jìn)一步增大；同樣10月13日、14日、18日入庫流量分別為1 280、1 000、4 390 m3/s，造成了10月～11月雙溫躍層的出現(xiàn)。

兩次雙溫躍層持續(xù)的時(shí)間不同，前一次持續(xù)時(shí)間長(約40 d)，后一次持續(xù)時(shí)間短(約20 d)，其原因可能是8月～10月氣溫、太陽輻射還處于下降初期，水庫散熱速率慢，且該時(shí)期屬流域的主汛期，入庫流量較大，中層同溫層穩(wěn)定存在，只能靠表層水體緩慢冷卻來破壞雙溫躍層結(jié)構(gòu)；而10月中下旬后入庫流量大幅度減少，中部同溫層無法維持強(qiáng)烈的湍流作用，加上水庫散熱速率增快，故水庫雙溫躍層結(jié)構(gòu)持續(xù)時(shí)間縮短。

4 結(jié) 論

本文以海南省松濤水庫為例，采用DYRESM構(gòu)建垂向一維水溫模型，得出以下結(jié)論：

(1)松濤水庫為典型的分層型，3月～翌年1月為水庫分層期，6、7月分層最明顯，垂向最大溫差達(dá)13.90 ℃，8月份進(jìn)入降溫期，水庫水體出現(xiàn)翻滾，分層強(qiáng)度開始減弱，1月底分層現(xiàn)象完全消失，2月為混合期。

(2)滯溫層水體受環(huán)境影響較小，全年保持19.0 ℃左右；庫表水溫與入流水溫、氣溫的關(guān)系為庫表水溫>入流水溫>氣溫，與其他地區(qū)水溫滯后于氣溫的規(guī)律不同；模擬期內(nèi)由于入庫流量突然增大，導(dǎo)致水庫出現(xiàn)兩次雙溫躍層結(jié)構(gòu)。