粒子群算法
PSO——连续
import numpy as np import matplotlib.pyplot as plt from pylab import * matplotlib.use('TkAgg') from pylab import * mpl.rcParams['font.sans-serif'] = ['SimHei'] # 指定默认字体 mpl.rcParams['axes.unicode_minus'] = False # 解决保存图像是负号'-'显示为方块的问题 from tqdm import tqdm # 进度条设置 import random class PSO(object): def __init__(self, M,T,lb,ub,l,n): self.w = 0.6 # 惯性权重 self.c1 = self.c2 = 2 self.M= M # 粒子群数量 self.T = T self.L = l self.N = n self.lb = lb self.ub = ub self.lv = -10 self.uv = 10 def init_x(self): x = np.random.uniform(self.lb, self.ub, (self.M, self.L, self.N)) return x def init_v(self): v = np.random.uniform(0,1,(self.M,self.L,self.N)) return v def fitness(self,x): """ (x1-90)**2+(x2-90)**2 :param x: :return: """ result =(x[0]-90) ** 2 + (x[1]-90) ** 2 return result def best(self,x): l = [] for m in range(self.M): for i in range(self.L): result = self.fitness(list(x[m, i])) l.append(result) return l.index(min(l)),l def updata_v(self,X,V,pbest,gbest): r1 = np.random.random((self.M, self.L, self.N)) # (粒子个数,1) r2 = np.random.random((self.M, self.L, self.N)) V = self.w * V + self.c1 * r1 * (pbest - X) + self.c2 * r2 * (gbest - X) # 直接对照公式写就好了 # 防止越界处理 V[V < self.lv] = np.random.random() * (self.uv - self.lv) + self.lv V[V > self.uv] = np.random.random() * (self.uv - self.lv) + self.lv return V def update_x(self, X, V): """ 根据公式更新粒子的位置 :param X: 粒子当前的位置矩阵，维度是 size*10 :param V: 粒子当前的速度举着，维度是 size*10 return 更新后的种群 """ X1= X.copy() for i in range(self.M): # 遍历每一个粒子 for j in range(self.L): # 遍历粒子中的每一个元素 X1[i,j] = X[i,j]+V[i,j] # 对当前个体进行限制，放在超重 return X1 def update_pbest(self, X, fitness, pbest, pbestfitness): """ 更新个体最优 :param X: 当前种群 :param fitness: 当前每个粒子的适应度 :param pbest: 更新前的个体最优解 :param pbestfitness: 更新前的个体最优适应度 :param m: 粒子数量 :return: 更新后的个体最优解、个体最优适应度 """ for i in range(self.M): if fitness[i] < pbestfitness[i]: pbest[i] = X[i] pbestfitness[i] = fitness[i] return pbest, pbestfitness def update_gbest(self, pbest, pbestfitness, gbest, gbestfitness): """ 更新全局最优解 :param pbest: 粒子群 :param pbestfitness: 个体适应度(个体最优) :param gbest: 全局最优解 :param gbestfitness: 全局最优解 :param m: 粒子数量 :return: gbest全局最优，g对应的全局最优解 """ for i in range(self.M): if pbestfitness[i] < gbestfitness: gbest = pbest[i] gbestfitness = pbestfitness[i] return gbest, gbestfitness def main(self): op_list = [] x = self.init_x() v = self.init_v() num,pbestfitness = self.best(x) gbestfitness = max(pbestfitness) op_list.append(gbestfitness) pbest = x gbest = x[num] for i in range(self.T): v = self.updata_v(x,v,pbest,gbest) x = self.update_x(x,v) num1,fitness = self.best(x) pbest, pbestfitness = self.update_pbest(x,fitness,pbest,pbestfitness) gbest, gbestfitness = self.update_gbest( pbest, pbestfitness, gbest, gbestfitness) op_list.append(gbestfitness) print("最优适应度是：%.5f" % op_list[-1]) print("最优解是", gbest) plt.plot(op_list, label='PSO') plt.xlabel('Iterations', size=13) plt.ylabel('Fitness', size=13) plt.legend() plt.show() if __name__ == '__main__': pso = PSO(30, 1000,-100,100,1,2) # M,T,lb,ub,l,n pso.main()
PSO-离散
import numpy as np import matplotlib.pyplot as plt from tqdm import tqdm # 进度条设置 import matplotlib.pyplot as plt import matplotlib; matplotlib.use('TkAgg') from pylab import * mpl.rcParams['font.sans-serif'] = ['SimHei'] # 指定默认字体 mpl.rcParams['axes.unicode_minus'] = False # 解决保存图像是负号'-'显示为方块的问题 #未修改 class BackPack(object): def __init__(self): self.Weight = [95, 75, 23, 73, 50, 22, 6, 57, 89, 98] # 物品体积 self.Value = [89, 59, 19, 43, 100, 72, 44, 16, 7, 64] # 物品价值 self.N = 50 # 群体粒子个数 self.D = 10 # 粒子维数 self.T = 200# 最大迭代次数 self.c1 = 2 # 学习因子1 self.c2 = 2 # 学习因子2 self.w = 1.5 # 惯性因子，一般取1 self.V_max = 1 # 速度最大值 self.V_min = -1 # 速度最小值 self.Weight_max = 300 # 背包容量 self.afa = 100# 惩罚系数 # 初始化种群 def init_x(self): """ :return: 随机生成的种群(二维list) """ X = np.random.choice([0, 1], size=(self.N, self.D)) return X # 初始化速度 def init_v(self): """ :return: 速度 """ V = np.random.random(size=(self.N, self.D)) ##10维 return V # 适应度值 def fitness_func(self, X): """ :param X: 种群 :return: """ """计算粒子的的适应度值，也就是目标函数值，X 的维度是 size * 10 """ size = np.shape(X)[0] # 种群个数 result = np.zeros((size, 1)) for i in range(size): # 遍历每一个粒子 fit = sum(X[i] * self.Value) # 物品价值 TotalSize = sum(X[i] * self.Weight) # 物品重量 l=list(X[i]) #xiandinglegeshu if TotalSize <= self.Weight_max and l.count(1)<=1: # 如果小于限制重量 fit = fit result[i] = fit else: if TotalSize - self.Weight_max>0: fit = fit - self.afa * (TotalSize - self.Weight_max) # 对fit进行惩罚 else: fit = fit +self.afa * (TotalSize - self.Weight_max) # 对fit进行惩罚 # fit=0 result[i] = fit return result # 我们要求result越大 # 速度更新公式 def velocity_update(self, V, X, pbest, gbest): """ 根据速度更新公式更新每个粒子的速度 种群size=100 :param V: 粒子当前的速度矩阵，size*10 的矩阵 :param X: 粒子当前的位置矩阵，size*10 的矩阵 :param pbest: 每个粒子历史最优位置，size*10 的矩阵 :param gbest: 种群历史最优位置，1*10 的矩阵 """ r1 = np.random.random((self.N, 1)) # (粒子个数,1) r2 = np.random.random((self.N, 1)) V = self.w * V + self.c1 * r1 * (pbest - X) + self.c2 * r2 * (gbest - X) # 直接对照公式写就好了 # 防止越界处理 V[V < self.V_min] = np.random.random() * (self.V_max - self.V_min) + self.V_min V[V > self.V_max] = self.V_max - np.random.random() * (self.V_max - self.V_min) # V[V < self.V_min] =self.V_min # V[V > self.V_max] =self.V_max return V # 位置更新公式 def position_update(self, X, V): """ 根据公式更新粒子的位置 :param X: 粒子当前的位置矩阵，维度是 size*10 :param V: 粒子当前的速度举着，维度是 size*10 return 更新后的种群 """ for i in range(self.N): # 遍历每一个粒子 # 修改速度为sigmoid形式 V[i, :] = 1. / (1 + np.exp(-np.array(V[i, :]))) for j in range(self.D): # 遍历粒子中的每一个元素 rand = np.random.random() # 生成 0-1之间的随机数 if V[i, j] > rand: X[i, j] = 1 else: X[i, j] = 0 # 对当前个体进行限制，放在超重 while np.sum(X[i] * self.Weight) > self.Weight_max: # 如果当前粒子超重 X[i] = np.random.choice([0, 1], size=(1, self.D)) return X def update_pbest(self, X, fitness, pbest, pbestfitness, m): """ 更新个体最优 :param X: 当前种群 :param fitness: 当前每个粒子的适应度 :param pbest: 更新前的个体最优解 :param pbestfitness: 更新前的个体最优适应度 :param m: 粒子数量 :return: 更新后的个体最优解、个体最优适应度 """ for i in range(m): if fitness[i] > pbestfitness[i]: pbest[i] = X[i] pbestfitness[i] = fitness[i] return pbest, pbestfitness def update_gbest(self, pbest, pbestfitness, gbest, gbestfitness, m): """ 更新全局最优解 :param pbest: 粒子群 :param pbestfitness: 个体适应度(个体最优) :param gbest: 全局最优解 :param gbestfitness: 全局最优解 :param m: 粒子数量 :return: gbest全局最优，g对应的全局最优解 """ for i in range(m): if pbestfitness[i] > gbestfitness: gbest = pbest[i] gbestfitness = pbestfitness[i] return gbest, gbestfitness def main(self): fitneess_value_list = [] # 记录每次迭代过程中的种群适应度值变化 x = self.init_x() # 初始化x v = self.init_v() # 初始化v # 计算种群各个粒子的初始适应度值 p_fitness = self.fitness_func(x) # 计算种群的初始最优适应度值 g_fitness = p_fitness.max() # 讲添加到记录中 fitneess_value_list.append(g_fitness) # 初始的个体最优位置和种群最优位置 pbest = x gbest = x[p_fitness.argmax()] # # 接下来就是不断迭代了 for i in tqdm(range(self.T)): pbest = pbest.copy() # 必须加，不然会被篡改值，造成结果错 p_fitness = p_fitness.copy() # 必须加，不然会被篡改值，造成结果错 gbest = gbest.copy() # 必须加，不然会被篡改值，造成结果错 g_fitness = g_fitness.copy() # 必须加，不然会被篡改值，造成结果错 v = self.velocity_update(v, x, pbest, gbest) # 更新速度 x = self.position_update(x, v) # 更新位置 p_fitness2 = self.fitness_func(x) # 计算子代各个粒子的适应度 # 更新每个粒子的历史最优位置 pbest, p_fitness = self.update_pbest(x, p_fitness2, pbest, p_fitness, self.N) # 更新群体的最优位置 gbest, g_fitness = self.update_gbest(pbest, p_fitness, gbest, g_fitness, self.N) # 记录最优迭代结果 fitneess_value_list.append(g_fitness) print("最优适应度是：%.5f" % fitneess_value_list[-1]) print("最优解是", gbest) print('最优解对应的体积', np.sum(gbest * self.Weight)) print('最优解对应的价值', np.sum(gbest * self.Value)) plt.plot(fitneess_value_list, label='迭代曲线') plt.xlabel('迭代次数') plt.ylabel('适应度') plt.legend() plt.show() if __name__ == '__main__': pso = BackPack() pso.main()
