趙 敏,孫 勇
(貴州大學(xué) 喀斯特環(huán)境與地質(zhì)災(zāi)害重點(diǎn)實(shí)驗(yàn)室,貴陽(yáng) 550002)
大壩是極其重要的基礎(chǔ)設(shè)施,是綜合防洪工程體系和水力發(fā)電工程的重要組成部分,往往規(guī)模較大,其價(jià)值體現(xiàn)在防洪、發(fā)電、灌溉、航運(yùn)、旅游、水產(chǎn)養(yǎng)殖等多個(gè)方面[1]。然而,一旦大壩發(fā)生潰決,不僅會(huì)造成社會(huì)經(jīng)濟(jì)損失,更者會(huì)導(dǎo)致人員傷亡、生態(tài)環(huán)境問(wèn)題等。國(guó)內(nèi)外潰壩事故均有發(fā)生,我國(guó)1975年的板橋事件造成河南省29個(gè)縣市、1 100萬(wàn)人受災(zāi),2.6萬(wàn)人死亡,直接經(jīng)濟(jì)損失近百億元[2]。1979年印度曼朱二號(hào)水庫(kù)潰壩事件,造成5 000~10 000人的死亡,經(jīng)濟(jì)損失巨大。
大壩潰口形成過(guò)程和泄洪過(guò)程線(可細(xì)分為潰口特征參數(shù)、水庫(kù)蓄水量和下泄流量)的分析對(duì)已發(fā)洪水的危害評(píng)估至關(guān)重要[3]。分析大壩破壞以及泄洪水位線有許多模型,其中DAMBRK和HEC-1應(yīng)用較廣,Singh和Snorrason采用這兩個(gè)模型對(duì)20個(gè)潰壩案例進(jìn)行了關(guān)鍵潰口參數(shù)的敏感度研究[4]。Petrascheck和Sydler提出了一個(gè)數(shù)學(xué)模型來(lái)模擬洪水的外泄過(guò)程,研究潰口泄量、洪水位和洪水抵達(dá)時(shí)間對(duì)潰口寬度及潰口形成時(shí)間的變化的敏感度。該研究表明,可靠的估計(jì)大壩附近及下游區(qū)域洪水情況,準(zhǔn)確地預(yù)測(cè)潰口參數(shù)是必要的[5]。
對(duì)潰口參數(shù)的評(píng)估有多種方式,Wahl總結(jié)為3種:類比法、基于物理原理的模型和預(yù)測(cè)方程[6]。類比法的前提是大壩間具有相似的幾何和水文情況;基于物理原理的模型運(yùn)用水文學(xué)、泥沙運(yùn)動(dòng)學(xué)模擬潰決過(guò)程,往往因?yàn)槿狈?duì)潰決機(jī)理的了解而過(guò)于簡(jiǎn)化;而預(yù)測(cè)方程法在今天依然使用最為廣泛。預(yù)測(cè)方程參數(shù)一般是通過(guò)使用數(shù)據(jù)庫(kù)中對(duì)有限實(shí)例數(shù)據(jù)回歸分析所得最小擬合誤差作為唯一指標(biāo)來(lái)得到的,而過(guò)度擬合常常導(dǎo)致公式十分復(fù)雜并影響預(yù)測(cè)方程的準(zhǔn)確性。
本文運(yùn)用貝葉斯推論的方法,對(duì)ICOLD統(tǒng)計(jì)的101個(gè)案例數(shù)據(jù)進(jìn)行分析,選擇最優(yōu)模型,從而提出可以平衡擬合能力和復(fù)雜性的預(yù)測(cè)方程。
潰口參數(shù)包括潰口的幾何參數(shù)、形成時(shí)間參數(shù)和水文參數(shù)。將潰壩理想化為一個(gè)不規(guī)則四邊形,見(jiàn)圖1。幾何參數(shù)即包括潰口深度(潰口高度)Hb(即壩頂至潰口反彎點(diǎn)的距離)、潰口頂部寬度Bt、潰口底部寬度Bb、潰口平均寬度Bave、潰口側(cè)坡因子Z。只要知道其中的任意3個(gè)參數(shù)就足以描述完整潰口的幾何情況。
圖1 理想化大壩潰口參數(shù)圖
潰決時(shí)間Tf,Singh和Snorrason定義其為潰裂開(kāi)始到完成的持續(xù)時(shí)間[4]。潰口開(kāi)始時(shí)間起始于水流首次流過(guò)大壩,此時(shí)流出水量較小,隨著出流量和侵蝕迅速增加,潰口形成,潰口開(kāi)始時(shí)間終止。潰口形成時(shí)間指大壩上游面首次潰決到潰口完全形成的時(shí)間。兩個(gè)參數(shù)較難區(qū)分。洪峰泄量Qp,取決于水庫(kù)的水位記錄或?qū)?cè)坡的測(cè)量。
基于潰壩的破壞機(jī)理,本文中涉及的參數(shù)只有4個(gè):大壩寬度Hd、潰口水頭Hw(潰口深度Hb)、大壩水庫(kù)的容積Vd、潰口容積Vw。
預(yù)測(cè)方程可以有加的形式和積的形式。
考慮到測(cè)量誤差的存在,輸出值以預(yù)測(cè)值和誤差項(xiàng)的總和的形式:
(1)
(2)
積的形式:
(3)
為了更方便地進(jìn)行回歸分析,對(duì)式(3)兩邊取對(duì)數(shù)從而轉(zhuǎn)化為加的形式:
lny=β0lnx0+β1lnx1+…+βNxlnxNx
(4)
近10年來(lái),貝葉斯推論被逐漸廣泛運(yùn)用在工程中。為避免過(guò)度擬合導(dǎo)致公式過(guò)于復(fù)雜和影響預(yù)測(cè)方程的準(zhǔn)確性,Yan等人提出了基于貝葉斯推論的模型選擇框架[7]。該方法有助于在一系列候選模型中挑選出一個(gè)具有最合理的輸入輸出關(guān)系的模型。貝葉斯模型選擇框架與現(xiàn)有的預(yù)測(cè)方程參數(shù)只考慮最小擬合誤差不同,旨在通過(guò)衡量擬合準(zhǔn)確度(goodness-of-fit)和魯棒性(robustness)兩者,從而挑選出最合適的模型。模型合理性按下式計(jì)算:
j=1,2,…,NM
(5)
(6)
(7)
(8)
(9)
式中:p(Mj|U)為用戶對(duì)模型Mj的最初合理性的判斷;p(D|Mj,U)為給定D時(shí)計(jì)算模型Mj的一個(gè)量;θj為在模型Mj中定義的參數(shù)向量,包含著未知擬合參數(shù)αj和σε;p(θj|Mj,U)為用戶指定的先驗(yàn)概率密度函數(shù);p(D|θj,Mj,U)為模型Mj準(zhǔn)確度的似然函數(shù),又寫為p(D|θj,Mj) ;Jg(α|D,Mj)為在擬合系數(shù)α下的準(zhǔn)確度函數(shù);Nj為模型Mj中不確定參數(shù)數(shù)量。
式(7)p(D|Mj,U)在實(shí)際計(jì)算中,經(jīng)常用拉普拉斯的漸進(jìn)展開(kāi)來(lái)近似計(jì)算:
(10)
奧克漢因子被用于評(píng)估模型Mj魯棒性,如果模型參數(shù)對(duì)測(cè)量值很敏感并非是好事,這意味著擬合系數(shù)一個(gè)小的變化就會(huì)引起預(yù)測(cè)值出現(xiàn)一個(gè)對(duì)應(yīng)的很大的改變。一個(gè)模型的奧克漢因子越高,魯棒性就越大,模型就越強(qiáng)健,對(duì)測(cè)量誤差的敏感度就越低[8]。
(11)
(12)
(13)
(14)
預(yù)測(cè)方程的廣義公式可以寫成如下形式:
(15)
式中:yi為向量y=[Hb,Bt,Bave,Qp,Tf]預(yù)測(cè)的破壞參數(shù);xj為向量x=[Hd,Hw,Vd,Vw]中控制變量的選擇;Pi為模型Mj控制變量的集合。
4個(gè)潰口參數(shù)通過(guò)排列組合形成15個(gè)模型種類(表1)。
表1 模型種類列表
見(jiàn)表2。
從表2可以看出,Bave與Hw有較大的相關(guān)性(ρ=0.61),與Vd和Vw次之(ρ大于0.5),說(shuō)明潰口的平均寬度與潰口水頭的關(guān)系較大。同時(shí),潰口水頭Hw相比較其他3個(gè)參數(shù),潰口頂部寬度Bt與其相關(guān)性最大(ρ=0.65);Hb與Hd(ρ=0.92)和Hw(ρ=0.97)都有很大相關(guān)性,說(shuō)明潰口深度與大壩寬度和潰口水頭有密切關(guān)系;洪峰泄量Qp則是與潰口水頭Hw和潰口容積Vw有相對(duì)較大且相同的相關(guān)性(ρ=0.66);對(duì)于潰決時(shí)間Tf而言,大壩水庫(kù)容積Vd和潰口容積Vw與其相關(guān)性較大(ρ=0.51)。
表2 相關(guān)系數(shù)表
4個(gè)控制變量與潰口參數(shù)之間的關(guān)系亦可通過(guò)散點(diǎn)圖得出。圖2以Hb為例,Vd、Vw、Hd、Hw的數(shù)據(jù)進(jìn)行了標(biāo)準(zhǔn)化處理。
圖2中可以看出,Hd和Hw與Hb呈明顯的正相關(guān)關(guān)系,同樣驗(yàn)證了潰口深度與大壩寬度和潰口水頭有密切關(guān)系。
每個(gè)候選模型的合理性和奧克漢因子值都可以通過(guò)基于貝葉斯推論的模型選擇框架獲得。圖3為Hb加的形式下的預(yù)測(cè)結(jié)果。圖3中同時(shí)包含只代表模型準(zhǔn)確度的大小的確定系數(shù)R2和奧克漢因子的對(duì)數(shù)形式,以便做出比較。
圖2 潰口深度Hb與各變量關(guān)系圖
圖3 Hb的模型分析圖
從圖3中可以看出,模型15具有最大的R2值。然而模型15的奧克漢因子對(duì)數(shù)值最低,這表明雖然該模型有很高的準(zhǔn)確性,擬合誤差最小,但是魯棒性最低。也就是說(shuō),模型15對(duì)測(cè)量值的變化很敏感,不利于預(yù)測(cè)。模型5不僅有最大R2值,而且?jiàn)W克漢因子對(duì)數(shù)值較小,可以作為最優(yōu)的模型選擇。同理,可以找出其他潰口參數(shù)的最優(yōu)模型。
潰口參數(shù)相關(guān)性最大的變量和最優(yōu)模型的R2值的關(guān)系見(jiàn)表3。從表3中可以看出,具有最大準(zhǔn)確度,即R2值最大的模型并不是最理想的。
表3 R2值與相關(guān)變量關(guān)系
貝葉斯推論近10年在工程中的應(yīng)用比較廣泛,也體現(xiàn)出其優(yōu)越性。預(yù)測(cè)方程參數(shù)通常使用最小擬合誤差作為唯一指標(biāo),模型選擇過(guò)程中只考慮準(zhǔn)確度,過(guò)度擬合不僅使公式過(guò)度復(fù)雜,也會(huì)導(dǎo)致預(yù)測(cè)結(jié)果不理想。本文運(yùn)用貝葉斯推論分析了水壩潰口參數(shù)與變量之間的關(guān)系,驗(yàn)證了模型選擇時(shí)不能只考慮準(zhǔn)確度,導(dǎo)致過(guò)度擬合,還應(yīng)該考慮模型的魯棒性,最后根據(jù)最優(yōu)模型得出的預(yù)測(cè)方程對(duì)測(cè)量值變化較不敏感,更具有可預(yù)測(cè)性。