优化算法的python实现——粒子群算法（PSO）

用处：用于优化连续非线性方程
特点：收敛速度快，参数少，算法简单
问题：可能陷入局部最优

我认为在中文互联网上找到的对于粒子群算法讲解得最详细的
粒子群优化算法(Particle Swarm Optimization, PSO)的详细解读

根据算法的定义，自己写出来的 Demo 如下, 优化的例子是
$$\underset{x,y}{\max }-(x-1{)}^2-(y-2{)}^2$$

2023.6.13 版本（eval()函数 & 动态变量定义与调用）

把函数外置了！这样就能不用在类里改方程了。主要用到的就是 eval()函数 和 动态变量名定义与调用。虽然这样，但是还不是很清楚为什么用 exec 的时候会报错orz

#import sympy
import numpy as np
class PSODemo:
    def __init__(self,N,dim,weight,learn,itermax,fun):
        self.num = N #   num of initial paricles
        self.dim = dim#   dim of variable
        self.minweight = weight[0]
        self.maxweight = weight[1]
        self.person = learn[0] #personal_learn
        self.group = learn[1] #group_learn
        self.itermax = itermax
        self.fun = fun
        
#    objective function
    def Function_Value(self,state):
        names = locals()
        for i in range(self.dim):
            names['x'+str(i+1)] = state[...,i]
        value = eval(self.fun)
        return value
    
#     intertia weight w
    def In_weight(self,iter):
        w = self.maxweight-(self.maxweight-self.minweight)*iter/self.itermax
        return w
        
#     update velocity and state
    def Update_VeloState(self,w,old_velo,old_state,old_pbest,old_gbest):
        r = np.random.rand(1,2)
        velo = w*old_velo+self.person*r[0,0]*(old_pbest-old_state)+self.group*r[0,1]*(np.ones((self.num,2))*old_gbest-old_state)
        state = old_state+velo
        return velo,state
    
#     update pbest and gbest
    def Update_Best(self,old_state,state,pbest,gbest,value_pbest,value_gbest):
        value = self.Function_Value(state)
        for i in range(0,self.num-1):
            if value[i] > value_pbest[i]:
                pbest[i,...] = state[i,...]
                value_pbest[i] = value[i]
        if max(value) > value_gbest:
            value_gbest = max(value)
            value_li = value.tolist()
            index = value_li.index(max(value_li))
            gbest = state[index,...]
            value_gbest = value[index]
        return pbest,gbest,value_pbest,value_gbest
    
#       loop
    def Loop(self):
        new_state = np.random.rand(self.num,self.dim)
        new_pbest = new_state
        new_value_pbest = self.Function_Value(new_state)
        value_li = new_value_pbest.tolist()
        index = value_li.index(max(value_li))
        new_gbest = new_state[index,...]
        new_value_gbest = new_value_pbest[index]
        new_velo = np.ones((self.num,self.dim))
        for idx in range(1,self.itermax):
            old_state = new_state
            old_pbest = new_pbest
            old_value_pbest = new_value_pbest
            old_gbest = new_gbest
            old_value_gbest = new_value_gbest
            old_velo = new_velo
            w = self.In_weight(idx)
            new_velo,new_state =  self.Update_VeloState(w,old_velo,old_state,old_pbest,old_gbest)
            new_pbest,new_gbest,new_value_pbest,new_value_gbest = self.Update_Best(old_state,new_state,old_pbest,old_gbest,old_value_pbest,old_value_gbest)
        return new_gbest,new_value_gbest 


if __name__ == '__main__':
    num = 40
    dim = 2
    weight = [0.9,0.1]# [maxweight,minweight]
    learn = [1.8,2]
    itermax = 1000
    fun = '-(x1-1)**2-(x2-2)**2' 
    # objective function [maximize], the [nth] variable is defined as [xn]
    PSO = PSODemo(num,dim,weight,learn,itermax,fun)
    
    new_gbest,new_value_gbest = PSO.Loop()
    print('the optimal solution = ',new_gbest)
    print('the optimal value = ',new_value_gbest)

2023.6.11 版本

很简单的凸问题。代码有个问题，优化目标函数需要每次都在 class 里修改，还在思考代码在这方面的优化（自己对 python 实在不是很熟悉，如果有大佬经过，指点一下就好了（跪

import sympy
import numpy as np
class PSODemo:
    def __init__(self,N,dim,weight,learn):
        self.num = N # num of initial paticles
        self.dim = dim # dim of variables
        self.minweight = weight[0]
        self.maxweight = weight[1]
        self.person = learn[0] #personal_learn
        self.group = learn[1] #group_learn
        
#    objective function
    def Function_Value(self,state):
#    max -(x-1)^2-(y-2)^2
#        x,y = sympy.symbols('x y')
        value = -(state[...,0]-1)**2-(state[...,1]-2)**2
        return value
    
#     intertia weight w
    def In_weight(self,iter,intermax):
        w = self.maxweight-(self.maxweight-self.minweight)*iter/itermax
        return w
        
#     update velocity and state
    def Update_VeloState(self,w,old_velo,old_state,old_pbest,old_gbest):
        r = np.random.rand(1,2)
        velo = w*old_velo+self.person*r[0,0]*(old_pbest-old_state)+self.group*r[0,1]*(np.ones((self.num,2))*old_gbest-old_state)
        state = old_state+velo
        return velo,state
    
#     update pbest and gbest
    def Update_Best(self,old_state,state,pbest,gbest,value_pbest,value_gbest):
        value = self.Function_Value(state)
        for i in range(0,self.num-1):
            if value[i] > value_pbest[i]:
                pbest[i,...] = state[i,...]
                value_pbest[i] = value[i]
        if max(value) > value_gbest:
            value_gbest = max(value)
            value_li = value.tolist()
            index = value_li.index(max(value_li))
            gbest = state[index,...]
            value_gbest = value[index]
        return pbest,gbest,value_pbest,value_gbest

停止条件是设置的最大迭代次数

num = 40
dim = 2
weight = [0.9,0.1]
learn = [1.6,2]
PSO = PSODemo(num,dim,weight,learn)

if __name__ == '__main__':
    itermax = 1000
    # num of interation
    new_state = np.random.rand(num,dim)
    new_pbest = new_state
    new_value_pbest = PSO.Function_Value(new_state)
    value_li = new_value_pbest.tolist()
    index = value_li.index(max(value_li))
    new_gbest = new_state[index,...]
    new_value_gbest = new_value_pbest[index]
    new_velo = np.ones((num,dim))
    for idx in range(1,itermax):
        old_state = new_state
        old_pbest = new_pbest
        old_value_pbest = new_value_pbest
        old_gbest = new_gbest
        old_value_gbest = new_value_gbest
        old_velo = new_velo
        w = PSO.In_weight(idx,itermax)
        new_velo,new_state =  PSO.Update_VeloState(w,old_velo,old_state,old_pbest,old_gbest)
        new_pbest,new_gbest,new_value_pbest,new_value_gbest = PSO.Update_Best(old_state,new_state,old_pbest,old_gbest,old_value_pbest,old_value_gbest)
    print(new_gbest,new_value_gbest)

输出的结果为
[1.00000008 1.99999998] -6.9334143407325055e-15
试过几次之后感觉跟参数的选择有很大的关系，也不懂是不是我写得代码问题。希望有大佬能指导指导。