文章详情

短信预约-IT技能 免费直播动态提醒

请输入下面的图形验证码

提交验证

短信预约提醒成功

模拟退火算法(Python)

2023-08-31 05:42

关注

一、模拟退火算法

1、模拟退火算法的定义

模拟退火算法是一种现代优化算法。基于蒙特卡洛迭代求解方法的随机寻优算法,模拟退火算法于1983 年成功地应用到组合优化领域。因固体物理退火过程与组合优化问题存在着相似性,模拟退火算法对固体物质的退火过程进行一定程度的模拟,来获得问题的最优解。

2、模拟退火算法的特点

优点

① 全局搜索能力强,统计上可以保证找到全局最优

缺点

① 找到最优解所耗费的时间较长,尤其是使用标准的Metropolis准则时

3、模拟退火算法的主要步骤

模拟退火算法本质是两层循环,外层循环控制温度由高向低变化;内层循环中,温度固定,对旧解添加随机扰动得到新解,并按一定规则接受新解。

① 解的编码

② 确定初始解

③ 邻域新解的生成

④ 确定能量函数

⑤ Metropolis准则

⑥ 搜索最优解

二、模拟退火算法基本模型

(以求解 x^{3}cosx(-1.57\leq x\leq 20.18) 的最小值为例)

 1、解的编码

由于该解的表示简单,采用实数编码即可。

2、确定初始解

根据定义域范围,随机生成初始解

#[start,end]为定义域def initialization(start, end):    return random.uniform(start, end)

3、邻域新解的生成

① 判断生成的解是否在定义域范围内

def in_range(x, start, end):    return True if start <= x <= end else False

② 生成邻域新解

def generate_new(x, start, end):    while True:        #采用高斯分布生成新解        upper_bound = end - x        lower_bound = start - x        sigma = max(upper_bound, lower_bound) / 3        new_x = random.gauss(x, sigma)        #判断是否在定义域内,在则返回;否则重复生成        if in_range(new_x, start, end):            return new_x

4、确定能量函数

能量越低,被接受的概率越高。因此,当求解最大值目标函数时,可以以目标函数的倒数作为能量函数;当求解最小值目标函数时,可以直接以目标函数作为能量函数

def e(x):    return pow(x, 3) * math.cos(x)

5、新解接受准则

新解接受准则的三条原则

① 在固定温度下,接受使目标函数下降的候选解的概率要大于使目标函数上升的候选解概率;

② 随着温度的下降,接受使目标函数上升的解的概率要逐渐减小;

③ 当温度趋于零时,只能接受目标函数下降的解。

Metropolis准则

新解接受常采用标准的Metropolis准则

在这里插入图片描述

 若\Delta E< 0,接受新解作为当前解;否则,按照概率判断是否接受新解。

def metropolis(e, new_e, t):    if new_e <= e:        return True    else:        p = math.exp((e - new_e) / t)        return True if random.random() < p else False

6、搜索最优解

def search(all_x):    best_e = 0xffff    best_x = all_x[0]    for x in all_x:        ex = e(x)        if ex < best_e:            best_e = ex            best_x = x    return best_x

7、主函数

#t0为初始温度,t_final为终止温度,alpha为冷却系数,inner_iter为内层迭代次数,[start,end]为定义域def sa(t0, t_final, alpha, inner_iter, start, end):    #记录每次内层迭代生成的解    all_x = []    #生成初始解并加入解集    init = initialization(start, end)    all_x.append(init)    #设置初始温度    t = t0    #外层循环    while t > t_final:        x = search(all_x)        ex = e(x)        #内层循环        for i in range(inner_iter):            new_x = generate_new(x, start, end)            new_ex = e(new_x)            if metropolis(ex, new_ex, t):                x = new_x                ex = new_ex        all_x.append(x)        t = alpha * t    return all_x

8、问题求解

all_x = sa(3000, 10, 0.98, 200, -1.57, 20.18)print(round(search(all_x), 3), round(e(search(all_x)), 3))iteration = len(all_x)all_x = np.array(all_x)all_ex = all_x ** 3 * np.cos(all_x)plt.xlabel("Iteration")plt.ylabel("f(x)")plt.plot(range(iteration), all_ex)plt.show()

 输出结果如下

894 -3945.849

 三、模拟退火算法求解TSP模型

以下图为例,求从0号城市出发访问每一座城市并回到0号城市的最短回路

1、基本模型

① 解的编码 

求解TSP模型采用自然数编码,0节点代表初始城市,1-7节点表示其他7座城市。如 1,7,6,3,2,5,4 表示从0出发依次经过城市1、7、6、3、2、5、4最后回到0

② 确定初始解

随机生成一个由数字1-7各出现一次的数列。

③ 邻域新解的生成

采用两点逆转的方法,即对数列设置两个点,然后数列在这两点中发生逆转,改变数列。逆转操作如下:假设数列为A=1,4,5,2,3,6,7,随机生成两个点2和6,那么A1=4,5,2,3,6,逆转后A1=6,3,2,5,4,A=1,6,3,2,5,4,7。

④ 确定能量函数

TSP模型以最短路径为目标函数,因此直接以路径长度为能量函数即可。

⑤ 新解接受准则

采用标准的Metropolis准则。

⑥ 搜索最优解

根据已生成的所有解的能量,搜索能量最低的解,即为模拟退火算法的最优解

2、代码求解

① 确定初始解

#numbers为需要经过的其他城市的数量def initialization(numbers):    path_random = np.random.choice(list(range(1, numbers + 1)), replace=False, size=numbers)    path_random = path_random.tolist()    return path_random

② 邻域新解的生成

def generate_new(path):    numbers = len(path)    #随机生成两个不重复的点    positions = np.random.choice(list(range(numbers)), replace=False, size=2)    lower_position = min(positions[0], positions[1])    upper_position = max(positions[0], positions[1])    #将数列中段逆转    mid_reversed = path[lower_position:upper_position + 1]    mid_reversed.reverse()    #拼接生成新的数列    new_path = path[:lower_position]    new_path.extend(mid_reversed)    new_path.extend(path[upper_position + 1:])    return new_path

 ③ 确定能量函数

#length_mat为各个节点之间的最短距离矩阵def e(path, length_mat):    numbers = len(path) - 1    length = length_mat[0][path[0]] + length_mat[path[-1]][0]    for i in range(numbers):        length += length_mat[path[i]][path[i + 1]]    return length

 ④ Metropolis准则

def metropolis(e, new_e, t):    if new_e <= e:        return True    else:        p = math.exp((e - new_e) / t)        return True if random.random() < p else False

 ⑤ 搜索最优解

def search(all_path, length_mat):    best_e = 0xffff    best_path = all_path[0]    for path in all_path:        ex = e(path, length_mat)        if ex < best_e:            best_e = ex            best_path = path    return best_path

⑥ 模拟退火算法主函数

#t0为初始温度,t_final为终止温度,alpha为冷却系数,inner_iter为内层迭代次数#city_numbers为所需要经过的其他城市,length_mat为各个节点之间的最短距离矩阵def sa(t0, t_final, alpha, inner_iter, city_numbers, length_mat):    all_path = []    all_ex = []    init = initialization(city_numbers)    all_path.append(init)    all_ex.append(e(init, length_mat))    t = t0    while t > t_final:        path = search(all_path, length_mat)        ex = e(path, length_mat)        for i in range(inner_iter):            new_path = generate_new(path)            new_ex = e(new_path, length_mat)            if metropolis(ex, new_ex, t):                path = new_path                ex = new_ex        all_path.append(path)        all_ex.append(ex)        t = alpha * t    return all_path, all_ex

⑦ 主程序

先通过Floyd算法求每个节点的最短路径

graph = nx.Graph()graph.add_nodes_from(range(0, 8))edges = [(0, 1, 3), (0, 4, 2),         (1, 4, 3), (1, 5, 4), (1, 7, 2),         (2, 3, 5), (2, 4, 1), (2, 5, 4), (2, 6, 3),         (3, 6, 1),         (4, 5, 5),         (5, 6, 2), (5, 7, 4),         (6, 7, 4)]graph.add_weighted_edges_from(edges)shortest_length = dict(nx.shortest_path_length(graph, weight='weight'))length_mat = []for i in range(0, 8):    i_j = []    for j in range(0, 8):        i_j.append(shortest_length[i][j])    length_mat.append(i_j)

(图论详见(23条消息) 图论(Python networkx)_Zengwh_02的博客-CSDN博客

再调用模拟退火算法求最优解,并可视化

all_path, all_ex = sa(3000, pow(10, -1), 0.98, 200, 7, length_mat)print(search(all_path, length_mat), round(e(search(all_path, length_mat), length_mat)))iteration = len(all_path)all_path = np.array(all_path)all_ex = np.array(all_ex)plt.xlabel("Iteration")plt.ylabel("Length")plt.plot(range(iteration), all_ex)plt.show()

 3、输出结果

[1, 7, 5, 6, 3, 2, 4] 19

来源地址:https://blog.csdn.net/weixin_58427214/article/details/125901431

阅读原文内容投诉

免责声明:

① 本站未注明“稿件来源”的信息均来自网络整理。其文字、图片和音视频稿件的所属权归原作者所有。本站收集整理出于非商业性的教育和科研之目的,并不意味着本站赞同其观点或证实其内容的真实性。仅作为临时的测试数据,供内部测试之用。本站并未授权任何人以任何方式主动获取本站任何信息。

② 本站未注明“稿件来源”的临时测试数据将在测试完成后最终做删除处理。有问题或投稿请发送至: 邮箱/279061341@qq.com QQ/279061341

软考中级精品资料免费领

  • 历年真题答案解析
  • 备考技巧名师总结
  • 高频考点精准押题
  • 2024年上半年信息系统项目管理师第二批次真题及答案解析(完整版)

    难度     813人已做
    查看
  • 【考后总结】2024年5月26日信息系统项目管理师第2批次考情分析

    难度     354人已做
    查看
  • 【考后总结】2024年5月25日信息系统项目管理师第1批次考情分析

    难度     318人已做
    查看
  • 2024年上半年软考高项第一、二批次真题考点汇总(完整版)

    难度     435人已做
    查看
  • 2024年上半年系统架构设计师考试综合知识真题

    难度     224人已做
    查看

相关文章

发现更多好内容

猜你喜欢

AI推送时光机
位置:首页-资讯-后端开发
咦!没有更多了?去看看其它编程学习网 内容吧
首页课程
资料下载
问答资讯