粒子群算法

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()