算法改进"/>
模拟退火算法改进
import numpy as np
import matplotlib.pyplot as plt
import math
import random
from scipy.stats import norm
from mpl_toolkits.mplot3d import Axes3D
# 目标函数
def Function(x, y):
return -20 * np.exp(-0.2*np.sqrt(0.5*(x*x+y*y)))\
-np.exp(0.5*(np.cos(2*np.pi*x)+np.cos(2*np.pi*y)))+20+np.e
# 初始化状态
def initN(low, high):
'''
随机生成一组状态取能量最低的状态
:param low:
:param high:
:return:
'''
x = random.uniform(low, high)
y = random.uniform(low, high)
z = Function(x, y)
for i in range(20):
x1 = random.uniform(low, high)
y1 = random.uniform(low, high)
z1 = Function(x1, y1)
if z1 < z:
x = x1
y = y1
z = z1
return x, y
# 绘制目标函数
def figureFuc(low, high):
X = np.linspace(low, high, 500)
Y = np.linspace(low, high, 500)
XX, YY = np.meshgrid(X, Y)
Z = -20 * np.exp(-0.2*np.sqrt(0.5*(XX*XX+YY*YY)))\
-np.exp(0.5*(np.cos(2*np.pi*XX)+np.cos(2*np.pi*YY)))+20+np.e
fig = plt.figure()
ax = Axes3D(fig)
ax.plot_surface(XX, YY, Z, cmap=plt.get_cmap('rainbow'))
plt.show()
# Metropolis准则接受新解
def Metropolis(detaF, T): # detaF为f(n+1) - f(x)
if detaF < 0:
return 1
else:
pTk = math.exp(-detaF/T)
if pTk > random.random():
return 1
else:
return 0
# 利用正态分布产生新解
def get_normal_random_number(x,y,scale): # 正态分布
'''
:param x: 均值
:param y: 均值
:param scale: 方差
:return:
'''
fx = np.random.normal(x, scale)
x = norm.ppf(fx)
fy = np.random.normal(y, scale)
y = norm.ppf(fy)
return x, y
# 利用均匀分布产生新解
def get_uniform_random_number(low, high):
'''
:param low:
:param high:
:return:
'''
x = np.random.uniform(low, high)
y = np.random.uniform(low, high)
return x, y
# 冷却函数
def descT(T, k):
# return T/np.log(1 + k)
return 0.9*T
# 主函数
def startMain():
# 初始化
low = -5
high = 5
T = 10000
Tmin = 10
k = 1
# figureFuc(low, high) # 画图
#x = random.uniform(low, high)
#y = random.uniform(low, high)
x, y = initN(low, high)
z = Function(x, y)
min_value = z
record_value = [] # 用数组记录被接受的新解并绘图,方便分析
while(T > Tmin and k <= 1000):
x1, y1 = get_normal_random_number(x, y, 2) # 利用正态分布产生新解
# x1, y1 = get_uniform_random_number(low, high) # 利用随机分布产生新解
if x1 < low or x1 > high or y1 < low or y1 > high: # 新解不在定义域内时跳出本次循环
break
z1 = Function(x1, y1) # 计算新解的目标函数值
deltaE = z1 - z
min_value = min(min_value, z1)
if Metropolis(deltaE, T) == 1: # 接受按照Metropolis准则接受新解
x = x1
y = y1
z = z1
record_value.append(z)
if deltaE > 0:
T = descT(T, k)
else:
k += 1
print('迭代到组后的解:', z)
print('记录下的最优解:', min_value)
# 打印解的变化曲线
x=[i+1 for i in range(len(record_value))]
plt.plot(x, record_value)
plt.show()
if __name__ == "__main__":
startMain()
更多推荐
模拟退火算法改进
发布评论