求解非線性方程組的HHHO算法及工程應(yīng)用

洪麗啦,莫愿斌*,2,鮑冬雪

(1. 廣西民族大學人工智能學院,廣西 南寧530006;2. 廣西民族大學廣西混雜計算與集成電路設(shè)計分析重點實驗室,廣西 南寧530006)

1 引言

現(xiàn)實生活中很多問題的求解,最終都可以歸結(jié)為求解復雜非線性方程組,由于非線性方程組的求解具有一定的復雜性,如何有效的獲得方程組的精確解及全部解成為了一個研究目標,在群智能優(yōu)化算法未被廣泛地應(yīng)用于求解優(yōu)化問題前,線性方程組的求解大多依賴于傳統(tǒng)的求解方法,例如,牛頓法、迭代法、最速下降法、單純型法等,但是這些方法對目標函數(shù)的初始值具有很強的依賴性和敏感性。與傳統(tǒng)的求解方法相比,群智能優(yōu)化算法能夠在全局范圍內(nèi)對最優(yōu)值進行自適應(yīng)、自動地向個體最優(yōu)值學習追蹤目標解,不依賴于初始值,且對目標函數(shù)是否可導不做限制。因此,用群智能優(yōu)化算法求解非線性方程組成為了一個研究方向,目前已經(jīng)有很多學者把智能優(yōu)化算法應(yīng)用于非線性方程組的求解,如文獻[1]的和聲搜索算法,文獻[2]的蝴蝶優(yōu)化算法、文獻[3]的混合布谷鳥算法、文獻[4]的差分進化算法,文獻[5]的BFGS差分進化算法等,雖然智能優(yōu)化算法計算速度快,但是存在求解精度低,求解根的個數(shù)不完全,收斂速度慢等問題,針對以上存在的問題,本文提出一種改進的混合哈里斯鷹優(yōu)化算法用于求解非線性方程組,首先,在全局階段,采用二次插值這種二次擬合的方式,增強了算法的全局尋優(yōu)能力,加快了搜索速度,其次,當算法迭代到后期時容易陷入局部最優(yōu),根據(jù)差分進化算法具有良好的全局探索能力,將差分進化中的變異、選擇操作引入到局部搜索階段,使算法增強擺脫陷入局部最優(yōu)的能力,最后,通過把HHHO算法用于解決幾何約束問題,結(jié)果說明HHHO算法在求解精度具有一定的優(yōu)勢。

2 非線性方程組

非線性方程組的數(shù)學模型如下

其中,m代表變量的個數(shù),n代表求解方程的個數(shù),可以將方程組的求解轉(zhuǎn)換成目標函數(shù),形式如下:

求解A(X)的根等價于求解非線性方程組的各個解。

3 哈里斯鷹優(yōu)化算法

哈里斯鷹優(yōu)化算法[6](Harris hawks optimization,HHO)是于2019年Heidari 等人通過模擬哈里斯鷹的捕食行為而提出的一種啟發(fā)式智能優(yōu)化算法,哈里斯鷹圍攻獵物可分為兩個階段,勘探階段和開發(fā)階段,勘探階段主要是對獵物進行等待和監(jiān)視,開發(fā)階段對獵物進行抓捕,在開發(fā)階段根據(jù)不同策略對獵物進行圍攻,每個階段詳細的過程如下。

3.1 勘探階段

在勘探階段,哈里斯鷹主要通過兩種概率相等的策略對獵物進行監(jiān)視。

(1)

(2)

其中,Xrand(t)代表當前種群中的一個隨機位置,Xrabbit(t)代表獵物的位置,即當前最優(yōu)解,Xm(t)代表當前種群位置的平均值,r1,r2,r3,r4,q均為0到1的隨機數(shù),UB和LB分別為種群的上下界,N為種群數(shù)量。

3.2 勘探階段向開發(fā)階段的轉(zhuǎn)換

勘探階段轉(zhuǎn)成開發(fā)階段是通過獵物的逃跑能E來控制,當能量的絕對值大于等于1,則進行全局搜索,當能量的絕對值小于等于1,則進行局部搜索。

(3)

其中,E0為獵物的初始狀態(tài),E0∈(-1,1),t為當前迭代次數(shù),T為最大迭代次數(shù)。

3.3 開發(fā)階段

1)軟圍困

當r≥0.5且|E|≥0.5時,此時獵物有足夠的能量并且想通過隨機的跳躍逃跑,哈里斯鷹將其輕輕包圍,進行突襲,其數(shù)學公式如下

X(t+1)=ΔX(t)-E|JXrabbit(t)-X(t)|

(4)

ΔX(t)=Xrabbit(t)-X(t)

(5)

J=2(1-r5)

(6)

其中,ΔX(t)是獵物的位置向量與當前迭代t中的個體向量的差值,J表示獵物在逃離哈里斯鷹圍困的過程中的跳躍程度,r5∈(0,1)。

2)硬圍困

當r≥0.5且|E|<0.5時,此時獵物已經(jīng)非常疲憊,沒有足夠的能量逃跑,哈里斯鷹將其進行硬圍困,其數(shù)學公式如下

X(t+1)=Xrabbit(t)-E|ΔX(t)|

(7)

3)軟性包圍和漸進式快速俯沖

當r<0.5且|E|≥0.5時,獵物有足夠的能量逃脫包圍圈,此時哈里斯鷹根據(jù)獵物(尤其是兔子)的欺騙性動作變換抓捕的位置和方向,對獵物形成一個包圍圈并進行快速俯沖圍困,其數(shù)學公式如下

Y=Xrabbit(t)-E|JXrabbit(t)-X(t)|

(8)

Z=Y+S×LF(D)

(9)

(10)

(11)

(12)

其中,D為種群的維度,S為1×D的隨機向量,LF為萊維飛行函數(shù),μ和v是0到1的隨機數(shù),β取值為1.5。

4)使用漸進的快速俯沖進行猛烈的圍攻

當r<0.5且|E|<0.5時,獵物沒有足夠的能量逃跑,此時老鷹通過縮小與獵物的平均位置距離,形成一個很小的包圍圈來抓捕獵物,其數(shù)學公式如下

Y=Xrabbit(t)-E|JXrabbit(t)-Xm(t)|

(13)

Z=Y+S×LF(D)

(14)

(15)

4 混合哈里斯鷹優(yōu)化算法

4.1 引入二次插值

在勘探階段,由于哈里斯鷹的初始位置是在解空間中隨機產(chǎn)生的,且哈里斯鷹尋找獵物是通過從解空間中隨機的同伴來獲取信息,對目標的信息量獲取較少,因此不容易找到目標,導致尋優(yōu)能力降低,為此在勘探階段引入二次插值算子[7]式(16),二次插值是一種曲線擬合方法,通過在解空間中構(gòu)造低次多項式的解來不斷的替換逼近原目標函數(shù),以此來更快的求得目標函數(shù)的解。

在每一次迭代的勘探階段,在種群中隨機的選擇兩個哈里斯鷹的位置作為二次插值的區(qū)間,根據(jù)式(16)求出二次曲面的極小值作為當前最優(yōu)值,把得到的最優(yōu)值和種群的最優(yōu)值進行比較,選擇較好的作為當前迭代的最優(yōu)個體,通過這種策略,可以使種群更快的向最優(yōu)個體靠攏,以此來提高算法的全局搜索能力。

(16)

(17)

4.2 結(jié)合差分變異機制

經(jīng)過追逐奔跑,獵物的逃跑能量漸漸衰弱,處于被圍困狀態(tài),此時哈里斯鷹聚集在獵物(最優(yōu)個體)周圍,但如果當前獵物的位置為局部極值,那么哈里斯鷹就容易陷入局部最優(yōu),即出現(xiàn)“早熟”現(xiàn)象,不利于種群向更優(yōu)的個體進化,導致算法在解空間的求解精度降低。因此需要采取策略來改變哈里斯鷹的過度聚集,即需要增強種群的多樣性,使部分哈里斯鷹朝著解空間其它地方搜索,以此來增強哈里斯鷹掙脫陷入局部最優(yōu)的能力,從而提高算法的求解精度。

根據(jù)早熟機制[8],當種群滿足早熟條件,說明種群多樣性變差,此時需要采取一定的手段來改善種群的多樣性。差分進化[9]是一種具有較強的全局搜索能力的智能優(yōu)化算法,其主要是依靠種群之間的競爭與合作行為來完成信息的共享與交流,首先,對父代種群進行變異、交叉操作得到子代,其次,根據(jù)貪婪選擇機制,選擇較優(yōu)的個體作為當代的新一代群體。差分進化有五種變異策略,采用了DE/best/1,具體的形式如下

DE/best/1:Vi(t)=Xbest(t)+F(Xr1(t)-Xr2(t))

(18)

其中,Xbest(t)為本次迭代的最優(yōu)個體,Xr1(t) 和Xr2(t) 是當前迭代種群中的不同的個體,F是變異系數(shù),其值設(shè)為0.4,通過該變異式子加強種群間的信息共享,通過個體之間的交叉組合,可以增強哈里斯鷹種群的多樣性,使陷入早熟的哈里斯鷹有機會朝著其它解空間繼續(xù)搜索其它潛在解,然后在再根據(jù)選擇操作保留下優(yōu)質(zhì)個體,選擇操作如下

(19)

4.3 HHHO算法步驟

1)初始化相應(yīng)的參數(shù)并給種群賦初值;

2)根據(jù)式(3)計算出獵物的逃跑能量E;

3)計算種群的適應(yīng)度函數(shù)、均值,找出種群的最優(yōu)個體;

4)若|E|≥1,則采用式(1)-(2)對種群進行全局搜索。

5)若|E|<1,根據(jù)以下條件對種群進行局部開發(fā),并用式(16)-式(17)進行二次插值。

a)若r≥0.5且|E|>0.5時,采用式(4)-(6)對種群進行局部開發(fā);

b)若r≥0.5且|E|<0.5時,采用式(7)對種群進行局部開發(fā);

c)若r<0.5且|E|≥0.5時,采用式(8)-(12)對種群進行局部開發(fā);

d)若r<0.5且|E|<0.5時,采用式(13)-(15)對種群進行局部開發(fā);

6)根據(jù)早熟機制,若算法陷入局部最優(yōu),則采用式(18)-式(19)對進行變異、選擇操作。

7)判斷是否達到最大迭代次數(shù),若沒有達到最大迭代次數(shù),則反復執(zhí)行3)-6),直到滿足最大迭代條件為止。

4.4 算法偽碼

Begin

Initialize algorithm parameters and generate initialization population

X={Xij,i=1,2…N,j=1,2…N}

Calculate escape energyE

While t

for i=1 to N do

The fitness value of population is calculated to find out the optimal individual

if(|E|>1)do

Update the position according to equations (1)and (2)

else if (|E|<1)do

if(r≥0.5&&|E|>0.5)do

Update the position according to equation (4)

if(r≥0.5&&|E|<0.5)do

Update the position according to equation (7)

if (r≥0.5&&|E|≥0.5)do

Update the position according to equations (8)-(12)

if(r<0.5&&|E|<0.5)do

Update the position according to equations (13)-(15)

end if

end for

t=t+1;

end while

end

5 仿真與結(jié)果分析

5.1 測試函數(shù)

為了測試HHHO算法的性能,本文選取了10個基本測試函數(shù)對算法進行測試,測試函數(shù)如表1所示。

表1 標準函數(shù)測試

5.2 測試結(jié)果

為了驗證混合哈里斯鷹優(yōu)化算法性能的有效性和可行性,本文將HHHO算法和HHO,WOA,GWO,ALO四種算法同時對同一個測試函數(shù)進行測試分析,種群規(guī)模設(shè)置為30,最大迭代次數(shù)為300,為了驗證算法的穩(wěn)定性和魯棒性,將算法獨立運行30次,計算出每種算法測試的最優(yōu)值、最差值、平均值,實驗數(shù)據(jù)如表2-表3所示,為了更直觀的觀察算法的收斂效果,畫出了效果收斂圖如圖1-圖10所示。

圖1 函數(shù)f1的收斂圖

圖2 函數(shù)f2的收斂圖

圖4 函數(shù)f4的收斂圖

圖5 函數(shù)f5的收斂圖

圖6 函數(shù)f6的收斂圖

圖7 函數(shù)f7的收斂圖

圖8 函數(shù)f8的收斂圖

圖9 函數(shù)f1的收斂圖

圖10 函數(shù)f10的收斂圖

表2 實驗數(shù)據(jù)(1)

表3 實驗數(shù)據(jù)(2)

5.3 測試結(jié)果分析

本文將HHHO算法和WOA、GWO、ALO、HHO四種算法進行比較,首先,由表2-表3可知,HHHO算法對f1,f3,f4,f7,f8,f10能搜索到理論最優(yōu)值0,且平均值和標準差都是0,可以看出HHHO具有較高的穩(wěn)定性,而對于函數(shù)f5,f6,f9雖然在300次迭代中不能搜尋到最優(yōu)值,但是和其它四個算法相比而言,求解精度是最高的,而對f2函數(shù)HHHO和HHO具有相同的求解精度。其次,從圖1、圖2、圖7、圖8的收斂曲線可以直觀的看出HHHO在前期的收斂速度較快,證明了二次插值提高了HHO在高維空間的搜索能力。

6 算法用于求解非線性方程組

采用5個非線性方程組進行測試,種群規(guī)模設(shè)置為30,迭代次數(shù)設(shè)置為800次,每個方程組獨立運行10次取平均值,實驗結(jié)果如表4-表8所示。

表4 方程組1結(jié)果比較

表5 方程組2結(jié)果比較

表6 方程組3結(jié)果比較

6.1 非線性方程組

其中,-1≤x1,x2≤1,理論值x*=(1,1,1),x*=(-1,-1,-1)

其中,-1≤x1,x2≤1,理論值x*=(0,1)

其中,-1.732≤x1,x2,x3≤1.732,理論值x*=(1,1,1)

-100≤x1,x2≤100,理論值為x*=(10,1)和x*=(-7.5,1)

其中,-2≤x1,x2≤2,該方程有10個解,分別是:(-1.81089,-0.34909),(-1.81089,0.34909),(-1.50222,-0.40908),(-1.50222,0.40908),(-1.79130,0.30193),(-1.791330,-0.30193),(-0.94724,0.78502),(-0.94724,-0.78502),(-0.21306,-1.25685)．

6.2 測試結(jié)果

6.3 測試結(jié)果分析

由表4-6可知,HHHO算法相比于文獻[10]在種群規(guī)模為100時具有更高的求解精度,且能夠搜索到理論最優(yōu)值,由表4可以看出,文獻[11]只能搜索到一個解,而HHHO能夠求到兩個精確解,由表7可知,HHHO算法能夠獲得方程組的理論解,但是方程組的值比文獻[11]略差了點,由表8可知,HHHO算法能夠求得方程的10個根,并且與理論值接近,雖然文獻[11]也能求到這個值,但是種群規(guī)模卻比HHHO的種群規(guī)模大得多。綜上,HHHO算法在求解精度及求解速度方面都具有一定的優(yōu)勢。

表7 方程組4結(jié)果比較

表8 方程組5結(jié)果

7 算法在工程上的應(yīng)用

7.1 幾何約束問題

圖11 幾何約束模型

x-y-100=0

(20)

(21)

通過仿真,求得方程組(21)的解如表9所示。

表9 仿真求得方程組(21)的解

表10 文獻[14]和HHHO兩算法求解比較

7.2 三角函數(shù)超越方程

三角函數(shù)超越方程[14]的零點求解對許多機械問題具有重要的求解意義,由于其龐大的數(shù)值解成為了研究的一個熱門和難點,本文采用HHHO算法來求解三角函數(shù)超越方程,其表達式如式(22)所示,其中,x1∈[-3π/2,3π/2],x2∈[-π,2π]。

(22)

通過把HHHO算法應(yīng)用在幾何約束問題和三角函數(shù)超越方程中,從實驗數(shù)據(jù)可以看出,相比于文獻[13],HHHO算法求解出來的值更接近于理論值,相比于文獻[14]采用的同倫法求解三角函數(shù)超越方程,HHHO算法求出了18個近似的解,并且解能精確到10位數(shù),相對于文獻[14],HHHO算法的求解精度更高,因此可以看出HHHO算法具有較高的求解性能,對工程應(yīng)用具有一定的實際意義。

8 結(jié)語

針對HHO算法在前期搜索速度慢,且由于在進化后期隨著種群多樣性減少,導致算法容易陷入局部最優(yōu)等缺陷,提出一種基于二次插值和差分進化優(yōu)化算法的HHHO算法,通過10個測試函數(shù)和5個非線性方程組的測試結(jié)果比較分析,驗證了HHHO算法在求解精度、避免算法陷于局部最優(yōu)、收斂速度快等方面都具有良好的表現(xiàn),證明了算法的全局探索能力和局部開發(fā)能力得到了改進。最后,將HHHO算法應(yīng)用于求解幾何約束問題和三角函數(shù)超越方程,并且得到了滿意的求解結(jié)果。