代码
# -*- coding: utf-8 -*-
"""
Created on Wed Jul 20 17:14:56 2022
@author: 83989
"""
import math
import random
import numpy as np
import matplotlib.pyplot as plt
import pylab as mpl
mpl.rcParams['font.sans-serif'] = ['SimHei']
class PSO:
def __init__(self, dimension, time, size, low, up, v_low, v_high):
# 初始化
self.dimension = dimension # 变量个数
self.time = time # 迭代的代数
self.size = size # 种群大小
self.bound = [] # 变量的约束范围
self.bound.append(low)
self.bound.append(up)
self.v_low = v_low
self.v_high = v_high
self.x = np.zeros((self.size, self.dimension)) # 所有粒子的位置
self.v = np.zeros((self.size, self.dimension)) # 所有粒子的速度
self.p_best = np.zeros((self.size, self.dimension)) # 每个粒子最优的位置
self.g_best = np.zeros((1, self.dimension))[0] # 全局最优的位置
# 初始化第0代初始全局最优解
temp = -1000000
for i in range(self.size):
for j in range(self.dimension):
self.x[i][j] = random.uniform(self.bound[0][j], self.bound[1][j])
self.v[i][j] = random.uniform(self.v_low, self.v_high)
self.p_best[i] = self.x[i] # 储存最优的个体
fit = self.fitness(self.p_best[i])
# 做出修改
if fit > temp:
self.g_best = self.p_best[i]
temp = fit
def fitness(self, x):
"""
个体适应值计算
"""
x1 = x[0]
x2 = x[1]
x3 = x[2]
x4 = x[3]
x5 = x[4]
y = math.floor((x2 * np.exp(x1) + x3 * np.sin(x2) + x4 + x5) * 100) / 100
# print(y)
return y
def update(self, size):
c1 = 2.0 # 学习因子
c2 = 2.0
w = 0.8 # 自身权重因子
for i in range(size):
# 更新速度(核心公式)
self.v[i] = w * self.v[i] + c1 * random.uniform(0, 1) * (
self.p_best[i] - self.x[i]) + c2 * random.uniform(0, 1) * (self.g_best - self.x[i])
# 速度限制
for j in range(self.dimension):
if self.v[i][j] < self.v_low:
self.v[i][j] = self.v_low
if self.v[i][j] > self.v_high:
self.v[i][j] = self.v_high
# 更新位置
self.x[i] = self.x[i] + self.v[i]
# 位置限制
for j in range(self.dimension):
if self.x[i][j] < self.bound[0][j]:
self.x[i][j] = self.bound[0][j]
if self.x[i][j] > self.bound[1][j]:
self.x[i][j] = self.bound[1][j]
# 更新p_best和g_best
if self.fitness(self.x[i]) > self.fitness(self.p_best[i]):
self.p_best[i] = self.x[i]
if self.fitness(self.x[i]) > self.fitness(self.g_best):
self.g_best = self.x[i]
def pso(self):
best = []
self.final_best = np.array([1, 2, 3, 4, 5])
for gen in range(self.time):
self.update(self.size)
if self.fitness(self.g_best) > self.fitness(self.final_best):
self.final_best = self.g_best.copy()
print('当前最佳位置：{}'.format(self.final_best))
temp = self.fitness(self.final_best)
print('当前的最佳适应度：{}'.format(temp))
best.append(temp)
t = [i for i in range(self.time)]
plt.figure()
plt.grid(axis='both')
plt.plot(t, best, color='red', marker='.', ms=10)
plt.rcParams['axes.unicode_minus'] = False
plt.margins(0)
plt.xlabel(u"迭代次数") # X轴标签
plt.ylabel(u"适应度") # Y轴标签
plt.title(u"迭代过程") # 标题
plt.show()
if __name__ == '__main__':
time = 50
size = 100
dimension = 5
v_low = -1
v_high = 1
low = [1, 1, 1, 1, 1]
up = [25, 25, 25, 25, 25]
pso = PSO(dimension, time, size, low, up, v_low, v_high)
pso.pso()
特点
收敛速度快但容易陷入局部最优解
遗传算法 大概思路代码
from __future__ import division
import numpy as np
import matplotlib.pyplot as plt
def f(x):
return 5*np.sin(6*x) + 6*np.cos(5*x)
class GeneticAlgorithm(object):
"""
遗传算法求多极小函数最优值
"""
def __init__(self, cross_rate, mutation_rate, n_population, n_iterations, DNA_size):
self.cross_rate = cross_rate #交配的可能性大小
self.mutate_rate = mutation_rate #基因突变率
self.n_population = n_population #种群大小
self.n_iterations = n_iterations #迭代次数
self.DNA_size = 6 # DNA的长度
self.x_bounder = [-3, 6] #搜索范围
def init_population(self):
"""
初始化种群
:return:
"""
population = np.random.randint(low=0, high=2, size=(self.n_population, self.DNA_size)).astype(np.int8)
return population
def transformDNA(self, population):
"""
编码：十进制转化为二进制
:param population:
:return:
"""
population_decimal = ( (population.dot(np.power(2, np.arange(self.DNA_size)[::-1])) / np.power(2, self.DNA_size) - 0.5) *
(self.x_bounder[1] - self.x_bounder[0]) + 0.5 * (self.x_bounder[0] + self.x_bounder[1]) )
return population_decimal
def fitness(self, population):
"""
适应性函数
:param population: 基因个体
:return: 与当前最优的解
"""
transform_population = self.transformDNA(population)
fitness_score = f(transform_population)
return fitness_score - fitness_score.min()
def select(self, population, fitness_score):
"""
根据适应值进行选择
:param population:
:param fitness_score:
:return:
"""
fitness_score = fitness_score + 1e-4
fitness_score = np.array([1/i for i in fitness_score]) #求最小值，需要反一下，值越小，越适应
idx = np.random.choice(np.arange(self.n_population), size=self.n_population, replace=True, p=fitness_score/fitness_score.sum())
return population[idx]
def create_child(self, parent, pop):
"""
交叉
:param parent:
:param pop:
:return:
"""
if np.random.rand() < self.cross_rate:
index = np.random.randint(0, self.n_population, size=1)
cross_points = np.random.randint(0, 2, self.DNA_size).astype(np.bool)
parent[cross_points] = pop[index, cross_points]
return parent
def mutate_child(self, child):
"""
变异
:param child:
:return:
"""
for i in range(self.DNA_size):
if np.random.rand() < self.mutate_rate:
child[i] = 1
else:
child[i] = 0
return child
# 进化
def evolution(self):
population = self.init_population()
NUM = []
for i in range(self.n_iterations):
fitness_score = self.fitness(population)
best_person = population[np.argmin(fitness_score)]
if (i+1)%100 == 0:
NUM.append(f(self.transformDNA(best_person)))
print(u'第%-4d次进化后, 最优结果在: %s处取得, 对应的最小值为: %f' % (i, self.transformDNA(
best_person),
f(self.transformDNA(
best_person))))
population = self.select(population, fitness_score)
population_copy = population.copy()
for parent in population:
child = self.create_child(parent, population_copy)
child = self.mutate_child(child)
parent[:] = child
population = population
x_num = [i for i in range(len(NUM))]
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
plt.figure(1)
plt.title('遗传算法')
plt.ylabel('适应度（最小值）')
plt.xlabel('进化次数（*100）')
plt.plot(x_num, NUM, 'g')
plt.show()
self.best_person = best_person
return self.transformDNA(best_person)
def main():
ga = GeneticAlgorithm(cross_rate=0.9, mutation_rate=0.1, n_population=200, n_iterations=2001, DNA_size=8)
best = ga.evolution()
# 绘图
plt.figure(2)
x = np.linspace(start=ga.x_bounder[0], stop=ga.x_bounder[1], num=200)
plt.plot(x, f(x))
plt.scatter(ga.transformDNA(ga.best_person), f(ga.transformDNA(ga.best_person)), s=200, lw=0, c='red', alpha=0.5)
ax = plt.gca()
ax.spines['right'].set_color('none')
ax.spines['top'].set_color('none')
ax.xaxis.set_ticks_position('bottom')
ax.yaxis.set_ticks_position('left')
plt.show()
return best
if __name__ == '__main__':
Bests = []
for _ in range(1):
Bests.append(f(main()))
Avg = np.mean(Bests)
Best = np.min(Bests)
Worst = np.max(Bests)
Var = np.var(Bests)
print("结果:", Bests)
print("nAVG: ", Avg,
"nBest: ", Best,
"nWorst: ", Worst,
"nVar: ", Var)
特点
算法能在较短的时间当中得到最优值，且种群的概念使得算法有并行的模式，多的个体同时搜索最优值。但算法的实现比较麻烦，设计到问题的编码和解码，且参数较多，不同的参数选择选择也没有理论的好坏，只能通过多次尝试。同时由于遗传算法的迭代眼里，使得初始解对最终解的影响比较大。
模拟退火算法 大概思路代码
import matplotlib.pyplot as plt
import numpy as np
import math
import time
np.random.seed(int(time.time()))
num = 300
K = 0.1
R = 0.9 #控制温度降低快慢
T_max = 30 #初始温度
T_min = 0.1 #下限温度
def Func(x):
return 5*np.sin(6*x) + 6*np.cos(5*x)
def main():
x = np.random.uniform(0,2*math.pi)
Best_A =Func(x)
Best_array = []
T_array = []
T = T_max
while(T >T_min):
for i in range(num):
x_temp = np.random.uniform(0,2*math.pi)
current = Func(x_temp)
dE = Best_A - current
if dE>=0:
Best_A = current
x = x_temp
else:
if math.exp(dE/T) > np.random.uniform(0,1):
Best_A = current
x = x_temp
T = R*T
T_array.append(T)
Best_array.append(Best_A)
Plot(Best_array,T_array)
# print("最小值 ： ",Best_A)
return Best_A,x
def Plot(num,T_array):
plt.figure(1)
x_num = [i for i in range(len(num))]
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
plt.subplot(211)
plt.title('模拟退火算法')
plt.ylabel('目标函数')
plt.plot(x_num,num,'g')
plt.subplot(212)
plt.ylabel('温度')
plt.plot(x_num,T_array,'r')
plt.show()
if __name__ == '__main__':
N = 1
Bests = []
xs = []
for _ in range(N):
y,x = main()
Bests.append(y)
xs.append(x)
Avg = np.mean(Bests)
Best = np.min(Bests)
Worst = np.max(Bests)
Var = np.var(Bests)
print("结果:",Bests,"坐标：",xs)
print("nAVG: ",Avg,
"nBest: ",Best,
"nWorst: ",Worst,
"nVar: ",Var)
plt.figure(2)
x = np.linspace(start=0, stop=2*math.pi, num=200)
plt.plot(x, Func(x))
plt.scatter(xs, Bests, s=200, lw=0, c='red', alpha=0.5)
ax = plt.gca()
ax.spines['right'].set_color('none')
ax.spines['top'].set_color('none')
ax.xaxis.set_ticks_position('bottom')
ax.yaxis.set_ticks_position('left')
plt.show()
特点
模拟退火算法计算过程比较简单，且鲁棒性强，多次实验都能达到较好的结果，且结果方差小，统计结果都有较好的表现。同时也发现了模拟退火算法的一些缺点，如果温度下降较快，很容易得到的不是全局最优解。为了得到全局最优解，就得使温度下降的足够慢，这样导致执行时间比较长。
