孫博聞，楊晰淯，暴 柱，劉曉波，劉 暢，高學(xué)平

(1.天津大學(xué) 水利工程仿真與安全國家重點(diǎn)實(shí)驗(yàn)室，天津 300350；2.水利部海河水利委員會(huì)引灤工程管理局 海河流域?yàn)春铀|(zhì)監(jiān)測中心，河北 唐山 064309；3.中國水利水電科學(xué)研究院 水生態(tài)環(huán)境研究所，北京 100038)

1 研究背景

深水湖庫通常存在季節(jié)性溫度分層現(xiàn)象，水溫不僅影響水體的運(yùn)動(dòng)規(guī)律，也影響著水體中生化反應(yīng)以及水生生物的新陳代謝過程，包括溶解氧飽和度[1]、微生物和藻類生長[2]、沉積物向水體釋放化學(xué)成分[3-5]以及魚類對(duì)棲息地適宜性[2，6]；由溫度分層引起溶解氧等環(huán)境指標(biāo)的分層還會(huì)誘發(fā)水環(huán)境水生態(tài)問題[7-8]，進(jìn)而威脅供水安全及下游生態(tài)系統(tǒng)健康。在氣候變化與水庫調(diào)控的共同影響下，上述變化更加趨于復(fù)雜[9-10]，因此研究水溫變化對(duì)湖庫水質(zhì)管理與生態(tài)安全十分必要。

當(dāng)前對(duì)水庫水溫的研究大多關(guān)注其長期演變規(guī)律[11-13]，中短期預(yù)報(bào)方面的研究相對(duì)較少。事實(shí)上，在未來10天這一時(shí)間尺度上預(yù)報(bào)湖庫水溫變化對(duì)管理者具有重要意義[14]，如Weber等[15]利用水溫-深度變化數(shù)據(jù)動(dòng)態(tài)確定Grosse Dhuenn水庫取水深度，在優(yōu)化下泄水溫同時(shí)還抑制了水庫缺氧的發(fā)展；Baracchini等[16]在日內(nèi)瓦湖的觀測發(fā)現(xiàn)，強(qiáng)風(fēng)驅(qū)動(dòng)產(chǎn)生的劇烈涌升流使部分湖區(qū)的上層水溫在2天內(nèi)驟降10 ℃，直接影響下游的供水與河流生態(tài)。此外水溫垂向分布還決定了湖庫熱分層強(qiáng)度，進(jìn)而決定湖庫的分層或混合狀態(tài)[17]。秋季湖庫由分層狀態(tài)過渡到混合狀態(tài)時(shí)，庫底積累的營養(yǎng)物質(zhì)和重金屬在整個(gè)水體混合，此時(shí)冷雨等氣象條件誘發(fā)的翻庫過程則會(huì)使垂向混合在短時(shí)間內(nèi)發(fā)生[18]，這類由剖面水溫驟變引發(fā)的水動(dòng)力過程也會(huì)對(duì)湖庫生態(tài)安全產(chǎn)生較大沖擊。通過上述分析可見，中短期湖庫水溫預(yù)報(bào)能夠在調(diào)控湖庫下泄水溫這類常規(guī)管理中發(fā)揮作用，還能在應(yīng)對(duì)極端天氣(如風(fēng)暴)等誘發(fā)的翻庫、涌升流這類威脅供水與生態(tài)安全的應(yīng)急管理工作中提供決策支撐。

數(shù)學(xué)模型是湖庫水環(huán)境管理與決策分析的重要工具，初始條件與模型參數(shù)的確定是提升模型預(yù)報(bào)準(zhǔn)確性面臨的挑戰(zhàn)[19-20]，隨著實(shí)時(shí)觀測數(shù)據(jù)的獲取和傳輸?shù)玫斤w速發(fā)展，使得將原型觀測數(shù)據(jù)引入水環(huán)境預(yù)報(bào)，以提高預(yù)報(bào)精度成為可能[21]。數(shù)據(jù)同化融合觀測數(shù)據(jù)與模型計(jì)算結(jié)果，綜合考慮模型結(jié)構(gòu)誤差、邊界條件誤差和觀測誤差，對(duì)模型狀態(tài)變量進(jìn)行最優(yōu)估計(jì)，實(shí)現(xiàn)了模型狀態(tài)和參數(shù)不斷更新[19，22]，目前已應(yīng)用于提高河流磷遷移估計(jì)和湖泊富營養(yǎng)化等方面的模擬精度[23-24]。集合卡爾曼濾波是應(yīng)用最為廣泛的數(shù)據(jù)同化算法，該方法能夠提供狀態(tài)量的均值及其相應(yīng)的誤差協(xié)方差，同時(shí)由于引入了集合的思想且不需要伴隨或線性算子，因此具有適用于非線性系統(tǒng)、程序設(shè)計(jì)相對(duì)簡單等特點(diǎn)[25]。

本文以引灤入津工程源頭水庫—大黑汀水庫為研究對(duì)象，建立大黑汀水庫的二維水動(dòng)力模型，基于集合卡爾曼濾波算法構(gòu)建可綜合考慮模型參數(shù)、邊界條件以及觀測數(shù)據(jù)不確定性的湖庫水溫?cái)?shù)據(jù)同化系統(tǒng)。系統(tǒng)利用未來水庫調(diào)度數(shù)據(jù)與氣象數(shù)據(jù)等作為預(yù)報(bào)條件，能夠進(jìn)行未來1～10 d內(nèi)的湖庫水溫預(yù)報(bào)，高精度的水溫中短期預(yù)報(bào)方法可以為深水湖庫供水與生態(tài)安全提供理論與技術(shù)支撐。

2 研究方法

2.1 CE-QUAL-W2模型CE-QUAL-W2是一個(gè)考慮縱向和垂向的二維水動(dòng)力水質(zhì)模型(簡稱W2模型)[26]，模型假設(shè)水體橫向平均，適宜模擬河道型水庫的水動(dòng)力與水質(zhì)變化，模型主要控制方程如下。

連續(xù)性方程：

(1)

X方向動(dòng)量方程：

(2)

Z方向動(dòng)量方程：

(3)

狀態(tài)方程：

ρ=f(Tw，ΦTDS，ΦSS)

(4)

水面線方程：

(5)

質(zhì)量/熱量守恒方程：

(6)

式中：U與W為河道縱向與垂向流速，m/s；q為單寬流量，m3/s；B為網(wǎng)格寬度，m；g為重力加速度，m2/s；α為河底線與水平面夾角， °；η為水面線高程，m；Bη為水面寬度，m；ρ為密度，kg/m3；τxx與τxz與分別縱向流速的垂直梯度在x方向和z方向所產(chǎn)生的單位面積的剪切應(yīng)力，Pa；P為大氣壓，kPa；f(Tw，ΦTDS，Φss)為考慮水溫、總?cè)芙夤腆w、無機(jī)懸浮物的密度函數(shù)；h為沿河道底坡法線方向的水深，m；Φ為各組分濃度，g/m3；Dx為縱向擴(kuò)散系數(shù)，m/s2；Dz為垂向擴(kuò)散系數(shù)，m/s2；qΦ為入流組分濃度，g/(m3·s)；SΦ為源匯項(xiàng)，g/(m3·s)。

圖1 基于數(shù)據(jù)同化的湖庫水溫中短期預(yù)報(bào)工作流程

(7)

(8)

(9)

(10)

(11)

(12)

(13)

2.3 湖庫水溫中短期預(yù)報(bào)流程基于數(shù)據(jù)同化的湖庫水溫中短期預(yù)報(bào)包括四個(gè)階段(圖1)：(1)建模階段，根據(jù)獲取的地形、氣象等數(shù)據(jù)建立W2數(shù)學(xué)模型；(2)參數(shù)優(yōu)化階段，對(duì)影響W2水溫模擬的敏感參數(shù)增加噪聲擾動(dòng)，運(yùn)行N組W2模型，將水溫模擬值與觀測值輸入EnKF算法，得到關(guān)鍵參數(shù)的優(yōu)化結(jié)果，此階段模擬時(shí)長不少于20 d；(3)變量同化階段，在得到的關(guān)鍵參數(shù)后驗(yàn)分布基礎(chǔ)上，為邊界條件增加噪聲擾動(dòng)以體現(xiàn)其不確定性，再次運(yùn)行N組W2模型，將水溫模擬值與觀測值二次輸入EnKF算法，得到水溫同化結(jié)果；(4)預(yù)報(bào)階段，將同化結(jié)果作為本階段的湖庫水溫初始條件與參數(shù)組合，在獲取未來氣象數(shù)據(jù)、水庫入流條件以及調(diào)度計(jì)劃作為預(yù)報(bào)邊界條件后，運(yùn)行CE-QUAL-W2模型實(shí)現(xiàn)對(duì)湖庫水溫在未來1～10 d的中短期預(yù)報(bào)。

3 研究實(shí)例

3.1 研究區(qū)域與數(shù)據(jù)來源大黑汀水庫位于河北省唐山市遷西縣灤河干流上，是引灤樞紐工程體系中的骨干水庫工程，通過承接其上游來水與潘家口水庫泄水為天津、唐山等地供水。水庫總庫容約3.37億m3，為大Ⅱ型水庫，庫區(qū)回水長度約23 km，具有河道型水庫狹長的地形特征(圖2)。壩前水深可達(dá)25 m，全庫表現(xiàn)出典型的季節(jié)性水溫分層特征(圖2(d))。水庫溫度分層直接影響著水體理化過程，是目前大黑汀水庫供水所面臨的難題之一[29]。對(duì)大壩下游而言，大黑汀水庫下泄水溫是控制河道水溫的主要因素，進(jìn)而影響著下游生態(tài)系統(tǒng)及耕種活動(dòng)；對(duì)庫區(qū)而言，水溫分層誘發(fā)的庫底缺氧環(huán)境使得沉積物中的污染物發(fā)生厭氧分解，加速水質(zhì)惡化并威脅供水安全[30]。本文水溫觀測數(shù)據(jù)包括長期定點(diǎn)觀測與定期走航觀測(圖2)，定點(diǎn)觀測位于大黑汀水庫壩前500 m，在該斷面1、7和13 m等深度布置自動(dòng)水溫記錄儀，觀測頻率30 min/組；氣象數(shù)據(jù)源自國家氣象中心以及布設(shè)于水庫的在線氣象站；入庫水溫、流量及水庫調(diào)度數(shù)據(jù)由水利部海河水利委員會(huì)引灤工程管理局提供。

圖2 大黑汀水庫概況與2019年水溫觀測數(shù)據(jù)

3.2 大黑汀水庫水動(dòng)力模型本文將大黑汀水庫縱向網(wǎng)格間距(DLX)設(shè)為300 m，垂向?qū)雍穸?H)設(shè)為1 m，模型底部高程設(shè)為103 m，頂部高程為133 m，包含虛擬網(wǎng)格共劃分69個(gè)單元(segment)、最多33層(layer)，考慮到大黑汀水庫水面波動(dòng)較小，因此設(shè)水面坡度(slope)為0。W2模型中影響水溫計(jì)算的關(guān)鍵參數(shù)及取值列于表1。

表1 模型關(guān)鍵參數(shù)取值

對(duì)于垂向擴(kuò)散系數(shù)W2模型要求用戶選定計(jì)算模式后，由模型內(nèi)部計(jì)算得出，本文選定的計(jì)算模式為W2N。率定過程中發(fā)現(xiàn)，AX和DX對(duì)水溫結(jié)構(gòu)影響較小，因此選取默認(rèn)值1，這與李娟等[31]的研究結(jié)果一致。采用W2模型中的TERM算法計(jì)算水面熱量交換，該算法考慮了日照短波與長波輻射的直接入射和反射、水面-大氣反向熱輻射、空氣-水面熱傳導(dǎo)以及蒸發(fā)潛熱幾方面的影響。率定結(jié)果表明模型模擬的水動(dòng)力過程與實(shí)測值符合較好，可用于模型計(jì)算。模型參數(shù)率定與驗(yàn)證過程參照姚嘉偉[32]研究成果，本文不再贅述。

3.3 數(shù)據(jù)同化方案設(shè)計(jì)合理的同化方案與參數(shù)設(shè)置是保證同化效果的前提，EnKF算法利用集合思想解決了實(shí)際應(yīng)用中背景協(xié)方差矩陣估計(jì)和預(yù)報(bào)困難的問題，集合數(shù)越多，集合均值協(xié)方差與真值協(xié)方差越接近，但計(jì)算負(fù)荷會(huì)顯著增加，集合數(shù)過少，則可能導(dǎo)致模型誤差協(xié)方差估計(jì)錯(cuò)誤[33]。為此本文設(shè)置7組集合數(shù)(3、10、20、50、100、200和400)進(jìn)行數(shù)值試驗(yàn)，尋找可兼顧模擬精度與效率的集合數(shù)。數(shù)據(jù)同化通常包括只同化狀態(tài)變量以及同時(shí)同化狀態(tài)變量與模型參數(shù)兩種模式，根據(jù)徐興亞等[24]等的研究，后者表現(xiàn)通常優(yōu)于前者，因此本文也采用該模式。參考相關(guān)研究，選擇底部熱交換系數(shù)(CBHE)、太陽輻射吸收系數(shù)(BETA)和風(fēng)遮蔽系數(shù)(WSC)這三個(gè)對(duì)W2模型水溫模擬最敏感參數(shù)[34-35]，作為參數(shù)優(yōu)化階段考慮的模型參數(shù)。雖然數(shù)據(jù)同化方法能夠在理論上實(shí)現(xiàn)模型模擬結(jié)果和觀測數(shù)據(jù)的最佳融合，但這種“最佳”與模型的模擬誤差協(xié)方差和觀測誤差協(xié)方差的精度密切相關(guān)，因此選擇合適的模擬誤差和觀測誤差對(duì)提高同化精度同樣至關(guān)重要。參考李港等[36]研究，將觀測誤差與模擬誤差分別設(shè)置為1%、10%、20%、30%，篩選觀測誤差和模擬誤差的最優(yōu)組合。

在此基礎(chǔ)上，選擇W2模型的入流流量以及氣象數(shù)據(jù)中的太陽輻射熱量和風(fēng)速作為增加噪聲擾動(dòng)的邊界條件，以壩前測點(diǎn)1、7和13 m水深處的水溫觀測為同化數(shù)據(jù)(觀測頻率12 h/次)，對(duì)W2模型的參數(shù)和水溫進(jìn)行數(shù)據(jù)同化。首先將4月28日至5月9日的數(shù)據(jù)進(jìn)行同化模擬，再自5月10日起開展基于數(shù)據(jù)同化結(jié)果和未來氣象數(shù)據(jù)、水庫入流條件以及調(diào)度計(jì)劃作為預(yù)報(bào)邊界條件的水溫預(yù)報(bào)，每次預(yù)報(bào)未來1～10 d的水溫，直至10月27日(預(yù)報(bào)流程詳見圖1)。采用均方根誤差(RMSE)、一致性系數(shù)(IOA)和百分比偏差系數(shù)(PBIAS)評(píng)價(jià)數(shù)據(jù)同化系統(tǒng)的表現(xiàn)。IOA用來表征模擬值和觀測值的一致性，IOA取值在0～1之間，越接近1表示二者之間的一致程度越高。PBIAS結(jié)果以百分比定量化給出模擬值比觀測值整體被低估或高估的平均趨勢(shì)[37]，PBIAS為0表示模擬效果最佳，PBIAS>0表示模擬值傾向于低估，反之PBIAS<0表示模型傾向于高估。IOA和PBIAS的計(jì)算公式分別為：

(14)

(15)

4 結(jié)果與分析

4.1 數(shù)據(jù)同化方案確定對(duì)壩前斷面處不同集合數(shù)下的數(shù)據(jù)同化水溫模擬結(jié)果與實(shí)測水溫進(jìn)行對(duì)比(圖3)，可以看出隨著集合數(shù)增大RMSE逐漸減小，當(dāng)集合數(shù)達(dá)到100后RMSE減小幅度大幅降低，同化結(jié)果趨于穩(wěn)定。將觀測誤差與模擬誤差分別設(shè)置為1%、10%、20%、30%共4種情況，圖3展示了不同的誤差組合下的同化結(jié)果?？梢钥闯?，當(dāng)觀測誤差取1%，模型誤差取10%、20%時(shí)，數(shù)據(jù)同化結(jié)果的RMSE均較小。綜上，本文的數(shù)據(jù)同化方案為：集合數(shù)取100，模擬誤差與觀測誤差分別取10%和1%。

圖3 不同同化方案模擬結(jié)果對(duì)比

4.2 數(shù)據(jù)同化效果圖4為有無數(shù)據(jù)同化情況下，壩前500 m斷面1、7和13 m水深處W2模型模擬結(jié)果與數(shù)據(jù)同化系統(tǒng)模擬結(jié)果相比實(shí)測水溫的偏差值?？梢钥闯?，在經(jīng)過EnKF算法對(duì)模型參數(shù)和狀態(tài)變量同時(shí)進(jìn)行更新后，數(shù)據(jù)同化系統(tǒng)在1、7和13 m水深處的平均模擬偏差分別為1.01、0.31和0.21 ℃，較無數(shù)據(jù)同化的模擬結(jié)果分別提升了68.8%、51.6%和41.2%，模擬精度顯著提升。本文參考Zhang等[38]提出的方法，以垂向溫度梯度1.5 ℃/m為閾值識(shí)別溫躍層，當(dāng)某日所有觀測水溫剖面全部出現(xiàn)溫躍層時(shí)，認(rèn)定進(jìn)入分層期，否則為混合期。針對(duì)分層期(6月5日至10月3日)、混合期和全時(shí)段分別評(píng)價(jià)W2模型與數(shù)據(jù)同化系統(tǒng)的模擬結(jié)果，結(jié)果列于表2?？梢钥闯?，無論是分層期還是混合期，W2模型的均方根誤差(RMSE)明顯大于同化系統(tǒng)，最大可達(dá)到1.37 ℃；W2模型的一致性系數(shù)(IOA)相比1的偏差值也均大于同化系統(tǒng)；W2模型的百分比偏差系數(shù)(PBIAS)同樣大于同化系統(tǒng)?？梢姅?shù)據(jù)同化方法可以有效的提升模型的模擬效果，為提高湖庫水溫的中短期預(yù)報(bào)精度奠定了基礎(chǔ)。

圖4 W2模型模擬結(jié)果與數(shù)據(jù)同化系統(tǒng)模擬結(jié)果相比實(shí)測水溫的偏差值(2019年)

表2 W2模型與數(shù)據(jù)同化系統(tǒng)的模擬結(jié)果評(píng)價(jià)

4.3 數(shù)據(jù)同化系統(tǒng)預(yù)報(bào)性能數(shù)據(jù)同化系統(tǒng)可為湖庫水溫預(yù)報(bào)提供可靠的模型參數(shù)與初始狀態(tài)，在此基礎(chǔ)上利用未來氣象數(shù)據(jù)、水庫入流條件以及調(diào)度計(jì)劃作為預(yù)報(bào)邊界條件，在考慮邊界條件的不確定性前提下，運(yùn)行W2模型進(jìn)行湖庫水溫預(yù)報(bào)。為分析數(shù)據(jù)同化系統(tǒng)的水溫預(yù)報(bào)性能，對(duì)分層期、混合期與全時(shí)段下各水深處未來1～10 d的水溫預(yù)報(bào)結(jié)果進(jìn)行評(píng)價(jià)，結(jié)果如圖5所示?？梢钥闯?，在1與13 m處分層期的預(yù)報(bào)性能優(yōu)于混合期，7 m處兩階段預(yù)報(bào)性能各有優(yōu)劣。從RMSE指標(biāo)來看，隨著預(yù)報(bào)時(shí)段由1 d延長至10 d，預(yù)報(bào)誤差由0.22～0.35 ℃增大至0.77～1.09 ℃，可見隨著預(yù)報(bào)時(shí)段延長，不確定因素的累計(jì)效應(yīng)對(duì)預(yù)報(bào)會(huì)產(chǎn)生一定影響。從IOA指標(biāo)來看，分層期7 d預(yù)報(bào)期內(nèi)各深度的IOA值均超過0.65，具有較好的一致性；預(yù)報(bào)期延長至10 d后，7 m處的IOA由0.78降至0.57，降幅遠(yuǎn)超1與13 m深度。近年來引灤工程常采用潘家口水庫大流量調(diào)度抑制大黑汀水庫缺氧發(fā)展[39]，該現(xiàn)象可能與此種調(diào)度有關(guān)。從PBIAS指標(biāo)來看，1和13 m處預(yù)報(bào)均傾向于高估，7 m處則傾向于低估，但高/低估均在±3%以內(nèi)。從各統(tǒng)計(jì)指標(biāo)結(jié)果可以看出，構(gòu)建的數(shù)據(jù)同化系統(tǒng)能夠較好的預(yù)報(bào)水庫1～10 d內(nèi)水溫變化。

圖5 不同預(yù)測時(shí)段水溫預(yù)報(bào)效果對(duì)比

以10 d預(yù)報(bào)期為例，在整個(gè)模擬周期內(nèi)選擇湖庫水溫分層發(fā)展、穩(wěn)定和減弱的三個(gè)典型時(shí)段進(jìn)一步分析，對(duì)應(yīng)時(shí)段分別5月10—19日、8月2—11日和10月18—27日，預(yù)報(bào)結(jié)果如圖6所示，圖6中虛線為95%預(yù)報(bào)置信區(qū)間。由圖6(a)可以看出，熱分層形成階段受氣象條件的影響，庫區(qū)水溫整體緩慢上升、表底溫差逐漸增大，1與13 m的溫差由4.3上升至8.6 ℃。沿水深方向來看，1 m處的預(yù)報(bào)誤差始終保持在1.3 ℃以內(nèi)，7和13 m處預(yù)報(bào)誤差也始終保持在0.8和0.3 ℃以內(nèi)。7和13 m處水溫受氣象條件影響較小，受入流條件和水庫取水等的影響較大，5月18日預(yù)報(bào)得到的水溫上升現(xiàn)象即由水庫取水所致，這也表明同化預(yù)報(bào)系統(tǒng)能夠較好的預(yù)報(bào)湖庫調(diào)度對(duì)水庫水溫的影響。由圖6(b)可以看出，此時(shí)水庫熱分層已完全形成，入庫水溫和氣溫均維持在全年的較高水平，整體而言7和13 m的預(yù)報(bào)精度依然高于1m。數(shù)據(jù)同化系統(tǒng)預(yù)報(bào)出了8月6日受短期強(qiáng)降雨和降溫引起的水庫在1和7 m處1.6～2.3 ℃的降溫過程，13 m處水溫波動(dòng)不大。由圖6(c)可以看出，水庫水溫整體逐漸降低且1和7 m處水溫已基本相等，此時(shí)較弱的熱分層使10月24日時(shí)在大風(fēng)降溫驅(qū)動(dòng)下，水庫可能產(chǎn)生較強(qiáng)的垂向摻混，發(fā)生“翻庫”等危害供水安全事件的風(fēng)險(xiǎn)較大，在運(yùn)行管理中值得關(guān)注。整體而言，數(shù)據(jù)同化系統(tǒng)在受氣象條件及水庫調(diào)度等內(nèi)外部因素影響下，能夠在10 d的預(yù)報(bào)期內(nèi)維持較高準(zhǔn)確性，為湖庫調(diào)度管理提供依據(jù)。

圖6 2019年3個(gè)典型時(shí)段未來10 d水庫水溫預(yù)報(bào)結(jié)果

5 討論

5.1 觀測數(shù)據(jù)密度的影響觀測數(shù)據(jù)在數(shù)據(jù)同化過程中扮演著關(guān)鍵作用，觀測數(shù)據(jù)的密度與同化模型的預(yù)報(bào)精度息息相關(guān)。但較高的觀測頻次往往需要投入更多的設(shè)備成本，此外由于設(shè)備、通訊等軟硬件限制，并非所有觀測點(diǎn)都能保證有時(shí)間序列上的完整數(shù)據(jù)，因此有必要分析觀測數(shù)據(jù)密度對(duì)數(shù)據(jù)同化效果的影響。在對(duì)本文預(yù)報(bào)期內(nèi)9個(gè)時(shí)段進(jìn)行10 d的同化預(yù)測中，分別設(shè)置3個(gè)深度的觀測數(shù)據(jù)密度為12、24和36 h/個(gè)，圖7為不同深度水溫預(yù)報(bào)結(jié)果(RMSE)。

圖7 不同觀測數(shù)據(jù)密度下各深度預(yù)測RMSE結(jié)果對(duì)比

從圖7可以看出，在本文的數(shù)據(jù)同化模型中，當(dāng)預(yù)報(bào)點(diǎn)位與W2模型模擬時(shí)間一致的條件下，觀測數(shù)據(jù)密度從12 h/個(gè)降低至36 h/個(gè)時(shí)，水溫預(yù)報(bào)的精度逐漸降低。但觀測數(shù)據(jù)密度并非越高越好，崔凱鵬等[40]的研究表明，過高的觀測時(shí)間密度可能會(huì)產(chǎn)生信息冗余，Thomas等[14]的研究也發(fā)現(xiàn)由于設(shè)備故障產(chǎn)生的數(shù)據(jù)損失并未對(duì)同化表現(xiàn)產(chǎn)生持續(xù)影響。因此在兼顧效率和模型預(yù)測準(zhǔn)確性的同時(shí)，確定合適的觀測數(shù)據(jù)密度，并在模擬期內(nèi)明確哪些時(shí)段的觀測數(shù)據(jù)能帶來最大價(jià)值，對(duì)選取觀測資料有著重要意義。

5.2 同化方案設(shè)定的影響基于EnKF算法的同化方案設(shè)定中，集合數(shù)(N)、觀測誤差、模擬誤差等的選取對(duì)模型精度有較大影響，對(duì)于上述指標(biāo)選取均需要進(jìn)行大量的模擬測試。集合數(shù)方面，Houtekamer等[41]的研究表明N在100這一量級(jí)對(duì)于數(shù)據(jù)同化研究是足夠的；Evensen[42]的研究表明模型誤差與N-1/2相關(guān)，當(dāng)該值更接近于0后會(huì)有利于系統(tǒng)的誤差控制。相關(guān)研究中N的取值通常在20～500之間[14，23，33]，本文根據(jù)調(diào)整集合數(shù)后模型的RMSE表現(xiàn)確定N取值為100，結(jié)果表明該值能夠較好地兼顧同化系統(tǒng)的模擬精度與計(jì)算效率。觀測誤差和模擬誤差也會(huì)對(duì)模型預(yù)報(bào)精度有較大的影響，通常認(rèn)為觀測數(shù)據(jù)精度較高且更接近真實(shí)值，因此觀測誤差取值往往小于模擬誤差[23，43]。從本文16組誤差組合下的模擬結(jié)果看，當(dāng)模擬誤差一定時(shí)，RMSE隨觀測誤差的增大而增大；當(dāng)觀測誤差一定時(shí)，RMSE隨模擬誤差的增大而增大的趨勢(shì)并不明顯，劉卓等[23]認(rèn)為觀測誤差對(duì)同化精度的影響較模擬誤差大是產(chǎn)生這一現(xiàn)象的主要原因，因此設(shè)置較小的觀測誤差是合理的。

6 結(jié)論與展望

本文基于集合卡爾曼濾波算法與CE-QUAL-W2模型，構(gòu)建可綜合考慮模型參數(shù)、邊界條件以及觀測數(shù)據(jù)不確定性的湖庫水溫?cái)?shù)據(jù)同化系統(tǒng)，將其應(yīng)用于大黑汀水庫進(jìn)行1～10 d的中短期水溫預(yù)報(bào)，主要結(jié)論包括：(1)同化方案設(shè)計(jì)對(duì)數(shù)據(jù)同化系統(tǒng)具有重要影響，當(dāng)集合數(shù)為100、模擬誤差和觀測誤差分別為10%和1%時(shí)，本文構(gòu)建的同化系統(tǒng)能夠兼顧較高的計(jì)算效率與模擬精度。(2)采用雙重假設(shè)方法，經(jīng)過兩次集合卡爾曼濾波算法對(duì)W2模型的水溫敏感參數(shù)和水溫狀態(tài)變量進(jìn)行同步更新后，數(shù)據(jù)同化系統(tǒng)在1、7和13 m水深處的水溫模擬準(zhǔn)確性較無數(shù)據(jù)同化的模擬結(jié)果分別提升了68.8%、51.6%和41.2%，模擬精度顯著提升。(3)利用數(shù)據(jù)同化系統(tǒng)進(jìn)行水溫預(yù)報(bào)時(shí)，隨著預(yù)報(bào)期由1 d延長至10 d，不同水深的預(yù)報(bào)誤差由0.22～0.35 ℃增大至0.77～1.09 ℃。無論水庫處于分層期或混合期，數(shù)據(jù)同化系統(tǒng)均能夠在預(yù)報(bào)期內(nèi)的氣象條件及水庫調(diào)度等內(nèi)外部因素驅(qū)動(dòng)下維持較高準(zhǔn)確性，為水庫中短期調(diào)度管理提供理論支撐。此外本文是針對(duì)W2模型開展數(shù)據(jù)同化研究的首個(gè)案例，所構(gòu)建的同化系統(tǒng)對(duì)于其他應(yīng)用W2模型開展水動(dòng)力水質(zhì)模擬與預(yù)報(bào)的研究具有參考價(jià)值。

本項(xiàng)研究還能在以下四方面進(jìn)行提升：(1)不確定性分析方面。本文的預(yù)報(bào)效果除了與模型參數(shù)和邊界條件有關(guān)，還會(huì)受觀測數(shù)據(jù)密度與同化方案的影響，因此系統(tǒng)分析上述要素對(duì)預(yù)報(bào)結(jié)果的不確定性很有必要。此外本文只選取了3個(gè)對(duì)W2模型水溫模擬最敏感的參數(shù)進(jìn)行數(shù)據(jù)同化，這也可能造成一定的不確定性低估。(2)觀測數(shù)據(jù)獲取方面。本文采用同一點(diǎn)位、不同深度的原位觀測數(shù)據(jù)進(jìn)行同化系統(tǒng)構(gòu)建，但在對(duì)空間異質(zhì)性較強(qiáng)的湖庫進(jìn)行預(yù)報(bào)時(shí)，增大觀測數(shù)據(jù)的空間密度十分必要。布設(shè)多個(gè)觀測站點(diǎn)，提升觀測數(shù)據(jù)的空間密度較容易實(shí)現(xiàn)，但也意味著需要投入成倍的觀測成本。遙感數(shù)據(jù)具有較強(qiáng)的空間代表性，利用該方法獲取觀測數(shù)據(jù)已有研究在開展[33，36]，多源觀測數(shù)據(jù)的融合將會(huì)進(jìn)一步提升數(shù)據(jù)同化系統(tǒng)的預(yù)報(bào)性能。(3)從信息融合角度而言，數(shù)據(jù)同化與機(jī)器學(xué)習(xí)在數(shù)學(xué)原理上采用的均為貝葉斯估計(jì)策略[44]，融合機(jī)器學(xué)習(xí)和與機(jī)理模型的水溫預(yù)測研究已經(jīng)逐漸開展[45]，數(shù)據(jù)同化方法可能再次扮演重要角色。(4)數(shù)字孿生是當(dāng)前水利工程領(lǐng)域的重要研究方向，建設(shè)數(shù)字孿生水利工程就是為了服務(wù)預(yù)報(bào)、預(yù)警、預(yù)演、預(yù)案四項(xiàng)基本功能，其中預(yù)報(bào)功能要實(shí)現(xiàn)對(duì)水安全要素的高精度短期預(yù)報(bào)和中期預(yù)測，本文也可為開展數(shù)字孿生水利工程的中短期預(yù)報(bào)研究提供思路和方法。