張 田,孫銘悅,何 理
(大連海事大學(xué) 交通運(yùn)輸工程學(xué)院,遼寧 大連 116026)
抗風(fēng)設(shè)計是橋梁設(shè)計的重要組成部分,橋址區(qū)設(shè)計風(fēng)速的推算是橋梁抗風(fēng)設(shè)計的基礎(chǔ),是保障橋梁抗風(fēng)安全首要解決的問題.目前,推算橋址區(qū)風(fēng)速的一般方法是按實(shí)測資料序列計算統(tǒng)計參數(shù),包括均值、標(biāo)準(zhǔn)差σ、離差系數(shù)Cv、偏差系數(shù)Cs,由以上參數(shù)確定一條理論頻率分布曲線(極值-I型曲線或Pearson-III型曲線)[1-4],再計算不同重現(xiàn)期下的設(shè)計風(fēng)速值.然而,由于實(shí)測資料觀測年限與總體相比相差很遠(yuǎn),而且經(jīng)驗(yàn)頻率點(diǎn)分布不是一條光滑連續(xù)曲線,直接由統(tǒng)計參數(shù)得到的理論頻率分布曲線與經(jīng)驗(yàn)頻率點(diǎn)偏離較大.因此,借鑒水文學(xué)中推算設(shè)計流量的方法[5],即適線法,通過調(diào)整統(tǒng)計參數(shù)值,以一定的準(zhǔn)則選配出一條與經(jīng)驗(yàn)頻率點(diǎn)擬合最好的理論頻率分布曲線,進(jìn)而推算設(shè)計風(fēng)速.
由于氣象記錄中的極大風(fēng)速值是觀測期內(nèi)的極大值,工程中可能遇到的極端風(fēng)速情況不能簡單地取氣象資料中觀測到的極大值,而需要利用最大風(fēng)速的頻率分布來推算橋位處的風(fēng)速極值.常用的風(fēng)速頻率分布曲線形式有極值-I型和Pearson-III型[6-9].
極值-I型分布函數(shù)為
其超過保證率的函數(shù)為
式中,x為隨機(jī)變量;xp為頻率p時的隨機(jī)變量取值;α和μ分別為分布的尺度參數(shù)和位置參數(shù),只要利用已知的樣本序列x1,x2,…,xn-1,xn,合理估計出參數(shù)α和μ的值,則F(x)就可確定.
重現(xiàn)期為T年(累積頻率為p=1/T)的風(fēng)速值可表示為
對于參數(shù)α和μ的估計,可以采用矩法.由數(shù)學(xué)期望和方差的計算公式,進(jìn)一步推導(dǎo)得
一般采用有限樣本容量的均值和標(biāo)準(zhǔn)差作為理論值Ex和σ的近似估計.
設(shè)風(fēng)速值為隨機(jī)變量x,則Pearson-III型曲線的密度函數(shù)為
式中,Γ(α)為α的伽馬函數(shù),α、β和x0為曲線的參數(shù),可以用樣本序列的3個統(tǒng)計參數(shù)即均值、離差系數(shù)Cv和偏差系數(shù)Cs來表示
式中,
式中,n為樣本數(shù),σ為樣本推算總體的標(biāo)準(zhǔn)差.
計算時不計算α、β和x0,只計算、Cv和Cs,以及給定工程設(shè)計的標(biāo)準(zhǔn),即重現(xiàn)期T的設(shè)計要求,則累積頻率p=1/T,于是有
對式(8)進(jìn)行積分計算,就可以得到頻率分布曲線縱坐標(biāo)值,即累積頻率p對應(yīng)的隨機(jī)變量值為
式中,Φ為離均系數(shù),是頻率和偏差系數(shù)的函數(shù),
適線法為選定最優(yōu)理論頻率曲線的方法,即通過不斷優(yōu)化調(diào)整統(tǒng)計參數(shù),繪制多條理論頻率曲線,選擇其中與實(shí)測資料吻合最好的曲線作為最終的理論頻率曲線.適線之前需要繪制經(jīng)驗(yàn)頻率點(diǎn)圖,其步驟為:將已知的風(fēng)速序列不論年序按數(shù)值大小遞減次序排列x1≥x2≥…≥xi≥…≥xn-1≥xn;按Weibull公式,即式(10),計算各風(fēng)速值的經(jīng)驗(yàn)頻率,并在海森機(jī)率格紙上繪出經(jīng)驗(yàn)頻率點(diǎn)群.
式中,P為風(fēng)速系列中第i項(xiàng)隨機(jī)變量的經(jīng)驗(yàn)頻率;n為風(fēng)速系列的容量;i為計算隨機(jī)變量的序號.
根據(jù)概率論與數(shù)理統(tǒng)計相關(guān)計算方法,適線準(zhǔn)則可以采用離(殘)差平方和最小準(zhǔn)則(OLS準(zhǔn)則)和離(殘)差絕對值和最小準(zhǔn)則(ABS準(zhǔn)則).
①由OLS準(zhǔn)則得目標(biāo)函數(shù)為
②由ABS準(zhǔn)則得目標(biāo)函數(shù)為
對極值-I型分布曲線,由上述兩種不同的適線準(zhǔn)則給出的目標(biāo)函數(shù),將式(3)代入式(11)、式(12)得
以O(shè)LS準(zhǔn)則為例給出計算過程,若Q1取最小值,則令和?Q1?α=0,則有
式(15)是關(guān)于α和μ的非線性方程組,可采用迭代法求解.令
設(shè)迭代的初始值為 1α和1μ,將式(4)計算結(jié)果作為初始值.將 1φ和 2φ在(1α,1μ)的鄰域作Taylor展開,略去高次項(xiàng)后的迭代方程組為
直到滿足下式時迭代終止,即
式中,1ε和 2ε為預(yù)先給定的允許誤差.
對Pearson-III型分布曲線,將式(9)代入式(11)、式(12)得
以O(shè)LS準(zhǔn)則為例,若Q1取最小值,則令和,有
由于經(jīng)驗(yàn)頻率pi為已知,所以給定一個Cs,可以查表求出與pi相對應(yīng)的Φi,由式(22)按上述迭代求解方法可解出此時使Q1最小的和Cv值.為此,將Cs從0開始,按步長Δ增長(取Δ為0.001),從而得到不同組值,每一組值都對應(yīng)一條Pearson-III型理論頻率曲線,也對應(yīng)一個按式(20)計算的Q1值,最后選用使Q1最小的一組作為計算結(jié)果.
對極值-I型或Pearson-III型分布曲線的參數(shù),雖然可以由上述的迭代方法求解,但從推導(dǎo)過程可以看出,該方法需要推求偏導(dǎo)數(shù)矩陣,且涉及到矩陣的求逆過程,求解復(fù)雜,收斂速度較慢.可以采用模擬退火算法[10-13],該方法對優(yōu)化目標(biāo)要求較低,一般只需要知道其數(shù)值關(guān)系,不受函數(shù)連續(xù)性、光滑性的限制,避免了大量求導(dǎo)數(shù)計算.
采用模擬退火算法對于一般無約束優(yōu)化求最小值問題的求解過程[14-17]如下.
(1)給定起始溫度T0、終止溫度Tfinal和初始解x0.
(2)在x0的鄰域內(nèi)模擬產(chǎn)生隨機(jī)擾動Δx.
(3)計算擾動引起的目標(biāo)函數(shù)值的變化ΔE.
(4)若ΔE≤0,接受新值x→x0+Δx,轉(zhuǎn)向步驟(5);否則,若exp(ΔE/T)>rand(0,1),rand(0,1)表示0~1內(nèi)的隨機(jī)數(shù),也接受新值x→x0+Δx,轉(zhuǎn)向步驟(5);否則就拒絕x→x0+Δx,轉(zhuǎn)向步驟(6).
(5)按退火策略進(jìn)行降溫T→update(T).
采用指數(shù)退溫方式,設(shè)tk、tk-1分別為第k和k-1次迭代時的溫度,t0為初始溫度,η為參數(shù),tk=η·tk-1.
(6)是否滿足T>Tfinal?若滿足,則輸出求解結(jié)果,結(jié)束迭代計算;否則,轉(zhuǎn)向步驟(2).
作為實(shí)例,以貴州省興仁縣氣象站1966-1998年的實(shí)測風(fēng)速資料按不同的推斷方法推算氣象站10 m高度處的設(shè)計風(fēng)速值,經(jīng)高度換算后的實(shí)測資料[7]列于表1,其中1975年、1992年的資料缺失.
表1 興仁縣10 m高度,年最大10 min平均風(fēng)速 Tab.1 annual maximum 10 min average wind speed data at 10 m altitude in Xingren County
根據(jù)表1的觀測資料可以計算經(jīng)驗(yàn)頻率點(diǎn)群圖,同時計算各統(tǒng)計參數(shù)的初始值,即α0、μ0,、Cv0、Cs0.為了使理論頻率曲線與經(jīng)驗(yàn)頻率點(diǎn)群吻合得更好,采用文中提出的模擬退火算法,按2種適線準(zhǔn)則逐步調(diào)整統(tǒng)計參數(shù)以獲得最優(yōu)解,計算的統(tǒng)計參數(shù)和目標(biāo)函數(shù)結(jié)果列于表2.
表2 不同方法獲得的參數(shù)及目標(biāo)函數(shù)值 Tab.2 parameters and objective function values obtained by different methods
從表2可以發(fā)現(xiàn),無論采用那種準(zhǔn)則,應(yīng)用模擬退火算法計算出的統(tǒng)計參數(shù)均優(yōu)于直接計算方法.利用上述方法求出的統(tǒng)計參數(shù),來擬合理論頻率曲線,見圖1,進(jìn)而求出不同重現(xiàn)期的設(shè)計風(fēng)速列于表3.分布曲線的吻合程度可以通過直接觀察分布曲線與經(jīng)驗(yàn)頻率點(diǎn)群圖像或判定系數(shù)檢驗(yàn)、誤差平方和檢驗(yàn)等方法確定.在圖1中將理論頻率分布曲線與經(jīng)驗(yàn)頻率點(diǎn)群分布比較,可以直觀地看出Pearson-III型曲線作為理論頻率分布曲線與經(jīng)驗(yàn)頻率點(diǎn)群吻合得更好.
圖1 經(jīng)驗(yàn)頻率點(diǎn)及擬合頻率 Fig.1 empirical frequency point and fitting frequency
表3 不同重現(xiàn)期時采用模擬退火算法推算的設(shè)計風(fēng)速 Tab.3 design wind speed calculated by simulated annealing algorithm in different return periods
由表3和圖1可以發(fā)現(xiàn),若采用Pearson-III型曲線作為理論頻率分布曲線,運(yùn)用模擬退火算法計算出的統(tǒng)計參數(shù),所得的不同重現(xiàn)期時的設(shè)計風(fēng)速差別不大,可見模擬退火算法計算統(tǒng)計參數(shù)是有效可行的.
(1)提出適線法,通過不斷調(diào)整統(tǒng)計參數(shù),使理論頻率曲線與經(jīng)驗(yàn)頻率點(diǎn)群盡可能吻合良好,同時給出相應(yīng)的適線準(zhǔn)則及采用模擬退火算法優(yōu)化統(tǒng)計參數(shù)的步驟.
(2)通過實(shí)例計算,得到:①基于經(jīng)驗(yàn)頻率點(diǎn)群,采用適線法優(yōu)化統(tǒng)計參數(shù)比由觀測數(shù)據(jù)直接計算統(tǒng)計參數(shù)能使理論頻率曲線與經(jīng)驗(yàn)頻率點(diǎn)群擬合得更好;②給出適線準(zhǔn)則:離(殘)差平方和最小準(zhǔn)則(OLS準(zhǔn)則)和離(殘)差絕對值和最小準(zhǔn)則(ABS準(zhǔn)則),兩種目標(biāo)函數(shù)都能很好用于優(yōu)化統(tǒng)計參數(shù);③采用模擬退火算法,可避免求解統(tǒng)計參數(shù)時大量的求導(dǎo)運(yùn)算,迭代過程簡單,對兩種適線準(zhǔn)則,都能使目標(biāo)函數(shù)快速收斂.