改進的入侵雜草優化算法在泰斯模型中的應用

邱云翔，劉國東，王 亮，李世鈺

(1.四川大學 水利水電學院，成都 610065；2.四川大學 水力學與山區河流開發保護國家重點實驗室，成都 610065)

0 引 言

在地下水資源開發利用的過程中，導水系數T和貯水系數S等水文地質參數是非常重要的基礎數據。為獲得符合實際情況的T和S，需對非穩定流抽水試驗數據進行分析。常用的分析方法有配線法和直線圖解法。數學模型的優化方法因其能高效模擬配線分析和Theis求解也成為反演T和S的分析方法，主要有高斯—牛頓迭代算法[1]，非線性最小二乘法[2]，遺傳算法(GA)[3]，梯度法[4]，粒子群算法(PSO)[5]，GA優化BP神經網絡法[6]等。上述方法均能反演出滿足一定精度要求的T和S，但對T和S的初始取值范圍取較大區間的情形研究較少。2006年，伊朗德黑蘭大學的C Lucas和A R Mehrabian首次提出入侵雜草優化算法(Invasive Weed Optimization,IWO)[7]，該算法具備廣度搜索的優勢，目前已被應用于函數優化，模式識別，自動控制，數據挖掘，決策分析等諸多領域[8-10]。基于此，本文對IWO算法作適當改進，探討其在T和S的初始取值范圍取較大區間時反演T和S的有效性與穩定性。

1 IWO算法

IWO算法的實現步驟如下[7,8]：

(1)確定初始種群大小N0、最大迭代次數itmax、問題的求解維數dim、最大種群數量pmax、最大種子數zmax、最小種子數zmin、非線性調制指數n、初始標準差σinitial、最終標準差σfinal、自變量初始搜索空間xini，根據N0、dim、xini隨機產生初始解。

(2)依據適應度函數計算初始解(父代種群)的適應度值，確定父代種群的最優適應度值和最劣適應度值。

(3)在具備初始解(父代種群)的基礎上依據其適應度值計算繁殖種子的數量，種子數量與其適應度值呈線性關系，適應度值最優個體繁殖出zmax個種子，適應度值最劣個體繁殖出zmin個種子。

(4)父代種群中所有個體種子產生完成后，根據當前所處迭代次數產生標準差σiter，σiter按下式計算：

(1)

(5)種子依據σiter，按照正態分布N(0,σ2)隨機擴散分布在父代個體的周圍形成子代種群。

(6)子代種群形成后，判斷是否達到pmax。若父代種群與子代種群的個體數量之和未達到pmax，則重復(2)～(5)，若超過pmax，將父代和子代個體按適應度值優劣排序，保留適應度值優秀的個體，使種群數量控制在pmax以內。

(7)判斷是否達到itmax，若未達到，則重復(2)～(6)，反之，將適應度值最優個體的信息輸出作為結果。

2 改進的IWO算法運用Theis公式

2.1 算法改進

IWO算法在廣度搜索空間中運用Theis公式時，搜索空間的下限值一般設置為零，為實現算法從廣度搜索逐漸過渡到深度搜索，算法參數中的σinitial和σfinal須分別設置成較大值與較小值。較大的σinitial常使算法在迭代初期擴散分布時產生負解，違背了水文地質參數的非負性，導致算法無法完成深度搜索，難以獲得接近實際的T和S。為此，本文嘗試引入兩個新的算法參數：初始擴散指數ninitial，最終擴散指數nfinal，以擴散指數niter的計算公式代替σiter的計算公式變IWO算法σiter的靜態設置為動態調整。改進的IWO算法步驟如下：

(1) 確定初始種群大小N0、最大迭代次數itmax、問題的求解維數dim、最大種群數量pmax、最大種子數zmax、最小種子數zmin、初始擴散指數ninitial、最終擴散指數nfinal、自變量初始搜索空間xini，根據N0、dim、xini隨機產生初始解。

(2)～(3)同IWO算法。

(4)父代種群中所有個體種子產生完成后，根據當前所處迭代次數產生擴散指數niter，niter按下式計算：

(2)

(5)種子以父代個體除以niter的商作為σiter，按照正態分布N(0,σ2)隨機擴散分布在父代個體的周圍形成子代種群。

(6)～(7)同IWO算法。

2.2 模型設計

2.2.1 求解公式

根據非穩定流抽水試驗資料，選擇Theis公式進行求解，其公式為：

(3)

(4)

(5)

式中：s為抽水影響范圍內任一點任一時刻的水位降深，m；Q為抽水井的流量，m3/s；T為導水系數，m2/s；t為自抽水開始到計算時刻的時間，s；r為計算點到抽水井的距離，m；S為含水層的貯水系數。

從式(4)可以看出，計算W(u)時，需要進行廣義積分，這將極大增加算法運行時間。故本文采用Rajesh Srivasiava[11]給出的W(u)近似表達式計算，其相對誤差小于0.001%，足以滿足生產上允許相對誤差2%左右的要求[12]。該近似表達式為：

u≤1:W(u)=-Inu+a0+a1u+a2u2+a3u3+a4u4+a5u5

(6)

(7)

式中，a0=-0.577 22；a1=0.999 99；a2=-0.249 91；a3=0.055 19；a4=-0.009 76；a5=0.001 08；b0=0.267 77；b1=8.634 76；b2=18.059 02；b3=8.573 33；c0=3.958 50；c1=21.099 65；c2=25.632 96；c3=9.573 32。

2.2.2 適應度函數

以算法搜索的T和S結合求解公式計算各時刻模擬降深值，以各時刻實測降深值與模擬降深值的離差平方和的均值為適應度值[5]：

(8)

T≥0,S≥0

(9)

式中：J為觀測時刻的個數；si為第i個觀測時刻的實測降深值。

2.3 模型測試

2.3.1 測試方法

2.3.2 參數設置

測試采用統一的參數設置標準：導水系數搜索空間下限值Tmin=0，導水系數搜索空間上限值Tmax取隨機值的500倍，貯水系數搜索空間下限值Smin=0，貯水系數搜索空間上限值Smax=0.5。經調試，改進的IWO算法參數取值如下：N0=5；pmax=30；zmax=4；zmin=0；ninitial=4；nfinal=100；itmax=200。GA參數取值如下：種群規模(M=100)；交叉概率(Pc=0.9)；變異概率(Pm=0.7)；遺傳代數(G=400)。

2.3.3 結果分析

圖1展示隨機生成的1 000組T和S，圖2和圖3分別展示改進的IWO算法和GA所計算的1 000組T和S的相對誤差值(大量數據點在坐標原點附近聚集)，表2對改進的IWO算法和GA的測試結果進行了分析。

圖2 1 000組相對誤差值(改進的IWO算法)Fig.2 1 000 relative errors (The modified IWO)

圖3 1 000組相對誤差值(GA)Fig.3 1000 relative errors (GA)

測試指標方 法改進的IWO算法GA相對誤差最小值1.38×10-72.32×10-10相對誤差最大值37.9597.81相對誤差合格率/%92.573.2單組運算時間/s1.093.53

從上述理論測試的表達結果可以看出，當T和S的初始取值范圍取較大區間時，改進的IWO算法合格率較GA更高，單組運算時間更短，相對誤差波動區間更小，說明改進的IWO算法在廣度搜索空間中運用Theis公式是有效的。

3 實例應用

3.1 數據來源

為開展實例應用，選取3個不同尺度的非穩定流抽水試驗。依次為文獻[12]、[13]和[14]中抽水10 min以后的實測降深等數據。

3.2 應用開展

為進一步展現改進的IWO算法在廣度搜索空間中運用Theis公式的優越性，將搜索空間進行擴展：Tmin=0，Tmax取配線法確定結果的1 000倍，Smin=0，Smax=0.5。為判別算法在廣度搜索空間中運用Theis公式的穩定性，將算法在每組實例應用中均運行100次。仍新增GA與改進的IWO算法在同等搜索空間下反演T和S的對比。經調試，改進的IWO算法參數設置如表3，GA參數設置如表4。

表3 改進的IWO算法參數設置Tab.3 Parameters of the modified IWO

表4 GA參數設置Tab.4 Parameters of the GA

3.3 結果分析

改進的IWO算法與其他方法求解結果的對比見表5與圖4，與GA處理上述3組抽水試驗的成果分別見圖5、6。

表5 不同方法的求解結果對比Tab.5 Results of different methods

圖4 82組實測降深的相對誤差值Fig.4 82 relative errors of drawdowns

從表5可以看出，當T與S較小時，本文方法的適應度值(fitness)略低于常用配線法及直線圖解法，各觀測時刻實測降深的平均相對誤差(MRE)與常用方法接近，當T與S較大時，fitness和MRE均較常用方法明顯降低。從圖4可以看出，在文獻[13]中，本文方法的多數觀測時刻實測降深的相對誤差(RE)介于最小值與最大值之間，在文獻[12]和[14](尤其文獻[14])中，多數RE的最小值由本文方法獲得。綜上所述，改進的IWO算法可以在廣度搜索空間中反演T和S，并且當T和S較大時其求解結果較常用方法更符合實際情況。

圖5 300個適應度值(改進的IWO算法)Fig.5 300 fitness values (The modified IWO)

圖6 300個適應度值(GA)Fig.6 300 fitness values (GA)

從與GA的對比還可看出，當T和S的初始取值范圍取較大區間時，GA需要更大的種群規模及更多的迭代次數(更長的運行時間)來完成相應的求解，且在100次的運算中適應度值存在明顯波動，尤其當T和S較小時(文獻[12])，這種波動更為顯著。改進的IWO算法在廣度搜索空間中反演T和S時較GA運行時間更短，所求適應度值更低，并且100次的運算結果更穩定。

4 結 語

改進的IWO算法能在廣度搜索空間中高效模擬配線分析和Theis求解，反演出接近實際的T和S，降低了預估T和S初值范圍的難度。算法編程較易實現，更換算法程序中抽水試驗的數據部分和對算法參數略作調整可處理不同尺度的非穩定流抽水試驗。在理論測試中，改進的IWO算法較GA更快更多地搜尋出滿足一定精度要求的參數。實例應用中，當T和S較大時，改進的IWO算法比常用配線法和直線圖解法求解精度更高，在廣度搜索空間的條件下，改進的IWO算法較GA更優更穩定。基于上述特點，改進的IWO算法可作為在廣度搜索空間中反演T和S的一種新方法。

隨著搜索空間的擴展和對參數求解精度要求的提高，如何進一步改善算法搜索原理，提高算法運行效率需要進一步的研究。如何將其他算法引入本算法并拓展應用需要進一步的嘗試。

□

[1] 齊學斌．非穩定流抽水試驗參數計算的迭代算法及計算機模擬[J]．水利學報，1995，(7)：67-71，66．

[2] 朱春龍．非穩定流抽水試驗參數計算的優化算法[J]．水科學進展，1999，10(1)：75-77．

[3] 霍小虎，黃國如．遺傳算法在水文地質參數確定中的應用[J]．地下水，2001，23(4)：195-197．

[4] 劉立才，陳鴻漢，張達政．梯度法在水文地質參數估值中的應用[J]．水文地質工程地質，2003，(3)：39-41．

[5] 郭建青，李 彥，王洪勝，等．粒子群優化算法在確定含水層參數中的應用[J]．中國農村水利水電，2008,(4)：4-7．

[6] 葉 咸，許 模，廖曉超，等．遺傳算法優化BP神經網絡在求解水文地質參數中的應用[J]．水電能源科學，2013，31(12)：55-57，65．

[7] A R Mehrabian，C．Lucas．A novel numerical optimization algorithm inspired from weed colonization[J]．Ecological Informatics，2006，1(4)：355-366．

[8] 彭 斌，胡常安，趙榮珍．基于混合雜草算法的神經網絡優化策略[J]．振動、測試與診斷，2013，33(4)：634-639，725．

[9] 劉 燕，焦永昌，張亞明，等．混合入侵雜草算法用于陣列天線波束賦形[J]．西安電子科技大學學報(自然科學版)，2014,41(4)：94-99．

[10] 周金虎．基于入侵雜草算法的數據挖掘聚類算法研究[D]．蘭州：蘭州理工大學，2014．

[11] Rajesh Srivastava．Implications of using approximate expressions for well function[J]．Journal of Irrigation and Drainage Engineering．1995，121(6)：459-462．

[12] 薛禹群，吳吉春．地下水動力學[M]．3版. 北京：地質出版社，2010．

[13] 曹萬金．地下水資源計算與評價[M]．北京：水利電力出版社，1987:97-109．

[14] U S Department of the Interior(USDI)．Groundwater Manual[M]．Bureau of Reclamation．Supt．of Documents，U S Govt．Printing Office，Washington，D C，1977：119．