摘要：巴基斯坦卡洛特水電站采用靜庫容反推計算的入庫流量存在“鋸齒”問題。為此，量化分析了卡洛特水庫水位觀測誤差對入庫流量的影響，并利用3點線性平滑法、5點線性平滑法、7點線性平滑法、五點三次函數(shù)擬合平滑法，對卡洛特水電站2023年8月和9月入庫流量中的“鋸齒”波動進(jìn)行平滑。結(jié)果表明：卡洛特水電站正常蓄水位461 m時，0.01 m的水位誤差造成的流量誤差為15 m3/s，超過允許入庫流量預(yù)報誤差（13 m3/s），庫水位觀測誤差對流量預(yù)報工作存在不利影響。五點三次平滑算法可以有效消除反推入庫流量數(shù)據(jù)中的“鋸齒”干擾成分，更全面地保留原入庫流量的水文特征，平滑結(jié)果較理想。研究成果可為調(diào)節(jié)庫容小、入庫流量預(yù)報精度要求高的中小水庫反推入庫流量提供參考。

關(guān)鍵詞：入庫流量； 靜庫容反推計算； 流量平滑修正； 五點三次平滑法； 卡洛特水庫； 巴基斯坦

中圖法分類號：TV72

文獻(xiàn)標(biāo)志碼：A

DOI：10.15974/j.cnki.slsdkb.2024.11.003

文章編號：1006-0081（2024）11-0021-05

0 引 言

巴基斯坦卡洛特電站工程設(shè)計初期考慮到水庫庫尾與上一梯級規(guī)劃的阿扎德帕坦（Azad Pattan）水電站壩址相接，采用水量平衡法計算入庫洪水，發(fā)現(xiàn)流量存在嚴(yán)重的“鋸齒”現(xiàn)象，增加了用其校核預(yù)報成果的難度。

入庫流量計算的準(zhǔn)確性問題難以從測量手段解決，閔要武、張俊等［1-2］提出“動庫容”計算的方式，在三峽水庫進(jìn)行了模擬并取得了一定的效果。然而，多數(shù)水電站不具備設(shè)立多個水庫水位站的條件。目前采用平滑算法對入庫流量結(jié)果進(jìn)行平滑修正是處理此類問題的常用手段［3-4］，主要的平滑算法有 3 點、5 點、7 點線性平滑，五點三次函數(shù)擬合平滑，最小二乘曲線［5］，樣條插值平滑［6］，以及一些其他的濾波算法［7］ 。本文以卡洛特水電站實際計算入庫流量為例，對比分析了 3 點、5 點、7點線性平滑法，以及五點三次函數(shù)擬合平滑法對入庫流量的修正效果。

1 工程概況

卡洛特水庫總庫容1.88億m3，正常蓄水位以下庫容1.52億m3，死庫容1.03億m3，調(diào)節(jié)庫容0.49億m3，具有日調(diào)節(jié)性能?？逄厮娬镜纳虡I(yè)模式為容量電價機(jī)制，根據(jù)巴基斯坦中央電力采購局對卡洛特水電站的發(fā)電指標(biāo)及水資源利用考核要求，當(dāng)水電站入庫流量不大于水電站最大發(fā)電引用流量時，電站每小時均需要充分利用水資源。根據(jù)水庫運(yùn)行模式，當(dāng)卡洛特水庫入庫流量為624.2～1 248.4 m3/s時，水庫運(yùn)行處于“徑流模式”，購電協(xié)議要求水位需維持在正常蓄水位461 m運(yùn)行，即保持出入庫平衡，這對水庫入庫流量預(yù)報精度提出了極高的挑戰(zhàn)。

2 入庫流量誤差分析

基于水量平衡的反推入庫流量誤差來源主要有庫容曲線誤差、水位觀測誤差和閘門啟閉引起的泄流量變化［8］，以上誤差導(dǎo)致反推入庫流量偏離實際值，流量過程線呈現(xiàn)“鋸齒”狀?？逄厮娬緦儆诤拥佬退畮欤?021年11月下閘蓄水，蓄水時間短，水庫淤積對庫容曲線影響有限。本文定量分析水位觀測誤差對反推入庫流量的影響。

水庫反推入庫實測流量，都是根據(jù)庫水位和總出庫流量來計算。根據(jù)水庫水量平衡方程，計算公式為

Q1+Q22Δt-q1+q22Δt=V2-V1（1）

式中：Q1，q1為時段初入庫、出庫流量，m3/s；Q2，q2為時段末入庫、出庫流量，m3/s；V1，V2為時段初、末的水庫蓄水量，m3，由庫水位查水庫庫容曲線獲得；Δt為時間步長，s。

對式（1）取時段平均值，有：

QΔt-qΔt=V2-V1（2）

式中：Q為時段平均入庫流量，m3/s；q為時段平均總出庫流量，m3/s；V1=V（H1），V2=V（H2），H1和H2分別為時段初始和末了時刻的庫水位，m，V1和V2為相應(yīng)時刻的庫容，m3，由庫水位查水庫庫容曲線獲得。

根據(jù)庫水位查算庫容，再由式（2）推算時段平均入庫流量，誤差存在放大現(xiàn)象，其放大的倍數(shù)與水庫水面面積成正比，與時段長成反比。設(shè)時段初始和末了時刻的庫水位差為ΔH=H2-H1，H1和H2間水庫橫斷面的平均面積為A（單位為km2），則式（2）可表示為

Q=ΔH×A3600×Δt+q（3）

當(dāng)時段始末庫水位存在誤差（單位為m）時，就有相應(yīng)的流量誤差ΔQ：

ΔQ=1000××A3.6×Δt（4）

對于現(xiàn)行的大多數(shù)水庫實時洪水預(yù)報系統(tǒng)，時段間隔為1 h，則式（4）的誤差ΔQ與水面面積呈線性比例關(guān)系，卡洛特水電站流量誤差與水位誤差關(guān)系見表1。

卡洛特水電站正常蓄水位461 m，電站運(yùn)行水位長期控制在460～461 m。庫水位為461 m時，0.01 m的水位誤差，其流量差值為15 m3/s；0.02 m的水位誤差，其流量差值為31 m3/s；當(dāng)水位誤差為0.05 m時，其流量差值就達(dá)到了77 m3/s?？逄厮娬?0月至次年2月多年平均流量約270 m3/s，按照5%的預(yù)報誤差，允許入庫流量預(yù)報誤差往往只有13 m3/s左右。因此，庫水位誤差變化0.01 m帶來的流量差值就足以影響預(yù)報精度評定結(jié)果。

3 研究方法

3.1 滑動平均法

滑動平均法在數(shù)學(xué)和統(tǒng)計學(xué)中是一種常用的數(shù)據(jù)平滑方法。通過對原始數(shù)據(jù)進(jìn)行平滑處理，可以減少數(shù)據(jù)中的噪聲和波動，使數(shù)據(jù)更加平穩(wěn)和可靠?；驹硎菍r間序列資料進(jìn)行逐項推移計算，可按照3點、5點或7點等長度逐步遞進(jìn)，依次計算各長度項數(shù)的時序平均值，以弱化隨機(jī)波動的影響［9］。計算方法見式（5）。

Yk=∑pi=-qwiXk+ik=q+1，q+2，…，N-p（5）

式中：Y為平滑后的時間序列；X為原始時間序列；n為滑動窗口長度；N為時間序列的長度。wi為加權(quán)系數(shù)；p，q 為小于n 的任一正整數(shù)，且p+q+1=n。

3.2 五點三次函數(shù)擬合平滑法

滑動平均法［9］應(yīng)用領(lǐng)域比較廣泛，五點三次函數(shù)擬合平滑法在最小二乘法的基礎(chǔ)上加入了數(shù)據(jù)平滑技術(shù)，使得該算法更加注重數(shù)據(jù)的細(xì)節(jié)處理［10-11］，且算法簡單、運(yùn)算速度快。國內(nèi)水庫的相關(guān)研究［12-14］顯示，五點三次函數(shù)擬合平滑法在平滑入庫流量的問題中取得了理想的結(jié)果。

設(shè)在2n+1個等距為h的節(jié)點x-n，x-n+1，…，x0， …，xn-1，xn上，數(shù)據(jù)分別為Y-n，Y-n+1，…，Y0，…，Yn-1，Yn，作交換t=（x-x0）/h，則原節(jié)點變?yōu)閠-n=-n， …，t0=0，…，tn=n，建立m次多項式，見式（6），系數(shù)由最小乘法確定，方差和見式（7）。

Yt=a0+a1t+a2t2+…+amtm（6）

F（a0，a1，…，am）=∑ni=-n∑mj=0ajtji-Yi2（7）

為使方差和最小，令：

Fak=0（k=0，1，…，m） （8）

可得方程組：

∑ni=-nYjtki=∑mj=0aj∑ni=-ntk+1i（9）

當(dāng)n取2（即5個節(jié)點）、m取3時，將計算所得的a0，a1，a2代入式（7），同時令t=-2，-1，0，1，2，得到五點三次平滑公式。設(shè)卡洛特水電站的入庫流量序列為Q1，Q2，…，Qt，…，Qn-1，Qn（ngt;5），將Qn代入五點三次平滑公式中的Yn，則得到五點三次平滑計算公式如下。當(dāng)節(jié)點個數(shù)超過5時，除在序列兩端4個節(jié)點分別使用式（10）、式（11）、式（13）及式（14）進(jìn)行平滑，其余節(jié)點均用式（12）進(jìn)行平滑。

Q1=170（69Q1+4Q2+6Q3+4Q4-Q5）（10）

Q2=130（2Q1+27Q2+12Q3-8Q4+2Q5）（11）

Qt=135（-3Qt-2+12Qt-1+17Qt+12Qt+1-3Qt+2）（2lt;tlt;n-1）（12）

Qn-1=1352Qn-4-8Qn-3+12Qn-2+27Qn-1+2Qn（13）

Qn=170-Qn-4+4Qn-3-6Qn-2+4Qn-1+69Qn（14）

式中：Q1為t=1時經(jīng)五點三次平滑后的改進(jìn)流量，m3/s；Qt為當(dāng)2lt;tlt;n-1時，平滑后的改進(jìn)流量，m3/s；Qn為t=n時平滑后的改進(jìn)流量，m3/s。

4 平滑效果評價指標(biāo)

本文引入平滑度和相關(guān)系數(shù)指標(biāo)評價平滑后的入庫流量。

（1） 平滑度由平滑后數(shù)據(jù)和原始數(shù)據(jù)差分平方和構(gòu)成。平滑度范圍為0～1，且越接近1曲線越光滑。

Esm=1-∑n-1i=1（Qde，i+1-Qde，i）2∑n-1i=1（Qori，i+1-Qori，i）2（15）

式中：Esm為平滑度，其值越大表示曲線越平滑；Qori，i，Qori，i+1 分別為第i，i+1時段平滑前的入庫流量，m3/s；Qde，i，Qde，i+1分別為第i，i+1時段平滑后的入庫流量，m3/s。

（2） 相關(guān)系數(shù)可以反映兩個隨機(jī)變量相關(guān)程度，記為P。P越接近于1，數(shù)據(jù)相關(guān)性越好。相關(guān)系數(shù)計算公式如式（16）所示。

P=1n-1∑ni=1［fi-uf］［gi-ug］σfσg（16）

式中：uf，ug分別為f（i）和g（i）的均值；σf，σg分別為f（i）和g（i）的標(biāo)準(zhǔn)差；f（i），g（i）分別為i時段平滑前、后的入庫流量。

5 平滑結(jié)果對比分析

本文根據(jù)卡洛特水電站機(jī)組過流特性，結(jié)合電站“徑流模式”流量考核區(qū)間，選取2023年8月和9月的逐時入庫流量作為對象，對比目前常采用的 3 點、5 點、7 點線性平滑法，五點三次函數(shù)擬合平滑法的平滑度和相關(guān)系數(shù)（表2）。

從表2可以看出：4種方法得到的8月和9月入庫流量平滑度從高到低依次是7點線性平滑法、5點線性平滑法、3點線性平滑法、五點三次函數(shù)擬合平滑法。在相關(guān)系數(shù)指標(biāo)上，4種方法平滑后的相關(guān)系數(shù)均在0.89以上；平滑度最高的7點線性平滑法在9月的相關(guān)系數(shù)反而是最低的，說明存在過度平滑現(xiàn)象；五點三次函數(shù)平滑在研究時段內(nèi)相關(guān)系數(shù)指標(biāo)最高，均在0.95以上，說明較好地保留了原始數(shù)據(jù)的特征。

選擇卡洛特水電站“徑流模式”上限流量的時間段8月5～9日、下限流量時間段8月15～19日和枯水期流量時間段9月22～25日的流量平滑過程放大觀察（圖1～3）。

可以看出，所有方法都有平滑效果，基本能準(zhǔn)確反映對應(yīng)洪峰時間，各方法的洪量相等，但過度平滑帶來的“削峰填谷”效應(yīng)在電站運(yùn)行中會帶來一定風(fēng)險。在洪峰流量過度削減情況下，電站調(diào)度若不能及時響應(yīng)可能會發(fā)生庫水位超限的情況。而低流量的“虛高”效應(yīng)也會導(dǎo)致庫水位急劇下降，導(dǎo)致水頭減小，影響發(fā)電效益。

卡洛特水電站在“徑流模式”下水庫需保持出入庫平衡，水位控制在461 m，過度平滑會使得電站運(yùn)行違反規(guī)定，引起購電方罰款。 五點三次平滑算法可以將反推入庫流量數(shù)據(jù)中“鋸齒”現(xiàn)象的干擾成分消除，保持?jǐn)?shù)據(jù)的變化趨勢，可以認(rèn)為該方法的平滑結(jié)果較理想。

6 結(jié) 論

本文根據(jù)卡洛特水電站特性，分析了庫水位的觀測誤差引起的反推入庫流量誤差，并基于3點線性平滑法、5點線性平滑法、7點線性平滑法、五點三次函數(shù)擬合平滑法，對卡洛特水電站2023年8月和9月入庫流量中的“鋸齒”問題進(jìn)行平滑，同時采用平滑度、相關(guān)系數(shù)指標(biāo)對其平滑效果進(jìn)行評價，得到結(jié)論如下。

（1） 卡洛特水電站正常蓄水位461 m時，0.01 m的水位誤差，其流量差值為15 m3/s；流域10月至次年2月的多年平均流量約270 m3/s，按照5%的預(yù)報誤差，允許入庫流量預(yù)報誤差往往只有13 m3/s左右，庫水位觀測誤差對預(yù)報工作極為不利。

（2） 卡洛特水電站8月和9月入庫流量平滑度從高到低依次是7點線性平滑法、5點線性平滑法、3點線性平滑法、五點三次函數(shù)擬合平滑法；五點三次平滑法相關(guān)系數(shù)在0.95以上，表現(xiàn)最佳。綜合分析認(rèn)為，五點三次平滑算法不僅可以消除反推入庫流量數(shù)據(jù)中“鋸齒”現(xiàn)象的干擾成分，還可以更全面地保證原入庫流量的水文特征，其平滑結(jié)果較理想，在調(diào)節(jié)庫容小、入庫流量預(yù)報精度要求高的中小水庫具有適用價值。

參考文獻(xiàn)：

［1］ 閔要武，王俊，陳力.三峽水庫入庫流量計算及調(diào)洪演算方法探討［J］.人民長江，2011，42（6）：49-52.

［2］ 張俊，閔要武，陳力，等.基于動庫容曲線的水庫調(diào)洪高水位查算方法研究［J］.水文，2016，36（5）：63-68.

［3］ 唐海華，陳森林，趙云發(fā)，等.三峽水庫入庫流量計算方法研究［J］.中國農(nóng)村水利水電，2008（4）：26-28.

［4］ 陳建，李允軍，王建平，等.基于梯度自適應(yīng)優(yōu)化的入庫流量數(shù)據(jù)平滑算法研究［J］.人民長江，2021，53（3）：92-97.

［5］ 張俊，廖勝利，程春田，等.基于最小二乘曲線的水電站入流平滑處理［J］.水電能源科學(xué)，2011，29（4）：55-56，154.

［6］ 楊海星，黎浩，曹凈.基于樣條函數(shù)的監(jiān)測數(shù)據(jù)平滑算法與Matlab實現(xiàn)［J］.長春工程學(xué)院學(xué)報（自然科學(xué)版），2018（3）：92-97，120.

［7］ 楊金標(biāo)，舒凱，張后來，等.水庫水量平衡計算中消除水位跳變影響的過濾算法［J］.人民長江，2019，50（8）：98-102.

［8］ 徐楊，劉亞新，汪濤，等.基于LSTM的三峽水庫短期上游水位預(yù)測方法研究［J］.水利水電快報，2022，43（10）：13-18.

［9］ 李征，房宏才，柯熙政，等.滑動平均法在MEMS陀螺信號趨勢項提取中的應(yīng)用［J］.電子測量與儀器學(xué)報，2019，33（7）：43-49.

［10］ 蔡天凈，唐瀚.Savitzky-Golay平滑濾波器的最小二乘擬合原理綜述［J］.數(shù)字通信，2011（1）：63-68，82.

［11］ 李匡，柳俊，張子平，等.基于自適應(yīng)權(quán)重的滑動平均入庫流量修正方法［J］.中國防汛抗旱，2021，31（5）：20-24.

［12］ 武煒，陳標(biāo)，吳劍鋒，等.基于五點三次平滑算法的入庫流量反推研究［J］.水利水電技術(shù)，2013，44（12）：100-102.

［13］ 金錚，張行南，夏達(dá)忠，等.基于多種平滑算法的三峽水庫反推入庫流量消噪［J］.中國農(nóng)村水利水電，2021（3）：48-53.

［14］ 李凱凱.張峰水庫入庫流量推算方法與應(yīng)用研究［J］.山西水土保持科技，2021（2）：38-40，48.

（編輯：江 文）

Smoothing correction method for inflow into Karot Hydropower Station

XU Zhi，LIU Tianyong，HE Song

（Jingjiang Bureau of Hydrology and Water Resources Survey，Bureau of Hydrology of Changjiang Water Resources Commission，Jingzhou 434000，China）

Abstract：

The \"jagged\" issue exists in the inflow into the Karot Hydropower Station calculated by reversing the static reservoir capacity.We quantitatively analyzed the impact of the Karot Reservoir′s water level measurement error on the inflow and smoothed the jagged fluctuations in the inflow for August and September in 2023 by using the 3-point linear smoothing method，5-point linear smoothing method，7-point linear smoothing method，and the cubical smoothing algorithm with 5-point approximation.The results showed that when the normal storage level of Karot Hydropower Station was 461 meters，a water level error of 0.01 meters caused a flow error of 15 m3/s，exceeding the allowable inflow forecast error （13 m3/s）.The reservoir water level measurement error had an adverse effect on the flow forecasting work.The quintic smoothing algorithm can effectively eliminate the \"jagged\" interference components in the reverse-calculated inflow data，comprehensively retain the hydrological characteristics of the original inflow，and the smoothing results were relatively more ideal.The research results can provide a reference for the small and medium reservoirs with small adjustable capacities and high requirements for inflow forecasting accuracy in reversing the inflow calculation.

Key words：

inflow into reservoir； calculation by reversing static reservoir capacity； flow smoothing correction； cubical smoothing algorithm with 5-point approximation； Karot Reservoir； Pakistan