【十分钟速成系列】--模拟退火算法（SA）（内含代码）

Zachary_coding

已于 2024-01-23 19:02:49 修改

阅读量3k

点赞数 39

文章标签：模拟退火算法算法机器学习

于 2024-01-22 15:37:57 首次发布

本文链接：https://blog.csdn.net/zzylalalala/article/details/135744432

版权

1.算法背景

模拟退火算法(Simulated Annealing,SA)最早的思想是由N. Metropolis等人于1953年提出。1983年,S. Kirkpatrick等成功地将退火思想引入到组合优化领域。它是基于Monte-Carlo 迭代求解策略的一种随机寻优算法，其出发点是基于物理中固体物质的退火过程与一般组合优化问题之间的相似性。

模拟退火算法来源于固体退火原理，将固体加温至充分高，再让其徐徐冷却，加温时，固体内部粒子随温升变为无序状，内能增大，而徐徐冷却时粒子渐趋有序，在每个温度都达到平衡态，最后在常温时达到基态，内能减为最小。

首先物体刚开始处于非晶体状态（左图）。我们将固体加温至充分高，固体内部粒子随温升变为无序状，能量增大，可以自由运动和重新排列（中图）。再让其缓慢冷却，粒子渐趋有序，在每个温度都达到平衡态，最后完全冷却时能量减为最小，物体形成晶体形态，这就是退火（右图）。

2.启发式算法

2.1启发式算法通俗理解：

启发式算法（Heuristic Algorithm）是一种求解最优化问题的通用方法，启发式搜索算法蕴含着许多人生哲学，它虽不是数学方法，但其思想更接近于人类解决问题的思想和一些人生道理。最后，用网上一段描述各种搜索算法的例子来作为总结:

为了找出地球上最高的山，一群有志的兔子们开始想办法。

(1) 兔子朝着比现在高的地方跳去。他们找到了不远处的最高山峰。但是这座山不一定是珠穆朗玛峰。这就是爬山法，它不能保证局部最优值就是全局最优值。

(2) 兔子喝醉了。他随机地跳了很长时间。这期间，它可能走向高处,也可能踏入平地。但是,他渐渐清醒了并朝他踏过的最方向跳去。这就是模拟退火。

(3) 兔子们知道一个兔的力量是渺小的。他们互相转告着，哪里的山已经找过，并组找过的每一座山他们都留下一只兔子做记号。他们制定了下一步去哪里寻找的策略。这就是禁忌搜索。

(4) 兔子们吃了失忆药片，并被发射到太空,然后随机落到了地球上的某些地方。他们不知道自己的使命是什么。但是,如果你过几年就杀死一部分海拔低的兔子，多产的兔子们自己就会找到珠穆朗玛峰。这就是遗传算法。

2.2元启发式算法概念：

元启发式算法是相对于最优化算法提出来的，一个问题的最优化算法可以求得该问题的最优解，而元启发式算法是一个基于直观或经验构造的算法，它可以在可接受的花费（指计算时间和空间）下给出问题的一个可行解，并且该可行解与最优解的偏离程度不一定可以事先预计。

2.3 爬山算法：

爬山算法：从当前的节点开始，和周围的邻居节点的值进行比较。如果当前节点是最大的，那么返回当前节点，作为最大值(即山峰最高点);反之就用最高的邻居节点来，替换当前节点，从而实现向山峰的高处攀爬的目的。如此循环直到达到最高点。
优点：避免遍历，通过启发选择部分节点，从而达到提高效率的目的。
缺点：因为不是全面搜索，所以结果可能不是最佳。

爬山算法本质上是沿着梯度方向进行寻优

如上图所示：假设C点为当前解，爬山算法搜索到A点这个局部最优解就会停止搜索，因为在A点无论向那个方向小幅度移动都不能得到更优的解。

3.模拟退火算法

3.1 算法步骤：

3.2 直观理解（链接）：

不同温度参数情况下，解空间寻找目标函数的全局最优解

3.3 步骤解释：

模拟退火算法包含两个部分即Metropolis算法和退火过程,，分别对应内循环和外循环。

外循环就是退火过程，将固体达到较高的温度（初始温度T（0）），然后按照降温系数alpha使温度按照一定的比例下降，当达到终止温度Tf时，冷却结束，即退火过程结束。
Metropolis算法是内循环，即在每次温度下，迭代L次，寻找在该温度下能量的最小值（即最优解）。下图中所示即为在一次温度下，跌代L次，固体能量发生的变化。在该温度下，整个迭代过程中温度不发生变化，能量发生变化，当前一个状态x(n)的能量大于后一个状态x(n+1)的能量时，状态x(n)的解没有状态x(n+1)的解好，所以接受状态x(n+1)。即上面步骤中的3.（1）

但是如果下一状态的能量比前一个状态的能量高时，该不该接受下一状态呢？在这里设置一个接受概率P,即如果下一状态的能量比前一个状态的能量高，则接受下一状态的概率为P.即上面步骤的3.（2）。

3.2 实例

如图中B、C、E为局部最优解，

假设开始状态在A，多次迭代之后更新到B的局部最优解，这时发现更新到B时，能力比A要低，则说明接近最优解了，因此百分百转移，状态到达B后，发现下一步能量上升了，如果是梯度下降则是不允许继续向前的，而这里会以一定的概率跳出这个坑，这各概率和当前的状态、能量等都有关系。所以说这个概率的设计是很重要的，下面从数学方面进行解释。

据Metropolis准则，粒子在温度T时趋于平衡的概率为 exp(-ΔE/(kT)) ，其中 E 为温度 T 时的内能，ΔE为其改变数, k 为 Boltzmann 常数。 Metropolis 准则常表示为

也就是说：

$\bullet$ 如果能量减小了，那么这种转移就被接受（概率为1）

$\bullet$ 如果能量增大了，就说明系统偏离全局最优值位置更远了，此时算法不会立刻将其抛弃，而是进行概率操作：首先在区间【0,1】产生一个均匀分布的随机数ϵ ，如果ϵ <P，则此种转移接受，否则拒绝转移，进入下一步，往复循环。其中P以能量的变化量和T进行决定概率P的大小，所以这个值是动态的。

总结：

若 f( Y(i+1) ) <= f( Y(i) ) (即移动后得到更优解)，则总是接受该移动；
若 f( Y(i+1) ) > f( Y(i) ) (即移动后的解比当前解要差)，则以一定的概率接受移动，而且这个概率随着时间推移逐渐降低（逐渐降低才能趋向稳定）相当于上图中，从 B 移向 BC 之间的小波峰时，每次右移（即接受一个更糟糕值）的概率在逐渐降低。如果这个坡特别长，那么很有可能最终我们并不会翻过这个坡。如果它不太长，这很有可能会翻过它，这取决于衰减 t 值的设定。

3.4 流程图：

4.应用场景

模拟退火算法作为一种通用的随机搜索算法，现已广泛用于VLSI设计、图像识别和神经网计算机的研究。模拟退火算法的应用如下：

模拟退火算法在VLSI设计中的应用
利用模拟退火算法进行VLSI(Very Large Scale Integration,超大规模集成电路)的最优设计，是目前模拟退火算法最成功的应用实例之一。用模拟退火算法几乎可以很好地完成所有优化的VLSI设计工作。如全局布线、布板、布局和逻辑最小化等等。
模拟退火算法作为一种通用的随机搜索算法，现已广泛用于VLSI设计、图像识别和神经网计算机的研究。模拟退火算法的应用如下：
模拟退火算法在VLSI设计中的应用
利用模拟退火算法进行VLSI(Very Large Scale Integration,超大规模集成电路)的最优设计，是目前模拟退火算法最成功的应用实例之一。用模拟退火算法几乎可以很好地完成所有优化的VLSI设计工作。如全局布线、布板、布局和逻辑最小化等等。

5.示例代码:(python）

用模拟退火算法解决旅行商问题：

import random
import math
import numpy as np
import matplotlib.pyplot as plt


class SA(object):
    def __init__(self, num_city, data):
        self.T0 = 4000
        self.Tend = 1e-3
        self.rate = 0.9995
        self.num_city = num_city
        self.scores = []
        self.location = data
        # fruits中存每一个个体是下标的list
        self.fires = []
        self.dis_mat = self.compute_dis_mat(num_city, data)
        self.fire = self.greedy_init(self.dis_mat,100,num_city)
        # 显示初始化后的路径
        init_pathlen = 1. / self.compute_pathlen(self.fire, self.dis_mat)
        init_best = self.location[self.fire]
        # 存储存储每个温度下的最终路径，画出收敛图
        self.iter_x = [0]
        self.iter_y = [1. / init_pathlen]
    def greedy_init(self, dis_mat, num_total, num_city):
        start_index = 0
        result = []
        for i in range(num_total):
            rest = [x for x in range(0, num_city)]
            # 所有起始点都已经生成了
            if start_index >= num_city:
                start_index = np.random.randint(0, num_city)
                result.append(result[start_index].copy())
                continue
            current = start_index
            rest.remove(current)
            # 找到一条最近邻路径
            result_one = [current]
            while len(rest) != 0:
                tmp_min = math.inf
                tmp_choose = -1
                for x in rest:
                    if dis_mat[current][x] < tmp_min:
                        tmp_min = dis_mat[current][x]
                        tmp_choose = x

                current = tmp_choose
                result_one.append(tmp_choose)
                rest.remove(tmp_choose)
            result.append(result_one)
            start_index += 1
        pathlens = self.compute_paths(result)
        sortindex = np.argsort(pathlens)
        index = sortindex[0]
        return result[index]

    # 初始化一条随机路径
    def random_init(self, num_city):
        tmp = [x for x in range(num_city)]
        random.shuffle(tmp)
        return tmp

    # 计算不同城市之间的距离
    def compute_dis_mat(self, num_city, location):
        dis_mat = np.zeros((num_city, num_city))
        for i in range(num_city):
            for j in range(num_city):
                if i == j:
                    dis_mat[i][j] = np.inf
                    continue
                a = location[i]
                b = location[j]
                tmp = np.sqrt(sum([(x[0] - x[1]) ** 2 for x in zip(a, b)]))
                dis_mat[i][j] = tmp
        return dis_mat

    # 计算路径长度
    def compute_pathlen(self, path, dis_mat):
        a = path[0]
        b = path[-1]
        result = dis_mat[a][b]
        for i in range(len(path) - 1):
            a = path[i]
            b = path[i + 1]
            result += dis_mat[a][b]
        return result

    # 计算一个温度下产生的一个群体的长度
    def compute_paths(self, paths):
        result = []
        for one in paths:
            length = self.compute_pathlen(one, self.dis_mat)
            result.append(length)
        return result

    # 产生一个新的解：随机交换两个元素的位置
    def get_new_fire(self, fire):
        fire = fire.copy()
        t = [x for x in range(len(fire))]
        a, b = np.random.choice(t, 2)
        fire[a:b] = fire[a:b][::-1]
        return fire

    # 退火策略，根据温度变化有一定概率接受差的解
    def eval_fire(self, raw, get, temp):
        len1 = self.compute_pathlen(raw, self.dis_mat)
        len2 = self.compute_pathlen(get, self.dis_mat)
        dc = len2 - len1
        p = max(1e-1, np.exp(-dc / temp))
        if len2 < len1:
            return get, len2
        elif np.random.rand() <= p:
            return get, len2
        else:
            return raw, len1

    # 模拟退火总流程
    def sa(self):
        count = 0
        # 记录最优解
        best_path = self.fire
        best_length = self.compute_pathlen(self.fire, self.dis_mat)

        while self.T0 > self.Tend:
            count += 1
            # 产生在这个温度下的随机解
            tmp_new = self.get_new_fire(self.fire.copy())
            # 根据温度判断是否选择这个解
            self.fire, file_len = self.eval_fire(best_path, tmp_new, self.T0)
            # 更新最优解
            if file_len < best_length:
                best_length = file_len
                best_path = self.fire
            # 降低温度
            self.T0 *= self.rate
            # 记录路径收敛曲线
            self.iter_x.append(count)
            self.iter_y.append(best_length)
            print(count, best_length)
        return best_length, best_path

    def run(self):
        best_length, best_path = self.sa()
        return self.location[best_path], best_length


# 读取数据
def read_tsp(path):
    lines = open(path, 'r').readlines()
    assert 'NODE_COORD_SECTION\n' in lines
    index = lines.index('NODE_COORD_SECTION\n')
    data = lines[index + 1:-1]
    tmp = []
    for line in data:
        line = line.strip().split(' ')
        if line[0] == 'EOF':
            continue
        tmpline = []
        for x in line:
            if x == '':
                continue
            else:
                tmpline.append(float(x))
        if tmpline == []:
            continue
        tmp.append(tmpline)
    data = tmp
    return data


data = read_tsp('data/st70.tsp')

data = np.array(data)
data = data[:, 1:]
show_data = np.vstack([data, data[0]])
Best, Best_path = math.inf, None

model = SA(num_city=data.shape[0], data=data.copy())
path, path_len = model.run()
print(path_len)
if path_len < Best:
    Best = path_len
    Best_path = path
# 加上一行因为会回到起点
Best_path = np.vstack([Best_path, Best_path[0]])
fig, axs = plt.subplots(2, 1, sharex=False, sharey=False)
axs[0].scatter(Best_path[:, 0], Best_path[:,1])
Best_path = np.vstack([Best_path, Best_path[0]])
axs[0].plot(Best_path[:, 0], Best_path[:, 1])
axs[0].set_title('规划结果')
iterations = model.iter_x
best_record = model.iter_y
axs[1].plot(iterations, best_record)
axs[1].set_title('收敛曲线')
plt.show()

参考博客：

通俗理解启发式算法 - 知乎 (zhihu.com)

模拟退火算法详解 - 知乎 (zhihu.com) 、模拟退火算法详细讲解（含实例python代码）_退火算法代码详解-CSDN博客

Zachary_coding

关注

39
点赞
踩
63

收藏

觉得还不错? 一键收藏
0
评论
【十分钟速成系列】--模拟退火算法（SA）（内含代码）

下图中所示即为在一次温度下，跌代L次，固体能量发生的变化。在该温度下，整个迭代过程中温度不发生变化，能量发生变化，当前一个状态x(n)的能量大于后一个状态x(n+1)的能量时，状态x(n)的解没有状态x(n+1)的解好，所以接受状态x(n+1)。假设开始状态在A，多次迭代之后更新到B的局部最优解，这时发现更新到B时，能力比A要低，则说明接近最优解了，因此百分百转移，状态到达B后，发现下一步能量上升了，如果是梯度下降则是不允许继续向前的，而这里会以一定的概率跳出这个坑，这各概率和当前的状态、能量等都有关系。
复制链接

扫一扫