模拟退火算法（附简单案例及详细matlab源码）

作者：非妃是公主
专栏：《智能优化算法》
博客地址：https://blog.csdn.net/myf_666
个性签：顺境不惰，逆境不馁，以心制境，万事可成。——曾国藩

专栏推荐

专栏名称	专栏地址
软件工程	专栏——软件工程
计算机图形学	专栏——计算机图形学
操作系统	专栏——操作系统
软件测试	专栏——软件测试
机器学习	专栏——机器学习
数据库	专栏——数据库
算法	专栏——算法

序

对于一个问题，解决他的第一步是建立模型，利用数学的方法，将问题抽象成一个数学模型（可能是一个函数，一个算法，一个方程……）。第二步，就是对这个问题进行求解，也叫做模型求解。

但对于模型求解来讲，不同的模型有着不同大小的难度。传统的算法，比如梯度下降，很可能会陷入局部最优解等问题，对于更加全局性的最优解搜索，智能优化算法有着很好的效果。而模拟退火算法，就是智能优化算法中的一种，有着很好的效果。

模拟退火算法（Simulated Annealing，SA）最早由Metropolis等人于 1953 年提出。1983 年Kirkpatrick等人第一次使用模拟退火算法求解组合优化问题后，它就发表在了 Science 上。¹直到今天，它依然被广泛使用，这篇文章将详细介绍模拟退火算法的基本原理，以及matlab的代码实现。

---

一、概论

模拟退火算法其实就是一个类似于仿生学的算法，模仿的就是物理退火的过程。我们炼钢的时候，如果我们急速冷凝，这时候的状态是不稳定的，原子间杂乱无章的排序，能量很高。而如果我们让钢水慢慢冷凝，很缓慢的降温，那么这个时候的状态就是很稳定的，各个分子都趋向于自己能量最低的位置。

而模拟退火算法，恰恰就是利用了物理退火这一过程的原理，求解一个优化目标（目标函数）的最小值。

---

二、物理退火

如果想说清楚模拟退火，必然绕不过物理退火！

1. 加温过程

增强粒子热运动，使其偏离平衡位置。当温度足够高时，固体将溶解为液体，从而消除系统原先可能存在的非均匀态，使随后进行的冷却过程以某一平衡态为起点。溶解过程于系统的能量增大过程相联系，系统能量也随温度的升高而增大。

---

2. 等温过程

通过物理学知识得知，对于与周围环境交换热量而温度不变的封闭系统，系统状态的自发变化总是朝着自由能降低的方向进行；当自由能达到最低时，系统达到平衡态。

---

3. 降温过程

随着温度降低，分子热运动减弱，趋向有序，系统能量逐渐降低，从而得到了低能量的晶体结构。

---

三、模拟退火原理

而模拟退火算法，就是要通过如上3个部分的操作，获得低能量的晶体结构（最优解）。

物理退火与模拟退火中的各个状态对应如下：

物理退火	模拟退火
粒子状态	每个状态对应一个（可行）解
能量最低态	最优解
溶解过程	设定初温（设定参数T的值）
等温过程	一个温度下，多次采样（Metropolis采样过程）
冷却	控制参数下降
能量	目标函数

主要思想，在搜索区间进行随机游走（通过生成随机数实现），再利用Metropolis抽样准则，使随机游走逐渐收敛于局部最优解。在这里，温度是一个重要的控制参数，这个参数的大小控制了随机过程向局部或全局最优解移动的快慢。具体的过程如下。

首先，定义一个$p$：
$$p=e^{-\frac{E_2-E_1}{T}}$$

这个公式表示，系统从$E_1$变化到$E_2$，其概率为$p$，就是上述公式左边的部分。

如果$E_2<E_1$，那么证明系统向更低的方向转移了，我们无条件接受此状态。

否则，以上述的概率，接受这个较坏的结果$E_2$。注意，这里的原则，不同于贪心，每次选择最好的，而是有一定的概率接受坏一些的。（这个概率受到能量差和温度两个因素的影响）

$$p(1\to 2)=\{\begin{array}{ll}1& E_2<E_1\\ e^{-\frac{E_2-E_1}{T}}& E_2>E_1\end{array}$$

这样通过一定的迭代次数，我们就会找到一个比较好的解了。

---

四、模拟退火的优点

总结一下模拟退火的优点如下：

• 以一定的概率接受恶化解。看似在接受恶化解，其实是寻找了全局最优解。
• 引进算法控制参数——温度。温度的引进，使得算法更加智能、灵活。在前期，进行跳跃搜索，更容易跳出局部最优解，提高搜索的全局性。在后期，进行幅度较小的搜索，更容易收敛。
• 对目标函数要求少。模拟退火算法其实是一种搜索或者说枚举的方法，我们经过大量的实验，最终得到最优解（或者一定程度上的较优解）。由于这种做法，因此对目标函数的要求很低。不需要依赖什么，直接猜就好。

---

五、算法具体流程

学习的过程是广度优先，层层深入的，前面大致了解了模拟退火算法的原理及思想，下面我们看一看模拟退火算法的具体算法流程。

---

1. 整体流程

1. 初始化温度$T_0$，初始解状态$X_0$（算法迭代起点）、每个$T$值的迭代次数$L$（也叫做Markov 链长度，就是在同一个温度下，我们猜测的次数）。
2. 对$k=1,\dots ,L$做第（3）至第（6）步；
3. 产生新解$X^{\prime }$；
4. 计算增量$\Delta E=E(X^{\prime })-E(X)$，其中$E(X)$为评价函数（越低越好）；
5. 若 $\Delta E＜0$，则接受 $X^{\prime }$作为新的当前解，否则以概率$e^{\frac{-\Delta E}{T}}$接受 $X^{\prime }$作为新的当前解；
6. 如果满足终止条件，则输出当前解为最优解，结束程序；
7. $T$逐渐减小，且$T\to 0$，然后转第 2 步。

流程图如下：

---

2. 细节处理

Ⅰ. 状态产生函数

主要就是通过在邻域内随机进行选择产生的。

---

Ⅱ. 初始温度

初始温度对于算法的影响较大，而且效果的好坏与求解问题的求解空间有关。初温越大，获得高质量解的可能性越高，但是相应的时间也就越长，具体可以在此处进行取舍。

一种很好的做法是，可以均匀抽样一组状态，以各状态目标值的方差作为初温。

• 这样，如果方差越大，证明该问题越有可能陷入到局部最优解（多峰型优化），因此初温较高，防止陷入。

• 如果方差很小，说明该问题不容易陷入局部最优(单峰型优化)，因此，初温较低，可以快速收敛。

---

Ⅲ. 退温函数

退温函数即温度更新函数，用于在外循环中修改温度值。目前，最常用的温度更新函数为指数退温函数，即$T(n+1)=K\times T(n)$，其中$0＜K＜1$，$K$为一个非常接近于1的常数。

---

Ⅳ. Markov 链长度

Markov 链长度是在等温条件下进行迭代优化的次数，其选取原则是在衰减参数 $T$ 的衰减函数已选定的前提下，还要产生随机数的次数，一般 L 取100～1000。

---

Ⅴ. 算法终止准则

算法停止的条件。常用的有，温度降低到一定的阈值结束，迭代一定的次数后结束，最优值连续保持不变（或者变化值$<\delta$）时停止搜索。

---

六、仿真实例：模拟退火求解函数最小值

1. 题目

计算函数$f(x)={\sum }_{i=1}^n{x}_i^2(-20<=x_i<=20)$的最小值，其中个体$x$的维数为$n=10$。

---

2. 分析

这时一个平方和函数，只有一个极小值$x=(0,0,...,0)$，理论最小值$f(0,0,...,0)=0$。

Markov链长度初始化为$L=200$，衰减参数设置为0.998，补偿因子为$S=0.01$，初始温度$T=100$，容差为$YZ=1\times 10^{-8}$（用于判断结束条件，如果最优值连续变化小于$YZ$，那么就终止搜索，输出这个较优解，其它参数上边都已经提到过了）；随机产生初始解，并计算目标函数值。

---

3. matlab求解

%%%%%%%%%%%%%%%%%%%%%%模拟退火算法解决函数极值%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%初始化%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all;                      %清除所有变量
close all;                      %清图
clc;                            %清屏
D=10;                           %变量维数 
Xs=20;                          %上限                                
Xx=-20;                         %下限
%%%%%%%%%%%%%%%%%%%%%%%%%%%冷却表参数%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
L = 200;                        %马可夫链长度
K = 0.998;                      %衰减参数
S = 0.01;                       %步长因子
T=100;                          %初始温度
YZ = 1e-8;                      %容差
P = 0;                          %Metropolis过程中总接受点
%%%%%%%%%%%%%%%%%%%%%%%%%%随机选点 初值设定%%%%%%%%%%%%%%%%%%%%%%%%%
PreX = rand(D,1)*(Xs-Xx)+Xx;
PreBestX = PreX;
PreX = rand(D,1)*(Xs-Xx)+Xx;
BestX = PreX;
%%%%%%%%%%%每迭代一次退火一次(降温), 直到满足迭代条件为止%%%%%%%%%%%%
deta=abs(func1(BestX)-func1(PreBestX));
while (deta > YZ) && (T>0.001)
    T=K*T; 
    %%%%%%%%%%%%%%%%%%%%%在当前温度T下迭代次数%%%%%%%%%%%%%%%%%%%%%%
    for i=1:L  
        %%%%%%%%%%%%%%%%%在此点附近随机选下一点%%%%%%%%%%%%%%%%%%%%%
        NextX = (1 - S) * PreX + S * (rand(D,1) *(Xs-Xx)+Xx);
        %%%%%%%%%%%%%%%%%边界条件处理%%%%%%%%%%%%%%%%%%%%%%%%%%
        for ii=1:D
            if NextX(ii)>Xs || NextX(ii)<Xx
                NextX(ii)=rand *(Xs-Xx)+Xx;
            end
        end            
        %%%%%%%%%%%%%%%%%%%%%%%是否全局最优解%%%%%%%%%%%%%%%%%%%%%%
        if (func1(BestX) > func1(NextX))
            %%%%%%%%%%%%%%%%%%保留上一个最优解%%%%%%%%%%%%%%%%%%%%%
            PreBestX = BestX;
            %%%%%%%%%%%%%%%%%%%此为新的最优解%%%%%%%%%%%%%%%%%%%%%%
            BestX=NextX;
        end
        %%%%%%%%%%%%%%%%%%%%%%%% Metropolis过程%%%%%%%%%%%%%%%%%%%
        if( func1(PreX) - func1(NextX) > 0 )
            %%%%%%%%%%%%%%%%%%%%%%%接受新解%%%%%%%%%%%%%%%%%%%%%%%%
            PreX=NextX;
            P=P+1;
        else
            changer = -1*(func1(NextX)-func1(PreX))/ T ;
            p1=exp(changer);
            %%%%%%%%%%%%%%%%%%%%%%%%接受较差的解%%%%%%%%%%%%%%%%%%%%
            if p1 > rand        
                PreX=NextX;
                P=P+1;         
            end
        end
        trace(P+1) = func1(BestX);
    end
    deta = abs(func1(BestX) - func1(PreBestX)); 
end
disp('最小值在点:');
BestX
disp( '最小值为:');
func1(BestX)
figure
plot(trace(2:end))
xlabel('迭代次数')
ylabel('目标函数值')
title('最优解变化曲线')

其中，目标函数（适应度函数）fun1的定义如下：

%%%%%%%%%%%%%%%%%%%%%%%%%适应度函数%%%%%%%%%%%%%%%%%%%%%%%
function result=func1(x)
summ=sum(x.^2);
result=summ;

---

4. 求解结果及分析

求解的适应度变化曲线如下：

可以看到，随着迭代次数的增加，目标函数值下降速度还是很快的，说明算法收敛速度较好。

最终搜索结果如下：

可以看到，我们在理论最优解 0 附近找到了近似最优解 0.0056，求解效果还是不错的。

继续运行代码，求解结果如下：

可以看到求解过程以及最终的结果都不太相同，因此，对于模拟退火这种启发式智能优化算法，求解结果是有一定随机性的。但是，会随着迭代次数的增加，越发趋于稳定。

---

七、模拟退火的一些改进方向

• 增加记忆功能
• 增加升温或重升温过程。
• 对每一当前状态，采用多次搜索策略，以概率接受区域内的最优状态，而不是标准模拟退火算法的单次比较方式。
• 与其他搜索机制的算法（如遗传算法、免疫算法等）相结合。可以综合其他方法的优点，提高运行效率和求解质量。²

---

the end……

模拟退火算法到这里就要结束啦~~到此既是缘分，欢迎您的点赞、评论、收藏！关注我，不迷路，我们下期再见！！

😘😘😘 我是Cherries，一位计算机科班在校大学生，写博客用来记录自己平时的所思所想！
💞💞💞 内容繁杂，又才疏学浅，难免存在错误，欢迎各位大佬的批评指正！
👋👋👋 我们相互交流，共同进步！

注：本文由非妃是公主发布于https://blog.csdn.net/myf_666，转载请务必标明原文链接：https://blog.csdn.net/myf_666/article/details/129061318

---

1. Kirkpatrick S,Gelatt C,Vecchi M.Optimization by simulated Anealing.Science,1983(220):671-680. ↩︎

2. 包子阳,智能优化算法及其matlab实例.电子工业出版社. ↩︎