基于粒子群算法的新安江模型參數(shù)的優(yōu)化

溫海燕

(河北省唐秦水文水資源勘測局,河北 唐山 063000)

正確的來水量預(yù)報是防汛抗旱和水利水電工程設(shè)計、施工、調(diào)度、管理等的重要依據(jù)。水文預(yù)報模型具有操作簡單、計算速度快的特點(diǎn),通常能在很短的時間內(nèi)得到較為準(zhǔn)確的結(jié)果。特別是伴隨著計算機(jī)技術(shù)的飛速發(fā)展,水文預(yù)報模型日益普及。但運(yùn)用水文模型預(yù)報的精度如何,則很大程度上取決于模型的參數(shù),因此,參數(shù)優(yōu)化非常重要。

1 基本資料

選定灤河流域唐秦地區(qū)大黑汀水庫(灤河干流)、桃林口水庫(一級支流青龍河)至灤縣水文站(灤河出口控制站)區(qū)間為研究區(qū)域,流域面積3550 km2。根據(jù)研究區(qū)域現(xiàn)有觀測站網(wǎng)自測報系統(tǒng)站網(wǎng)規(guī)劃,將研究區(qū)按照泰森多邊形法則劃分為 13個流域單元,分別為七道河、上營、遷西、東荒峪、太平寨、肖家營子、冷口、桑園、遷安、夏官營、桃林口、盧龍和灤縣。分別統(tǒng)計研究區(qū)域1979～2002年歷年用水量資料、各單元逐日降水量資料、灤縣站蒸發(fā)量資料及大黑汀水庫、桃林口水庫出庫和灤縣水文站年徑流量系列,作為進(jìn)行模型參數(shù)優(yōu)化的基本資料。

2 模型選擇

研究區(qū)域地處華北平原灤河流域的燕山迎風(fēng)區(qū),大部分地區(qū)呈半濕潤、半干旱狀態(tài),故預(yù)報模型選擇目前廣泛應(yīng)用于水文預(yù)報工作,適用于半濕潤半干旱地區(qū)的新安江模型[1],[2]。新安江模型主要包含以下參數(shù)[3],[4],[5]：

1)蒸散發(fā)參數(shù)K,WUM,WLM,C;K為蒸散發(fā)能力和蒸發(fā)皿實測值之比;WUM為上層土壤水容量,對蒸散發(fā)計算影響顯著;WLM為下層土壤水容量,常取60～90 mm;C為深層蒸散發(fā)系數(shù),與流域上深根植物的覆蓋面積比例有關(guān)。

2)產(chǎn)流參數(shù)IM,W M,B,FC,KF,BF,KG,KI;IM為不透水面積比例;WM為張力水容量,分為上層、下層與深層3層,分別為WUM,WLM,WDM;B為張力水蓄水容量曲線的方次,與流域面積有關(guān);FC為穩(wěn)定下 滲率,唯一地取決于水力傳導(dǎo)度k;F為滲透系數(shù),反映土壤含水量對下滲的影響;BF為反映流域上下滲率分布的不均勻性。

3 粒子群算法

粒子群優(yōu)化算法(Particle Swarm Optimization,PSO)是由 Kennedy和 Eberhart于 1995年提出的一種模擬生物群落行為的進(jìn)化算法[6]。該算法本質(zhì)上是一種多代理算法,是一種基于迭代的優(yōu)化工具,近年來廣泛應(yīng)用于函數(shù)優(yōu)化、車間調(diào)度等問題[7]。

3.1 PSO的基本原理

假設(shè)在一個D維的目標(biāo)搜索空間中,有m個粒子組成一個群落,其中第i個粒子表示為一個D維向量=(xi1,xi2,…,xiD),i=1,2,…,m,第i個粒子的“飛翔”速度也是一個D維向量,記為(vi1,vi2,…,viD)。每個粒子都有一個由被優(yōu)化的函數(shù)決定的適應(yīng)值(Fitness Value),并且知道自己到目前為止發(fā)現(xiàn)的最好位置(pbest)和現(xiàn)在的位置每個粒子的位置就是一個潛在的解。將?代入一個目標(biāo)函數(shù)就可以計算出適應(yīng)值,根據(jù)適應(yīng)值的大小衡量?的優(yōu)劣。這個可以看作是粒子自己的飛行經(jīng)驗。除此之外,每個粒子還知道到目前為止整個群體中所有粒子發(fā)現(xiàn)的最好位置(gbest是在 pbest中的最好值)。這個可以看作是粒子同伴的經(jīng)驗。每個粒子使用下列信息改變自己的當(dāng)前位置：1)當(dāng)前位置;2)當(dāng)前速度;3)當(dāng)前位置與自己最好位置之間的距離;4)當(dāng)前位置與群體最好位置之間的距離。

PSO是一種基于迭代的優(yōu)化工具。對于第k次迭代,每個粒子的變化：

式中,i=1,2,…,m,m是該群體中粒子的總數(shù);vkid 為第k次迭代粒子i飛行速度矢量的第d維分量;xkid為第k次迭代粒子i位置矢量的第d維分量;pid為粒子i個體最好位置pbest的第d維分量;pgd為群體最好位置gbest的第d維分量;c1,c2為權(quán)重因子;rand()為隨機(jī)函數(shù),產(chǎn)生 [0,1]的隨機(jī)數(shù)。

3.2 引入慣性權(quán)重w的粒子群算法計算原理

式(1)是根據(jù)粒子以前的速度和粒子現(xiàn)在位置與自己曾經(jīng)最好的位置之間的距離及當(dāng)前位置和群體最好位置之間的距離三者的和,來計算粒子新的速度。Shi與Eberhart的研究發(fā)現(xiàn),式(1)的第一部分vid由于具有隨機(jī)性且其本身缺乏記憶能力,有擴(kuò)大搜索空間、探索新的搜索區(qū)域的趨勢。因此,具有全局優(yōu)化的能力。在考慮實際優(yōu)化問題時,往往希望先采用全局搜索,使搜索空間快速收斂于某一區(qū)域,然后采用局部精細(xì)搜索以獲得高精度的解。

對PSO全局和局部之搜索能力的平衡主要被慣性權(quán)重控制[8],因此,在式(1)的vid前乘以慣性權(quán)重w,用來控制粒子以前速度對當(dāng)前速度的影響,它將影響粒子的全局和局部搜索能力,即

w較大則算法具有較強(qiáng)的全局搜索能力,w較小則算法傾向于局部搜索。選擇一個合適的w可以平衡全局和局部搜索能力,這樣可以以最少的迭代次數(shù)找到最優(yōu)解。

3.3 PSO算法求解最優(yōu)化問題

此處將以圖形的方式說明PSO算法求解最優(yōu)化問題。對一個最優(yōu)化問題而言,就是要找到全局最小(或最大)點(diǎn)。如圖1所示,為某一二維空間中最小化問題,全局最優(yōu)解點(diǎn)在o處。圖中標(biāo)明了第i個粒子從位置i1更新到位置i2的過程。其中v1是此時刻全局極值gbest引起的第i個粒子向gbest方向的速度;v2是此時刻第i個粒子的個體極值pbesti引起的第i個粒子向pbesti方向的速度;v3是此時刻第i個粒子自身具有的速度。在速度v1,v2,v3的共同作用下,最終粒子將以速度v到達(dá)新的位置i2。下一時刻,粒子i將從i2位置出發(fā),以同樣的方式向最優(yōu)點(diǎn)o靠近。如此迭代,所有粒子都將逐步向最優(yōu)點(diǎn)靠近,最終找到最優(yōu)點(diǎn)。

圖1 PSO算法求解最優(yōu)問題過程示意圖

4 粒子群算法優(yōu)化垂直混合產(chǎn)流模型參數(shù)

4.1 優(yōu)選參數(shù)的確定

選出模型中 8個敏感參數(shù),通過粒子群算法對其進(jìn)行優(yōu)選。優(yōu)選參數(shù)的物理意義及取值范圍見表 1。

表1 優(yōu)選參數(shù)的物理意義及取值范圍

4.2 基于粒子群優(yōu)化算法的新安江模型參數(shù)的優(yōu)化

4.2.1 目標(biāo)函數(shù)

利用日模型優(yōu)化參數(shù)時,其目標(biāo)函數(shù)為

式中,Q(i)為歷年實測年徑流量,Q′(i)為歷年模擬年徑流量。

4.2.2 粒子群算法參數(shù)優(yōu)化步驟

步驟1：初始化粒子,將種群中的每一個粒子隨機(jī)地分布在空間中;隨機(jī)地給每個粒子一個初始的速度vi,根據(jù)表 1中給出的各參數(shù)的范圍,得出每個參數(shù)的vmax,vi∈ (-vmax,vmax);根據(jù)每個粒子,調(diào)用新安江模型,計算各年徑流量,按式(4)計算目標(biāo)函數(shù)值,并將各粒子計算出的目標(biāo)函數(shù)值設(shè)為其相應(yīng)的個體極值pbest(即個體極值點(diǎn)的適應(yīng)度值),從個體極值中選出一個最好的作為全局極值gbest(即全局極值點(diǎn)的適應(yīng)度值),并將gbest設(shè)置為該最好粒子的當(dāng)前位置;給定學(xué)習(xí)因子c1和c2。

步驟2：根據(jù)迭代次數(shù)計算慣性權(quán)重w,通過式(3)計算粒子的速度,判斷速度vi是否超出其取值范圍(-vmax,vmax),若超出范圍,則將其設(shè)為邊界值。

步驟3：根據(jù)粒子速度vi和式(2)計算粒子的位置,若粒子的任何一維超出范圍,則返回到步驟 2,依靠隨機(jī)數(shù)重新計算粒子的速度和位置,直至粒子滿足范圍要求。

步驟4：對每個粒子,根據(jù)其在解空間的位置,調(diào)用垂直混合產(chǎn)流模型,進(jìn)行預(yù)報,通過式(4)計算目標(biāo)函數(shù)值,作為其適應(yīng)度。

步驟5：對每個粒子,比較其適應(yīng)度值和粒子所經(jīng)歷過的最優(yōu)解pbest,若前者優(yōu)于后者,則用當(dāng)前值更新粒子記憶的pbest;否則,保持pbest不變。

步驟6：比較當(dāng)前群體最優(yōu)解和gbest標(biāo)識的粒子所找到的最優(yōu)解,若前者優(yōu)于后者,則更新gbest;否則,保持gbest不變。

步驟 7：當(dāng)結(jié)束條件(如達(dá)到最大迭代次數(shù) )滿足時,退出計算,否則返回到步驟 2。

4.3 計算成果及分析

4.3.1 計算成果

選用灤河流域總出口斷面控制站灤縣水文站1979～2002年實測年徑流量資料,采用人工試錯法和粒子群算法率定模型參數(shù)。本次率定粒子個數(shù)分別采用了10個和 20個,迭代次數(shù)分別采用了 50次、100次、200次和400次。計算結(jié)果顯示,不管是10個粒子還是 20個粒子,當(dāng)?shù)螖?shù)增加到 400次時,其目標(biāo)函數(shù)值不再變化,因此認(rèn)定迭代 200次時,已經(jīng)達(dá)到最優(yōu)。參數(shù)最終率定成果見表 2,年徑流模擬成果見表3。

表2 1979～2002年徑流參數(shù)優(yōu)化結(jié)果

表3 1979～2002年徑流量模擬結(jié)果

4.3.2 成果分析

從表3可以看出,用 PSO算法進(jìn)行模型參數(shù)率定,并運(yùn)用此參數(shù)值對年徑流量進(jìn)行模擬,其絕對誤差和相對誤差均比人工試錯法低。若以模擬年徑流量與實測徑流量相對誤差在±20%內(nèi)為合格,則人工試錯法合格率為 41.7%,PSO算法的合格率為 58.3%,精度較人工試錯法有所提高。

若以1979～2002年平均徑流量為界,將這24年分為大水年和小水年,則大水年 PSO優(yōu)化和人工試錯法模擬徑流量的合格率分別為 80%和70%;小水年 PSO優(yōu)化和人工試錯法模擬徑流量的合格率分別為37.5%和 18.8%?？梢?利用 PSO優(yōu)化計算結(jié)果明顯優(yōu)于傳統(tǒng)的人工試錯法。

5 結(jié)束語

粒子群算法是一種全新的優(yōu)化算法,文中采用灤河流域大黑汀—桃林口—灤縣區(qū)間長系列的降水、流量等實際資料,對區(qū)間新安江預(yù)報模型的參數(shù)進(jìn)行了優(yōu)化,并比較了傳統(tǒng)的人工試錯法和粒子群算法對區(qū)間年徑流量的模擬成果。從結(jié)果可以看出,粒子群算法優(yōu)化參數(shù)具有模擬精度高、調(diào)算速度快的特點(diǎn),肯定了該算法在新安江模型參數(shù)優(yōu)化中的適應(yīng)性。

[1]包為民,張建云.水文預(yù)報 [M].第 4版.北京：中國水利水電出版社,2009.

[2]芮孝芳.產(chǎn)匯流理論 [M].北京：水利電力出版社,1995.

[3]王佩蘭,趙人俊.新安江模型(三水源)參數(shù)的客觀優(yōu)選方法 [J].河海大學(xué)學(xué)報,1989,17(4)：65-68.

[4]趙人俊,王佩蘭,胡鳳彬.新安江模型的根據(jù)及模型參數(shù)與自然條件的關(guān)系 [J].河海大學(xué)學(xué)報,1992,20(1)：52-59.

[5]楊小輝,郭菲菲.率定新安江模型(三水源)產(chǎn)匯流參數(shù) [J].貴州水力發(fā)電,2000,14(1)：7-10.

[6]Kennedy,J.and Eberhart,R.Particle Swarm Optimization[C].//IEEE International Conference on Neural Networks(Perth,Australia),IEEE Service Center,Piscataway,N J,1995,IV：1942-1948.

[7]張利彪,周春光,劉小華,等.粒子群算法在求解優(yōu)化問題中的應(yīng)用 [J].吉林大學(xué)學(xué)報,2005,23(4)：385-389.

[8]郭文忠,陳國龍.粒子群優(yōu)化算法的研究進(jìn)展 [J].福建電腦,2005(4)：7-8.