1.什么是模拟退火算法(简介)

模拟退火算法是一种基于概率的全局寻优方法,用来在一个大的搜寻空间内找寻命题的最优解。它是基于Monte-Carlo迭代求解策略的一种随机寻优算法。

2.模拟退火算法核心思想

模拟退火算法的(Simulated Annealing,SA)是一种基于概率的全局寻优方法,已在理论上被证明以概率1 收敛于全局最优解。模拟退火算法模拟物理退火过程,从某一较高初温出发,随着温度的不断下降,以一定概率突跳在全局进行寻优,并最终趋于全局最优,搜索过程中趋于零概率的突跳特性可有效避免算法陷入局部最优。模拟退火算法依赖现有求解规则,是一种对已有规则进行改造的算法,它的解与初始值无关;其核心思想是以1概率接受较优解,以较小概率接受劣质解(Metropolis准则)

3.模拟退火的热力学原理

根据热力学的原理,在温度为T时,出现能量差为dE的降温的概率为P(dE),表示为:

P(dE) = exp( dE/(kT) )

其中k是一个常数,exp表示自然指数,且dE<0(温度总是降低的)。这条公式指明了:

1) 温度越高,出现一次能量差为dE的降温的概率就越大。

2) 温度越低,则出现降温的概率就越小。又由于dE总是小于0(不然怎么叫退火),因此dE/kT < 0 ,exp(dE/kT)取值是(0,1),那么P(dE)的函数取值范围是(0,1) 。

随着温度T的降低,P(dE)会逐渐降低。我们将一次向较差解的移动看做一次温度跳变过程,我们以概率P(dE)来接受这样的移动。也就是说,在用固体退火模拟组合优化问题,将内能E模拟为目标函数值 f,温度T演化成控制参数 t,即得到解组合优化问题的模拟退火演算法:

由初始解 i 和控制参数初值 t 开始,对当前解重复“产生新解→计算目标函数差→接受或丢弃”的迭代,并逐步衰减 t 值,算法终止时的当前解即为所得近似最优解。

4.模拟退火算法伪代码

1 代码
2
3 * f(y): 表示当前状态的目标函数值
4 * f(x): 表示新的可能状态的目标函数值
5 * T: 当前的温度
6 * T_min: 温度的下限,当温度低于该值时,停止搜索
7 * y: 当前的状态
8 * x: 新的可能的状态
9 * rand(0,1): 生成一个[0,1]区间内均匀分布的随机数
10 * exp(-dE/T): Metropolis准则,dE为目标函数的增量,即f(x)-f(y)
11 * 当dE<0时,即新状态的目标函数值小于当前状态,接受新状态
12 * 当dE>=0时,以一定的概率接受新状态,概率为exp(-dE/T)
13
14 while (T > T_min)
15 {
16     GetNewState(x); //产生新的状态x
17     dE = f(x) - f(y);
18     if (dE < 0)
19     {
20         y = x; //接受新的状态x
21     }
22     else
23     {
24         if (exp(-dE/T) > rand(0,1))
25         {
26             y = x; //以一定概率接受新的状态x,概率为exp(-dE/T)
27         }
28     }
29     T = T * r; //(r<1)是降温系数,随着迭代次数增加,温度逐渐降低
30     //r取值通常在(0.8,0.99)之间,r越小,降温越快,搜索过程越快,但可能会错过全局最优解
31     //r取大一些,降温过程慢一些,搜索时间相应增加,但有可能会得到更好的解。
32 }

5.使用模拟退火算法解决34个城市的TSP问题

核心算子代码:

# 计算总距离
def dis_cal(path, dist_mat):
    """
    适应度函数,计算目标函数值.
    :param dist_mat: 城市的距离矩阵
    :param path: 一个解
    :return: 目标函数值,即总距离
    """
    distance = 0
    for i in range(len(path) - 1):
        distance += dist_mat[path[i]][path[i + 1]]
    distance += dist_mat[path[-1]][path[0]]  # 从最后一个城市回到出发城市
    return distance


# 随机交换两个位置的方式,返回新的路径
def swap_operate(city_num, cur_path):
    path = cur_path.copy()
    u = np.random.randint(city_num)
    while True:
        v = np.random.randint(city_num)
        if u != v: break

    path[u], path[v] = path[v], path[u]
    return path


def gen_dist_mat(original_cities):
    D = original_cities[['经度', '纬度']].values * math.pi / 180
    city_cnt = len(original_cities)
    dist_mat = np.zeros((city_cnt, city_cnt))
    for i in range(city_cnt):
        for j in range(city_cnt):
            if i == j:
                # 相同城市不允许访问
                dist_mat[i][j] = 1000000
            else:
                # 单位:km
                dist_mat[i][j] = 6378.14 * math.acos(
                    math.cos(D[i][1]) * math.cos(D[j][1]) * math.cos(D[i][0] - D[j][0]) +
                    math.sin(D[i][1]) * math.sin(D[j][1]))
    return dist_mat

6.结果显示和图形绘制

# 结果输出
初始路径:
[18, 21, 23, 28, 14, 13, 30, 5, 33, 11, 4, 26, 1, 19, 24, 22, 10, 8, 9, 27, 16, 20, 17, 6, 32, 31, 0, 15, 25, 12, 2, 29, 7, 3]
优化路径:
[19, 22, 21, 23, 24, 25, 26, 27, 28, 31, 30, 32, 33, 0, 1, 29, 3, 4, 6, 5, 17, 16, 14, 15, 20, 10, 9, 12, 13, 11, 8, 2, 7, 18]

初始路径:

优化路径(红色部分)和初始路径(蓝色部分):

7.相关阅读

1.模拟退火算法求解TSP问题-python实现_tsp退火算法python-CSDN博客

2.干货 | 用模拟退火(SA, Simulated Annealing)算法解决旅行商问题

3.模拟退火算法求解旅行商问题(python实现)

4.模拟退火算法求解TSP问题-python实现

相关视频会在后续发布,欢迎关注我的bilibili:无形忍者的个人空间-无形忍者个人主页-哔哩哔哩视频

Logo

开放原子开发者工作坊旨在鼓励更多人参与开源活动,与志同道合的开发者们相互交流开发经验、分享开发心得、获取前沿技术趋势。工作坊有多种形式的开发者活动,如meetup、训练营等,主打技术交流,干货满满,真诚地邀请各位开发者共同参与!

更多推荐