以Ras函數(shù)(Rastrigin's Function)為目標函數(shù),求其在x1,x2∈[-5,5]上的最小值。這個函數(shù)對模擬退火、進化計算等算法具有很強的欺騙性,因為它有非常多的局部最小值點和局部最大值點,很容易使算法陷入局部最優(yōu),而不能得到全局最優(yōu)解。如下圖所示,該函數(shù)只在(0,0)處存在全局最小值0。
Python代碼實現(xiàn)
import numpy as np
import matplotlib.pyplot as plt
# 目標函數(shù)定義
def ras(x):
y = 20 + x[0] ** 2 + x[1] ** 2 - 10 * (np.cos(2 * np.pi * x[0]) + np.cos(2 * np.pi * x[1]))
return y
# 參數(shù)初始化
w = 1.0
c1 = 1.49445
c2 = 1.49445
maxgen = 200 # 進化次數(shù)
sizepop = 20 # 種群規(guī)模
# 粒子速度和位置的范圍
Vmax = 1
Vmin = -1
popmax = 5
popmin = -5
# 產(chǎn)生初始粒子和速度
pop = 5 * np.random.uniform(-1, 1, (2, sizepop))
v = np.random.uniform(-1, 1, (2, sizepop))
fitness = ras(pop) # 計算適應度
i = np.argmin(fitness) # 找最好的個體
gbest = pop # 記錄個體最優(yōu)位置
zbest = pop[:, i] # 記錄群體最優(yōu)位置
fitnessgbest = fitness # 個體最佳適應度值
fitnesszbest = fitness[i] # 全局最佳適應度值
# 迭代尋優(yōu)
t = 0
record = np.zeros(maxgen)
while t < maxgen:
# 速度更新
v = w * v + c1 * np.random.random() * (gbest - pop) + c2 * np.random.random() * (zbest.reshape(2, 1) - pop)
v[v > Vmax] = Vmax # 限制速度
v[v < Vmin] = Vmin
# 位置更新
pop = pop + 0.5 * v
pop[pop > popmax] = popmax # 限制位置
pop[pop < popmin] = popmin
'''
# 自適應變異
p = np.random.random() # 隨機生成一個0~1內(nèi)的數(shù)
if p > 0.8: # 如果這個數(shù)落在變異概率區(qū)間內(nèi),則進行變異處理
k = np.random.randint(0,2) # 在[0,2)之間隨機選一個整數(shù)
pop[:,k] = np.random.random() # 在選定的位置進行變異
'''
# 計算適應度值
fitness = ras(pop)
# 個體最優(yōu)位置更新
index = fitness < fitnessgbest
fitnessgbest[index] = fitness[index]
gbest[:, index] = pop[:, index]
# 群體最優(yōu)更新
j = np.argmin(fitness)
if fitness[j] < fitnesszbest:
zbest = pop[:, j]
fitnesszbest = fitness[j]
record[t] = fitnesszbest # 記錄群體最優(yōu)位置的變化
t = t + 1
# 結(jié)果分析
print(zbest)
plt.plot(record, 'b-')
plt.xlabel('generation')
plt.ylabel('fitness')
plt.title('fitness curve')
plt.show()
結(jié)果為
[0.99699579 0.00148844]
可以知道求解的點非最小值,算法陷入了局部最小值。
刪除自適應變異部分的注釋,運行后結(jié)果如下,可以看出收斂到全局最優(yōu)解。
[0.00022989 0.00014612]
Matlab有個自帶的粒子群優(yōu)化函數(shù)particleswarm也可以使用。本例的代碼如下:
y = @(x) 20 + x(1).^2 + x(2).^2 - 10*(cos(2*pi*x(1))+cos(2*pi*x(2)));
rng default
options = optimoptions('particleswarm','SwarmSize',200,'HybridFcn',@fmincon,'MaxIterations',200, 'Display','iter');
lb = [-5 -5]; % 這是變量的下限
ub = [5 5]; % 這是變量的上限
[x,fval,exitflag,output] = particleswarm(y,length(lb),lb,ub,options);
結(jié)果如下
-
MATLAB仿真
+關(guān)注
關(guān)注
4文章
175瀏覽量
19876 -
粒子群算法
+關(guān)注
關(guān)注
0文章
63瀏覽量
13001 -
python
+關(guān)注
關(guān)注
54文章
4759瀏覽量
84297
發(fā)布評論請先 登錄
相關(guān)推薦
評論