一篇很久之前写的文章了,不过过了这么多年,模拟退火由于其效果和实现简单的优势,依然是智能算法中一个比较热门的算法,故老博客搬运而来。

模拟退火算法从外部来看就是一个优化问题的解析器,我们给他传递初始解和产生新解的方法,它就能不断产生新解,并比较最终返回一个近似最优解.由于数学建模对算法的时间限制不严,而模拟退火又较易于实现,因此它也是数学建模里较常用的一种智能算法.

模拟退火算法求解最大值 python 模拟退火算法应用实例_ci

快速使用

在介绍具体算法前,我们完全可以在短时间内使用上模拟退火.

例1:求min(x^2+y^2),x,y∈R.首先,我们提供一个初始解.文档main.m.1x=[2,2];

其次,构造出一个评价函数(或称目标函数).文档OptFun.m.

3function =(x)
y=x(1)^2+x(2)^2;
end

接着,构造一个能够不断根据旧解产生新解的函数.这里我们根据旧解以正态随机函数的形式产生新解.文档Arrise.m.1


function X=Arrise(x)
X(1)=normrnd(x(1),2);
X(2)=normrnd(x(2),2);
end

复制EzSA.m到文档夹.

最后,调用现成的模拟退火函数EzSA.文档main.m.1[ x,res ]=EzSA(x,@myFirstSA,@Arrise)

如果你看到一个进度条,那么恭喜你,你已经会使用模拟退火算法了!

让我们看看结果:

模拟退火算法求解最大值 python 模拟退火算法应用实例_ci_02

图像记载了我们之前尝试的解值,可以看出在数次迭代后数值处于稳定状态,表示这次模拟退火算法成功了.

同时,x 返回较优解,res返回较优值.1

2
3
4x =
-0.0035 -0.0027
res =
1.9562e-05

让我们总结一下模拟退火函数的使用步骤:提供或初始化一个初始解.

构造出一个评价函数(或称目标函数),该函数接收解,并返回一个数值(视值越小解越优).

构造一个能够不断根据旧解产生新解的函数(注意,这个函数的设计优劣直接影响到模拟退火效果的好坏).

调用现成的模拟退火函数EzSA(初始解,评价函数句柄,产生新解函数句柄).

一段时间后模拟退火算法结束,返回较优解和解值

例2:旅行商问题(TSP)

现有五个城市,彼此间距离如图所示,现在旅行商需要经过所有城市一次并回到出发点.我们需要为他规划最短路线.

模拟退火算法求解最大值 python 模拟退火算法应用实例_迭代_03

首先,以邻接矩阵存储图并提供初始解.文档main.m.1


global n
global graph
n=5;
graph=[0,7,6,1,3;7,0,3,7,8;6,3,0,12,11;1,7,12,0,2;3,8,11,2,0];
city=1:5;

其次是评价函数,设city为五个城市的访问顺序.文档computerTour.m.1


function len=computerTour(city) %计算路线总长度,每个城市只计算和下家城市之间的距离。
global n
global graph
len=0;
for i=1:n-1
len=len+graph(city(i),city(i+1));
end
len=len+graph(city(n),city(1));
end

接着,根据旧解产生新解的函数.文档perturbTour.m.1


function city=perturbTour(city)

%随机置换两个不同的城市的坐标

%产生随机扰动

global n
p1=randi([1,n]);
p2=randi([1,n]);
tmp=city(p1);
city(p1)=city(p2);
city(p2)=tmp;
end

最后,调用模拟退火函数(与第一步写在同一文档),并运行.文档main.m1[city,res]=EzSA(city,@computerTour,@perturbTour)

结果:


city =
4 1 3 2 5
res =
20

原理讲解

源代码:


36function [X resEnd]=EzSA(X,ObjFun,ArriseNew,iter,zero)
[ra,co]=size(X);
RES=[ObjFun(X)]; %每次迭代后的结果
temperature=100*co; %初始温度
if nargin==3
zero=1e-2;
iter=5e2; %内部蒙特卡洛循环迭代次数
end
if nargin==4
zero=1e-2;
end
h=waitbar(0,'SAing....');
while temperature>zero %停止迭代温度
for i=1:iter %多次迭代扰动,一种蒙特卡洛方法,温度降低之前多次实验
preRes=ObjFun(X); %目标函数计算结果
tmpX=ArriseNew(X); %产生随机扰动
newRes=ObjFun(tmpX); %计算新结果
delta_e=newRes-preRes; %新老结果的差值,相当于能量
if delta_e<0 %新结果好于旧结果,用新路线代替旧路线
X=tmpX;
else %温度越低,越不太可能接受新解;新老距离差值越大,越不太可能接受新解
if exp(-delta_e/temperature)>rand() %以概率选择是否接受新解 p=exp(-ΔE/T)
X=tmpX; %可能得到较差的解
end
end
end
RES=[RES ObjFun(X)];
temperature=temperature*0.99; %温度不断下降
waitbar((log(temperature/(100*co))/log(0.99))/(log(zero/(100*co))/log(0.99)),h,sprintf('Now Temperature:%.2f',temperature));
end
close(h)
plot(RES);
resEnd=RES(end);
end

结合代码再看开头的流程图.

初始化,计算初始解的解值,设置初始温度.

模拟退火结构上就是两重循环,外部循环检查温度并降温,内部不断地产生新解并与旧解比较.

若新解优于旧解则新解无条件被旧解替代.

否则,有一定概率(exp(-ΔE/T))新解取代旧解.注意这个环节正是模拟退火能跳脱局部最优解,取得全局最优解的关键.

由此,我们可以得知影响模拟退火效果的主要因素有:终止温度.一般上,终止温度越低,取得解越优.

内部迭代次数.一般上,内部迭代次数越多,取得解越优.

产生新解函数.

总结

模拟退火是对热力学退火过程的模拟,使算法在多项式时间内能给出一个近似最优解.由于MATLAB自带的模拟退火工具箱调用复杂且执行效果不理想,本文给出了较简单的函数原型和调用方法.该算法也包含以下优缺点(个人见解):

优点:相较于一般的蒙特卡洛算法,有更少的尝试次数,同时实现上并不比蒙特卡洛花更多时间.

相较于遗传算法等大型智能算法,模拟退火实现简单,并能返回较满意的结果.

目标函数可以自己定制,相较于普通的规划解析器,模拟退火能适用于更广的范围(NPC问题,甚至给给神经网络做优化).

对于离散型的变量有更优秀的效果.

缺点:内部本质上还是蒙特卡洛算法,新解与旧解本质上无关联.

相较于遗传算法,模拟退火难以控制算法的运行时间,EzSA的后面两个可选参数就是内部迭代次数和0度温度.而迭代次数给少了效果不理想,给多了有会增加等待时间.

对连续型的规划问题效果并不好.