赞
踩
作者简介:本人擅长运筹优化建模及算法设计,包括各类车辆路径问题、生产车间调度、二三维装箱问题,熟悉CPLEX和gurobi求解器
微信公众号:运筹优化与学习
若有运筹优化建模及算法定制需求,欢迎联系我们私聊沟通
粒子群算法(Particle Swarm Optimization, PSO)的思想源于对鸟群觅食行为的研究,其核心思想是通过群体中个体之间的协作和信息共享来寻找最优解。
相较于遗传算法,粒子群算法具有收敛速度快、参数少、算法简单易实现的优点(对高维度优化问题,比遗传算法更快收敛于最优解),但是同样存在陷入局部最优解的问题。
粒子群算法通过设计粒子来模拟鸟群中的鸟,粒子具有两个属性:速度和位置,其中速度表示粒子下一步迭代时移动的方向和距离,位置代表移动的方向。每个粒子在搜索空间中单独的搜寻最优解,并将其记为当前个体极值,并将个体极值与整个粒子群里的其他粒子共享,找到最优的那个个体极值作为整个粒子群的当前全局最优解,粒子群中的所有粒子根据自己找到的当前个体极值和整个粒子群共享的当前全局最优解来调整自己的速度和位置。
假设在 D D D维搜索空间 S S S中,有 N N N个粒子,每个粒子代表一个解,则第 i i i个粒子的速度、位置表达式如下:
速度 V i = ( v 1 , v 2 , ⋯ , v D ) V_{i}=(v_{1},v_{2},\cdots,v_{D}) Vi=(v1,v2,⋯,vD)
位置 X i = ( x 1 , x 2 , ⋯ , x D ) X_{i}=(x_{1},x_{2},\cdots,x_{D}) Xi=(x1,x2,⋯,xD)
在每一次的迭代中,粒子通过跟踪两个“极值”(pbest,gbest)来更新自己的速度和位置,更新表达式如下:
V i = ω × V i + c 1 × r a n d ( ) × ( p b e s t i − X i ) + c 2 × r a n d ( ) × ( g b e s t i − X i ) V_{i}=\omega\times V_{i}+c_{1}\times rand()\times(pbest_{i}−X_{i})+c_{2}\times rand()\times(gbest_{i}−X_{i}) Vi=ω×Vi+c1×rand()×(pbesti−Xi)+c2×rand()×(gbesti−Xi)
X i = X i + V i X_{i} = X_{i} + V_{i} Xi=Xi+Vi
式中, ω \omega ω表示惯性因子,现有研究表明采用动态惯性因子能够获得更好的寻优结果,目前最常用的为线性递减权值(Linear Decreasing Weight,LDW)策略。 c 1 c_1 c1为个体学习因子, c 2 c_2 c2为群体学习因子。pbest表示当前粒子搜索最好解,gbest表示群体中搜索最好解。
第一部分 V i V_{i} Vi和 X i X_{i} Xi称为惯性部分,表示上次迭代中粒子的速度、位置的影响;
第二部分 c 1 × r a n d ( ) × ( p b e s t i − X i ) c_{1}\times rand()\times(pbest_{i}−X_{i}) c1×rand()×(pbesti−Xi)称为自身认知部分,可理解为粒子当前位置与自身历史最优位置之间的距离和方向,表示粒子的动作来源于自己经验的部分;
第三部分 c 2 × r a n d ( ) × ( g b e s t i − X i ) c_{2}\times rand()\times(gbest_{i}−X_{i}) c2×rand()×(gbesti−Xi)称为群体认知部分,可理解为粒子当前位置与群体历史最优位置之间的距离和方向,反映了粒子间的协同合作和知识共享。
粒子群优化算法的步骤如下:
class particle_single(): ''' position -> 粒子位置 velocity -> 粒子速度 best_p -> 全局最优解 obj_function -> 目标函数 constraints -> 约束 obj_value -> 目标函数值 ''' def __init__(self,obj_func,attribute_number,constr=[],vmax=np.array(np.nan),l_bound=np.array(np.nan),u_bound=np.array(np.nan),integer=np.array(np.nan),position=np.array(np.nan),velocity=np.array(np.nan)): self.obj_function = obj_func if type(constr)!=list: constr = [constr] self.constraints = constr self.att = attribute_number if np.all(np.isnan(position))==False: self.position = position else: try: if attribute_number!=1: init_pos = [] for i in range(self.att): if integer[i]==False: init_pos.append(random.uniform(l_bound[i],u_bound[i])) else: init_pos.append(random.randint(l_bound[i],u_bound[i])) self.position = np.array(init_pos) else: if integer==False: self.position= random.uniform(l_bound,u_bound) else: self.position= random.randint(l_bound,u_bound) except: print('We need lower and upper bound for init position') self.obj_value = self.calc_obj_value() if np.all(np.isnan(velocity))==False: self.velocity = velocity else: try: if attribute_number!=1: self.velocity = np.array([random.uniform(-vmax[i],vmax[i]) for i in range(self.att)]) else: self.velocity = random.uniform(-vmax,vmax) except: print('we need an vmax for init velocity') self.best_p=np.nan
class particle_multi(particle_single): ''' position -> 粒子位置 velocity -> 粒子速度 best_p -> best particle in the past gets just set by first_set_best obj_functions -> 目标函数 constraints -> 约束条件 obj_values -> 目标函数值 S -> 当前解支配解集 n -> 支配解数量 distance -> 拥挤距离 rank -> 支配排序 ''' def __init__(self,obj_func,attribute_number,constr=[],vmax=np.array(np.nan),l_bound=np.array(np.nan),u_bound=np.array(np.nan),integer=np.array(np.nan),position=np.array(np.nan),velocity=np.array(np.nan)): self.obj_functions = obj_func self.constraints=constr self.att = attribute_number if np.all(np.isnan(position))==False: self.position=position else: try: if attribute_number!=1: init_pos = [] for i in range(self.att): if integer[i]==False: init_pos.append(random.uniform(l_bound[i],u_bound[i])) else: init_pos.append(random.randint(l_bound[i],u_bound[i])) self.position = np.array(init_pos) else: if integer==False: self.position= random.uniform(l_bound,u_bound) else: self.position= random.randint(l_bound,u_bound) except: print('We need lower and upper bound for init position') self.obj_values =self.calc_obj_value() if np.all(np.isnan(velocity))==False: self.velocity = velocity else: try: if attribute_number!=1: self.velocity = np.array([random.uniform(-vmax[i],vmax[i]) for i in range(self.att)]) else: self.velocity = random.uniform(-vmax,vmax) except: print('we need an vmax for init velocity') self.best_p=np.nan self.S = [] self.n = np.nan self.rank = np.nan self.distance = np.nan
def runner(self,Iterations, time_termination): t0 = time.time() for i in range(Iterations): # 迭代次数上限 if time_termination != -1 and time.time()-t0 > time_termination: # 求解时间上限 break if self.multi: self.comp_swarm=[] self.comp_swarm.append(self.g_best) for part in self.swarm: # 更新速度 r1 = random.random() r2 = random.random() new_v = self.v_weight*part.velocity + self.c_param*r1*(part.best_p.position-part.position) + self.s_param*r2*(self.g_best.position-part.position) # control vmax new_v = np.array([new_v[i] if new_v[i]>-self.vmax[i] else -self.vmax[i] for i in range(len(new_v))]) new_v = np.array([new_v[i] if new_v[i]<self.vmax[i] else self.vmax[i] for i in range(len(new_v))]) # 更新位置 new_p = part.position + new_v for i in range(len(new_p)): if self.integer[i]: new_p[i] = int(new_p[i]) # 检验边界条件 new_p = np.array([new_p[i] if new_p[i]>self.l_bound[i] else self.u_bound[i] - abs(self.l_bound[i]-new_p[i])%(self.u_bound[i]-self.l_bound[i]) for i in range(len(new_p))]) new_p = np.array([new_p[i] if new_p[i]<self.u_bound[i] else self.l_bound[i] + abs(self.u_bound[i]-new_p[i])%(self.u_bound[i]-self.l_bound[i]) for i in range(len(new_p))]) part.set_velocity(new_v) part.set_position(new_p) if self.multi: self.comp_swarm.append(part) self.comp_swarm.append(part.best_p) if self.multi: self.non_dom_sort() # 非支配排序 self.g_best = copy.deepcopy(self.comp_swarm[0]) j=1 new_swarm = self.comp_swarm[1:-1:2] self.swarm = copy.deepcopy(new_swarm) j+=1 for part in self.swarm: if self.multi: part.best_p = copy.deepcopy(self.comp_swarm[j]) j+=2 if part.compare(self.g_best): self.g_best = copy.deepcopy(part) # 更新当前最优解 part.compare_p_best()
若有运筹优化建模及算法定制需求,欢迎联系我们私聊沟通
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。